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