]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWG1/AliPerformanceEff.cxx
Production point stored in AliFastGlauber.
[u/mrichter/AliRoot.git] / PWG1 / AliPerformanceEff.cxx
CommitLineData
777a0ba6 1//------------------------------------------------------------------------------\r
2// Implementation of AliPerformanceEff class. It keeps information from \r
3// comparison of reconstructed and MC particle tracks. In addtion, \r
4// it keeps selection cuts used during comparison. The comparison \r
5// information is stored in the ROOT histograms. Analysis of these \r
6// histograms can be done by using Analyse() class function. The result of \r
7// the analysis (histograms/graphs) are stored in the folder which is \r
8// a data member of AliPerformanceEff.\r
9// \r
10// Author: J.Otwinowski 04/02/2008 \r
11//------------------------------------------------------------------------------\r
12\r
13/*\r
14 \r
15 // after running comparison task, read the file, and get component\r
16 gROOT->LoadMacro("$ALICE_ROOT/PWG1/Macros/LoadMyLibs.C");\r
17 LoadMyLibs();\r
18 TFile f("Output.root");\r
19 //AliPerformanceEff * compObj = (AliPerformanceEff*)f.Get("AliPerformanceEff");\r
20 AliPerformanceEff * compObj = (AliPerformanceEff*)cOutput->FindObject("AliPerformanceEff");\r
21\r
22 // Analyse comparison data\r
23 compObj->Analyse();\r
24\r
25 // the output histograms/graphs will be stored in the folder "folderEff" \r
26 compObj->GetAnalysisFolder()->ls("*");\r
27\r
28 // user can save whole comparison object (or only folder with anlysed histograms) \r
29 // in the seperate output file (e.g.)\r
30 TFile fout("Analysed_Eff.root","recreate");\r
31 compObj->Write(); // compObj->GetAnalysisFolder()->Write();\r
32 fout.Close();\r
33\r
34*/\r
35\r
36#include <TAxis.h>\r
37#include <TH1D.h>\r
38\r
39// \r
40#include "AliESDtrack.h"\r
41#include "AliRecInfoCuts.h" \r
42#include "AliMCInfoCuts.h" \r
43#include "AliLog.h" \r
44#include "AliESDVertex.h" \r
45#include "AliExternalTrackParam.h" \r
46#include "AliTracker.h" \r
47#include "AliESDEvent.h" \r
48#include "AliMCEvent.h" \r
49#include "AliMCParticle.h" \r
50#include "AliHeader.h" \r
51#include "AliGenEventHeader.h" \r
52#include "AliStack.h" \r
53#include "AliPerformanceEff.h" \r
54\r
55using namespace std;\r
56\r
57ClassImp(AliPerformanceEff)\r
58\r
59//_____________________________________________________________________________\r
60AliPerformanceEff::AliPerformanceEff():\r
61 AliPerformanceObject("AliPerformanceEff"),\r
62\r
63 // histograms\r
64 \r
65 fEffHisto(0),\r
66\r
67 // Cuts \r
68 fCutsRC(0), \r
69 fCutsMC(0),\r
70\r
71 // histogram folder \r
72 fAnalysisFolder(0)\r
73{\r
74 // default consttructor \r
75 Init();\r
76}\r
77\r
78//_____________________________________________________________________________\r
79AliPerformanceEff::AliPerformanceEff(Char_t* name="AliPerformanceEff",Char_t*title="AliPerformanceEff",Int_t analysisMode=0, Bool_t hptGenerator=kFALSE):\r
80 AliPerformanceObject(name,title),\r
81\r
82 // histograms\r
83 fEffHisto(0),\r
84\r
85 // Cuts \r
86 fCutsRC(0), \r
87 fCutsMC(0),\r
88\r
89 // histogram folder \r
90 fAnalysisFolder(0)\r
91{\r
92 // named constructor\r
93 //\r
94 SetAnalysisMode(analysisMode);\r
95 SetHptGenerator(hptGenerator);\r
96\r
97 Init();\r
98}\r
99\r
100\r
101//_____________________________________________________________________________\r
102AliPerformanceEff::~AliPerformanceEff()\r
103{\r
104// destructor\r
105\r
106 if(fEffHisto) delete fEffHisto; fEffHisto=0;\r
107 if(fAnalysisFolder) delete fAnalysisFolder; fAnalysisFolder=0;\r
108}\r
109\r
110//_____________________________________________________________________________\r
111void AliPerformanceEff::Init()\r
112{\r
113 // Init histograms\r
114 //\r
115 // set pt bins\r
116 Int_t nPtBins = 50;\r
117 Double_t ptMin = 1.e-2, ptMax = 10.;\r
118\r
119 Double_t *binsPt = 0;\r
120 if (IsHptGenerator()) { \r
121 nPtBins = 100; ptMax = 100.;\r
122 binsPt = CreateLogAxis(nPtBins,ptMin,ptMax);\r
123 } else {\r
124 binsPt = CreateLogAxis(nPtBins,ptMin,ptMax);\r
125 }\r
126\r
127 /*\r
128 Int_t nPtBins = 31;\r
129 Double_t binsPt[32] = {0.,0.05,0.1,0.15,0.2,0.25,0.3,0.35,0.4,0.45,0.5,0.55,0.6,0.7,0.8,0.9,1.0,1.2,1.4,1.6,1.8,2.0,2.25,2.5,2.75,3.,3.5,4.,5.,6.,8.,10.};\r
130 Double_t ptMin = 0., ptMax = 10.;\r
131\r
132 if(IsHptGenerator() == kTRUE) {\r
133 nPtBins = 100;\r
134 ptMin = 0.; ptMax = 100.;\r
135 }\r
136 */\r
137\r
138 //mceta:mcphi:mcpt:pid:recStatus:findable\r
139 Int_t binsEffHisto[6]={30,144,nPtBins,5,2,2};\r
140 Double_t minEffHisto[6]={-1.5,0.,ptMin,0.,0.,0.};\r
141 Double_t maxEffHisto[6]={ 1.5,2.*TMath::Pi(), ptMax,5.,2.,2.};\r
142\r
143 fEffHisto = new THnSparseF("fEffHisto","mceta:mcphi:mcpt:pid:recStatus:findable",6,binsEffHisto,minEffHisto,maxEffHisto);\r
144 fEffHisto->SetBinEdges(2,binsPt);\r
145\r
146 fEffHisto->GetAxis(0)->SetTitle("#eta_{mc}");\r
147 fEffHisto->GetAxis(1)->SetTitle("#phi_{mc} (rad)");\r
148 fEffHisto->GetAxis(2)->SetTitle("p_{Tmc} (GeV/c)");\r
149 fEffHisto->GetAxis(3)->SetTitle("pid");\r
150 fEffHisto->GetAxis(4)->SetTitle("recStatus");\r
151 fEffHisto->GetAxis(5)->SetTitle("findable");\r
152 fEffHisto->Sumw2();\r
153\r
154 // init cuts\r
155 if(!fCutsMC) \r
156 AliDebug(AliLog::kError, "ERROR: Cannot find AliMCInfoCuts object");\r
157 if(!fCutsRC) \r
158 AliDebug(AliLog::kError, "ERROR: Cannot find AliRecInfoCuts object");\r
159\r
160 // init folder\r
161 fAnalysisFolder = CreateFolder("folderEff","Analysis Efficiency Folder");\r
162}\r
163\r
164//_____________________________________________________________________________\r
165void AliPerformanceEff::ProcessTPC(AliMCEvent* const mcEvent, AliESDEvent *const esdEvent)\r
166{\r
167 // Fill TPC only efficiency comparison information \r
168 Int_t *labelsRec = new Int_t[esdEvent->GetNumberOfTracks()];\r
169 if(!labelsRec) \r
170 AliDebug(AliLog::kError, "Cannot create labelsRec");\r
171\r
172 Int_t *labelsAllRec = new Int_t[esdEvent->GetNumberOfTracks()];\r
173 if(!labelsAllRec) \r
174 AliDebug(AliLog::kError, "Cannot create labelsAllRec");\r
175\r
176 // loop over rec. tracks\r
177 AliESDtrack *track=0;\r
178 for (Int_t iTrack = 0; iTrack < esdEvent->GetNumberOfTracks(); iTrack++) \r
179 { \r
180 track = esdEvent->GetTrack(iTrack);\r
181 if(!track) continue;\r
182 if(track->Charge()==0) continue;\r
183 Int_t label = TMath::Abs(track->GetLabel()); \r
184 if(label == 0) continue;\r
185\r
186 labelsAllRec[iTrack]=label;\r
187\r
188 // TPC only\r
189 if(IsRecTPC(track) != 0) \r
190 labelsRec[iTrack]=label;\r
191 }\r
192\r
193 // \r
194 // MC histograms for efficiency studies\r
195 //\r
196 \r
197 AliStack *stack = mcEvent->Stack();\r
198 if (!stack) {\r
199 AliDebug(AliLog::kError, "Stack not available");\r
200 return;\r
201 }\r
202\r
203 //Int_t nPart = stack->GetNtrack();\r
204 Int_t nPart = stack->GetNprimary();\r
205 for (Int_t iMc = 0; iMc < nPart; ++iMc) \r
206 {\r
207 TParticle* particle = stack->Particle(iMc);\r
208 if (!particle) continue;\r
209 if (particle->GetPDG()->Charge() == 0.0) continue;\r
210 \r
211 // physical primary\r
212 //Bool_t prim = stack->IsPhysicalPrimary(iMc);\r
213\r
214 Bool_t findable = kFALSE;\r
215 for(Int_t iRec=0; iRec<esdEvent->GetNumberOfTracks(); ++iRec) \r
216 {\r
217 // check findable\r
218 if(iMc > 0 && iMc == labelsAllRec[iRec]) \r
219 {\r
220 findable = IsFindable(mcEvent,iMc);\r
221 break;\r
222 }\r
223 } \r
224\r
225 Bool_t recStatus = kFALSE;\r
226 for(Int_t iRec=0; iRec<esdEvent->GetNumberOfTracks(); ++iRec) \r
227 {\r
228 // check reconstructed\r
229 if(iMc > 0 && iMc == labelsRec[iRec]) \r
230 {\r
231 recStatus = kTRUE;\r
232 break;\r
233 }\r
234 }\r
235\r
236 // Only 5 charged particle species (e,mu,pi,K,p)\r
237 if (fCutsMC->IsPdgParticle(TMath::Abs(particle->GetPdgCode())) == kFALSE) continue; \r
238\r
239 // transform Pdg to Pid\r
240 Int_t pid = TransformToPID(particle);\r
241\r
242 Float_t mceta = particle->Eta();\r
243 Float_t mcphi = particle->Phi();\r
244 if(mcphi<0) mcphi += 2.*TMath::Pi();\r
245 Float_t mcpt = particle->Pt();\r
246\r
247 // Fill histograms\r
248 Double_t vEffHisto[6] = { mceta, mcphi, mcpt, pid, recStatus, findable}; \r
249 fEffHisto->Fill(vEffHisto);\r
250 }\r
251\r
252 if(labelsRec) delete [] labelsRec; labelsRec = 0;\r
253 if(labelsAllRec) delete [] labelsAllRec; labelsAllRec = 0;\r
254}\r
255\r
256//_____________________________________________________________________________\r
257void AliPerformanceEff::ProcessTPCITS(AliMCEvent* const mcEvent, AliESDEvent *const esdEvent)\r
258{\r
259 // Fill efficiency comparison information\r
260 Int_t *labelsRec = new Int_t[esdEvent->GetNumberOfTracks()];\r
261 if(!labelsRec) \r
262 AliDebug(AliLog::kError, "Cannot create labelsRec");\r
263\r
264 Int_t *labelsAllRec = new Int_t[esdEvent->GetNumberOfTracks()];\r
265 if(!labelsAllRec) \r
266 AliDebug(AliLog::kError, "Cannot create labelsAllRec");\r
267\r
268 // loop over rec. tracks\r
269 AliESDtrack *track=0;\r
270 for (Int_t iTrack = 0; iTrack < esdEvent->GetNumberOfTracks(); iTrack++) \r
271 { \r
272 track = esdEvent->GetTrack(iTrack);\r
273 if(!track) continue;\r
274 if(track->Charge()==0) continue;\r
275 Int_t label = TMath::Abs(track->GetLabel()); \r
276 if(label == 0) continue;\r
277\r
278 labelsAllRec[iTrack]=label;\r
279\r
280 // iTPC+ITS\r
281 if(IsRecTPCITS(track) != 0) \r
282 labelsRec[iTrack]=label;\r
283 }\r
284\r
285 // \r
286 // MC histograms for efficiency studies\r
287 //\r
288 \r
289 AliStack *stack = mcEvent->Stack();\r
290 if (!stack) {\r
291 AliDebug(AliLog::kError, "Stack not available");\r
292 return;\r
293 }\r
294\r
295 //Int_t nPart = stack->GetNtrack();\r
296 Int_t nPart = stack->GetNprimary();\r
297 for (Int_t iMc = 0; iMc < nPart; ++iMc) \r
298 {\r
299 TParticle* particle = stack->Particle(iMc);\r
300 if (!particle) continue;\r
301 if (particle->GetPDG()->Charge() == 0.0) continue;\r
302 \r
303 // physical primary\r
304 //Bool_t prim = stack->IsPhysicalPrimary(iMc);\r
305\r
306 Bool_t findable = kFALSE;\r
307 for(Int_t iRec=0; iRec<esdEvent->GetNumberOfTracks(); ++iRec) \r
308 {\r
309 // check findable\r
310 if(iMc > 0 && iMc == labelsAllRec[iRec]) \r
311 {\r
312 findable = IsFindable(mcEvent,iMc);\r
313 break;\r
314 }\r
315 } \r
316\r
317 Bool_t recStatus = kFALSE;\r
318 for(Int_t iRec=0; iRec<esdEvent->GetNumberOfTracks(); ++iRec) \r
319 {\r
320 // check reconstructed\r
321 if(iMc > 0 && iMc == labelsRec[iRec]) \r
322 {\r
323 recStatus = kTRUE;\r
324 break;\r
325 }\r
326 }\r
327\r
328 // Only 5 charged particle species (e,mu,pi,K,p)\r
329 if (fCutsMC->IsPdgParticle(TMath::Abs(particle->GetPdgCode())) == kFALSE) continue; \r
330\r
331 // transform Pdg to Pid\r
332 Int_t pid = TransformToPID(particle);\r
333\r
334 Float_t mceta = particle->Eta();\r
335 Float_t mcphi = particle->Phi();\r
336 if(mcphi<0) mcphi += 2.*TMath::Pi();\r
337 Float_t mcpt = particle->Pt();\r
338\r
339 // Fill histograms\r
340 Double_t vEffHisto[6] = { mceta, mcphi, mcpt, pid, recStatus, findable}; \r
341 fEffHisto->Fill(vEffHisto);\r
342 }\r
343\r
344 if(labelsRec) delete [] labelsRec; labelsRec = 0;\r
345 if(labelsAllRec) delete [] labelsAllRec; labelsAllRec = 0;\r
346}\r
347\r
348//_____________________________________________________________________________\r
349void AliPerformanceEff::ProcessConstrained(AliMCEvent* const mcEvent, AliESDEvent *const esdEvent)\r
350{\r
351 // Process comparison information \r
352 Int_t *labelsRec = new Int_t[esdEvent->GetNumberOfTracks()];\r
353 if(!labelsRec) \r
354 AliDebug(AliLog::kError, "Cannot create labelsRec");\r
355\r
356 Int_t *labelsAllRec = new Int_t[esdEvent->GetNumberOfTracks()];\r
357 if(!labelsAllRec) \r
358 AliDebug(AliLog::kError, "Cannot create labelsAllRec");\r
359\r
360 // loop over rec. tracks\r
361 AliESDtrack *track=0;\r
362 for (Int_t iTrack = 0; iTrack < esdEvent->GetNumberOfTracks(); iTrack++) \r
363 { \r
364 track = esdEvent->GetTrack(iTrack);\r
365 if(!track) continue;\r
366 if(track->Charge()==0) continue;\r
367 Int_t label = TMath::Abs(track->GetLabel()); \r
368 if(label == 0) continue;\r
369\r
370 labelsAllRec[iTrack]=label;\r
371\r
372 // Constrained\r
373 if(IsRecConstrained(track) != 0) \r
374 labelsRec[iTrack]=label;\r
375\r
376 }\r
377\r
378 // \r
379 // MC histograms for efficiency studies\r
380 //\r
381 \r
382 AliStack *stack = mcEvent->Stack();\r
383 if (!stack) {\r
384 AliDebug(AliLog::kError, "Stack not available");\r
385 return;\r
386 }\r
387\r
388 //Int_t nPart = stack->GetNtrack();\r
389 Int_t nPart = stack->GetNprimary();\r
390 for (Int_t iMc = 0; iMc < nPart; ++iMc) \r
391 {\r
392 TParticle* particle = stack->Particle(iMc);\r
393 if (!particle) continue;\r
394 if (particle->GetPDG()->Charge() == 0.0) continue;\r
395 \r
396 // physical primary\r
397 //Bool_t prim = stack->IsPhysicalPrimary(iMc);\r
398\r
399 Bool_t findable = kFALSE;\r
400 for(Int_t iRec=0; iRec<esdEvent->GetNumberOfTracks(); ++iRec) \r
401 {\r
402 // check findable\r
403 if(iMc > 0 && iMc == labelsAllRec[iRec]) \r
404 {\r
405 findable = IsFindable(mcEvent,iMc);\r
406 break;\r
407 }\r
408 } \r
409\r
410 Bool_t recStatus = kFALSE;\r
411 for(Int_t iRec=0; iRec<esdEvent->GetNumberOfTracks(); ++iRec) \r
412 {\r
413 // check reconstructed\r
414 if(iMc > 0 && iMc == labelsRec[iRec]) \r
415 {\r
416 recStatus = kTRUE;\r
417 break;\r
418 }\r
419 }\r
420\r
421 // Only 5 charged particle species (e,mu,pi,K,p)\r
422 if (fCutsMC->IsPdgParticle(TMath::Abs(particle->GetPdgCode())) == kFALSE) continue; \r
423\r
424 // transform Pdg to Pid\r
425 Int_t pid = TransformToPID(particle);\r
426\r
427 Float_t mceta = particle->Eta();\r
428 Float_t mcphi = particle->Phi();\r
429 if(mcphi<0) mcphi += 2.*TMath::Pi();\r
430 Float_t mcpt = particle->Pt();\r
431\r
432 // Fill histograms\r
433 Double_t vEffHisto[6] = { mceta, mcphi, mcpt, pid, recStatus, findable}; \r
434 fEffHisto->Fill(vEffHisto);\r
435 }\r
436\r
437 if(labelsRec) delete [] labelsRec; labelsRec = 0;\r
438 if(labelsAllRec) delete [] labelsAllRec; labelsAllRec = 0;\r
439}\r
440\r
441//_____________________________________________________________________________\r
442void AliPerformanceEff::Exec(AliMCEvent* const mcEvent, AliESDEvent* const esdEvent, const Bool_t bUseMC)\r
443{\r
444 // Process comparison information \r
445 //\r
446 if(!esdEvent) \r
447 {\r
448 AliDebug(AliLog::kError, "esdEvent not available");\r
449 return;\r
450 }\r
451 AliHeader* header = 0;\r
452 AliGenEventHeader* genHeader = 0;\r
453 AliStack* stack = 0;\r
454 TArrayF vtxMC(3);\r
455 \r
456 if(bUseMC)\r
457 {\r
458 if(!mcEvent) {\r
459 AliDebug(AliLog::kError, "mcEvent not available");\r
460 return;\r
461 }\r
462 // get MC event header\r
463 header = mcEvent->Header();\r
464 if (!header) {\r
465 AliDebug(AliLog::kError, "Header not available");\r
466 return;\r
467 }\r
468 // MC particle stack\r
469 stack = mcEvent->Stack();\r
470 if (!stack) {\r
471 AliDebug(AliLog::kError, "Stack not available");\r
472 return;\r
473 }\r
474 // get MC vertex\r
475 genHeader = header->GenEventHeader();\r
476 if (!genHeader) {\r
477 AliDebug(AliLog::kError, "Could not retrieve genHeader from Header");\r
478 return;\r
479 }\r
480 genHeader->PrimaryVertex(vtxMC);\r
481\r
482 } // end bUseMC\r
483\r
484 //\r
485 // Process events\r
486 //\r
487\r
488 if(GetAnalysisMode() == 0) ProcessTPC(mcEvent,esdEvent);\r
489 else if(GetAnalysisMode() == 1) ProcessTPCITS(mcEvent,esdEvent);\r
490 else if(GetAnalysisMode() == 2) ProcessConstrained(mcEvent,esdEvent);\r
491 else {\r
492 printf("ERROR: AnalysisMode %d \n",fAnalysisMode);\r
493 return;\r
494 }\r
495}\r
496\r
497//_____________________________________________________________________________\r
498Int_t AliPerformanceEff::TransformToPID(TParticle *particle) \r
499{\r
500// transform Pdg to Pid\r
501// Pdg convension is different for hadrons and leptons \r
502// (e.g. K+/K- = 321/-321; e+/e- = -11/11 ) \r
503\r
504 Int_t pid = -1;\r
505 if( TMath::Abs(particle->GetPdgCode())==fCutsMC->GetEM() ) pid = 0; \r
506 if( TMath::Abs(particle->GetPdgCode())==fCutsMC->GetMuM() ) pid = 1; \r
507 if( TMath::Abs(particle->GetPdgCode())==fCutsMC->GetPiP() ) pid = 2; \r
508 if( TMath::Abs(particle->GetPdgCode())==fCutsMC->GetKP() ) pid = 3; \r
509 if( TMath::Abs(particle->GetPdgCode())==fCutsMC->GetProt() ) pid = 4; \r
510\r
511return pid;\r
512}\r
513\r
514//_____________________________________________________________________________\r
515Bool_t AliPerformanceEff::IsFindable(AliMCEvent *mcEvent, Int_t label) \r
516{\r
517if(!mcEvent) return kFALSE;\r
518if(label==0) return kFALSE;\r
519\r
520 AliMCParticle *mcParticle = mcEvent->GetTrack(label);\r
521 if(!mcParticle) return kFALSE;\r
522\r
523 Int_t counter; \r
524 Float_t tpcTrackLength = mcParticle->GetTPCTrackLength(AliTracker::GetBz(),0.1,counter,3.0); \r
525 //printf("tpcTrackLength %f \n", tpcTrackLength);\r
526\r
527return (tpcTrackLength>fCutsMC->GetMinTrackLength()); \r
528}\r
529\r
530//_____________________________________________________________________________\r
531Bool_t AliPerformanceEff::IsRecTPC(AliESDtrack *esdTrack) \r
532{\r
533if(!esdTrack) return kFALSE;\r
534\r
535 const AliExternalTrackParam *track = esdTrack->GetTPCInnerParam();\r
536 if(!track) return kFALSE;\r
537\r
538 Float_t dca[2], cov[3]; // dca_xy, dca_z, sigma_xy, sigma_xy_z, sigma_z\r
539 esdTrack->GetImpactParametersTPC(dca,cov);\r
540 \r
541 Int_t label = TMath::Abs(esdTrack->GetLabel()); \r
542 if(label == 0) return kFALSE;\r
543\r
544 Bool_t recStatus = kFALSE;\r
545 if(esdTrack->GetTPCNcls()>fCutsRC->GetMinNClustersTPC() && \r
546 TMath::Abs(dca[0])<fCutsRC->GetMaxDCAToVertexXY() && \r
547 TMath::Abs(dca[1])<fCutsRC->GetMaxDCAToVertexZ())\r
548 {\r
549 recStatus = kTRUE;\r
550 }\r
551\r
552return recStatus;\r
553}\r
554\r
555//_____________________________________________________________________________\r
556Bool_t AliPerformanceEff::IsRecTPCITS(AliESDtrack *esdTrack) \r
557{\r
558if(!esdTrack) return kFALSE;\r
559\r
560 Float_t dca[2], cov[3]; // dca_xy, dca_z, sigma_xy, sigma_xy_z, sigma_z\r
561 esdTrack->GetImpactParameters(dca,cov);\r
562\r
563 Int_t label = TMath::Abs(esdTrack->GetLabel()); \r
564 if(label == 0) return kFALSE;\r
565\r
566 Bool_t recStatus = kFALSE;\r
567\r
568 if ((esdTrack->GetStatus()&AliESDtrack::kTPCrefit)==0) return kFALSE; // TPC refit\r
569 if (esdTrack->GetTPCNcls()<fCutsRC->GetMinNClustersTPC()) return kFALSE; // min. nb. TPC clusters\r
570 Int_t clusterITS[200];\r
571 if(esdTrack->GetITSclusters(clusterITS)<fCutsRC->GetMinNClustersITS()) return kFALSE; // min. nb. ITS clusters\r
572\r
573 if(TMath::Abs(dca[0])<fCutsRC->GetMaxDCAToVertexXY() && \r
574 TMath::Abs(dca[1])<fCutsRC->GetMaxDCAToVertexZ())\r
575 {\r
576 recStatus = kTRUE;\r
577 }\r
578\r
579return recStatus;\r
580}\r
581\r
582//_____________________________________________________________________________\r
583Bool_t AliPerformanceEff::IsRecConstrained(AliESDtrack *esdTrack) \r
584{\r
585 if(!esdTrack) return kFALSE;\r
586\r
587 const AliExternalTrackParam * track = esdTrack->GetConstrainedParam();\r
588 if(!track) return kFALSE;\r
589\r
590 Float_t dca[2], cov[3]; // dca_xy, dca_z, sigma_xy, sigma_xy_z, sigma_z\r
591 esdTrack->GetImpactParameters(dca,cov);\r
592\r
593 Int_t label = TMath::Abs(esdTrack->GetLabel()); \r
594 if(label == 0) return kFALSE;\r
595\r
596 Bool_t recStatus = kFALSE;\r
597\r
598 if ((esdTrack->GetStatus()&AliESDtrack::kTPCrefit)==0) return kFALSE; // TPC refit\r
599 if (esdTrack->GetTPCNcls()<fCutsRC->GetMinNClustersTPC()) return kFALSE; // min. nb. TPC clusters\r
600 Int_t clusterITS[200];\r
601 if(esdTrack->GetITSclusters(clusterITS)<fCutsRC->GetMinNClustersITS()) return kFALSE; // min. nb. ITS clusters\r
602\r
603 if(TMath::Abs(dca[0])<fCutsRC->GetMaxDCAToVertexXY() && \r
604 TMath::Abs(dca[1])<fCutsRC->GetMaxDCAToVertexZ())\r
605 {\r
606 recStatus = kTRUE;\r
607 }\r
608\r
609return recStatus;\r
610}\r
611\r
612//_____________________________________________________________________________\r
613Long64_t AliPerformanceEff::Merge(TCollection* const list) \r
614{\r
615 // Merge list of objects (needed by PROOF)\r
616\r
617 if (!list)\r
618 return 0;\r
619\r
620 if (list->IsEmpty())\r
621 return 1;\r
622\r
623 TIterator* iter = list->MakeIterator();\r
624 TObject* obj = 0;\r
625\r
626 // collection of generated histograms\r
627\r
628 Int_t count=0;\r
629 while((obj = iter->Next()) != 0) \r
630 {\r
631 AliPerformanceEff* entry = dynamic_cast<AliPerformanceEff*>(obj);\r
632 if (entry == 0) continue; \r
633 \r
634 fEffHisto->Add(entry->fEffHisto);\r
635 count++;\r
636 }\r
637\r
638return count;\r
639}\r
640 \r
641//_____________________________________________________________________________\r
642void AliPerformanceEff::Analyse() \r
643{\r
644 // Analyse comparison information and store output histograms\r
645 // in the folder "folderEff" \r
646 //\r
647 TH1::AddDirectory(kFALSE);\r
648 TObjArray *aFolderObj = new TObjArray;\r
649 char title[256];\r
650\r
651 //\r
652 // efficiency vs pt\r
653 //\r
654 fEffHisto->GetAxis(0)->SetRangeUser(-0.9,0.9); // eta range\r
655 fEffHisto->GetAxis(2)->SetRangeUser(0.1,10.); // pt range\r
656\r
657 // rec efficiency vs pt\r
658 TH1D *ptAll = fEffHisto->Projection(2);\r
659\r
660 fEffHisto->GetAxis(4)->SetRangeUser(1.,1.); // reconstructed \r
661 TH1D *ptRec = fEffHisto->Projection(2);\r
662 TH1D *ptRec_c = (TH1D*)ptRec->Clone();\r
663 ptRec_c->Divide(ptRec,ptAll,1,1,"B");\r
664 ptRec_c->SetName("ptRecEff");\r
665\r
666 ptRec_c->GetXaxis()->SetTitle(fEffHisto->GetAxis(2)->GetTitle());\r
667 ptRec_c->GetYaxis()->SetTitle("efficiency");\r
668 sprintf(title,"%s vs %s","rec. efficiency",fEffHisto->GetAxis(2)->GetTitle());\r
669 ptRec_c->SetTitle(title);\r
670\r
671 ptRec_c->SetBit(TH1::kLogX);\r
672 aFolderObj->Add(ptRec_c);\r
673\r
674 // rec efficiency vs pid vs pt\r
675\r
676 fEffHisto->GetAxis(4)->SetRangeUser(0.,1.); \r
677 fEffHisto->GetAxis(3)->SetRangeUser(2.,2.); // pions\r
678\r
679 TH1D *ptAllPi = fEffHisto->Projection(2);\r
680\r
681 fEffHisto->GetAxis(4)->SetRangeUser(1.,1.); // reconstructed\r
682 TH1D *ptRecPi = fEffHisto->Projection(2);\r
683 TH1D *ptRecPi_c = (TH1D*)ptRecPi->Clone();\r
684 ptRecPi_c->Divide(ptRecPi,ptAllPi,1,1,"B");\r
685 ptRecPi_c->SetName("ptRecEffPi");\r
686\r
687 ptRecPi_c->GetXaxis()->SetTitle(fEffHisto->GetAxis(2)->GetTitle());\r
688 ptRecPi_c->GetYaxis()->SetTitle("efficiency");\r
689 sprintf(title,"%s vs %s","rec. efficiency (pions)",fEffHisto->GetAxis(2)->GetTitle());\r
690 ptRecPi_c->SetTitle(title);\r
691\r
692 ptRecPi_c->SetBit(TH1::kLogX);\r
693 aFolderObj->Add(ptRecPi_c);\r
694\r
695 fEffHisto->GetAxis(4)->SetRangeUser(0.,1.); \r
696 fEffHisto->GetAxis(3)->SetRangeUser(3.,3.); // kaons\r
697 TH1D *ptAllK = fEffHisto->Projection(2);\r
698\r
699 fEffHisto->GetAxis(4)->SetRangeUser(1.,1.); // reconstructed\r
700 TH1D *ptRecK = fEffHisto->Projection(2);\r
701\r
702 TH1D *ptRecK_c = (TH1D*)ptRecK->Clone();\r
703 ptRecK_c->Divide(ptRecK,ptAllK,1,1,"B");\r
704 ptRecK_c->SetName("ptRecEffK");\r
705\r
706 ptRecK_c->GetXaxis()->SetTitle(fEffHisto->GetAxis(2)->GetTitle());\r
707 ptRecK_c->GetYaxis()->SetTitle("efficiency");\r
708 sprintf(title,"%s vs %s","rec. efficiency (kaons)",fEffHisto->GetAxis(2)->GetTitle());\r
709 ptRecK_c->SetTitle(title);\r
710\r
711\r
712 ptRecK_c->SetBit(TH1::kLogX);\r
713 aFolderObj->Add(ptRecK_c);\r
714\r
715 fEffHisto->GetAxis(4)->SetRangeUser(0.,1.); \r
716 fEffHisto->GetAxis(3)->SetRangeUser(4.,4.); // protons\r
717 TH1D *ptAllP = fEffHisto->Projection(2);\r
718\r
719 fEffHisto->GetAxis(4)->SetRangeUser(1.,1.); // reconstructed\r
720 TH1D *ptRecP = fEffHisto->Projection(2);\r
721 TH1D *ptRecP_c = (TH1D*)ptRecP->Clone();\r
722 ptRecP_c->Divide(ptRecP,ptAllP,1,1,"B");\r
723 ptRecP_c->SetName("ptRecEffP");\r
724\r
725 ptRecP_c->GetXaxis()->SetTitle(fEffHisto->GetAxis(2)->GetTitle());\r
726 ptRecP_c->GetYaxis()->SetTitle("efficiency");\r
727 sprintf(title,"%s vs %s","rec. efficiency (protons)",fEffHisto->GetAxis(2)->GetTitle());\r
728 ptRecP_c->SetTitle(title);\r
729\r
730 ptRecP_c->SetBit(TH1::kLogX);\r
731 aFolderObj->Add(ptRecP_c);\r
732\r
733 // findable efficiency vs pt\r
734\r
735 fEffHisto->GetAxis(3)->SetRangeUser(0.,4.); \r
736 fEffHisto->GetAxis(4)->SetRangeUser(0.,1.); \r
737 fEffHisto->GetAxis(5)->SetRangeUser(1.,1.); // findable\r
738 TH1D *ptAllF = fEffHisto->Projection(2);\r
739\r
740 fEffHisto->GetAxis(4)->SetRangeUser(1.,1.);\r
741 fEffHisto->GetAxis(5)->SetRangeUser(1.,1.);\r
742\r
743 TH1D *ptRecF = fEffHisto->Projection(2); // rec findable\r
744 TH1D *ptRecF_c = (TH1D*)ptRecF->Clone();\r
745 ptRecF_c->Divide(ptRecF,ptAllF,1,1,"B");\r
746 ptRecF_c->SetName("ptRecEffF");\r
747\r
748 ptRecF_c->GetXaxis()->SetTitle(fEffHisto->GetAxis(2)->GetTitle());\r
749 ptRecF_c->GetYaxis()->SetTitle("efficiency");\r
750 sprintf(title,"%s vs %s","rec. efficiency (findable)",fEffHisto->GetAxis(2)->GetTitle());\r
751 ptRecF_c->SetTitle(title);\r
752\r
753 ptRecF_c->SetBit(TH1::kLogX);\r
754 aFolderObj->Add(ptRecF_c);\r
755\r
756 //\r
757 // efficiency vs eta\r
758 //\r
759\r
760 fEffHisto->GetAxis(0)->SetRangeUser(-1.5,1.5); // eta range\r
761 fEffHisto->GetAxis(2)->SetRangeUser(0.2,10.); // pt range\r
762 fEffHisto->GetAxis(4)->SetRangeUser(0.,1.); // all\r
763 fEffHisto->GetAxis(5)->SetRangeUser(0.,1.); // all\r
764\r
765 TH1D *etaAll = fEffHisto->Projection(0);\r
766\r
767 fEffHisto->GetAxis(4)->SetRangeUser(1.,1.); // reconstructed \r
768 TH1D *etaRec = fEffHisto->Projection(0);\r
769 TH1D *etaRec_c = (TH1D*)etaRec->Clone();\r
770 etaRec_c->Divide(etaRec,etaAll,1,1,"B");\r
771 etaRec_c->SetName("etaRecEff");\r
772\r
773 etaRec_c->GetXaxis()->SetTitle(fEffHisto->GetAxis(0)->GetTitle());\r
774 etaRec_c->GetYaxis()->SetTitle("efficiency");\r
775 sprintf(title,"%s vs %s","rec. efficiency",fEffHisto->GetAxis(0)->GetTitle());\r
776 etaRec_c->SetTitle(title);\r
777\r
778 aFolderObj->Add(etaRec_c);\r
779\r
780 // rec efficiency vs pid vs eta\r
781 fEffHisto->GetAxis(4)->SetRangeUser(0.,1.); \r
782 fEffHisto->GetAxis(3)->SetRangeUser(2.,2.); // pions\r
783\r
784 TH1D *etaAllPi = fEffHisto->Projection(0);\r
785\r
786 fEffHisto->GetAxis(4)->SetRangeUser(1.,1.); // reconstructed\r
787 TH1D *etaRecPi = fEffHisto->Projection(0);\r
788 TH1D *etaRecPi_c = (TH1D*)etaRecPi->Clone();\r
789 etaRecPi_c->Divide(etaRecPi,etaAllPi,1,1,"B");\r
790 etaRecPi_c->SetName("etaRecEffPi");\r
791\r
792 etaRecPi_c->GetXaxis()->SetTitle(fEffHisto->GetAxis(0)->GetTitle());\r
793 etaRecPi_c->GetYaxis()->SetTitle("efficiency");\r
794 sprintf(title,"%s vs %s","rec. efficiency (pions)",fEffHisto->GetAxis(0)->GetTitle());\r
795 etaRecPi_c->SetTitle(title);\r
796\r
797 aFolderObj->Add(etaRecPi_c);\r
798\r
799 fEffHisto->GetAxis(4)->SetRangeUser(0.,1.); \r
800 fEffHisto->GetAxis(3)->SetRangeUser(3.,3.); // kaons\r
801 TH1D *etaAllK = fEffHisto->Projection(0);\r
802\r
803 fEffHisto->GetAxis(4)->SetRangeUser(1.,1.); // reconstructed\r
804 TH1D *etaRecK = fEffHisto->Projection(0);\r
805\r
806 TH1D *etaRecK_c = (TH1D*)etaRecK->Clone();\r
807 etaRecK_c->Divide(etaRecK,etaAllK,1,1,"B");\r
808 etaRecK_c->SetName("etaRecEffK");\r
809\r
810 etaRecK_c->GetXaxis()->SetTitle(fEffHisto->GetAxis(0)->GetTitle());\r
811 etaRecK_c->GetYaxis()->SetTitle("efficiency");\r
812 sprintf(title,"%s vs %s","rec. efficiency (kaons)",fEffHisto->GetAxis(0)->GetTitle());\r
813 etaRecK_c->SetTitle(title);\r
814\r
815\r
816 aFolderObj->Add(etaRecK_c);\r
817\r
818 fEffHisto->GetAxis(4)->SetRangeUser(0.,1.); \r
819 fEffHisto->GetAxis(3)->SetRangeUser(4.,4.); // protons\r
820 TH1D *etaAllP = fEffHisto->Projection(0);\r
821\r
822 fEffHisto->GetAxis(4)->SetRangeUser(1.,1.); // reconstructed\r
823 TH1D *etaRecP = fEffHisto->Projection(0);\r
824 TH1D *etaRecP_c = (TH1D*)etaRecP->Clone();\r
825 etaRecP_c->Divide(etaRecP,etaAllP,1,1,"B");\r
826 etaRecP_c->SetName("etaRecEffP");\r
827\r
828 etaRecP_c->GetXaxis()->SetTitle(fEffHisto->GetAxis(0)->GetTitle());\r
829 etaRecP_c->GetYaxis()->SetTitle("efficiency");\r
830 sprintf(title,"%s vs %s","rec. efficiency (protons)",fEffHisto->GetAxis(0)->GetTitle());\r
831 etaRecP_c->SetTitle(title);\r
832\r
833 aFolderObj->Add(etaRecP_c);\r
834\r
835 // findable efficiency vs eta\r
836\r
837 fEffHisto->GetAxis(3)->SetRangeUser(0.,4.); \r
838 fEffHisto->GetAxis(4)->SetRangeUser(0.,1.); \r
839 fEffHisto->GetAxis(5)->SetRangeUser(1.,1.); // findable\r
840 TH1D *etaAllF = fEffHisto->Projection(0);\r
841\r
842 fEffHisto->GetAxis(4)->SetRangeUser(1.,1.);\r
843 fEffHisto->GetAxis(5)->SetRangeUser(1.,1.);\r
844\r
845 TH1D *etaRecF = fEffHisto->Projection(0); // rec findable\r
846 TH1D *etaRecF_c = (TH1D*)etaRecF->Clone();\r
847 etaRecF_c->Divide(etaRecF,etaAllF,1,1,"B");\r
848 etaRecF_c->SetName("etaRecEffF");\r
849\r
850 etaRecF_c->GetXaxis()->SetTitle(fEffHisto->GetAxis(0)->GetTitle());\r
851 etaRecF_c->GetYaxis()->SetTitle("efficiency");\r
852 sprintf(title,"%s vs %s","rec. efficiency (findable)",fEffHisto->GetAxis(0)->GetTitle());\r
853 etaRecF_c->SetTitle(title);\r
854\r
855 aFolderObj->Add(etaRecF_c);\r
856\r
857 //\r
858 // efficiency vs phi\r
859 //\r
860\r
861 fEffHisto->GetAxis(0)->SetRangeUser(-0.9,0.9); // eta range\r
862 fEffHisto->GetAxis(2)->SetRangeUser(0.2,10.); // pt range\r
863 fEffHisto->GetAxis(4)->SetRangeUser(0.,1.); // all\r
864 fEffHisto->GetAxis(5)->SetRangeUser(0.,1.); // all\r
865\r
866 TH1D *phiAll = fEffHisto->Projection(1);\r
867\r
868 fEffHisto->GetAxis(4)->SetRangeUser(1.,1.); // reconstructed \r
869 TH1D *phiRec = fEffHisto->Projection(1);\r
870 TH1D *phiRec_c = (TH1D*)phiRec->Clone();\r
871 phiRec_c->Divide(phiRec,phiAll,1,1,"B");\r
872 phiRec_c->SetName("phiRecEff");\r
873\r
874 phiRec_c->GetXaxis()->SetTitle(fEffHisto->GetAxis(1)->GetTitle());\r
875 phiRec_c->GetYaxis()->SetTitle("efficiency");\r
876 sprintf(title,"%s vs %s","rec. efficiency",fEffHisto->GetAxis(1)->GetTitle());\r
877 phiRec_c->SetTitle(title);\r
878\r
879 aFolderObj->Add(phiRec_c);\r
880\r
881 // rec efficiency vs pid vs phi\r
882 fEffHisto->GetAxis(4)->SetRangeUser(0.,1.); \r
883 fEffHisto->GetAxis(3)->SetRangeUser(2.,2.); // pions\r
884\r
885 TH1D *phiAllPi = fEffHisto->Projection(1);\r
886\r
887 fEffHisto->GetAxis(4)->SetRangeUser(1.,1.); // reconstructed\r
888 TH1D *phiRecPi = fEffHisto->Projection(1);\r
889 TH1D *phiRecPi_c = (TH1D*)phiRecPi->Clone();\r
890 phiRecPi_c->Divide(phiRecPi,phiAllPi,1,1,"B");\r
891 phiRecPi_c->SetName("phiRecEffPi");\r
892\r
893 phiRecPi_c->GetXaxis()->SetTitle(fEffHisto->GetAxis(1)->GetTitle());\r
894 phiRecPi_c->GetYaxis()->SetTitle("efficiency");\r
895 sprintf(title,"%s vs %s","rec. efficiency (pions)",fEffHisto->GetAxis(1)->GetTitle());\r
896 phiRecPi_c->SetTitle(title);\r
897\r
898 aFolderObj->Add(phiRecPi_c);\r
899\r
900 fEffHisto->GetAxis(4)->SetRangeUser(0.,1.); \r
901 fEffHisto->GetAxis(3)->SetRangeUser(3.,3.); // kaons\r
902 TH1D *phiAllK = fEffHisto->Projection(1);\r
903\r
904 fEffHisto->GetAxis(4)->SetRangeUser(1.,1.); // reconstructed\r
905 TH1D *phiRecK = fEffHisto->Projection(1);\r
906\r
907 TH1D *phiRecK_c = (TH1D*)phiRecK->Clone();\r
908 phiRecK_c->Divide(phiRecK,phiAllK,1,1,"B");\r
909 phiRecK_c->SetName("phiRecEffK");\r
910\r
911 phiRecK_c->GetXaxis()->SetTitle(fEffHisto->GetAxis(1)->GetTitle());\r
912 phiRecK_c->GetYaxis()->SetTitle("efficiency");\r
913 sprintf(title,"%s vs %s","rec. efficiency (kaons)",fEffHisto->GetAxis(1)->GetTitle());\r
914 phiRecK_c->SetTitle(title);\r
915\r
916\r
917 aFolderObj->Add(phiRecK_c);\r
918\r
919 fEffHisto->GetAxis(4)->SetRangeUser(0.,1.); \r
920 fEffHisto->GetAxis(3)->SetRangeUser(4.,4.); // protons\r
921 TH1D *phiAllP = fEffHisto->Projection(1);\r
922\r
923 fEffHisto->GetAxis(4)->SetRangeUser(1.,1.); // reconstructed\r
924 TH1D *phiRecP = fEffHisto->Projection(1);\r
925 TH1D *phiRecP_c = (TH1D*)phiRecP->Clone();\r
926 phiRecP_c->Divide(phiRecP,phiAllP,1,1,"B");\r
927 phiRecP_c->SetName("phiRecEffP");\r
928\r
929 phiRecP_c->GetXaxis()->SetTitle(fEffHisto->GetAxis(1)->GetTitle());\r
930 phiRecP_c->GetYaxis()->SetTitle("efficiency");\r
931 sprintf(title,"%s vs %s","rec. efficiency (protons)",fEffHisto->GetAxis(1)->GetTitle());\r
932 phiRecP_c->SetTitle(title);\r
933\r
934 aFolderObj->Add(phiRecP_c);\r
935\r
936 // findable efficiency vs phi\r
937\r
938 fEffHisto->GetAxis(3)->SetRangeUser(0.,4.); \r
939 fEffHisto->GetAxis(4)->SetRangeUser(0.,1.); \r
940 fEffHisto->GetAxis(5)->SetRangeUser(1.,1.); // findable\r
941 TH1D *phiAllF = fEffHisto->Projection(1);\r
942\r
943 fEffHisto->GetAxis(4)->SetRangeUser(1.,1.);\r
944 fEffHisto->GetAxis(5)->SetRangeUser(1.,1.);\r
945\r
946 TH1D *phiRecF = fEffHisto->Projection(1); // rec findable\r
947 TH1D *phiRecF_c = (TH1D*)phiRecF->Clone();\r
948 phiRecF_c->Divide(phiRecF,phiAllF,1,1,"B");\r
949 phiRecF_c->SetName("phiRecEffF");\r
950\r
951 phiRecF_c->GetXaxis()->SetTitle(fEffHisto->GetAxis(1)->GetTitle());\r
952 phiRecF_c->GetYaxis()->SetTitle("efficiency");\r
953 sprintf(title,"%s vs %s","rec. efficiency (findable)",fEffHisto->GetAxis(1)->GetTitle());\r
954 phiRecF_c->SetTitle(title);\r
955\r
956 aFolderObj->Add(phiRecF_c);\r
957\r
958 // export objects to analysis folder\r
959 fAnalysisFolder = ExportToFolder(aFolderObj);\r
960\r
961 // delete only TObjArray\r
962 if(aFolderObj) delete aFolderObj;\r
963}\r
964\r
965//_____________________________________________________________________________\r
966TFolder* AliPerformanceEff::ExportToFolder(TObjArray * array) \r
967{\r
968 // recreate folder avery time and export objects to new one\r
969 //\r
970 AliPerformanceEff * comp=this;\r
971 TFolder *folder = comp->GetAnalysisFolder();\r
972\r
973 TString name, title;\r
974 TFolder *newFolder = 0;\r
975 Int_t i = 0;\r
976 Int_t size = array->GetSize();\r
977\r
978 if(folder) { \r
979 // get name and title from old folder\r
980 name = folder->GetName(); \r
981 title = folder->GetTitle(); \r
982\r
983 // delete old one\r
984 delete folder;\r
985\r
986 // create new one\r
987 newFolder = CreateFolder(name.Data(),title.Data());\r
988 newFolder->SetOwner();\r
989\r
990 // add objects to folder\r
991 while(i < size) {\r
992 newFolder->Add(array->At(i));\r
993 i++;\r
994 }\r
995 }\r
996\r
997return newFolder;\r
998}\r
999\r
1000\r
1001//_____________________________________________________________________________\r
1002TFolder* AliPerformanceEff::CreateFolder(TString name,TString title) { \r
1003// create folder for analysed histograms\r
1004//\r
1005TFolder *folder = 0;\r
1006 folder = new TFolder(name.Data(),title.Data());\r
1007\r
1008 return folder;\r
1009}\r