AliPhysics  b752f14 (b752f14)
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
AliV0ReaderV1.h
Go to the documentation of this file.
1 #ifndef ALIV0READERV1_H
2 #define ALIV0READERV1_H
3 
4 #include "AliAnalysisTaskSE.h"
5 #include "AliAODv0.h"
6 #include "AliESDv0.h"
8 #include "AliConvEventCuts.h"
9 #include "AliExternalTrackParam.h"
10 #include "TObject.h"
11 #include "AliMCEvent.h"
12 #include "AliESDEvent.h"
13 #include "AliKFParticle.h"
14 #include "TParticle.h"
15 #include <iterator>
16 #include <vector>
17 #include "AliESDpid.h"
18 #include "TF1.h"
19 #include "TRandom3.h"
20 #include "AliAnalysisManager.h"
21 
23 class TRandom3;
24 class AliStack;
25 class TList;
27 class TString;
28 class TClonesArray;
29 class TH1F;
30 class TH2F;
32 
33 using namespace std;
34 
36 
37  public:
38 
39  class iterator : public std::iterator<std::bidirectional_iterator_tag, AliConversionPhotonBase> {
40  public:
41  enum Direction_t {
42  kForwardDirection = 0,
43  kBackwardDirection = 1
44  };
45  iterator(const AliV0ReaderV1 *reader, Direction_t dir, int position);
46  iterator(const iterator &other);
47  iterator &operator=(const iterator &other);
48  virtual ~iterator() {}
49 
50  bool operator!=(iterator &other) const;
51  iterator operator++(int);
52  iterator &operator++();
53  iterator operator--(int);
54  iterator &operator--();
55  AliConversionPhotonBase *operator*();
56 
57  private:
61  };
62 
63  AliV0ReaderV1(const char *name="V0ReaderV1");
64  virtual ~AliV0ReaderV1(); //virtual destructor
65 
66  void UserCreateOutputObjects();
67  virtual Bool_t Notify();
68  virtual void UserExec(Option_t *option);
69  virtual void Terminate(Option_t *);
70  virtual void Init();
71 
72  Bool_t ProcessEvent (AliVEvent *inputEvent, AliMCEvent *mcEvent=NULL);
73  Bool_t IsEventSelected() {return fEventIsSelected;}
74 
75  // Return Reconstructed Gammas
76  TClonesArray* GetReconstructedGammas() const {return fConversionGammas;}
77  Int_t GetNReconstructedGammas() const {if(fConversionGammas){return fConversionGammas->GetEntriesFast();} else{ return 0;}}
78  AliConversionPhotonBase *operator[](int index) const;
79 
80  AliConversionPhotonCuts* GetConversionCuts() {return fConversionCuts;}
81  AliConvEventCuts* GetEventCuts() {return fEventCuts;}
82  TList* GetCutHistograms() {if(fConversionCuts) {return fConversionCuts->GetCutHistograms();}
83  return NULL;}
84  TList* GetEventCutHistograms() {if(fEventCuts) {return fEventCuts->GetCutHistograms();}
85  return NULL;}
86  TString GetCurrentFileName() {return fCurrentFileName;}
87  // Set Options
88  void SetAddv0sInESDFilter(Bool_t addv0s) {kAddv0sInESDFilter = addv0s;}
89  void CountTracks();
90  void CountTPCoutTracks();
91  void SetConversionCuts(const TString cut);
92  void SetConversionCuts(AliConversionPhotonCuts *cuts) {fConversionCuts=cuts; return;}
93  void SetEventCuts(const TString cut);
94  void SetEventCuts(AliConvEventCuts *cuts) {fEventCuts=cuts; return;}
95 
96  void SetUseOwnXYZCalculation(Bool_t flag) {fUseOwnXYZCalculation=flag; return;}
97  void SetUseConstructGamma(Bool_t flag) {fUseConstructGamma=flag; return;}
98  void SetUseAODConversionPhoton(Bool_t b) {if(b){ cout<<"Setting Outputformat to AliAODConversionPhoton "<<endl;}
99  else { cout<<"Setting Outputformat to AliKFConversionPhoton "<<endl;}
100  kUseAODConversionPhoton=b; return;}
101 
102  void SetCreateAODs(Bool_t k=kTRUE) {fCreateAOD=k; return;}
103  void SetDeltaAODFilename(TString s) {fDeltaAODFilename=s; return;}
104  void SetDeltaAODBranchName(TString string) {fDeltaAODBranchName = string;
105  fRelabelAODs = kTRUE;
106  AliInfo(Form("Set DeltaAOD BranchName to: %s",fDeltaAODBranchName.Data()));
107  AliInfo(Form("Relabeling of AODs has automatically been switched ON!"));
108  return;}
109 
110  void RelabelAODs(Bool_t relabel=kTRUE) {fRelabelAODs=relabel; return;}
111  Bool_t AreAODsRelabeled() {return fRelabelAODs;}
112  Int_t IsReaderPerformingRelabeling() {return fPreviousV0ReaderPerformsAODRelabeling;}
113  void RelabelAODPhotonCandidates(AliAODConversionPhoton *PhotonCandidate);
114  void SetPeriodName(TString name) {fPeriodName = name;
115  AliInfo(Form("Set PeriodName to: %s",fPeriodName.Data()));
116  return;}
117  TString GetPeriodName() {return fPeriodName;}
118  Int_t GetPtHardFromFile() {return fPtHardBin;}
119  Int_t GetNumberOfPrimaryTracks() {return fNumberOfPrimaryTracks;}
120  Int_t GetNumberOfTPCoutTracks() {return fNumberOfTPCoutTracks;}
121  void SetUseMassToZero (Bool_t b) {if(b){ cout<<"enable set mass to zero for AliAODConversionPhoton"<<endl;}
122  else { cout<<"disable set mass to zero for AliAODConversionPhoton "<<endl;}
123  fUseMassToZero=b; return;}
124 
125  void SetProduceV0FindingEfficiency(Bool_t b) {fProduceV0findingEffi = b;
126  if(b) AliInfo("Enabled V0finding Efficiency");
127  return;}
128 
129  Bool_t GetProduceV0FindingEfficiency() {return fProduceV0findingEffi;}
130  TList* GetV0FindingEfficiencyHistograms() {return fHistograms;}
131  void SetProduceImpactParamHistograms(Bool_t b) {fProduceImpactParamHistograms = b;
132  if(b) AliInfo("Producing additional impact parameter histograms");
133  return;}
134 
135  Bool_t GetProduceImpactParamHistograms() {return fProduceImpactParamHistograms;}
136  TList* GetImpactParamHistograms() {return fImpactParamHistograms;}
137 
138  Bool_t ParticleIsConvertedPhoton(AliStack *MCStack, TParticle *particle, Double_t etaMax, Double_t rMax, Double_t zMax);
139  void CreatePureMCHistosForV0FinderEffiESD();
140  void FillRecMCHistosForV0FinderEffiESD( AliESDv0* currentV0);
141  void FillImpactParamHistograms(AliVTrack *ptrack, AliVTrack* ntrack, AliESDv0 *fCurrentV0, AliKFConversionPhoton *fCurrentMotherKF);
142  Bool_t CheckVectorOnly(vector<Int_t> &vec, Int_t tobechecked);
143  Bool_t CheckVectorForDoubleCount(vector<Int_t> &vec, Int_t tobechecked);
144  void SetImprovedPsiPair(Int_t p) {fImprovedPsiPair=p;return;}
145  Int_t GetImprovedPsiPair() {return fImprovedPsiPair;}
146 
147  iterator begin() const {return iterator(this, iterator::kForwardDirection, 0);}
148  iterator end() const {return iterator(this, iterator::kForwardDirection, GetNReconstructedGammas());}
149  iterator rbegin() const {return iterator(this, iterator::kBackwardDirection, GetNReconstructedGammas() -1); }
150  iterator rend() const {return iterator(this, iterator::kBackwardDirection, -1);}
151 
152  protected:
153  // Reconstruct Gammas
154  Bool_t ProcessESDV0s();
155  AliKFConversionPhoton* ReconstructV0(AliESDv0* fCurrentV0,Int_t currentV0Index);
156  void FillAODOutput();
157  void FindDeltaAODBranchName();
158  Bool_t GetAODConversionGammas();
159 
160  // Getter Functions
161  const AliExternalTrackParam* GetExternalTrackParam(AliESDv0 *fCurrentV0, Int_t &tracklabel, Int_t charge);
162  const AliExternalTrackParam* GetExternalTrackParamP(AliESDv0 *fCurrentV0, Int_t &tracklabel) {return GetExternalTrackParam(fCurrentV0,tracklabel,1);}
163  const AliExternalTrackParam* GetExternalTrackParamN(AliESDv0 *fCurrentV0, Int_t &tracklabel) {return GetExternalTrackParam(fCurrentV0,tracklabel,-1);}
164  AliKFParticle* GetPositiveKFParticle(AliAODv0 *fCurrentV0, Int_t fTrackLabel[2]);
165  AliKFParticle* GetNegativeKFParticle(AliAODv0 *fCurrentV0, Int_t fTrackLabel[2]);
166  AliKFParticle* GetPositiveKFParticle(AliESDv0 *fCurrentV0, Int_t fTrackLabel[2]);
167  AliKFParticle* GetNegativeKFParticle(AliESDv0 *fCurrentV0, Int_t fTrackLabel[2]);
168 
169  Bool_t GetConversionPoint(const AliExternalTrackParam *pparam, const AliExternalTrackParam *nparam, Double_t convpos[3], Double_t dca[2]);
170  Bool_t GetHelixCenter(const AliExternalTrackParam *track, Double_t center[2]);
171  Double_t GetPsiPair(const AliESDv0* v0, const AliExternalTrackParam *positiveparam, const AliExternalTrackParam *negativeparam, const Double_t convpos[3]) const;
172  Bool_t kAddv0sInESDFilter; // Add PCM v0s to AOD created in ESD filter
173  TBits *fPCMv0BitField; // Pointer to bitfield of PCM v0s
174  AliConversionPhotonCuts *fConversionCuts; // Pointer to the ConversionCut Selection
175  AliConvEventCuts *fEventCuts; // Pointer to the ConversionCut Selection
176  TClonesArray *fConversionGammas; // TClonesArray holding the reconstructed photons
177  Bool_t fUseImprovedVertex; // set flag to improve primary vertex estimation by adding photons
178  Bool_t fUseOwnXYZCalculation; //flag that determines if we use our own calculation of xyz (markus)
179  Bool_t fUseConstructGamma; //flag that determines if we use ConstructGamma method from AliKF
180  Bool_t kUseAODConversionPhoton; // set flag to use AOD instead of KF output format for photons
181  Bool_t fCreateAOD; // set flag for AOD creation
182  TString fDeltaAODBranchName; // File where Gamma Conv AOD is located, if not in default AOD
183  TString fDeltaAODFilename; // set filename for delta/satellite aod
185  Int_t fPreviousV0ReaderPerformsAODRelabeling; // 0->not set, meaning V0Reader has not yet determined if it should do AODRelabeling, 1-> V0Reader perfomrs relabeling, 2-> previous V0Reader in list perfomrs relabeling
187  Int_t fNumberOfPrimaryTracks; // Number of Primary Tracks in AOD or ESD
188  Int_t fNumberOfTPCoutTracks; // Number of TPC Tracks with TPCout flag
190  Int_t fPtHardBin; // ptHard bin from file
191  Bool_t fUseMassToZero; // switch on setting the mass to 0 for AODConversionPhotons
192  Bool_t fProduceV0findingEffi; // enable histograms for V0finding efficiency
193  Bool_t fProduceImpactParamHistograms; // enable histograms of impact parameters
194  Float_t fCurrentInvMassPair; // Invariant mass of the pair
195  Int_t fImprovedPsiPair; // enables the calculation of PsiPair after the precise calculation of R and use of the proper function for propagation
196  TList *fHistograms; // list of histograms for V0 finding efficiency
197  TList *fImpactParamHistograms; // list of histograms of impact parameters
198  TH2F *fHistoMCGammaPtvsR; // histogram with all converted gammas vs Pt and R (eta < 0.9)
199  TH2F *fHistoMCGammaPtvsPhi; // histogram with all converted gammas vs Pt and Phi (eta < 0.9)
200  TH2F *fHistoMCGammaPtvsEta; // histogram with all converted gammas vs Pt and Eta
201  TH2F *fHistoMCGammaRvsPhi; // histogram with all converted gammas vs R and Phi (eta < 0.9)
202  TH2F *fHistoMCGammaRvsEta; // histogram with all converted gammas vs R and Eta
203  TH2F *fHistoMCGammaPhivsEta; // histogram with all converted gammas vs Phi and Eta
204  TH2F *fHistoRecMCGammaPtvsR; // histogram with all reconstructed converted gammas vs Pt and R (eta < 0.9)
205  TH2F *fHistoRecMCGammaPtvsPhi; // histogram with all reconstructed converted gammas vs Pt and Phi (eta < 0.9)
206  TH2F *fHistoRecMCGammaPtvsEta; // histogram with all reconstructed converted gammas vs Pt and Eta
207  TH2F *fHistoRecMCGammaRvsPhi; // histogram with all reconstructed converted gammas vs R and Phi (eta < 0.9)
208  TH2F *fHistoRecMCGammaRvsEta; // histogram with all reconstructed converted gammas vs R and Eta
209  TH2F *fHistoRecMCGammaPhivsEta; // histogram with all reconstructed converted gammas vs Phi and Eta
210  TH1F *fHistoRecMCGammaMultiPt; // histogram with all at least double counted photons vs Pt (eta < 0.9)
211  TH2F *fHistoRecMCGammaMultiPtvsEta; // histogram with all at least double counted photons vs Pt vs Eta
212  TH1F *fHistoRecMCGammaMultiR; // histogram with all at least double counted photons vs R (eta < 0.9)
213  TH1F *fHistoRecMCGammaMultiPhi; // histogram with all at least double counted photons vs Phi (eta < 0.9)
214  TH1F *fHistoPosTrackImpactParamZ; //impact parameter z of positive track of V0
226  TH2F *fHistoImpactParamZvsR; // conversion point z vs conversion radius
227  TH2F *fHistoImpactParamZvsR2; // after cuts
228  TH1F *fHistoPt;
229  TH1F *fHistoPt2; // Pt after Impact parameter and causality cuts
231  TH1F *fHistoDCAzPhoton2; // photon dca after impact parameter and causality cuts
232  TH1F *fHistoR; // conversion radius
233  TH1F *fHistoRrecalc; // recalculated conversion radius
234  TH1F *fHistoRviaAlpha; // conversion radius
235  TH1F *fHistoRviaAlphaRecalc; // recalculated conversion radius
236  TH1F *fHistoRdiff; // difference in R between conflict cluster and conversion radius
237  TH1F *fHistoImpactParameterStudy; // info about which cut rejected how many V0s
238  TTree *fImpactParamTree; // tree with y, pt and conversion radius
239 
240  vector<Int_t> fVectorFoundGammas; // vector with found MC labels of gammas
241  TString fCurrentFileName; // current file name
242  Bool_t fMCFileChecked; // vector with MC file names which are broken
243 
244  private:
245  AliV0ReaderV1(AliV0ReaderV1 &original);
246  AliV0ReaderV1 &operator=(const AliV0ReaderV1 &ref);
247 
248  ClassDef(AliV0ReaderV1, 16)
249 
250 };
251 
253  if(fConversionCuts != NULL){
254  delete fConversionCuts;
255  fConversionCuts=NULL;
256  }
257  if(fConversionCuts == NULL){
258  fConversionCuts=new AliConversionPhotonCuts("V0ReaderCuts","V0ReaderCuts");
259  fConversionCuts->InitializeCutsFromCutString(cut.Data());
260  }
261 }
262 
263 inline void AliV0ReaderV1::SetEventCuts(const TString cut){
264  if(fEventCuts != NULL){
265  delete fEventCuts;
266  fEventCuts=NULL;
267  }
268  if(fEventCuts == NULL){
269  fEventCuts=new AliConvEventCuts("V0ReaderEventCuts","V0ReaderEventCuts");
270  fEventCuts->InitializeCutsFromCutString(cut.Data());
271  }
272 }
273 
274 #endif
Int_t IsReaderPerformingRelabeling()
Int_t charge
Direction_t fDirection
Iterator in forward direction.
Definition: AliV0ReaderV1.h:60
void SetCreateAODs(Bool_t k=kTRUE)
TH1F * fHistoNegTrackImpactParamZ
TH2F * fHistoPosTrackImpactParamZvsPt
double Double_t
Definition: External.C:58
Bool_t fUseConstructGamma
TH2F * fHistoRecMCGammaPhivsEta
Bool_t fMCFileChecked
Definition: External.C:236
Int_t GetNReconstructedGammas() const
Definition: AliV0ReaderV1.h:77
void SetConversionCuts(const TString cut)
TH2F * fHistoImpactParamZvsR
TH1F * fHistoPosTrackImpactParamX
TH1F * fHistoRviaAlphaRecalc
Bool_t AreAODsRelabeled()
const AliV0ReaderV1 * fkData
V0 reader used to iterate over.
Definition: AliV0ReaderV1.h:58
Bool_t IsEventSelected()
Definition: AliV0ReaderV1.h:73
TString fDeltaAODFilename
TH1F * fHistoRecMCGammaMultiPhi
Int_t fNumberOfPrimaryTracks
Int_t GetNumberOfPrimaryTracks()
Float_t fCurrentInvMassPair
TH2F * fHistoRecMCGammaPtvsR
TList * GetEventCutHistograms()
Definition: AliV0ReaderV1.h:84
void SetUseOwnXYZCalculation(Bool_t flag)
Definition: AliV0ReaderV1.h:96
TString GetCurrentFileName()
Definition: AliV0ReaderV1.h:86
TString GetPeriodName()
TH1F * fHistoPosTrackImpactParamY
TH1F * fHistoDCAzPhoton
TList * GetCutHistograms()
Definition: AliV0ReaderV1.h:82
void SetProduceV0FindingEfficiency(Bool_t b)
TH1F * fHistoDCAzPhoton2
TH1F * fHistoNegTrackImpactParamX
TH2F * fHistoMCGammaPhivsEta
AliConvEventCuts * fEventCuts
void SetImprovedPsiPair(Int_t p)
const AliExternalTrackParam * GetExternalTrackParamN(AliESDv0 *fCurrentV0, Int_t &tracklabel)
Int_t GetImprovedPsiPair()
TH2F * fHistoNegTrackImpactParamXvsPt
void SetConversionCuts(AliConversionPhotonCuts *cuts)
Definition: AliV0ReaderV1.h:92
TH2F * fHistoRecMCGammaPtvsEta
Bool_t GetProduceV0FindingEfficiency()
Bool_t kUseAODConversionPhoton
TList * fHistograms
Bool_t GetProduceImpactParamHistograms()
iterator begin() const
Int_t fNumberOfTPCoutTracks
iterator rbegin() const
int fCurrentIndex
Index of the current element.
Definition: AliV0ReaderV1.h:59
AliConversionPhotonCuts * fConversionCuts
int Int_t
Definition: External.C:63
iterator end() const
TH1F * fHistoRrecalc
TH1F * fHistoPosTrackImpactParamZ
Bool_t fProduceImpactParamHistograms
float Float_t
Definition: External.C:68
void SetEventCuts(AliConvEventCuts *cuts)
Definition: AliV0ReaderV1.h:94
Int_t GetNumberOfTPCoutTracks()
AliConversionPhotonCuts * GetConversionCuts()
Definition: AliV0ReaderV1.h:80
TH1F * fHistoRecMCGammaMultiR
TString fCurrentFileName
Class handling all kinds of selection cuts for Gamma Conversion analysis.
TH1F * fHistoRviaAlpha
void SetUseAODConversionPhoton(Bool_t b)
Definition: AliV0ReaderV1.h:98
TH1F * fHistoImpactParameterStudy
TH2F * fHistoNegTrackImpactParamZvsPt
TString fPeriodName
TBits * fPCMv0BitField
Bool_t fUseImprovedVertex
void SetUseConstructGamma(Bool_t flag)
Definition: AliV0ReaderV1.h:97
TH2F * fHistoRecMCGammaPtvsPhi
TString fDeltaAODBranchName
TH2F * fHistoPosTrackImpactParamXvsPt
TH1F * fHistoNegTrackImpactParamY
TH2F * fHistoMCGammaRvsEta
Int_t fPreviousV0ReaderPerformsAODRelabeling
TList * GetImpactParamHistograms()
vector< Int_t > fVectorFoundGammas
Bool_t fUseOwnXYZCalculation
TClonesArray * GetReconstructedGammas() const
Definition: AliV0ReaderV1.h:76
Bool_t fEventIsSelected
Int_t fImprovedPsiPair
Int_t GetPtHardFromFile()
void SetProduceImpactParamHistograms(Bool_t b)
Class handling all kinds of selection cuts for Gamma Conversion analysis.
TH2F * fHistoPosTrackImpactParamYvsPt
TH1F * fHistoRdiff
TH2F * fHistoRecMCGammaRvsPhi
Bool_t fProduceV0findingEffi
TH2F * fHistoNegTrackImpactParamYvsPt
TH1F * fHistoRecMCGammaMultiPt
AliConvEventCuts * GetEventCuts()
Definition: AliV0ReaderV1.h:81
TList * fImpactParamHistograms
const char Option_t
Definition: External.C:48
TClonesArray * fConversionGammas
Bool_t fRelabelAODs
const AliExternalTrackParam * GetExternalTrackParamP(AliESDv0 *fCurrentV0, Int_t &tracklabel)
void SetDeltaAODFilename(TString s)
TH2F * fHistoMCGammaRvsPhi
void SetDeltaAODBranchName(TString string)
TH2F * fHistoMCGammaPtvsPhi
bool Bool_t
Definition: External.C:53
TH2F * fHistoMCGammaPtvsR
void SetAddv0sInESDFilter(Bool_t addv0s)
Definition: AliV0ReaderV1.h:88
Bool_t fUseMassToZero
TTree * fImpactParamTree
void SetPeriodName(TString name)
void SetUseMassToZero(Bool_t b)
void SetEventCuts(const TString cut)
TH2F * fHistoMCGammaPtvsEta
void RelabelAODs(Bool_t relabel=kTRUE)
TH2F * fHistoImpactParamZvsR2
iterator rend() const
Bool_t kAddv0sInESDFilter
TH2F * fHistoRecMCGammaRvsEta
TList * GetV0FindingEfficiencyHistograms()
TDirectoryFile * dir
TH2F * fHistoRecMCGammaMultiPtvsEta