bugs corrected
[u/mrichter/AliRoot.git] / PWG4 / PartCorr / AliV0Reader.cxx
CommitLineData
9fdb2372 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
a6f44de5 80 fChi2CutConversion(0.),\r
81 fChi2CutMeson(0.),\r
9fdb2372 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
a6f44de5 130 fChi2CutConversion(original.fChi2CutConversion),\r
131 fChi2CutMeson(original.fChi2CutMeson),\r
9fdb2372 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
190\r
191\r
192Bool_t AliV0Reader::NextV0(){\r
193 //see header file for documentation\r
194\r
195 Bool_t iResult=kFALSE;\r
196 while(fCurrentV0IndexNumber<fESDEvent->GetNumberOfV0s()){\r
197 fCurrentV0 = fESDEvent->GetV0(fCurrentV0IndexNumber);\r
198 \r
199 //checks if on the fly mode is set\r
200 if ( !fCurrentV0->GetOnFlyStatus() ){\r
201 fCurrentV0IndexNumber++;\r
202 fHistograms->FillHistogram("V0MassDebugCut1",GetMotherCandidateMass());\r
203 continue;\r
204 }\r
205\r
206 if(fESDEvent->GetPrimaryVertex()->GetNContributors()<=0) {//checks if we have a vertex\r
207 fCurrentV0IndexNumber++;\r
208 fHistograms->FillHistogram("V0MassDebugCut2",GetMotherCandidateMass());\r
209 continue;\r
210 }\r
211\r
212 if(CheckPIDProbability(fPIDProbabilityCutNegativeParticle,fPIDProbabilityCutPositiveParticle)==kFALSE){\r
213 fCurrentV0IndexNumber++;\r
214 fHistograms->FillHistogram("V0MassDebugCut3",GetMotherCandidateMass());\r
215 continue;\r
216 }\r
217\r
218\r
219 fCurrentV0->GetXYZ(fCurrentXValue,fCurrentYValue,fCurrentZValue);\r
220 \r
221 if(GetXYRadius()>fMaxR){ // cuts on distance from collision point\r
222 fCurrentV0IndexNumber++;\r
223 fHistograms->FillHistogram("V0MassDebugCut4",GetMotherCandidateMass());\r
224 continue;\r
225 }\r
226\r
227 UpdateV0Information();\r
228 \r
229 if(fUseKFParticle){\r
230 if(fCurrentMotherKFCandidate->GetNDF()<=0){\r
231 fCurrentV0IndexNumber++;\r
232 fHistograms->FillHistogram("V0MassDebugCut5",GetMotherCandidateMass());\r
233 continue;\r
234 }\r
235 Double_t chi2V0 = fCurrentMotherKFCandidate->GetChi2()/fCurrentMotherKFCandidate->GetNDF();\r
a6f44de5 236 if(chi2V0 > fChi2CutConversion || chi2V0 <=0){\r
9fdb2372 237 fCurrentV0IndexNumber++;\r
238 fHistograms->FillHistogram("V0MassDebugCut6",GetMotherCandidateMass());\r
239 continue;\r
240 }\r
241 \r
242 if(TMath::Abs(fMotherCandidateLorentzVector->Eta())> fEtaCut){\r
243 fCurrentV0IndexNumber++;\r
244 fHistograms->FillHistogram("V0MassDebugCut7",GetMotherCandidateMass());\r
245 continue;\r
246 }\r
247 \r
248 if(fMotherCandidateLorentzVector->Pt()<fPtCut){\r
249 fCurrentV0IndexNumber++;\r
250 fHistograms->FillHistogram("V0MassDebugCut8",GetMotherCandidateMass());\r
251 continue;\r
252 }\r
253 \r
254 }\r
255 else if(fUseESDTrack){\r
256 //TODO\r
257 }\r
258\r
259 iResult=kTRUE;//means we have a v0 who survived all the cuts applied\r
260\r
261 fCurrentV0IndexNumber++;\r
262 \r
263 break;\r
264 }\r
265 return iResult; \r
266}\r
267\r
268void AliV0Reader::UpdateV0Information(){\r
269 //see header file for documentation\r
270 \r
271 if(fCurrentNegativeKFParticle != NULL){\r
272 delete fCurrentNegativeKFParticle;\r
273 }\r
274 fCurrentNegativeKFParticle = new AliKFParticle(*(fCurrentV0->GetParamN()),fNegativeTrackPID);\r
275 \r
276 if(fCurrentPositiveKFParticle != NULL){\r
277 delete fCurrentPositiveKFParticle;\r
278 }\r
279 fCurrentPositiveKFParticle = new AliKFParticle(*(fCurrentV0->GetParamP()),fPositiveTrackPID);\r
280 \r
281 if(fCurrentMotherKFCandidate != NULL){\r
282 delete fCurrentMotherKFCandidate;\r
283 }\r
284 fCurrentMotherKFCandidate = new AliKFParticle(*fCurrentNegativeKFParticle,*fCurrentPositiveKFParticle);\r
285\r
286 fCurrentNegativeESDTrack = fESDEvent->GetTrack(fCurrentV0->GetNindex());\r
287\r
288 fCurrentPositiveESDTrack = fESDEvent->GetTrack(fCurrentV0->GetPindex());\r
289\r
290 if(fPositiveTrackPID==-11 && fNegativeTrackPID==11){\r
291 fCurrentMotherKFCandidate->SetMassConstraint(0,fNSigmaMass);\r
292 }\r
293\r
294 if(fUseImprovedVertex == kTRUE){\r
295 AliKFVertex primaryVertexImproved(*GetPrimaryVertex());\r
296 primaryVertexImproved+=*fCurrentMotherKFCandidate;\r
297 fCurrentMotherKFCandidate->SetProductionVertex(primaryVertexImproved);\r
298 }\r
299\r
300 fCurrentMotherKFCandidate->GetMass(fMotherCandidateKFMass,fMotherCandidateKFWidth);\r
301\r
302\r
303 if(fNegativeTrackLorentzVector != NULL){\r
304 delete fNegativeTrackLorentzVector;\r
305 }\r
306 if(fUseKFParticle){\r
307 fNegativeTrackLorentzVector = new TLorentzVector(fCurrentNegativeKFParticle->Px(),fCurrentNegativeKFParticle->Py(),fCurrentNegativeKFParticle->Pz());\r
308 }\r
309 else if(fUseESDTrack){\r
310 fNegativeTrackLorentzVector = new TLorentzVector(fCurrentNegativeESDTrack->Px(),fCurrentNegativeESDTrack->Py(),fCurrentNegativeESDTrack->Pz());\r
311 }\r
312\r
313 if(fPositiveTrackLorentzVector != NULL){\r
314 delete fPositiveTrackLorentzVector;\r
315 }\r
316 if(fUseKFParticle){\r
317 fPositiveTrackLorentzVector = new TLorentzVector(fCurrentPositiveKFParticle->Px(),fCurrentPositiveKFParticle->Py(),fCurrentPositiveKFParticle->Pz());\r
318 }\r
319 else if(fUseESDTrack){\r
320 fPositiveTrackLorentzVector = new TLorentzVector(fCurrentPositiveESDTrack->Px(),fCurrentPositiveESDTrack->Py(),fCurrentPositiveESDTrack->Pz());\r
321 }\r
322\r
323 if(fMotherCandidateLorentzVector != NULL){\r
324 delete fMotherCandidateLorentzVector;\r
325 }\r
326 if(fUseKFParticle){\r
327 fMotherCandidateLorentzVector = new TLorentzVector(*fNegativeTrackLorentzVector + *fPositiveTrackLorentzVector);\r
328 }\r
329 else if(fUseESDTrack){\r
330 fMotherCandidateLorentzVector = new TLorentzVector(*fNegativeTrackLorentzVector + *fPositiveTrackLorentzVector);\r
331 }\r
332\r
333 if(fPositiveTrackPID==-11 && fNegativeTrackPID==11){\r
334 fMotherCandidateLorentzVector->SetXYZM(fMotherCandidateLorentzVector->Px() ,fMotherCandidateLorentzVector->Py(),fMotherCandidateLorentzVector->Pz(),0.); \r
335 }\r
336 \r
a6f44de5 337 if(fDoMC == kTRUE){\r
9fdb2372 338 fNegativeMCParticle = fMCStack->Particle(TMath::Abs(fESDEvent->GetTrack(fCurrentV0->GetNindex())->GetLabel()));\r
339 fPositiveMCParticle = fMCStack->Particle(TMath::Abs(fESDEvent->GetTrack(fCurrentV0->GetPindex())->GetLabel()));\r
340 }\r
a6f44de5 341 fCurrentEventGoodV0s.push_back(*fCurrentMotherKFCandidate);\r
9fdb2372 342}\r
343\r
344Bool_t AliV0Reader::HasSameMCMother(){\r
345 //see header file for documentation\r
346\r
347 Bool_t iResult = kFALSE;\r
348 if(fDoMC == kTRUE){\r
349 if(fNegativeMCParticle != NULL && fPositiveMCParticle != NULL){\r
350 if(fNegativeMCParticle->GetMother(0) == fPositiveMCParticle->GetMother(0))\r
351 fMotherMCParticle = fMCStack->Particle(fPositiveMCParticle->GetMother(0));\r
352 iResult = kTRUE;\r
353 }\r
354 }\r
355 return iResult;\r
356}\r
357\r
358Bool_t AliV0Reader::CheckPIDProbability(Double_t negProbCut, Double_t posProbCut){\r
359 //see header file for documentation\r
360\r
361 Bool_t iResult=kFALSE;\r
362\r
363 Double_t *posProbArray = new Double_t[10];\r
364 Double_t *negProbArray = new Double_t[10];\r
365 AliESDtrack* negTrack = fESDEvent->GetTrack(fCurrentV0->GetNindex());\r
366 AliESDtrack* posTrack = fESDEvent->GetTrack(fCurrentV0->GetPindex());\r
367 \r
368 negTrack->GetTPCpid(negProbArray);\r
369 posTrack->GetTPCpid(posProbArray);\r
370\r
371 if(negProbArray!=NULL && posProbArray!=NULL){\r
372 if(negProbArray[GetSpeciesIndex(-1)]>=negProbCut && posProbArray[GetSpeciesIndex(1)]>=posProbCut){\r
373 iResult=kTRUE;\r
374 }\r
375 }\r
376 delete [] posProbArray;\r
377 delete [] negProbArray;\r
378 return iResult;\r
379}\r
380\r
381void AliV0Reader::UpdateEventByEventData(){\r
382 //see header file for documentation\r
383\r
384 if(fCurrentEventGoodV0s.size() >0 ){\r
385 fPreviousEventGoodV0s.clear();\r
386 fPreviousEventGoodV0s = fCurrentEventGoodV0s;\r
387 }\r
388 fCurrentEventGoodV0s.clear();\r
389 \r
390 fCurrentV0IndexNumber=0;\r
391}\r
392\r
393Double_t AliV0Reader::GetNegativeTrackPhi() const{\r
394 //see header file for documentation\r
395\r
396 Double_t offset=0;\r
397 if(fNegativeTrackLorentzVector->Phi()> TMath::Pi()){\r
398 offset = -2*TMath::Pi();\r
399 }\r
400 return fNegativeTrackLorentzVector->Phi()+offset;\r
401}\r
402\r
403Double_t AliV0Reader::GetPositiveTrackPhi() const{\r
404 //see header file for documentation\r
405\r
406 Double_t offset=0;\r
407 if(fPositiveTrackLorentzVector->Phi()> TMath::Pi()){\r
408 offset = -2*TMath::Pi();\r
409 }\r
410 return fPositiveTrackLorentzVector->Phi()+offset;\r
411}\r
412\r
413Double_t AliV0Reader::GetMotherCandidatePhi() const{\r
414 //see header file for documentation\r
415\r
416 Double_t offset=0;\r
417 if(fMotherCandidateLorentzVector->Phi()> TMath::Pi()){\r
418 offset = -2*TMath::Pi();\r
419 }\r
420 return fMotherCandidateLorentzVector->Phi()+offset;\r
421}\r
422\r
423Int_t AliV0Reader::GetSpeciesIndex(Int_t chargeOfTrack){\r
424 //see header file for documentation\r
425\r
426 Int_t iResult = 10; // Unknown particle\r
427\r
428 if(chargeOfTrack==-1){ //negative track\r
429 switch(abs(fNegativeTrackPID)){\r
430 case 11: //electron\r
431 iResult = 0;\r
432 break;\r
433 case 13: //muon\r
434 iResult = 1;\r
435 break;\r
436 case 211: //pion\r
437 iResult = 2;\r
438 break;\r
439 case 321: //kaon\r
440 iResult = 3;\r
441 break;\r
442 case 2212: //proton\r
443 iResult = 4;\r
444 break;\r
445 case 22: //photon\r
446 iResult = 5;\r
447 break;\r
448 case 111: //pi0\r
449 iResult = 6;\r
450 break;\r
451 case 2112: //neutron\r
452 iResult = 7;\r
453 break;\r
454 case 311: //K0\r
455 iResult = 8;\r
456 break;\r
457 \r
458 //Put in here for kSPECIES::kEleCon ????\r
459 }\r
460 }\r
461 else if(chargeOfTrack==1){ //positive track\r
462 switch(abs(fPositiveTrackPID)){\r
463 case 11: //electron\r
464 iResult = 0;\r
465 break;\r
466 case 13: //muon\r
467 iResult = 1;\r
468 break;\r
469 case 211: //pion\r
470 iResult = 2;\r
471 break;\r
472 case 321: //kaon\r
473 iResult = 3;\r
474 break;\r
475 case 2212: //proton\r
476 iResult = 4;\r
477 break;\r
478 case 22: //photon\r
479 iResult = 5;\r
480 break;\r
481 case 111: //pi0\r
482 iResult = 6;\r
483 break;\r
484 case 2112: //neutron\r
485 iResult = 7;\r
486 break;\r
487 case 311: //K0\r
488 iResult = 8;\r
489 break;\r
490\r
491 //Put in here for kSPECIES::kEleCon ????\r
492 }\r
493 }\r
494 else{\r
495 //Wrong parameter.. Print warning\r
496 }\r
497 return iResult;\r
498}\r