Major fixes in the digits and raw-data reco. Now both are equivalent. To be propagate...
[u/mrichter/AliRoot.git] / TOF / AliTOFT0maker.cxx
CommitLineData
8f589502 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 **************************************************************************/
6a4e212e 15/* $Id: AliTOFT0maker.cxx,v 1.8 2010/01/19 16:32:20 noferini Exp $ */
8f589502 16
17/////////////////////////////////////////////////////////////////////////////
18// //
19// This class contains the basic functions for the time zero //
6a4e212e 20// evaluation with TOF detector informations. //
21// Use case in an analysis task: //
22// //
23// Create the object in the task constructor (fTOFmaker is a private var) //
24// fTOFmaker = new AliTOFT0maker(); //
03bd764d 25// fTOFmaker->SetTimeResolution(130.0); // if you want set the TOF res //
6a4e212e 26// 115 ps is the TOF default resolution value //
27// //
28// Use the RemakePID method in the task::Exec //
29// Double_t* calcolot0; //
30// calcolot0=fTOFmaker->RemakePID(fESD); //
31// //calcolot0[0] = calculated event time //
32// //calcolot0[1] = event time time resolution //
33// //calcolot0[2] = average event time for the current fill //
03bd764d 34// //calcolot0[3] = tracks at TOF //
35// //calcolot0[4] = calculated event time (only TOF) //
36// //calcolot0[5] = event time time resolution (only TOF) //
37// //calcolot0[6] = sigma t0 fill //
38// //calcolot0[7] = tracks at TOF really used in tht algorithm //
6a4e212e 39// //
40// Let consider that: //
41// - the PIF is automatically recalculated with the event time subtrction //
8f589502 42// //
43/////////////////////////////////////////////////////////////////////////////
44
536031f2 45#include "AliTOFT0v1.h"
46#include "AliTOFT0maker.h"
47#include "AliTOFcalibHisto.h"
48#include "AliPID.h"
49#include "AliESDpid.h"
03bd764d 50#include "AliESDEvent.h"
51#include "TFile.h"
52#include "TH1F.h"
536031f2 53
54ClassImp(AliTOFT0maker)
55
56//____________________________________________________________________________
03bd764d 57 AliTOFT0maker::AliTOFT0maker():
58 fCalib(new AliTOFcalibHisto()),
59 fnT0(0),
60 fiT0(0),
61 fNoTOFT0(0),
62 fESDswitch(0),
63 fTimeResolution(115),
64 fT0sigma(1000),
65 fHmapChannel(0),
66 fKmask(0)
536031f2 67{
6a4e212e 68 // ctr
6a4e212e 69 fCalculated[0] = 0;
70 fCalculated[1] = 0;
71 fCalculated[2] = 0;
03bd764d 72 fCalculated[3] = 0;
6a4e212e 73
03bd764d 74 // fCalib->SetCalibParFileName("./AliTOFcalibPar.LHC10b.7000GeV.20100405.root");
536031f2 75 fCalib->LoadCalibPar();
76
77 if(AliPID::ParticleMass(0) == 0) new AliPID();
03bd764d 78
79 SetESDdata();
536031f2 80}
81//____________________________________________________________________________
8f589502 82AliTOFT0maker::AliTOFT0maker(const AliTOFT0maker & t) :
03bd764d 83 TObject(),
8f589502 84 fCalib(t.fCalib),
03bd764d 85 fnT0(t.fnT0),
86 fiT0(t.fiT0),
87 fNoTOFT0(t.fNoTOFT0),
8f589502 88 fESDswitch(t.fESDswitch),
89 fTimeResolution(t.fTimeResolution),
03bd764d 90 fT0sigma(t.fT0sigma),
91 fHmapChannel(t.fHmapChannel),
92 fKmask(t.fKmask)
8f589502 93{
03bd764d 94 // copy ctr
8f589502 95}
96
97//____________________________________________________________________________
98AliTOFT0maker& AliTOFT0maker::operator=(const AliTOFT0maker &t)
99{
6a4e212e 100 //
8f589502 101 // assign. operator
102 //
103
104 if (this == &t)
105 return *this;
106 fCalib = t.fCalib;
107 fESDswitch = t.fESDswitch;
108 fTimeResolution = t.fTimeResolution;
109 fT0sigma = t.fT0sigma;
110
111 return *this;
112}
113//____________________________________________________________________________
536031f2 114AliTOFT0maker::~AliTOFT0maker()
115{
116 // dtor
117 if(fCalib) delete fCalib;
118}
119//____________________________________________________________________________
120Double_t* AliTOFT0maker::RemakePID(AliESDEvent *esd,Double_t t0time,Double_t t0sigma){
8f589502 121 //
122 // Remake TOF PID probabilities
123 //
124
8f589502 125 Double_t *t0tof;
536031f2 126
03bd764d 127 if(fKmask) ApplyMask(esd);
128
536031f2 129 AliTOFT0v1* t0maker=new AliTOFT0v1(esd);
03bd764d 130// t0maker->SetCalib(fCalib);
536031f2 131 t0maker->SetTimeResolution(fTimeResolution*1e-12);
132
133 if(! fESDswitch){
536031f2 134 TakeTimeRawCorrection(esd);
135 }
03bd764d 136
137 t0tof=t0maker->DefineT0("all");
536031f2 138
8f589502 139 Float_t lT0Current=0.;
536031f2 140 fT0sigma=1000;
141
03bd764d 142// Int_t nrun = esd->GetRunNumber();
143 Double_t t0fill = GetT0Fill();
144 t0time += t0fill;
8f589502 145
03bd764d 146 Float_t sigmaFill = (t0fill - Int_t(t0fill))*1000;
147 if(sigmaFill < 0) sigmaFill += 1000;
8f589502 148
03bd764d 149 fCalculated[0]=-1000*t0tof[0]; // best t0
150 fCalculated[1]=1000*t0tof[1]; // sigma best t0
151 fCalculated[2] = t0fill; //t0 fill
152 fCalculated[3] = t0tof[2]; // n TOF tracks
153 fCalculated[4]=-1000*t0tof[0]; // TOF t0
154 fCalculated[5]=1000*t0tof[1]; // TOF t0 sigma
155 fCalculated[6]=sigmaFill; // sigma t0 fill
156 fCalculated[7] = t0tof[3]; // n TOF tracks used for T0
157
158 if(fCalculated[1] < sigmaFill){
159 if(fnT0 < 10){
160 fT0fill[fiT0] = fCalculated[0];
161 fT0sigmaTOF[fiT0] = fCalculated[1];
162 fiT0++;
163 fnT0++;
164 }
165 else if(TMath::Abs(fCalculated[0] - t0fill) < 500){
166 fT0fill[fiT0] = fCalculated[0];
167 fT0sigmaTOF[fiT0] = fCalculated[1];
168 fiT0++;
169 fnT0++;
170 }
171
172 // printf("%i - %i) %f\n",fiT0,fnT0,t0fill);
173 }
174 if(fnT0==10) fiT0=0;
175
176 if(fiT0 > fgkNmaxT0step-1) fiT0=0;
177
178 if(fnT0 < 100){
179 t0time -= t0fill;
180 sigmaFill=200;
181 t0fill=0;
182 fCalculated[2] = t0fill; //t0 fill
183 }
184
185 if(fCalculated[1] < sigmaFill && TMath::Abs(fCalculated[0] - t0fill) < 500){
6a4e212e 186 fT0sigma=fCalculated[1];
187 lT0Current=fCalculated[0];
536031f2 188 }
03bd764d 189 else{
190 fCalculated[4] = t0fill;
191 fCalculated[5] = sigmaFill;
192 }
193
194 if(fCalculated[1] < 1 || fT0sigma > sigmaFill){
195 fT0sigma =1000;
196 fCalculated[4] = t0fill;
197 fCalculated[5] = sigmaFill;
198 }
536031f2 199
200 if(t0sigma < 1000){
201 if(fT0sigma < 1000){
202 Double_t w1 = 1./t0sigma/t0sigma;
6a4e212e 203 Double_t w2 = 1./fCalculated[1]/fCalculated[1];
536031f2 204
205 Double_t wtot = w1+w2;
206
6a4e212e 207 lT0Current = (w1*t0time + w2*fCalculated[0]) / wtot;
536031f2 208 fT0sigma = TMath::Sqrt(1./wtot);
209 }
210 else{
8f589502 211 lT0Current=t0time;
536031f2 212 fT0sigma=t0sigma;
213 }
214 }
215
03bd764d 216 if(fT0sigma < sigmaFill && TMath::Abs(lT0Current - t0fill) < 500){
217 fCalculated[1]=fT0sigma;
218 fCalculated[0]=lT0Current;
219 }
220
221 if(fT0sigma >= 1000 || fNoTOFT0){
8f589502 222 lT0Current = t0fill;
03bd764d 223 fT0sigma = sigmaFill;
8f589502 224
6a4e212e 225 fCalculated[0] = t0fill;
03bd764d 226 fCalculated[1] = sigmaFill;
536031f2 227 }
228
536031f2 229
03bd764d 230
231 RemakeTOFpid(esd,lT0Current);
232
6a4e212e 233 return fCalculated;
536031f2 234}
235//____________________________________________________________________________
8f589502 236void AliTOFT0maker::TakeTimeRawCorrection(AliESDEvent * const esd){
237 //
238 // Take raw corrections for time measurements
239 //
03bd764d 240
536031f2 241 Int_t ntracks = esd->GetNumberOfTracks();
03bd764d 242
536031f2 243 while (ntracks--) {
244 AliESDtrack *t=esd->GetTrack(ntracks);
245
246 if ((t->GetStatus()&AliESDtrack::kTOFout)==0) continue;
247
248 Double_t time=t->GetTOFsignalRaw();
03bd764d 249
536031f2 250 Double_t tot = t->GetTOFsignalToT();
251 Int_t chan = t->GetTOFCalChannel();
252 Double_t corr = fCalib->GetFullCorrection(chan,tot) - fCalib->GetCorrection(AliTOFcalibHisto::kTimeSlewingCorr,chan,0);
253 time -= corr*1000;
254
536031f2 255 t->SetTOFsignal(time);
03bd764d 256
536031f2 257 }
258}
259//____________________________________________________________________________
260void AliTOFT0maker::RemakeTOFpid(AliESDEvent *esd,Float_t timezero){
8f589502 261 //
262 // Recalculate TOF PID probabilities
263 //
264
10d100d4 265 AliESDpid pidESD;
266 pidESD.GetTOFResponse().SetTimeResolution(TMath::Sqrt(fT0sigma*fT0sigma + fTimeResolution*fTimeResolution));
267 pidESD.MakePID(esd,kFALSE,timezero);
536031f2 268
269}
270//____________________________________________________________________________
03bd764d 271Double_t AliTOFT0maker::GetT0Fill() const {
8f589502 272 //
273 // Return T0 of filling
274 //
03bd764d 275
276 Double_t t0=0.200;
277
278 Int_t n=fnT0;
8f589502 279
03bd764d 280 if(n >10 && n <= 20) n = 10;
281 else if(n > 20){
282 n -= 10;
283 }
536031f2 284
03bd764d 285 if(n > fgkNmaxT0step) n = fgkNmaxT0step;
286
287 if(n>1){
288 Double_t lT0av=0;
289 Double_t lT0sigmaav=0;
290 Double_t lT0avErr=0;
291 for(Int_t i=0;i<n;i++){
292 lT0av+=fT0fill[i];
293 lT0sigmaav += fT0sigmaTOF[fiT0];
294 lT0avErr+=fT0fill[i]*fT0fill[i];
295 }
296 lT0avErr -= lT0av*lT0av/n;
297 lT0av /= n;
298 lT0sigmaav /= n;
299 lT0avErr = TMath::Sqrt(TMath::Max(lT0avErr/(n-1) - lT0sigmaav*lT0sigmaav,0.00001));
300
301
302 if(lT0avErr > 300) lT0avErr = 300;
303
304 lT0av = Int_t(lT0av) + lT0avErr/1000.;
305
306 return lT0av;
307 }
308
309
536031f2 310 return t0;
311}
03bd764d 312//____________________________________________________________________________
313void AliTOFT0maker::LoadChannelMap(char *filename){
314 // Load the histo with the channel off map
315 TFile *f= new TFile(filename);
316 if(!f){
317 printf("Cannot open the channel map file (%s)\n",filename);
318 return;
319 }
320
321 fHmapChannel = (TH1F *) f->Get("hChEnabled");
322
323 if(!fHmapChannel){
324 printf("Cannot laod the channel map histo (from %s)\n",filename);
325 return;
326 }
327
328}
329//____________________________________________________________________________
330void AliTOFT0maker::ApplyMask(AliESDEvent * const esd){
331 // Switch off the disable channel
332 if(!fHmapChannel){
333 printf("Channel Map is not available\n");
334 return;
335 }
336
337 Int_t ntracks = esd->GetNumberOfTracks();
338
339 while (ntracks--) {
340 AliESDtrack *t=esd->GetTrack(ntracks);
341
342 if ((t->GetStatus()&AliESDtrack::kTOFout)==0) continue;
343
344 Int_t chan = t->GetTOFCalChannel();
345
346 if(fHmapChannel->GetBinContent(chan) < 0.01){
347 t->ResetStatus(AliESDtrack::kTOFout);
348 }
349 }
350}
351
352//____________________________________________________________________________