]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWG0/dNdPt/AlidNdPtCutAnalysisPbPb.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
coverity
[u/mrichter/AliRoot.git] / PWG0 / dNdPt / AlidNdPtCutAnalysisPbPb.cxx
CommitLineData
dac20c1b 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// AlidNdPtCutAnalysisPbPb 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
10360c58 33#include "AliInputEventHandler.h" \r
34#include "AliCentrality.h" \r
dac20c1b 35#include "AliStack.h" \r
36#include "AliESDEvent.h" \r
37#include "AliMCEvent.h" \r
38#include "AliESDtrackCuts.h" \r
39#include "AliLog.h" \r
40#include "AliTracker.h" \r
41\r
42#include "AlidNdPtEventCuts.h"\r
43#include "AlidNdPtAcceptanceCuts.h"\r
44#include "AlidNdPtBackgroundCuts.h"\r
45#include "AlidNdPtAnalysis.h"\r
46#include "AliPhysicsSelection.h"\r
47\r
48#include "AliPWG0Helper.h"\r
49#include "AlidNdPtHelper.h"\r
50#include "AlidNdPtCutAnalysisPbPb.h"\r
51\r
52using namespace std;\r
53\r
54ClassImp(AlidNdPtCutAnalysisPbPb)\r
55\r
56//_____________________________________________________________________________\r
57 AlidNdPtCutAnalysisPbPb::AlidNdPtCutAnalysisPbPb(): AlidNdPt(),\r
58 fAnalysisFolder(0),\r
59 fEventCount(0),\r
60 fRecEventHist(0),\r
61 fMCEventHist(0),\r
62 fRecMCEventHist(0),\r
10360c58 63 fRecMCTrackHist(0),\r
64 fCentralityEstimator(0)\r
dac20c1b 65{\r
66 // default constructor\r
67 Init();\r
68}\r
69\r
70//_____________________________________________________________________________\r
71AlidNdPtCutAnalysisPbPb::AlidNdPtCutAnalysisPbPb(Char_t* name, Char_t* title): AlidNdPt(name,title),\r
72 fAnalysisFolder(0),\r
73 fEventCount(0),\r
74 fRecEventHist(0),\r
75 fMCEventHist(0),\r
76 fRecMCEventHist(0),\r
10360c58 77 fRecMCTrackHist(0),\r
78 fCentralityEstimator(0)\r
dac20c1b 79{\r
80 // constructor\r
81 Init();\r
82}\r
83\r
84//_____________________________________________________________________________\r
85AlidNdPtCutAnalysisPbPb::~AlidNdPtCutAnalysisPbPb() {\r
86 // \r
87 if(fEventCount) delete fEventCount; fEventCount=0;\r
88 if(fRecEventHist) delete fRecEventHist; fRecEventHist=0;\r
89 if(fMCEventHist) delete fMCEventHist; fMCEventHist=0;\r
90 if(fRecMCEventHist) delete fRecMCEventHist; fRecMCEventHist=0;\r
91 if(fRecMCTrackHist) delete fRecMCTrackHist; fRecMCTrackHist=0;\r
92\r
93 if(fAnalysisFolder) delete fAnalysisFolder; fAnalysisFolder=0;\r
94}\r
95\r
96//_____________________________________________________________________________\r
97void AlidNdPtCutAnalysisPbPb::Init(){\r
98 //\r
99 // Init histograms\r
100 //\r
10360c58 101 //const Int_t ptNbins = 58; \r
102 //const Double_t ptMin = 0.; \r
103 //const Double_t ptMax = 20.; \r
104 /*\r
dffe2c06 105 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,18.0, 20.};\r
10360c58 106 */\r
107\r
108 // set pt bins\r
109 const Int_t ptNbins = 50;\r
110 const Double_t ptMin = 1.e-2, ptMax = 50.;\r
111 Double_t *binsPt = CreateLogAxis(ptNbins,ptMin,ptMax);\r
112\r
113 // centrality bins\r
114 const Int_t centrNbins = 3;\r
115 const Double_t centrMin = 0, centrMax = 1;\r
116 Double_t binsCentr[centrNbins+1] = {0.0, 0.2, 0.5, 1.0};\r
dac20c1b 117\r
118 // \r
119 Int_t binsEventCount[2]={2,2};\r
120 Double_t minEventCount[2]={0,0}; \r
121 Double_t maxEventCount[2]={2,2}; \r
122 fEventCount = new THnSparseF("fEventCount","trig vs trig+vertex",2,binsEventCount,minEventCount,maxEventCount);\r
123 fEventCount->GetAxis(0)->SetTitle("trig");\r
124 fEventCount->GetAxis(1)->SetTitle("trig+vert");\r
125 fEventCount->Sumw2();\r
126\r
127 //Xv:Yv:Zv:ResZv:Mult\r
10360c58 128 Double_t kFact = 0.1;\r
dac20c1b 129\r
130 Int_t binsRecEventHist[5]={80,80,100,80,150};\r
10360c58 131 Double_t minRecEventHist[5]={-1.*kFact,-1.*kFact,-35.,0.,0.}; \r
132 Double_t maxRecEventHist[5]={1.*kFact,1.*kFact,35.,10.,150.}; \r
dac20c1b 133 fRecEventHist = new THnSparseF("fRecEventHist","Xv:Yv:Zv:ResZv:Mult",5,binsRecEventHist,minRecEventHist,maxRecEventHist);\r
134 fRecEventHist->GetAxis(0)->SetTitle("Xv (cm)");\r
135 fRecEventHist->GetAxis(1)->SetTitle("Yv (cm)");\r
136 fRecEventHist->GetAxis(2)->SetTitle("Zv (cm)");\r
137 fRecEventHist->GetAxis(3)->SetTitle("ResZv (cm)");\r
138 fRecEventHist->GetAxis(4)->SetTitle("Mult");\r
139 fRecEventHist->Sumw2();\r
140\r
141 //Xv:Yv:Zv\r
142 Int_t binsMCEventHist[3]={80,80,100};\r
143 Double_t minMCEventHist[3]={-0.1,-0.1,-35.}; \r
144 Double_t maxMCEventHist[3]={0.1,0.1,35.}; \r
145 fMCEventHist = new THnSparseF("fMCEventHist","mcXv:mcYv:mcZv",3,binsMCEventHist,minMCEventHist,maxMCEventHist);\r
146 fMCEventHist->GetAxis(0)->SetTitle("mcXv (cm)");\r
147 fMCEventHist->GetAxis(1)->SetTitle("mcYv (cm)");\r
148 fMCEventHist->GetAxis(2)->SetTitle("mcZv (cm)");\r
149 fMCEventHist->Sumw2();\r
150\r
151 //Xv-mcXv:Yv-mcYv:Zv-mcZv:Mult\r
152 Int_t binsRecMCEventHist[4]={100,100,100,150};\r
153 Double_t minRecMCEventHist[4]={-1.0*kFact,-1.0*kFact,-1.0*kFact,0.}; \r
154 Double_t maxRecMCEventHist[4]={1.0*kFact,1.0*kFact,1.0*kFact,150.}; \r
155 fRecMCEventHist = new THnSparseF("fRecMCEventHist","Xv-mcXv:Yv-mcYv:Zv-mcZv:Mult",4,binsRecMCEventHist,minRecMCEventHist,maxRecMCEventHist);\r
156 fRecMCEventHist->GetAxis(0)->SetTitle("Xv-mcXv (cm)");\r
157 fRecMCEventHist->GetAxis(1)->SetTitle("Yv-mcYv (cm)");\r
158 fRecMCEventHist->GetAxis(2)->SetTitle("Zv-mcZv (cm)");\r
159 fRecMCEventHist->GetAxis(3)->SetTitle("Mult");\r
160 fRecMCEventHist->Sumw2();\r
161\r
162 //\r
163 // THnSparse track histograms\r
164 //\r
10360c58 165 //\r
166 //\r
167 //nCrossRows:chi2PerClust:nCrossRows/nFindableClust:fracSharedClust:DCAy:DCAz:eta:phi:pt:hasStrangeMother:isFromMaterial:isPrim:charge:centr\r
168 Int_t binsRecMCTrackHist[14]= {160, 10, 20, 20, 50, 50, 20, 90, ptNbins, 2, 2, 2, 3, centrNbins};\r
169 Double_t minRecMCTrackHist[14]={0., 0., 0., 0.,-0.5,-0.5,-1.0, 0., ptMin, 0., 0., 0.,-1., centrMin};\r
98af3b87 170 Double_t maxRecMCTrackHist[14]={160.,10.,1., 1., 0.5, 0.5, 1.0, 2.*TMath::Pi(), ptMax, 2., 2., 2., 2., centrMax};\r
10360c58 171\r
172 fRecMCTrackHist = new THnSparseF("fRecMCTrackHist","nCrossRows:chi2PerClust:nCrossRows/nFindableClust:fracSharedClust:DCAy:DCAz:eta:phi:pt:hasStrangeMother:isFromMaterial:isPrim:charge:centr",14,binsRecMCTrackHist,minRecMCTrackHist,maxRecMCTrackHist);\r
173 fRecMCTrackHist->SetBinEdges(8,binsPt);\r
174 fRecMCTrackHist->SetBinEdges(13,binsCentr);\r
175\r
176 fRecMCTrackHist->GetAxis(0)->SetTitle("nCrossRows");\r
177 fRecMCTrackHist->GetAxis(1)->SetTitle("chi2PerClust");\r
178 fRecMCTrackHist->GetAxis(2)->SetTitle("nCrossRows/nFindableClust");\r
179 fRecMCTrackHist->GetAxis(3)->SetTitle("fracSharedClust");\r
180 fRecMCTrackHist->GetAxis(4)->SetTitle("DCAy (cm)");\r
181 fRecMCTrackHist->GetAxis(5)->SetTitle("DCAz (cm)");\r
182 fRecMCTrackHist->GetAxis(6)->SetTitle("#eta");\r
183 fRecMCTrackHist->GetAxis(7)->SetTitle("#phi (rad)");\r
184 fRecMCTrackHist->GetAxis(8)->SetTitle("p_{T} (GeV/c)");\r
185 fRecMCTrackHist->GetAxis(9)->SetTitle("hasStrangeMother");\r
186 fRecMCTrackHist->GetAxis(10)->SetTitle("isFromMaterial");\r
187 fRecMCTrackHist->GetAxis(11)->SetTitle("isPrim");\r
188 fRecMCTrackHist->GetAxis(12)->SetTitle("charge");\r
189 fRecMCTrackHist->GetAxis(13)->SetTitle("centrality");\r
190 fRecMCTrackHist->Sumw2();\r
191\r
192\r
193 //nClust:chi2PerClust:nClust/nFindableClust:DCAy:DCAz:eta:phi:pt:hasStrangeMother:isFromMaterial:isPrim:charge:centr\r
194 /*\r
195 Int_t binsRecMCTrackHist[13]={160,20,40,50,50,20,90,ptNbins,2,2,2,3,centrNbins};\r
196 Double_t minRecMCTrackHist[13]={0., 0., 0., -0.5,-0.5,-1.0, 0., ptMin, 0., 0., 0., -1., centrMin};\r
197 Double_t maxRecMCTrackHist[13]={160.,10.,1.2, 0.5,0.5,1.0, 2.*TMath::Pi(), ptMax, 2., 2., 2., 2.,centrMax};\r
dac20c1b 198\r
6a01c05f 199 fRecMCTrackHist = new THnSparseF("fRecMCTrackHist","nClust:chi2PerClust:nClust/nFindableClust:DCAy:DCAz:eta:phi:pt:hasStrangeMother:isFromMaterial:isPrim:charge",12,binsRecMCTrackHist,minRecMCTrackHist,maxRecMCTrackHist);\r
dac20c1b 200 fRecMCTrackHist->SetBinEdges(7,binsPt);\r
10360c58 201 fRecMCTrackHist->SetBinEdges(12,binsCentr);\r
dac20c1b 202\r
203 fRecMCTrackHist->GetAxis(0)->SetTitle("nClust");\r
204 fRecMCTrackHist->GetAxis(1)->SetTitle("chi2PerClust");\r
205 fRecMCTrackHist->GetAxis(2)->SetTitle("nClust/nFindableClust");\r
206 fRecMCTrackHist->GetAxis(3)->SetTitle("DCAy (cm)");\r
207 fRecMCTrackHist->GetAxis(4)->SetTitle("DCAz (cm)");\r
208 fRecMCTrackHist->GetAxis(5)->SetTitle("#eta");\r
209 fRecMCTrackHist->GetAxis(6)->SetTitle("#phi (rad)");\r
210 fRecMCTrackHist->GetAxis(7)->SetTitle("p_{T} (GeV/c)");\r
e84b59cb 211 fRecMCTrackHist->GetAxis(8)->SetTitle("hasStrangeMother");\r
6a01c05f 212 fRecMCTrackHist->GetAxis(9)->SetTitle("isFromMaterial");\r
213 fRecMCTrackHist->GetAxis(10)->SetTitle("isPrim");\r
214 fRecMCTrackHist->GetAxis(11)->SetTitle("charge");\r
10360c58 215 fRecMCTrackHist->GetAxis(12)->SetTitle("centrality");\r
dac20c1b 216 fRecMCTrackHist->Sumw2();\r
10360c58 217 */\r
dac20c1b 218\r
219 // init output folder\r
220 fAnalysisFolder = CreateFolder("folderdNdPt","Analysis dNdPt Folder");\r
221\r
222}\r
223\r
224//_____________________________________________________________________________\r
225void AlidNdPtCutAnalysisPbPb::Process(AliESDEvent *const esdEvent, AliMCEvent * const mcEvent)\r
226{\r
227 //\r
228 // Process real and/or simulated events\r
229 //\r
230 if(!esdEvent) {\r
231 AliDebug(AliLog::kError, "esdEvent not available");\r
232 return;\r
233 }\r
234\r
235 // get selection cuts\r
236 AlidNdPtEventCuts *evtCuts = GetEventCuts(); \r
237 AlidNdPtAcceptanceCuts *accCuts = GetAcceptanceCuts(); \r
238 AliESDtrackCuts *esdTrackCuts = GetTrackCuts(); \r
239\r
240 if(!evtCuts || !accCuts || !esdTrackCuts) {\r
241 AliDebug(AliLog::kError, "cuts not available");\r
242 return;\r
243 }\r
244\r
245 // trigger selection\r
10360c58 246\r
dac20c1b 247 Bool_t isEventTriggered = kTRUE;\r
10360c58 248 AliPhysicsSelection *physicsSelection = NULL;\r
249 AliTriggerAnalysis* triggerAnalysis = NULL;\r
250\r
251 // \r
252 AliInputEventHandler* inputHandler = (AliInputEventHandler*) AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler();\r
253 if (!inputHandler)\r
254 {\r
255 Printf("ERROR: Could not receive input handler");\r
256 return;\r
257 }\r
dac20c1b 258\r
259 if(evtCuts->IsTriggerRequired()) \r
260 {\r
10360c58 261 // always MB\r
262 isEventTriggered = inputHandler->IsEventSelected() & AliVEvent::kMB;\r
263\r
264 physicsSelection = static_cast<AliPhysicsSelection*> (inputHandler->GetEventSelection());\r
265 if(!physicsSelection) return;\r
266 //SetPhysicsTriggerSelection(physicsSelection);\r
267\r
268 if (isEventTriggered && (GetTrigger() == AliTriggerAnalysis::kV0AND)) {\r
269 // set trigger (V0AND)\r
270 triggerAnalysis = physicsSelection->GetTriggerAnalysis();\r
271 if(!triggerAnalysis) return;\r
272 isEventTriggered = triggerAnalysis->IsOfflineTriggerFired(esdEvent, GetTrigger());\r
dac20c1b 273 }\r
10360c58 274 }\r
dac20c1b 275\r
10360c58 276 // centrality determination\r
277 Float_t centralityF = -1;\r
278 AliCentrality *esdCentrality = esdEvent->GetCentrality();\r
279 centralityF = esdCentrality->GetCentralityPercentile(fCentralityEstimator.Data());\r
dac20c1b 280\r
281 // use MC information\r
282 AliHeader* header = 0;\r
283 AliGenEventHeader* genHeader = 0;\r
284 AliStack* stack = 0;\r
285 TArrayF vtxMC(3);\r
286\r
287 Int_t multMCTrueTracks = 0;\r
288 if(IsUseMCInfo())\r
289 {\r
290 //\r
291 if(!mcEvent) {\r
292 AliDebug(AliLog::kError, "mcEvent not available");\r
293 return;\r
294 }\r
295 // get MC event header\r
296 header = mcEvent->Header();\r
297 if (!header) {\r
298 AliDebug(AliLog::kError, "Header not available");\r
299 return;\r
300 }\r
301 // MC particle stack\r
302 stack = mcEvent->Stack();\r
303 if (!stack) {\r
304 AliDebug(AliLog::kError, "Stack not available");\r
305 return;\r
306 }\r
307\r
308 // get MC vertex\r
309 genHeader = header->GenEventHeader();\r
310 if (!genHeader) {\r
311 AliDebug(AliLog::kError, "Could not retrieve genHeader from Header");\r
312 return;\r
313 }\r
314 genHeader->PrimaryVertex(vtxMC);\r
315\r
316 Double_t vMCEventHist[3]={vtxMC[0],vtxMC[1],vtxMC[2]};\r
317 fMCEventHist->Fill(vMCEventHist);\r
318\r
319 // multipliticy of all MC primary tracks\r
320 // in Zv, pt and eta ranges)\r
321 multMCTrueTracks = AlidNdPtHelper::GetMCTrueTrackMult(mcEvent,evtCuts,accCuts);\r
322\r
323 } // end bUseMC\r
324\r
325 // get reconstructed vertex \r
326 const AliESDVertex* vtxESD = 0; \r
62e3b4b6 327 if(GetAnalysisMode() == AlidNdPtHelper::kTPC) {\r
dffe2c06 328 vtxESD = esdEvent->GetPrimaryVertexTPC();\r
62e3b4b6 329 }\r
330 else if(GetAnalysisMode() == AlidNdPtHelper::kTPCITS) {\r
331 vtxESD = esdEvent->GetPrimaryVertexTracks();\r
332 }\r
333 else {\r
82292e7f 334 return;\r
dac20c1b 335 }\r
62e3b4b6 336\r
337 if(!vtxESD) return;\r
338\r
82292e7f 339 Bool_t isEventOK = evtCuts->AcceptEvent(esdEvent,mcEvent,vtxESD); \r
dac20c1b 340 //printf("isEventOK %d, isEventTriggered %d \n",isEventOK, isEventTriggered);\r
341 //printf("GetAnalysisMode() %d \n",GetAnalysisMode());\r
342\r
343 TObjArray *allChargedTracks=0;\r
344 Int_t multAll=0;\r
345 \r
346 //\r
347 // event counter\r
348 // \r
349 //printf("isEventOK %d, isEventTriggered %d \n",isEventOK,isEventTriggered);\r
350\r
351 Bool_t isTrigAndVertex = isEventTriggered && isEventOK;\r
352 Double_t vEventCount[2] = { isEventTriggered, isTrigAndVertex};\r
353 fEventCount->Fill(vEventCount);\r
354\r
355\r
356 // check event cuts\r
357 if(isEventOK && isEventTriggered)\r
358 {\r
359 // get all charged tracks\r
360 allChargedTracks = AlidNdPtHelper::GetAllChargedTracks(esdEvent,GetAnalysisMode());\r
361 if(!allChargedTracks) return;\r
362\r
363 Int_t entries = allChargedTracks->GetEntries();\r
364 for(Int_t i=0; i<entries;++i) \r
365 {\r
366 AliESDtrack *track = (AliESDtrack*)allChargedTracks->At(i);\r
367 if(!track) continue;\r
368\r
82292e7f 369 if(!esdTrackCuts->AcceptTrack(track)) continue;\r
10360c58 370 FillHistograms(track, stack, centralityF);\r
dac20c1b 371 multAll++;\r
372 }\r
373\r
374 Double_t vRecEventHist[5] = {vtxESD->GetXv(),vtxESD->GetYv(),vtxESD->GetZv(),vtxESD->GetZRes(),multAll};\r
375 fRecEventHist->Fill(vRecEventHist);\r
376\r
377 if(IsUseMCInfo()) {\r
378 Double_t vRecMCEventHist[5] = {vtxESD->GetXv()-vtxMC[0],vtxESD->GetYv()-vtxMC[1],vtxESD->GetZv()-vtxMC[2],multAll};\r
379 fRecMCEventHist->Fill(vRecMCEventHist);\r
380 }\r
381 }\r
382\r
383 if(allChargedTracks) delete allChargedTracks; allChargedTracks = 0;\r
384\r
385}\r
386\r
387//_____________________________________________________________________________\r
10360c58 388void AlidNdPtCutAnalysisPbPb::FillHistograms(AliESDtrack *const esdTrack, AliStack *const stack, Float_t centralityF) const\r
dac20c1b 389{\r
390 //\r
391 // Fill ESD track and MC histograms \r
392 //\r
393 if(!esdTrack) return;\r
394 if(esdTrack->Charge() == 0.) return;\r
395\r
396 Float_t pt = esdTrack->Pt();\r
397 Float_t eta = esdTrack->Eta();\r
398 Float_t phi = esdTrack->Phi();\r
d1c6e451 399\r
400 Int_t nClust = 0;\r
401 if(GetAnalysisMode() == AlidNdPtHelper::kTPC) { \r
402 nClust = esdTrack->GetTPCNclsIter1();\r
403 } else {\r
404 nClust = esdTrack->GetTPCclusters(0);\r
405 }\r
dac20c1b 406\r
407 Float_t chi2PerCluster = 0.;\r
d1c6e451 408 if(GetAnalysisMode() == AlidNdPtHelper::kTPC) { \r
409 if(nClust>0.) chi2PerCluster = esdTrack->GetTPCchi2Iter1()/Float_t(nClust);\r
410 } else {\r
411 chi2PerCluster = esdTrack->GetTPCchi2()/Float_t(nClust);\r
412 }\r
dac20c1b 413\r
414 Float_t clustPerFindClust = 0.;\r
d1c6e451 415 Int_t nFindableClust = esdTrack->GetTPCNclsF();\r
dac20c1b 416 if(nFindableClust>0.) clustPerFindClust = Float_t(nClust)/nFindableClust;\r
417\r
418 Float_t b[2], bCov[3];\r
419 esdTrack->GetImpactParameters(b,bCov);\r
420\r
10360c58 421 // \r
422 Float_t nCrossedRowsTPC = esdTrack->GetTPCClusterInfo(2,1);\r
423\r
424 Float_t ratioCrossedRowsOverFindableClustersTPC = 1.0;\r
425 if (esdTrack->GetTPCNclsF()>0) {\r
426 ratioCrossedRowsOverFindableClustersTPC = esdTrack->GetTPCClusterInfo(2,1)/esdTrack->GetTPCNclsF();\r
427 }\r
428\r
429 //\r
430 Int_t nClustersTPCShared = esdTrack->GetTPCnclsS();\r
431 Float_t fracClustersTPCShared = -1.;\r
432 fracClustersTPCShared = Float_t(nClustersTPCShared)/Float_t(nClust);\r
433\r
434 // \r
dac20c1b 435 // kink idx\r
436 Int_t kinkIdx = 0;\r
437 //if(esdTrack->GetKinkIndex(0) > 0.) isKink = kTRUE;\r
438 if(esdTrack->GetKinkIndex(0) > 0) kinkIdx = 1; // kink daughter\r
439 else if(esdTrack->GetKinkIndex(0) < 0) kinkIdx = -1; // kink mother\r
440 else kinkIdx = 0; // not kink\r
441\r
442 //printf("esdTrack->GetKinkIndex(0) %d \n", esdTrack->GetKinkIndex(0));\r
443 //printf("esdTrack->GetKinkIndex(1) %d \n", esdTrack->GetKinkIndex(1));\r
444 //printf("esdTrack->GetKinkIndex(2) %d \n", esdTrack->GetKinkIndex(2));\r
445 //printf("kinkIdx %d \n", kinkIdx);\r
446\r
447 //\r
448 // Fill rec vs MC information\r
449 //\r
450 Bool_t isPrim = kTRUE;\r
dffe2c06 451 Bool_t hasStrangeMother = kFALSE;\r
6a01c05f 452 Bool_t isFromMaterial = kFALSE;\r
dac20c1b 453\r
454 if(IsUseMCInfo()) {\r
455 if(!stack) return;\r
456 Int_t label = TMath::Abs(esdTrack->GetLabel()); \r
457 TParticle* particle = stack->Particle(label);\r
458 if(!particle) return;\r
459 if(particle->GetPDG() && particle->GetPDG()->Charge()==0.) return;\r
460 isPrim = stack->IsPhysicalPrimary(label);\r
461\r
e84b59cb 462 // check whether has stange mother\r
463 //\r
464 Int_t motherPdg = -1; \r
465 TParticle* mother = 0; \r
466 \r
467 Int_t motherLabel = particle->GetMother(0); \r
468 if(motherLabel>0) mother = stack->Particle(motherLabel); \r
469 if(mother) motherPdg = TMath::Abs(mother->GetPdgCode()); // take abs for visualisation only \r
470 Int_t mech = particle->GetUniqueID(); // production mechanism \r
dffe2c06 471\r
472 if( (motherPdg == 3122) || (motherPdg == -3122) || (motherPdg == 310)) // lambda, antilambda, k0s\r
473 {\r
474 if( (mech == 4) || (mech == 5) ) hasStrangeMother = kTRUE;\r
6a01c05f 475 } \r
476 else {\r
477 //if(isPrim==0 && mech == 13) \r
478 //printf("mech %d \n", mech);\r
479 if(!isPrim) isFromMaterial = kTRUE; \r
dffe2c06 480 }\r
dac20c1b 481 }\r
482 \r
483 // fill histo\r
6a01c05f 484 Int_t charge = esdTrack->Charge();\r
10360c58 485 //Double_t vRecMCTrackHist[13] = { nClust,chi2PerCluster,clustPerFindClust,b[0],b[1],eta,phi,pt,hasStrangeMother,isFromMaterial,isPrim,charge, centralityF }; \r
486 Double_t vRecMCTrackHist[14] = { nCrossedRowsTPC, chi2PerCluster, ratioCrossedRowsOverFindableClustersTPC, fracClustersTPCShared , b[0], b[1], eta, phi, pt, hasStrangeMother, isFromMaterial, isPrim, charge, centralityF }; \r
dac20c1b 487 fRecMCTrackHist->Fill(vRecMCTrackHist);\r
488}\r
489\r
490//_____________________________________________________________________________\r
491Long64_t AlidNdPtCutAnalysisPbPb::Merge(TCollection* const list) \r
492{\r
493 // Merge list of objects (needed by PROOF)\r
494\r
495 if (!list)\r
496 return 0;\r
497\r
498 if (list->IsEmpty())\r
499 return 1;\r
500\r
501 TIterator* iter = list->MakeIterator();\r
502 TObject* obj = 0;\r
503\r
504 //TList *collPhysSelection = new TList;\r
505\r
506 // collection of generated histograms\r
507 Int_t count=0;\r
508 while((obj = iter->Next()) != 0) {\r
509 AlidNdPtCutAnalysisPbPb* entry = dynamic_cast<AlidNdPtCutAnalysisPbPb*>(obj);\r
510 if (entry == 0) continue; \r
511 \r
512 // event histo\r
513 fEventCount->Add(entry->fEventCount);\r
514 fRecEventHist->Add(entry->fRecEventHist);\r
515 fRecMCEventHist->Add(entry->fRecMCEventHist);\r
516 fMCEventHist->Add(entry->fMCEventHist);\r
517\r
518 // track histo\r
519 fRecMCTrackHist->Add(entry->fRecMCTrackHist);\r
520\r
521 // physics selection\r
522 //collPhysSelection->Add(entry->GetPhysicsTriggerSelection());\r
523 \r
524 count++;\r
525 }\r
526\r
527 //AliPhysicsSelection *trigSelection = GetPhysicsTriggerSelection();\r
528 //trigSelection->Merge(collPhysSelection);\r
529\r
530 //if(collPhysSelection) delete collPhysSelection;\r
531\r
532return count;\r
533}\r
534\r
535//_____________________________________________________________________________\r
536void AlidNdPtCutAnalysisPbPb::Analyse() \r
537{\r
538 //\r
539 // Analyse histograms\r
540 //\r
541 TH1::AddDirectory(kFALSE);\r
542 TObjArray *aFolderObj = new TObjArray;\r
543 TH1D *h1D = 0; \r
544 TH2D *h2D = 0; \r
545\r
546\r
547 //\r
548 // get cuts\r
549 //\r
550 AlidNdPtEventCuts *evtCuts = GetEventCuts(); \r
551 AlidNdPtAcceptanceCuts *accCuts = GetAcceptanceCuts(); \r
552 AliESDtrackCuts *esdTrackCuts = GetTrackCuts(); \r
553\r
554 if(!evtCuts || !accCuts || !esdTrackCuts) {\r
555 Error("AlidNdPtCutAnalysisPbPb::Analyse()", "cuts not available");\r
556 return;\r
557 }\r
558\r
559 //\r
560 // set min and max values\r
561 //\r
562 Double_t minPt = accCuts->GetMinPt();\r
563 Double_t maxPt = accCuts->GetMaxPt();\r
564 Double_t minEta = accCuts->GetMinEta();\r
565 Double_t maxEta = accCuts->GetMaxEta()-0.00001;\r
566\r
567 Double_t maxDCAr = accCuts->GetMaxDCAr();\r
568\r
569 //\r
570 // Event counters\r
571 //\r
572 h2D = (TH2D*)fEventCount->Projection(0,1);\r
573 h2D->SetName("trig_vs_trigANDvertex");\r
574 aFolderObj->Add(h2D);\r
575\r
576 fEventCount->GetAxis(0)->SetRange(2,2); // triggered\r
577 h1D = (TH1D*)fEventCount->Projection(1);\r
578 h1D->SetTitle("rec. vertex for triggered events");\r
579 h1D->SetName("trigANDvertex");\r
580 aFolderObj->Add(h1D);\r
581\r
582 //\r
583 // Create rec. event histograms\r
584 //\r
585 h1D = (TH1D *)fRecEventHist->Projection(0);\r
586 h1D->SetName("rec_xv");\r
587 aFolderObj->Add(h1D);\r
588\r
589 h1D = (TH1D *)fRecEventHist->Projection(1);\r
590 h1D->SetName("rec_yv");\r
591 aFolderObj->Add(h1D);\r
592\r
593 h1D = (TH1D *)fRecEventHist->Projection(2);\r
594 h1D->SetName("rec_zv");\r
595 aFolderObj->Add(h1D);\r
596\r
597 h2D = (TH2D *)fRecEventHist->Projection(3,4);\r
598 h2D->SetName("rec_resZv_vs_Mult");\r
599 aFolderObj->Add(h2D);\r
600\r
601 h2D = (TH2D *)fRecEventHist->Projection(0,1);\r
602 h2D->SetName("rec_xv_vs_yv");\r
603 aFolderObj->Add(h2D);\r
604\r
605 h2D = (TH2D *)fRecEventHist->Projection(0,2);\r
606 h2D->SetName("rec_xv_vs_zv");\r
607 aFolderObj->Add(h2D);\r
608\r
609 h2D = (TH2D *)fRecEventHist->Projection(3,4);\r
610 h2D->SetName("rec_resZv_vs_Mult");\r
611 aFolderObj->Add(h2D);\r
612\r
613 //\r
614 // MC available\r
615 //\r
616 if(IsUseMCInfo()) {\r
617\r
618 //\r
619 // Create mc event histograms\r
620 //\r
621 h2D = (TH2D *)fMCEventHist->Projection(0,1);\r
622 h2D->SetName("mc_xv_vs_yv");\r
623 aFolderObj->Add(h2D);\r
624\r
625 h2D = (TH2D *)fMCEventHist->Projection(0,2);\r
626 h2D->SetName("mc_xv_vs_zv");\r
627 aFolderObj->Add(h2D);\r
628\r
629 //\r
630 // Create rec-mc event histograms\r
631 //\r
632 h2D = (TH2D *)fRecMCEventHist->Projection(0,3);\r
633 h2D->SetName("rec_mc_deltaXv_vs_mult");\r
634 aFolderObj->Add(h2D);\r
635\r
636 h2D = (TH2D *)fRecMCEventHist->Projection(1,3);\r
637 h2D->SetName("rec_mc_deltaYv_vs_mult");\r
638 aFolderObj->Add(h2D);\r
639\r
640 h2D = (TH2D *)fRecMCEventHist->Projection(2,3);\r
641 h2D->SetName("rec_mc_deltaZv_vs_mult");\r
642 aFolderObj->Add(h2D);\r
643\r
644 } // end use MC info \r
645\r
646\r
647\r
648 //\r
649 // Create rec-mc track track histograms \r
650 //\r
651\r
652 // DCA cuts\r
653 fRecMCTrackHist->GetAxis(3)->SetRangeUser(-maxDCAr,maxDCAr);\r
654 fRecMCTrackHist->GetAxis(4)->SetRangeUser(-maxDCAr,maxDCAr);\r
655\r
656 h2D = (TH2D *)fRecMCTrackHist->Projection(7,5);\r
657 h2D->SetName("pt_vs_eta");\r
658 aFolderObj->Add(h2D);\r
659\r
660 fRecMCTrackHist->GetAxis(7)->SetRangeUser(minPt,maxPt); \r
661\r
662 h2D = (TH2D *)fRecMCTrackHist->Projection(0,5);\r
663 h2D->SetName("nClust_vs_eta");\r
664 aFolderObj->Add(h2D);\r
665\r
666 h2D = (TH2D *)fRecMCTrackHist->Projection(1,5);\r
667 h2D->SetName("chi2PerClust_vs_eta");\r
668 aFolderObj->Add(h2D);\r
669\r
670 h2D = (TH2D *)fRecMCTrackHist->Projection(2,5);\r
671 h2D->SetName("ratio_nClust_nFindableClust_vs_eta");\r
672 aFolderObj->Add(h2D);\r
673\r
674 h2D = (TH2D *)fRecMCTrackHist->Projection(5,6);\r
675 h2D->SetName("eta_vs_phi");\r
676 aFolderObj->Add(h2D);\r
677\r
678 //\r
679 fRecMCTrackHist->GetAxis(5)->SetRangeUser(minEta,maxEta); \r
680\r
681 h2D = (TH2D *)fRecMCTrackHist->Projection(0,6);\r
682 h2D->SetName("nClust_vs_phi");\r
683 aFolderObj->Add(h2D);\r
684\r
685 h2D = (TH2D *)fRecMCTrackHist->Projection(1,6);\r
686 h2D->SetName("chi2PerClust_vs_phi");\r
687 aFolderObj->Add(h2D);\r
688\r
689 h2D = (TH2D *)fRecMCTrackHist->Projection(2,6);\r
690 h2D->SetName("ratio_nClust_nFindableClust_vs_phi");\r
691 aFolderObj->Add(h2D);\r
692\r
693 //\r
694 fRecMCTrackHist->GetAxis(7)->SetRangeUser(0.0,maxPt); \r
695\r
696 h2D = (TH2D *)fRecMCTrackHist->Projection(0,7);\r
697 h2D->SetName("nClust_vs_pt");\r
698 aFolderObj->Add(h2D);\r
699\r
700 h2D = (TH2D *)fRecMCTrackHist->Projection(1,7);\r
701 h2D->SetName("chi2PerClust_vs_pt");\r
702 aFolderObj->Add(h2D);\r
703\r
704 h2D = (TH2D *)fRecMCTrackHist->Projection(2,7);\r
705 h2D->SetName("ratio_nClust_nFindableClust_vs_pt");\r
706 aFolderObj->Add(h2D);\r
707\r
708 h2D = (TH2D *)fRecMCTrackHist->Projection(6,7);\r
709 h2D->SetName("phi_vs_pt");\r
710 aFolderObj->Add(h2D);\r
711\r
712\r
713 // fiducial volume\r
714 fRecMCTrackHist->GetAxis(5)->SetRangeUser(minEta,maxEta); \r
715 fRecMCTrackHist->GetAxis(7)->SetRangeUser(minPt,maxPt); \r
716\r
717 // DCA cuts\r
718 fRecMCTrackHist->GetAxis(3)->SetRangeUser(-maxDCAr,maxDCAr);\r
719 fRecMCTrackHist->GetAxis(4)->SetRangeUser(-maxDCAr,maxDCAr);\r
720\r
721 h2D = (TH2D *)fRecMCTrackHist->Projection(0,1);\r
722 h2D->SetName("nClust_vs_chi2PerClust");\r
723 aFolderObj->Add(h2D);\r
724\r
725 h2D = (TH2D *)fRecMCTrackHist->Projection(0,2);\r
726 h2D->SetName("nClust_vs_ratio_nClust_nFindableClust");\r
727 aFolderObj->Add(h2D);\r
728\r
729 //\r
730 // DCAy cuts\r
731 //\r
732 fRecMCTrackHist->GetAxis(0)->SetRange(50,160); // nClust/track > 50\r
733 fRecMCTrackHist->GetAxis(1)->SetRangeUser(0.,3.9999); // chi2/cluster < 4.0\r
734 fRecMCTrackHist->GetAxis(3)->SetRange(1,fRecMCTrackHist->GetAxis(3)->GetNbins());\r
735 //fRecMCTrackHist->GetAxis(4)->SetRangeUser(-1.0,1.0);\r
736 fRecMCTrackHist->GetAxis(4)->SetRange(1,fRecMCTrackHist->GetAxis(4)->GetNbins());\r
737\r
738 // sec\r
739 fRecMCTrackHist->GetAxis(9)->SetRange(1,1);\r
740 h1D = (TH1D *)fRecMCTrackHist->Projection(3);\r
741 h1D->SetName("dcay_sec");\r
742 aFolderObj->Add(h1D);\r
743\r
744 // prim\r
745 fRecMCTrackHist->GetAxis(9)->SetRange(2,2);\r
746 h1D = (TH1D *)fRecMCTrackHist->Projection(3);\r
747 h1D->SetName("dcay_prim");\r
748 aFolderObj->Add(h1D);\r
749\r
750 // DCAz cuts\r
751 //fRecMCTrackHist->GetAxis(3)->SetRangeUser(-1.0,1.0);\r
752 fRecMCTrackHist->GetAxis(4)->SetRange(1,fRecMCTrackHist->GetAxis(4)->GetNbins());\r
753\r
754 // sec\r
755 fRecMCTrackHist->GetAxis(9)->SetRange(1,1);\r
756 h1D = (TH1D *)fRecMCTrackHist->Projection(4);\r
757 h1D->SetName("dcaz_sec");\r
758 aFolderObj->Add(h1D);\r
759\r
760 // prim\r
761 fRecMCTrackHist->GetAxis(9)->SetRange(2,2);\r
762 h1D = (TH1D *)fRecMCTrackHist->Projection(4);\r
763 h1D->SetName("dcaz_prim");\r
764 aFolderObj->Add(h1D);\r
765\r
766\r
767 // export objects to analysis folder\r
768 fAnalysisFolder = ExportToFolder(aFolderObj);\r
769\r
770 // delete only TObjArray\r
771 if(aFolderObj) delete aFolderObj;\r
772}\r
773\r
774//_____________________________________________________________________________\r
775TFolder* AlidNdPtCutAnalysisPbPb::ExportToFolder(TObjArray * const array) \r
776{\r
777 // recreate folder avery time and export objects to new one\r
778 //\r
779 AlidNdPtCutAnalysisPbPb * comp=this;\r
780 TFolder *folder = comp->GetAnalysisFolder();\r
781\r
782 TString name, title;\r
783 TFolder *newFolder = 0;\r
784 Int_t i = 0;\r
785 Int_t size = array->GetSize();\r
786\r
787 if(folder) { \r
788 // get name and title from old folder\r
789 name = folder->GetName(); \r
790 title = folder->GetTitle(); \r
791\r
792 // delete old one\r
793 delete folder;\r
794\r
795 // create new one\r
796 newFolder = CreateFolder(name.Data(),title.Data());\r
797 newFolder->SetOwner();\r
798\r
799 // add objects to folder\r
800 while(i < size) {\r
801 newFolder->Add(array->At(i));\r
802 i++;\r
803 }\r
804 }\r
805\r
806return newFolder;\r
807}\r
808\r
809//_____________________________________________________________________________\r
810TFolder* AlidNdPtCutAnalysisPbPb::CreateFolder(TString name,TString title) { \r
811// create folder for analysed histograms\r
812//\r
813TFolder *folder = 0;\r
814 folder = new TFolder(name.Data(),title.Data());\r
815\r
816 return folder;\r
817}\r