]> git.uio.no Git - u/mrichter/AliRoot.git/blame - STEER/STEERBase/AliTOFPIDResponse.cxx
Merge branch 'master' into TPCdev
[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 //
2ec09a5e 185 Double_t expt[AliPID::kSPECIESC];
fc9b31a7 186 track->GetIntegratedTimes(expt,AliPID::kSPECIESC);
67376d1d 187 if (type<=AliPID::kProton) return expt[type];
188 else {
2ec09a5e
RS
189 if (expt[type]<1.E-1) {
190 Double_t p = track->P();
191 Double_t massZ = AliPID::ParticleMassZ(type);
192 return expt[0]/p*massZ*TMath::Sqrt(1.+p*p/massZ/massZ);
193 } else return expt[type];
67376d1d 194 }
195}
196//_________________________________________________________________________
6c68754e 197Int_t AliTOFPIDResponse::GetMomBin(Float_t p) const{
f858b00e 198 //
199 // Returns the momentum bin index
200 //
201
6c68754e 202 Int_t i=0;
203 while(p > fPCutMin[i] && i < fNmomBins) i++;
204 if(i > 0) i--;
10d100d4 205
6c68754e 206 return i;
207}
208//_________________________________________________________________________
209void AliTOFPIDResponse::SetMomBoundary(){
f858b00e 210 //
211 // Set boundaries for momentum bins
212 //
213
6c68754e 214 fPCutMin[0] = 0.3;
215 fPCutMin[1] = 0.5;
216 fPCutMin[2] = 0.6;
217 fPCutMin[3] = 0.7;
218 fPCutMin[4] = 0.8;
219 fPCutMin[5] = 0.9;
220 fPCutMin[6] = 1;
221 fPCutMin[7] = 1.2;
222 fPCutMin[8] = 1.5;
223 fPCutMin[9] = 2;
224 fPCutMin[10] = 3;
225}
f858b00e 226//_________________________________________________________________________
7170298c 227Float_t AliTOFPIDResponse::GetStartTime(Float_t mom) const {
f858b00e 228 //
229 // Returns event_time value as estimated by TOF combinatorial algorithm
230 //
231
232 Int_t ibin = GetMomBin(mom);
233 return GetT0bin(ibin);
234
235}
236//_________________________________________________________________________
7170298c 237Float_t AliTOFPIDResponse::GetStartTimeRes(Float_t mom) const {
f858b00e 238 //
239 // Returns event_time resolution as estimated by TOF combinatorial algorithm
240 //
241
242 Int_t ibin = GetMomBin(mom);
243 return GetT0binRes(ibin);
244
245}
d4d15b21 246//_________________________________________________________________________
247Int_t AliTOFPIDResponse::GetStartTimeMask(Float_t mom) const {
248 //
249 // Returns event_time mask
250 //
251
252 Int_t ibin = GetMomBin(mom);
253 return GetT0binMask(ibin);
254
255}
a2c30af1 256//_________________________________________________________________________
c53e310b 257Double_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 258{
c53e310b 259
260 // To add mismatch
261 Float_t mismAdd = addmism*0.01;
262 if(pt>1.0) mismAdd /= pt;
263
264 if(mismAdd > 0.01){ // apply additional mismatch
265 if(gRandom->Rndm() < mismAdd){
266 return GetMismatchRandomValue(eta)-time;
267 }
268 }
269
a2c30af1 270 if(fTOFtailResponse)
271 return fTOFtailResponse->GetRandom();
272 else
273 return 0.0;
274}
3cb11031 275//_________________________________________________________________________
276Double_t AliTOFPIDResponse::GetMismatchRandomValue(Float_t eta) // generate a random value for mismatched tracks (for MC analyses)
277{
278 if(!fHmismTOF){
279 TFile *fmism = new TFile("$ALICE_ROOT/TOF/data/TOFmismatchDistr.root");
280 if(fmism) fHmismTOF = (TH1F *) fmism->Get("TOFmismDistr");
281 if(!fHmismTOF){
282 printf("I cannot retrive TOF mismatch histos... skipped!");
283 return -10000.;
284 }
4fc0f532 285 fHmismTOF->Scale(TMath::Sqrt(2*TMath::Pi())/(fHmismTOF->Integral(1,fHmismTOF->GetNbinsX()) * fHmismTOF->GetBinWidth(1)));
3cb11031 286
287 TFile *fchDist = new TFile("$ALICE_ROOT/TOF/data/TOFchannelDist.root");
288 if(fchDist) fHchannelTOFdistr = (TH1D *) fchDist->Get("hTOFchanDist");
289 if(!fHchannelTOFdistr){
290 printf("I cannot retrive TOF channel distance distribution... skipped!");
291 return -10000.;
292 }
293 }
294
295 Float_t etaAbs = TMath::Abs(eta);
296 Int_t channel = Int_t(4334.09 - 4758.36 * etaAbs -1989.71 * etaAbs*etaAbs + 1957.62*etaAbs*etaAbs*etaAbs);
297 if(channel < 1 || etaAbs > 1) channel = 1;
298 Float_t distIP = fHchannelTOFdistr->GetBinContent(channel);
299
300 return fHmismTOF->GetRandom() + distIP*3.35655419905265973e+01;
301}
c53e310b 302//_________________________________________________________________________
303Int_t AliTOFPIDResponse::GetTOFchannel(AliVParticle *trk) const{
304 Float_t etaAbs = TMath::Abs(trk->Eta());
305 Int_t channel = Int_t(4334.09 - 4758.36 * etaAbs -1989.71 * etaAbs*etaAbs + 1957.62*etaAbs*etaAbs*etaAbs);
306 if(channel < 1 || etaAbs > 1) channel = 1;
307
308 return channel;
309}