]> git.uio.no Git - u/mrichter/AliRoot.git/blame - STEER/STEERBase/AliTOFPIDResponse.cxx
Fix for Coverity defect 21184.
[u/mrichter/AliRoot.git] / STEER / STEERBase / AliTOFPIDResponse.cxx
CommitLineData
10d100d4 1/**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
6 * *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
15
16//-----------------------------------------------------------------//
17// //
18// Implementation of the TOF PID class //
19// Origin: Iouri Belikov, CERN, Jouri.Belikov@cern.ch //
20// //
21//-----------------------------------------------------------------//
22
23#include "TMath.h"
24#include "AliLog.h"
a2c30af1 25#include "TF1.h"
3cb11031 26#include "TH1F.h"
27#include "TH1D.h"
28#include "TFile.h"
c53e310b 29#include "TRandom.h"
10d100d4 30
31#include "AliTOFPIDResponse.h"
32
33ClassImp(AliTOFPIDResponse)
34
a2c30af1 35TF1 *AliTOFPIDResponse::fTOFtailResponse = NULL; // function to generate a TOF tail
3cb11031 36TH1F *AliTOFPIDResponse::fHmismTOF = NULL; // TOF mismatch distribution
37TH1D *AliTOFPIDResponse::fHchannelTOFdistr=NULL; // TOF channel distance distribution
a2c30af1 38
10d100d4 39//_________________________________________________________________________
40AliTOFPIDResponse::AliTOFPIDResponse():
41 fSigma(0),
42 fPmax(0), // zero at 0.5 GeV/c for pp
43 fTime0(0)
44{
679de35e 45 fPar[0] = 0.008;
46 fPar[1] = 0.008;
47 fPar[2] = 0.002;
48 fPar[3] = 40.0;
95ad1018 49
a2c30af1 50 if(!fTOFtailResponse){
c5fb644a 51 fTOFtailResponse = new TF1("fTOFtail","[0]*TMath::Exp(-(x-[1])*(x-[1])/2/[2]/[2])* (x < [1]+[3]*[2]) + (x > [1]+[3]*[2])*[0]*TMath::Exp(-(x-[1]-[3]*[2]*0.5)*[3]/[2] * 0.0111)*0.018",-1000,1000);
a2c30af1 52 fTOFtailResponse->SetParameter(0,1);
c5fb644a 53 fTOFtailResponse->SetParameter(1,-26);
a2c30af1 54 fTOFtailResponse->SetParameter(2,1);
c5fb644a 55 fTOFtailResponse->SetParameter(3,0.89);
a2c30af1 56 fTOFtailResponse->SetNpx(10000);
57 }
58
59
6c68754e 60 // Reset T0 info
61 ResetT0info();
62 SetMomBoundary();
10d100d4 63}
64//_________________________________________________________________________
65AliTOFPIDResponse::AliTOFPIDResponse(Double_t *param):
66 fSigma(param[0]),
67 fPmax(0), // zero at 0.5 GeV/c for pp
68 fTime0(0)
69{
70 //
71 // The main constructor
72 //
73 //
74
75 //fPmax=TMath::Exp(-0.5*3*3)/fSigma; // ~3 sigma at 0.5 GeV/c for PbPb
6c68754e 76
679de35e 77 fPar[0] = 0.008;
78 fPar[1] = 0.008;
79 fPar[2] = 0.002;
80 fPar[3] = 40.0;
95ad1018 81
c5fb644a 82 if(!fTOFtailResponse){
83 fTOFtailResponse = new TF1("fTOFtail","[0]*TMath::Exp(-(x-[1])*(x-[1])/2/[2]/[2])* (x < [1]+[3]*[2]) + (x > [1]+[3]*[2])*[0]*TMath::Exp(-(x-[1]-[3]*[2]*0.5)*[3]/[2] * 0.0111)*0.018",-1000,1000);
84 fTOFtailResponse->SetParameter(0,1);
85 fTOFtailResponse->SetParameter(1,-26);
86 fTOFtailResponse->SetParameter(2,1);
87 fTOFtailResponse->SetParameter(3,0.89);
88 fTOFtailResponse->SetNpx(10000);
89 }
90
6c68754e 91 // Reset T0 info
92 ResetT0info();
93 SetMomBoundary();
10d100d4 94}
95//_________________________________________________________________________
b3f687a1 96void AliTOFPIDResponse::SetTOFtail(Float_t tail){
97 if(!fTOFtailResponse){
98 fTOFtailResponse = new TF1("fTOFtail","[0]*TMath::Exp(-(x-[1])*(x-[1])/2/[2]/[2])* (x < [1]+[3]*[2]) + (x > [1]+[3]*[2])*[0]*TMath::Exp(-(x-[1]-[3]*[2]*0.5)*[3]/[2] * 0.0111)*0.018",-1000,1000);
99 fTOFtailResponse->SetParameter(0,1);
100 fTOFtailResponse->SetParameter(1,-26);
101 fTOFtailResponse->SetParameter(2,1);
102 fTOFtailResponse->SetParameter(3,tail);
103 fTOFtailResponse->SetNpx(10000);
104 }
105 else{
106 fTOFtailResponse->SetParameter(3,tail);
107 }
108}
109//_________________________________________________________________________
10d100d4 110Double_t
4fc0f532 111AliTOFPIDResponse::GetMismatchProbability(Double_t time,Double_t eta) const {
112 if(!fHmismTOF){
113 TFile *fmism = new TFile("$ALICE_ROOT/TOF/data/TOFmismatchDistr.root");
114 if(fmism) fHmismTOF = (TH1F *) fmism->Get("TOFmismDistr");
115 if(!fHmismTOF){
116 printf("I cannot retrive TOF mismatch histos... skipped!");
117 return 1E-4;
118 }
119 fHmismTOF->Scale(TMath::Sqrt(2*TMath::Pi())/(fHmismTOF->Integral(1,fHmismTOF->GetNbinsX()) * fHmismTOF->GetBinWidth(1)));
120
121 TFile *fchDist = new TFile("$ALICE_ROOT/TOF/data/TOFchannelDist.root");
122 if(fchDist) fHchannelTOFdistr = (TH1D *) fchDist->Get("hTOFchanDist");
123 if(!fHchannelTOFdistr){
124 printf("I cannot retrive TOF channel distance distribution... skipped!");
125 return 1E-4;
126 }
127 }
10d100d4 128
4fc0f532 129 Float_t etaAbs = TMath::Abs(eta);
130 Int_t channel = Int_t(4334.09 - 4758.36 * etaAbs -1989.71 * etaAbs*etaAbs + 1957.62*etaAbs*etaAbs*etaAbs);
131 if(channel < 1 || etaAbs > 1) channel = 1;
132 Float_t distIP = fHchannelTOFdistr->GetBinContent(channel);
133
134 Double_t mismWeight = fHmismTOF->Interpolate(time - distIP*3.35655419905265973e+01);
10d100d4 135
4fc0f532 136 return mismWeight;
10d100d4 137}
138//_________________________________________________________________________
139Double_t AliTOFPIDResponse::GetExpectedSigma(Float_t mom, Float_t time, Float_t mass) const {
140 //
141 // Return the expected sigma of the PID signal for the specified
67376d1d 142 // particle mass/Z.
10d100d4 143 // If the operation is not possible, return a negative value.
144 //
145
95ad1018 146 Double_t dpp=fPar[0] + fPar[1]*mom + fPar[2]*mass/mom; //mean relative pt resolution;
10d100d4 147
148
149 Double_t sigma = dpp*time/(1.+ mom*mom/(mass*mass));
6c68754e 150
151 Int_t index = GetMomBin(mom);
152
153 Double_t t0res = fT0resolution[index];
154
95ad1018 155 return TMath::Sqrt(sigma*sigma + fPar[3]*fPar[3]/mom/mom + fSigma*fSigma + t0res*t0res);
10d100d4 156
10d100d4 157}
6c68754e 158//_________________________________________________________________________
67376d1d 159Double_t AliTOFPIDResponse::GetExpectedSigma(Float_t mom, Float_t time, AliPID::EParticleType type) const {
160 //
161 // Return the expected sigma of the PID signal for the specified
162 // particle type.
163 // If the operation is not possible, return a negative value.
164 //
165
166 Double_t mass = AliPID::ParticleMassZ(type);
167 Double_t dpp=fPar[0] + fPar[1]*mom + fPar[2]*mass/mom; //mean relative pt resolution;
168
169
170 Double_t sigma = dpp*time/(1.+ mom*mom/(mass*mass));
171
172 Int_t index = GetMomBin(mom);
173
174 Double_t t0res = fT0resolution[index];
175
176 return TMath::Sqrt(sigma*sigma + fPar[3]*fPar[3]/mom/mom + fSigma*fSigma + t0res*t0res);
177
178}
179//_________________________________________________________________________
180Double_t AliTOFPIDResponse::GetExpectedSignal(const AliVTrack* track,AliPID::EParticleType type) const {
181 //
182 // Return the expected signal of the PID signal for the particle type
183 // If the operation is not possible, return a negative value.
184 //
185 Double_t expt[5];
186 track->GetIntegratedTimes(expt);
187 if (type<=AliPID::kProton) return expt[type];
188 else {
189 Double_t p = track->P();
190 Double_t massZ = AliPID::ParticleMassZ(type);
191 return expt[0]/p*massZ*TMath::Sqrt(1.+p*p/massZ/massZ);
192 }
193}
194//_________________________________________________________________________
6c68754e 195Int_t AliTOFPIDResponse::GetMomBin(Float_t p) const{
f858b00e 196 //
197 // Returns the momentum bin index
198 //
199
6c68754e 200 Int_t i=0;
201 while(p > fPCutMin[i] && i < fNmomBins) i++;
202 if(i > 0) i--;
10d100d4 203
6c68754e 204 return i;
205}
206//_________________________________________________________________________
207void AliTOFPIDResponse::SetMomBoundary(){
f858b00e 208 //
209 // Set boundaries for momentum bins
210 //
211
6c68754e 212 fPCutMin[0] = 0.3;
213 fPCutMin[1] = 0.5;
214 fPCutMin[2] = 0.6;
215 fPCutMin[3] = 0.7;
216 fPCutMin[4] = 0.8;
217 fPCutMin[5] = 0.9;
218 fPCutMin[6] = 1;
219 fPCutMin[7] = 1.2;
220 fPCutMin[8] = 1.5;
221 fPCutMin[9] = 2;
222 fPCutMin[10] = 3;
223}
f858b00e 224//_________________________________________________________________________
7170298c 225Float_t AliTOFPIDResponse::GetStartTime(Float_t mom) const {
f858b00e 226 //
227 // Returns event_time value as estimated by TOF combinatorial algorithm
228 //
229
230 Int_t ibin = GetMomBin(mom);
231 return GetT0bin(ibin);
232
233}
234//_________________________________________________________________________
7170298c 235Float_t AliTOFPIDResponse::GetStartTimeRes(Float_t mom) const {
f858b00e 236 //
237 // Returns event_time resolution as estimated by TOF combinatorial algorithm
238 //
239
240 Int_t ibin = GetMomBin(mom);
241 return GetT0binRes(ibin);
242
243}
d4d15b21 244//_________________________________________________________________________
245Int_t AliTOFPIDResponse::GetStartTimeMask(Float_t mom) const {
246 //
247 // Returns event_time mask
248 //
249
250 Int_t ibin = GetMomBin(mom);
251 return GetT0binMask(ibin);
252
253}
a2c30af1 254//_________________________________________________________________________
c53e310b 255Double_t AliTOFPIDResponse::GetTailRandomValue(Float_t pt,Float_t eta,Float_t time,Float_t addmism) // generate a random value to add a tail to TOF time (for MC analyses)
a2c30af1 256{
c53e310b 257
258 // To add mismatch
259 Float_t mismAdd = addmism*0.01;
260 if(pt>1.0) mismAdd /= pt;
261
262 if(mismAdd > 0.01){ // apply additional mismatch
263 if(gRandom->Rndm() < mismAdd){
264 return GetMismatchRandomValue(eta)-time;
265 }
266 }
267
a2c30af1 268 if(fTOFtailResponse)
269 return fTOFtailResponse->GetRandom();
270 else
271 return 0.0;
272}
3cb11031 273//_________________________________________________________________________
274Double_t AliTOFPIDResponse::GetMismatchRandomValue(Float_t eta) // generate a random value for mismatched tracks (for MC analyses)
275{
276 if(!fHmismTOF){
277 TFile *fmism = new TFile("$ALICE_ROOT/TOF/data/TOFmismatchDistr.root");
278 if(fmism) fHmismTOF = (TH1F *) fmism->Get("TOFmismDistr");
279 if(!fHmismTOF){
280 printf("I cannot retrive TOF mismatch histos... skipped!");
281 return -10000.;
282 }
4fc0f532 283 fHmismTOF->Scale(TMath::Sqrt(2*TMath::Pi())/(fHmismTOF->Integral(1,fHmismTOF->GetNbinsX()) * fHmismTOF->GetBinWidth(1)));
3cb11031 284
285 TFile *fchDist = new TFile("$ALICE_ROOT/TOF/data/TOFchannelDist.root");
286 if(fchDist) fHchannelTOFdistr = (TH1D *) fchDist->Get("hTOFchanDist");
287 if(!fHchannelTOFdistr){
288 printf("I cannot retrive TOF channel distance distribution... skipped!");
289 return -10000.;
290 }
291 }
292
293 Float_t etaAbs = TMath::Abs(eta);
294 Int_t channel = Int_t(4334.09 - 4758.36 * etaAbs -1989.71 * etaAbs*etaAbs + 1957.62*etaAbs*etaAbs*etaAbs);
295 if(channel < 1 || etaAbs > 1) channel = 1;
296 Float_t distIP = fHchannelTOFdistr->GetBinContent(channel);
297
298 return fHmismTOF->GetRandom() + distIP*3.35655419905265973e+01;
299}
c53e310b 300//_________________________________________________________________________
301Int_t AliTOFPIDResponse::GetTOFchannel(AliVParticle *trk) const{
302 Float_t etaAbs = TMath::Abs(trk->Eta());
303 Int_t channel = Int_t(4334.09 - 4758.36 * etaAbs -1989.71 * etaAbs*etaAbs + 1957.62*etaAbs*etaAbs*etaAbs);
304 if(channel < 1 || etaAbs > 1) channel = 1;
305
306 return channel;
307}