Line data Source code
1 : #ifndef ALITPCPIDRESPONSE_H
2 : #define ALITPCPIDRESPONSE_H
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
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
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
17 : //-------------------------------------------------------
18 : #include <Rtypes.h>
19 :
20 : #include <TNamed.h>
21 : #include <TVectorF.h>
22 : #include <TObjArray.h>
23 : #include <TF1.h>
24 : #include <TString.h>
25 :
26 : #include "AliPID.h"
27 : #include "AliVTrack.h"
28 :
29 : class TH2D;
30 : class TSpline3;
31 : class AliOADBContainer;
32 :
33 : class AliTPCPIDResponse: public TNamed {
34 : public:
35 : AliTPCPIDResponse();
36 : //TODO Remove? AliTPCPIDResponse(const Double_t *param);
37 : AliTPCPIDResponse(const AliTPCPIDResponse&);
38 : AliTPCPIDResponse& operator=(const AliTPCPIDResponse&);
39 : virtual ~AliTPCPIDResponse();
40 :
41 : enum EChamberStatus {
42 : kChamberOff=0,
43 : kChamberHighGain=1,
44 : kChamberLowGain=2,
45 : kChamberInvalid=3
46 : };
47 :
48 : enum ETPCgainScenario {
49 : kDefault= 0,
50 : kALLhigh = 1,
51 : kOROChigh = 2,
52 : kGainScenarioInvalid = 3
53 : };
54 :
55 : static const Int_t fgkNumberOfParticleSpecies=AliPID::kSPECIESC;
56 : static const Int_t fgkNumberOfGainScenarios=3;
57 : static const Int_t fgkNumberOfdEdxSourceScenarios=3;
58 :
59 : enum ETPCdEdxSource {
60 : kdEdxDefault=0, // use combined dEdx from IROC+OROC (assumes ideal detector)
61 : kdEdxOROC=1, // use only OROC
62 : kdEdxHybrid=2, // Use IROC+OROC dEdx only where IROCS are good (high gain), otherwise fall back to OROC only
63 : kdEdxInvalid=3 //invalid
64 : };
65 :
66 : enum ETPCdEdxType {
67 : kdEdxTrack=0,
68 : kdEdxInfo=1
69 : };
70 :
71 : void SetSigma(Float_t res0, Float_t resN2);
72 : void SetBetheBlochParameters(Double_t kp1,
73 : Double_t kp2,
74 : Double_t kp3,
75 : Double_t kp4,
76 : Double_t kp5
77 : );
78 : //Better prevent user from setting fMIP != 50. because fMIP set fix to 50 for much other code:
79 0 : void SetMip(Float_t mip) { fMIP = mip; } // Set overall normalisation; mean dE/dx for MIP
80 : Double_t Bethe(Double_t bg) const;
81 0 : void SetUseDatabase(Bool_t useDatabase) { fUseDatabase = useDatabase;}
82 0 : Bool_t GetUseDatabase() const { return fUseDatabase;}
83 :
84 0 : void SetResponseFunction(AliPID::EParticleType type, TObject * const o) { fResponseFunctions.AddAt(o,(Int_t)type); }
85 268 : const TObject * GetResponseFunction(AliPID::EParticleType type) const { return fResponseFunctions.At((Int_t)type); }
86 0 : void SetVoltage(Int_t n, Float_t v) {fVoltageMap[n]=v;}
87 0 : void SetVoltageMap(const TVectorF& a) {fVoltageMap=a;} //resets ownership, ~ will not delete contents
88 0 : Float_t GetVoltage(Int_t n) const {return fVoltageMap[n];}
89 0 : void SetLowGainIROCthreshold(Float_t v) {fLowGainIROCthreshold=v;}
90 0 : void SetBadIROCthreshold(Float_t v) {fBadIROCthreshhold=v;}
91 0 : void SetLowGainOROCthreshold(Float_t v) {fLowGainOROCthreshold=v;}
92 0 : void SetBadOROCthreshold(Float_t v) {fBadOROCthreshhold=v;}
93 0 : void SetMaxBadLengthFraction(Float_t f) {fMaxBadLengthFraction=f;}
94 :
95 0 : void SetMagField(Double_t mf) { fMagField=mf; }
96 :
97 0 : const TH2D* GetEtaCorrMap() const { return fhEtaCorr; };
98 : Bool_t SetEtaCorrMap(TH2D* hMap);
99 :
100 : Double_t GetTrackTanTheta(const AliVTrack *track) const;
101 :
102 : Double_t GetEtaCorrection(const AliVTrack *track, AliPID::EParticleType species, ETPCdEdxSource dedxSource = kdEdxDefault) const;
103 :
104 : Double_t GetEtaCorrectedTrackdEdx(const AliVTrack *track, AliPID::EParticleType species, ETPCdEdxSource dedxSource = kdEdxDefault) const;
105 :
106 0 : const TH2D* GetSigmaPar1Map() const { return fhEtaSigmaPar1; };
107 0 : Double_t GetSigmaPar0() const { return fSigmaPar0; };
108 : Bool_t SetSigmaParams(TH2D* hSigmaPar1Map, Double_t sigmaPar0);
109 :
110 : Double_t GetSigmaPar1(const AliVTrack *track, AliPID::EParticleType species, ETPCdEdxSource dedxSource = kdEdxDefault) const;
111 :
112 :
113 0 : const TF1* GetMultiplicityCorrectionFunction() const { return fCorrFuncMultiplicity; };
114 : void SetParameterMultiplicityCorrection(Int_t parIndex, Double_t parValue)
115 0 : { if (fCorrFuncMultiplicity) fCorrFuncMultiplicity->SetParameter(parIndex, parValue); };
116 :
117 0 : const TF1* GetMultiplicityCorrectionFunctionTanTheta() const { return fCorrFuncMultiplicityTanTheta; };
118 : void SetParameterMultiplicityCorrectionTanTheta(Int_t parIndex, Double_t parValue)
119 0 : { if (fCorrFuncMultiplicityTanTheta) fCorrFuncMultiplicityTanTheta->SetParameter(parIndex, parValue); };
120 :
121 0 : const TF1* GetMultiplicitySigmaCorrectionFunction() const { return fCorrFuncSigmaMultiplicity; };
122 : void SetParameterMultiplicitySigmaCorrection(Int_t parIndex, Double_t parValue)
123 0 : { if (fCorrFuncSigmaMultiplicity) fCorrFuncSigmaMultiplicity->SetParameter(parIndex, parValue); };
124 :
125 : void ResetMultiplicityCorrectionFunctions();
126 :
127 0 : void SetCurrentEventMultiplicity(Int_t value) { fCurrentEventMultiplicity = value; };
128 0 : Int_t GetCurrentEventMultiplicity() const { return fCurrentEventMultiplicity; };
129 :
130 : Double_t GetMultiplicityCorrection(const AliVTrack *track, AliPID::EParticleType species, ETPCdEdxSource dedxSource = kdEdxDefault) const;
131 :
132 : Double_t GetMultiplicitySigmaCorrection(const AliVTrack *track, AliPID::EParticleType species, ETPCdEdxSource dedxSource = kdEdxDefault) const;
133 :
134 : Double_t GetMultiplicityCorrectedTrackdEdx(const AliVTrack *track, AliPID::EParticleType species, ETPCdEdxSource dedxSource = kdEdxDefault) const;
135 :
136 : Double_t GetEtaAndMultiplicityCorrectedTrackdEdx(const AliVTrack *track, AliPID::EParticleType species,
137 : ETPCdEdxSource dedxSource = kdEdxDefault) const;
138 :
139 : // Fast functions for expert use only
140 : Double_t GetEtaCorrectionFast(const AliVTrack *track, Double_t dEdxSplines) const;
141 :
142 : Double_t GetMultiplicityCorrectionFast(const AliVTrack *track, Double_t dEdxExpected, Int_t multiplicity) const;
143 :
144 : Double_t GetMultiplicitySigmaCorrectionFast(Double_t dEdxExpected, Int_t multiplicity) const;
145 :
146 : Double_t GetSigmaPar1Fast(const AliVTrack *track, AliPID::EParticleType species,
147 : Double_t dEdx, const TSpline3* responseFunction) const;
148 :
149 : //NEW
150 : void SetSigma(Float_t res0, Float_t resN2, ETPCgainScenario gainScenario );
151 : Double_t GetExpectedSignal( const AliVTrack* track,
152 : AliPID::EParticleType species,
153 : ETPCdEdxSource dedxSource = kdEdxDefault,
154 : Bool_t correctEta = kFALSE,
155 : Bool_t correctMultiplicity = kFALSE) const;
156 : Double_t GetExpectedSigma( const AliVTrack* track,
157 : AliPID::EParticleType species,
158 : ETPCdEdxSource dedxSource = kdEdxDefault,
159 : Bool_t correctEta = kFALSE,
160 : Bool_t correctMultiplicity = kFALSE) const;
161 : Float_t GetNumberOfSigmas( const AliVTrack* track,
162 : AliPID::EParticleType species,
163 : ETPCdEdxSource dedxSource = kdEdxDefault,
164 : Bool_t correctEta = kFALSE,
165 : Bool_t correctMultiplicity = kFALSE) const;
166 :
167 : Float_t GetSignalDelta( const AliVTrack* track,
168 : AliPID::EParticleType species,
169 : ETPCdEdxSource dedxSource = kdEdxDefault,
170 : Bool_t correctEta = kFALSE,
171 : Bool_t correctMultiplicity = kFALSE,
172 : Bool_t ratio = kFALSE) const;
173 :
174 : void SetResponseFunction(TObject* o,
175 : AliPID::EParticleType type,
176 : ETPCgainScenario gainScenario);
177 : void Print(Option_t* option="") const;
178 : TSpline3* GetResponseFunction( AliPID::EParticleType species,
179 : ETPCgainScenario gainScenario ) const;
180 : TSpline3* GetResponseFunction( const AliVTrack* track,
181 : AliPID::EParticleType species,
182 : ETPCdEdxSource dedxSource = kdEdxDefault) const;
183 : Bool_t ResponseFunctiondEdxN(const AliVTrack* track,
184 : AliPID::EParticleType species,
185 : ETPCdEdxSource dedxSource,
186 : Double_t& dEdx, Int_t& nPoints, ETPCgainScenario& gainScenario, TSpline3** responseFunction) const;
187 : Bool_t sectorNumbersInOut(Double_t* trackPositionInner,
188 : Double_t* trackPositionOuter,
189 : Float_t& phiIn, Float_t& phiOut,
190 : Int_t& in, Int_t& out ) const;
191 : AliTPCPIDResponse::EChamberStatus TrackStatus(const AliVTrack* track, Int_t layer) const;
192 : Float_t MaxClusterRadius(const AliVTrack* track) const;
193 : Bool_t TrackApex(const AliVTrack* track, Float_t magField, Double_t position[3]) const;
194 0 : static const char* GainScenarioName(Int_t n) {return fgkGainScenarioName[(n>fgkNumberOfGainScenarios)?fgkNumberOfGainScenarios:n];}
195 : Int_t ResponseFunctionIndex( AliPID::EParticleType species,
196 : ETPCgainScenario gainScenario ) const;
197 : void ResetSplines();
198 :
199 : //OLD
200 : Double_t GetExpectedSignal(Float_t mom,
201 : AliPID::EParticleType n=AliPID::kKaon) const;
202 : Double_t GetExpectedSigma(Float_t mom, Int_t nPoints,
203 : AliPID::EParticleType n=AliPID::kKaon) const;
204 : Float_t GetNumberOfSigmas(Float_t mom,
205 : Float_t dEdx,
206 : Int_t nPoints,
207 : AliPID::EParticleType n=AliPID::kKaon) const {
208 : //
209 : // Deprecated function (for backward compatibility). Please use
210 : // GetNumberOfSigmas(const AliVTrack *track, AliPID::EParticleType species, ETPCdEdxSource dedxSource,
211 : // Bool_t correctEta, Bool_t correctMultiplicity)
212 : // instead!TODO
213 : //
214 :
215 0 : Double_t bethe=GetExpectedSignal(mom,n);
216 0 : Double_t sigma=GetExpectedSigma(mom,nPoints,n);
217 0 : return (dEdx-bethe)/sigma;
218 : }
219 :
220 0 : Double_t GetMIP() const { return fMIP;}
221 0 : Float_t GetRes0() const { return fRes0[0]; }
222 0 : Float_t GetResN2() const { return fResN2[0]; }
223 0 : Float_t GetRes0(ETPCgainScenario s) const { return fRes0[s]; }
224 0 : Float_t GetResN2(ETPCgainScenario s) const { return fResN2[s]; }
225 :
226 : Bool_t RegisterSpline(const char * name, Int_t index);
227 : Double_t EvaldEdxSpline(Double_t bg,Int_t entry);
228 0 : static Double_t SEvaldEdx(Double_t bg,Int_t entry){ return (fgInstance!=0)? fgInstance->EvaldEdxSpline(bg,entry):0;};
229 :
230 : //===| dEdx type functions |==================================================
231 : void SetdEdxType(ETPCdEdxType dEdxType, Int_t dEdxChargeType=0, Int_t dEdxWeightType=0, Double_t dEdxIROCweight=1., Double_t dEdxOROCmedWeight=1., Double_t dEdxOROClongWeight=1.) {
232 0 : fdEdxType=dEdxType; fdEdxChargeType=dEdxChargeType; fdEdxWeightType=dEdxWeightType; fIROCweight=dEdxIROCweight; fOROCmedWeight=dEdxOROCmedWeight; fOROClongWeight=dEdxOROClongWeight; }
233 0 : ETPCdEdxType GetdEdxType() const { return fdEdxType; }
234 0 : Int_t GetdEdxChargeType() const { return fdEdxChargeType; }
235 0 : Int_t GetdEdxWeightType() const { return fdEdxWeightType; }
236 0 : Double_t GetIROCweight() const { return fIROCweight; }
237 0 : Double_t GetOROCmedWeight() const { return fOROCmedWeight; }
238 0 : Double_t GetOROClongWeight() const { return fOROClongWeight; }
239 :
240 : Double_t GetTrackdEdx(const AliVTrack* track) const;
241 :
242 : //===| Initialisation |=======================================================
243 : Bool_t InitFromOADB(const Int_t run, const char* pass, const char* oadbFile="$ALICE_PHYSICS/OADB/COMMON/PID/data/TPCPIDResponseOADB.root", Bool_t initMultiplicityCorrection=kTRUE);
244 :
245 : Bool_t SetSplinesFromArray (const TObjArray* arrSplines);
246 : Bool_t SetMultiplicityCorrectionFromString(const TString& multiplicityData);
247 : Bool_t SetdEdxTypeFromString (const TString& dEdxTypeSet);
248 : Bool_t SetdEdxResolutionFromString (const TString& dEdxTypeSet);
249 :
250 : //===| Helpers |==============================================================
251 : static TString GetChecksum(const TObject* obj);
252 : static TObjArray* GetMultiplicityCorrectionArrayFromString(const TString& corrections);
253 :
254 : protected:
255 : Double_t GetExpectedSignal(const AliVTrack* track,
256 : AliPID::EParticleType species,
257 : Double_t dEdx,
258 : const TSpline3* responseFunction,
259 : Bool_t correctEta,
260 : Bool_t correctMultiplicity) const;
261 :
262 : Double_t GetExpectedSigma(const AliVTrack* track,
263 : AliPID::EParticleType species,
264 : ETPCgainScenario gainScenario,
265 : Double_t dEdx,
266 : Int_t nPoints,
267 : const TSpline3* responseFunction,
268 : Bool_t correctEta,
269 : Bool_t correctMultiplicity) const;
270 : //
271 : // function for numberical debugging 0 registed splines can be used in the TFormula and tree visualizations
272 : //
273 : private:
274 : Float_t fMIP; // dEdx for MIP
275 : Float_t fRes0[fgkNumberOfGainScenarios]; // relative dEdx resolution rel sigma = fRes0*sqrt(1+fResN2/npoint)
276 : Float_t fResN2[fgkNumberOfGainScenarios]; // relative Npoint dependence rel sigma = fRes0*sqrt(1+fResN2/npoint)
277 :
278 : Double_t fKp1; // Parameters
279 : Double_t fKp2; // of
280 : Double_t fKp3; // the ALEPH
281 : Double_t fKp4; // Bethe-Bloch
282 : Double_t fKp5; // formula
283 :
284 : Bool_t fUseDatabase; // flag if fine-tuned database-response or simple ALEPH BB should be used
285 :
286 : TObjArray fResponseFunctions; //! ObjArray of response functions individually for each particle
287 : AliOADBContainer* fOADBContainer; //! OADB container with response functions
288 : TVectorF fVoltageMap; //!stores a map of voltages wrt nominal for all chambers
289 : Float_t fLowGainIROCthreshold; //voltage threshold below which the IROC is considered low gain
290 : Float_t fBadIROCthreshhold; //voltage threshold for bad IROCS
291 : Float_t fLowGainOROCthreshold; //voltage threshold below which the OROC is considered low gain
292 : Float_t fBadOROCthreshhold; //voltage threshold for bad OROCS
293 : Float_t fMaxBadLengthFraction; //the maximum allowed fraction of track length in a bad sector.
294 :
295 : Int_t sectorNumber(Double_t phi) const;
296 :
297 : Double_t fMagField; //! Magnetic field
298 :
299 : static const char* fgkGainScenarioName[fgkNumberOfGainScenarios+1];
300 :
301 : TH2D* fhEtaCorr; //! Map for TPC eta correction
302 : TH2D* fhEtaSigmaPar1; //! Map for parameter 1 of the IROCdEdx sigma parametrisation
303 :
304 : Double_t fSigmaPar0; // Parameter 0 of the dEdx sigma parametrisation
305 :
306 : Int_t fCurrentEventMultiplicity; // Multiplicity of the current event
307 : TF1* fCorrFuncMultiplicity; //! Function to correct for the multiplicity dependence of the TPC dEdx
308 : TF1* fCorrFuncMultiplicityTanTheta; //! Function to correct the additional tanTheta dependence of the multiplicity dependence of the TPC dEdx
309 : TF1* fCorrFuncSigmaMultiplicity; //! Function to correct for the multiplicity dependence of the TPC dEdx resolution
310 :
311 : // dEdx type information
312 : ETPCdEdxType fdEdxType; // source for dEdx information to use
313 : Int_t fdEdxChargeType; // charge type to use for dEdx calculation from AliTPCdEdxInfo
314 : Int_t fdEdxWeightType; // tracklet weighting type to use for dEdx calculation from AliTPCdEdxInfo
315 : Double_t fIROCweight; // IROC weight to use for dEdx calculation from AliTPCdEdxInfo
316 : Double_t fOROCmedWeight; // OROC medium pad size weight to use for dEdx calculation from AliTPCdEdxInfo
317 : Double_t fOROClongWeight; // OROC long pad size weight to use for dEdx calculation from AliTPCdEdxInfo
318 :
319 : //
320 : //
321 : static AliTPCPIDResponse* fgInstance; //! Instance of this class (singleton implementation)
322 : TObjArray fSplineArray; //array of registered splines
323 176 : ClassDef(AliTPCPIDResponse,6) // TPC PID class
324 : };
325 :
326 :
327 : #endif
328 :
329 :
|
