]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGHF/hfe/AliHFEInclusiveSpectrumQA.cxx
Fix of sigmaZ for crossing tracklets from Alex
[u/mrichter/AliRoot.git] / PWGHF / hfe / AliHFEInclusiveSpectrumQA.cxx
CommitLineData
959ea9d8 1
2/**************************************************************************
3* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4* *
5* Author: The ALICE Off-line Project. *
6* Contributors are mentioned in the code where appropriate. *
7* *
8* Permission to use, copy, modify and distribute this software and its *
9* documentation strictly for non-commercial purposes is hereby granted *
10* without fee, provided that the above copyright notice appears in all *
11* copies and that both the copyright notice and this permission notice *
12* appear in the supporting documentation. The authors make no claims *
13* about the suitability of this software for any purpose. It is *
14* provided "as is" without express or implied warranty. *
15**************************************************************************/
16//
17// Class for spectrum correction
18// Subtraction of hadronic background, Unfolding of the data and
19// Renormalization done here
20// The following containers have to be set:
21// - Correction framework container for real data
22// - Correction framework container for MC (Efficiency Map)
23// - Correction framework container for background coming from data
24// - Correction framework container for background coming from MC
25//
26// Author:
27// Raphaelle Bailhache <R.Bailhache@gsi.de>
28//
29
30#include <TArrayD.h>
31#include <TH1.h>
32#include <TList.h>
33#include <TObjArray.h>
34#include <TObject.h>
35#include <TROOT.h>
36#include <TCanvas.h>
37#include <TLegend.h>
38#include <TStyle.h>
39#include <TMath.h>
40#include <TAxis.h>
41#include <TGraphErrors.h>
42#include <TFile.h>
43#include <TPad.h>
44#include <TH2D.h>
45#include <TF1.h>
46
47#include "AliPID.h"
48#include "AliCFContainer.h"
49#include "AliCFDataGrid.h"
50#include "AliCFEffGrid.h"
51#include "AliCFGridSparse.h"
52#include "AliCFUnfolding.h"
53#include "AliLog.h"
54
55#include "AliHFEInclusiveSpectrumQA.h"
56#include "AliHFECorrectSpectrumBase.h"
57#include "AliHFEcuts.h"
58#include "AliHFEcontainer.h"
59#include "AliHFEtools.h"
60
61ClassImp(AliHFEInclusiveSpectrumQA)
62
63const Char_t *AliHFEInclusiveSpectrumQA::fgkNameCanvas[AliHFEInclusiveSpectrumQA::kNTypeEfficiency] = {
64 "V0Efficiency",
65 "MCEfficiency",
66 "ParametrizedEfficiency"
67};
68
afb48e1d 69const Char_t *AliHFEInclusiveSpectrumQA::fgkNameCanvasND[AliHFEInclusiveSpectrumQA::kNTypeEfficiency] = {
70 "V0EfficiencyND",
71 "MCEfficiencyND",
72 "ParametrizedEfficiencyND"
73};
74
959ea9d8 75//____________________________________________________________
76AliHFEInclusiveSpectrumQA::AliHFEInclusiveSpectrumQA():
77 TNamed(),
afb48e1d 78 fPtMax(10.0),
959ea9d8 79 fListOfResult(),
80 fWriteToFile(kTRUE)
81{
82 //
83 // Default constructor
84 //
85
86 fListOfResult = new TObjArray(kNResults);
87 fListOfResult->SetName("ListOfResults");
88
89
90}
91//____________________________________________________________
92AliHFEInclusiveSpectrumQA::AliHFEInclusiveSpectrumQA(const char *name):
93 TNamed(name, ""),
afb48e1d 94 fPtMax(10.0),
959ea9d8 95 fListOfResult(),
96 fWriteToFile(kTRUE)
97{
98 //
99 // Default constructor
100 //
101
102 fListOfResult = new TObjArray(kNResults);
103 fListOfResult->SetName("ListOfResults");
104
105
106}
107
108//____________________________________________________________
109AliHFEInclusiveSpectrumQA::~AliHFEInclusiveSpectrumQA(){
110 //
111 // Destructor
112 //
113 if(fListOfResult) delete fListOfResult;
114
115}
116//____________________________________________________________
117void AliHFEInclusiveSpectrumQA::AddResultAt(TObject *obj,Int_t index)
118{
119 //
120 // Init what we need for the correction:
121 //
122
123 if(fListOfResult) fListOfResult->AddAt(obj,index);
124
125}
126//____________________________________________________________
127TObject *AliHFEInclusiveSpectrumQA::GetResult(Int_t index)
128{
129 //
130 // Get result
131 //
132
133 if(fListOfResult) return fListOfResult->UncheckedAt(index);
134 else return 0x0;
135
136}
137//____________________________________________________________
138void AliHFEInclusiveSpectrumQA::DrawProjections() const
139{
140 //
141 // get spectrum for beauty 2nd method
142 //
143 //
144 AliCFContainer *data = (AliCFContainer *) fListOfResult->UncheckedAt(kDataProjection);
145 THnSparseF *correlation = (THnSparseF *) fListOfResult->UncheckedAt(kCMProjection);
146 if(!data || !correlation) return;
147
148 Int_t ndimcont = data->GetNVar();
149 Int_t ndimcor = correlation->GetNdimensions();
150 Int_t charge = 3;
151 Int_t centrality = 5;
152 Int_t eta = 1;
153
154 TCanvas * canvas = new TCanvas("Projections","Projections",1000,700);
155 Int_t n = 0;
156 if(charge < ndimcont) n++;
157 if(centrality < ndimcont) n++;
158 if(eta < ndimcont) n++;
159 canvas->Divide(2,n);
160 Int_t counter = 1;
161
162 if(charge < ndimcont) {
163
164 canvas->cd(counter);
165 TH1 *checkcharge = (TH1 *) data->Project(data->GetNStep()-1,charge);
166 checkcharge->Draw();
167 counter++;
168 canvas->cd(counter);
169 TH2F* projectioncharge = (TH2F *) correlation->Projection(charge,charge+((Int_t)(ndimcor/2.)));
170 projectioncharge->Draw("colz");
171 counter++;
172
173 }
174
175 if(centrality < ndimcont) {
176 canvas->cd(counter);
177 TH1 *checkcentrality = (TH1 *) data->Project(data->GetNStep()-1,centrality);
178 checkcentrality->Draw();
179 counter++;
180 canvas->cd(counter);
181 TH2F *projectioncentrality = (TH2F *) correlation->Projection(centrality,centrality+((Int_t)(ndimcor/2.)));
182 projectioncentrality->Draw("colz");
183 counter++;
184 }
185
186 if(eta < ndimcont) {
187 canvas->cd(counter);
188 TH1 *checketa = (TH1 *) data->Project(data->GetNStep()-1,eta);
189 checketa->Draw();
190 counter++;
191 canvas->cd(counter);
192 TH2D* projectioneta = (TH2D *) correlation->Projection(eta,eta+((Int_t)(ndimcor/2.)));
193 projectioneta->Draw("colz");
194 }
afb48e1d 195
196 if(fWriteToFile) canvas->SaveAs("Projections.png");
197
959ea9d8 198
199
200}
201//____________________________________________________________
202void AliHFEInclusiveSpectrumQA::DrawSubtractContamination() const
203{
204 //
205 // get spectrum for beauty 2nd method
206 //
207 //
208 TH1D *measuredTH1Daftersubstraction = (TH1D *) fListOfResult->UncheckedAt(kAfterSC);
209 TH1D *measuredTH1Dbeforesubstraction = (TH1D *) fListOfResult->UncheckedAt(kBeforeSC);
210 if(!measuredTH1Daftersubstraction || !measuredTH1Dbeforesubstraction) return;
211
212 SetStyle();
213
214 TCanvas * cbackgroundsubtraction = new TCanvas("backgroundsubtraction","backgroundsubtraction",1000,700);
215 cbackgroundsubtraction->Divide(2,1);
216 cbackgroundsubtraction->cd(1);
217 gPad->SetLogy();
218 gPad->SetTicks();
219 measuredTH1Daftersubstraction->SetStats(0);
220 measuredTH1Daftersubstraction->SetTitle("");
221 measuredTH1Daftersubstraction->GetYaxis()->SetTitleOffset(1.5);
222 measuredTH1Daftersubstraction->GetYaxis()->SetTitle("dN/dp_{T} [(GeV/c)^{-1}]");
223 measuredTH1Daftersubstraction->GetXaxis()->SetTitle("p^{rec}_{T} [GeV/c]");
224 measuredTH1Daftersubstraction->GetXaxis()->SetRangeUser(0.0,fPtMax);
225 measuredTH1Daftersubstraction->SetMarkerStyle(25);
226 measuredTH1Daftersubstraction->SetMarkerColor(kBlack);
227 measuredTH1Daftersubstraction->SetLineColor(kBlack);
228 measuredTH1Dbeforesubstraction->SetStats(0);
229 measuredTH1Dbeforesubstraction->SetTitle("");
230 measuredTH1Dbeforesubstraction->GetYaxis()->SetTitle("dN/dp_{T} [(GeV/c)^{-1}]");
231 measuredTH1Dbeforesubstraction->GetXaxis()->SetTitle("p^{rec}_{T} [GeV/c]");
232 measuredTH1Dbeforesubstraction->GetXaxis()->SetRangeUser(0.0,fPtMax);
233 measuredTH1Dbeforesubstraction->SetMarkerStyle(24);
234 measuredTH1Dbeforesubstraction->SetMarkerColor(kBlue);
235 measuredTH1Dbeforesubstraction->SetLineColor(kBlue);
236 measuredTH1Daftersubstraction->Draw();
237 measuredTH1Dbeforesubstraction->Draw("same");
238 TLegend *legsubstraction = new TLegend(0.4,0.6,0.89,0.89);
239 legsubstraction->AddEntry(measuredTH1Dbeforesubstraction,"With hadron contamination","p");
240 legsubstraction->AddEntry(measuredTH1Daftersubstraction,"Without hadron contamination ","p");
241 legsubstraction->SetFillStyle(0);
242 legsubstraction->SetLineStyle(0);
243 legsubstraction->SetLineColor(0);
244 legsubstraction->Draw("same");
245 cbackgroundsubtraction->cd(2);
246 gPad->SetLogy();
247 gPad->SetTicks();
248 TH1D* ratiomeasuredcontamination = (TH1D*)measuredTH1Dbeforesubstraction->Clone();
249 ratiomeasuredcontamination->SetName("ratiomeasuredcontamination");
250 ratiomeasuredcontamination->SetTitle("");
251 ratiomeasuredcontamination->GetYaxis()->SetTitleOffset(1.5);
252 ratiomeasuredcontamination->GetYaxis()->SetTitle("(with contamination - without contamination) / with contamination");
253 ratiomeasuredcontamination->GetXaxis()->SetTitle("p^{rec}_{T} [GeV/c]");
254 ratiomeasuredcontamination->GetYaxis()->SetRangeUser(0.8,1.2);
255 ratiomeasuredcontamination->GetXaxis()->SetRangeUser(0.0,fPtMax);
256 ratiomeasuredcontamination->Sumw2();
257 ratiomeasuredcontamination->Add(measuredTH1Daftersubstraction,-1.0);
258 ratiomeasuredcontamination->Divide(measuredTH1Dbeforesubstraction);
259 ratiomeasuredcontamination->SetStats(0);
260 ratiomeasuredcontamination->SetMarkerStyle(26);
261 ratiomeasuredcontamination->SetMarkerColor(kBlack);
262 ratiomeasuredcontamination->SetLineColor(kBlack);
263 for(Int_t k=0; k < ratiomeasuredcontamination->GetNbinsX(); k++){
264 ratiomeasuredcontamination->SetBinError(k+1,0.0);
265 }
266 ratiomeasuredcontamination->Draw("P");
267 if(fWriteToFile) cbackgroundsubtraction->SaveAs("BackgroundSubtracted.png");
268
269}
270
afb48e1d 271//____________________________________________________________
272void AliHFEInclusiveSpectrumQA::DrawSubtractContaminationND() const
273{
274 //
275 // subtract the hadron contamination
276 //
277 //
278 AliCFDataGrid *afterE = (AliCFDataGrid *) fListOfResult->UncheckedAt(kAfterSCND);
279 AliCFDataGrid *beforeE = (AliCFDataGrid *) fListOfResult->UncheckedAt(kBeforeSCND);
280 AliCFDataGrid *contamination = (AliCFDataGrid *) fListOfResult->UncheckedAt(kHadronContaminationND);
281
282 if(!afterE || !beforeE || !contamination) return;
283
284 SetStyle();
285
286 TCanvas * cbackgroundsubtractionND = new TCanvas("backgroundsubtractionND","backgroundsubtractionND",1000,700);
287 cbackgroundsubtractionND->Divide(3,1);
288 cbackgroundsubtractionND->cd(1);
289 gPad->SetLogz();
290 gPad->SetTicks();
291 TH2D *measuredTH2Dbeforesubstraction = (TH2D *) beforeE->Project(0,1);
292 measuredTH2Dbeforesubstraction->SetStats(0);
293 measuredTH2Dbeforesubstraction->SetTitle("Before contamination");
294 measuredTH2Dbeforesubstraction->GetZaxis()->SetTitleOffset(1.5);
295 measuredTH2Dbeforesubstraction->GetZaxis()->SetTitle("dN/dp_{T} [(GeV/c)^{-1}]");
296 measuredTH2Dbeforesubstraction->GetXaxis()->SetTitle("p^{rec}_{T} [GeV/c]");
297 measuredTH2Dbeforesubstraction->GetXaxis()->SetRangeUser(0.0,fPtMax);
298 measuredTH2Dbeforesubstraction->Draw("lego");
299 cbackgroundsubtractionND->cd(2);
300 gPad->SetLogz();
301 gPad->SetTicks();
302 TH2D *measuredTH2Daftersubstraction = (TH2D *) afterE->Project(0,1);
303 measuredTH2Daftersubstraction->SetStats(0);
304 measuredTH2Daftersubstraction->SetTitle("After contamination");
305 measuredTH2Daftersubstraction->GetZaxis()->SetTitleOffset(1.5);
306 measuredTH2Daftersubstraction->GetZaxis()->SetTitle("dN/dp_{T} [(GeV/c)^{-1}]");
307 measuredTH2Daftersubstraction->GetXaxis()->SetTitle("p^{rec}_{T} [GeV/c]");
308 measuredTH2Daftersubstraction->GetXaxis()->SetRangeUser(0.0,fPtMax);
309 measuredTH2Daftersubstraction->Draw("lego");
310 cbackgroundsubtractionND->cd(3);
311 gPad->SetLogz();
312 gPad->SetTicks();
313 TH2D *measuredsubstraction = (TH2D *) contamination->Project(0,1);
314 measuredsubstraction->SetStats(0);
315 measuredsubstraction->SetTitle("Contamination");
316 measuredsubstraction->GetZaxis()->SetTitleOffset(1.5);
317 measuredsubstraction->GetZaxis()->SetTitle("dN/dp_{T} [(GeV/c)^{-1}]");
318 measuredsubstraction->GetXaxis()->SetTitle("p^{rec}_{T} [GeV/c]");
319 measuredsubstraction->GetXaxis()->SetRangeUser(0.0,fPtMax);
320 measuredsubstraction->SetMarkerStyle(25);
321 measuredsubstraction->SetMarkerColor(kBlack);
322 measuredsubstraction->SetLineColor(kBlack);
323 measuredsubstraction->Draw("colz");
324 if(fWriteToFile) cbackgroundsubtractionND->SaveAs("BackgroundSubtractedND.png");
325
326}
327
63bdf450 328//____________________________________________________________
329void AliHFEInclusiveSpectrumQA::DrawSubtractPhotonicBackground() const
330{
331 //
332 // get spectrum
333 //
334
335 TH1D *measuredTH1Dafterphotonicsubstraction = (TH1D *) fListOfResult->UncheckedAt(kAfterSPB);
336 TH1D *measuredTH1Dbeforephotonicsubstraction = (TH1D *) fListOfResult->UncheckedAt(kBeforeSPB);
337 if(!measuredTH1Dafterphotonicsubstraction || !measuredTH1Dbeforephotonicsubstraction) return;
338
339 SetStyle();
340
341 TCanvas * cphotonic = new TCanvas("Photonic Subtraction","Photonic Subtraction",1000,700);
342 cphotonic->Divide(2,1);
343 cphotonic->cd(1);
344 gPad->SetLogy();
345 gPad->SetTicks();
346 measuredTH1Dafterphotonicsubstraction->SetStats(0);
347 measuredTH1Dafterphotonicsubstraction->SetTitle("");
348 measuredTH1Dafterphotonicsubstraction->GetYaxis()->SetTitleOffset(1.5);
349 measuredTH1Dafterphotonicsubstraction->GetYaxis()->SetTitle("dN/dp_{T} [(GeV/c)^{-1}]");
350 measuredTH1Dafterphotonicsubstraction->GetXaxis()->SetTitle("p^{rec}_{T} [GeV/c]");
351 measuredTH1Dafterphotonicsubstraction->GetXaxis()->SetRangeUser(0.0,fPtMax);
352 measuredTH1Dafterphotonicsubstraction->SetMarkerStyle(25);
353 measuredTH1Dafterphotonicsubstraction->SetMarkerColor(kBlack);
354 measuredTH1Dafterphotonicsubstraction->SetLineColor(kBlack);
355 measuredTH1Dbeforephotonicsubstraction->SetStats(0);
356 measuredTH1Dbeforephotonicsubstraction->SetTitle("");
357 measuredTH1Dbeforephotonicsubstraction->GetYaxis()->SetTitle("dN/dp_{T} [(GeV/c)^{-1}]");
358 measuredTH1Dbeforephotonicsubstraction->GetXaxis()->SetTitle("p^{rec}_{T} [GeV/c]");
359 measuredTH1Dbeforephotonicsubstraction->GetXaxis()->SetRangeUser(0.0,fPtMax);
360 measuredTH1Dbeforephotonicsubstraction->SetMarkerStyle(24);
361 measuredTH1Dbeforephotonicsubstraction->SetMarkerColor(kBlue);
362 measuredTH1Dbeforephotonicsubstraction->SetLineColor(kBlue);
363 measuredTH1Dafterphotonicsubstraction->Draw();
364 measuredTH1Dbeforephotonicsubstraction->Draw("same");
365 TLegend *legsubstraction = new TLegend(0.4,0.6,0.89,0.89);
366 legsubstraction->AddEntry(measuredTH1Dbeforephotonicsubstraction,"With photonic background","p");
367 legsubstraction->AddEntry(measuredTH1Dafterphotonicsubstraction,"Without photonic background","p");
368 legsubstraction->SetFillStyle(0);
369 legsubstraction->SetLineStyle(0);
370 legsubstraction->SetLineColor(0);
371 legsubstraction->Draw("same");
372 cphotonic->cd(2);
373 gPad->SetLogy();
374 gPad->SetTicks();
375 TH1D* ratiomeasuredphotonic = (TH1D*)measuredTH1Dbeforephotonicsubstraction->Clone();
376 ratiomeasuredphotonic->SetName("ratiomeasuredphotonic");
377 ratiomeasuredphotonic->SetTitle("");
378 ratiomeasuredphotonic->GetYaxis()->SetTitleOffset(1.5);
379 ratiomeasuredphotonic->GetYaxis()->SetTitle("(with photonic background - without photonic background) / with photonic background");
380 ratiomeasuredphotonic->GetXaxis()->SetTitle("p^{rec}_{T} [GeV/c]");
381 ratiomeasuredphotonic->GetYaxis()->SetRangeUser(0.8,1.2);
382 ratiomeasuredphotonic->GetXaxis()->SetRangeUser(0.0,fPtMax);
383 ratiomeasuredphotonic->Sumw2();
384 ratiomeasuredphotonic->Add(measuredTH1Dafterphotonicsubstraction,-1.0);
385 ratiomeasuredphotonic->Divide(measuredTH1Dbeforephotonicsubstraction);
386 ratiomeasuredphotonic->SetStats(0);
387 ratiomeasuredphotonic->SetMarkerStyle(26);
388 ratiomeasuredphotonic->SetMarkerColor(kBlack);
389 ratiomeasuredphotonic->SetLineColor(kBlack);
390 for(Int_t k=0; k < ratiomeasuredphotonic->GetNbinsX(); k++){
391 ratiomeasuredphotonic->SetBinError(k+1,0.0);
392 }
393 ratiomeasuredphotonic->Draw("P");
394 if(fWriteToFile) cphotonic->SaveAs("PhotonicSubtracted.png");
395
396}
959ea9d8 397
398//____________________________________________________________
399void AliHFEInclusiveSpectrumQA::DrawCorrectWithEfficiency(Int_t typeeff) const
400{
401 //
402 // Correct the spectrum for efficiency and unfolding
403 // with both method and compare
404 //
405
406 TH1D *afterE = 0x0;
407 TH1D *beforeE = 0x0;
408 TH1D *efficiencyDproj = 0x0;
409 TF1 *efficiencyparametrized = 0x0;
410
411 if(typeeff== kV0) {
412 afterE = (TH1D *) fListOfResult->UncheckedAt(kAfterV0);
413 beforeE = (TH1D *) fListOfResult->UncheckedAt(kBeforeV0);
414 efficiencyDproj = (TH1D *) fListOfResult->UncheckedAt(kV0Efficiency);
415 }
416 if(typeeff== kMC) {
417 afterE = (TH1D *) fListOfResult->UncheckedAt(kAfterMCE);
418 beforeE = (TH1D *) fListOfResult->UncheckedAt(kBeforeMCE);
419 efficiencyDproj = (TH1D *) fListOfResult->UncheckedAt(kMCEfficiency);
420 }
421 if(typeeff== kParametrized) {
422 afterE = (TH1D *) fListOfResult->UncheckedAt(kAfterPE);
423 beforeE = (TH1D *) fListOfResult->UncheckedAt(kBeforePE);
424 efficiencyparametrized = (TF1 *) fListOfResult->UncheckedAt(kPEfficiency);
425 }
426
32c709ed 427 if(!afterE || !beforeE) return;
428
afb48e1d 429 if((typeeff==kV0 || typeeff==kMC) && (!efficiencyDproj)) return;
430 if(typeeff==kParametrized && (!efficiencyparametrized)) return;
959ea9d8 431
63bdf450 432 SetStyle();
959ea9d8 433
434 TCanvas * cEfficiency = new TCanvas(AliHFEInclusiveSpectrumQA::fgkNameCanvas[typeeff],AliHFEInclusiveSpectrumQA::fgkNameCanvas[typeeff],1000,700);
435 cEfficiency->Divide(2,1);
436 cEfficiency->cd(1);
437 gPad->SetLogy();
438 gPad->SetTicks();
439 afterE->SetStats(0);
440 afterE->SetTitle("");
441 afterE->GetYaxis()->SetTitleOffset(1.5);
442 afterE->GetYaxis()->SetTitle("dN/dp_{T} [(GeV/c)^{-1}]");
443 afterE->GetXaxis()->SetTitle("p^{rec}_{T} [GeV/c]");
444 afterE->GetXaxis()->SetRangeUser(0.0,fPtMax);
445 afterE->SetMarkerStyle(25);
446 afterE->SetMarkerColor(kBlack);
447 afterE->SetLineColor(kBlack);
448 beforeE->SetStats(0);
449 beforeE->SetTitle("");
450 beforeE->GetYaxis()->SetTitle("dN/dp_{T} [(GeV/c)^{-1}]");
451 beforeE->GetXaxis()->SetTitle("p_{T} [GeV/c]");
452 beforeE->GetXaxis()->SetRangeUser(0.0,fPtMax);
453 beforeE->SetMarkerStyle(24);
454 beforeE->SetMarkerColor(kBlue);
455 beforeE->SetLineColor(kBlue);
456 afterE->Draw();
457 beforeE->Draw("same");
458 TLegend *legefficiency = new TLegend(0.4,0.6,0.89,0.89);
459 legefficiency->AddEntry(beforeE,"Before Efficiency correction","p");
460 legefficiency->AddEntry(afterE,"After Efficiency correction","p");
461 legefficiency->SetFillStyle(0);
462 legefficiency->SetLineStyle(0);
463 legefficiency->SetLineColor(0);
464 legefficiency->Draw("same");
465 cEfficiency->cd(2);
466 gPad->SetTicks();
467 if((typeeff==kV0 || typeeff==kMC)) {
91e50e2b 468 if(efficiencyDproj) {
469 efficiencyDproj->SetTitle("");
470 efficiencyDproj->SetStats(0);
471 efficiencyDproj->GetYaxis()->SetTitleOffset(1.5);
472 efficiencyDproj->GetYaxis()->SetRangeUser(0.0,1.0);
473 efficiencyDproj->GetYaxis()->SetTitle("Efficiency");
474 efficiencyDproj->GetXaxis()->SetTitle("p^{rec}_{T} [GeV/c]");
475 efficiencyDproj->GetXaxis()->SetRangeUser(0.0,fPtMax);
476 efficiencyDproj->SetMarkerStyle(25);
477 efficiencyDproj->Draw();
478 }
959ea9d8 479 }
480 if(typeeff==kParametrized) {
afb48e1d 481 if(efficiencyparametrized) {
482 efficiencyparametrized->GetYaxis()->SetTitleOffset(1.5);
483 efficiencyparametrized->GetYaxis()->SetRangeUser(0.0,1.0);
484 efficiencyparametrized->GetYaxis()->SetTitle("Efficiency");
485 efficiencyparametrized->GetXaxis()->SetTitle("p^{rec}_{T} [GeV/c]");
486 efficiencyparametrized->GetXaxis()->SetRangeUser(0.0,fPtMax);
487 efficiencyparametrized->Draw();
488 }
959ea9d8 489 }
490
491 if(fWriteToFile) {
492 if(typeeff==kV0) cEfficiency->SaveAs("EfficiencyV0.png");
493 if(typeeff==kMC) cEfficiency->SaveAs("EfficiencyMC.png");
494 if(typeeff==kParametrized) cEfficiency->SaveAs("EfficiencyParametrized.png");
495 }
496
afb48e1d 497}
498//____________________________________________________________
499void AliHFEInclusiveSpectrumQA::DrawCorrectWithEfficiencyND(Int_t typeeff) const
500{
501 //
502 // Correct the spectrum for efficiency and unfolding
503 // with both method and compare
504 //
505
506 AliCFDataGrid *afterE = 0x0;
507 AliCFDataGrid *beforeE = 0x0;
508 AliCFEffGrid *efficiencyND = 0x0;
509 TF1 *efficiencyparametrized = 0x0;
510
511 if(typeeff== kV0) {
512 afterE = (AliCFDataGrid *) fListOfResult->UncheckedAt(kAfterV0ND);
513 beforeE = (AliCFDataGrid *) fListOfResult->UncheckedAt(kBeforeV0ND);
514 efficiencyND = (AliCFEffGrid *) fListOfResult->UncheckedAt(kV0EfficiencyND);
515 }
516 if(typeeff== kMC) {
517 afterE = (AliCFDataGrid *) fListOfResult->UncheckedAt(kAfterMCEND);
518 beforeE = (AliCFDataGrid *) fListOfResult->UncheckedAt(kBeforeMCEND);
519 efficiencyND = (AliCFEffGrid *) fListOfResult->UncheckedAt(kMCEfficiencyND);
520 }
521 if(typeeff== kParametrized) {
522 afterE = (AliCFDataGrid *) fListOfResult->UncheckedAt(kAfterPEND);
523 beforeE = (AliCFDataGrid *) fListOfResult->UncheckedAt(kBeforePEND);
524 efficiencyparametrized = (TF1 *) fListOfResult->UncheckedAt(kPEfficiencyND);
525 }
526
527 if(!afterE || !beforeE) return;
528
529 if((typeeff==kV0 || typeeff==kMC) && (!efficiencyND)) return;
530 if(typeeff==kParametrized && (!efficiencyparametrized)) return;
531
532 SetStyle();
533
534 TCanvas * cEfficiency = new TCanvas(AliHFEInclusiveSpectrumQA::fgkNameCanvasND[typeeff],AliHFEInclusiveSpectrumQA::fgkNameCanvasND[typeeff],1000,700);
535 cEfficiency->Divide(3,1);
536 cEfficiency->cd(1);
537 gPad->SetLogz();
538 gPad->SetTicks();
539 TH2D *b2D = (TH2D *) beforeE->Project(0,1);
540 b2D->SetStats(0);
541 b2D->SetTitle("Before efficiency correction");
542 b2D->GetZaxis()->SetTitleOffset(1.5);
543 b2D->GetZaxis()->SetTitle("dN/dp_{T} [(GeV/c)^{-1}]");
544 b2D->GetXaxis()->SetTitle("p^{rec}_{T} [GeV/c]");
545 b2D->GetXaxis()->SetRangeUser(0.0,fPtMax);
546 b2D->Draw("lego");
547 cEfficiency->cd(2);
548 gPad->SetLogz();
549 gPad->SetTicks();
550 TH2D *a2D = (TH2D *) afterE->Project(0,1);
551 a2D->SetStats(0);
552 a2D->SetTitle("After efficiency correction");
553 a2D->GetZaxis()->SetTitleOffset(1.5);
554 a2D->GetZaxis()->SetTitle("dN/dp_{T} [(GeV/c)^{-1}]");
555 a2D->GetXaxis()->SetTitle("p^{rec}_{T} [GeV/c]");
556 a2D->GetXaxis()->SetRangeUser(0.0,fPtMax);
557 a2D->Draw("lego");
558 cEfficiency->cd(3);
559 gPad->SetTicks();
560 if((typeeff==kV0 || typeeff==kMC)) {
561 if(efficiencyND) {
562 THnSparseF *gride = (THnSparseF *) efficiencyND->GetGrid();
563 TH2D *e2D = (TH2D *) gride->Projection(1,0);
564 e2D->SetStats(0);
565 e2D->SetTitle("");
566 e2D->SetStats(0);
567 e2D->GetZaxis()->SetTitleOffset(1.5);
568 e2D->GetZaxis()->SetRangeUser(0.0,1.0);
569 e2D->GetZaxis()->SetTitle("Efficiency");
570 e2D->GetXaxis()->SetTitle("p^{rec}_{T} [GeV/c]");
571 e2D->GetXaxis()->SetRangeUser(0.0,fPtMax);
572 e2D->Draw("lego");
573 }
574 }
575 if(typeeff==kParametrized) {
576 if(efficiencyparametrized) {
577 efficiencyparametrized->GetYaxis()->SetTitleOffset(1.5);
578 efficiencyparametrized->GetYaxis()->SetRangeUser(0.0,1.0);
579 efficiencyparametrized->GetYaxis()->SetTitle("Efficiency");
580 efficiencyparametrized->GetXaxis()->SetTitle("p^{rec}_{T} [GeV/c]");
581 efficiencyparametrized->GetXaxis()->SetRangeUser(0.0,fPtMax);
582 efficiencyparametrized->Draw();
583 }
584 }
585
586 if(fWriteToFile) {
587 if(typeeff==kV0) cEfficiency->SaveAs("EfficiencyV0ND.png");
588 if(typeeff==kMC) cEfficiency->SaveAs("EfficiencyMCND.png");
589 if(typeeff==kParametrized) cEfficiency->SaveAs("EfficiencyParametrizedND.png");
590 }
591
959ea9d8 592}
593
594//____________________________________________________________
595void AliHFEInclusiveSpectrumQA::DrawUnfolding() const
596{
597 //
598 // Draw unfolding
599 //
600 TH1D *measuredspectrumD = (TH1D *) fListOfResult->UncheckedAt(kBeforeU);
601 TH1D *residualspectrumD = (TH1D *) fListOfResult->UncheckedAt(kResidualU);
602 TH1D *efficiencyDproj = (TH1D *) fListOfResult->UncheckedAt(kUEfficiency);
603 THnSparseF *correlation = (THnSparseF *) fListOfResult->UncheckedAt(kCMProjection);
604
605 if(!measuredspectrumD || !residualspectrumD || !efficiencyDproj || !correlation) return;
606
607 Int_t ndimcor = (Int_t) correlation->GetNdimensions()/2.;
608
609 SetStyle();
610
611 TCanvas * cunfolding = new TCanvas("unfolding","unfolding",1000,700);
612 cunfolding->Divide(2,2);
613 cunfolding->cd(1);
614 gPad->SetLogy();
615 gPad->SetTicks();
616 residualspectrumD->GetYaxis()->SetTitle("dN/dp_{T} [(GeV/c)^{-1}]");
617 residualspectrumD->GetXaxis()->SetTitle("p^{rec}_{T} [GeV/c]");
618 residualspectrumD->GetXaxis()->SetRangeUser(0.0,fPtMax);
619 residualspectrumD->SetStats(0);
620 residualspectrumD->SetTitle("");
621 residualspectrumD->GetYaxis()->SetTitleOffset(1.5);
622 residualspectrumD->SetMarkerStyle(26);
623 residualspectrumD->SetMarkerColor(kBlue);
624 residualspectrumD->SetLineColor(kBlue);
625 residualspectrumD->Sumw2();
626 residualspectrumD->Draw("P");
627 measuredspectrumD->SetStats(0);
628 measuredspectrumD->SetTitle("");
629 measuredspectrumD->GetYaxis()->SetTitleOffset(1.5);
630 measuredspectrumD->SetMarkerStyle(25);
631 measuredspectrumD->SetMarkerColor(kBlack);
632 measuredspectrumD->SetLineColor(kBlack);
633 measuredspectrumD->Draw("same");
634 TLegend *legres = new TLegend(0.4,0.6,0.89,0.89);
635 legres->AddEntry(residualspectrumD,"Residual","p");
636 legres->AddEntry(measuredspectrumD,"Measured","p");
637 legres->SetFillStyle(0);
638 legres->SetLineStyle(0);
639 legres->SetLineColor(0);
640 legres->Draw("same");
641 cunfolding->cd(2);
642 gPad->SetTicks();
643 TH1D* ratioresidual = (TH1D*)residualspectrumD->Clone();
644 ratioresidual->SetName("ratioresidual");
645 ratioresidual->SetTitle("");
646 ratioresidual->GetYaxis()->SetRangeUser(0.6,1.4);
647 ratioresidual->GetYaxis()->SetTitle("Folded/Measured");
648 ratioresidual->GetXaxis()->SetTitle("p_{T} [GeV/c]");
649 ratioresidual->Divide(measuredspectrumD);
650 ratioresidual->SetStats(0);
651 ratioresidual->Draw();
652 cunfolding->cd(3);
653 gPad->SetTicks();
654 efficiencyDproj->GetYaxis()->SetTitle("Efficiency");
655 efficiencyDproj->GetXaxis()->SetTitle("p^{MC}_{T} [GeV/c]");
656 efficiencyDproj->GetXaxis()->SetRangeUser(0.0,fPtMax);
657 efficiencyDproj->GetYaxis()->SetRangeUser(0.0,1.0);
658 efficiencyDproj->SetStats(0);
659 efficiencyDproj->SetTitle("");
660 efficiencyDproj->GetYaxis()->SetTitleOffset(1.5);
661 efficiencyDproj->SetMarkerStyle(26);
662 efficiencyDproj->SetMarkerColor(kBlue);
663 efficiencyDproj->SetLineColor(kBlue);
664 efficiencyDproj->Sumw2();
665 efficiencyDproj->Draw("P");
666 cunfolding->cd(4);
667 TH2F *projectioncorr = (TH2F *) correlation->Projection(0,ndimcor);
668 projectioncorr->GetYaxis()->SetTitle("p^{ESD}_{T} [GeV/c]");
669 projectioncorr->GetXaxis()->SetTitle("p^{MC}_{T} [GeV/c]");
670 projectioncorr->GetXaxis()->SetRangeUser(0.0,fPtMax);
671 projectioncorr->GetYaxis()->SetRangeUser(0.0,fPtMax);
672 projectioncorr->SetStats(0);
673 projectioncorr->SetTitle("");
674 projectioncorr->Draw("colz");
675
676 if(fWriteToFile){
677 cunfolding->SaveAs("Unfolding.png");
678 }
679
680}
681//____________________________________________________________
682void AliHFEInclusiveSpectrumQA::DrawResult()
683{
684 //
685 // Draw Results
686 //
687 TGraphErrors* correctedspectrumD = (TGraphErrors *) fListOfResult->UncheckedAt(kFinalResultUnfolded);
688 TGraphErrors* alltogetherspectrumD = (TGraphErrors *) fListOfResult->UncheckedAt( kFinalResultDirectEfficiency);
689 if(!correctedspectrumD || !alltogetherspectrumD) return;
690
afb48e1d 691
959ea9d8 692 correctedspectrumD->SetTitle("");
693 correctedspectrumD->GetYaxis()->SetTitleOffset(1.5);
694 correctedspectrumD->GetYaxis()->SetRangeUser(0.000001,100.0);
695 correctedspectrumD->GetXaxis()->SetRangeUser(0.0,fPtMax);
696 correctedspectrumD->SetMarkerStyle(26);
697 correctedspectrumD->SetMarkerColor(kBlue);
698 correctedspectrumD->SetLineColor(kBlue);
959ea9d8 699 alltogetherspectrumD->SetTitle("");
700 alltogetherspectrumD->GetYaxis()->SetTitleOffset(1.5);
701 alltogetherspectrumD->GetYaxis()->SetRangeUser(0.000000001,1.0);
702 alltogetherspectrumD->SetMarkerStyle(25);
703 alltogetherspectrumD->SetMarkerColor(kBlack);
704 alltogetherspectrumD->SetLineColor(kBlack);
959ea9d8 705 TLegend *legcorrected = new TLegend(0.4,0.6,0.89,0.89);
706 legcorrected->AddEntry(correctedspectrumD,"Unfolded","p");
707 legcorrected->AddEntry(alltogetherspectrumD,"Direct corrected","p");
708 legcorrected->SetFillStyle(0);
709 legcorrected->SetLineStyle(0);
710 legcorrected->SetLineColor(0);
afb48e1d 711
712
959ea9d8 713 TH1D* ratiocorrected = DivideSpectra(correctedspectrumD,alltogetherspectrumD);
afb48e1d 714 if(ratiocorrected) {
715 ratiocorrected->SetName("ratiocorrected");
716 ratiocorrected->SetTitle("");
717 ratiocorrected->GetYaxis()->SetTitleOffset(1.5);
718 ratiocorrected->GetYaxis()->SetTitle("Unfolded/DirectCorrected");
719 ratiocorrected->GetXaxis()->SetTitle("p_{T} [GeV/c]");
720 ratiocorrected->GetXaxis()->SetRangeUser(0.0,fPtMax);
721 ratiocorrected->GetYaxis()->SetRangeUser(0.4,1.4);
722 ratiocorrected->SetStats(0);
723 }
724
725
726 TCanvas * ccorrected = new TCanvas("corrected","corrected",1000,700);
727 if(ratiocorrected) ccorrected->Divide(2,1);
728 SetStyle();
729 ccorrected->cd(1);
730 gPad->SetLogy();
731 gPad->SetTicks();
732 correctedspectrumD->Draw("AP");
733 alltogetherspectrumD->Draw("P");
734 legcorrected->Draw("same");
735 if(ratiocorrected) {
736 ccorrected->cd(2);
737 gPad->SetTicks();
738 ratiocorrected->Draw();
739 }
959ea9d8 740 if(fWriteToFile)ccorrected->SaveAs("CorrectedResults.png");
741
742}
743//____________________________________________________________
744TH1D *AliHFEInclusiveSpectrumQA::DivideSpectra(TGraphErrors *ga, TGraphErrors *gb)
745{
746 //
747 // Divide Spectra
748 //
749
750 TH1D *afterE = (TH1D *) fListOfResult->UncheckedAt(kAfterMCE);
751 if(!afterE) return 0x0;
752
753 TH1D *histoB = (TH1D*) afterE->Clone();
754 histoB->Sumw2();
755 histoB->SetName("ratio");
756 TH1D *histoa = (TH1D*) afterE->Clone();
757 histoa->Sumw2();
758 histoa->SetName("a");
759 TH1D *histob = (TH1D*) afterE->Clone();
760 histob->Sumw2();
761 histob->SetName("b");
762
763 double xa,ya,xb,yb,eya,eyb;
764 Int_t npointsa = ga->GetN();
765 Int_t npointsb = gb->GetN();
766 if(npointsa != npointsb) {
767 printf("Problem the two spectra have not the same number of points\n");
768 return 0x0;
769 }
770 for(Int_t k = 0; k < npointsa; k++){
771 ga->GetPoint(k,xa,ya);
772 gb->GetPoint(k,xb,yb);
773 //
774 Double_t centerhisto = histoa->GetBinCenter(k+1);
775 //
776 if((TMath::Abs(xa-xb) > 0.0001) || (TMath::Abs(xa-centerhisto) > 0.0001)) {
777 printf("Problem not the same x axis\n");
778 return 0x0;
779 }
780 histoa->SetBinContent(k+1,ya);
781 histob->SetBinContent(k+1,yb);
782 //
783 eya = ga->GetErrorY(k);
784 eyb = gb->GetErrorY(k);
785 //
786 histoa->SetBinError(k+1,eya);
787 histob->SetBinError(k+1,eyb);
788
789 }
790
791 histoB->Sumw2();
792 histoB->Divide(histoa,histob,1.0,1.0,"B");
793
794 return histoB;
795
796}
797//__________________________________________
798void AliHFEInclusiveSpectrumQA::SetStyle() const
799{
800 //
801 // Set style
802 //
803
804 gStyle->SetPalette(1);
805 gStyle->SetOptStat(1111);
806 gStyle->SetPadBorderMode(0);
807 gStyle->SetCanvasColor(10);
808 gStyle->SetPadLeftMargin(0.13);
809 gStyle->SetPadRightMargin(0.13);
810
811}