]>
Commit | Line | Data |
---|---|---|
a65a7e70 | 1 | void CorrectForEfficiencypptrd(const char *filedata,const char *fileMC,const char *NameList,Bool_t withtrd=kTRUE, Float_t CorrType=0, Int_t conf=2, Bool_t smoothing=kFALSE, Bool_t hadroncontaminationsubtracted=kTRUE, Bool_t UnsetCorrelatedErrors=kFALSE); |
2 | void CorrectForEfficiencyBeauty(const char *filedata,const char *fileMC, const char *fileMCbg); | |
3 | TList *GetResults(const Char_t *testfile,const Char_t *plus=""); | |
4 | TList *GetQA(const Char_t *testfile,const Char_t *plus=""); | |
5 | TObject* GetSpectrum(AliCFContainer *c, Int_t step); | |
6 | THnSparseF* GetHadronEffbyIPcut(TList *hfeqa); | |
7 | void CorrectFromTheWidth(TH1D *h1); | |
8 | ||
9 | void CorrectForEfficiencypptrd(const char *filedata,const char *fileMC,const char *NameList,Bool_t withtrd, Float_t CorrType, Int_t conf, Bool_t smoothing, Bool_t hadroncontaminationsubtracted, Bool_t UnsetCorrelatedErrors) { | |
10 | ||
11 | /////////////////////////////////////////////////////////////////////////////////////////////////////// | |
12 | // Will produce finalSpectrum.root with results inside | |
13 | // TGraphErrors UnfoldingCorrectedSpectrum -> unfolding procedure (doesn't work in 2D with trunk) | |
14 | // TGraphErrors AlltogetherSpectrum -> corrected spectrum with other method | |
15 | // TH1D RatioUnfoldingAlltogetherSpectrum | |
16 | // THnSparseF UnfoldingCorrectedNotNormalizedSpectrum -> unfolding procedure not yet normalized | |
17 | // AliCFDataGrid AlltogetherCorrectedNotNormalizedSpectrum -> other method not yet normalized | |
18 | // | |
19 | // NameList - possible extension for the full name of the data and MC TLists (e.g. "_PID2" for HFE_Results_PID2) | |
20 | // CorrType - 0 use fixed TPC PID 0.5 | |
21 | // 1 use TPC PID efficiency from V0s | |
22 | // | |
23 | // Conf - TRD efficiency | |
24 | // double TRDeffs[5] = {0.7,0.8,0.9,0.75,0.85}; | |
25 | // conf = 2, TRDeff = 0.8... | |
26 | // (Only relevant for withtrd = kTRUE) | |
27 | // | |
28 | // f/ smoothing kFALSE | |
29 | // | |
30 | // g/ Flag for hadron contamination subtraction | |
31 | // kTRUE means we subtract the hadron contamination | |
32 | // kFALSE not | |
33 | // | |
34 | /////////////////////////////////////////////////////////////////////////////////////////////////////// | |
35 | ||
36 | double TRDeffs[5] = {0.7,0.8,0.9,0.75,0.85}; | |
37 | ||
38 | gStyle->SetPalette(1); | |
39 | gStyle->SetOptStat(1111); | |
40 | gStyle->SetPadBorderMode(0); | |
41 | gStyle->SetCanvasColor(10); | |
42 | gStyle->SetPadLeftMargin(0.13); | |
43 | gStyle->SetPadRightMargin(0.13); | |
44 | ||
45 | if (CorrType>1){ | |
46 | printf("Wrong argument for the type of correction, it is maybe missing, or it has to be <=2!!! \n"); | |
47 | return; | |
48 | } | |
49 | ||
50 | ///////////////////// | |
51 | // Take the stuff | |
52 | ///////////////////// | |
53 | ||
54 | printf("Get results for data \n"); | |
55 | TList *resultsdata = GetResults(filedata,NameList); | |
56 | if(!resultsdata){ | |
57 | printf("No output objects for data: Calculation will terminate here\n"); | |
58 | return; | |
59 | } | |
60 | ||
61 | /////////////////////////////////////////////////////////// | |
62 | // Event container for the normalization to CINB1 | |
63 | // Take the number of events as function of z and number of contributors | |
64 | // For those without primary vertex from tracks, take the fraction in the z range from MC | |
65 | // For 10cm: LHC10f6a (88.0756%), LHC10f6 (86.7461%) | |
66 | //////////////////////////////////////////////////////////// | |
67 | resultsdata->Print(); | |
68 | ||
69 | AliCFContainer *containerdata = (AliCFContainer *) resultsdata->FindObject("eventContainer"); | |
70 | if(!containerdata) { | |
71 | printf("No container \n"); | |
72 | return; | |
73 | } | |
74 | AliCFDataGrid *dataGrid = (AliCFDataGrid *) GetSpectrum(containerdata,AliHFEcuts::kEventStepZRange); | |
75 | //dataGrid->PrintBinLimits(); | |
76 | dataGrid->PrintNBins(); | |
77 | ||
78 | // | |
79 | TH2D *spectrum_zrange = (TH2D *) dataGrid->Project(0,3); | |
80 | TAxis *axis_0 = spectrum_zrange->GetXaxis(); | |
81 | TAxis *axis_1 = spectrum_zrange->GetYaxis(); | |
82 | printf("Number of bins in x %d and in y %d\n",axis_0->GetNbins(),axis_1->GetNbins()); | |
83 | TH1D *spectrum1Da = (TH1D *) spectrum_zrange->ProjectionX("bin0",1,1); | |
84 | TH1D *spectrum1Db = (TH1D *) spectrum_zrange->ProjectionX("bin>0",2,12); | |
85 | TH1D *spectrum1Dc = (TH1D *) spectrum_zrange->ProjectionX("total"); | |
86 | Double_t nbinbin0 = spectrum1Da->Integral(); | |
87 | Double_t nbinnobin0 = spectrum1Db->Integral(); | |
88 | Double_t nbintotal = spectrum1Dc->Integral(); | |
89 | ||
90 | //TCanvas *c1teste = new TCanvas("eventcontainertest","eventcontainertest",800,800); | |
91 | //c1teste->cd(1); | |
92 | //spectrum_zrange->Draw("colz"); | |
93 | ||
94 | ||
95 | Float_t numberofentries = nbinnobin0+nbinbin0*0.880756; | |
96 | printf("Number in bin0 %f, number out %f, total %f, normalized %f\n",nbinbin0,nbinnobin0,nbintotal,numberofentries); | |
97 | ||
98 | ////////////////////////////// | |
99 | // Take more stuff | |
100 | /////////////////////////////// | |
101 | ||
102 | printf("Get results for MC \n"); | |
103 | TList *resultsmc = GetResults(fileMC,NameList); | |
104 | if(!resultsmc){ | |
105 | printf("No output objects for mc: Calculation will terminate here\n"); | |
106 | return; | |
107 | } | |
108 | ||
109 | AliHFEcontainer *datahfecontainer = (AliHFEcontainer *) resultsdata->FindObject("trackContainer"); | |
110 | if(!datahfecontainer) { | |
111 | printf("No container for data \n"); | |
112 | return; | |
113 | } | |
114 | AliCFContainer *sumcontainer = datahfecontainer->GetCFContainer("recTrackContReco"); | |
115 | ||
116 | //////////////////////////////////////////// | |
117 | //Plot raw spectrum for TPC TOF scenario | |
118 | //////////////////////////////////////////// | |
119 | AliCFDataGrid *dataGrida = 0x0; | |
120 | if(!withtrd) dataGrida = (AliCFDataGrid *) GetSpectrum(sumcontainer,AliHFEcuts::kStepHFEcutsTRD + 2); | |
121 | else dataGrida = (AliCFDataGrid *) GetSpectrum(sumcontainer,AliHFEcuts::kStepHFEcutsTRD + 3); | |
122 | ||
123 | TH1D *spectrumpt = (TH1D *) dataGrida->Project(0); | |
124 | CorrectFromTheWidth(spectrumpt); | |
125 | ||
126 | TH1D *total = new TH1D("total","",1,0.3,10.3); | |
127 | total->SetMaximum(1.0e+07); | |
128 | total->SetMinimum(1000); | |
129 | total->SetXTitle("p_{T} [GeV/c]"); | |
130 | total->SetYTitle("dN/dp_{T}[GeV/c]^{-1}"); | |
131 | total->SetTitleOffset(1.5,"Y"); | |
132 | //total->Scale(0.9); | |
133 | total->GetXaxis()->SetRangeUser(0.38,10.3); | |
134 | ||
135 | ||
136 | TCanvas *c1test = new TCanvas("rawspectrum","rawspectrum",800,800); | |
137 | c1test->cd(1); | |
138 | gPad->SetLeftMargin(0.13); | |
139 | gPad->SetLogy(); | |
140 | gPad->SetTicks(); | |
141 | //gPad->SetGridx(); | |
142 | //gPad->SetGridy(); | |
143 | gPad->SetFillColor(10); | |
144 | gPad->SetFrameFillColor(10); | |
145 | gStyle->SetOptStat(0); | |
146 | gStyle->SetOptFit(1111); | |
147 | gStyle->SetOptTitle(0); | |
148 | total->SetStats(0); | |
149 | spectrumpt->SetStats(0); | |
150 | total->Draw(); | |
151 | spectrumpt->Draw("same"); | |
152 | TPaveText* t1=new TPaveText(0.49,0.45,0.79,0.52,"NDC"); | |
153 | t1->SetFillStyle(0); | |
154 | t1->SetBorderSize(0); | |
155 | t1->AddText(0.,0.,"pp, #sqrt{s} = 7 TeV"); | |
156 | t1->SetTextColor(kRed); | |
157 | //t1->SetTextSize(20); | |
158 | t1->SetTextFont(42); | |
159 | t1->Draw(); | |
160 | TPaveText* t2=new TPaveText(0.49,0.35,0.79,0.42,"NDC"); | |
161 | t2->SetFillStyle(0); | |
162 | t2->SetBorderSize(0); | |
163 | t2->AddText(0.,0.,"L = 2.6 nb^{-1}"); | |
164 | t2->SetTextColor(kRed); | |
165 | //t1->SetTextSize(20); | |
166 | t2->SetTextFont(42); | |
167 | t2->Draw(); | |
168 | TPaveText* t3=new TPaveText(0.49,0.35,0.79,0.42,"NDC"); | |
169 | t3->SetFillStyle(0); | |
170 | t3->SetBorderSize(0); | |
171 | t3->AddText(0.,0.,"|#eta| < 0.8"); | |
172 | t3->SetTextColor(kRed); | |
173 | //t1->SetTextSize(20); | |
174 | t3->SetTextFont(42); | |
175 | t3->Draw(); | |
176 | ||
177 | ///////////////////////////// | |
178 | // Check number of events | |
179 | ///////////////////////////// | |
180 | ||
181 | ||
182 | Int_t numberOfEventsafterallcuts = (Int_t) datahfecontainer->GetNumberOfEvents(); | |
183 | ||
184 | printf("Number of events not corrected %d\n",numberOfEventsafterallcuts); | |
185 | printf("Number of events corrected %f\n",numberofentries); | |
186 | ||
187 | AliHFEcontainer *mchfecontainer = (AliHFEcontainer *) resultsmc->FindObject("trackContainer"); | |
188 | if(!mchfecontainer) { | |
189 | printf("No mc container \n"); | |
190 | return; | |
191 | } | |
192 | ||
193 | printf("Find the container V0\n"); | |
194 | AliHFEcontainer *containerhfeV0 = (AliHFEcontainer *) resultsdata->FindObject("containerV0"); | |
195 | if(!containerhfeV0) { | |
196 | printf("No hfe container \n"); | |
197 | return; | |
198 | } | |
199 | ||
200 | ////////////// | |
201 | // Correct | |
202 | ///////////// | |
203 | ||
204 | AliHFEspectrum *spectrum = new AliHFEspectrum("HFEspectrum"); | |
205 | spectrum->SetNbDimensions(1); | |
206 | // If you want to correct positive (0) or negative (1) separately | |
207 | //spectrum->SetChargeChoosen(0); | |
208 | spectrum->SetSmoothing(smoothing); | |
209 | spectrum->SetUnSetCorrelatedErrors(UnsetCorrelatedErrors); | |
210 | spectrum->SetDebugLevel(1); | |
211 | spectrum->SetNumberOfEvents((Int_t)numberofentries); | |
212 | spectrum->SetDumpToFile(kTRUE); | |
213 | // True step in MC (events in range +- 10 cm and no pile up) | |
214 | spectrum->SetMCTruthStep(AliHFEcuts::kStepMCGeneratedZOutNoPileUpCentralityFine); | |
215 | ||
216 | // Step where we correct from MC (tracking + TOF) | |
217 | spectrum->SetMCEffStep(AliHFEcuts::kNcutStepsMCTrack + AliHFEcuts::kStepHFEcutsTRD + 1); | |
218 | printf("!!! Step we correct with MC %d\n!!!",AliHFEcuts::kNcutStepsMCTrack + AliHFEcuts::kStepHFEcutsTRD + 1); | |
219 | ||
220 | ||
221 | // Step from data we correct for | |
222 | //////////****************************////////////////////////////// | |
223 | if(withtrd) spectrum->SetStepToCorrect(AliHFEcuts::kStepHFEcutsTRD + 3); | |
224 | else spectrum->SetStepToCorrect(AliHFEcuts::kStepHFEcutsTRD + 2); | |
225 | //////////****************************////////////////////////////// | |
226 | ||
227 | // Steps to be corrected with V0 | |
228 | if (CorrType == 1){ | |
229 | printf("\n This time we correct with PID efficiency from V0s !!! \n"); | |
230 | if(withtrd) { | |
231 | spectrum->SetStepBeforeCutsV0(AliHFEcuts::kStepHFEcutsTRD + 2); | |
232 | spectrum->SetStepAfterCutsV0(AliHFEcuts::kStepHFEcutsTRD + 3); | |
233 | } | |
234 | else { | |
235 | spectrum->SetStepBeforeCutsV0(AliHFEcuts::kStepHFEcutsTRD + 1); | |
236 | spectrum->SetStepAfterCutsV0(AliHFEcuts::kStepHFEcutsTRD + 2); | |
237 | } | |
238 | } | |
239 | ||
240 | if (CorrType < 1){ | |
241 | ||
242 | Double_t effTPC = 0.5, | |
243 | effTRD = TRDeffs[conf-1]; | |
244 | printf("effTRD = %f\n", effTRD); | |
245 | // Place for a correction as an efficiency parametrized as function of pt of pt,eta,phi or p=pt | |
246 | TF1 *hEfficiency = new TF1("efficiency", "[0]", 0., 20.); | |
247 | if(withtrd) hEfficiency->SetParameter(0, effTPC*effTRD); | |
248 | else hEfficiency->SetParameter(0, effTPC); | |
249 | spectrum->SetEfficiencyFunction(hEfficiency); | |
250 | } | |
251 | ||
252 | // Give everything (data container, mc container and V0 data container) | |
253 | if (CorrType == 1){ | |
254 | printf("\n This time we correct with PID efficiency from V0s !!! \n"); | |
255 | spectrum->Init(datahfecontainer,mchfecontainer,containerhfeV0); | |
256 | } | |
257 | else { | |
258 | printf("\n This time we correct with PID efficiency from MC or fixed value !!! \n"); | |
259 | spectrum->Init(datahfecontainer,mchfecontainer); | |
260 | } | |
261 | ||
262 | // kTRUE means subtract hadron background, kFALSE means do not subtract hadron background | |
263 | spectrum->Correct(hadroncontaminationsubtracted); | |
264 | ||
265 | ||
266 | } | |
267 | ||
268 | //_____________________________________________________________________ | |
269 | void CorrectForEfficiencyBeauty(const char *filedata,const char *fileMC, const char *fileMCbg) { | |
270 | ||
271 | /////////////////////////////////////////////////////////////////////////////////////////////////////// | |
272 | // Will produce finalSpectrum.root with results inside | |
273 | // TGraphErrors UnfoldingCorrectedSpectrum -> unfolding procedure (doesn't work in 2D with trunk) | |
274 | // TGraphErrors AlltogetherSpectrum -> corrected spectrum with other method | |
275 | // TH1D RatioUnfoldingAlltogetherSpectrum | |
276 | // THnSparseF UnfoldingCorrectedNotNormalizedSpectrum -> unfolding procedure not yet normalized | |
277 | // AliCFDataGrid AlltogetherCorrectedNotNormalizedSpectrum -> other method not yet normalized | |
278 | /////////////////////////////////////////////////////////////////////////////////////////////////////// | |
279 | ||
280 | gStyle->SetPalette(1); | |
281 | gStyle->SetOptStat(1111); | |
282 | gStyle->SetPadBorderMode(0); | |
283 | gStyle->SetCanvasColor(10); | |
284 | gStyle->SetPadLeftMargin(0.13); | |
285 | gStyle->SetPadRightMargin(0.13); | |
286 | ||
287 | ///////////////////// | |
288 | // Take the stuff | |
289 | ///////////////////// | |
290 | ||
291 | TList *resultsdata = GetResults(filedata,"_PID2"); | |
292 | //TList *resultsdata = GetResults(filedata); | |
293 | if(!resultsdata){ | |
294 | printf("No output objects for data: Calculation will terminate here\n"); | |
295 | return; | |
296 | } | |
297 | /////////////////////////////////////////////////////////// | |
298 | // Event container for the normalization to CINB1 | |
299 | // Take the number of events as function of z and number of contributors | |
300 | // For those without primary vertex from tracks, take the fraction in the z range from MC | |
301 | // For 10cm: LHC10f6a (88.0756%), LHC10f6 (86.7461%) | |
302 | //////////////////////////////////////////////////////////// | |
303 | ||
304 | AliCFContainer *containerdata = (AliCFContainer *) resultsdata->FindObject("eventContainer"); | |
305 | if(!containerdata) { | |
306 | printf("No container \n"); | |
307 | return; | |
308 | } | |
309 | AliCFDataGrid *dataGrid = (AliCFDataGrid *) GetSpectrum(containerdata,AliHFEcuts::kEventStepZRange); | |
310 | TH2D *spectrum_zrange = (TH2D *) dataGrid->Project(0,2); | |
311 | TH1D *spectrum1Da = (TH1D *) spectrum_zrange->ProjectionX("bin0",1,1); | |
312 | TH1D *spectrum1Db = (TH1D *) spectrum_zrange->ProjectionX("bin>0",2,12); | |
313 | TH1D *spectrum1Dc = (TH1D *) spectrum_zrange->ProjectionX("total"); | |
314 | Double_t nbinbin0 = spectrum1Da->Integral(); | |
315 | Double_t nbinnobin0 = spectrum1Db->Integral(); | |
316 | Double_t nbintotal = spectrum1Dc->Integral(); | |
317 | ||
318 | Float_t numberofentries = nbinnobin0+nbinbin0*0.880756; | |
319 | printf("Number in bin0 %f, number out %f, total %f, normalized %f\n",nbinbin0,nbinnobin0,nbintotal,numberofentries); | |
320 | ||
321 | ||
322 | AliHFEcontainer *datahfecontainer = (AliHFEcontainer *) resultsdata->FindObject("trackContainer"); | |
323 | if(!datahfecontainer) { | |
324 | printf("No container for data \n"); | |
325 | return; | |
326 | } | |
327 | ||
328 | ||
329 | ////////////////////////////// | |
330 | // Check MC # of events | |
331 | /////////////////////////////// | |
332 | ||
333 | TList *resultsmc = GetResults(fileMC,"_PID2"); | |
334 | if(!resultsmc){ | |
335 | printf("No output objects for mc: Calculation will terminate here\n"); | |
336 | return; | |
337 | } | |
338 | ||
339 | AliCFContainer *containermc = (AliCFContainer *) resultsmc->FindObject("eventContainer"); | |
340 | if(!containermc) { | |
341 | printf("No container \n"); | |
342 | return; | |
343 | } | |
344 | AliCFDataGrid *mcGrid = (AliCFDataGrid *) GetSpectrum(containermc,AliHFEcuts::kEventStepZRange); | |
345 | TH2D *spectrum_zrangemc = (TH2D *) mcGrid->Project(0,2); | |
346 | TH1D *spectrum1Damc = (TH1D *) spectrum_zrangemc->ProjectionX("bin0",1,1); | |
347 | TH1D *spectrum1Dbmc = (TH1D *) spectrum_zrangemc->ProjectionX("bin>0",2,12); | |
348 | TH1D *spectrum1Dcmc = (TH1D *) spectrum_zrangemc->ProjectionX("total"); | |
349 | Double_t nbinbin0mc = spectrum1Damc->Integral(); | |
350 | Double_t nbinnobin0mc = spectrum1Dbmc->Integral(); | |
351 | Double_t nbintotalmc = spectrum1Dcmc->Integral(); | |
352 | ||
353 | Float_t numberofentriesmc = nbinnobin0mc+nbinbin0mc*0.880756; | |
354 | printf("MC: Number in bin0 %f, number out %f, total %f, normalized %f\n",nbinbin0mc,nbinnobin0mc,nbintotalmc,numberofentriesmc); | |
355 | ||
356 | ||
357 | ||
358 | ///////////////////////////// | |
359 | // Check number of events | |
360 | ///////////////////////////// | |
361 | ||
362 | Int_t numberOfEventsafterallcuts = (Int_t) datahfecontainer->GetNumberOfEvents(); | |
363 | ||
364 | printf("Number of events not corrected %d\n",numberOfEventsafterallcuts); | |
365 | printf("Number of events corrected %f\n",numberofentries); | |
366 | ||
367 | AliHFEcontainer *mchfecontainer = (AliHFEcontainer *) resultsmc->FindObject("trackContainer"); | |
368 | if(!mchfecontainer) { | |
369 | printf("No mc container \n"); | |
370 | return; | |
371 | } | |
372 | ||
373 | printf("Find the container V0\n"); | |
374 | AliHFEcontainer *containerhfeV0 = (AliHFEcontainer *) resultsdata->FindObject("containerV0"); | |
375 | if(!containerhfeV0) { | |
376 | printf("No hfe container \n"); | |
377 | return; | |
378 | } | |
379 | ||
380 | // nonHFE backgrounds | |
381 | TList *resultsmcbg = GetResults(fileMCbg,"_PID2"); | |
382 | if(!resultsmcbg){ | |
383 | printf("No output objects for mc: Calculation will terminate here\n"); return; | |
384 | } | |
385 | AliHFEcontainer *mcbghfecontainer = (AliHFEcontainer *) resultsmcbg->FindObject("trackContainer"); | |
386 | if(!mcbghfecontainer) { | |
387 | printf("No mc container \n"); | |
388 | return; | |
389 | } | |
390 | ||
391 | AliCFContainer *containermcbg = (AliCFContainer *) resultsmcbg->FindObject("eventContainer"); | |
392 | if(!containermcbg) { | |
393 | printf("No container \n"); | |
394 | return; | |
395 | } | |
396 | ||
397 | AliCFDataGrid *mcbgGrid = (AliCFDataGrid *) GetSpectrum(containermcbg,AliHFEcuts::kEventStepZRange); | |
398 | TH2D *spectrum_zrangemcbg = (TH2D *) mcbgGrid->Project(0,2); | |
399 | TH1D *spectrum1Damcbg = (TH1D *) spectrum_zrangemcbg->ProjectionX("bin0",1,1); | |
400 | TH1D *spectrum1Dbmcbg = (TH1D *) spectrum_zrangemcbg->ProjectionX("bin>0",2,12); | |
401 | TH1D *spectrum1Dcmcbg = (TH1D *) spectrum_zrangemcbg->ProjectionX("total"); | |
402 | Double_t nbinbin0mcbg = spectrum1Damcbg->Integral(); | |
403 | Double_t nbinnobin0mcbg = spectrum1Dbmcbg->Integral(); | |
404 | Double_t nbintotalmcbg = spectrum1Dcmcbg->Integral(); | |
405 | ||
406 | Float_t numberofentriesmcbg = nbinnobin0mcbg+nbinbin0mcbg*0.880756; | |
407 | printf("MCbg: Number in bin0 %f, number out %f, total %f, normalized %f\n",nbinbin0mcbg,nbinnobin0mcbg,nbintotalmcbg,numberofentriesmcbg); | |
408 | ||
409 | ||
410 | // hadron contamination after IP cuts | |
411 | TList *hfeqa = GetQA(fileMC,"_PID2"); | |
412 | THnSparseF* hsHadronEffbyIPcut = (THnSparseF* )GetHadronEffbyIPcut(hfeqa); | |
413 | ||
414 | ////////////// | |
415 | // Correct | |
416 | ///////////// | |
417 | ||
418 | AliHFEspectrum *spectrum = new AliHFEspectrum("HFEspectrum"); | |
419 | ||
420 | spectrum->SetBeautyAnalysis(); | |
421 | // Enable background subtraction | |
422 | spectrum->EnableIPanaHadronBgSubtract(); | |
423 | spectrum->EnableIPanaCharmBgSubtract(); | |
424 | spectrum->EnableIPanaConversionBgSubtract(); | |
425 | spectrum->EnableIPanaNonHFEBgSubtract(); | |
426 | //spectrum->SetNonHFEBackground2ndMethod(); | |
427 | ||
428 | spectrum->SetNbDimensions(1); | |
429 | // If you want to correct positive (0) or negative (1) separately | |
430 | //spectrum->SetChargeChoosen(0); | |
431 | spectrum->SetDebugLevel(1); | |
432 | spectrum->SetNumberOfEvents((Int_t)numberofentries); | |
433 | spectrum->SetNumberOfMCEvents((Int_t)numberofentriesmc); | |
434 | spectrum->SetNumberOfMC2Events((Int_t)numberofentriesmcbg); | |
435 | spectrum->SetDumpToFile(smoothing); | |
436 | // True step in MC (events in range +- 10 cm and no pile up) | |
437 | spectrum->SetMCTruthStep(AliHFEcuts::kStepMCGeneratedZOutNoPileUpCentralityFine); | |
438 | // Step where we correct from MC (tracking + TOF) | |
439 | spectrum->SetMCEffStep(AliHFEcuts::kNcutStepsMCTrack + AliHFEcuts::kStepHFEcutsTRD + 1); | |
440 | // Step from data we correct for | |
441 | spectrum->SetStepToCorrect(AliHFEcuts::kStepHFEcutsTRD + 3); // +3 = tracking + TOF + TPC + IPcut | |
442 | // Steps to be corrected with V0 (TPC PID only) | |
443 | spectrum->SetStepBeforeCutsV0(AliHFEcuts::kStepHFEcutsTRD + 1); | |
444 | spectrum->SetStepAfterCutsV0(AliHFEcuts::kStepHFEcutsTRD + 2); | |
445 | ||
446 | spectrum->SetHadronEffbyIPcut(hsHadronEffbyIPcut); | |
447 | // Give everything (data container, mc container and V0 data container) | |
448 | spectrum->Init(datahfecontainer,mchfecontainer,containerhfeV0,mcbghfecontainer); | |
449 | ||
450 | // kTRUE means subtract hadron background, kFALSE means do not subtract hadron background | |
451 | spectrum->CorrectBeauty(kTRUE); | |
452 | ||
453 | } | |
454 | ||
455 | //_____________________________________________________________________ | |
456 | TList *GetResults(const Char_t *testfile,const Char_t *plus){ | |
457 | // | |
458 | // read output | |
459 | // | |
460 | TFile *f = TFile::Open(testfile); | |
461 | if(!f || f->IsZombie()){ | |
462 | printf("File not readable\n"); | |
463 | return 0x0; | |
464 | } | |
465 | TString name(plus); | |
466 | printf("Name of TList %s\n",(const char*)name); | |
467 | TList *l = dynamic_cast<TList *>(f->Get((const char*)name)); | |
468 | if(!l){ | |
469 | printf("Results list was not found\n"); | |
470 | f->Close(); delete f; | |
471 | return 0x0; | |
472 | } | |
473 | TList *returnlist = dynamic_cast<TList *>(l->Clone()); | |
474 | f->Close(); delete f; | |
475 | return returnlist; | |
476 | } | |
477 | ||
478 | //_____________________________________________________________________ | |
479 | TList *GetQA(const Char_t *testfile,const Char_t *plus){ | |
480 | // | |
481 | // read output | |
482 | // | |
483 | TFile *f = TFile::Open(testfile); | |
484 | if(!f || f->IsZombie()){ | |
485 | printf("File not readable\n"); | |
486 | return 0x0; | |
487 | } | |
488 | TString name("HFE_QA"); | |
489 | name += plus; | |
490 | printf("Name of TList %s\n",(const char*)name); | |
491 | TList *l = dynamic_cast<TList *>(f->Get((const char*)name)); | |
492 | if(!l){ | |
493 | printf("QA list was not found\n"); | |
494 | f->Close(); delete f; | |
495 | return 0x0; | |
496 | } | |
497 | TList *returnlist = dynamic_cast<TList *>(l->Clone()); | |
498 | f->Close(); delete f; | |
499 | return returnlist; | |
500 | } | |
501 | ||
502 | //_____________________________________________________________________ | |
503 | THnSparseF* GetHadronEffbyIPcut(TList *hfeqa){ | |
504 | ||
505 | ||
506 | // get hadron reduction factors due to the IP cuts | |
507 | TList* tl=(TList*)hfeqa->FindObject("list_TaskQA"); | |
508 | TH1F* hbefore=(TH1F*)tl->FindObject("hadronsBeforeIPcut"); | |
509 | TH1F* hafter=(TH1F*)tl->FindObject("hadronsAfterIPcut"); | |
510 | TH1F* hreduction= (TH1F*)hafter->Clone("hreduction"); | |
511 | hbefore->Sumw2(); | |
512 | hafter->Sumw2(); | |
513 | hreduction->Sumw2(); | |
514 | hreduction->Divide(hbefore); | |
515 | ||
516 | Double_t* binEdges[0]; | |
517 | Int_t hnBin = hreduction->GetNbinsX(); | |
518 | Int_t nBin[1] = {hnBin}; | |
519 | ||
520 | for(int i=0; i<nBin[0]; i++){ | |
521 | binEdges[0][i] = hreduction->GetBinLowEdge(i+1); | |
522 | } | |
523 | binEdges[0][nBin[0]] = hreduction->GetBinLowEdge(nBin[0]) + hreduction->GetBinWidth(nBin[0]); | |
524 | ||
525 | THnSparseF* hsreduction = new THnSparseF("hadroncontamin", "hadroncontamin; pt[GeV/c]", 1, nBin); | |
526 | hsreduction->SetBinEdges(0, binEdges[0]); | |
527 | ||
528 | Double_t dataE[1]; | |
529 | Double_t yval; | |
530 | for(int i=0; i<nBin[0]; i++){ | |
531 | dataE[0]=hreduction->GetBinCenter(i+1); | |
532 | yval=hreduction->GetBinContent(i+1); | |
533 | hsreduction->Fill(dataE,yval); | |
534 | } | |
535 | ||
536 | Int_t* ibins; | |
537 | ibins = new Int_t[nBin[0] + 1]; | |
538 | hsreduction->SetBinError(ibins,0); | |
539 | ||
540 | ||
541 | return hsreduction; | |
542 | } | |
543 | ||
544 | //_________________________________________________________________________ | |
545 | TObject* GetSpectrum(AliCFContainer *c, Int_t step) { | |
546 | AliCFDataGrid* data = new AliCFDataGrid("data","",*c,step); | |
547 | //data->SetMeasured(step); | |
548 | return data; | |
549 | } | |
550 | //___________________________________________________________________________ | |
551 | void CorrectFromTheWidth(TH1D *h1) { | |
552 | // | |
553 | // Correct from the width of the bins --> dN/dp_{T} (GeV/c)^{-1} | |
554 | // | |
555 | ||
556 | TAxis *axis = h1->GetXaxis(); | |
557 | Int_t nbinX = h1->GetNbinsX(); | |
558 | ||
559 | for(Int_t i = 1; i <= nbinX; i++) { | |
560 | ||
561 | Double_t width = axis->GetBinWidth(i); | |
562 | Double_t content = h1->GetBinContent(i); | |
563 | Double_t error = h1->GetBinError(i); | |
564 | h1->SetBinContent(i,content/width); | |
565 | h1->SetBinError(i,error/width); | |
566 | } | |
567 | ||
568 | } |