- introduction of gain scenarios (e.g. OROC only - for homogeneous gain in 11h)
[u/mrichter/AliRoot.git] / STEER / STEERBase / AliTPCPIDResponse.h
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 //-------------------------------------------------------
14 #include <Rtypes.h>
15
16 #include <TNamed.h>
17 #include <TVectorF.h>
18 #include <TObjArray.h>
19
20 #include "AliPID.h"
21
22 class AliVTrack;
23 class TSpline3;
24
25 class AliTPCPIDResponse: public TNamed {
26 public:
27   AliTPCPIDResponse();
28   AliTPCPIDResponse(const Double_t *param);
29   AliTPCPIDResponse(const AliTPCPIDResponse&);
30   AliTPCPIDResponse& operator=(const AliTPCPIDResponse&);
31   virtual ~AliTPCPIDResponse() {}
32
33   enum EChamberStatus {
34     kChamberOff=0,
35     kChamberHighGain=1,
36     kChamberLowGain=2,
37     kChamberInvalid=3
38   };
39   
40   enum ETPCgainScenario {
41     kDefault= 0,
42     kALLhigh = 1,
43     kOROChigh = 2,
44     kGainScenarioInvalid = 3
45   };
46
47   static const Int_t fgkNumberOfParticleSpecies=AliPID::kSPECIESC;
48   static const Int_t fgkNumberOfGainScenarios=3;
49   static const Int_t fgkNumberOfdEdxSourceScenarios=3;
50
51   enum ETPCdEdxSource {
52     kdEdxDefault=0,        // use combined dEdx from IROC+OROC (assumes ideal detector)
53     kdEdxOROC=1,       // use only OROC
54     kdEdxHybrid=2,   // Use IROC+OROC dEdx only where IROCS are good (high gain), otherwise fall back to OROC only
55     kdEdxInvalid=3     //invalid
56   };
57
58   void SetSigma(Float_t res0, Float_t resN2);
59   void SetBetheBlochParameters(Double_t kp1,
60                                Double_t kp2,
61                                Double_t kp3,
62                                Double_t kp4,
63                                Double_t kp5
64                                );
65   void SetMip(Float_t mip) { fMIP = mip; } // Set overall normalisation; mean dE/dx for MIP
66   Double_t Bethe(Double_t bg) const;
67   void SetUseDatabase(Bool_t useDatabase) { fUseDatabase = useDatabase;}
68   Bool_t GetUseDatabase() const { return fUseDatabase;}
69   
70   void SetResponseFunction(AliPID::EParticleType type, TObject * const o) { fResponseFunctions.AddAt(o,(Int_t)type); }
71   const TObject * GetResponseFunction(AliPID::EParticleType type) { return fResponseFunctions.At((Int_t)type); }
72   void SetVoltage(Int_t n, Float_t v) {fVoltageMap[n]=v;}
73   void SetVoltageMap(const TVectorF& a) {fVoltageMap=a;} //resets ownership, ~ will not delete contents
74   Float_t GetVoltage(Int_t n) const {return fVoltageMap[n];}
75   void SetLowGainIROCthreshold(Float_t v) {fLowGainIROCthreshold=v;}
76   void SetBadIROCthreshold(Float_t v) {fBadIROCthreshhold=v;}
77   void SetLowGainOROCthreshold(Float_t v) {fLowGainOROCthreshold=v;}
78   void SetBadOROCthreshold(Float_t v) {fBadOROCthreshhold=v;}
79   void SetMaxBadLengthFraction(Float_t f) {fMaxBadLengthFraction=f;}
80
81   void SetMagField(Double_t mf) { fMagField=mf; }
82   
83   //NEW
84   void SetSigma(Float_t res0, Float_t resN2, ETPCgainScenario gainScenario );
85   Double_t GetExpectedSignal( Double_t momentum,
86                               AliPID::EParticleType species,
87                               const TSpline3* responseFunction ) const;
88   Double_t GetExpectedSignal( const AliVTrack* track,
89                               AliPID::EParticleType species,
90                               ETPCdEdxSource dedxSource );
91   Double_t GetExpectedSigma( const AliVTrack* track, 
92                              AliPID::EParticleType species,
93                              ETPCdEdxSource dedxSource );
94   Double_t GetExpectedSigma( Double_t mom,
95                              Int_t nPoints,
96                              AliPID::EParticleType species,
97                              ETPCgainScenario gainScenario,
98                              const TSpline3* responseFunction) const;
99   Float_t GetNumberOfSigmas( const AliVTrack* track,
100                              AliPID::EParticleType species,
101                              ETPCdEdxSource dedxSource );
102
103   void SetResponseFunction(TObject* o,
104                            AliPID::EParticleType type,
105                            ETPCgainScenario gainScenario);
106   void Print(Option_t* option="") const;
107   TSpline3* GetResponseFunction( AliPID::EParticleType species,
108                                  ETPCgainScenario gainScenario ) const;
109   TSpline3* GetResponseFunction( const AliVTrack* track,
110                                  AliPID::EParticleType species,
111                                  ETPCdEdxSource dedxSource );
112   Bool_t ResponseFunctiondEdxN(const AliVTrack* track, 
113                                AliPID::EParticleType species,
114                                ETPCdEdxSource dedxSource);
115   Bool_t sectorNumbersInOut(const AliVTrack* track, 
116                             Double_t innerRadius, Double_t outerRadius, 
117                             Float_t& phiIn, Float_t& phiOut, 
118                             Int_t& in, Int_t& out ) const;
119   AliTPCPIDResponse::EChamberStatus TrackStatus(const AliVTrack* track, Int_t layer) const;
120   Float_t MaxClusterRadius(const AliVTrack* track) const;
121   Bool_t TrackApex(const AliVTrack* track, Float_t magField, Double_t position[3]) const;
122   static const char* GainScenarioName(Int_t n) {return fgkGainScenarioName[(n>fgkNumberOfdEdxSourceScenarios)?fgkNumberOfdEdxSourceScenarios+1:n];}
123   Int_t ResponseFunctionIndex( AliPID::EParticleType species,
124                                ETPCgainScenario gainScenario ) const;
125   void ResetSplines();
126
127   void InvalidateCurrentValues();
128   TSpline3* GetCurrentResponseFunction() const {return fCurrentResponseFunction;}
129   Double_t GetCurrentdEdx() const {return fCurrentdEdx;}
130   Int_t GetCurrentNPoints() const {return fCurrentNPoints;}
131   ETPCgainScenario GetCurrentGainScenario() const {return fCurrentGainScenario;}
132
133   //OLD
134   Double_t GetExpectedSignal(const Float_t mom,
135                      AliPID::EParticleType n=AliPID::kKaon) const;
136   Double_t GetExpectedSigma(const Float_t mom, const Int_t nPoints,
137                             AliPID::EParticleType n=AliPID::kKaon) const;
138   Float_t  GetNumberOfSigmas(const Float_t mom, 
139                              const Float_t dEdx, 
140                                                const Int_t nPoints,
141                              AliPID::EParticleType n=AliPID::kKaon) const {
142
143     Double_t bethe=GetExpectedSignal(mom,n);
144     Double_t sigma=GetExpectedSigma(mom,nPoints,n);
145     return (dEdx-bethe)/sigma;
146   }
147
148   Double_t GetMIP() const { return fMIP;} 
149   Float_t  GetRes0()  const { return fRes0[0];  }
150   Float_t  GetResN2() const { return fResN2[0]; }
151   Float_t  GetRes0(ETPCgainScenario s)  const { return fRes0[s];  }
152   Float_t  GetResN2(ETPCgainScenario s) const { return fResN2[s]; }
153
154 private:
155   Float_t fMIP;          // dEdx for MIP
156   Float_t fRes0[fgkNumberOfGainScenarios];  // relative dEdx resolution  rel sigma = fRes0*sqrt(1+fResN2/npoint)
157   Float_t fResN2[fgkNumberOfGainScenarios]; // relative Npoint dependence rel  sigma = fRes0*sqrt(1+fResN2/npoint)
158
159   Double_t fKp1;   // Parameters
160   Double_t fKp2;   //    of
161   Double_t fKp3;   // the ALEPH
162   Double_t fKp4;   // Bethe-Bloch
163   Double_t fKp5;   // formula
164
165   Bool_t fUseDatabase; // flag if fine-tuned database-response or simple ALEPH BB should be used
166   
167   TObjArray fResponseFunctions; //! ObjArray of response functions individually for each particle
168   TVectorF fVoltageMap; //!stores a map of voltages wrt nominal for all chambers
169   Float_t fLowGainIROCthreshold;  //voltage threshold below which the IROC is considered low gain
170   Float_t fBadIROCthreshhold;     //voltage threshold for bad IROCS
171   Float_t fLowGainOROCthreshold;  //voltage threshold below which the OROC is considered low gain
172   Float_t fBadOROCthreshhold;     //voltage threshold for bad OROCS
173   Float_t fMaxBadLengthFraction;  //the maximum allowed fraction of track length in a bad sector.
174
175   TSpline3* fCurrentResponseFunction;      //!response function for current track
176   Double_t fCurrentdEdx;                  //!dEdx for currently processed track
177   Int_t fCurrentNPoints;                  //!number of points used for dEdx calculation for current track
178   ETPCgainScenario fCurrentGainScenario;  //!gain scenario used for current track
179   Int_t sectorNumber(Double_t phi) const;
180
181   Double_t fMagField;  //! Magnetic field
182
183   static const char* fgkGainScenarioName[fgkNumberOfGainScenarios+1];
184
185   ClassDef(AliTPCPIDResponse,4)   // TPC PID class
186 };
187
188 #endif
189
190