]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWG4/GammaConv/AliV0Reader.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Bugfix, fixed memory leaks.
[u/mrichter/AliRoot.git] / PWG4 / GammaConv / AliV0Reader.cxx
CommitLineData
a0b94e5c 1/**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: Ana Marin, Kathrin Koch, Kenneth Aamodt *
5 * Version 1.0 *
6 * *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
15
16////////////////////////////////////////////////
17//---------------------------------------------
18// Class used to do analysis on conversion pairs
19//---------------------------------------------
20////////////////////////////////////////////////
21
22// --- ROOT system ---
23#include <TMath.h>
24
25//---- ANALYSIS system ----
26#include "AliV0Reader.h"
27#include "AliAnalysisManager.h"
28#include "AliESDInputHandler.h"
29#include "AliESDtrack.h"
30#include "AliMCEvent.h"
31#include "AliKFVertex.h"
32
33#include "AliStack.h"
34#include "AliMCEventHandler.h"
9640a3d1 35#include "AliTPCpidESD.h"
a0b94e5c 36
37class iostream;
38class AliESDv0;
39class TFormula;
40
41using namespace std;
42
43ClassImp(AliV0Reader)
44
45
46
47AliV0Reader::AliV0Reader() :
48 TObject(),
49 fMCStack(NULL),
50 fMCTruth(NULL),
51 fMCEvent(NULL), // for CF
52 fChain(NULL),
53 fESDHandler(NULL),
54 fESDEvent(NULL),
4a6157dc 55 fCFManager(NULL),
9640a3d1 56 fTPCpid(NULL),
a0b94e5c 57 fHistograms(NULL),
58 fCurrentV0IndexNumber(0),
59 fCurrentV0(NULL),
60 fCurrentNegativeKFParticle(NULL),
61 fCurrentPositiveKFParticle(NULL),
62 fCurrentMotherKFCandidate(NULL),
63 fCurrentNegativeESDTrack(NULL),
64 fCurrentPositiveESDTrack(NULL),
65 fNegativeTrackLorentzVector(NULL),
66 fPositiveTrackLorentzVector(NULL),
67 fMotherCandidateLorentzVector(NULL),
68 fCurrentXValue(0),
69 fCurrentYValue(0),
70 fCurrentZValue(0),
71 fPositiveTrackPID(0),
72 fNegativeTrackPID(0),
73 fNegativeMCParticle(NULL),
74 fPositiveMCParticle(NULL),
75 fMotherMCParticle(NULL),
76 fMotherCandidateKFMass(0),
77 fMotherCandidateKFWidth(0),
78 fUseKFParticle(kTRUE),
79 fUseESDTrack(kFALSE),
80 fDoMC(kFALSE),
81 fMaxR(10000),// 100 meter(outside of ALICE)
82 fEtaCut(0.),
83 fPtCut(0.),
84 fMaxZ(0.),
85 fLineCutZRSlope(0.),
86 fLineCutZValue(0.),
87 fChi2CutConversion(0.),
88 fChi2CutMeson(0.),
89 fPIDProbabilityCutNegativeParticle(0),
90 fPIDProbabilityCutPositiveParticle(0),
9640a3d1 91 fDodEdxSigmaCut(kFALSE),
92 fPIDnSigmaAboveElectronLine(100),
93 fPIDnSigmaBelowElectronLine(-100),
94 fPIDnSigmaAbovePionLine(-100),
95 fPIDMinPnSigmaAbovePionLine(100),
a0b94e5c 96 fXVertexCut(0.),
97 fYVertexCut(0.),
98 fZVertexCut(0.),
99 fNSigmaMass(0.),
100 fUseImprovedVertex(kFALSE),
101 fUseOwnXYZCalculation(kFALSE),
1e7846f4 102 fDoCF(kFALSE),
a0b94e5c 103 fCurrentEventGoodV0s(),
104 fPreviousEventGoodV0s()
105{
9640a3d1 106 fTPCpid = new AliTPCpidESD;
a0b94e5c 107}
108
109
110AliV0Reader::AliV0Reader(const AliV0Reader & original) :
111 TObject(original),
112 fMCStack(original.fMCStack),
113 fMCTruth(original.fMCTruth),
114 fMCEvent(original.fMCEvent), // for CF
115 fChain(original.fChain),
116 fESDHandler(original.fESDHandler),
117 fESDEvent(original.fESDEvent),
4a6157dc 118 fCFManager(original.fCFManager),
9640a3d1 119 fTPCpid(original.fTPCpid),
a0b94e5c 120 fHistograms(original.fHistograms),
121 fCurrentV0IndexNumber(original.fCurrentV0IndexNumber),
122 fCurrentV0(original.fCurrentV0),
123 fCurrentNegativeKFParticle(original.fCurrentNegativeKFParticle),
124 fCurrentPositiveKFParticle(original.fCurrentPositiveKFParticle),
125 fCurrentMotherKFCandidate(original.fCurrentMotherKFCandidate),
126 fCurrentNegativeESDTrack(original.fCurrentNegativeESDTrack),
127 fCurrentPositiveESDTrack(original.fCurrentPositiveESDTrack),
128 fNegativeTrackLorentzVector(original.fNegativeTrackLorentzVector),
129 fPositiveTrackLorentzVector(original.fPositiveTrackLorentzVector),
130 fMotherCandidateLorentzVector(original.fMotherCandidateLorentzVector),
131 fCurrentXValue(original.fCurrentXValue),
132 fCurrentYValue(original.fCurrentYValue),
133 fCurrentZValue(original.fCurrentZValue),
134 fPositiveTrackPID(original.fPositiveTrackPID),
135 fNegativeTrackPID(original.fNegativeTrackPID),
136 fNegativeMCParticle(original.fNegativeMCParticle),
137 fPositiveMCParticle(original.fPositiveMCParticle),
138 fMotherMCParticle(original.fMotherMCParticle),
139 fMotherCandidateKFMass(original.fMotherCandidateKFMass),
140 fMotherCandidateKFWidth(original.fMotherCandidateKFWidth),
141 fUseKFParticle(kTRUE),
142 fUseESDTrack(kFALSE),
143 fDoMC(kFALSE),
144 fMaxR(original.fMaxR),
145 fEtaCut(original.fEtaCut),
146 fPtCut(original.fPtCut),
147 fMaxZ(original.fMaxZ),
148 fLineCutZRSlope(original.fLineCutZRSlope),
149 fLineCutZValue(original.fLineCutZValue),
150 fChi2CutConversion(original.fChi2CutConversion),
151 fChi2CutMeson(original.fChi2CutMeson),
152 fPIDProbabilityCutNegativeParticle(original.fPIDProbabilityCutNegativeParticle),
153 fPIDProbabilityCutPositiveParticle(original.fPIDProbabilityCutPositiveParticle),
9640a3d1 154 fDodEdxSigmaCut(original.fDodEdxSigmaCut),
155 fPIDnSigmaAboveElectronLine(original.fPIDnSigmaAboveElectronLine),
156 fPIDnSigmaBelowElectronLine(original.fPIDnSigmaBelowElectronLine),
157 fPIDnSigmaAbovePionLine(original.fPIDnSigmaAbovePionLine),
158 fPIDMinPnSigmaAbovePionLine(original.fPIDMinPnSigmaAbovePionLine),
a0b94e5c 159 fXVertexCut(original.fXVertexCut),
160 fYVertexCut(original.fYVertexCut),
161 fZVertexCut(original.fZVertexCut),
162 fNSigmaMass(original.fNSigmaMass),
163 fUseImprovedVertex(original.fUseImprovedVertex),
164 fUseOwnXYZCalculation(original.fUseOwnXYZCalculation),
1e7846f4 165 fDoCF(original.fDoCF),
a0b94e5c 166 fCurrentEventGoodV0s(original.fCurrentEventGoodV0s),
167 fPreviousEventGoodV0s(original.fPreviousEventGoodV0s)
168{
169
170}
171
172
173AliV0Reader & AliV0Reader::operator = (const AliV0Reader & /*source*/)
174{
175 // assignment operator
176 return *this;
177}
9640a3d1 178AliV0Reader::~AliV0Reader()
179{
180 if(fTPCpid){
181 delete fTPCpid;
182 }
183}
a0b94e5c 184
185void AliV0Reader::Initialize(){
186 //see header file for documentation
187
188 // Get the input handler from the manager
189 fESDHandler = (AliESDInputHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
190 if(fESDHandler == NULL){
191 //print warning here
192 }
193
194 // Get pointer to esd event from input handler
195 fESDEvent = fESDHandler->GetEvent();
196 if(fESDEvent == NULL){
197 //print warning here
198 }
199
200 //Get pointer to MCTruth
201 fMCTruth = (AliMCEventHandler*)((AliAnalysisManager::GetAnalysisManager())->GetMCtruthEventHandler());
202 if(fMCTruth == NULL){
203 //print warning here
204 }
205
206 //Get pointer to the mc stack
207 fMCStack = fMCTruth->MCEvent()->Stack();
208 if(fMCStack == NULL){
209 //print warning here
210 }
211
212
213 // for CF
214 //Get pointer to the mc event
1e7846f4 215 if(fDoCF){
216 fMCEvent = fMCTruth->MCEvent();
217 if(fMCEvent == NULL){
218 //print warning here
219 fDoCF = kFALSE;
220 }
221 }
a0b94e5c 222
223 AliKFParticle::SetField(fESDEvent->GetMagneticField());
224
225}
226
227AliESDv0* AliV0Reader::GetV0(Int_t index){
228 //see header file for documentation
229 fCurrentV0 = fESDEvent->GetV0(index);
230 UpdateV0Information();
231 return fCurrentV0;
232}
233
234Bool_t AliV0Reader::CheckForPrimaryVertex(){
235 return fESDEvent->GetPrimaryVertex()->GetNContributors()>0;
236}
237
238
239
240Bool_t AliV0Reader::NextV0(){
241 //see header file for documentation
242
243 Bool_t iResult=kFALSE;
244 while(fCurrentV0IndexNumber<fESDEvent->GetNumberOfV0s()){
245 fCurrentV0 = fESDEvent->GetV0(fCurrentV0IndexNumber);
246
247 // moved it up here so that the correction framework can access pt and eta information
248 if(UpdateV0Information() == kFALSE){
249 fCurrentV0IndexNumber++;
250 continue;
251 }
252
253 Double_t containerInput[3];
1e7846f4 254 if(fDoCF){
255 containerInput[0] = GetMotherCandidatePt();
256 containerInput[1] = GetMotherCandidateEta();
257 containerInput[2] = GetMotherCandidateMass();
258 }
a0b94e5c 259
260 //checks if on the fly mode is set
261 if ( !fCurrentV0->GetOnFlyStatus() ){
262 if(fHistograms != NULL){
263 fHistograms->FillHistogram("ESD_CutGetOnFly_InvMass",GetMotherCandidateMass());
264 }
265 fCurrentV0IndexNumber++;
266 continue;
267 }
1e7846f4 268 if(fDoCF){
269 fCFManager->GetParticleContainer()->Fill(containerInput,kStepGetOnFly); // for CF
270 }
271
a0b94e5c 272 //checks if we have a prim vertex
273 if(fESDEvent->GetPrimaryVertex()->GetNContributors()<=0) {
274 if(fHistograms != NULL){
275 fHistograms->FillHistogram("ESD_CutNContributors_InvMass",GetMotherCandidateMass());
276 }
277 fCurrentV0IndexNumber++;
278 continue;
279 }
1e7846f4 280 if(fDoCF){
281 fCFManager->GetParticleContainer()->Fill(containerInput,kStepNContributors); // for CF
282 }
a0b94e5c 283
284 //Check the pid probability
285 if(CheckPIDProbability(fPIDProbabilityCutNegativeParticle,fPIDProbabilityCutPositiveParticle)==kFALSE){
286 if(fHistograms != NULL){
287 fHistograms->FillHistogram("ESD_CutPIDProb_InvMass",GetMotherCandidateMass());
288 }
289 fCurrentV0IndexNumber++;
290 continue;
291 }
1e7846f4 292 if(fDoCF){
293 fCFManager->GetParticleContainer()->Fill(containerInput,kStepTPCPID); // for CF
294 }
a0b94e5c 295
296
4a6157dc 297 if(fUseOwnXYZCalculation == kFALSE){
298 fCurrentV0->GetXYZ(fCurrentXValue,fCurrentYValue,fCurrentZValue);
299 }
300 else{
301 Double_t convpos[2];
302 convpos[0]=0;
303 convpos[1]=0;
304 GetConvPosXY(GetPositiveESDTrack(),GetNegativeESDTrack(),GetMagneticField(),convpos);
305 fCurrentXValue = convpos[0];
306 fCurrentYValue = convpos[1];
307 fCurrentZValue = GetConvPosZ(GetPositiveESDTrack(),GetNegativeESDTrack(),GetMagneticField());
308 }
a0b94e5c 309
310 if(GetXYRadius()>fMaxR){ // cuts on distance from collision point
311 if(fHistograms != NULL){
312 fHistograms->FillHistogram("ESD_CutR_InvMass",GetMotherCandidateMass());
313 }
314 fCurrentV0IndexNumber++;
315 continue;
1e7846f4 316 }
317 if(fDoCF){
318 fCFManager->GetParticleContainer()->Fill(containerInput,kStepR); // for CF
319 }
a0b94e5c 320
321
322 if((TMath::Abs(fCurrentZValue)*fLineCutZRSlope)-fLineCutZValue > GetXYRadius() ){ // cuts out regions where we do not reconstruct
323 if(fHistograms != NULL){
324 fHistograms->FillHistogram("ESD_CutLine_InvMass",GetMotherCandidateMass());
325 }
326 fCurrentV0IndexNumber++;
327 continue;
1e7846f4 328 }
329 if(fDoCF){
330 fCFManager->GetParticleContainer()->Fill(containerInput,kStepLine); // for CF
331 }
a0b94e5c 332
333 if(TMath::Abs(fCurrentZValue) > fMaxZ ){ // cuts out regions where we do not reconstruct
334 if(fHistograms != NULL){
335 fHistograms->FillHistogram("ESD_CutZ_InvMass",GetMotherCandidateMass());
336 }
337 fCurrentV0IndexNumber++;
338 continue;
1e7846f4 339 }
340 if(fDoCF){
341 fCFManager->GetParticleContainer()->Fill(containerInput,kStepZ); // for CF
342 }
a0b94e5c 343
344 /* Moved further up so corr framework can work
345 if(UpdateV0Information() == kFALSE){
346 fCurrentV0IndexNumber++;
347 continue;
348 }
349 */
9640a3d1 350
a0b94e5c 351
352 if(fUseKFParticle){
353 if(fCurrentMotherKFCandidate->GetNDF()<=0){
354 if(fHistograms != NULL){
355 fHistograms->FillHistogram("ESD_CutNDF_InvMass",GetMotherCandidateMass());
356 }
357 fCurrentV0IndexNumber++;
358 continue;
359 }
1e7846f4 360 if(fDoCF){
361 fCFManager->GetParticleContainer()->Fill(containerInput,kStepNDF); // for CF
362 }
a0b94e5c 363
364 Double_t chi2V0 = fCurrentMotherKFCandidate->GetChi2()/fCurrentMotherKFCandidate->GetNDF();
365 if(chi2V0 > fChi2CutConversion || chi2V0 <=0){
366 if(fHistograms != NULL){
367 fHistograms->FillHistogram("ESD_CutChi2_InvMass",GetMotherCandidateMass());
368 }
369 fCurrentV0IndexNumber++;
370 continue;
371 }
1e7846f4 372 if(fDoCF){
373 fCFManager->GetParticleContainer()->Fill(containerInput,kStepChi2); // for CF
374 }
a0b94e5c 375
376 if(TMath::Abs(fMotherCandidateLorentzVector->Eta())> fEtaCut){
377 if(fHistograms != NULL){
378 fHistograms->FillHistogram("ESD_CutEta_InvMass",GetMotherCandidateMass());
379 }
380 fCurrentV0IndexNumber++;
381 continue;
382 }
1e7846f4 383 if(fDoCF){
384 fCFManager->GetParticleContainer()->Fill(containerInput,kStepEta); // for CF
385 }
a0b94e5c 386
387 if(fMotherCandidateLorentzVector->Pt()<fPtCut){
388 if(fHistograms != NULL){
389 fHistograms->FillHistogram("ESD_CutPt_InvMass",GetMotherCandidateMass());
390 }
391 fCurrentV0IndexNumber++;
392 continue;
393 }
1e7846f4 394 if(fDoCF){
395 fCFManager->GetParticleContainer()->Fill(containerInput,kStepPt); // for CF
396 }
a0b94e5c 397
398 }
399 else if(fUseESDTrack){
400 //TODO
401 }
9640a3d1 402
403 if(fHistograms != NULL){
404 fHistograms->FillHistogram("ESD_GoodV0s_InvMass",GetMotherCandidateMass());
405 }
406
a0b94e5c 407 fCurrentEventGoodV0s.push_back(*fCurrentMotherKFCandidate);
408
409 iResult=kTRUE;//means we have a v0 who survived all the cuts applied
410
411 fCurrentV0IndexNumber++;
412
413 break;
414 }
415 return iResult;
416}
417
418Bool_t AliV0Reader::UpdateV0Information(){
419 //see header file for documentation
420
421 Bool_t iResult=kTRUE; // for taking out not refitted, kinks and like sign tracks
422
423 Bool_t switchTracks = kFALSE;
424
425 fCurrentNegativeESDTrack = fESDEvent->GetTrack(fCurrentV0->GetNindex());
426 fCurrentPositiveESDTrack = fESDEvent->GetTrack(fCurrentV0->GetPindex());
427
428 if(fCurrentNegativeESDTrack->GetSign() == fCurrentPositiveESDTrack->GetSign()){ // avoid like sign
429 iResult=kFALSE;
430 if(fHistograms != NULL){
431 fHistograms->FillHistogram("ESD_CutLikeSign_InvMass",GetMotherCandidateMass());
432 }
433 }
434
435 if(fCurrentPositiveESDTrack->GetSign() == -1 && fCurrentNegativeESDTrack->GetSign() == 1){ // switch wrong signed tracks
436 fCurrentNegativeESDTrack = fESDEvent->GetTrack(fCurrentV0->GetPindex());
437 fCurrentPositiveESDTrack = fESDEvent->GetTrack(fCurrentV0->GetNindex());
438 switchTracks = kTRUE;
439 }
440
441 if( !(fCurrentNegativeESDTrack->GetStatus() & AliESDtrack::kTPCrefit) ||
442 !(fCurrentPositiveESDTrack->GetStatus() & AliESDtrack::kTPCrefit) ){
443 // if( !(fCurrentNegativeESDTrack->GetStatus() & AliESDtrack::kITSrefit) ||
444 // !(fCurrentPositiveESDTrack->GetStatus() & AliESDtrack::kITSrefit) ){
445
446 iResult=kFALSE;
447 if(fHistograms != NULL){
448 fHistograms->FillHistogram("ESD_CutRefit_InvMass",GetMotherCandidateMass());
449 }
450 }
451
452 if( fCurrentNegativeESDTrack->GetKinkIndex(0) > 0 ||
453 fCurrentPositiveESDTrack->GetKinkIndex(0) > 0) {
454
455 iResult=kFALSE;
456 if(fHistograms != NULL){
457 fHistograms->FillHistogram("ESD_CutKink_InvMass",GetMotherCandidateMass());
458 }
459 }
9640a3d1 460
461 if(fDodEdxSigmaCut == kTRUE){
462
463 if( fTPCpid->GetNumberOfSigmas(fCurrentPositiveESDTrack,AliPID::kElectron)<fPIDnSigmaBelowElectronLine ||
464 fTPCpid->GetNumberOfSigmas(fCurrentPositiveESDTrack,AliPID::kElectron)>fPIDnSigmaAboveElectronLine ||
465 fTPCpid->GetNumberOfSigmas(fCurrentNegativeESDTrack,AliPID::kElectron)<fPIDnSigmaBelowElectronLine ||
466 fTPCpid->GetNumberOfSigmas(fCurrentNegativeESDTrack,AliPID::kElectron)>fPIDnSigmaAboveElectronLine ){
467 iResult=kFALSE;
468 if(fHistograms != NULL){
469 fHistograms->FillHistogram("ESD_CutdEdxSigmaElectronLine_InvMass",GetMotherCandidateMass());
470 }
471 }
472 if( fCurrentPositiveESDTrack->P()>fPIDMinPnSigmaAbovePionLine){
473 if(fTPCpid->GetNumberOfSigmas(fCurrentPositiveESDTrack,AliPID::kElectron)>fPIDnSigmaBelowElectronLine &&
474 fTPCpid->GetNumberOfSigmas(fCurrentPositiveESDTrack,AliPID::kElectron)<fPIDnSigmaAboveElectronLine&&
475 fTPCpid->GetNumberOfSigmas(fCurrentPositiveESDTrack,AliPID::kPion)<fPIDnSigmaAbovePionLine){
476 iResult=kFALSE;
477 if(fHistograms != NULL){
478 fHistograms->FillHistogram("ESD_CutdEdxSigmaPionLine_InvMass",GetMotherCandidateMass());
479 }
480 }
481 }
482
483 if( fCurrentNegativeESDTrack->P()>fPIDMinPnSigmaAbovePionLine){
484 if(fTPCpid->GetNumberOfSigmas(fCurrentNegativeESDTrack,AliPID::kElectron)>fPIDnSigmaBelowElectronLine &&
485 fTPCpid->GetNumberOfSigmas(fCurrentNegativeESDTrack,AliPID::kElectron)<fPIDnSigmaAboveElectronLine&&
486 fTPCpid->GetNumberOfSigmas(fCurrentNegativeESDTrack,AliPID::kPion)<fPIDnSigmaAbovePionLine){
487 iResult=kFALSE;
488 if(fHistograms != NULL){
489 fHistograms->FillHistogram("ESD_CutdEdxSigmaPionLine_InvMass",GetMotherCandidateMass());
490 }
491 }
492 }
493 }
494
495
496
a0b94e5c 497
498 if(fCurrentNegativeKFParticle != NULL){
499 delete fCurrentNegativeKFParticle;
500 }
501 if(switchTracks == kFALSE){
502 fCurrentNegativeKFParticle = new AliKFParticle(*(fCurrentV0->GetParamN()),fNegativeTrackPID);
503 }
504 else{
505 fCurrentNegativeKFParticle = new AliKFParticle(*(fCurrentV0->GetParamP()),fNegativeTrackPID);
506 }
507
508 if(fCurrentPositiveKFParticle != NULL){
509 delete fCurrentPositiveKFParticle;
510 }
511 if(switchTracks == kFALSE){
512 fCurrentPositiveKFParticle = new AliKFParticle(*(fCurrentV0->GetParamP()),fPositiveTrackPID);
513 }
514 else{
515 fCurrentPositiveKFParticle = new AliKFParticle(*(fCurrentV0->GetParamN()),fPositiveTrackPID);
516 }
517
518 if(fCurrentMotherKFCandidate != NULL){
519 delete fCurrentMotherKFCandidate;
520 }
521 fCurrentMotherKFCandidate = new AliKFParticle(*fCurrentNegativeKFParticle,*fCurrentPositiveKFParticle);
522
523
524 if(fPositiveTrackPID==-11 && fNegativeTrackPID==11){
525 fCurrentMotherKFCandidate->SetMassConstraint(0,fNSigmaMass);
526 }
527
528 if(fUseImprovedVertex == kTRUE){
529 AliKFVertex primaryVertexImproved(*GetPrimaryVertex());
530 primaryVertexImproved+=*fCurrentMotherKFCandidate;
531 fCurrentMotherKFCandidate->SetProductionVertex(primaryVertexImproved);
532 }
533
534 fCurrentMotherKFCandidate->GetMass(fMotherCandidateKFMass,fMotherCandidateKFWidth);
535
536
537 if(fNegativeTrackLorentzVector != NULL){
538 delete fNegativeTrackLorentzVector;
539 }
540 if(fUseKFParticle){
541 fNegativeTrackLorentzVector = new TLorentzVector(fCurrentNegativeKFParticle->Px(),fCurrentNegativeKFParticle->Py(),fCurrentNegativeKFParticle->Pz());
542 }
543 else if(fUseESDTrack){
544 fNegativeTrackLorentzVector = new TLorentzVector(fCurrentNegativeESDTrack->Px(),fCurrentNegativeESDTrack->Py(),fCurrentNegativeESDTrack->Pz());
545 }
546
547 if(fPositiveTrackLorentzVector != NULL){
548 delete fPositiveTrackLorentzVector;
549 }
550 if(fUseKFParticle){
551 fPositiveTrackLorentzVector = new TLorentzVector(fCurrentPositiveKFParticle->Px(),fCurrentPositiveKFParticle->Py(),fCurrentPositiveKFParticle->Pz());
552 }
553 else if(fUseESDTrack){
554 fPositiveTrackLorentzVector = new TLorentzVector(fCurrentPositiveESDTrack->Px(),fCurrentPositiveESDTrack->Py(),fCurrentPositiveESDTrack->Pz());
555 }
556
557 if(fMotherCandidateLorentzVector != NULL){
558 delete fMotherCandidateLorentzVector;
559 }
560 if(fUseKFParticle){
561 fMotherCandidateLorentzVector = new TLorentzVector(*fNegativeTrackLorentzVector + *fPositiveTrackLorentzVector);
562 }
563 else if(fUseESDTrack){
564 fMotherCandidateLorentzVector = new TLorentzVector(*fNegativeTrackLorentzVector + *fPositiveTrackLorentzVector);
565 }
566
567 if(fPositiveTrackPID==-11 && fNegativeTrackPID==11){
568 fMotherCandidateLorentzVector->SetXYZM(fMotherCandidateLorentzVector->Px() ,fMotherCandidateLorentzVector->Py(),fMotherCandidateLorentzVector->Pz(),0.);
569 }
570
571
572 if(fDoMC == kTRUE){
573 fMotherMCParticle= NULL;
574 fNegativeMCParticle = fMCStack->Particle(TMath::Abs(fESDEvent->GetTrack(fCurrentV0->GetNindex())->GetLabel()));
575 fPositiveMCParticle = fMCStack->Particle(TMath::Abs(fESDEvent->GetTrack(fCurrentV0->GetPindex())->GetLabel()));
576 if(fPositiveMCParticle->GetMother(0)>-1){
577 fMotherMCParticle = fMCStack->Particle(fPositiveMCParticle->GetMother(0));
578 }
579 }
580
581 // if(iResult==kTRUE){
582 // fCurrentEventGoodV0s.push_back(*fCurrentMotherKFCandidate); // moved it to NextV0() after all the cuts are applied
583 // }
584
585
586 // for CF
587 Double_t containerInput[3];
1e7846f4 588 if(fDoCF){
589 containerInput[0] = GetMotherCandidatePt();
590 containerInput[1] = GetMotherCandidateEta();
591 containerInput[2] = GetMotherCandidateMass();
592
593 fCFManager->GetParticleContainer()->Fill(containerInput,kStepLikeSign); // for CF
594 fCFManager->GetParticleContainer()->Fill(containerInput,kStepTPCRefit); // for CF
595 fCFManager->GetParticleContainer()->Fill(containerInput,kStepKinks); // for CF
596 }
a0b94e5c 597 return iResult;
598}
599
600
601
602Bool_t AliV0Reader::HasSameMCMother(){
603 //see header file for documentation
604
605 Bool_t iResult = kFALSE;
606 if(fDoMC == kTRUE){
607 if(fNegativeMCParticle != NULL && fPositiveMCParticle != NULL){
608 if(fNegativeMCParticle->GetMother(0) == fPositiveMCParticle->GetMother(0))
609 if(fMotherMCParticle){
610 iResult = kTRUE;
611 }
612 }
613 }
614 return iResult;
615}
616
617Bool_t AliV0Reader::CheckPIDProbability(Double_t negProbCut, Double_t posProbCut){
618 //see header file for documentation
619
620 Bool_t iResult=kFALSE;
621
622 Double_t *posProbArray = new Double_t[10];
623 Double_t *negProbArray = new Double_t[10];
624 AliESDtrack* negTrack = fESDEvent->GetTrack(fCurrentV0->GetNindex());
625 AliESDtrack* posTrack = fESDEvent->GetTrack(fCurrentV0->GetPindex());
626
627 negTrack->GetTPCpid(negProbArray);
628 posTrack->GetTPCpid(posProbArray);
629
4a6157dc 630 // if(negProbArray != NULL && posProbArray != NULL){ // this is not allowed anymore for some reason(RC19)
631 if(negProbArray && posProbArray){
a0b94e5c 632 if(negProbArray[GetSpeciesIndex(-1)]>=negProbCut && posProbArray[GetSpeciesIndex(1)]>=posProbCut){
633 iResult=kTRUE;
634 }
635 }
636 delete [] posProbArray;
637 delete [] negProbArray;
638 return iResult;
639}
640
641void AliV0Reader::GetPIDProbability(Double_t &negPIDProb,Double_t & posPIDProb){
4a6157dc 642 // see header file for documentation
643
a0b94e5c 644 Double_t *posProbArray = new Double_t[10];
645 Double_t *negProbArray = new Double_t[10];
646 AliESDtrack* negTrack = fESDEvent->GetTrack(fCurrentV0->GetNindex());
647 AliESDtrack* posTrack = fESDEvent->GetTrack(fCurrentV0->GetPindex());
648
649 negTrack->GetTPCpid(negProbArray);
650 posTrack->GetTPCpid(posProbArray);
651
4a6157dc 652 // if(negProbArray!=NULL && posProbArray!=NULL){ // this is not allowed anymore for some reason(RC19)
653 if(negProbArray && posProbArray){
a0b94e5c 654 negPIDProb = negProbArray[GetSpeciesIndex(-1)];
655 posPIDProb = posProbArray[GetSpeciesIndex(1)];
656 }
657 delete [] posProbArray;
658 delete [] negProbArray;
659}
660
661void AliV0Reader::UpdateEventByEventData(){
662 //see header file for documentation
663
664 if(fCurrentEventGoodV0s.size() >0 ){
665 // fPreviousEventGoodV0s.clear();
666 // fPreviousEventGoodV0s = fCurrentEventGoodV0s;
667 if(fPreviousEventGoodV0s.size()>19){
668 for(UInt_t nCurrent=0;nCurrent<fCurrentEventGoodV0s.size();nCurrent++){
669 fPreviousEventGoodV0s.erase(fPreviousEventGoodV0s.begin());
670 fPreviousEventGoodV0s.push_back(fCurrentEventGoodV0s.at(nCurrent));
671 }
672 }
673 else{
674 for(UInt_t nCurrent=0;nCurrent<fCurrentEventGoodV0s.size();nCurrent++){
675 if(fPreviousEventGoodV0s.size()<20){
676 fPreviousEventGoodV0s.push_back(fCurrentEventGoodV0s.at(nCurrent));
677 }
678 else{
679 fPreviousEventGoodV0s.erase(fPreviousEventGoodV0s.begin());
680 fPreviousEventGoodV0s.push_back(fCurrentEventGoodV0s.at(nCurrent));
681 }
682 }
683 }
684 }
685 fCurrentEventGoodV0s.clear();
686
687 fCurrentV0IndexNumber=0;
688}
689
690
691Double_t AliV0Reader::GetNegativeTrackPhi() const{
692 //see header file for documentation
693
694 Double_t offset=0;
695 if(fNegativeTrackLorentzVector->Phi()> TMath::Pi()){
696 offset = -2*TMath::Pi();
697 }
698 return fNegativeTrackLorentzVector->Phi()+offset;
699}
700
701Double_t AliV0Reader::GetPositiveTrackPhi() const{
702 //see header file for documentation
703
704 Double_t offset=0;
705 if(fPositiveTrackLorentzVector->Phi()> TMath::Pi()){
706 offset = -2*TMath::Pi();
707 }
708 return fPositiveTrackLorentzVector->Phi()+offset;
709}
710
711Double_t AliV0Reader::GetMotherCandidatePhi() const{
712 //see header file for documentation
713
714 Double_t offset=0;
715 if(fMotherCandidateLorentzVector->Phi()> TMath::Pi()){
716 offset = -2*TMath::Pi();
717 }
718 return fMotherCandidateLorentzVector->Phi()+offset;
719}
720
721
722Double_t AliV0Reader::GetMotherCandidateRapidity() const{
723 //see header file for documentation
724
725 Double_t rapidity=0;
726 if(fMotherCandidateLorentzVector->Energy() - fMotherCandidateLorentzVector->Pz() == 0 || fMotherCandidateLorentzVector->Energy() + fMotherCandidateLorentzVector->Pz() == 0) rapidity=0;
727 else rapidity = 0.5*(TMath::Log((fMotherCandidateLorentzVector->Energy() + fMotherCandidateLorentzVector->Pz()) / (fMotherCandidateLorentzVector->Energy()-fMotherCandidateLorentzVector->Pz())));
728 return rapidity;
729
730}
731
732
733
734
735
736Int_t AliV0Reader::GetSpeciesIndex(Int_t chargeOfTrack){
737 //see header file for documentation
738
739 Int_t iResult = 10; // Unknown particle
740
741 if(chargeOfTrack==-1){ //negative track
742 switch(abs(fNegativeTrackPID)){
743 case 11: //electron
744 iResult = 0;
745 break;
746 case 13: //muon
747 iResult = 1;
748 break;
749 case 211: //pion
750 iResult = 2;
751 break;
752 case 321: //kaon
753 iResult = 3;
754 break;
755 case 2212: //proton
756 iResult = 4;
757 break;
758 case 22: //photon
759 iResult = 5;
760 break;
761 case 111: //pi0
762 iResult = 6;
763 break;
764 case 2112: //neutron
765 iResult = 7;
766 break;
767 case 311: //K0
768 iResult = 8;
769 break;
770
771 //Put in here for kSPECIES::kEleCon ????
772 }
773 }
774 else if(chargeOfTrack==1){ //positive track
775 switch(abs(fPositiveTrackPID)){
776 case 11: //electron
777 iResult = 0;
778 break;
779 case 13: //muon
780 iResult = 1;
781 break;
782 case 211: //pion
783 iResult = 2;
784 break;
785 case 321: //kaon
786 iResult = 3;
787 break;
788 case 2212: //proton
789 iResult = 4;
790 break;
791 case 22: //photon
792 iResult = 5;
793 break;
794 case 111: //pi0
795 iResult = 6;
796 break;
797 case 2112: //neutron
798 iResult = 7;
799 break;
800 case 311: //K0
801 iResult = 8;
802 break;
803
804 //Put in here for kSPECIES::kEleCon ????
805 }
806 }
807 else{
808 //Wrong parameter.. Print warning
809 }
810 return iResult;
811}
812
4a6157dc 813Bool_t AliV0Reader::GetHelixCenter(AliESDtrack* track, Double_t b,Int_t charge, Double_t center[2]){
a0b94e5c 814 // see header file for documentation
815
816 Double_t pi = 3.14159265358979323846;
817
818 Double_t helix[6];
819 track->GetHelixParameters(helix,b);
820
821 Double_t xpos = helix[5];
822 Double_t ypos = helix[0];
823 Double_t radius = TMath::Abs(1./helix[4]);
824 Double_t phi = helix[2];
825
826 if(phi < 0){
827 phi = phi + 2*pi;
828 }
829
830 phi -= pi/2.;
831 Double_t xpoint = radius * TMath::Cos(phi);
832 Double_t ypoint = radius * TMath::Sin(phi);
833
834 if(charge > 0){
835 xpoint = - xpoint;
836 ypoint = - ypoint;
837 }
838
839 if(charge < 0){
840 xpoint = xpoint;
841 ypoint = ypoint;
842 }
843 center[0] = xpos + xpoint;
844 center[1] = ypos + ypoint;
845
846 return 1;
847}
848
4a6157dc 849Bool_t AliV0Reader::GetConvPosXY(AliESDtrack* ptrack, AliESDtrack* ntrack, Double_t b, Double_t convpos[2]){
a0b94e5c 850 //see header file for documentation
851
852 Double_t helixcenterpos[2];
853 GetHelixCenter(ptrack,b,ptrack->Charge(),helixcenterpos);
854
855 Double_t helixcenterneg[2];
856 GetHelixCenter(ntrack,b,ntrack->Charge(),helixcenterneg);
857
858 Double_t poshelix[6];
859 ptrack->GetHelixParameters(poshelix,b);
860 Double_t posradius = TMath::Abs(1./poshelix[4]);
861
862 Double_t neghelix[6];
863 ntrack->GetHelixParameters(neghelix,b);
864 Double_t negradius = TMath::Abs(1./neghelix[4]);
865
866 Double_t xpos = helixcenterpos[0];
867 Double_t ypos = helixcenterpos[1];
868 Double_t xneg = helixcenterneg[0];
869 Double_t yneg = helixcenterneg[1];
870
871 convpos[0] = (xpos*negradius + xneg*posradius)/(negradius+posradius);
872 convpos[1] = (ypos*negradius+ yneg*posradius)/(negradius+posradius);
873
874 return 1;
875}
876
877
878
4a6157dc 879Double_t AliV0Reader::GetConvPosZ(AliESDtrack* ptrack,AliESDtrack* ntrack, Double_t b){
a0b94e5c 880 //see header file for documentation
881
882 Double_t helixpos[6];
883 ptrack->GetHelixParameters(helixpos,b);
884
885 Double_t helixneg[6];
886 ntrack->GetHelixParameters(helixneg,b);
887
888 Double_t negtrackradius = TMath::Abs(1./helixneg[4]);
889 Double_t postrackradius = TMath::Abs(1./helixpos[4]);
890
891 Double_t pi = 3.14159265358979323846;
892
893 Double_t convpos[2];
894 GetConvPosXY(ptrack,ntrack,b,convpos);
895
896 Double_t convposx = convpos[0];
897 Double_t convposy = convpos[1];
898
899 Double_t helixcenterpos[2];
900 GetHelixCenter(ptrack,b,ptrack->Charge(),helixcenterpos);
901
902 Double_t helixcenterneg[2];
903 GetHelixCenter(ntrack,b,ntrack->Charge(),helixcenterneg);
904
905 Double_t xpos = helixcenterpos[0];
906 Double_t ypos = helixcenterpos[1];
907 Double_t xneg = helixcenterneg[0];
908 Double_t yneg = helixcenterneg[1];
909
4a6157dc 910 Double_t deltaXPos = convposx - xpos;
911 Double_t deltaYPos = convposy - ypos;
a0b94e5c 912
4a6157dc 913 Double_t deltaXNeg = convposx - xneg;
914 Double_t deltaYNeg = convposy - yneg;
a0b94e5c 915
4a6157dc 916 Double_t alphaPos = pi + TMath::ATan2(-deltaYPos,-deltaXPos);
917 Double_t alphaNeg = pi + TMath::ATan2(-deltaYNeg,-deltaXNeg);
a0b94e5c 918
4a6157dc 919 Double_t vertexXNeg = xneg + TMath::Abs(negtrackradius)*
920 TMath::Cos(alphaNeg);
921 Double_t vertexYNeg = yneg + TMath::Abs(negtrackradius)*
922 TMath::Sin(alphaNeg);
a0b94e5c 923
4a6157dc 924 Double_t vertexXPos = xpos + TMath::Abs(postrackradius)*
925 TMath::Cos(alphaPos);
926 Double_t vertexYPos = ypos + TMath::Abs(postrackradius)*
927 TMath::Sin(alphaPos);
a0b94e5c 928
929 Double_t x0neg = helixneg[5];
930 Double_t y0neg = helixneg[0];
931
932 Double_t x0pos = helixpos[5];
933 Double_t y0pos = helixpos[0];
934
4a6157dc 935 Double_t dNeg = TMath::Sqrt((vertexXNeg - x0neg)*(vertexXNeg - x0neg)
936 +(vertexYNeg - y0neg)*(vertexYNeg - y0neg));
a0b94e5c 937
4a6157dc 938 Double_t dPos = TMath::Sqrt((vertexXPos - x0pos)*(vertexXPos - x0pos)
939 +(vertexYPos - y0pos)*(vertexYPos - y0pos));
a0b94e5c 940
4a6157dc 941 Double_t rNeg = TMath::Sqrt(negtrackradius*negtrackradius -
942 dNeg*dNeg/4.);
a0b94e5c 943
4a6157dc 944 Double_t rPos = TMath::Sqrt(postrackradius*postrackradius -
945 dPos*dPos/4.);
a0b94e5c 946
4a6157dc 947 Double_t deltabetaNeg = 2*(pi + TMath::ATan2(-dNeg/2.,-rNeg));
948 Double_t deltabetaPos = 2*(pi + TMath::ATan2(-dPos/2.,-rPos));
a0b94e5c 949
4a6157dc 950 Double_t deltaUNeg = negtrackradius*deltabetaNeg;
951 Double_t deltaUPos = postrackradius*deltabetaPos;
a0b94e5c 952
4a6157dc 953 Double_t zphaseNeg = ntrack->GetZ() + deltaUNeg * ntrack->GetTgl();
954 Double_t zphasePos = ptrack->GetZ() + deltaUPos * ptrack->GetTgl();
a0b94e5c 955
956 Double_t convposz =
4a6157dc 957 (zphasePos*negtrackradius+zphaseNeg*postrackradius)/(negtrackradius+postrackradius);
a0b94e5c 958
959 return convposz;
960}