File indexing completed on 2025-02-23 09:21:47
0001
0002 void BinLogX(TH1* h)
0003 {
0004 TAxis *axis = h->GetXaxis();
0005 int bins = axis->GetNbins();
0006
0007 Axis_t from = axis->GetXmin();
0008 Axis_t to = axis->GetXmax();
0009 Axis_t width = (to - from) / bins;
0010 Axis_t *new_bins = new Axis_t[bins + 1];
0011
0012 for (int i = 0; i <= bins; i++) {
0013 new_bins[i] = TMath::Power(10, from + i * width);
0014 }
0015 axis->Set(bins, new_bins);
0016 delete []new_bins;
0017 }
0018
0019
0020 void plot(){
0021
0022 double dens = 1000;
0023
0024 double R =50;
0025
0026
0027 TFile *fin = new TFile("AuNP.root");
0028 TH1F *h1neve = (TH1F*)fin->Get("h1Events");
0029 TH1F *h1Edep = (TH1F*)fin->Get("h1Edep");
0030 TH1F *h1secnp_cha = (TH1F*)fin->Get("h1SecEnergyNP_charged");
0031 TH1F *h1secnp_nut = (TH1F*)fin->Get("h1SecEnergyNP_nutral");
0032 TH1F *h1secnpsurf_cha= (TH1F*)fin->Get("h1SecEnergyNPSurf_charged");
0033 TH1F *h1secnpsurf_nut= (TH1F*)fin->Get("h1SecEnergyNPSurf_nutral");
0034 TH1F *h1sec_cha = (TH1F*)fin->Get("h1Sec_charged");
0035 TH1F *h1sec_nut = (TH1F*)fin->Get("h1Sec_nutral");
0036 TH1F *h1chem_0 = (TH1F*)fin->Get("h1Chem_0");
0037 TH1F *h1chem_1 = (TH1F*)fin->Get("h1Chem_1");
0038 TH1F *h1chem_2 = (TH1F*)fin->Get("h1Chem_2");
0039 TH1F *h1chem_3 = (TH1F*)fin->Get("h1Chem_3");
0040 TH1F *h1chem_4 = (TH1F*)fin->Get("h1Chem_4");
0041 TH1F *h1chem_5 = (TH1F*)fin->Get("h1Chem_5");
0042 TH1F *h1chem_6 = (TH1F*)fin->Get("h1Chem_6");
0043 TH1F *h1chem_7 = (TH1F*)fin->Get("h1Chem_7");
0044
0045 TH2F *h2Edep = (TH2F*)fin->Get("h2Edep");
0046 TH2F *h2sec2_cha = (TH2F*)fin->Get("h2SecEnergyAbs_charged");
0047 TH2F *h2sec2_nut = (TH2F*)fin->Get("h2SecEnergyAbs_nutral");
0048
0049 double neve = h1neve->GetBinContent(1);
0050 h1Edep ->Scale(1./neve);
0051 h1secnp_cha ->Scale(1./neve);
0052 h1secnp_nut ->Scale(1./neve);
0053 h1secnpsurf_cha->Scale(1./neve);
0054 h1secnpsurf_nut->Scale(1./neve);
0055 h1sec_cha ->Scale(1./neve);
0056 h1sec_nut ->Scale(1./neve);
0057 h1chem_0 ->Scale(1./neve);
0058 h1chem_1 ->Scale(1./neve);
0059 h1chem_2 ->Scale(1./neve);
0060 h1chem_3 ->Scale(1./neve);
0061 h1chem_4 ->Scale(1./neve);
0062 h1chem_5 ->Scale(1./neve);
0063 h1chem_6 ->Scale(1./neve);
0064 h1chem_7 ->Scale(1./neve);
0065 h2Edep ->Scale(1./neve);
0066 h2sec2_cha ->Scale(1./neve);
0067 h2sec2_nut ->Scale(1./neve);
0068
0069 double val_cha=0;
0070 double err_cha=0;
0071 double val_nut=0;
0072 double err_nut=0;
0073
0074 int NR = h1sec_cha -> GetNbinsX();
0075 TH1D *h1sec_tot = new TH1D("h1Sec_tot","h1Sec_tot",NR,1,5);
0076 BinLogX(h1sec_tot);
0077 for(int i=0;i<NR;i++){
0078 val_cha = h1sec_cha->GetBinContent(i+1);
0079 err_cha = h1sec_cha->GetBinError (i+1);
0080 val_nut = h1sec_nut->GetBinContent(i+1);
0081 err_nut = h1sec_nut->GetBinError (i+1);
0082
0083
0084
0085
0086 h1sec_cha -> SetBinContent(i+1,val_cha);
0087 h1sec_cha -> SetBinError (i+1,err_cha);
0088 h1sec_nut -> SetBinContent(i+1,val_nut);
0089 h1sec_nut -> SetBinError (i+1,err_nut);
0090 h1sec_tot -> SetBinContent(i+1,val_cha+val_nut);
0091 h1sec_tot -> SetBinError (i+1,err_cha+err_nut);
0092 }
0093
0094 TH1D *h1chem_tot = new TH1D("h1Chem_tot","h1Chem_tot",NR,1,5);
0095 BinLogX(h1chem_tot);
0096 for(int i=0;i<NR;i++){
0097 double val_0 = h1chem_0->GetBinContent(i+1)/h1chem_0->GetBinWidth(i+1);
0098 double val_1 = h1chem_1->GetBinContent(i+1)/h1chem_1->GetBinWidth(i+1);
0099 double val_2 = h1chem_2->GetBinContent(i+1)/h1chem_2->GetBinWidth(i+1);
0100 double val_3 = h1chem_3->GetBinContent(i+1)/h1chem_3->GetBinWidth(i+1);
0101 double val_4 = h1chem_4->GetBinContent(i+1)/h1chem_4->GetBinWidth(i+1);
0102 double val_5 = h1chem_5->GetBinContent(i+1)/h1chem_5->GetBinWidth(i+1);
0103 double val_6 = h1chem_6->GetBinContent(i+1)/h1chem_6->GetBinWidth(i+1);
0104 double val_7 = h1chem_7->GetBinContent(i+1)/h1chem_7->GetBinWidth(i+1);
0105 double err_0 = h1chem_0->GetBinError (i+1)/h1chem_0->GetBinWidth(i+1);
0106 double err_1 = h1chem_1->GetBinError (i+1)/h1chem_1->GetBinWidth(i+1);
0107 double err_2 = h1chem_2->GetBinError (i+1)/h1chem_2->GetBinWidth(i+1);
0108 double err_3 = h1chem_3->GetBinError (i+1)/h1chem_3->GetBinWidth(i+1);
0109 double err_4 = h1chem_4->GetBinError (i+1)/h1chem_4->GetBinWidth(i+1);
0110 double err_5 = h1chem_5->GetBinError (i+1)/h1chem_5->GetBinWidth(i+1);
0111 double err_6 = h1chem_6->GetBinError (i+1)/h1chem_6->GetBinWidth(i+1);
0112 double err_7 = h1chem_7->GetBinError (i+1)/h1chem_7->GetBinWidth(i+1);
0113 h1chem_tot -> SetBinContent(i+1,val_0+val_1+val_2+val_3+val_4+val_5+val_6+val_7);
0114 h1chem_tot -> SetBinError (i+1,err_0+err_1+err_2+err_3+err_4+err_5+err_6+err_7);
0115 }
0116
0117
0118 NR = h1secnp_cha -> GetNbinsX();
0119 TH1D *h1secnp_tot = new TH1D("h1SecEnergyNP" ,"Energy Spectra in NP" ,NR,0,6);
0120 BinLogX(h1secnp_tot);
0121 NR = h1secnpsurf_cha -> GetNbinsX();
0122 TH1D *h1secnpsurf_tot = new TH1D("h1SecEnergyNPSurf","Energy Spectra at NP surface",NR,0,6);
0123 BinLogX(h1secnpsurf_tot);
0124 for(int i=0;i<NR;i++){
0125 val_cha = h1secnp_cha->GetBinContent(i+1);
0126 err_cha = h1secnp_cha->GetBinError (i+1);
0127 val_nut = h1secnp_nut->GetBinContent(i+1);
0128 err_nut = h1secnp_nut->GetBinError (i+1);
0129
0130
0131
0132
0133 h1secnp_cha -> SetBinContent(i+1,val_cha);
0134 h1secnp_cha -> SetBinError (i+1,err_cha);
0135 h1secnp_nut -> SetBinContent(i+1,val_nut);
0136 h1secnp_nut -> SetBinError (i+1,err_nut);
0137 h1secnp_tot -> SetBinContent(i+1,val_cha+val_nut);
0138 h1secnp_tot -> SetBinError (i+1,err_cha+err_nut);
0139
0140
0141 val_cha = h1secnpsurf_cha->GetBinContent(i+1);
0142 err_cha = h1secnpsurf_cha->GetBinError (i+1);
0143 val_nut = h1secnpsurf_nut->GetBinContent(i+1);
0144 err_nut = h1secnpsurf_nut->GetBinError (i+1);
0145
0146
0147
0148
0149 h1secnpsurf_cha -> SetBinContent(i+1,val_cha);
0150 h1secnpsurf_cha -> SetBinError (i+1,err_cha);
0151 h1secnpsurf_nut -> SetBinContent(i+1,val_nut);
0152 h1secnpsurf_nut -> SetBinError (i+1,err_nut);
0153 h1secnpsurf_tot -> SetBinContent(i+1,val_cha+val_nut);
0154 h1secnpsurf_tot -> SetBinError (i+1,err_cha+err_nut);
0155
0156 }
0157
0158 NR = h1Edep -> GetNbinsX();
0159 for(int i=0;i<NR;i++){
0160 double rmax = h1Edep->GetXaxis()->GetBinUpEdge (i+1)*1.E-9;
0161 double rmin = h1Edep->GetXaxis()->GetBinLowEdge(i+1)*1.E-9;
0162 double vol = 4./3.*TMath::Pi()*(pow(rmax,3)-pow(rmin,3));
0163 double mass = dens*vol;
0164 h1Edep -> SetBinContent(i+1,h1Edep->GetBinContent(i+1)/mass);
0165 h1Edep -> SetBinError (i+1,h1Edep->GetBinError (i+1)/mass);
0166 }
0167
0168
0169 int Nazm =h2Edep->GetXaxis()->GetNbins()-1;
0170 NR =h2Edep->GetYaxis()->GetNbins()-1;
0171 double Rmax =h2Edep->GetYaxis()->GetBinUpEdge(NR);
0172 TH2D *h2Edep_pol = new TH2D("h2Edep_pol","h2Edep_pol",Nazm,0,360,NR,0,Rmax);
0173 for(int i=0;i<Nazm;i++){
0174 for(int j=0;j<NR;j++){
0175 double height = 10*1.E-9;
0176 double rmax = h2Edep_pol->GetYaxis()->GetBinUpEdge (j+1)*1.E-9;
0177 double rmin = h2Edep_pol->GetYaxis()->GetBinLowEdge (j+1)*1.E-9;
0178 double vol = (TMath::Pi()*pow(rmax,2)-TMath::Pi()*pow(rmin,2))*height;
0179 double mass = dens*vol;
0180 h2Edep_pol->SetBinContent(i+1,j+1,h2Edep->GetBinContent(i+1,j+1)/mass);
0181 }
0182 }
0183 int NRX = h2sec2_cha->GetXaxis()->GetNbins();
0184 int NRY = h2sec2_cha->GetYaxis()->GetNbins();
0185 TH2D *h2sec2_tot = (TH2D*)h2sec2_cha->Clone("h2sec2_tot");
0186 for(int i=0;i<NRX;i++){
0187 for(int j=0;j<NRY;j++){
0188 val_cha = h2sec2_cha->GetBinContent(i+1,j+1);
0189 val_nut = h2sec2_nut->GetBinContent(i+1,j+1);
0190 h2sec2_tot->SetBinContent(i+1,j+1,val_cha+val_nut);
0191 }
0192 }
0193
0194
0195
0196
0197
0198
0199
0200 TCanvas *cedep = new TCanvas("cedep","cedep",700,700);
0201 cedep->cd(1)->SetTheta(90);
0202 cedep->cd(1)->SetPhi(0);
0203 cedep->cd(1)->SetLogz();
0204 h2Edep_pol->RebinX(6);
0205 h2Edep_pol->RebinY(10);
0206
0207 h2Edep_pol->SetAxisRange(0,150,"Y");
0208 h2Edep_pol->SetTitle("");
0209 h2Edep_pol->SetXTitle("Angle [deg]");
0210 h2Edep_pol->Draw("surf1POL");
0211
0212 TCanvas *cedepsec = new TCanvas("cedepsec","cedepsec",1800,700);
0213 cedepsec->Divide(2,1);
0214 cedepsec->cd(1)->SetLogy();
0215 cedepsec->cd(1)->SetLogx();
0216 h1Edep->SetYTitle("Energy Deposite [Gy/event]");
0217 h1Edep->SetXTitle("Distance from NP center[nm]");
0218 h1Edep->SetMaximum(100);
0219 h1Edep->SetMinimum(0.0000000001);
0220 h1Edep->SetAxisRange(0,1000000);
0221 h1Edep->Draw();
0222 cedepsec->cd(2)->SetLogy();
0223 cedepsec->cd(2)->SetLogx();
0224 h1sec_tot->SetLineColor(kBlack);
0225 h1sec_cha->SetLineColor(kRed);
0226 h1sec_nut->SetLineColor(kBlue);
0227 h1sec_tot->SetMarkerColor(kBlack);
0228 h1sec_cha->SetMarkerColor(kRed);
0229 h1sec_nut->SetMarkerColor(kBlue);
0230 h1sec_tot->SetAxisRange(0,1000000);
0231 h1sec_tot->SetMaximum(10);
0232 h1sec_tot->SetMinimum(0.000001);
0233 h1sec_tot->SetTitle ("Number of Generated Secondaries");
0234 h1sec_tot->SetYTitle("N/nm/event");
0235 h1sec_tot->SetXTitle("Distance from NP center[nm]");
0236
0237
0238 h1sec_tot->Draw();
0239
0240
0241
0242
0243
0244
0245
0246
0247
0248
0249
0250
0251
0252 TCanvas *csec = new TCanvas("csec","csec",2400,700);
0253 csec->Divide(3,1);
0254 csec->cd(1)->SetLogx();
0255 csec->cd(1)->SetLogy();
0256 h1secnp_tot->SetLineColor(kBlack);
0257 h1secnp_cha->SetLineColor(kRed);
0258 h1secnp_nut->SetLineColor(kBlue);
0259 h1secnp_tot->SetMarkerColor(kBlack);
0260 h1secnp_cha->SetMarkerColor(kRed);
0261 h1secnp_nut->SetMarkerColor(kBlue);
0262 h1secnp_tot->RebinX(10);
0263 h1secnp_tot->Scale(1./10.);
0264 h1secnp_tot->SetMaximum(20);
0265 h1secnp_tot->SetMinimum(0.000000001);
0266 h1secnp_tot->SetAxisRange(1,1.E6);
0267 h1secnp_tot->SetTitle ("Secondary Energy produced in NP");
0268 h1secnp_tot->SetYTitle("N/event");
0269 h1secnp_tot->SetXTitle("Secandary Energy [eV]");
0270
0271
0272 h1secnp_tot->DrawCopy();
0273 csec->cd(2)->SetLogx();
0274 csec->cd(2)->SetLogy();
0275 h1secnpsurf_tot->SetLineColor(kBlack);
0276 h1secnpsurf_cha->SetLineColor(kRed);
0277 h1secnpsurf_nut->SetLineColor(kBlue);
0278 h1secnpsurf_tot->SetMarkerColor(kBlack);
0279 h1secnpsurf_cha->SetMarkerColor(kRed);
0280 h1secnpsurf_nut->SetMarkerColor(kBlue);
0281 h1secnpsurf_tot->RebinX(10);
0282 h1secnpsurf_tot->Scale(1./10.);
0283 h1secnpsurf_tot->SetMaximum(20);
0284 h1secnpsurf_tot->SetMinimum(0.000000001);
0285 h1secnpsurf_tot->SetAxisRange(1,1.E6);
0286 h1secnpsurf_tot->SetTitle ("Secondary Energy at NP surface");
0287 h1secnpsurf_tot->SetYTitle("N/event");
0288 h1secnpsurf_tot->SetXTitle("Secandary Energy [eV]");
0289
0290
0291 h1secnpsurf_tot->DrawCopy();
0292 csec->cd(3)->SetLogy();
0293 h2sec2_tot->SetTitle ("Secondary Energy vs production point");
0294 h2sec2_tot->SetYTitle("Secondary Energy [eV]");
0295 h2sec2_tot->SetXTitle("Distance from NP center [nm]");
0296 h2sec2_tot->GetXaxis()->SetRange(0,1000);
0297 h2sec2_tot->Draw("colz");
0298
0299 }
0300