]>
Commit | Line | Data |
---|---|---|
10d100d4 | 1 | #ifndef ALITPCPIDRESPONSE_H |
2 | #define ALITPCPIDRESPONSE_H | |
8c6a71ab | 3 | /* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. * |
4 | * See cxx source for full Copyright notice */ | |
5 | ||
6 | //------------------------------------------------------- | |
7 | // TPC PID class | |
8 | // A very naive design... Should be made better by the detector experts... | |
9 | // Origin: Iouri Belikov, CERN, Jouri.Belikov@cern.ch | |
aea7a46d | 10 | // With many additions and modifications suggested by |
11 | // Alexander Kalweit, GSI, alexander.philipp.kalweit@cern.ch | |
12 | // Dariusz Miskowiec, GSI, D.Miskowiec@gsi.de | |
f84b18dd | 13 | // ...and some modifications by |
14 | // Mikolaj Krzewicki, GSI, mikolaj.krzewicki@cern.ch | |
15 | // ...and some modifications plus eta correction functions by | |
16 | // Benjamin Hess, University of Tuebingen, bhess@cern.ch | |
8c6a71ab | 17 | //------------------------------------------------------- |
18 | #include <Rtypes.h> | |
644666df | 19 | |
20 | #include <TNamed.h> | |
21 | #include <TVectorF.h> | |
d2aa6df0 | 22 | #include <TObjArray.h> |
8c6a71ab | 23 | |
aea7a46d | 24 | #include "AliPID.h" |
f84b18dd | 25 | #include "AliVTrack.h" |
aea7a46d | 26 | |
f84b18dd | 27 | class TH2D; |
644666df | 28 | class TSpline3; |
29 | ||
30 | class AliTPCPIDResponse: public TNamed { | |
8c6a71ab | 31 | public: |
10d100d4 | 32 | AliTPCPIDResponse(); |
1b45e564 | 33 | //TODO Remove? AliTPCPIDResponse(const Double_t *param); |
644666df | 34 | AliTPCPIDResponse(const AliTPCPIDResponse&); |
35 | AliTPCPIDResponse& operator=(const AliTPCPIDResponse&); | |
f84b18dd | 36 | virtual ~AliTPCPIDResponse(); |
644666df | 37 | |
38 | enum EChamberStatus { | |
39 | kChamberOff=0, | |
40 | kChamberHighGain=1, | |
41 | kChamberLowGain=2, | |
42 | kChamberInvalid=3 | |
43 | }; | |
44 | ||
45 | enum ETPCgainScenario { | |
46 | kDefault= 0, | |
47 | kALLhigh = 1, | |
48 | kOROChigh = 2, | |
49 | kGainScenarioInvalid = 3 | |
50 | }; | |
51 | ||
52 | static const Int_t fgkNumberOfParticleSpecies=AliPID::kSPECIESC; | |
53 | static const Int_t fgkNumberOfGainScenarios=3; | |
54 | static const Int_t fgkNumberOfdEdxSourceScenarios=3; | |
55 | ||
56 | enum ETPCdEdxSource { | |
57 | kdEdxDefault=0, // use combined dEdx from IROC+OROC (assumes ideal detector) | |
58 | kdEdxOROC=1, // use only OROC | |
59 | kdEdxHybrid=2, // Use IROC+OROC dEdx only where IROCS are good (high gain), otherwise fall back to OROC only | |
60 | kdEdxInvalid=3 //invalid | |
61 | }; | |
62 | ||
10d100d4 | 63 | void SetSigma(Float_t res0, Float_t resN2); |
aea7a46d | 64 | void SetBetheBlochParameters(Double_t kp1, |
65 | Double_t kp2, | |
66 | Double_t kp3, | |
67 | Double_t kp4, | |
68 | Double_t kp5 | |
69 | ); | |
1b45e564 | 70 | //Better prevent user from setting fMIP != 50. because fMIP set fix to 50 for much other code: |
10d100d4 | 71 | void SetMip(Float_t mip) { fMIP = mip; } // Set overall normalisation; mean dE/dx for MIP |
aea7a46d | 72 | Double_t Bethe(Double_t bg) const; |
d2aa6df0 | 73 | void SetUseDatabase(Bool_t useDatabase) { fUseDatabase = useDatabase;} |
644666df | 74 | Bool_t GetUseDatabase() const { return fUseDatabase;} |
d2aa6df0 | 75 | |
76 | void SetResponseFunction(AliPID::EParticleType type, TObject * const o) { fResponseFunctions.AddAt(o,(Int_t)type); } | |
deae51a8 | 77 | const TObject * GetResponseFunction(AliPID::EParticleType type) { return fResponseFunctions.At((Int_t)type); } |
644666df | 78 | void SetVoltage(Int_t n, Float_t v) {fVoltageMap[n]=v;} |
79 | void SetVoltageMap(const TVectorF& a) {fVoltageMap=a;} //resets ownership, ~ will not delete contents | |
80 | Float_t GetVoltage(Int_t n) const {return fVoltageMap[n];} | |
81 | void SetLowGainIROCthreshold(Float_t v) {fLowGainIROCthreshold=v;} | |
82 | void SetBadIROCthreshold(Float_t v) {fBadIROCthreshhold=v;} | |
83 | void SetLowGainOROCthreshold(Float_t v) {fLowGainOROCthreshold=v;} | |
84 | void SetBadOROCthreshold(Float_t v) {fBadOROCthreshhold=v;} | |
85 | void SetMaxBadLengthFraction(Float_t f) {fMaxBadLengthFraction=f;} | |
86 | ||
87 | void SetMagField(Double_t mf) { fMagField=mf; } | |
d2aa6df0 | 88 | |
f84b18dd | 89 | const TH2D* GetEtaCorrMap() const { return fhEtaCorr; }; |
90 | Bool_t SetEtaCorrMap(TH2D* hMap); | |
91 | ||
92 | Double_t GetEtaCorrection(const AliVTrack *track, AliPID::EParticleType species, ETPCdEdxSource dedxSource = kdEdxDefault) const; | |
93 | ||
94 | Double_t GetEtaCorrectedTrackdEdx(const AliVTrack *track, AliPID::EParticleType species, ETPCdEdxSource dedxSource = kdEdxDefault) const; | |
95 | ||
96 | const TH2D* GetSigmaPar1Map() const { return fhEtaSigmaPar1; }; | |
97 | Double_t GetSigmaPar0() const { return fSigmaPar0; }; | |
98 | Bool_t SetSigmaParams(TH2D* hSigmaPar1Map, Double_t sigmaPar0); | |
99 | ||
100 | Double_t GetSigmaPar1(const AliVTrack *track, AliPID::EParticleType species, ETPCdEdxSource dedxSource = kdEdxDefault) const; | |
101 | ||
644666df | 102 | //NEW |
103 | void SetSigma(Float_t res0, Float_t resN2, ETPCgainScenario gainScenario ); | |
644666df | 104 | Double_t GetExpectedSignal( const AliVTrack* track, |
105 | AliPID::EParticleType species, | |
f84b18dd | 106 | ETPCdEdxSource dedxSource = kdEdxDefault, |
107 | Bool_t correctEta = kFALSE) const;//TODO: In future, default kTRUE | |
644666df | 108 | Double_t GetExpectedSigma( const AliVTrack* track, |
109 | AliPID::EParticleType species, | |
f84b18dd | 110 | ETPCdEdxSource dedxSource = kdEdxDefault, |
111 | Bool_t correctEta = kFALSE) const;//TODO: In future, default kTRUE | |
644666df | 112 | Float_t GetNumberOfSigmas( const AliVTrack* track, |
113 | AliPID::EParticleType species, | |
f84b18dd | 114 | ETPCdEdxSource dedxSource = kdEdxDefault, |
115 | Bool_t correctEta = kFALSE) const;//TODO: In future, default kTRUE | |
644666df | 116 | |
117 | void SetResponseFunction(TObject* o, | |
118 | AliPID::EParticleType type, | |
119 | ETPCgainScenario gainScenario); | |
120 | void Print(Option_t* option="") const; | |
121 | TSpline3* GetResponseFunction( AliPID::EParticleType species, | |
122 | ETPCgainScenario gainScenario ) const; | |
123 | TSpline3* GetResponseFunction( const AliVTrack* track, | |
124 | AliPID::EParticleType species, | |
f84b18dd | 125 | ETPCdEdxSource dedxSource = kdEdxDefault) const; |
644666df | 126 | Bool_t ResponseFunctiondEdxN(const AliVTrack* track, |
127 | AliPID::EParticleType species, | |
f84b18dd | 128 | ETPCdEdxSource dedxSource, |
129 | Double_t& dEdx, Int_t& nPoints, ETPCgainScenario& gainScenario, TSpline3** responseFunction) const; | |
f85a3764 | 130 | Bool_t sectorNumbersInOut(Double_t* trackPositionInner, |
131 | Double_t* trackPositionOuter, | |
644666df | 132 | Float_t& phiIn, Float_t& phiOut, |
133 | Int_t& in, Int_t& out ) const; | |
134 | AliTPCPIDResponse::EChamberStatus TrackStatus(const AliVTrack* track, Int_t layer) const; | |
135 | Float_t MaxClusterRadius(const AliVTrack* track) const; | |
136 | Bool_t TrackApex(const AliVTrack* track, Float_t magField, Double_t position[3]) const; | |
137 | static const char* GainScenarioName(Int_t n) {return fgkGainScenarioName[(n>fgkNumberOfdEdxSourceScenarios)?fgkNumberOfdEdxSourceScenarios+1:n];} | |
138 | Int_t ResponseFunctionIndex( AliPID::EParticleType species, | |
139 | ETPCgainScenario gainScenario ) const; | |
140 | void ResetSplines(); | |
141 | ||
644666df | 142 | //OLD |
10d100d4 | 143 | Double_t GetExpectedSignal(const Float_t mom, |
aea7a46d | 144 | AliPID::EParticleType n=AliPID::kKaon) const; |
10d100d4 | 145 | Double_t GetExpectedSigma(const Float_t mom, const Int_t nPoints, |
644666df | 146 | AliPID::EParticleType n=AliPID::kKaon) const; |
147 | Float_t GetNumberOfSigmas(const Float_t mom, | |
148 | const Float_t dEdx, | |
149 | const Int_t nPoints, | |
150 | AliPID::EParticleType n=AliPID::kKaon) const { | |
f84b18dd | 151 | // |
152 | // Deprecated function (for backward compatibility). Please use | |
153 | // GetNumberOfSigmas(const AliVTrack *track, AliPID::EParticleType species, ETPCdEdxSource dedxSource ) | |
154 | // instead!TODO | |
155 | // | |
156 | ||
10d100d4 | 157 | Double_t bethe=GetExpectedSignal(mom,n); |
158 | Double_t sigma=GetExpectedSigma(mom,nPoints,n); | |
159 | return (dEdx-bethe)/sigma; | |
160 | } | |
aea7a46d | 161 | |
10d100d4 | 162 | Double_t GetMIP() const { return fMIP;} |
644666df | 163 | Float_t GetRes0() const { return fRes0[0]; } |
164 | Float_t GetResN2() const { return fResN2[0]; } | |
165 | Float_t GetRes0(ETPCgainScenario s) const { return fRes0[s]; } | |
166 | Float_t GetResN2(ETPCgainScenario s) const { return fResN2[s]; } | |
b981edcd | 167 | |
f84b18dd | 168 | protected: |
169 | Double_t GetExpectedSignal(const AliVTrack* track, | |
170 | AliPID::EParticleType species, | |
171 | Double_t dEdx, | |
172 | const TSpline3* responseFunction, | |
173 | Bool_t correctEta) const; | |
174 | ||
175 | Double_t GetExpectedSigma(const AliVTrack* track, | |
176 | AliPID::EParticleType species, | |
177 | ETPCgainScenario gainScenario, | |
178 | Double_t dEdx, | |
179 | Int_t nPoints, | |
180 | const TSpline3* responseFunction, | |
181 | Bool_t correctEta) const; | |
182 | ||
183 | Double_t GetEtaCorrection(const AliVTrack *track, Double_t dEdxSplines) const; | |
184 | ||
185 | Double_t GetSigmaPar1(const AliVTrack *track, AliPID::EParticleType species, | |
186 | Double_t dEdx, const TSpline3* responseFunction) const; | |
187 | ||
8c6a71ab | 188 | private: |
10d100d4 | 189 | Float_t fMIP; // dEdx for MIP |
644666df | 190 | Float_t fRes0[fgkNumberOfGainScenarios]; // relative dEdx resolution rel sigma = fRes0*sqrt(1+fResN2/npoint) |
191 | Float_t fResN2[fgkNumberOfGainScenarios]; // relative Npoint dependence rel sigma = fRes0*sqrt(1+fResN2/npoint) | |
aea7a46d | 192 | |
193 | Double_t fKp1; // Parameters | |
194 | Double_t fKp2; // of | |
195 | Double_t fKp3; // the ALEPH | |
196 | Double_t fKp4; // Bethe-Bloch | |
197 | Double_t fKp5; // formula | |
198 | ||
d2aa6df0 | 199 | Bool_t fUseDatabase; // flag if fine-tuned database-response or simple ALEPH BB should be used |
644666df | 200 | |
d2aa6df0 | 201 | TObjArray fResponseFunctions; //! ObjArray of response functions individually for each particle |
644666df | 202 | TVectorF fVoltageMap; //!stores a map of voltages wrt nominal for all chambers |
203 | Float_t fLowGainIROCthreshold; //voltage threshold below which the IROC is considered low gain | |
204 | Float_t fBadIROCthreshhold; //voltage threshold for bad IROCS | |
205 | Float_t fLowGainOROCthreshold; //voltage threshold below which the OROC is considered low gain | |
206 | Float_t fBadOROCthreshhold; //voltage threshold for bad OROCS | |
207 | Float_t fMaxBadLengthFraction; //the maximum allowed fraction of track length in a bad sector. | |
208 | ||
644666df | 209 | Int_t sectorNumber(Double_t phi) const; |
210 | ||
211 | Double_t fMagField; //! Magnetic field | |
212 | ||
213 | static const char* fgkGainScenarioName[fgkNumberOfGainScenarios+1]; | |
d2aa6df0 | 214 | |
f84b18dd | 215 | TH2D* fhEtaCorr; //! Map for TPC eta correction |
216 | TH2D* fhEtaSigmaPar1; //! Map for parameter 1 of the dEdx sigma parametrisation | |
217 | ||
218 | Double_t fSigmaPar0; // Parameter 0 of the dEdx sigma parametrisation | |
219 | ||
220 | ClassDef(AliTPCPIDResponse,5) // TPC PID class | |
8c6a71ab | 221 | }; |
222 | ||
223 | #endif | |
224 | ||
225 |