3 /* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * See cxx source for full Copyright notice */
6 ////////////////////////////////////////////////
7 //---------------------------------------------
8 // Class used to do analysis on conversion pairs
9 //---------------------------------------------
10 ////////////////////////////////////////////////
12 // --- ROOT system ---
14 #include "AliMCEvent.h" // for CF
16 #include "AliESDEvent.h"
17 #include "AliKFParticle.h"
18 #include "TParticle.h"
19 #include "AliGammaConversionHistograms.h"
21 #include "AliCFManager.h"
22 #include "AliGammaConversionBGHandler.h"
23 #include "AliESDpid.h"
32 //--- AliRoot system ---
35 class AliMCEvent; // for CF
37 class AliMCEventHandler;
38 class AliESDInputHandler;
43 class AliCFManager; // for CF
44 class AliCFContainer; // for CF
45 //class AliESDpid; // for dEdx cut based on nSigma to particle lines
46 class AliESDtrackCuts;
50 class AliV0Reader : public TObject {
58 kStepReconstructable = 1,
63 kStepdEdxElectronselection = 6,
64 kStepdEdxPionrejection = 7,
65 kStepNContributors = 8,
79 AliV0Reader(); //constructor
80 AliV0Reader(const AliV0Reader & g); //copy constructor
81 AliV0Reader & operator = (const AliV0Reader & g); //assignment operator
82 // virtual ~AliV0Reader() {;} //virtual destructor
83 virtual ~AliV0Reader(); //virtual destructor
85 *Initialize the reader
88 void SetInputAndMCEvent(AliVEvent* esd, AliMCEvent* mc) ;
91 virtual void SetInputEvent(AliVEvent* const input) {fESDEvent = dynamic_cast<AliESDEvent*>(input);}
92 virtual void SetMC(AliMCEvent* const mc) {fMCEvent = mc;}
96 void SetCFManager(AliCFManager * const io){fCFManager = io;};
97 AliCFManager *GetCFManager() const {return fCFManager;}
103 * Returns AliESDEvent
105 AliESDEvent* GetESDEvent() const{return fESDEvent;}
108 *Returns the number of v0s in the event, no cuts applied.
110 Int_t GetNumberOfV0s() const{return fESDEvent->GetNumberOfV0s();}
113 *Returns the number of contributors to the vertex
115 // Int_t GetNumberOfContributorsVtx() const{return fESDEvent->GetPrimaryVertex()->GetNContributors();}
116 Int_t GetNumberOfContributorsVtx();
119 * Check if there are any more good v0s left in the v0 stack
120 * if so, fCurrent v0 is set to this v0 and can be retrieved
121 * by GetCurrentV0 function.
122 * returns kFALSE if there is no more good v0s in the v0 stack
127 * Returns the v0 at the given index, no checks are done on the v0.
129 AliESDv0* GetV0(Int_t index);
132 * Returns the current v0
134 AliESDv0* GetCurrentV0() const{return fCurrentV0;}
137 * Returns the negative ESD track which belongs to fCurrentV0
139 // AliESDtrack* GetNegativeESDTrack(){return fESDEvent->GetTrack(fCurrentV0->GetNindex());}
140 AliESDtrack* GetNegativeESDTrack(){return fCurrentNegativeESDTrack;}
143 * Returns the positive ESD track which belongs to fCurrentV0
145 // AliESDtrack* GetPositiveESDTrack(){return fESDEvent->GetTrack(fCurrentV0->GetPindex());}
146 AliESDtrack* GetPositiveESDTrack(){return fCurrentPositiveESDTrack;}
149 * Returns the negative KF particle which belongs to fCurrentV0
151 AliKFParticle* GetNegativeKFParticle() const{return fCurrentNegativeKFParticle;}
154 * Returns the positive KF particle which belongs to fCurrentV0
156 AliKFParticle* GetPositiveKFParticle() const{return fCurrentPositiveKFParticle;}
159 * Returns the KFParticle object of the 2 tracks.
161 AliKFParticle* GetMotherCandidateKFCombination() const{return fCurrentMotherKFCandidate;}
164 * Checks the probablity that the PID of the particle is what we want it to be.
166 Bool_t CheckPIDProbability(Double_t negProbCut, Double_t posProbCut);
169 * Checks if the PID of the two particles are within our cuts.
171 void GetPIDProbability(Double_t &negPIDProb, Double_t &posPIDProb);
174 * Checks if the PID of the two particles are within our cuts.
176 void GetPIDProbabilityMuonPion(Double_t &negPIDProb, Double_t &posPIDProb);
179 *Get the negative MC TParticle from the stack
181 TParticle * GetNegativeMCParticle() const{return fNegativeMCParticle;}
184 *Get the positive MC TParticle from the stack
186 TParticle * GetPositiveMCParticle() const{return fPositiveMCParticle;}
189 *Get the mother MC TParticle from the stack
191 TParticle * GetMotherMCParticle() const{return fMotherMCParticle;}
194 * Flag to see if the v0 particles share the same mother
196 Bool_t HasSameMCMother();
200 *Get the PID of the MC mother particle
202 Int_t GetMotherMCParticlePDGCode() const{return fMotherMCParticle->GetPdgCode();}
207 AliStack* GetMCStack() const{return fMCStack;}
211 * Setup AliMCEventHandler
213 // AliMCEventHandler* GetMCTruth() const{return fMCTruth;} // for CF
219 AliMCEvent* GetMCEvent() const{return fMCEvent;} // for CF
223 *Get the magnetic field from the ESD event
225 Double_t GetMagneticField() const{return fESDEvent->GetMagneticField();}
228 *Get the primary vertex from the esd event
230 const AliESDVertex *GetPrimaryVertex() const {return fESDEvent->GetPrimaryVertex();}
233 * Set the PID of the negative track
235 void SetNegativeTrackPID(Int_t negTrackPID){fNegativeTrackPID=negTrackPID;}
238 * Set the PID of the positive track
240 void SetPositiveTrackPID(Int_t posTrackPID){fPositiveTrackPID=posTrackPID;}
243 * Set the flag to use the kfparticle class. Will also disable the use of esd tracks
245 void UseKFParticle(){fUseKFParticle = kTRUE; fUseESDTrack = kFALSE;}
248 * Set the flag to use the esd track class. Will also disable the use of kf particles
250 void UseESDTrack(){fUseESDTrack = kTRUE; fUseKFParticle = kFALSE;}
253 * Set the flag to use improved vertex or not
255 void SetUseImprovedVertex(Bool_t useImprovedVertex){fUseImprovedVertex=useImprovedVertex;}
258 * Return the number in the species array belonging to the negative or positive track pid.
260 Int_t GetSpeciesIndex(Int_t chargeOfTrack);
263 * Return the x coordinate of the v0
265 Double_t GetX() const{return fCurrentXValue;}
268 * Return the y coordinate of the v0
270 Double_t GetY() const{return fCurrentYValue;}
273 * Return the Z coordinate of the v0
275 Double_t GetZ() const{return fCurrentZValue;}
278 * Return the radius of the v0
280 Double_t GetXYRadius() const{return sqrt((Double_t)(fCurrentXValue*fCurrentXValue + fCurrentYValue*fCurrentYValue));}
283 * Get the opening angle between the two tracks
285 Double_t GetOpeningAngle(){return fNegativeTrackLorentzVector->Angle(fPositiveTrackLorentzVector->Vect());}
288 * Get the Cos Pointing angle between the two tracks
290 Double_t GetCosPointingAngle(){return fCurrentV0->GetV0CosineOfPointingAngle();}
293 * Get the DCA between the two tracks
295 Double_t GetDcaDaughters(){return fCurrentV0->GetDcaV0Daughters();}
298 * Get the Normalized DCA between the two tracks
300 Double_t GetNormDcaDistDaughters(){return fCurrentV0->GetDcaV0Daughters()/fCurrentV0->GetDistSigma();}
303 * Get the Likelihood for a Conversion
305 Double_t GetLikelihoodAP(){return fCurrentV0->GetLikelihoodAP(0,0);}
308 * Gets the Energy of the negative track.
310 Double_t GetNegativeTrackEnergy() const{return fCurrentNegativeKFParticle->E();}
313 * Gets the Energy of the positive track.
315 Double_t GetPositiveTrackEnergy() const{return fCurrentPositiveKFParticle->E();}
318 * Gets the Energy of the mother candidate.
320 Double_t GetMotherCandidateEnergy() const{return fCurrentMotherKFCandidate->E();}
323 * Gets the Pt of the negative track.
325 Double_t GetNegativeTrackPt() const{return fNegativeTrackLorentzVector->Pt();}
328 * Gets the Pt of the positive track.
330 Double_t GetPositiveTrackPt() const{return fPositiveTrackLorentzVector->Pt();}
334 * Gets the Pt of the mother candidate.
336 Double_t GetMotherCandidatePt() const{return fMotherCandidateLorentzVector->Pt();}
340 * Gets the P of the mother candidate.
342 Double_t GetMotherCandidateP() const{return fMotherCandidateLorentzVector->P();}
346 * Gets the Eta of the negative track.
348 Double_t GetNegativeTrackEta() const{return fNegativeTrackLorentzVector->Eta();}
350 * Gets the Eta of the positive track.
352 Double_t GetPositiveTrackEta() const{return fPositiveTrackLorentzVector->Eta();}
354 * Gets the Eta of the mother candidate.
356 Double_t GetMotherCandidateEta() const{return fMotherCandidateLorentzVector->Eta();}
359 * Gets the NDF of the mother candidate.
361 Double_t GetMotherCandidateNDF() const{return fCurrentMotherKFCandidate->GetNDF();}
364 * Gets the Chi2 of the mother candidate.
366 Double_t GetMotherCandidateChi2() const{return fCurrentMotherKFCandidate->GetChi2();}
369 * Gets the Mass of the mother candidate.
371 Double_t GetMotherCandidateMass() const{return fMotherCandidateKFMass;}
374 * Gets the Width of the mother candidate.
376 Double_t GetMotherCandidateWidth() const{return fMotherCandidateKFWidth;}
379 * Gets the Phi of the negative track.
381 Double_t GetNegativeTrackPhi() const;
384 * Gets the Phi of the positive track.
386 Double_t GetPositiveTrackPhi() const;
389 * Gets the Phi of the mother candidate.
391 Double_t GetMotherCandidatePhi() const;
394 * Gets the Rapidity of the mother candidate.
396 Double_t GetMotherCandidateRapidity() const;
400 * Gets the P of the negative track.
402 Double_t GetNegativeTrackP() const{return fNegativeTrackLorentzVector->P();}
405 * Gets the P of the positive track.
407 Double_t GetPositiveTrackP() const{return fPositiveTrackLorentzVector->P();}
410 * Gets the dE/dx in the TPC of the negative track.
412 Double_t GetNegativeTrackTPCdEdx() const{return fCurrentNegativeESDTrack->GetTPCsignal();}
415 * Gets the dE/dx in the TPC of the positive track.
417 Double_t GetPositiveTrackTPCdEdx() const{return fCurrentPositiveESDTrack->GetTPCsignal();}
420 * Gets the Number of the TPC clusters of the negative track.
422 Int_t GetNegativeTracknTPCClusters() const{return fCurrentNegativeESDTrack->GetNcls(1);}
425 * Gets the Number of the TPC clusters of the positive track.
427 Int_t GetPositiveTracknTPCClusters() const{return fCurrentPositiveESDTrack->GetNcls(1);}
430 * Get the TOFsignal for negative/positive track. RRnewTOF
432 Double_t GetNegativeTrackTOFsignal() const{return fCurrentNegativeESDTrack->GetTOFsignal();}
433 Double_t GetPositiveTrackTOFsignal() const{return fCurrentPositiveESDTrack->GetTOFsignal();}
436 * Gets the Number of the TPC findable clusters of the negative track.
438 Int_t GetNegativeTracknTPCFClusters() const{return fCurrentNegativeESDTrack->GetTPCNclsF();}
441 * Gets the Number of the TPC findable clusters of the positive track.
443 Int_t GetPositiveTracknTPCFClusters() const{return fCurrentPositiveESDTrack->GetTPCNclsF();}
446 * Gets the Number of the ITS clusters of the negative track.
448 Int_t GetNegativeTracknITSClusters() const{return fCurrentNegativeESDTrack->GetNcls(0);}
451 * Gets the Number of the ITS clusters of the positive track.
453 Int_t GetPositiveTracknITSClusters() const{return fCurrentPositiveESDTrack->GetNcls(0);}
456 * Gets the chi2 of the TPC negative track.
458 Double_t GetNegativeTrackTPCchi2() const{return fCurrentNegativeESDTrack->GetTPCchi2();}
461 * Gets the chi2 of the TPC the positive track.
463 Double_t GetPositiveTrackTPCchi2() const{return fCurrentPositiveESDTrack->GetTPCchi2();}
466 * Update data which need to be updated every event.
468 void UpdateEventByEventData();
471 * Gets the MaxRCut value.
473 Double_t GetMaxVertexZ() const{return fMaxVertexZ;}
476 * Gets the MaxRCut value.
478 Double_t GetMaxRCut() const{return fMaxR;}
481 * Gets the MinRCut value.
483 Double_t GetMinRCut() const{return fMinR;}
486 * Gets the Eta cut value.
488 Double_t GetEtaCut() const{return fEtaCut;}
491 * Gets the Rapidity Meson cut value.
493 Double_t GetRapidityMesonCut() const{return fRapidityMesonCut;}
496 * Gets the Pt cut value.
498 Double_t GetPtCut() const{return fPtCut;}
502 * Gets the MaxZCut value.
504 Double_t GetMaxZCut() const{return fMaxZ;}
508 * Gets the MinClsTPC value.
510 Double_t GetMinClsTPCCut() const{return fMinClsTPC;}
513 * Gets the MinClsTPC value.
515 Double_t GetMinClsTPCCutToF() const{return fMinClsTPCToF;}
520 * Gets the line cut values.
522 Double_t GetLineCutZRSlope() const{return fLineCutZRSlope;}
523 Double_t GetLineCutZValue() const{return fLineCutZValue;}
526 * Gets the Chi2 cut value for the conversions.
528 Double_t GetChi2CutConversion() const{return fChi2CutConversion;}
531 * Gets the Chi2 cut value for the mesons.
533 Double_t GetChi2CutMeson() const{return fChi2CutMeson;}
536 * Gets the alpha cut value for the mesons.
538 Double_t GetAlphaCutMeson() const{return fAlphaCutMeson;}
541 * Gets the Minimum alpha cut value for the mesons.
543 Double_t GetAlphaMinCutMeson() const{return fAlphaMinCutMeson;}
545 Double_t GetPositiveTrackLength() const{return fCurrentPositiveESDTrack->GetIntegratedLength();}
546 Double_t GetNegativeTrackLength() const{return fCurrentNegativeESDTrack->GetIntegratedLength();}
548 Double_t GetPositiveNTPCClusters() const{return fCurrentPositiveESDTrack->GetTPCNcls();}
549 Double_t GetNegativeNTPCClusters() const{return fCurrentNegativeESDTrack->GetTPCNcls();}
552 * Sets the MaxVertexZ value.
554 void SetMaxVertexZ(Double_t maxVertexZ){fMaxVertexZ=maxVertexZ;}
557 * Sets the MaxRCut value.
559 void SetMaxRCut(Double_t maxR){fMaxR=maxR;}
561 * Sets the MinRCut value.
563 void SetMinRCut(Double_t minR){fMinR=minR;}
566 * Sets the EtaCut value.
568 void SetEtaCut(Double_t etaCut){fEtaCut=etaCut;}
571 * Sets the Rapidity Meson Cut value.
573 void SetRapidityMesonCut(Double_t RapidityMesonCut){fRapidityMesonCut=RapidityMesonCut;}
576 * Sets the PtCut value.
578 void SetPtCut(Double_t ptCut){fPtCut=ptCut;}
581 * Sets the PtCut value.
583 void SetSinglePtCut(Double_t singleptCut){fSinglePtCut=singleptCut;}
587 * Sets the MaxZCut value.
589 void SetMaxZCut(Double_t maxZ){fMaxZ=maxZ;}
592 * Sets the MinClsTPC value.
594 void SetMinClsTPCCut(Double_t minClsTPC){fMinClsTPC=minClsTPC;}
597 * Sets the MinClsTPC value.
599 void SetMinClsTPCCutToF(Double_t minClsTPCToF){fMinClsTPCToF=minClsTPCToF;}
603 * Sets the LineCut values.
605 void SetLineCutZRSlope(Double_t LineCutZRSlope){fLineCutZRSlope=LineCutZRSlope;}
606 void SetLineCutZValue(Double_t LineCutZValue){fLineCutZValue=LineCutZValue;}
609 * Sets the Chi2Cut value for conversions.
611 void SetChi2CutConversion(Double_t chi2){fChi2CutConversion=chi2;}
614 * Sets the Chi2Cut for the mesons.
616 void SetChi2CutMeson(Double_t chi2){fChi2CutMeson=chi2;}
619 * Sets the AlphaCut for the mesons.
621 void SetAlphaCutMeson(Double_t alpha){fAlphaCutMeson=alpha;}
625 * Sets the AlphaCut for the mesons.
627 void SetAlphaMinCutMeson(Double_t alpha){fAlphaMinCutMeson=alpha;}
631 * Sets the XVertexCut value.
633 void SetXVertexCut(Double_t xVtx){fCurrentXValue=xVtx;}
636 * Sets the YVertexCut value.
638 void SetYVertexCut(Double_t yVtx){fCurrentYValue=yVtx;}
641 * Sets the ZVertexCut value.
643 void SetZVertexCut(Double_t zVtx){fCurrentZValue=zVtx;}
646 * Sets the PIDProbabilityCut value for track particles.
648 void SetPIDProbability(Double_t pidProb){fPIDProbabilityCutPositiveParticle=pidProb; fPIDProbabilityCutNegativeParticle=pidProb;}
651 * Sets the PIDProbability cut value for the negative track.
653 void SetPIDProbabilityNegativeParticle(Double_t pidProb){fPIDProbabilityCutNegativeParticle=pidProb;}
656 * Sets the PIDProbability cut value for the positive track.
658 void SetPIDProbabilityPositiveParticle(Double_t pidProb){fPIDProbabilityCutPositiveParticle=pidProb;}
661 * Sets the PIDnSigmaAboveElectron cut value for the tracks.
663 void SetPIDnSigmaAboveElectronLine(Double_t nSigmaAbove){fPIDnSigmaAboveElectronLine=nSigmaAbove;}
664 void SetTofPIDnSigmaAboveElectronLine(Double_t nTofSigmaAbove){fTofPIDnSigmaAboveElectronLine=nTofSigmaAbove;} // RRnewTOF
667 * Sets the PIDnSigmaBelowElectron cut value for the tracks.
669 void SetPIDnSigmaBelowElectronLine(Double_t nSigmaBelow){fPIDnSigmaBelowElectronLine=nSigmaBelow;}
670 void SetTofPIDnSigmaBelowElectronLine(Double_t nTofSigmaBelow){fTofPIDnSigmaBelowElectronLine=nTofSigmaBelow;} // RRnewTOF
673 * Sets the PIDnSigmaAbovePion cut value for the tracks.
675 void SetPIDnSigmaAbovePionLine(Double_t nSigmaAbovePion){fPIDnSigmaAbovePionLine=nSigmaAbovePion;}
678 * Sets the PIDnSigmaAbovePion cut value for the tracks.
680 void SetPIDnSigmaAbovePionLineHighPt(Double_t nSigmaAbovePionHighPt){fPIDnSigmaAbovePionLineHighPt=nSigmaAbovePionHighPt;}
683 * Sets the PIDMinPnSigmaAbovePion cut value for the tracks.
685 void SetPIDMinPnSigmaAbovePionLine(Double_t MinPnSigmaAbovePion){fPIDMinPnSigmaAbovePionLine=MinPnSigmaAbovePion;}
688 * Sets the PIDMinPnSigmaAbovePion cut value for the tracks.
690 void SetPIDMaxPnSigmaAbovePionLine(Double_t MaxPnSigmaAbovePion){fPIDMaxPnSigmaAbovePionLine=MaxPnSigmaAbovePion;}
693 * Sets the SigmaMassCut value.
695 void SetSigmaMass(Double_t sigmaMass){fNSigmaMass=sigmaMass;}
698 * Sets the flag to enable/disable the usage of MC information.
700 void SetDoMCTruth(Bool_t doMC){fDoMC = doMC;}
703 * Sets the flag to enable/disable the usage of MC information.
705 Bool_t GetDoMCTruth() const {return fDoMC;}
708 * Sets the flag to enable/disable the cut dedx N sigma
711 void SetDodEdxSigmaCut( Bool_t dodEdxSigmaCut){fDodEdxSigmaCut=dodEdxSigmaCut;}
712 void SetDoTOFsigmaCut( Bool_t doTOFsigmaCut){fDoTOFsigmaCut=doTOFsigmaCut;} //RRnewTOF
713 void SetDoPhotonAsymmetryCut( Bool_t doPhotonAsymmetryCut){fDoPhotonAsymmetryCut=doPhotonAsymmetryCut;}
715 void SetMinPPhotonAsymmetryCut(Double_t minPPhotonAsymmetryCut){fMinPPhotonAsymmetryCut=minPPhotonAsymmetryCut;}
716 void SetMinPhotonAsymmetry(Double_t minPhotonAsymmetry){fMinPhotonAsymmetry=minPhotonAsymmetry;}
718 * Sets the flag to enable/disable the cut dedx N sigma for Kaon Rejection at low p
720 void SetDoKaonRejectionLowP( Bool_t doKaonRejectionLowP){fDoKaonRejectionLowP=doKaonRejectionLowP;}
722 * Sets the flag to enable/disable the cut dedx N sigma for Proton Rejection at low p
724 void SetDoProtonRejectionLowP( Bool_t doProtonRejectionLowP){fDoProtonRejectionLowP=doProtonRejectionLowP;}
727 * Sets the flag to enable/disable the cut dedx N sigma for Pion Rejection at low p
729 void SetDoPionRejectionLowP( Bool_t doPionRejectionLowP){fDoPionRejectionLowP=doPionRejectionLowP;}
732 * Sets the PIDMinPnSigmaAroundKaon cut value for the tracks.
734 void SetPIDnSigmaAtLowPAroundKaonLine(Double_t nSigmaAtLowPAroundKaon){fPIDnSigmaAtLowPAroundKaonLine =nSigmaAtLowPAroundKaon;}
737 * Sets the PIDMinPnSigmaAroundProton cut value for the tracks.
739 void SetPIDnSigmaAtLowPAroundProtonLine(Double_t nSigmaAtLowPAroundProton){fPIDnSigmaAtLowPAroundProtonLine =nSigmaAtLowPAroundProton;}
742 * Sets the PIDMinPnSigmaAroundPion cut value for the tracks.
744 void SetPIDnSigmaAtLowPAroundPionLine(Double_t nSigmaAtLowPAroundPion){fPIDnSigmaAtLowPAroundPionLine =nSigmaAtLowPAroundPion;}
747 * Sets the PIDMinPnSigmaAbovePion cut value for the tracks.
749 void SetPIDMinPKaonRejectionLowP(Double_t PIDMinPKaonRejectionLowP ){fPIDMinPKaonRejectionLowP=PIDMinPKaonRejectionLowP;}
752 * Sets the PIDMinPnSigmaAbovePion cut value for the tracks.
754 void SetPIDMinPProtonRejectionLowP(Double_t PIDMinPProtonRejectionLowP ){fPIDMinPProtonRejectionLowP=PIDMinPProtonRejectionLowP;}
756 * Sets the PIDMinPnSigmaAbovePion cut value for the tracks.
758 void SetPIDMinPPionRejectionLowP(Double_t PIDMinPPionRejectionLowP ){fPIDMinPPionRejectionLowP=PIDMinPPionRejectionLowP;}
761 *Set if we want to use Gamma Selection based on Qt from Armenteros
763 void SetDoQtGammaSelection(Bool_t doQtGammaSelection){fDoQtGammaSelection=doQtGammaSelection;}
764 void SetDoHighPtQtGammaSelection(Bool_t doHighPtQtGammaSelection){fDoHighPtQtGammaSelection=doHighPtQtGammaSelection;} // RRnew
766 * Sets the MaxQtCut value.
768 void SetQtMax(Double_t qtMax){fQtMax=qtMax;}
769 void SetHighPtQtMax(Double_t qtMaxHighPt){fHighPtQtMax=qtMaxHighPt;} // RRnew
770 void SetPtBorderForQt(Double_t ptBorderForQt){fPtBorderForQt=ptBorderForQt;} // RRnew
773 * Updates the V0 information of the current V0.
775 Bool_t UpdateV0Information();
778 * Resets the V0 index.
780 void ResetV0IndexNumber(){fCurrentV0IndexNumber=0;}
784 * Returns number of good v0s in the event
786 Int_t GetNGoodV0s() const {return fNumberOfGoodV0s;}
789 * Sets the histograms.
791 void SetHistograms(AliGammaConversionHistograms * const histograms){fHistograms=histograms;}
794 * Check for primary vertex.
796 Bool_t CheckForPrimaryVertex();
799 * Check for primary vertex Z.
801 Bool_t CheckForPrimaryVertexZ();
804 * Gets a vector of good v0s.
806 TClonesArray* GetCurrentEventGoodV0s() const{return fCurrentEventGoodV0s;}
809 * Gets the vector of previous events v0s (for bacground analysis)
811 AliGammaConversionKFVector* GetBGGoodV0s(Int_t event) const;
812 // vector<AliKFParticle> GetPreviousEventGoodV0s() const{return fPreviousEventGoodV0s;}
814 void SetUseOwnXYZCalculation(Bool_t flag){fUseOwnXYZCalculation=flag;}
816 void SetUseConstructGamma(Bool_t flag){fUseConstructGamma=flag;}
818 Bool_t GetHelixCenter(AliESDtrack* track, Double_t b,Int_t charge, Double_t center[2]);
820 Bool_t GetConvPosXY(AliESDtrack* ptrack,AliESDtrack* ntrack, Double_t b, Double_t convpos[2]);
822 Double_t GetConvPosZ(AliESDtrack* ptrack,AliESDtrack* ntrack, Double_t b);
824 Bool_t GetArmenterosQtAlfa(const AliKFParticle * posKFparticle, const AliKFParticle * negKFparticle, const AliKFParticle * gamKFparticle,Double_t armenterosQtAlfa[2]);
826 void SetDoCF(Bool_t flag){fDoCF = flag;}
828 Bool_t CheckV0FinderStatus(Int_t index);
830 void SetOnFlyFlag(Bool_t flag){fUseOnFlyV0Finder = flag;}
832 Int_t GetNBGEvents(){return fBGEventHandler->GetNBGEvents();}
834 void SetCalculateBackground(Bool_t flag){fCalculateBackground=flag;}
836 AliGammaConversionBGHandler* GetBGHandler() const {return fBGEventHandler;}
838 Double_t GetVertexZ(){return fESDEvent->GetPrimaryVertex()->GetZ();}
840 Int_t GetMultiplicity(){return CountESDTracks();}
842 void SetESDtrackCuts(AliESDtrackCuts * const trackCuts){fEsdTrackCuts = trackCuts;}
844 void SetNEventsForBG(Int_t nev){fNEventsForBGCalculation = nev;}
846 Int_t CountESDTracks();
848 Int_t GetCurrentV0IndexNumber() const {return fCurrentV0IndexNumber;}
850 Bool_t CheckIfPi0IsMother(Int_t label);
851 Bool_t CheckIfEtaIsMother(Int_t label);
853 static void InitESDpid(Int_t type=0);
854 static void SetESDpid(AliESDpid * const pid) {fgESDpid=pid;}
855 static AliESDpid* GetESDpid() {return fgESDpid;}
857 void SetUseChargedTracksMultiplicityForBG(Bool_t flag){fUseChargedTrackMultiplicityForBG = flag;}
858 void SetIsHeavyIon(Int_t isHeavyIon) {fIsHeavyIon=isHeavyIon;}
859 Int_t GetIsHeavyIon() const { return fIsHeavyIon;}
861 Int_t GetPindex(Int_t i) {return fV0Pindex.at(i);}
862 Int_t GetNindex(Int_t i) {return fV0Nindex.at(i);}
864 void ResetNGoodV0s(){fNumberOfGoodV0s=0;}
865 Int_t GetFirstTPCRow(Double_t radius);
867 void SetUseCorrectedTPCClsInfo(Bool_t flag){fUseCorrectedTPCClsInfo = flag;}
868 Bool_t GetUseCorrectedTPCClsInfo() const {return fUseCorrectedTPCClsInfo;}
870 void SetUseMCPSmearing(Int_t useMCPSmearing) {fUseMCPSmearing=useMCPSmearing;}
871 void SetPBremSmearing(Double_t pBremSmearing){fPBremSmearing=pBremSmearing;}
872 void SetPSigSmearing(Double_t pSigSmearing){fPSigSmearing=pSigSmearing;}
873 void SetPSigSmearingCte(Double_t pSigSmearingCte){fPSigSmearingCte=pSigSmearingCte;}
874 void SmearKFParticle(AliKFParticle * kfParticle);
877 AliStack * fMCStack; // pointer to MonteCarlo particle stack
878 // AliMCEventHandler* fMCTruth; // for CF pointer to the MC object
879 AliMCEvent *fMCEvent; // for CF pointer to MC event
880 TChain * fChain; // pointer to the TChain
882 // AliESDInputHandler* fESDHandler; //! pointer to esd object
883 AliESDEvent *fESDEvent; //! pointer to esd object
887 AliCFManager *fCFManager; // pointer to the cf manager
888 // AliCFContainer *container;
890 // for dEdx cut based on nSigma to a particle line
891 //AliESDpid * fESDpid; // esd pid
893 AliGammaConversionHistograms *fHistograms; // pointer to histogram handling class
895 Int_t fCurrentV0IndexNumber; // the current v0 index number
896 AliESDv0 * fCurrentV0; //! pointer to the current v0
897 AliKFParticle * fCurrentNegativeKFParticle; //! pointer to the negative KF particle
898 AliKFParticle * fCurrentPositiveKFParticle; //! pointer to the positive KF particle
899 AliKFParticle * fCurrentMotherKFCandidate; //! pointer to the positive KF particle
901 AliESDtrack * fCurrentNegativeESDTrack; //! pointer to the negative ESD track
902 AliESDtrack * fCurrentPositiveESDTrack; //! pointer to the positive ESD track
904 TLorentzVector * fNegativeTrackLorentzVector; //! pointer to the negative Track Lorentz Vector
905 TLorentzVector * fPositiveTrackLorentzVector; //! pointer to the positive Track Lorentz Vector
906 TLorentzVector * fMotherCandidateLorentzVector; //! pointer to the mother candidate Track Lorentz Vector
908 Double_t fCurrentXValue; // current x value
909 Double_t fCurrentYValue; // current y value
910 Double_t fCurrentZValue; // current z value
912 Int_t fPositiveTrackPID; // positive track pid
913 Int_t fNegativeTrackPID; // negative track pid
915 TParticle *fNegativeMCParticle; //!
916 TParticle *fPositiveMCParticle; //!
917 TParticle *fMotherMCParticle; //!
919 Double_t fMotherCandidateKFMass; // mass of mother candidate KF particle
920 Double_t fMotherCandidateKFWidth; // width of mother candidate KF particle
922 Bool_t fUseKFParticle; // flag
923 Bool_t fUseESDTrack; // flag
924 Bool_t fDoMC; // flag
927 Double_t fMaxVertexZ; // max z vertex cut
929 Double_t fMaxR; //r cut
930 Double_t fMinR; //r cut
931 Double_t fEtaCut; //eta cut
932 Double_t fRapidityMesonCut; //rapidity for meson cut
933 Double_t fPtCut; // pt cut
934 Double_t fSinglePtCut; // pt cut for electron/positron
935 Double_t fMaxZ; //z cut
936 Double_t fMinClsTPC; // minimum clusters in the TPC
937 Double_t fMinClsTPCToF; // minimum clusters to findable clusters
938 Double_t fLineCutZRSlope; //linecut
939 Double_t fLineCutZValue; //linecut
940 Double_t fChi2CutConversion; //chi2cut
941 Double_t fChi2CutMeson; //chi2cut
942 Double_t fAlphaCutMeson; //alphacut
943 Double_t fAlphaMinCutMeson; //alphacut
944 Double_t fPIDProbabilityCutNegativeParticle; //pid cut
945 Double_t fPIDProbabilityCutPositiveParticle; //pid cut
946 Bool_t fDodEdxSigmaCut; // flag to use the dEdxCut based on sigmas
947 Bool_t fDoTOFsigmaCut; // flag to use TOF pid cut RRnewTOF
948 Double_t fPIDnSigmaAboveElectronLine; // sigma cut
949 Double_t fPIDnSigmaBelowElectronLine; // sigma cut
950 Double_t fTofPIDnSigmaAboveElectronLine; // sigma cut RRnewTOF
951 Double_t fTofPIDnSigmaBelowElectronLine; // sigma cut RRnewTOF
952 Double_t fPIDnSigmaAbovePionLine; // sigma cut
953 Double_t fPIDnSigmaAbovePionLineHighPt; // sigma cut
954 Double_t fPIDMinPnSigmaAbovePionLine; // sigma cut
955 Double_t fPIDMaxPnSigmaAbovePionLine; // sigma cut
956 Double_t fDoKaonRejectionLowP; // Kaon rejection at low p
957 Double_t fDoProtonRejectionLowP; // Proton rejection at low p
958 Double_t fDoPionRejectionLowP; // Pion rejection at low p
959 Double_t fPIDnSigmaAtLowPAroundKaonLine; // sigma cut
960 Double_t fPIDnSigmaAtLowPAroundProtonLine; // sigma cut
961 Double_t fPIDnSigmaAtLowPAroundPionLine; // sigma cut
962 Double_t fPIDMinPKaonRejectionLowP; // Momentum limit to apply kaon rejection
963 Double_t fPIDMinPProtonRejectionLowP; // Momentum limit to apply proton rejection
964 Double_t fPIDMinPPionRejectionLowP; // Momentum limit to apply proton rejection
965 Bool_t fDoQtGammaSelection; // Select gammas using qtMax
966 Bool_t fDoHighPtQtGammaSelection; // RRnew Select gammas using qtMax for high pT
967 Double_t fQtMax; // Maximum Qt from Armenteros to select Gammas
968 Double_t fHighPtQtMax; // RRnew Maximum Qt for High pT from Armenteros to select Gammas
969 Double_t fPtBorderForQt; // RRnew
970 Double_t fXVertexCut; //vertex cut
971 Double_t fYVertexCut; //vertex cut
972 Double_t fZVertexCut; // vertexcut
974 Double_t fNSigmaMass; //nsigma cut
976 Bool_t fUseImprovedVertex; //flag
978 Bool_t fUseOwnXYZCalculation; //flag that determines if we use our own calculation of xyz (markus)
980 Bool_t fUseConstructGamma; //flag that determines if we use ConstructGamma method from AliKF
984 Bool_t fUseOnFlyV0Finder; //flag
986 Bool_t fUpdateV0AlreadyCalled; //flag
988 TClonesArray* fCurrentEventGoodV0s; //vector of good v0s
990 vector<Int_t> fV0Pindex; // index of positive track belonging to a V0
991 vector<Int_t> fV0Nindex; // index of positive track belonging to a V0
992 // vector<AliKFParticle> fPreviousEventGoodV0s; // vector of good v0s from prevous events
994 Bool_t fCalculateBackground; //flag
995 AliGammaConversionBGHandler *fBGEventHandler; // background handler
996 Bool_t fBGEventInitialized; //flag
998 AliESDtrackCuts *fEsdTrackCuts; // track cuts
999 Int_t fNumberOfESDTracks; //track counter
1001 static AliESDpid* fgESDpid; // ESD pid object
1003 Int_t fNEventsForBGCalculation; // Number of events used for background calculation
1005 Bool_t fUseChargedTrackMultiplicityForBG; // flag
1006 Int_t fNumberOfGoodV0s; // number of good V0s
1007 Int_t fIsHeavyIon; // flag
1008 Bool_t fUseCorrectedTPCClsInfo;
1009 Int_t fUseMCPSmearing;
1010 Double_t fPBremSmearing;
1011 Double_t fPSigSmearing;
1012 Double_t fPSigSmearingCte;
1015 Bool_t fDoPhotonAsymmetryCut; // flag to use the PhotonAsymetryCut
1016 Double_t fMinPPhotonAsymmetryCut;
1017 Double_t fMinPhotonAsymmetry;
1018 ClassDef(AliV0Reader,22) // RRnewTOF
1021 inline void AliV0Reader::InitESDpid(Int_t type)
1024 // initialize PID parameters
1025 // type=0 is simulation
1028 if (!fgESDpid) fgESDpid=new AliESDpid;
1029 Double_t alephParameters[5];
1031 alephParameters[0] = 2.15898e+00/50.;
1032 alephParameters[1] = 1.75295e+01;
1033 alephParameters[2] = 3.40030e-09;
1034 alephParameters[3] = 1.96178e+00;
1035 alephParameters[4] = 3.91720e+00;
1036 fgESDpid->GetTOFResponse().SetTimeResolution(80.);
1040 alephParameters[0] = 0.0283086/0.97;
1041 alephParameters[1] = 2.63394e+01;
1042 alephParameters[2] = 5.04114e-11;
1043 alephParameters[3] = 2.12543e+00;
1044 alephParameters[4] = 4.88663e+00;
1045 fgESDpid->GetTOFResponse().SetTimeResolution(130.);
1046 fgESDpid->GetTPCResponse().SetMip(50.);
1049 fgESDpid->GetTPCResponse().SetBetheBlochParameters(
1050 alephParameters[0],alephParameters[1],alephParameters[2],
1051 alephParameters[3],alephParameters[4]);
1053 fgESDpid->GetTPCResponse().SetSigma(3.79301e-03, 2.21280e+04);