20582a531b4a164f67521970bea29107a20ffbe1
[u/mrichter/AliRoot.git] / PWG0 / dNdPt / AlidNdPtCutAnalysis.cxx
1 /**************************************************************************\r
2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *\r
3  *                                                                        *\r
4  * Author: The ALICE Off-line Project.                                    *\r
5  * Contributors are mentioned in the code where appropriate.              *\r
6  *                                                                        *\r
7  * Permission to use, copy, modify and distribute this software and its   *\r
8  * documentation strictly for non-commercial purposes is hereby granted   *\r
9  * without fee, provided that the above copyright notice appears in all   *\r
10  * copies and that both the copyright notice and this permission notice   *\r
11  * appear in the supporting documentation. The authors make no claims     *\r
12  * about the suitability of this software for any purpose. It is          *\r
13  * provided "as is" without express or implied warranty.                  *\r
14  **************************************************************************/\r
15 //------------------------------------------------------------------------------\r
16 // AlidNdPtCutAnalysis class. \r
17 //\r
18 // a. functionality:\r
19 // - fills generic cut histograms\r
20 // - generates cuts (selection criteria)\r
21 //\r
22 // b. data members:\r
23 // - generic cut histograms\r
24 // - control histograms\r
25 //\r
26 // Author: J.Otwinowski 04/11/2008 \r
27 //------------------------------------------------------------------------------\r
28 #include "TH1.h"\r
29 #include "TH2.h"\r
30 \r
31 #include "AliHeader.h"  \r
32 #include "AliGenEventHeader.h"  \r
33 #include "AliStack.h"  \r
34 #include "AliESDEvent.h"  \r
35 #include "AliMCEvent.h"  \r
36 #include "AliESDtrackCuts.h"  \r
37 #include "AliLog.h" \r
38 #include "AliTracker.h" \r
39 \r
40 #include "AlidNdPtEventCuts.h"\r
41 #include "AlidNdPtAcceptanceCuts.h"\r
42 #include "AlidNdPtBackgroundCuts.h"\r
43 #include "AliPhysicsSelection.h"\r
44 \r
45 #include "AliPWG0Helper.h"\r
46 #include "AlidNdPtHelper.h"\r
47 #include "AlidNdPtCutAnalysis.h"\r
48 \r
49 using namespace std;\r
50 \r
51 ClassImp(AlidNdPtCutAnalysis)\r
52 \r
53 //_____________________________________________________________________________\r
54   AlidNdPtCutAnalysis::AlidNdPtCutAnalysis(): AlidNdPt(),\r
55   fAnalysisFolder(0),\r
56   fEventCount(0),\r
57   fRecEventHist(0),\r
58   fMCEventHist(0),\r
59   fRecMCEventHist(0),\r
60   fRecMCTrackHist(0)\r
61 {\r
62   // default constructor\r
63   Init();\r
64 }\r
65 \r
66 //_____________________________________________________________________________\r
67 AlidNdPtCutAnalysis::AlidNdPtCutAnalysis(Char_t* name, Char_t* title): AlidNdPt(name,title),\r
68   fAnalysisFolder(0),\r
69   fEventCount(0),\r
70   fRecEventHist(0),\r
71   fMCEventHist(0),\r
72   fRecMCEventHist(0),\r
73   fRecMCTrackHist(0)\r
74 {\r
75   // constructor\r
76   Init();\r
77 }\r
78 \r
79 //_____________________________________________________________________________\r
80 AlidNdPtCutAnalysis::~AlidNdPtCutAnalysis() {\r
81   // \r
82   if(fEventCount) delete fEventCount; fEventCount=0;\r
83   if(fRecEventHist) delete fRecEventHist; fRecEventHist=0;\r
84   if(fMCEventHist) delete fMCEventHist; fMCEventHist=0;\r
85   if(fRecMCEventHist) delete fRecMCEventHist; fRecMCEventHist=0;\r
86   if(fRecMCTrackHist) delete fRecMCTrackHist; fRecMCTrackHist=0;\r
87 \r
88   if(fAnalysisFolder) delete fAnalysisFolder; fAnalysisFolder=0;\r
89 }\r
90 \r
91 //_____________________________________________________________________________\r
92 void AlidNdPtCutAnalysis::Init(){\r
93   //\r
94   // Init histograms\r
95   //\r
96   const Int_t ptNbins = 56; \r
97   const Double_t ptMin = 0.; \r
98   const Double_t ptMax = 16.; \r
99 \r
100   Double_t binsPt[ptNbins+1] = {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.65,0.7,0.75,0.8,0.85,0.9,0.95,1.0,1.1,1.2,1.3,1.4,1.5,1.6,1.7,1.8,1.9,2.0,2.2,2.4,2.6,2.8,3.0,3.2,3.4,3.6,3.8,4.0,4.5,5.0,5.5,6.0,6.5,7.0,7.5,8.0,9.0,10.0,11.0,12.0,13.0,14.0,15.0,16.0};\r
101 \r
102   // \r
103   Int_t binsEventCount[2]={2,2};\r
104   Double_t minEventCount[2]={0,0}; \r
105   Double_t maxEventCount[2]={2,2}; \r
106   fEventCount = new THnSparseF("fEventCount","trig vs trig+vertex",2,binsEventCount,minEventCount,maxEventCount);\r
107   fEventCount->GetAxis(0)->SetTitle("trig");\r
108   fEventCount->GetAxis(1)->SetTitle("trig+vert");\r
109   fEventCount->Sumw2();\r
110 \r
111   //Xv:Yv:Zv:ResZv:Mult\r
112   Double_t kFact = 1.0;\r
113 \r
114   Int_t binsRecEventHist[5]={80,80,100,80,150};\r
115   Double_t minRecEventHist[5]={-3.*kFact,-3.*kFact,-35.,0.,0.}; \r
116   Double_t maxRecEventHist[5]={3.*kFact,3.*kFact,35.,10.,150.}; \r
117   fRecEventHist = new THnSparseF("fRecEventHist","Xv:Yv:Zv:ResZv:Mult",5,binsRecEventHist,minRecEventHist,maxRecEventHist);\r
118   fRecEventHist->GetAxis(0)->SetTitle("Xv (cm)");\r
119   fRecEventHist->GetAxis(1)->SetTitle("Yv (cm)");\r
120   fRecEventHist->GetAxis(2)->SetTitle("Zv (cm)");\r
121   fRecEventHist->GetAxis(3)->SetTitle("ResZv (cm)");\r
122   fRecEventHist->GetAxis(4)->SetTitle("Mult");\r
123   fRecEventHist->Sumw2();\r
124 \r
125   //Xv:Yv:Zv\r
126   Int_t binsMCEventHist[3]={80,80,100};\r
127   Double_t minMCEventHist[3]={-0.1,-0.1,-35.}; \r
128   Double_t maxMCEventHist[3]={0.1,0.1,35.}; \r
129   fMCEventHist = new THnSparseF("fMCEventHist","mcXv:mcYv:mcZv",3,binsMCEventHist,minMCEventHist,maxMCEventHist);\r
130   fMCEventHist->GetAxis(0)->SetTitle("mcXv (cm)");\r
131   fMCEventHist->GetAxis(1)->SetTitle("mcYv (cm)");\r
132   fMCEventHist->GetAxis(2)->SetTitle("mcZv (cm)");\r
133   fMCEventHist->Sumw2();\r
134 \r
135   //Xv-mcXv:Yv-mcYv:Zv-mcZv:Mult\r
136   Int_t binsRecMCEventHist[4]={100,100,100,150};\r
137   Double_t minRecMCEventHist[4]={-1.0*kFact,-1.0*kFact,-1.0*kFact,0.}; \r
138   Double_t maxRecMCEventHist[4]={1.0*kFact,1.0*kFact,1.0*kFact,150.}; \r
139   fRecMCEventHist = new THnSparseF("fRecMCEventHist","Xv-mcXv:Yv-mcYv:Zv-mcZv:Mult",4,binsRecMCEventHist,minRecMCEventHist,maxRecMCEventHist);\r
140   fRecMCEventHist->GetAxis(0)->SetTitle("Xv-mcXv (cm)");\r
141   fRecMCEventHist->GetAxis(1)->SetTitle("Yv-mcYv (cm)");\r
142   fRecMCEventHist->GetAxis(2)->SetTitle("Zv-mcZv (cm)");\r
143   fRecMCEventHist->GetAxis(3)->SetTitle("Mult");\r
144   fRecMCEventHist->Sumw2();\r
145 \r
146   //\r
147   // THnSparse track histograms\r
148   //\r
149 \r
150   //nClust:chi2PerClust:nClust/nFindableClust:DCAy:DCAz:eta:phi:pt:isKink:isPrim:polarity\r
151   Int_t binsRecMCTrackHist[11]={160,80,80,80,80,30,90,ptNbins, 2, 2, 2};\r
152   Double_t minRecMCTrackHist[11]={0., 0., 0., -10.,-10.,-1.5, 0., ptMin, 0., 0., 0.};\r
153   Double_t maxRecMCTrackHist[11]={160.,10.,1.2, 10.,10.,1.5, 2.*TMath::Pi(), ptMax, 2.,2., 2.};\r
154 \r
155   fRecMCTrackHist = new THnSparseF("fRecMCTrackHist","nClust:chi2PerClust:nClust/nFindableClust:DCAy:DCAz:eta:phi:pt:isKink:isPrim:polarity",11,binsRecMCTrackHist,minRecMCTrackHist,maxRecMCTrackHist);\r
156   fRecMCTrackHist->SetBinEdges(7,binsPt);\r
157 \r
158   fRecMCTrackHist->GetAxis(0)->SetTitle("nClust");\r
159   fRecMCTrackHist->GetAxis(1)->SetTitle("chi2PerClust");\r
160   fRecMCTrackHist->GetAxis(2)->SetTitle("nClust/nFindableClust");\r
161   fRecMCTrackHist->GetAxis(3)->SetTitle("DCAy (cm)");\r
162   fRecMCTrackHist->GetAxis(4)->SetTitle("DCAz (cm)");\r
163   fRecMCTrackHist->GetAxis(5)->SetTitle("#eta");\r
164   fRecMCTrackHist->GetAxis(6)->SetTitle("#phi (rad)");\r
165   fRecMCTrackHist->GetAxis(7)->SetTitle("p_{T} (GeV/c)");\r
166   fRecMCTrackHist->GetAxis(8)->SetTitle("isKink");\r
167   fRecMCTrackHist->GetAxis(9)->SetTitle("isPrim");\r
168   fRecMCTrackHist->GetAxis(10)->SetTitle("polarity");\r
169   fRecMCTrackHist->Sumw2();\r
170 \r
171   // init output folder\r
172   fAnalysisFolder = CreateFolder("folderdNdPt","Analysis dNdPt Folder");\r
173 \r
174 }\r
175 \r
176 //_____________________________________________________________________________\r
177 void AlidNdPtCutAnalysis::Process(AliESDEvent *const esdEvent, AliMCEvent * const mcEvent)\r
178 {\r
179   //\r
180   // Process real and/or simulated events\r
181   //\r
182   if(!esdEvent) {\r
183     AliDebug(AliLog::kError, "esdEvent not available");\r
184     return;\r
185   }\r
186 \r
187   // get selection cuts\r
188   AlidNdPtEventCuts *evtCuts = GetEventCuts(); \r
189   AlidNdPtAcceptanceCuts *accCuts = GetAcceptanceCuts(); \r
190   AliESDtrackCuts *esdTrackCuts = GetTrackCuts(); \r
191 \r
192   if(!evtCuts || !accCuts  || !esdTrackCuts) {\r
193     AliDebug(AliLog::kError, "cuts not available");\r
194     return;\r
195   }\r
196 \r
197   // trigger selection\r
198   Bool_t isEventTriggered = kTRUE;\r
199   if(evtCuts->IsTriggerRequired())  \r
200   {\r
201     AliPhysicsSelection *trigSel = GetPhysicsTriggerSelection();\r
202     if(!trigSel) {\r
203       AliDebug(AliLog::kError, "cannot get trigSel");\r
204       return;\r
205     }\r
206 \r
207     if(IsUseMCInfo()) { \r
208       //static AliTriggerAnalysis* triggerAnalysis = new AliTriggerAnalysis;\r
209       //isEventTriggered = triggerAnalysis->IsTriggerFired(esdEvent, GetTrigger());\r
210       trigSel->SetAnalyzeMC();\r
211       isEventTriggered = trigSel->IsCollisionCandidate(esdEvent);\r
212     }\r
213     else {\r
214       isEventTriggered = trigSel->IsCollisionCandidate(esdEvent);\r
215     }\r
216   }\r
217 \r
218   // use MC information\r
219   AliHeader* header = 0;\r
220   AliGenEventHeader* genHeader = 0;\r
221   AliStack* stack = 0;\r
222   TArrayF vtxMC(3);\r
223   //AlidNdPtHelper::MCProcessType evtType = AlidNdPtHelper::kInvalidProcess;\r
224   AliPWG0Helper::MCProcessType evtType = AliPWG0Helper::kInvalidProcess;\r
225 \r
226   if(IsUseMCInfo())\r
227   {\r
228     if(!mcEvent) {\r
229       AliDebug(AliLog::kError, "mcEvent not available");\r
230       return;\r
231     }\r
232 \r
233     // get MC event header\r
234     header = mcEvent->Header();\r
235     if (!header) {\r
236       AliDebug(AliLog::kError, "Header not available");\r
237       return;\r
238     }\r
239 \r
240     // MC particle stack\r
241     stack = mcEvent->Stack();\r
242     if (!stack) {\r
243       AliDebug(AliLog::kError, "Stack not available");\r
244       return;\r
245     }\r
246 \r
247     // get event type (ND=0x1, DD=0x2, SD=0x4)\r
248     evtType = AliPWG0Helper::GetEventProcessType(header);\r
249     AliDebug(AliLog::kDebug+1, Form("Found process type %d", evtType));\r
250 \r
251     // get MC vertex\r
252     genHeader = header->GenEventHeader();\r
253     if (!genHeader) {\r
254       AliDebug(AliLog::kError, "Could not retrieve genHeader from Header");\r
255       return;\r
256     }\r
257     genHeader->PrimaryVertex(vtxMC);\r
258 \r
259     // Fill MC event histogram\r
260     Double_t vMCEventHist[3]={vtxMC[0],vtxMC[1],vtxMC[2]};\r
261     fMCEventHist->Fill(vMCEventHist);\r
262 \r
263   } // end bUseMC\r
264 \r
265   // get reconstructed vertex  \r
266   Bool_t bRedoTPCVertex = evtCuts->IsRedoTPCVertex();\r
267   Bool_t bUseConstraints = evtCuts->IsUseBeamSpotConstraint();\r
268   const AliESDVertex* vtxESD = AlidNdPtHelper::GetVertex(esdEvent,evtCuts,accCuts,esdTrackCuts,GetAnalysisMode(),kFALSE,bRedoTPCVertex,bUseConstraints); \r
269   Bool_t isRecVertex = AlidNdPtHelper::TestRecVertex(vtxESD, GetAnalysisMode(), kFALSE);\r
270   Bool_t isEventOK = evtCuts->AcceptEvent(esdEvent,mcEvent,vtxESD) && isRecVertex;\r
271 \r
272   TObjArray *allChargedTracks=0;\r
273   Int_t multAll=0;\r
274   \r
275   //\r
276   // event counter\r
277   // \r
278   //printf("isEventOK %d, isEventTriggered %d \n",isEventOK,isEventTriggered);\r
279 \r
280   Bool_t isTrigAndVertex = isEventTriggered && isEventOK;\r
281   Double_t vEventCount[2] = { isEventTriggered, isTrigAndVertex};\r
282   fEventCount->Fill(vEventCount);\r
283 \r
284   //\r
285   // cosmic background and splitted tracks\r
286   //\r
287   if(GetParticleMode() == AlidNdPtHelper::kBackgroundTrack) \r
288   {\r
289     AlidNdPtBackgroundCuts *backCuts = GetBackgroundCuts(); \r
290     if(!backCuts) return;\r
291 \r
292     for (Int_t iTrack = 0; iTrack < esdEvent->GetNumberOfTracks(); iTrack++) \r
293     {\r
294       AliESDtrack *track1 = esdEvent->GetTrack(iTrack);\r
295       if(!track1) continue; \r
296       if(track1->Charge()==0) continue; \r
297 \r
298       for (Int_t jTrack = iTrack+1; jTrack < esdEvent->GetNumberOfTracks(); jTrack++) \r
299       {\r
300         AliESDtrack *track2 = esdEvent->GetTrack(jTrack);\r
301         if(!track2) continue; \r
302         if(track2->Charge()==0) continue; \r
303 \r
304         //printf("track2->Charge() %d\n",track2->Charge());\r
305 \r
306         backCuts->IsBackgroundTrack(track1, track2);\r
307       }\r
308     }\r
309   }\r
310 \r
311   // check event cuts\r
312   if(isEventOK && isEventTriggered)\r
313   {\r
314     // get all charged tracks\r
315     allChargedTracks = AlidNdPtHelper::GetAllChargedTracks(esdEvent,GetAnalysisMode());\r
316     if(!allChargedTracks) return;\r
317 \r
318     Int_t entries = allChargedTracks->GetEntries();\r
319     for(Int_t i=0; i<entries;++i) \r
320     {\r
321       AliESDtrack *track = (AliESDtrack*)allChargedTracks->At(i);\r
322       if(!track) continue;\r
323 \r
324       //\r
325       Bool_t isOK = kFALSE;\r
326       Double_t x[3]; track->GetXYZ(x);\r
327       Double_t b[3]; AliTracker::GetBxByBz(x,b);\r
328 \r
329       //\r
330       // if TPC-ITS hybrid tracking (kTPCITSHybrid)\r
331       // replace track parameters with TPC-ony track parameters\r
332       //\r
333       if( GetAnalysisMode() == AlidNdPtHelper::kTPCITSHybrid ) \r
334       {\r
335         // Relate TPC-only tracks to SPD vertex\r
336         isOK = track->RelateToVertexTPCBxByBz(esdEvent->GetPrimaryVertexSPD(), b, kVeryBig);\r
337         if(!isOK) continue;\r
338 \r
339         // replace esd track parameters with TPCinner\r
340         AliExternalTrackParam  *tpcTrack  = new AliExternalTrackParam(*(track->GetTPCInnerParam()));\r
341         if (!tpcTrack) return;\r
342         track->Set(tpcTrack->GetX(),tpcTrack->GetAlpha(),tpcTrack->GetParameter(),tpcTrack->GetCovariance());\r
343 \r
344         if(tpcTrack) delete tpcTrack; \r
345       } \r
346 \r
347       //\r
348       if (GetAnalysisMode()==AlidNdPtHelper::kTPCSPDvtxUpdate) \r
349       {\r
350         //\r
351         // update track parameters\r
352         //\r
353         AliExternalTrackParam cParam;\r
354         isOK = track->RelateToVertexTPCBxByBz(esdEvent->GetPrimaryVertexSPD(), b, kVeryBig, &cParam);\r
355         if(!isOK) continue;\r
356 \r
357         track->Set(cParam.GetX(),cParam.GetAlpha(),cParam.GetParameter(),cParam.GetCovariance());\r
358       }\r
359 \r
360       FillHistograms(track, stack);\r
361       multAll++;\r
362     }\r
363 \r
364     Double_t vRecEventHist[5] = {vtxESD->GetXv(),vtxESD->GetYv(),vtxESD->GetZv(),vtxESD->GetZRes(),multAll};\r
365     fRecEventHist->Fill(vRecEventHist);\r
366 \r
367     if(IsUseMCInfo()) {\r
368       Double_t vRecMCEventHist[5] = {vtxESD->GetXv()-vtxMC[0],vtxESD->GetYv()-vtxMC[1],vtxESD->GetZv()-vtxMC[2],multAll};\r
369       fRecMCEventHist->Fill(vRecMCEventHist);\r
370     }\r
371   }\r
372 \r
373   if(allChargedTracks) delete allChargedTracks; allChargedTracks = 0;\r
374 \r
375 }\r
376 \r
377 //_____________________________________________________________________________\r
378 void AlidNdPtCutAnalysis::FillHistograms(AliESDtrack *const esdTrack, AliStack *const stack) const\r
379 {\r
380   //\r
381   // Fill ESD track and MC histograms \r
382   //\r
383   if(!esdTrack) return;\r
384   if(esdTrack->Charge() == 0.) return;\r
385 \r
386   Float_t pt = esdTrack->Pt();\r
387   Float_t eta = esdTrack->Eta();\r
388   Float_t phi = esdTrack->Phi();\r
389   //Int_t nClust = esdTrack->GetTPCclusters(0);\r
390   Int_t nClust = esdTrack->GetTPCNclsIter1();\r
391   Int_t nFindableClust = esdTrack->GetTPCNclsF();\r
392 \r
393   Float_t chi2PerCluster = 0.;\r
394   //if(nClust>0.) chi2PerCluster = esdTrack->GetTPCchi2()/Float_t(nClust);\r
395   if(nClust>0.) chi2PerCluster = esdTrack->GetTPCchi2Iter1()/Float_t(nClust);\r
396 \r
397   Float_t clustPerFindClust = 0.;\r
398   if(nFindableClust>0.) clustPerFindClust = Float_t(nClust)/nFindableClust;\r
399 \r
400   Float_t b[2], bCov[3];\r
401   esdTrack->GetImpactParameters(b,bCov);\r
402 \r
403   // kink daughter\r
404   Bool_t isKink = kFALSE;\r
405   if(esdTrack->GetKinkIndex(0)>0) isKink=kTRUE;\r
406 \r
407   //\r
408   // Fill rec vs MC information\r
409   //\r
410 \r
411   Bool_t isPrim = kTRUE;\r
412 \r
413   if(IsUseMCInfo()) {\r
414     if(!stack) return;\r
415     Int_t label = TMath::Abs(esdTrack->GetLabel()); \r
416     TParticle* particle = stack->Particle(label);\r
417     if(!particle) return;\r
418     if(particle->GetPDG() && particle->GetPDG()->Charge()==0.) return;\r
419     isPrim = stack->IsPhysicalPrimary(label);\r
420   }\r
421 \r
422   // fill histo\r
423   Int_t polarity = -2;\r
424   if (esdTrack->Charge() < 0.) polarity = 0; \r
425   else polarity = 1; \r
426   Double_t vRecMCTrackHist[11] = { nClust,chi2PerCluster,clustPerFindClust,b[0],b[1],eta,phi,pt,isKink,isPrim, polarity }; \r
427   fRecMCTrackHist->Fill(vRecMCTrackHist);\r
428 }\r
429 \r
430 //_____________________________________________________________________________\r
431 Long64_t AlidNdPtCutAnalysis::Merge(TCollection* const list) \r
432 {\r
433   // Merge list of objects (needed by PROOF)\r
434 \r
435   if (!list)\r
436   return 0;\r
437 \r
438   if (list->IsEmpty())\r
439   return 1;\r
440 \r
441   TIterator* iter = list->MakeIterator();\r
442   TObject* obj = 0;\r
443 \r
444   //TList *collPhysSelection = new TList;\r
445 \r
446   // collection of generated histograms\r
447   Int_t count=0;\r
448   while((obj = iter->Next()) != 0) {\r
449     AlidNdPtCutAnalysis* entry = dynamic_cast<AlidNdPtCutAnalysis*>(obj);\r
450     if (entry == 0) continue; \r
451   \r
452     // event histo\r
453     fEventCount->Add(entry->fEventCount);\r
454     fRecEventHist->Add(entry->fRecEventHist);\r
455     fRecMCEventHist->Add(entry->fRecMCEventHist);\r
456     fMCEventHist->Add(entry->fMCEventHist);\r
457 \r
458     // track histo\r
459     fRecMCTrackHist->Add(entry->fRecMCTrackHist);\r
460 \r
461     // physics selection\r
462     //collPhysSelection->Add(entry->GetPhysicsTriggerSelection());\r
463     \r
464   count++;\r
465   }\r
466 \r
467   //AliPhysicsSelection *trigSelection = GetPhysicsTriggerSelection();\r
468   //trigSelection->Merge(collPhysSelection);\r
469 \r
470   //if(collPhysSelection) delete collPhysSelection;\r
471 \r
472 return count;\r
473 }\r
474 \r
475 //_____________________________________________________________________________\r
476 void AlidNdPtCutAnalysis::Analyse() \r
477 {\r
478   //\r
479   // Analyse histograms\r
480   //\r
481   TH1::AddDirectory(kFALSE);\r
482   TObjArray *aFolderObj = new TObjArray;\r
483   TH1D *h1D = 0; \r
484   TH2D *h2D = 0; \r
485 \r
486 \r
487   //\r
488   // get cuts\r
489   //\r
490   AlidNdPtEventCuts *evtCuts = GetEventCuts(); \r
491   AlidNdPtAcceptanceCuts *accCuts = GetAcceptanceCuts(); \r
492   AliESDtrackCuts *esdTrackCuts = GetTrackCuts(); \r
493 \r
494   if(!evtCuts || !accCuts || !esdTrackCuts) {\r
495     Error("AlidNdPtCutAnalysis::Analyse()", "cuts not available");\r
496     return;\r
497   }\r
498 \r
499   //\r
500   // set min and max values\r
501   //\r
502   Double_t minPt = accCuts->GetMinPt();\r
503   Double_t maxPt = accCuts->GetMaxPt();\r
504   Double_t minEta = accCuts->GetMinEta();\r
505   Double_t maxEta = accCuts->GetMaxEta()-0.00001;\r
506 \r
507   Double_t maxDCAr = accCuts->GetMaxDCAr();\r
508 \r
509   //\r
510   // Event counters\r
511   //\r
512   h2D = (TH2D*)fEventCount->Projection(0,1);\r
513   h2D->SetName("trig_vs_trigANDvertex");\r
514   aFolderObj->Add(h2D);\r
515 \r
516   fEventCount->GetAxis(0)->SetRange(2,2); // triggered\r
517   h1D = (TH1D*)fEventCount->Projection(1);\r
518   h1D->SetTitle("rec. vertex for triggered events");\r
519   h1D->SetName("trigANDvertex");\r
520   aFolderObj->Add(h1D);\r
521 \r
522   //\r
523   // Create rec. event histograms\r
524   //\r
525   h1D = (TH1D *)fRecEventHist->Projection(0);\r
526   h1D->SetName("rec_xv");\r
527   aFolderObj->Add(h1D);\r
528 \r
529   h1D = (TH1D *)fRecEventHist->Projection(1);\r
530   h1D->SetName("rec_yv");\r
531   aFolderObj->Add(h1D);\r
532 \r
533   h1D = (TH1D *)fRecEventHist->Projection(2);\r
534   h1D->SetName("rec_zv");\r
535   aFolderObj->Add(h1D);\r
536 \r
537   h2D = (TH2D *)fRecEventHist->Projection(3,4);\r
538   h2D->SetName("rec_resZv_vs_Mult");\r
539   aFolderObj->Add(h2D);\r
540 \r
541   h2D = (TH2D *)fRecEventHist->Projection(0,1);\r
542   h2D->SetName("rec_xv_vs_yv");\r
543   aFolderObj->Add(h2D);\r
544 \r
545   h2D = (TH2D *)fRecEventHist->Projection(0,2);\r
546   h2D->SetName("rec_xv_vs_zv");\r
547   aFolderObj->Add(h2D);\r
548 \r
549   h2D = (TH2D *)fRecEventHist->Projection(3,4);\r
550   h2D->SetName("rec_resZv_vs_Mult");\r
551   aFolderObj->Add(h2D);\r
552 \r
553   //\r
554   // MC available\r
555   //\r
556   if(IsUseMCInfo()) {\r
557 \r
558   //\r
559   // Create mc event histograms\r
560   //\r
561   h2D = (TH2D *)fMCEventHist->Projection(0,1);\r
562   h2D->SetName("mc_xv_vs_yv");\r
563   aFolderObj->Add(h2D);\r
564 \r
565   h2D = (TH2D *)fMCEventHist->Projection(0,2);\r
566   h2D->SetName("mc_xv_vs_zv");\r
567   aFolderObj->Add(h2D);\r
568 \r
569   //\r
570   // Create rec-mc event histograms\r
571   //\r
572   h2D = (TH2D *)fRecMCEventHist->Projection(0,3);\r
573   h2D->SetName("rec_mc_deltaXv_vs_mult");\r
574   aFolderObj->Add(h2D);\r
575 \r
576   h2D = (TH2D *)fRecMCEventHist->Projection(1,3);\r
577   h2D->SetName("rec_mc_deltaYv_vs_mult");\r
578   aFolderObj->Add(h2D);\r
579 \r
580   h2D = (TH2D *)fRecMCEventHist->Projection(2,3);\r
581   h2D->SetName("rec_mc_deltaZv_vs_mult");\r
582   aFolderObj->Add(h2D);\r
583 \r
584   } // end use MC info \r
585 \r
586 \r
587 \r
588   //\r
589   // Create rec-mc track track histograms \r
590   //\r
591 \r
592   // DCA cuts\r
593   fRecMCTrackHist->GetAxis(3)->SetRangeUser(-maxDCAr,maxDCAr);\r
594   fRecMCTrackHist->GetAxis(4)->SetRangeUser(-maxDCAr,maxDCAr);\r
595 \r
596   h2D = (TH2D *)fRecMCTrackHist->Projection(7,5);\r
597   h2D->SetName("pt_vs_eta");\r
598   aFolderObj->Add(h2D);\r
599 \r
600   fRecMCTrackHist->GetAxis(7)->SetRangeUser(minPt,maxPt);  \r
601 \r
602   h2D = (TH2D *)fRecMCTrackHist->Projection(0,5);\r
603   h2D->SetName("nClust_vs_eta");\r
604   aFolderObj->Add(h2D);\r
605 \r
606   h2D = (TH2D *)fRecMCTrackHist->Projection(1,5);\r
607   h2D->SetName("chi2PerClust_vs_eta");\r
608   aFolderObj->Add(h2D);\r
609 \r
610   h2D = (TH2D *)fRecMCTrackHist->Projection(2,5);\r
611   h2D->SetName("ratio_nClust_nFindableClust_vs_eta");\r
612   aFolderObj->Add(h2D);\r
613 \r
614   h2D = (TH2D *)fRecMCTrackHist->Projection(5,6);\r
615   h2D->SetName("eta_vs_phi");\r
616   aFolderObj->Add(h2D);\r
617 \r
618   //\r
619   fRecMCTrackHist->GetAxis(5)->SetRangeUser(minEta,maxEta);  \r
620 \r
621   h2D = (TH2D *)fRecMCTrackHist->Projection(0,6);\r
622   h2D->SetName("nClust_vs_phi");\r
623   aFolderObj->Add(h2D);\r
624 \r
625   h2D = (TH2D *)fRecMCTrackHist->Projection(1,6);\r
626   h2D->SetName("chi2PerClust_vs_phi");\r
627   aFolderObj->Add(h2D);\r
628 \r
629   h2D = (TH2D *)fRecMCTrackHist->Projection(2,6);\r
630   h2D->SetName("ratio_nClust_nFindableClust_vs_phi");\r
631   aFolderObj->Add(h2D);\r
632 \r
633   //\r
634   fRecMCTrackHist->GetAxis(7)->SetRangeUser(0.0,maxPt);  \r
635 \r
636   h2D = (TH2D *)fRecMCTrackHist->Projection(0,7);\r
637   h2D->SetName("nClust_vs_pt");\r
638   aFolderObj->Add(h2D);\r
639 \r
640   h2D = (TH2D *)fRecMCTrackHist->Projection(1,7);\r
641   h2D->SetName("chi2PerClust_vs_pt");\r
642   aFolderObj->Add(h2D);\r
643 \r
644   h2D = (TH2D *)fRecMCTrackHist->Projection(2,7);\r
645   h2D->SetName("ratio_nClust_nFindableClust_vs_pt");\r
646   aFolderObj->Add(h2D);\r
647 \r
648   h2D = (TH2D *)fRecMCTrackHist->Projection(6,7);\r
649   h2D->SetName("phi_vs_pt");\r
650   aFolderObj->Add(h2D);\r
651 \r
652 \r
653   // fiducial volume\r
654   fRecMCTrackHist->GetAxis(5)->SetRangeUser(minEta,maxEta);  \r
655   fRecMCTrackHist->GetAxis(7)->SetRangeUser(minPt,maxPt);  \r
656 \r
657   // DCA cuts\r
658   fRecMCTrackHist->GetAxis(3)->SetRangeUser(-maxDCAr,maxDCAr);\r
659   fRecMCTrackHist->GetAxis(4)->SetRangeUser(-maxDCAr,maxDCAr);\r
660 \r
661   h2D = (TH2D *)fRecMCTrackHist->Projection(0,1);\r
662   h2D->SetName("nClust_vs_chi2PerClust");\r
663   aFolderObj->Add(h2D);\r
664 \r
665   h2D = (TH2D *)fRecMCTrackHist->Projection(0,2);\r
666   h2D->SetName("nClust_vs_ratio_nClust_nFindableClust");\r
667   aFolderObj->Add(h2D);\r
668 \r
669   //\r
670   // DCAy cuts\r
671   //\r
672   fRecMCTrackHist->GetAxis(0)->SetRange(50,160); // nClust/track > 50\r
673   fRecMCTrackHist->GetAxis(1)->SetRangeUser(0.,3.9999); // chi2/cluster < 4.0\r
674   fRecMCTrackHist->GetAxis(3)->SetRange(1,fRecMCTrackHist->GetAxis(3)->GetNbins());\r
675   //fRecMCTrackHist->GetAxis(4)->SetRangeUser(-1.0,1.0);\r
676   fRecMCTrackHist->GetAxis(4)->SetRange(1,fRecMCTrackHist->GetAxis(4)->GetNbins());\r
677 \r
678   // sec\r
679   fRecMCTrackHist->GetAxis(9)->SetRange(1,1);\r
680   h1D = (TH1D *)fRecMCTrackHist->Projection(3);\r
681   h1D->SetName("dcay_sec");\r
682   aFolderObj->Add(h1D);\r
683 \r
684   // prim\r
685   fRecMCTrackHist->GetAxis(9)->SetRange(2,2);\r
686   h1D = (TH1D *)fRecMCTrackHist->Projection(3);\r
687   h1D->SetName("dcay_prim");\r
688   aFolderObj->Add(h1D);\r
689 \r
690   // DCAz cuts\r
691   //fRecMCTrackHist->GetAxis(3)->SetRangeUser(-1.0,1.0);\r
692   fRecMCTrackHist->GetAxis(4)->SetRange(1,fRecMCTrackHist->GetAxis(4)->GetNbins());\r
693 \r
694   // sec\r
695   fRecMCTrackHist->GetAxis(9)->SetRange(1,1);\r
696   h1D = (TH1D *)fRecMCTrackHist->Projection(4);\r
697   h1D->SetName("dcaz_sec");\r
698   aFolderObj->Add(h1D);\r
699 \r
700   // prim\r
701   fRecMCTrackHist->GetAxis(9)->SetRange(2,2);\r
702   h1D = (TH1D *)fRecMCTrackHist->Projection(4);\r
703   h1D->SetName("dcaz_prim");\r
704   aFolderObj->Add(h1D);\r
705 \r
706 \r
707   // export objects to analysis folder\r
708   fAnalysisFolder = ExportToFolder(aFolderObj);\r
709 \r
710   // delete only TObjArray\r
711   if(aFolderObj) delete aFolderObj;\r
712 }\r
713 \r
714 //_____________________________________________________________________________\r
715 TFolder* AlidNdPtCutAnalysis::ExportToFolder(TObjArray * const array) \r
716 {\r
717   // recreate folder avery time and export objects to new one\r
718   //\r
719   AlidNdPtCutAnalysis * comp=this;\r
720   TFolder *folder = comp->GetAnalysisFolder();\r
721 \r
722   TString name, title;\r
723   TFolder *newFolder = 0;\r
724   Int_t i = 0;\r
725   Int_t size = array->GetSize();\r
726 \r
727   if(folder) { \r
728      // get name and title from old folder\r
729      name = folder->GetName();  \r
730      title = folder->GetTitle();  \r
731 \r
732          // delete old one\r
733      delete folder;\r
734 \r
735          // create new one\r
736      newFolder = CreateFolder(name.Data(),title.Data());\r
737      newFolder->SetOwner();\r
738 \r
739          // add objects to folder\r
740      while(i < size) {\r
741            newFolder->Add(array->At(i));\r
742            i++;\r
743          }\r
744   }\r
745 \r
746 return newFolder;\r
747 }\r
748 \r
749 //_____________________________________________________________________________\r
750 TFolder* AlidNdPtCutAnalysis::CreateFolder(TString name,TString title) { \r
751 // create folder for analysed histograms\r
752 //\r
753 TFolder *folder = 0;\r
754   folder = new TFolder(name.Data(),title.Data());\r
755 \r
756   return folder;\r
757 }\r