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