]> git.uio.no Git - u/mrichter/AliRoot.git/blob - PWGGA/GammaConv/AliDalitzElectronCuts.h
set magnetic field if not set
[u/mrichter/AliRoot.git] / PWGGA / GammaConv / AliDalitzElectronCuts.h
1 #ifndef ALIDALITZELECTRONCUTS_H
2 #define ALIDALITZELECTRONCUTS_H
3
4 // Class handling all kinds of selection cuts for electrons
5
6 // Authors: Svein Lindal, Daniel Lohner                                                                                         *
7
8
9 #include "AliAODpidUtil.h"
10 //#include "AliConversionPhotonBase.h"
11 //#include "AliAODConversionMother.h"
12 #include "AliAODTrack.h"
13 #include "AliESDtrack.h"
14 #include "AliVTrack.h"
15 #include "AliAODTrack.h"
16 #include "AliStack.h"
17 #include "AliAnalysisCuts.h"
18 #include "AliESDtrackCuts.h"
19 #include "TH1F.h"
20
21 class AliESDEvent;
22 class AliAODEvent;
23 class AliConversionPhotonBase;
24 class AliKFVertex;
25 class AliKFParticle;
26 class TH1F;
27 class TH2F;
28 class AliPIDResponse;
29 class AliAnalysisCuts;
30 class iostream;
31 class TList;
32 class AliAnalysisManager;
33
34
35 using namespace std;
36
37 class AliDalitzElectronCuts : public AliAnalysisCuts {
38         
39  public: 
40
41
42   enum cutIds {
43         kgoodId=0, 
44         kededxSigmaITSCut,
45         kededxSigmaTPCCut,
46         kpidedxSigmaTPCCut,
47         kpiMinMomdedxSigmaTPCCut,
48         kpiMaxMomdedxSigmaTPCCut,
49         kLowPRejectionSigmaCut,
50         kTOFelectronPID,
51         kclsITSCut,
52         kclsTPCCut,
53         ketaCut,
54         kPsiPair,
55         kRejectSharedElecGamma,
56         kBackgroundScheme,
57         kNumberOfRotations,
58         kptCut,
59         kDCACut,
60         kmassCut,
61         kWeights,
62         kNCuts
63   };
64
65
66  enum electronCuts {
67       kElectronIn=0,
68       kNoTracks,
69       kTrackCuts,
70       kdEdxCuts,
71       kElectronOut
72  };
73
74
75   Bool_t SetCutIds(TString cutString); 
76   Int_t fCuts[kNCuts];
77   Bool_t SetCut(cutIds cutID, Int_t cut);
78   Bool_t UpdateCutString();
79   static const char * fgkCutNames[kNCuts];
80
81
82   Bool_t InitializeCutsFromCutString(const TString analysisCutSelection); 
83   
84
85   AliDalitzElectronCuts(const char *name="ElectronCuts", const char * title="Electron Cuts");
86   virtual ~AliDalitzElectronCuts();                            //virtual destructor
87
88   virtual Bool_t IsSelected(TObject* /*obj*/){return kTRUE;}
89   virtual Bool_t IsSelected(TList* /*list*/) {return kTRUE;}
90
91   TString GetCutNumber();
92
93     // Cut Selection
94   Bool_t ElectronIsSelectedMC(Int_t labelParticle,AliStack *fMCStack);
95   Bool_t TrackIsSelected(AliESDtrack* lTrack);
96   Bool_t ElectronIsSelected(AliESDtrack* lTrack);
97   void InitAODpidUtil(Int_t type);
98   static AliDalitzElectronCuts * GetStandardCuts2010PbPb();
99   static AliDalitzElectronCuts * GetStandardCuts2010pp();
100   Bool_t InitPIDResponse();
101   
102   void SetPIDResponse(AliPIDResponse * pidResponse) {fPIDResponse = pidResponse;}
103   AliPIDResponse * GetPIDResponse() { return fPIDResponse;}
104   
105   void PrintCuts();
106
107   void InitCutHistograms(TString name="",Bool_t preCut = kTRUE,TString cutName="");
108   void SetFillCutHistograms(TString name="",Bool_t preCut = kTRUE,TString cutName=""){if(!fHistograms){InitCutHistograms(name,preCut,cutName);};}
109   TList *GetCutHistograms(){return fHistograms;}
110
111   static AliVTrack * GetTrack(AliVEvent * event, Int_t label);
112
113   
114
115   ///Cut functions
116   Bool_t dEdxCuts(AliVTrack * track);
117   Bool_t PIDProbabilityCut(AliConversionPhotonBase *photon, AliVEvent * event);
118   Bool_t RejectSharedElecGamma(TList *photons, Int_t indexEle);
119   Bool_t IsFromGammaConversion( Double_t psiPair, Double_t deltaPhi );
120   Bool_t MassCut(Double_t pi0CandidatePt,Double_t vphotonCandidateMass);
121
122   // Event Cuts
123
124   //Double_t GetPsiPair( const AliESDtrack *trackPos, const AliESDtrack *trackNeg );
125
126   Bool_t SetTPCdEdxCutPionLine(Int_t pidedxSigmaCut);
127   Bool_t SetTPCdEdxCutElectronLine(Int_t ededxSigmaCut);
128   Bool_t SetITSdEdxCutElectronLine(Int_t ededxSigmaCut);
129   Bool_t SetMinMomPiondEdxTPCCut(Int_t piMomdedxSigmaCut);
130   Bool_t SetMaxMomPiondEdxTPCCut(Int_t piMomdedxSigmaCut);
131   Bool_t SetITSClusterCut(Int_t clsITSCut);
132   Bool_t SetTPCClusterCut(Int_t clsTPCCut);
133   Bool_t SetEtaCut(Int_t etaCut);
134   Bool_t SetPtCut(Int_t ptCut);
135   Bool_t SetDCACut(Int_t dcaCut);
136   void SetEtaShift(Double_t etaShift){fEtaShift = etaShift;}
137   Bool_t SetMinMomPiondEdxCut(Int_t piMinMomdedxSigmaCut);
138   Bool_t SetMaxMomPiondEdxCut(Int_t piMaxMomdedxSigmaCut);
139   Bool_t SetLowPRejectionCuts(Int_t LowPRejectionSigmaCut);
140   Bool_t SetTOFElectronPIDCut(Int_t TOFelectronPID);
141   Bool_t SetPsiPairCut(Int_t psiCut);
142   Bool_t SetRejectSharedElecGamma(Int_t RCut);
143   Bool_t SetBackgroundScheme(Int_t BackgroundScheme);
144   Bool_t SetNumberOfRotations(Int_t NumberOfRotations);
145   Bool_t SetMassCut(Int_t massCut);
146   Bool_t SetDoWeights(Int_t opc);
147   
148   // Request Flags
149
150   Double_t GetEtaCut(){ return  fEtaCut;}
151   Double_t GetPsiPairCut(){ return fPsiPairCut; }
152   Double_t DoRejectSharedElecGamma(){ return fDoRejectSharedElecGamma;}
153   Double_t DoPsiPairCut(){return fDoPsiPairCut;}
154   Double_t GetNFindableClustersTPC(AliESDtrack* lTrack);
155   Bool_t   UseTrackMultiplicity(){ return fUseTrackMultiplicityForBG;}
156   Int_t    GetBKGMethod(){ return fBKGMethod; }
157   Int_t    NumberOfRotationEvents(){return fnumberOfRotationEventsForBG;}
158   Bool_t   DoMassCut(){return  fDoMassCut;}
159   Double_t GetMassCutLowPt(){return fMassCutLowPt;}
160   Double_t GetMassCutHighPt(){return fMassCutHighPt;}
161   Double_t GetPtMinMassCut(){return fMassCutPtMin;}
162   Bool_t   DoWeights(){return fDoWeights;}
163   
164
165   
166   protected:
167
168   TList *fHistograms;
169   AliPIDResponse *fPIDResponse;
170   AliESDtrackCuts *fesdTrackCuts;
171
172   Double_t fEtaCut; //eta cutç
173   Double_t fEtaShift;
174   Bool_t   fDoEtaCut;
175   Double_t fPtCut;
176   Double_t fRadiusCut; // radius cut
177   Double_t fPsiPairCut;
178   Double_t fDeltaPhiCutMin;
179   Double_t fDeltaPhiCutMax;
180   Int_t fMinClsTPC; // minimum clusters in the TPC
181   Double_t fMinClsTPCToF; // minimum clusters to findable clusters
182   Bool_t   fDodEdxSigmaITSCut; // flag to use the dEdxCut ITS based on sigmas
183   Bool_t   fDodEdxSigmaTPCCut; // flag to use the dEdxCut TPC based on sigmas
184   Bool_t   fDoTOFsigmaCut; // flag to use TOF pid cut RRnewTOF
185   Bool_t   fDoRejectSharedElecGamma; //Reject electrons from the gammas with Radius < RadiusCut
186   Bool_t   fDoPsiPairCut; // PsiPair Cut
187   Double_t fPIDnSigmaAboveElectronLineITS; // sigma cut
188   Double_t fPIDnSigmaBelowElectronLineITS; // sigma cut
189   Double_t fPIDnSigmaAboveElectronLineTPC; // sigma cut
190   Double_t fPIDnSigmaBelowElectronLineTPC; // sigma cut
191   Double_t fPIDnSigmaAbovePionLineTPC;
192   Double_t fPIDnSigmaAbovePionLineTPCHighPt;
193   Double_t fTofPIDnSigmaAboveElectronLine; // sigma cut RRnewTOF
194   Double_t fTofPIDnSigmaBelowElectronLine; // sigma cut RRnewTOF 
195   Double_t fPIDMinPnSigmaAbovePionLineTPC; // sigma cut
196   Double_t fPIDMaxPnSigmaAbovePionLineTPC; // sigma cut
197   Double_t fDoKaonRejectionLowP;   // Kaon rejection at low p
198   Double_t fDoProtonRejectionLowP; // Proton rejection at low p
199   Double_t fDoPionRejectionLowP;   // Pion rejection at low p
200   Double_t fPIDnSigmaAtLowPAroundKaonLine; // sigma cut
201   Double_t fPIDnSigmaAtLowPAroundProtonLine; // sigma cut
202   Double_t fPIDnSigmaAtLowPAroundPionLine; // sigma cut
203   Double_t fPIDMinPKaonRejectionLowP; // Momentum limit to apply kaon rejection
204   Double_t fPIDMinPProtonRejectionLowP; // Momentum limit to apply proton rejection
205   Double_t fPIDMinPPionRejectionLowP; // Momentum limit to apply proton rejection
206
207   Bool_t   fUseCorrectedTPCClsInfo; // flag to use corrected tpc cl info
208   Bool_t   fUseTOFpid; // flag to use tof pid
209   Bool_t   fRequireTOF; //flg to analyze only tracks with TOF signal
210   Bool_t   fUseTrackMultiplicityForBG; // use multiplicity
211   Int_t    fBKGMethod;
212   Int_t    fnumberOfRotationEventsForBG;
213   Bool_t   fDoMassCut;
214   Double_t fMassCutLowPt;
215   Double_t fMassCutHighPt;
216   Double_t fMassCutPtMin;
217   Bool_t   fDoWeights;
218
219
220   // Histograms
221   TObjString *fCutString; // cut number used for analysis
222   TH1F *hCutIndex; // bookkeeping for cuts
223   TH1F *hdEdxCuts;  // bookkeeping for dEdx cuts
224   TH2F *hITSdEdxbefore; // ITS dEdx before cuts
225   TH2F *hITSdEdxafter;
226   TH2F *hTPCdEdxbefore; // TPC dEdx before cuts
227   TH2F *hTPCdEdxafter; // TPC dEdx after cuts
228   TH2F *hTPCdEdxSignalbefore; //TPC dEdx signal before
229   TH2F *hTPCdEdxSignalafter; //TPC dEdx signal  after
230   TH2F *hTOFbefore; // TOF after cuts
231   TH2F *hTOFafter; // TOF after cuts
232   TH2F *hTrackDCAxyPtbefore;
233   TH2F *hTrackDCAxyPtafter;
234   TH2F *hTrackDCAzPtbefore;
235   TH2F *hTrackDCAzPtafter;
236   TH2F *hTrackNFindClsPtTPCbefore;
237   TH2F *hTrackNFindClsPtTPCafter;
238   
239   
240
241 private:
242
243   AliDalitzElectronCuts(const AliDalitzElectronCuts&); // not implemented
244   AliDalitzElectronCuts& operator=(const AliDalitzElectronCuts&); // not implemented
245
246
247   ClassDef(AliDalitzElectronCuts,2)
248 };
249
250
251 inline void AliDalitzElectronCuts::InitAODpidUtil(Int_t type) {
252   if (!fPIDResponse) fPIDResponse = new AliAODpidUtil();
253   Double_t alephParameters[5];
254   // simulation
255   alephParameters[0] = 2.15898e+00/50.;
256   alephParameters[1] = 1.75295e+01;
257   alephParameters[2] = 3.40030e-09;
258   alephParameters[3] = 1.96178e+00;
259   alephParameters[4] = 3.91720e+00;
260   fPIDResponse->GetTOFResponse().SetTimeResolution(80.);
261   
262   // data
263   if (type==1){
264     alephParameters[0] = 0.0283086/0.97;
265     alephParameters[1] = 2.63394e+01;
266     alephParameters[2] = 5.04114e-11;
267     alephParameters[3] = 2.12543e+00;
268     alephParameters[4] = 4.88663e+00;
269     fPIDResponse->GetTOFResponse().SetTimeResolution(130.);
270     fPIDResponse->GetTPCResponse().SetMip(50.);
271   }
272   
273   fPIDResponse->GetTPCResponse().SetBetheBlochParameters(
274     alephParameters[0],alephParameters[1],alephParameters[2],
275     alephParameters[3],alephParameters[4]);
276   
277   fPIDResponse->GetTPCResponse().SetSigma(3.79301e-03, 2.21280e+04);
278 }
279
280
281 #endif