Cleaning up files, added documentation. Added ntuple for v0s.
[u/mrichter/AliRoot.git] / PWG4 / GammaConv / AliV0Reader.cxx
CommitLineData
3c538586 1/**************************************************************************\r
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *\r
3 * *\r
4 * Author: Ana Marin, Kathrin Koch, Kenneth Aamodt *\r
5 * Version 1.0 *\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////////////////////////////////////////////////\r
17//--------------------------------------------- \r
18// Class used to do analysis on conversion pairs\r
19//---------------------------------------------\r
20////////////////////////////////////////////////\r
21\r
22// --- ROOT system ---\r
23#include <TMath.h>\r
24\r
25//---- ANALYSIS system ----\r
26#include "AliV0Reader.h"\r
27#include "AliAnalysisManager.h"\r
28#include "AliESDInputHandler.h"\r
29#include "AliMCEvent.h"\r
30#include "AliKFVertex.h"\r
31\r
32#include "AliStack.h"\r
33#include "AliMCEventHandler.h"\r
34\r
35\r
36class iostream;\r
37class AliESDv0;\r
38class TFormula;\r
39\r
40using namespace std;\r
41\r
42ClassImp(AliV0Reader)\r
43\r
44\r
45\r
46 AliV0Reader::AliV0Reader() :\r
47 TObject(),\r
48 fMCStack(NULL),\r
49 fMCTruth(NULL),\r
50 fChain(NULL),\r
51 fESDHandler(NULL),\r
52 fESDEvent(NULL),\r
53 fHistograms(NULL),\r
54 fCurrentV0IndexNumber(0),\r
55 fCurrentV0(NULL),\r
56 fCurrentNegativeKFParticle(NULL),\r
57 fCurrentPositiveKFParticle(NULL),\r
58 fCurrentMotherKFCandidate(NULL),\r
59 fCurrentNegativeESDTrack(NULL),\r
60 fCurrentPositiveESDTrack(NULL),\r
61 fNegativeTrackLorentzVector(NULL),\r
62 fPositiveTrackLorentzVector(NULL),\r
63 fMotherCandidateLorentzVector(NULL),\r
64 fCurrentXValue(0),\r
65 fCurrentYValue(0),\r
66 fCurrentZValue(0),\r
67 fPositiveTrackPID(0),\r
68 fNegativeTrackPID(0),\r
69 fNegativeMCParticle(NULL),\r
70 fPositiveMCParticle(NULL),\r
71 fMotherMCParticle(NULL),\r
72 fMotherCandidateKFMass(0),\r
73 fMotherCandidateKFWidth(0),\r
74 fUseKFParticle(kTRUE),\r
75 fUseESDTrack(kFALSE),\r
76 fDoMC(kFALSE),\r
77 fMaxR(10000),// 100 meter(outside of ALICE)\r
78 fEtaCut(0.),\r
79 fPtCut(0.),\r
80 fChi2CutConversion(0.),\r
81 fChi2CutMeson(0.),\r
82 fPIDProbabilityCutNegativeParticle(0),\r
83 fPIDProbabilityCutPositiveParticle(0),\r
84 fXVertexCut(0.),\r
85 fYVertexCut(0.),\r
86 fZVertexCut(0.),\r
87 fNSigmaMass(0.),\r
88 fUseImprovedVertex(kFALSE),\r
89 fCurrentEventGoodV0s(),\r
90 fPreviousEventGoodV0s()\r
91{\r
92\r
93}\r
94\r
95\r
96AliV0Reader::AliV0Reader(const AliV0Reader & original) :\r
97 TObject(original),\r
98 fMCStack(original.fMCStack),\r
99 fMCTruth(original.fMCTruth),\r
100 fChain(original.fChain),\r
101 fESDHandler(original.fESDHandler),\r
102 fESDEvent(original.fESDEvent),\r
103 fHistograms(original.fHistograms),\r
104 fCurrentV0IndexNumber(original.fCurrentV0IndexNumber),\r
105 fCurrentV0(original.fCurrentV0),\r
106 fCurrentNegativeKFParticle(original.fCurrentNegativeKFParticle),\r
107 fCurrentPositiveKFParticle(original.fCurrentPositiveKFParticle),\r
108 fCurrentMotherKFCandidate(original.fCurrentMotherKFCandidate),\r
109 fCurrentNegativeESDTrack(original.fCurrentNegativeESDTrack),\r
110 fCurrentPositiveESDTrack(original.fCurrentPositiveESDTrack),\r
111 fNegativeTrackLorentzVector(original.fNegativeTrackLorentzVector),\r
112 fPositiveTrackLorentzVector(original.fPositiveTrackLorentzVector),\r
113 fMotherCandidateLorentzVector(original.fMotherCandidateLorentzVector),\r
114 fCurrentXValue(original.fCurrentXValue),\r
115 fCurrentYValue(original.fCurrentYValue),\r
116 fCurrentZValue(original.fCurrentZValue),\r
117 fPositiveTrackPID(original.fPositiveTrackPID),\r
118 fNegativeTrackPID(original.fNegativeTrackPID),\r
119 fNegativeMCParticle(original.fNegativeMCParticle),\r
120 fPositiveMCParticle(original.fPositiveMCParticle),\r
121 fMotherMCParticle(original.fMotherMCParticle),\r
122 fMotherCandidateKFMass(original.fMotherCandidateKFMass),\r
123 fMotherCandidateKFWidth(original.fMotherCandidateKFWidth),\r
124 fUseKFParticle(kTRUE),\r
125 fUseESDTrack(kFALSE),\r
126 fDoMC(kFALSE),\r
127 fMaxR(original.fMaxR),\r
128 fEtaCut(original.fEtaCut),\r
129 fPtCut(original.fPtCut),\r
130 fChi2CutConversion(original.fChi2CutConversion),\r
131 fChi2CutMeson(original.fChi2CutMeson),\r
132 fPIDProbabilityCutNegativeParticle(original.fPIDProbabilityCutNegativeParticle),\r
133 fPIDProbabilityCutPositiveParticle(original.fPIDProbabilityCutPositiveParticle),\r
134 fXVertexCut(original.fXVertexCut),\r
135 fYVertexCut(original.fYVertexCut),\r
136 fZVertexCut(original.fZVertexCut),\r
137 fNSigmaMass(original.fNSigmaMass),\r
138 fUseImprovedVertex(original.fUseImprovedVertex),\r
139 fCurrentEventGoodV0s(original.fCurrentEventGoodV0s),\r
140 fPreviousEventGoodV0s(original.fPreviousEventGoodV0s)\r
141{\r
142\r
143}\r
144\r
145\r
146AliV0Reader & AliV0Reader::operator = (const AliV0Reader & /*source*/)\r
147{\r
148 // assignment operator\r
149 return *this;\r
150}\r
151\r
152void AliV0Reader::Initialize(){\r
153 //see header file for documentation\r
154\r
155 // Get the input handler from the manager\r
156 fESDHandler = (AliESDInputHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());\r
157 if(fESDHandler == NULL){\r
158 //print warning here\r
159 }\r
160 \r
161 // Get pointer to esd event from input handler\r
162 fESDEvent = fESDHandler->GetEvent();\r
163 if(fESDEvent == NULL){\r
164 //print warning here\r
165 }\r
166\r
167 //Get pointer to MCTruth\r
168 fMCTruth = (AliMCEventHandler*)((AliAnalysisManager::GetAnalysisManager())->GetMCtruthEventHandler());\r
169 if(fMCTruth == NULL){\r
170 //print warning here\r
171 }\r
172\r
173 //Get pointer to the mc stack\r
174 fMCStack = fMCTruth->MCEvent()->Stack();\r
175 if(fMCStack == NULL){\r
176 //print warning here\r
177 }\r
178\r
179 AliKFParticle::SetField(fESDEvent->GetMagneticField());\r
180\r
181}\r
182\r
183AliESDv0* AliV0Reader::GetV0(Int_t index){\r
184 //see header file for documentation\r
185\r
186 fCurrentV0 = fESDEvent->GetV0(index);\r
187 UpdateV0Information();\r
188 return fCurrentV0;\r
189}\r
a19c3402 190Bool_t AliV0Reader::CheckForPrimaryVertex(){\r
191 return fESDEvent->GetPrimaryVertex()->GetNContributors()>0;\r
192}\r
3c538586 193\r
194Bool_t AliV0Reader::NextV0(){\r
195 //see header file for documentation\r
196\r
197 Bool_t iResult=kFALSE;\r
198 while(fCurrentV0IndexNumber<fESDEvent->GetNumberOfV0s()){\r
199 fCurrentV0 = fESDEvent->GetV0(fCurrentV0IndexNumber);\r
200 \r
201 //checks if on the fly mode is set\r
202 if ( !fCurrentV0->GetOnFlyStatus() ){\r
203 fCurrentV0IndexNumber++;\r
a19c3402 204 if(fHistograms != NULL){\r
205 fHistograms->FillHistogram("V0MassDebugCut1",GetMotherCandidateMass());\r
206 }\r
3c538586 207 continue;\r
208 }\r
209\r
210 if(fESDEvent->GetPrimaryVertex()->GetNContributors()<=0) {//checks if we have a vertex\r
211 fCurrentV0IndexNumber++;\r
a19c3402 212 if(fHistograms != NULL){\r
213 fHistograms->FillHistogram("V0MassDebugCut2",GetMotherCandidateMass());\r
214 }\r
3c538586 215 continue;\r
216 }\r
217\r
218 if(CheckPIDProbability(fPIDProbabilityCutNegativeParticle,fPIDProbabilityCutPositiveParticle)==kFALSE){\r
219 fCurrentV0IndexNumber++;\r
a19c3402 220 if(fHistograms != NULL){\r
221 fHistograms->FillHistogram("V0MassDebugCut3",GetMotherCandidateMass());\r
222 }\r
3c538586 223 continue;\r
224 }\r
225\r
3c538586 226 fCurrentV0->GetXYZ(fCurrentXValue,fCurrentYValue,fCurrentZValue);\r
227 \r
228 if(GetXYRadius()>fMaxR){ // cuts on distance from collision point\r
229 fCurrentV0IndexNumber++;\r
a19c3402 230 if(fHistograms != NULL){\r
231 fHistograms->FillHistogram("V0MassDebugCut4",GetMotherCandidateMass());\r
232 }\r
3c538586 233 continue;\r
234 }\r
235\r
236 UpdateV0Information();\r
237 \r
238 if(fUseKFParticle){\r
239 if(fCurrentMotherKFCandidate->GetNDF()<=0){\r
240 fCurrentV0IndexNumber++;\r
a19c3402 241 if(fHistograms != NULL){\r
242 fHistograms->FillHistogram("V0MassDebugCut5",GetMotherCandidateMass());\r
243 }\r
3c538586 244 continue;\r
245 }\r
246 Double_t chi2V0 = fCurrentMotherKFCandidate->GetChi2()/fCurrentMotherKFCandidate->GetNDF();\r
247 if(chi2V0 > fChi2CutConversion || chi2V0 <=0){\r
248 fCurrentV0IndexNumber++;\r
a19c3402 249 if(fHistograms != NULL){\r
250 fHistograms->FillHistogram("V0MassDebugCut6",GetMotherCandidateMass());\r
251 }\r
3c538586 252 continue;\r
253 }\r
254 \r
255 if(TMath::Abs(fMotherCandidateLorentzVector->Eta())> fEtaCut){\r
256 fCurrentV0IndexNumber++;\r
a19c3402 257 if(fHistograms != NULL){\r
258 fHistograms->FillHistogram("V0MassDebugCut7",GetMotherCandidateMass());\r
259 }\r
3c538586 260 continue;\r
261 }\r
262 \r
263 if(fMotherCandidateLorentzVector->Pt()<fPtCut){\r
264 fCurrentV0IndexNumber++;\r
a19c3402 265 if(fHistograms != NULL){\r
266 fHistograms->FillHistogram("V0MassDebugCut8",GetMotherCandidateMass());\r
267 }\r
3c538586 268 continue;\r
269 }\r
3c538586 270 }\r
271 else if(fUseESDTrack){\r
272 //TODO\r
273 }\r
274\r
275 iResult=kTRUE;//means we have a v0 who survived all the cuts applied\r
276\r
277 fCurrentV0IndexNumber++;\r
278 \r
279 break;\r
280 }\r
281 return iResult; \r
282}\r
283\r
284void AliV0Reader::UpdateV0Information(){\r
285 //see header file for documentation\r
286 \r
287 if(fCurrentNegativeKFParticle != NULL){\r
288 delete fCurrentNegativeKFParticle;\r
289 }\r
290 fCurrentNegativeKFParticle = new AliKFParticle(*(fCurrentV0->GetParamN()),fNegativeTrackPID);\r
291 \r
292 if(fCurrentPositiveKFParticle != NULL){\r
293 delete fCurrentPositiveKFParticle;\r
294 }\r
295 fCurrentPositiveKFParticle = new AliKFParticle(*(fCurrentV0->GetParamP()),fPositiveTrackPID);\r
296 \r
297 if(fCurrentMotherKFCandidate != NULL){\r
298 delete fCurrentMotherKFCandidate;\r
299 }\r
300 fCurrentMotherKFCandidate = new AliKFParticle(*fCurrentNegativeKFParticle,*fCurrentPositiveKFParticle);\r
301\r
302 fCurrentNegativeESDTrack = fESDEvent->GetTrack(fCurrentV0->GetNindex());\r
303\r
304 fCurrentPositiveESDTrack = fESDEvent->GetTrack(fCurrentV0->GetPindex());\r
305\r
306 if(fPositiveTrackPID==-11 && fNegativeTrackPID==11){\r
307 fCurrentMotherKFCandidate->SetMassConstraint(0,fNSigmaMass);\r
308 }\r
309\r
310 if(fUseImprovedVertex == kTRUE){\r
311 AliKFVertex primaryVertexImproved(*GetPrimaryVertex());\r
312 primaryVertexImproved+=*fCurrentMotherKFCandidate;\r
313 fCurrentMotherKFCandidate->SetProductionVertex(primaryVertexImproved);\r
314 }\r
315\r
316 fCurrentMotherKFCandidate->GetMass(fMotherCandidateKFMass,fMotherCandidateKFWidth);\r
317\r
318\r
319 if(fNegativeTrackLorentzVector != NULL){\r
320 delete fNegativeTrackLorentzVector;\r
321 }\r
322 if(fUseKFParticle){\r
323 fNegativeTrackLorentzVector = new TLorentzVector(fCurrentNegativeKFParticle->Px(),fCurrentNegativeKFParticle->Py(),fCurrentNegativeKFParticle->Pz());\r
324 }\r
325 else if(fUseESDTrack){\r
326 fNegativeTrackLorentzVector = new TLorentzVector(fCurrentNegativeESDTrack->Px(),fCurrentNegativeESDTrack->Py(),fCurrentNegativeESDTrack->Pz());\r
327 }\r
328\r
329 if(fPositiveTrackLorentzVector != NULL){\r
330 delete fPositiveTrackLorentzVector;\r
331 }\r
332 if(fUseKFParticle){\r
333 fPositiveTrackLorentzVector = new TLorentzVector(fCurrentPositiveKFParticle->Px(),fCurrentPositiveKFParticle->Py(),fCurrentPositiveKFParticle->Pz());\r
334 }\r
335 else if(fUseESDTrack){\r
336 fPositiveTrackLorentzVector = new TLorentzVector(fCurrentPositiveESDTrack->Px(),fCurrentPositiveESDTrack->Py(),fCurrentPositiveESDTrack->Pz());\r
337 }\r
338\r
339 if(fMotherCandidateLorentzVector != NULL){\r
340 delete fMotherCandidateLorentzVector;\r
341 }\r
342 if(fUseKFParticle){\r
343 fMotherCandidateLorentzVector = new TLorentzVector(*fNegativeTrackLorentzVector + *fPositiveTrackLorentzVector);\r
344 }\r
345 else if(fUseESDTrack){\r
346 fMotherCandidateLorentzVector = new TLorentzVector(*fNegativeTrackLorentzVector + *fPositiveTrackLorentzVector);\r
347 }\r
348\r
349 if(fPositiveTrackPID==-11 && fNegativeTrackPID==11){\r
350 fMotherCandidateLorentzVector->SetXYZM(fMotherCandidateLorentzVector->Px() ,fMotherCandidateLorentzVector->Py(),fMotherCandidateLorentzVector->Pz(),0.); \r
351 }\r
352 \r
353 if(fDoMC == kTRUE){\r
a19c3402 354 fMotherMCParticle= NULL;\r
3c538586 355 fNegativeMCParticle = fMCStack->Particle(TMath::Abs(fESDEvent->GetTrack(fCurrentV0->GetNindex())->GetLabel()));\r
356 fPositiveMCParticle = fMCStack->Particle(TMath::Abs(fESDEvent->GetTrack(fCurrentV0->GetPindex())->GetLabel()));\r
a19c3402 357 if(fPositiveMCParticle->GetMother(0)>-1){\r
358 fMotherMCParticle = fMCStack->Particle(fPositiveMCParticle->GetMother(0));\r
359 }\r
3c538586 360 }\r
361 fCurrentEventGoodV0s.push_back(*fCurrentMotherKFCandidate);\r
362}\r
363\r
364Bool_t AliV0Reader::HasSameMCMother(){\r
365 //see header file for documentation\r
366\r
367 Bool_t iResult = kFALSE;\r
368 if(fDoMC == kTRUE){\r
369 if(fNegativeMCParticle != NULL && fPositiveMCParticle != NULL){\r
370 if(fNegativeMCParticle->GetMother(0) == fPositiveMCParticle->GetMother(0))\r
a19c3402 371 if(fMotherMCParticle){\r
372 iResult = kTRUE;\r
373 }\r
3c538586 374 }\r
375 }\r
376 return iResult;\r
377}\r
378\r
379Bool_t AliV0Reader::CheckPIDProbability(Double_t negProbCut, Double_t posProbCut){\r
380 //see header file for documentation\r
381\r
382 Bool_t iResult=kFALSE;\r
383\r
384 Double_t *posProbArray = new Double_t[10];\r
385 Double_t *negProbArray = new Double_t[10];\r
386 AliESDtrack* negTrack = fESDEvent->GetTrack(fCurrentV0->GetNindex());\r
387 AliESDtrack* posTrack = fESDEvent->GetTrack(fCurrentV0->GetPindex());\r
388 \r
389 negTrack->GetTPCpid(negProbArray);\r
390 posTrack->GetTPCpid(posProbArray);\r
391\r
392 if(negProbArray!=NULL && posProbArray!=NULL){\r
393 if(negProbArray[GetSpeciesIndex(-1)]>=negProbCut && posProbArray[GetSpeciesIndex(1)]>=posProbCut){\r
394 iResult=kTRUE;\r
395 }\r
396 }\r
397 delete [] posProbArray;\r
398 delete [] negProbArray;\r
399 return iResult;\r
400}\r
401\r
a19c3402 402void AliV0Reader::GetPIDProbability(Double_t &negPIDProb,Double_t & posPIDProb){\r
403\r
404 Double_t *posProbArray = new Double_t[10];\r
405 Double_t *negProbArray = new Double_t[10];\r
406 AliESDtrack* negTrack = fESDEvent->GetTrack(fCurrentV0->GetNindex());\r
407 AliESDtrack* posTrack = fESDEvent->GetTrack(fCurrentV0->GetPindex());\r
408 \r
409 negTrack->GetTPCpid(negProbArray);\r
410 posTrack->GetTPCpid(posProbArray);\r
411\r
412 if(negProbArray!=NULL && posProbArray!=NULL){\r
413 negPIDProb = negProbArray[GetSpeciesIndex(-1)];\r
414 posPIDProb = posProbArray[GetSpeciesIndex(1)];\r
415 }\r
416 delete [] posProbArray;\r
417 delete [] negProbArray;\r
418}\r
419\r
3c538586 420void AliV0Reader::UpdateEventByEventData(){\r
421 //see header file for documentation\r
422\r
423 if(fCurrentEventGoodV0s.size() >0 ){\r
424 fPreviousEventGoodV0s.clear();\r
425 fPreviousEventGoodV0s = fCurrentEventGoodV0s;\r
426 }\r
427 fCurrentEventGoodV0s.clear();\r
428 \r
429 fCurrentV0IndexNumber=0;\r
430}\r
431\r
432Double_t AliV0Reader::GetNegativeTrackPhi() const{\r
433 //see header file for documentation\r
434\r
435 Double_t offset=0;\r
436 if(fNegativeTrackLorentzVector->Phi()> TMath::Pi()){\r
437 offset = -2*TMath::Pi();\r
438 }\r
439 return fNegativeTrackLorentzVector->Phi()+offset;\r
440}\r
441\r
442Double_t AliV0Reader::GetPositiveTrackPhi() const{\r
443 //see header file for documentation\r
444\r
445 Double_t offset=0;\r
446 if(fPositiveTrackLorentzVector->Phi()> TMath::Pi()){\r
447 offset = -2*TMath::Pi();\r
448 }\r
449 return fPositiveTrackLorentzVector->Phi()+offset;\r
450}\r
451\r
452Double_t AliV0Reader::GetMotherCandidatePhi() const{\r
453 //see header file for documentation\r
454\r
455 Double_t offset=0;\r
456 if(fMotherCandidateLorentzVector->Phi()> TMath::Pi()){\r
457 offset = -2*TMath::Pi();\r
458 }\r
459 return fMotherCandidateLorentzVector->Phi()+offset;\r
460}\r
461\r
462Int_t AliV0Reader::GetSpeciesIndex(Int_t chargeOfTrack){\r
463 //see header file for documentation\r
464\r
465 Int_t iResult = 10; // Unknown particle\r
466\r
467 if(chargeOfTrack==-1){ //negative track\r
468 switch(abs(fNegativeTrackPID)){\r
469 case 11: //electron\r
470 iResult = 0;\r
471 break;\r
472 case 13: //muon\r
473 iResult = 1;\r
474 break;\r
475 case 211: //pion\r
476 iResult = 2;\r
477 break;\r
478 case 321: //kaon\r
479 iResult = 3;\r
480 break;\r
481 case 2212: //proton\r
482 iResult = 4;\r
483 break;\r
484 case 22: //photon\r
485 iResult = 5;\r
486 break;\r
487 case 111: //pi0\r
488 iResult = 6;\r
489 break;\r
490 case 2112: //neutron\r
491 iResult = 7;\r
492 break;\r
493 case 311: //K0\r
494 iResult = 8;\r
495 break;\r
496 \r
497 //Put in here for kSPECIES::kEleCon ????\r
498 }\r
499 }\r
500 else if(chargeOfTrack==1){ //positive track\r
501 switch(abs(fPositiveTrackPID)){\r
502 case 11: //electron\r
503 iResult = 0;\r
504 break;\r
505 case 13: //muon\r
506 iResult = 1;\r
507 break;\r
508 case 211: //pion\r
509 iResult = 2;\r
510 break;\r
511 case 321: //kaon\r
512 iResult = 3;\r
513 break;\r
514 case 2212: //proton\r
515 iResult = 4;\r
516 break;\r
517 case 22: //photon\r
518 iResult = 5;\r
519 break;\r
520 case 111: //pi0\r
521 iResult = 6;\r
522 break;\r
523 case 2112: //neutron\r
524 iResult = 7;\r
525 break;\r
526 case 311: //K0\r
527 iResult = 8;\r
528 break;\r
529\r
530 //Put in here for kSPECIES::kEleCon ????\r
531 }\r
532 }\r
533 else{\r
534 //Wrong parameter.. Print warning\r
535 }\r
536 return iResult;\r
537}\r