AliPhysics  v5-06-11-01 (156c7f3)
 All Classes Namespaces Files Functions Variables Enumerations Enumerator Macros
AliNeutralMesonSelection.h
Go to the documentation of this file.
1 #ifndef ALINEUTRALMESONSELECTION_H
2 #define ALINEUTRALMESONSELECTION_H
3 /* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4  * See cxx source for full Copyright notice */
5 
6 //_________________________________________________________________________
18 //_________________________________________________________________________
19 
20 // --- ROOT system ---
21 #include<TObject.h>
22 #include<TArrayD.h>
23 #include<TString.h>
24 
25 class TLorentzVector ;
26 class TList ;
27 class TH2F ;
28 
29 class AliNeutralMesonSelection : public TObject {
30 
31  public:
32  AliNeutralMesonSelection() ; // default ctor
33 
35  virtual ~AliNeutralMesonSelection() { ; }
36 
37  // General
38 
39  TList * GetCreateOutputObjects();
40 
41  void InitParameters();
42 
43  void Print(const Option_t * opt) const;
44 
47 
48  Bool_t SelectPair(TLorentzVector particlei, TLorentzVector particlej, Int_t calo) ;
49 
50  void SetParticle(TString particleName) ; // Do some default settings for "Pi0" or "Eta"
51  TString GetParticle() const { return fParticle ; }
52 
53  // Asymmetry selection
54 
55  Float_t GetAsymmetryCut() const { return fAsymmetryCut ; }
56  void SetAsymmetryCut(Float_t asy) { fAsymmetryCut = asy ; }
57 
60 
61  //Opening angle selection
62 
63  Double_t GetAngleMaxParam(Int_t i) const { return fAngleMaxParam.At(i) ; }
64  void SetAngleMaxParam(Int_t i, Double_t par) { fAngleMaxParam.AddAt(par,i) ; }
65 
66  void SetShiftMinAngleCut(Float_t a, Float_t b) { fShiftMinAngle[0] = a ;
67  fShiftMinAngle[1] = b ; }
68 
69  void SwitchOnAngleSelection() { fUseAngleCut = kTRUE ; }
70  void SwitchOffAngleSelection() { fUseAngleCut = kFALSE ; }
71 
72  Bool_t IsAngleInWindow(Float_t angle, Float_t e) const ;
73 
74  //Invariant mass selection
75 
76  Double_t GetInvMassMaxCut() const { return fInvMassMaxCut ; }
77  Double_t GetInvMassMinCut() const { return fInvMassMinCut ; }
78 
79  void SetInvMassCutRange(Double_t invmassmin, Double_t invmassmax)
80  { fInvMassMaxCut =invmassmax; fInvMassMinCut =invmassmin ; }
81 
82  void SetSideBandCutRanges( Double_t lmin, Double_t lmax,
83  Double_t rmin, Double_t rmax )
84  { fLeftBandMinCut = lmin ; fLeftBandMaxCut = lmax ;
85  fRightBandMinCut = rmin ; fRightBandMaxCut = rmax ; }
86 
87  void SetInvMassCutMaxParameters(Float_t a, Float_t b, Float_t c)
88  { fInvMassMaxCutParam[0] = a ;
89  fInvMassMaxCutParam[1] = b ;
90  fInvMassMaxCutParam[2] = c ; }
91 
92  Double_t GetMass() const { return fM ; }
93  void SetMass(Double_t m) { fM = m ; }
94 
95 
96  // Decay photon bit setting
97  static const Int_t fgkMaxNDecayBits = 8;
98 
99  enum decayTypes { kPi0 = 0, kEta = 1,
102 
103  UInt_t GetDecayBit() const { return fDecayBit ; }
104 
105  void SetDecayBit(Int_t &tag, UInt_t set) const {
106  // Set bit of type set (decayTypes) in tag
107  tag |= (1<<set) ;
108  }
109 
110  void SetDecayBit(Int_t &tag) const {
111  // Set bit of type set (decayTypes) in tag
112  tag |= (1<<fDecayBit) ;
113  }
114 
115  Bool_t CheckDecayBit(Int_t tag, UInt_t test) const {
116  // Check if in tag the bit test (decayTypes) is set.
117  if (tag & (1<<test) ) return kTRUE ;
118  else return kFALSE ;
119  }
120 
121  Bool_t CheckDecayBit(Int_t tag) const {
122  // Check if in tag the bit test (decayTypes) is set.
123  if (tag & (1<<fDecayBit) ) return kTRUE ;
124  else return kFALSE ;
125  }
126 
127  // Histograms setters and getters
128 
129  virtual void SetHistoERangeAndNBins(Float_t min, Float_t max, Int_t n) {
130  fHistoNEBins = n ;
131  fHistoEMax = max ;
132  fHistoEMin = min ;
133  }
134 
135  Int_t GetHistoNEBins() const { return fHistoNEBins ; }
136  Float_t GetHistoEMin() const { return fHistoEMin ; }
137  Float_t GetHistoEMax() const { return fHistoEMax ; }
138 
139  virtual void SetHistoAngleRangeAndNBins(Float_t min, Float_t max, Int_t n) {
140  fHistoNAngleBins = n ;
141  fHistoAngleMax = max ;
142  fHistoAngleMin = min ;
143  }
144 
145  Int_t GetHistoNAngleBins() const { return fHistoNAngleBins ; }
146  Float_t GetHistoAngleMin() const { return fHistoAngleMin ; }
147  Float_t GetHistoAngleMax() const { return fHistoAngleMax ; }
148 
149  virtual void SetHistoIMRangeAndNBins(Float_t min, Float_t max, Int_t n) {
150  fHistoNIMBins = n ;
151  fHistoIMMax = max ;
152  fHistoIMMin = min ;
153  }
154 
155  Int_t GetHistoNIMBins() const { return fHistoNIMBins ; }
156  Float_t GetHistoIMMin() const { return fHistoIMMin ; }
157  Float_t GetHistoIMMax() const { return fHistoIMMax ; }
158 
159  private:
160 
161  Float_t fAsymmetryCut ;
162 
164 
165  Double_t fM ;
166 
167  Double_t fInvMassMaxCut ;
168 
169  Double_t fInvMassMinCut ;
170 
171  Double_t fInvMassMaxCutParam[3];
172 
173  Double_t fLeftBandMinCut ;
174 
175  Double_t fLeftBandMaxCut ;
176 
177  Double_t fRightBandMinCut ;
178 
179  Double_t fRightBandMaxCut ;
180 
181  TArrayD fAngleMaxParam ;
182 
183  Bool_t fUseAngleCut ;
184 
185  Float_t fShiftMinAngle[2] ;
186 
188 
189  TString fParticle ;
190 
191  UInt_t fDecayBit;
192 
193  //Histograms
195 
197 
199 
201 
202 
204 
206 
208 
210 
212 
214 
216 
217 
218  //Histograms binning and range
219  Int_t fHistoNEBins ;
220 
221  Float_t fHistoEMax ;
222 
223  Float_t fHistoEMin ;
224 
226 
227  Float_t fHistoAngleMax ;
228 
229  Float_t fHistoAngleMin ;
230 
231 
232  Int_t fHistoNIMBins ;
233 
234  Float_t fHistoIMMax ;
235 
236  Float_t fHistoIMMin ;
237 
240 
243 
245  ClassDef(AliNeutralMesonSelection,8) ;
247 
248 } ;
249 
250 #endif //ALINEUTRALMESONSELECTION_H
251 
252 
253 
Double_t GetAngleMaxParam(Int_t i) const
Float_t fHistoAngleMin
Minimum value of angle histogram range.
Float_t fHistoAngleMax
Maximum value of angle histogram range.
Double_t fLeftBandMaxCut
Side Band selection, max left band cut.
Float_t fAsymmetryCut
Asymmetry cut.
virtual void SetHistoAngleRangeAndNBins(Float_t min, Float_t max, Int_t n)
void InitParameters()
Initialize the parameters of the analysis.
Bool_t AreNeutralMesonSelectionHistosKept() const
TH2F * fhInvMassPairAllCut
! Invariant mass of decay photons, all cuts.
AliNeutralMesonSelection()
Default constructor. Initialize parameters.
Double_t fRightBandMinCut
Side Band selection, min right band cut.
void SetAngleMaxParam(Int_t i, Double_t par)
void SetDecayBit(Int_t &tag) const
void Print(const Option_t *opt) const
Print some relevant parameters set for the analysis.
void SetInvMassCutRange(Double_t invmassmin, Double_t invmassmax)
void KeepNeutralMesonSelectionHistos(Bool_t keep)
Bool_t CheckDecayBit(Int_t tag, UInt_t test) const
Double_t fLeftBandMinCut
Side Band selection, min left band cut.
TH2F * fhAsymmetryOpeningAngleCut
! Asymmetry of decay photons, cut on opening angle.
virtual ~AliNeutralMesonSelection()
Virtual destructor.
TH2F * fhInvMassPairNoCut
! Invariant mass of decay photons, no cuts.
Bool_t fKeepNeutralMesonHistos
Keep neutral meson selection histograms.
Bool_t IsAngleInWindow(Float_t angle, Float_t e) const
Bool_t fUseAngleCut
Select pairs depending on their opening angle.
Double_t fRightBandMaxCut
Side Band selection, max right band cut.
void SetSideBandCutRanges(Double_t lmin, Double_t lmax, Double_t rmin, Double_t rmax)
Int_t fHistoNIMBins
Number of bins in Invariant Mass axis.
Float_t fHistoEMax
Maximum value of pi0 E histogram range.
Double_t fInvMassMaxCutParam[3]
Variable invariant mass max cut, for pi0 in EMCAL.
Double_t fInvMassMaxCut
Invariant Mass cut maximum.
virtual void SetHistoERangeAndNBins(Float_t min, Float_t max, Int_t n)
TArrayD fAngleMaxParam
Maximum opening angle selection parameters.
void SetParticle(TString particleName)
Set some default parameters for selection of pi0 or eta.
Float_t fShiftMinAngle[2]
Correction shift for min angle from true kinematic limit, resolution effects.
Bool_t SelectPair(TLorentzVector particlei, TLorentzVector particlej, Int_t calo)
TH2F * fhAnglePairAllCut
! Aperture angle of decay photons, all cuts.
Double_t fM
Mass of the neutral meson.
Int_t fHistoNEBins
Number of bins in pi0 E axis.
TH2F * fhAnglePairNoCut
! Aperture angle of decay photons, no cuts.
void SetDecayBit(Int_t &tag, UInt_t set) const
TH2F * fhAsymmetryNoCut
! Asymmetry of decay photons, no cuts.
Float_t fHistoEMin
Minimum value of pi0 E histogram range.
Float_t fHistoIMMin
Minimum value of Invariant Mass histogram range.
virtual void SetHistoIMRangeAndNBins(Float_t min, Float_t max, Int_t n)
TH2F * fhAnglePairOpeningAngleCut
! Aperture angle of decay photons, cut on opening angle.
UInt_t fDecayBit
Decay type flag, set while selecting, depending on fParticle and side range. See enum decayTypes for ...
Bool_t CheckDecayBit(Int_t tag) const
Float_t fHistoIMMax
Maximum value of Invariant Mass histogram range.
Class that contains methods to select candidate cluster pairs to neutral meson.
AliNeutralMesonSelection & operator=(const AliNeutralMesonSelection &nm)
Assignment operator not implemented.
TH2F * fhInvMassPairAsymmetryCut
! Invariant mass of decay photons, asymmetry cut.
Bool_t fUseAsymmetryCut
Use the asymmetry cut.
Double_t fInvMassMinCut
Invariant Masscut minimun.
Int_t fHistoNAngleBins
Number of bins in angle axis.
void SetInvMassCutMaxParameters(Float_t a, Float_t b, Float_t c)
TH2F * fhAsymmetryAllCut
! Asymmetry of decay photons, all cuts.
TH2F * fhAnglePairAsymmetryCut
! Aperture angle of decay photons, asymmetry cut.
TString fParticle
Meutral meson name (Pi0, Eta, +SideBand).
TH2F * fhInvMassPairOpeningAngleCut
! Invariant mass of decay photons, cut on opening angle.
void SetShiftMinAngleCut(Float_t a, Float_t b)