dna/microprox/plot.C

0001 // *********************************************************************
0002 // To execute this macro under ROOT after your simulation ended,
0003 //   1 - launch ROOT (usually type 'root' at your machine's prompt)
0004 //   2 - type '.X plot.C' at the ROOT session prompt
0005 //
0006 // Author: Sebastien Incerti
0007 // Date: March 2nd, 2019
0008 // The Geant4-DNA collaboration
0009 // *********************************************************************
0010
0011 void SetLeafAddress(TNtuple* ntuple, const char* name, void* address);
0012
0013 void plot()
0014 {
0015 gROOT->Reset();
0016 gStyle->SetPalette(1);
0017 gROOT->SetStyle("Plain");
0018 gStyle->SetOptStat(00000);
0019
0020 //********************
0021 Int_t nbRadius = 101;
0022 //********************
0023
0024 TCanvas *c1 = new TCanvas ("c1","",60,60,800,800);
0025 Int_t mycolor;
0026
0027 TFile f("t.root");
0028 mycolor=4;
0029
0030 TNtuple* ntuple;
0031 ntuple = (TNtuple*)f.Get("t");
0032
0033 bool rowWise = true;
0034 TBranch* eventBranch = ntuple->FindBranch("row_wise_branch");
0035 if ( ! eventBranch ) rowWise = false;
0036
0037 Double_t radius1,nofHits,nbEdep,edep,radius2,Einc;
0038 Int_t noRadius;
0039
0040 if ( ! rowWise ) {
0041   ntuple->SetBranchAddress("radius1",&radius1);
0042   ntuple->SetBranchAddress("noRadius",&noRadius);
0043   ntuple->SetBranchAddress("nbHits",&nofHits);
0044   ntuple->SetBranchAddress("nbScoredHits",&nbEdep);
0045   ntuple->SetBranchAddress("edep",&edep);
0046   ntuple->SetBranchAddress("radius2",&radius2);
0047   ntuple->SetBranchAddress("Einc",&Einc);
0048   }
0049 else {
0050   SetLeafAddress(ntuple, "radius1",&radius1);
0051   SetLeafAddress(ntuple, "noRadius",&noRadius);
0052   SetLeafAddress(ntuple, "nbHits",&nofHits);
0053   SetLeafAddress(ntuple, "nbScoredHits",&nbEdep);
0054   SetLeafAddress(ntuple, "edep",&edep);
0055   SetLeafAddress(ntuple, "radius2",&radius2);
0056   SetLeafAddress(ntuple, "Einc",&Einc);
0057 }
0058
0059 Int_t nentries = (Int_t)ntuple->GetEntries();
0060
0061 //
0062
0063 Double_t t[1000]; // 1000 is the max number of radius values
0064 Double_t population[1000];
0065 Double_t myRad[1000];
0066
0067 for (Int_t i=0; i<1000; i++)
0068 {
0069  t[i]=0;
0070  population[i]=0;
0071  myRad[i]=0;
0072 }
0073
0074 Int_t event = 0;
0075
0076 for (Int_t i=0; i<nentries; i++)
0077 {
0078   ntuple->GetEntry(i);
0079   t[noRadius] = t[noRadius] + edep;
0080   population[noRadius]=population[noRadius]+1;
0081   myRad[noRadius] = radius1;
0082 }
0083
0084 // Mean
0085
0086 for (Int_t j=1; j<nbRadius; j++)
0087 {
0088   t[j] = t[j]/population[j];
0089   t[j] = t[j]/(myRad[j+1]-myRad[j]);
0090   //cout << j << " " << myRad[j] << " " << myRad[j+1]
0091   // << " " << t[j] << " " << population[j] << endl;
0092 }
0093
0094 //
0095
0096  c1->cd(1);
0097
0098  TGraph* gr1 =new TGraph(nbRadius,myRad,t);
0099  gr1->SetMarkerColor(2);
0100  gr1->SetMarkerStyle(20);
0101  gr1->SetMarkerSize(1);
0102  gr1->SetLineColor(2);
0103  gr1->SetTitle("");
0104  gr1->GetXaxis()->SetLimits(0.1,100);
0105  gr1->GetYaxis()->SetLimits(0.,22);
0106  gr1->GetXaxis()->SetLabelSize(0.025);
0107  gr1->GetYaxis()->SetLabelSize(0.025);
0108  gr1->GetXaxis()->SetTitleSize(0.035);
0109  gr1->GetYaxis()->SetTitleSize(0.035);
0110  gr1->GetXaxis()->SetTitleOffset(1.4);
0111  gr1->GetYaxis()->SetTitleOffset(1.4);
0112  gr1->GetXaxis()->SetTitle("r (nm)");
0113  gr1->GetYaxis()->SetTitle("t (eV/nm)");
0114  gr1->Draw("");
0115  gPad->SetLogx();
0116
0117 }
0118
0119 void SetLeafAddress(TNtuple* ntuple, const char* name, void* address) {
0120   TLeaf* leaf = ntuple->FindLeaf(name);
0121   if ( ! leaf ) {
0122     std::cerr << "Error in <SetLeafAddress>: unknown leaf --> " << name << std::endl;
0123     return;
0124   }
0125   leaf->SetAddress(address);
0126 }