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