]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGGA/GammaConv/macros/Extract_Pi0_Characteristics.C
Transition PWG4 --> PWGGA
[u/mrichter/AliRoot.git] / PWGGA / GammaConv / macros / Extract_Pi0_Characteristics.C
CommitLineData
baf06eb6 1
2#include <fstream>
3#include <Riostream.h>
4/*
5 *
6 *
7 *
8 */
9
10extern TRandom *gRandom;
11extern TBenchmark *gBenchmark;
12extern TSystem *gSystem;
13
14void Extract_Pi0_Characteristics(const char *cutSelection="", const char *inputRootFile = "AnalysisResults",const char *path = "./",const char*outputDir = "./Output/",Bool_t makeMappingPlots=kTRUE){
15
16 gROOT->Reset();
17 gROOT->SetStyle("Plain");
18 gStyle->SetOptFit(0);
19 // rebinning
20 Int_t RebinMass = 2;
21
22 Double_t BGFit_range[2] = {0.17,0.3};
23
24 Double_t Pi0Mass = 0.135;
25 Double_t Pi0Mass_range[2] = {0.12,0.15};
26 Double_t Pi0Width = 0.003;
27 Double_t Pi0Width_range[2] = {0.001,0.007};
28 Double_t Pi0Slope = 0.007;
29 Double_t Pi0Slope_range[2] = {0.001,0.016};
30
31 Double_t Pi0FitRange[2] = {0.1,0.16};
32
33 // pt bins
34 Int_t NBinsPt = 31;
35 Double_t BinsPt[] = {0.0, 0.4, 0.8, 1.2, 1.6, 2.0, 2.4, 2.8, 3.2, 3.6, 4., 5., 6.,8.,10.,13.,15.,16.,17.,18.,19.,20.,21.,22.,23.,24.,25.,26.,27.,28.,29.,30.};
36
37 // input file
38 TString filename = Form("%s%s.root",path,inputRootFile);
39
40 TFile f(filename.Data());
41
42 TFile *outputFile;
43 if(makeMappingPlots==kTRUE){
44 outputFile = new TFile(Form("%sPi0Characteristics.root",outputDir),"UPDATE");
45 }
46 TDirectory *pwg4dir =(TDirectory*)f.Get(Form("PWG4_GammaConversion_%s",cutSelection));
47 TList *fHistosGammaConversion = (TList*)pwg4dir->Get(Form("histogramsAliGammaConversion_%s",cutSelection));
48 TList *fESDContainer = (TList*)fHistosGammaConversion->FindObject("ESD histograms");
49 TList *fMCContainer = (TList*)fHistosGammaConversion->FindObject("MC histograms");
50 TList *fBackgroundContainer = (TList*)fHistosGammaConversion->FindObject("Background histograms");
51
10e3319b 52 TH1F *ESD_Mother_InvMass =fESDContainer->FindObject("ESD_Mother_InvMass");
53 ESD_Mother_InvMass->SetName(Form("ESD_Mother_InvMass_%s",cutSelection));
54 ESD_Mother_InvMass->Write();
55
56 TH1F *ESD_Background_InvMass =fBackgroundContainer->FindObject("ESD_Background_InvMass");
57 ESD_Background_InvMass->SetName(Form("ESD_Background_InvMass_%s",cutSelection));
58 ESD_Background_InvMass->Write();
59
60 TH1F *ESD_dedx =fESDContainer->FindObject("ESD_ConvGamma_E_dEdxP");
61 ESD_dedx->SetName(Form("ESD_ConvGamma_E_dEdxP_%s",cutSelection));
62 ESD_dedx->Write();
63
64 TH2F *ESD_alfaqt = fESDContainer->FindObject("ESD_ConvGamma_alfa_qt");
65 ESD_alfaqt->SetName(Form("ESD_ConvGamma_alfa_qt_%s",cutSelection));
66 ESD_alfaqt->Write();
67
68
baf06eb6 69 TH2F *ESD_Mother_InvMass_vs_Pt = fESDContainer->FindObject("ESD_Mother_InvMass_vs_Pt");
70 ESD_Mother_InvMass_vs_Pt->Sumw2();
71 TH2F *ESD_Background_InvMass_vs_Pt = fBackgroundContainer->FindObject("ESD_Background_InvMass_vs_Pt");
72 ESD_Background_InvMass_vs_Pt->Sumw2();
73
74 //histos for comparison
75 TH1F *histoNormYield_Pi0 = new TH1F(Form("Norm_Yield_Pi0_%s",cutSelection),"",NBinsPt,BinsPt);
76 TH1F *histoRawYield_Pi0 = new TH1F(Form("Raw_Yield_Pi0_%s",cutSelection),"",NBinsPt,BinsPt);
77 TH1F *histoMass_Pi0 = new TH1F(Form("Mass_Pi0_%s",cutSelection),"",NBinsPt,BinsPt);
78 TH1F *histoFWHM_Pi0 = new TH1F(Form("FWHM_Pi0_%s",cutSelection),"",NBinsPt,BinsPt);
79 TH1F *histoSignificance_Pi0 = new TH1F(Form("Significance_Pi0_%s",cutSelection),"",NBinsPt,BinsPt);
80 TH1F *histoSB_Pi0 = new TH1F(Form("SB_Pi0_%s",cutSelection),"",NBinsPt,BinsPt);
81 TH1F *histoSBRatio_Pi0 = new TH1F(Form("SBRatio_Pi0_%s",cutSelection),"",NBinsPt,BinsPt);
82
83 // for yield correction in pt bins
84 TH1F *deltaPt = new TH1F("deltaPt","",NBinsPt,BinsPt);
85
86 TH1F * ESD_NumberOfContributorsVtx = fESDContainer->FindObject("ESD_NumberOfContributorsVtx");
87 ESD_NumberOfContributorsVtx->SetAxisRange(1.,100.);
88 Float_t nGoodEvents = ESD_NumberOfContributorsVtx->Integral();
89 cout <<"# good events: "<< nGoodEvents << endl;
90
91
92 // binning / max pt from analysis
93 Double_t maxPt = ESD_Mother_InvMass_vs_Pt->GetYaxis()->GetBinCenter(ESD_Mother_InvMass_vs_Pt->GetNbinsY()) + ESD_Mother_InvMass_vs_Pt->GetYaxis()->GetBinWidth(ESD_Mother_InvMass_vs_Pt->GetNbinsY())/2.;
94 Double_t binningPt = ESD_Mother_InvMass_vs_Pt->GetNbinsY();
95
96
97 // ---------------------------------------------------------------------------------------------
98
99 // extract yield, mass, significance and FWHM in pt bins
100
101 TH1D* Mapping_Reco_InvMass_PtBin[32]; // aray of histos for pt slices
102 TH1D* Mapping_Back_InvMass_PtBin[32];
103 TH1D* Mapping_Signal_InvMass_PtBin[32];
104 TH1D* Mapping_SBRatio_InvMass_PtBin[32];
105
106 TCanvas *canvas = new TCanvas("canvas","",200,10,600,600); // gives the page size
107 canvas->SetFillColor(0);
108
109 for(Int_t iPt = 2; iPt < NBinsPt+1; iPt++){
110 TString histonameReco = Form("Mapping_Reco_InvMass_in_Pt_Bin%s%02d",cutSelection ,iPt);
111 TString histonameBack = Form("Mapping_Back_InvMass_in_Pt_Bin%s%02d",cutSelection, iPt);
112
113 Mapping_Reco_InvMass_PtBin[iPt]=new TH1D(histonameReco.Data(),histonameReco.Data(),ESD_Mother_InvMass_vs_Pt->GetNbinsX(),0.,1.);
114 Mapping_Back_InvMass_PtBin[iPt]=new TH1D(histonameBack.Data(),histonameBack.Data(),ESD_Background_InvMass_vs_Pt->GetNbinsX(),0.,1.);
115
116 Int_t startBin = ESD_Mother_InvMass_vs_Pt->GetYaxis()->FindBin(BinsPt[iPt-1])-1;
117 Int_t endBin = ESD_Mother_InvMass_vs_Pt->GetYaxis()->FindBin(BinsPt[iPt])-1;
118
119 Double_t startpt = startBin*maxPt/binningPt;
120 Double_t endpt = endBin*maxPt/binningPt;
121
122 ESD_Mother_InvMass_vs_Pt->ProjectionX(histonameReco.Data(),startBin,endBin);
123 ESD_Background_InvMass_vs_Pt->ProjectionX(histonameBack.Data(),startBin,endBin);
124
125 Mapping_Reco_InvMass_PtBin[iPt]=(TH1D*)gDirectory->Get(histonameReco.Data());
126 Mapping_Reco_InvMass_PtBin[iPt]->Rebin(RebinMass);
127 Mapping_Back_InvMass_PtBin[iPt]=(TH1D*)gDirectory->Get(histonameBack.Data());
128 Mapping_Back_InvMass_PtBin[iPt]->Rebin(RebinMass);
129
130
131 // normalisation of background
132 TAxis *xaxis_reco = Mapping_Reco_InvMass_PtBin[iPt]->GetXaxis();
133 Int_t r_1 = xaxis_reco->FindBin(BGFit_range[0]);
134 Int_t r_2 = xaxis_reco->FindBin(BGFit_range[1]);
135 Double_t r = Mapping_Reco_InvMass_PtBin[iPt]->Integral(r_1,r_2); // Integral(75,125)
136 TAxis *xaxis_back = Mapping_Back_InvMass_PtBin[iPt]->GetXaxis();
137 Int_t b_1 = xaxis_back->FindBin(BGFit_range[0]);
138 Int_t b_2 = xaxis_back->FindBin(BGFit_range[1]);
139 Double_t b = Mapping_Back_InvMass_PtBin[iPt]->Integral(b_1,b_2);
140 Double_t norm = 1;
141 if(b != 0) norm = r/b;
142 Mapping_Back_InvMass_PtBin[iPt]->Sumw2();
143 Mapping_Back_InvMass_PtBin[iPt]->Scale(norm);
144
145 TString histonameSignal = Form("Mapping_Signal_InvMass_in_Pt_Bin%s%02d",cutSelection, iPt);
146 Mapping_Signal_InvMass_PtBin[iPt] = (TH1D*)Mapping_Reco_InvMass_PtBin[iPt]->Clone();
147 Mapping_Signal_InvMass_PtBin[iPt]->SetName(histonameSignal.Data());
148 Mapping_Signal_InvMass_PtBin[iPt]->SetTitle(histonameSignal.Data());
149 Mapping_Signal_InvMass_PtBin[iPt]->Sumw2();
150 Mapping_Signal_InvMass_PtBin[iPt]->Add(Mapping_Back_InvMass_PtBin[iPt],-1.);
151
152 /* TString histonameSBRatio= Form("Mapping_SBRatio_InvMass_in_Pt_Bin%02d", iPt);
153 Mapping_SBRatio_InvMass_PtBin[iPt]=->ProjectionX(histonameSBRatio.Data(),startBin,endBin);;
154 Mapping_SBRatio_InvMass_PtBin[iPt]->Add(Mapping_Reco_InvMass_PtBin[iPt],1);
155 Mapping_SBRatio_InvMass_PtBin[iPt]->Divide(Mapping_Back_InvMass_PtBin[iPt]);
156 Mapping_SBRatio_InvMass_PtBin[iPt]->Write();
157 */
158 Int_t bins = Mapping_Signal_InvMass_PtBin[iPt]->GetNbinsX();
159 for(Int_t ib = 1; ib < bins; ib++){
160 Double_t error = Mapping_Signal_InvMass_PtBin[iPt]->GetBinError(ib);
161 if(abs(error) < 0.1) error = Mapping_Signal_InvMass_PtBin[iPt]->GetBinError(ib) +1;
162 Mapping_Signal_InvMass_PtBin[iPt]->SetBinError(ib,error);
163 }
164
165 Double_t Pi0Amplitude = Mapping_Signal_InvMass_PtBin[iPt]->GetMaximum();
166 Double_t Pi0Amplitude_Min = Pi0Amplitude*80/100;
167 Double_t Pi0Amplitude_Max = Pi0Amplitude*100/100;
168
169 TAxis *xaxis = Mapping_Signal_InvMass_PtBin[iPt]->GetXaxis();
170 Int_t ilowmassPi0_bin = xaxis->FindBin(Pi0FitRange[0]);
171 Int_t ihighmassPi0_bin =xaxis->FindBin(Pi0FitRange[1]);
172 if(makeMappingPlots){
173 Mapping_Signal_InvMass_PtBin[iPt]->Write();
174 }
175 TF1 *fReco = new TF1(Form("Exp+Gauss PtBin[%02d]",iPt),"(x<[1])*([0]*(exp(-0.5*((x-[1])/[2])^2)+exp((x-[1])/[3])*(1-exp(-0.5*((x-[1])/[2])^2))))+(x>=[1])*([0]*exp(-0.5*((x-[1])/[2])^2))",Pi0FitRange[0],Pi0FitRange[1]);
176
177 // fit the peak around pi0 mass and draw
178 fReco->SetParameter(0,Pi0Amplitude);
179 fReco->SetParameter(1,Pi0Mass);
180 fReco->SetParameter(2,Pi0Width);
181 fReco->SetParameter(3,Pi0Slope);
182
183 fReco->SetParLimits(0,Pi0Amplitude_Min,Pi0Amplitude_Max);
184 fReco->SetParLimits(1,Pi0Mass_range[0],Pi0Mass_range[1]);
185 fReco->SetParLimits(2,Pi0Width_range[0],Pi0Width_range[1]);
186 fReco->SetParLimits(3,Pi0Slope_range[0],Pi0Slope_range[1]);
187
188 Mapping_Signal_InvMass_PtBin[iPt]->Fit(fReco,"RME");
189 // Mapping_Signal_InvMass_PtBin[iPt]->Draw();
190 fReco->SetLineColor(3);
191 fReco->SetLineWidth(0.7);
192 fReco->SetLineStyle(4);
193 fReco->Draw("same");
194 if(makeMappingPlots){
195 fReco->Write();
196 }
197 cout<<"Kenneth: Setting Bin content for bin: "<<iPt<<" to "<<fReco->GetParameter(1)<<endl;
198 histoMass_Pi0->SetBinContent(iPt, fReco->GetParameter(1));
199 cout<<"Kenneth: Setting Bin content for bin: "<<iPt<<" to "<<fReco->GetParError(1)<<endl;
200 histoMass_Pi0->SetBinError(iPt, fReco->GetParError(1));
201
202 TAxis *xaxis = Mapping_Signal_InvMass_PtBin[iPt]->GetXaxis();
203 Int_t ilowmassPi0_bin = xaxis->FindBin(Pi0FitRange[0]);
204 Int_t ihighmassPi0_bin =xaxis->FindBin(Pi0FitRange[1]);
205
206 Double_t Reco_error;
207 Double_t Background_error;
208
209 Double_t IntReco = Mapping_Reco_InvMass_PtBin[iPt]->IntegralAndError(ilowmassPi0_bin,ihighmassPi0_bin,Reco_error);
210 Double_t IntBackground = Mapping_Back_InvMass_PtBin[iPt]->IntegralAndError(ilowmassPi0_bin,ihighmassPi0_bin,Background_error);
211 if(IntBackground>0){
212 Double_t SB= IntReco/IntBackground;
213 Double_t SB_error = sqrt( pow(Reco_error/IntBackground,2) + pow(IntReco/IntBackground/IntBackground*Background_error,2));
214
215 histoSB_Pi0->SetBinContent(iPt,SB);
216 histoSB_Pi0->SetBinError(iPt, SB_error);
217 }
218 // canvas->Update();
219 // canvas->Clear();
220
221 Double_t parameter[4];
222 fReco->GetParameters(parameter);
223
224 // calculation of FWHM + error
225 Int_t Pi0Amplitude_fit = fReco->GetMaximum(Pi0FitRange[0],Pi0FitRange[1]);
226 Int_t Pi0Amplitude_bin = fReco->GetXaxis()->FindBin(Pi0Amplitude_fit);
227 Double_t MassFWHM1 = fReco->GetX(Pi0Amplitude_fit*0.5,0.1,0.135);
228 Double_t MassFWHM2 = fReco->GetX(Pi0Amplitude_fit*0.5,0.135,0.16);
229 Double_t FWHM = MassFWHM2 - MassFWHM1;
230
231 // + error
232 TF1 *fReco_EPlus = new TF1("Exp+Gauss","(x<[1])*([0]*(exp(-0.5*((x-[1])/[2])^2)+exp((x-[1])/[3])*(1-exp(-0.5*((x-[1])/[2])^2))))+(x>=[1])*([0]*exp(-0.5*((x-[1])/[2])^2))",Pi0FitRange[0],Pi0FitRange[1]);
233 fReco_EPlus->SetParameter(0,parameter[0] + fReco->GetParError(0));
234 fReco_EPlus->SetParameter(1,parameter[1] + fReco->GetParError(1));
235 fReco_EPlus->SetParameter(2,parameter[2] + fReco->GetParError(2));
236 fReco_EPlus->SetParameter(3,parameter[3] + fReco->GetParError(3));
237
238 Int_t Pi0Amplitude_EPlus_fit = fReco_EPlus->GetMaximum(Pi0FitRange[0],Pi0FitRange[1]);
239 Int_t Pi0Amplitude_EPlus_bin = fReco_EPlus->GetXaxis()->FindBin(Pi0Amplitude_EPlus_fit);
240
241 Double_t MassFWHM1_EPlus = fReco_EPlus->GetX(Pi0Amplitude_EPlus_fit*0.5,0.1,0.135);
242 Double_t MassFWHM2_EPlus = fReco_EPlus->GetX(Pi0Amplitude_EPlus_fit*0.5,0.135,0.16);
243 Double_t FWHM_EPlus = MassFWHM2_EPlus - MassFWHM1_EPlus;
244
245 // - error
246 TF1 *fReco_EMinus = new TF1("Exp+Gauss","(x<[1])*([0]*(exp(-0.5*((x-[1])/[2])^2)+exp((x-[1])/[3])*(1-exp(-0.5*((x-[1])/[2])^2))))+(x>=[1])*([0]*exp(-0.5*((x-[1])/[2])^2))",Pi0FitRange[0],Pi0FitRange[1]);
247 fReco_EMinus->SetParameter(0,parameter[0] - fReco->GetParError(0));
248 fReco_EMinus->SetParameter(1,parameter[1] - fReco->GetParError(1));
249 fReco_EMinus->SetParameter(2,parameter[2] - fReco->GetParError(2));
250 fReco_EMinus->SetParameter(3,parameter[3] - fReco->GetParError(3));
251
252 Int_t Pi0Amplitude_EMinus_fit = fReco_EMinus->GetMaximum(Pi0FitRange[0],Pi0FitRange[1]);
253 Int_t Pi0Amplitude_EMinus_bin = fReco_EMinus->GetXaxis()->FindBin(Pi0Amplitude_EMinus_fit);
254 Double_t MassFWHM1_EMinus = fReco_EMinus->GetX(Pi0Amplitude_EMinus_fit*0.5,0.1,0.135);
255 Double_t MassFWHM2_EMinus = fReco_EMinus->GetX(Pi0Amplitude_EMinus_fit*0.5,0.135,0.16);
256 Double_t FWHM_EMinus = MassFWHM2_EMinus - MassFWHM1_EMinus;
257
258 Double_t Error1 = TMath::Abs(FWHM-FWHM_EPlus);
259 Double_t Error2 = TMath::Abs(FWHM-FWHM_EMinus);
260
261 Double_t FWHM_error;
262 if(Error1>=Error2) FWHM_error = Error1;
263 if(Error1 < Error2) FWHM_error = Error2;
264
265 histoFWHM_Pi0->SetBinContent(iPt, FWHM);
266 histoFWHM_Pi0->SetBinError(iPt,FWHM_error);
267
268 Double_t Reco_error;
269 Double_t Background_error;
270
271 Double_t Reco = Mapping_Reco_InvMass_PtBin[iPt]->IntegralAndError(ilowmassPi0_bin,ihighmassPi0_bin,Reco_error);
272 Double_t Background = Mapping_Back_InvMass_PtBin[iPt]->IntegralAndError(ilowmassPi0_bin,ihighmassPi0_bin,Background_error);
273
274 if(Background>0){
275 Double_t Significance = Reco/sqrt(Background);
276 Double_t Significance_error = sqrt( pow(Reco_error/sqrt(Background),2) + pow( -0.5 * Reco/pow(Background,1.5) * Background_error,2) );
277
278 histoSignificance_Pi0->SetBinContent(iPt,Significance);
279 histoSignificance_Pi0->SetBinError(iPt, Significance_error);
280 }
281 histoNormYield_Pi0->SetBinContent(iPt,(Reco-Background)/nGoodEvents);
282 histoNormYield_Pi0->SetBinError(iPt,sqrt(Reco_error*Reco_error + Background_error*Background_error)/nGoodEvents);
283
284 histoRawYield_Pi0->SetBinContent(iPt,Reco-Background);
285 histoRawYield_Pi0->SetBinError(iPt,sqrt(Reco_error*Reco_error + Background_error*Background_error));
286
287
288 deltaPt->SetBinContent(iPt, endpt-startpt);
289 deltaPt->SetBinError(iPt,0);
290 if(iPt == 2){
291 deltaPt->SetBinContent(1,startpt-0);
292 deltaPt->SetBinError(1,0);
293 }
294 }
295
296 delete [] Mapping_Reco_InvMass_PtBin;
297 delete [] Mapping_Back_InvMass_PtBin;
298 delete [] Mapping_Signal_InvMass_PtBin;
299 delete [] Mapping_SBRatio_InvMass_PtBin;
300 // end of extract yield, mass .....
301
302 // write the created histos into the output root file
303
304 if(makeMappingPlots == kFALSE){// open the root file here so the mapping histograms is not included
305 outputFile = new TFile(Form("%sPi0Characteristics.root",outputDir),"UPDATE");
306 }
307
308
309 // ESD_Mother_InvMass_vs_Pt->ProjectionX(histonameReco.Data(),startBin,endBin);
310 // ESD_Background_InvMass_vs_Pt->ProjectionX(histonameBack.Data(),startBin,endBin);
311 /*
312 histoSBRatio_Pi0->Add(ESD_Mother_InvMass_vs_Pt,1);
313 histoSBRatio_Pi0->Divide(ESD_Background_InvMass_vs_Pt);
314 histoSBRatio_Pi0->Write();
315 */
316 histoNormYield_Pi0->SetXTitle("p_{t} (GeV/c)");
317 histoNormYield_Pi0->Divide(deltaPt);
318 histoNormYield_Pi0->SetYTitle("#pi^{0} raw yield/per event");
319 histoNormYield_Pi0->DrawCopy("e1");
320 histoNormYield_Pi0->Write();
321 // outputFile->Write();
322 canvas->Update();
323 canvas->Clear();
324
325 histoRawYield_Pi0->SetXTitle("p_{t} (GeV/c)");
326 histoRawYield_Pi0->Divide(deltaPt);
327 histoRawYield_Pi0->SetYTitle("#pi^{0} raw yield");
328 histoRawYield_Pi0->DrawCopy("e1");
329 histoRawYield_Pi0->Write();
330 // outputFile->Write();
331 canvas->Update();
332 canvas->Clear();
333
334 histoMass_Pi0->SetXTitle("p_{t} (GeV/c)");
335 // histoMass_Pi0->Divide(deltaPt);
336 histoMass_Pi0->SetYTitle("#pi^{0} mass");
337 histoMass_Pi0->DrawCopy("e1");
338 histoMass_Pi0->Write();
339 canvas->Update();
340 canvas->Clear();
341
342
343 histoSB_Pi0->SetXTitle("p_{t} (GeV/c)");
344 // histoMass_Pi0->Divide(deltaPt);
345 histoSB_Pi0->SetYTitle("S/B");
346 histoSB_Pi0->DrawCopy("e1");
347 histoSB_Pi0->Write();
348 canvas->Update();
349 canvas->Clear();
350
351
352 histoFWHM_Pi0->SetXTitle("p_{t} (GeV/c)");
353 // histoFWHM_Pi0->Divide(deltaPt);
354 histoFWHM_Pi0->SetYTitle("#pi^{0} FWHM/2.36");
355 histoFWHM_Pi0->DrawCopy("e1");
356 histoFWHM_Pi0->Write();
357 canvas->Update();
358 canvas->Clear();
359
360 histoSignificance_Pi0->SetXTitle("p_{t} (GeV/c)");
361 // histoSignificance_Pi0->Divide(deltaPt);
362 histoSignificance_Pi0->SetYTitle("Significance");
363 histoSignificance_Pi0->DrawCopy("e1");
364 histoSignificance_Pi0->Write();
365 canvas->Update();
366 canvas->Clear();
367
368 histoSBRatio_Pi0->Delete();
369 histoNormYield_Pi0->Delete();
370 histoRawYield_Pi0->Delete();
371 histoMass_Pi0->Delete();
372 histoFWHM_Pi0->Delete();
373 histoSignificance_Pi0->Delete();
374 histoSB_Pi0->Delete();
375 deltaPt->Delete();
376
377 canvas->Delete();
378
379 outputFile->Write();
380 outputFile->Close();
381}