]> git.uio.no Git - u/mrichter/AliRoot.git/blame - TOF/AliTOFcalib.cxx
Coding conventions (C.Oppedisano)
[u/mrichter/AliRoot.git] / TOF / AliTOFcalib.cxx
CommitLineData
6dc9348d 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
762446e0 16/*
17$Log$
28dd10b6 18Revision 1.3 2006/03/28 14:57:02 arcelli
19updates to handle new V5 geometry & some re-arrangements
20
d4ad0d6b 21Revision 1.2 2006/02/13 17:22:26 arcelli
22just Fixing Log info
23
762446e0 24Revision 1.1 2006/02/13 16:10:48 arcelli
25Add classes for TOF Calibration (C.Zampolli)
26
6dc9348d 27author: Chiara Zampolli, zampolli@bo.infn.it
762446e0 28*/
6dc9348d 29
30///////////////////////////////////////////////////////////////////////////////
31// //
32// class for TOF calibration //
33// //
34///////////////////////////////////////////////////////////////////////////////
35
36#include "AliTOFcalib.h"
37#include "AliRun.h"
38#include <TTask.h>
39#include <TFile.h>
40#include <TROOT.h>
41#include <TSystem.h>
42#include "AliTOF.h"
43#include "AliTOFcalibESD.h"
44#include "AliESD.h"
45#include <TObject.h>
46#include "TF1.h"
47#include "TH1F.h"
48#include "TH2F.h"
49#include "AliESDtrack.h"
50#include "AliTOFChannel.h"
51#include "AliTOFChSim.h"
d4ad0d6b 52#include "AliTOFGeometryV5.h"
6dc9348d 53#include "TClonesArray.h"
54#include "AliTOFCal.h"
55#include "TRandom.h"
56#include "AliTOFcluster.h"
57#include "TList.h"
58#include "AliCDBManager.h"
59#include "AliCDBMetaData.h"
60#include "AliCDBStorage.h"
61#include "AliCDBId.h"
62#include "AliCDBEntry.h"
63
64extern TROOT *gROOT;
65extern TStyle *gStyle;
66
67ClassImp(AliTOFcalib)
68
69const Int_t AliTOFcalib::fgkchannel = 5000;
6dc9348d 70//_______________________________________________________________________
71AliTOFcalib::AliTOFcalib():TTask("AliTOFcalib",""){
d4ad0d6b 72
6dc9348d 73 fArrayToT = 0x0;
74 fArrayTime = 0x0;
6dc9348d 75 fESDsel = 0x0;
d4ad0d6b 76 AliTOFGeometry *geom=new AliTOFGeometryV5();
77 AliInfo("V5 TOF Geometry is taken as the default");
78 fNSector = geom->NSectors();
79 fNPlate = geom->NPlates();
80 fNStripA = geom->NStripA();
81 fNStripB = geom->NStripB();
82 fNStripC = geom->NStripC();
83 fNpadZ = geom->NpadZ();
84 fNpadX = geom->NpadX();
85 fNChannels = fNSector*(2*(fNStripC+fNStripB)+fNStripA)*fNpadZ*fNpadX; //generalized version
86 fTOFCal = new AliTOFCal(geom);
87 fTOFSimCal = new AliTOFCal(geom);
88 fTOFCal->CreateArray();
89 fTOFSimCal->CreateArray();
90 delete geom;
6dc9348d 91}
d4ad0d6b 92//_______________________________________________________________________
93AliTOFcalib::AliTOFcalib(AliTOFGeometry *geom):TTask("AliTOFcalib",""){
6dc9348d 94 fArrayToT = 0x0;
95 fArrayTime = 0x0;
6dc9348d 96 fESDsel = 0x0;
d4ad0d6b 97 fNSector = geom->NSectors();
98 fNPlate = geom->NPlates();
99 fNStripA = geom->NStripA();
100 fNStripB = geom->NStripB();
101 fNStripC = geom->NStripC();
102 fNpadZ = geom->NpadZ();
103 fNpadX = geom->NpadX();
104 fNChannels = fNSector*(2*(fNStripC+fNStripB)+fNStripA)*fNpadZ*fNpadX; //generalized version
105 fTOFCal = new AliTOFCal(geom);
106 fTOFSimCal = new AliTOFCal(geom);
107 fTOFCal->CreateArray();
108 fTOFSimCal->CreateArray();
6dc9348d 109}
110//____________________________________________________________________________
111
112AliTOFcalib::AliTOFcalib(const AliTOFcalib & calib):TTask("AliTOFcalib","")
113{
114 fNSector = calib.fNSector;
115 fNPlate = calib.fNPlate;
116 fNStripA = calib.fNStripA;
117 fNStripB = calib.fNStripB;
118 fNStripC = calib.fNStripC;
119 fNpadZ = calib.fNpadZ;
120 fNpadX = calib.fNpadX;
d4ad0d6b 121 fNChannels = calib.fNChannels;
6dc9348d 122 fArrayToT = calib.fArrayToT;
123 fArrayTime = calib.fArrayTime;
6dc9348d 124 fTOFCal=calib.fTOFCal;
d4ad0d6b 125 fTOFSimCal = calib.fTOFSimCal;
6dc9348d 126}
127
128//____________________________________________________________________________
129
130AliTOFcalib::~AliTOFcalib()
131{
132 delete fArrayToT;
133 delete fArrayTime;
6dc9348d 134 delete fTOFCal;
d4ad0d6b 135 delete fTOFSimCal;
6dc9348d 136 delete fESDsel;
137}
6dc9348d 138//__________________________________________________________________________
139
140TF1* AliTOFcalib::SetFitFunctions(TH1F *histo){
141 TF1 * fpol[3];
142 const Int_t nbins = histo->GetNbinsX();
143 Float_t Delta = histo->GetBinWidth(1); //all the bins have the same width
144 Double_t max = histo->GetBinLowEdge(nbins)+Delta;
145 max = 15;
146 fpol[0]=new TF1("poly3","pol3",5,max);
147 fpol[1]=new TF1("poly4","pol4",5,max);
148 fpol[2]=new TF1("poly5","pol5",5,max);
149 char npoly[10];
150 Double_t chi[3]={1E6,1E6,1E6};
151 Int_t ndf[3]={-1,-1,-1};
152 Double_t Nchi[3]={1E6,1E6,1E6};
153 Double_t bestchi=1E6;
154 TF1 * fGold=0x0;
155 Int_t nonzero =0;
156 Int_t numberOfpar =0;
157 for (Int_t j=0; j<nbins; j++){
158 if (histo->GetBinContent(j)!=0) {
159 nonzero++;
160 }
161 }
162 Int_t norderfit = 0;
163 if (nonzero<=4) {
164 AliError(" Too few points in the histo. No fit performed.");
165 return 0x0;
166 }
167 else if (nonzero<=5) {
168 norderfit = 3;
169 AliInfo(" Only 3rd order polynomial fit possible.");
170 }
171 else if (nonzero<=6) {
172 norderfit = 4;
173 AliInfo(" Only 3rd and 4th order polynomial fit possible.");
174 }
175 else {
176 norderfit = 5;
177 AliInfo(" All 3rd, 4th and 5th order polynomial fit possible.");
178 }
179 for (Int_t ifun=norderfit-3;ifun<norderfit-2;ifun++){
180 sprintf(npoly,"poly%i",ifun+3);
181 histo->Fit(npoly, "ERN", " ", 5.,14.);
182 chi[ifun] = fpol[ifun]->GetChisquare();
183 ndf[ifun] = fpol[ifun]->GetNDF();
184 Nchi[ifun] = (Double_t)chi[ifun]/ndf[ifun];
185 if (Nchi[ifun]<bestchi) {
186 bestchi=Nchi[ifun];
187 fGold = fpol[ifun];
188 numberOfpar = fGold->GetNpar();
189 }
190 }
191 fGold=fpol[2]; //Gold fit function set to pol5 in any case
192 histo->Fit(fGold,"ER " ," ",5.,15.);
193 return fGold;
194}
195//____________________________________________________________________________
196
d4ad0d6b 197void AliTOFcalib::SelectESD(AliESD *event)
6dc9348d 198{
6dc9348d 199 Float_t LowerMomBound=0.8; // [GeV/c] default value Pb-Pb
200 Float_t UpperMomBound=1.8 ; // [GeV/c] default value Pb-Pb
201 Int_t ntrk =0;
202 Int_t ngoodtrkfinalToT = 0;
203 ntrk=event->GetNumberOfTracks();
204 fESDsel = new TObjArray(ntrk);
205 fESDsel->SetOwner();
206 TObjArray UCdatatemp(ntrk);
207 Int_t ngoodtrk = 0;
208 Int_t ngoodtrkfinal = 0;
209 Float_t mintime =1E6;
210 for (Int_t itrk=0; itrk<ntrk; itrk++) {
211 AliESDtrack *t=event->GetTrack(itrk);
212 //track selection: reconstrution to TOF:
213 if ((t->GetStatus()&AliESDtrack::kTOFout)==0) {
214 continue;
215 }
216 //IsStartedTimeIntegral
217 if ((t->GetStatus()&AliESDtrack::kTIME)==0) {
218 continue;
219 }
220 Double_t time=t->GetTOFsignal();
221 time*=1.E-3; // tof given in nanoseconds
222 if(time <= mintime)mintime=time;
223 Double_t mom=t->GetP();
224 if (!(mom<=UpperMomBound && mom>=LowerMomBound))continue;
225 UInt_t AssignedTOFcluster=t->GetTOFcluster();//index of the assigned TOF cluster, >0 ?
226 if(AssignedTOFcluster==0){ // not matched
227 continue;
228 }
6dc9348d 229 AliTOFcalibESD *unc = new AliTOFcalibESD;
230 unc->CopyFromAliESD(t);
6dc9348d 231 Double_t c1[15];
232 unc->GetExternalCovariance(c1);
233 UCdatatemp.Add(unc);
234 ngoodtrk++;
235 }
236 for (Int_t i = 0; i < ngoodtrk ; i ++){
237 AliTOFcalibESD *unc = (AliTOFcalibESD*)UCdatatemp.At(i);
238 if((unc->GetTOFsignal()-mintime*1.E3)<5.E3){
239 fESDsel->Add(unc);
240 ngoodtrkfinal++;
241 ngoodtrkfinalToT++;
242 }
243 }
244 fESDsel->Sort();
245}
246//_____________________________________________________________________________
247
248void AliTOFcalib::CombESDId(){
249 Float_t t0offset=0;
250 Float_t loffset=0;
251 Int_t ntracksinset=6;
252 Float_t exptof[6][3];
253 Float_t momentum[6]={0.,0.,0.,0.,0.,0.};
254 Int_t assparticle[6]={3,3,3,3,3,3};
255 Float_t massarray[3]={0.13957,0.493677,0.9382723};
256 Float_t timeofflight[6]={0.,0.,0.,0.,0.,0.};
257 Float_t beta[6]={0.,0.,0.,0.,0.,0.};
258 Float_t texp[6]={0.,0.,0.,0.,0.,0.};
259 Float_t sqMomError[6]={0.,0.,0.,0.,0.,0.};
260 Float_t sqTrackError[6]={0.,0.,0.,0.,0.,0.};
261 Float_t tracktoflen[6]={0.,0.,0.,0.,0.,0.};
262 Float_t TimeResolution = 0.90e-10; // 90 ps by default
263 Float_t timeresolutioninns=TimeResolution*(1.e+9); // convert in [ns]
264 Float_t timezero[6]={0.,0.,0.,0.,0.,0.};
265 Float_t weightedtimezero[6]={0.,0.,0.,0.,0.,0.};
266 Float_t besttimezero[6]={0.,0.,0.,0.,0.,0.};
267 Float_t bestchisquare[6]={0.,0.,0.,0.,0.,0.};
268 Float_t bestweightedtimezero[6]={0.,0.,0.,0.,0.,0.};
269 Float_t bestsqTrackError[6]={0.,0.,0.,0.,0.,0.};
270
271 Int_t nelements = fESDsel->GetEntries();
272 Int_t nset= (Int_t)(nelements/ntracksinset);
273 for (Int_t i=0; i< nset; i++) {
274
275 AliTOFcalibESD **unc=new AliTOFcalibESD*[ntracksinset];
276 for (Int_t itrk=0; itrk<ntracksinset; itrk++) {
277 Int_t index = itrk+i*ntracksinset;
278 AliTOFcalibESD *element=(AliTOFcalibESD*)fESDsel->At(index);
279 unc[itrk]=element;
280 }
281
282 for (Int_t j=0; j<ntracksinset; j++) {
283 AliTOFcalibESD *element=unc[j];
284 Double_t mom=element->GetP();
285 Double_t time=element->GetTOFsignal()*1.E-3; // in ns
286 Double_t exptime[10];
287 element->GetIntegratedTimes(exptime);
288 Double_t toflen=element->GetIntegratedLength()/100.; // in m
289 timeofflight[j]=time+t0offset;
290 tracktoflen[j]=toflen+loffset;
291 exptof[j][0]=exptime[2]*1.E-3+0.005;
292 exptof[j][1]=exptime[3]*1.E-3+0.005;
293 exptof[j][2]=exptime[4]*1.E-3+0.005;
294 momentum[j]=mom;
295 }
296 Float_t t0best=999.;
297 Float_t Et0best=999.;
298 Float_t chisquarebest=999.;
299 for (Int_t i1=0; i1<3;i1++) {
300 beta[0]=momentum[0]/sqrt(massarray[i1]*massarray[i1]+momentum[0]*momentum[0]);
301 texp[0]=exptof[0][i1];
302 for (Int_t i2=0; i2<3;i2++) {
303 beta[1]=momentum[1]/sqrt(massarray[i2]*massarray[i2]+momentum[1]*momentum[1]);
304 texp[1]=exptof[1][i2];
305 for (Int_t i3=0; i3<3;i3++) {
306 beta[2]=momentum[2]/sqrt(massarray[i3]*massarray[i3]+momentum[2]*momentum[2]);
307 texp[2]=exptof[2][i3];
308 for (Int_t i4=0; i4<3;i4++) {
309 beta[3]=momentum[3]/sqrt(massarray[i4]*massarray[i4]+momentum[3]*momentum[3]);
310 texp[3]=exptof[3][i4];
311
312 for (Int_t i5=0; i5<3;i5++) {
313 beta[4]=momentum[4]/sqrt(massarray[i5]*massarray[i5]+momentum[4]*momentum[4]);
314 texp[4]=exptof[4][i5];
315 for (Int_t i6=0; i6<3;i6++) {
316 beta[5]=momentum[5]/sqrt(massarray[i6]*massarray[i6]+momentum[5]*momentum[5]);
317 texp[5]=exptof[5][i6];
318
319 Float_t sumAllweights=0.;
320 Float_t meantzero=0.;
321 Float_t Emeantzero=0.;
322
323 for (Int_t itz=0; itz<ntracksinset;itz++) {
324 sqMomError[itz]=
325 ((1.-beta[itz]*beta[itz])*0.025)*
326 ((1.-beta[itz]*beta[itz])*0.025)*
327 (tracktoflen[itz]/
328 (0.299792*beta[itz]))*
329 (tracktoflen[itz]/
330 (0.299792*beta[itz]));
331 sqTrackError[itz]=
332 (timeresolutioninns*
333 timeresolutioninns
334 +sqMomError[itz]);
335
336 timezero[itz]=texp[itz]-timeofflight[itz];
337 weightedtimezero[itz]=timezero[itz]/sqTrackError[itz];
338 sumAllweights+=1./sqTrackError[itz];
339 meantzero+=weightedtimezero[itz];
340
341 } // end loop for (Int_t itz=0; itz<15;itz++)
342
343 meantzero=meantzero/sumAllweights; // it is given in [ns]
344 Emeantzero=sqrt(1./sumAllweights); // it is given in [ns]
345
346 // calculate chisquare
347
348 Float_t chisquare=0.;
349 for (Int_t icsq=0; icsq<ntracksinset;icsq++) {
350 chisquare+=(timezero[icsq]-meantzero)*(timezero[icsq]-meantzero)/sqTrackError[icsq];
351 } // end loop for (Int_t icsq=0; icsq<15;icsq++)
352 // cout << " chisquare " << chisquare << endl;
353
354 Int_t npion=0;
355 if(i1==0)npion++;
356 if(i2==0)npion++;
357 if(i3==0)npion++;
358 if(i4==0)npion++;
359 if(i5==0)npion++;
360 if(i6==0)npion++;
361
362 if(chisquare<=chisquarebest && ((Float_t) npion/ ((Float_t) ntracksinset)>0.3)){
363 // if(chisquare<=chisquarebest){
364
365 for(Int_t iqsq = 0; iqsq<ntracksinset; iqsq++) {
366 bestsqTrackError[iqsq]=sqTrackError[iqsq];
367 besttimezero[iqsq]=timezero[iqsq];
368 bestweightedtimezero[iqsq]=weightedtimezero[iqsq];
369 bestchisquare[iqsq]=(timezero[iqsq]-meantzero)*(timezero[iqsq]-meantzero)/sqTrackError[iqsq];
370 }
371
372 assparticle[0]=i1;
373 assparticle[1]=i2;
374 assparticle[2]=i3;
375 assparticle[3]=i4;
376 assparticle[4]=i5;
377 assparticle[5]=i6;
378
379 chisquarebest=chisquare;
380 t0best=meantzero;
381 Et0best=Emeantzero;
382 } // close if(dummychisquare<=chisquare)
383 } // end loop on i6
384 } // end loop on i5
385 } // end loop on i4
386 } // end loop on i3
387 } // end loop on i2
388 } // end loop on i1
389
390
391 Float_t confLevel=999;
392 if(chisquarebest<999.){
393 Double_t dblechisquare=(Double_t)chisquarebest;
394 confLevel=(Float_t)TMath::Prob(dblechisquare,ntracksinset-1);
395 }
396 // assume they are all pions for fake sets
397 if(confLevel<0.01 || confLevel==999. ){
398 for (Int_t itrk=0; itrk<ntracksinset; itrk++)assparticle[itrk]=0;
399 }
400 for (Int_t itrk=0; itrk<ntracksinset; itrk++) {
401 Int_t index = itrk+i*ntracksinset;
402 AliTOFcalibESD *element=(AliTOFcalibESD*)fESDsel->At(index);
403 element->SetCombID(assparticle[itrk]);
404 }
405 }
406}
407
408//_____________________________________________________________________________
409
410void AliTOFcalib::CalibrateESD(){
411 Int_t nelements = fESDsel->GetEntries();
d4ad0d6b 412 Int_t *number=new Int_t[fNChannels];
413 fArrayToT = new AliTOFArray(fNChannels);
414 fArrayTime = new AliTOFArray(fNChannels);
415 for (Int_t i=0; i<fNChannels; i++){
6dc9348d 416 number[i]=0;
417 fArrayToT->AddArray(i, new TArrayF(fgkchannel));
418 TArrayF * parrToT = fArrayToT->GetArray(i);
419 TArrayF & refaToT = * parrToT;
420 fArrayTime->AddArray(i, new TArrayF(fgkchannel));
421 TArrayF * parrTime = fArrayToT->GetArray(i);
422 TArrayF & refaTime = * parrTime;
423 for (Int_t j = 0;j<AliTOFcalib::fgkchannel;j++){
424 refaToT[j]=0.; //ToT[i][j]=j;
425 refaTime[j]=0.; //Time[i][j]=j;
426 }
427 }
428
429 for (Int_t i=0; i< nelements; i++) {
430 AliTOFcalibESD *element=(AliTOFcalibESD*)fESDsel->At(i);
431 Int_t ipid = element->GetCombID();
432 Double_t etime = 0; //expected time
433 Double_t expTime[10];
434 element->GetIntegratedTimes(expTime);
435 if (ipid == 0) etime = expTime[2]*1E-3; //ns
436 else if (ipid == 1) etime = expTime[3]*1E-3; //ns
437 else if (ipid == 2) etime = expTime[4]*1E-3; //ns
438 else AliError("No pid from combinatorial algo for this track");
439 Double_t mtime = (Double_t)element->GetTOFsignal()*1E-3; //measured time
440 Double_t mToT = (Double_t) element->GetToT(); //measured ToT, ns
6dc9348d 441 //select the correspondent channel with its simulated ToT spectrum
442 //summing up everything, index = 0 for all channels:
d4ad0d6b 443 Int_t index = element->GetTOFCalChannel();
6dc9348d 444 Int_t index2 = number[index];
445 TArrayF * parrToT = fArrayToT->GetArray(index);
446 TArrayF & refaToT = * parrToT;
447 refaToT[index2] = (Float_t)mToT;
448 TArrayF * parrTime = fArrayTime->GetArray(index);
449 TArrayF & refaTime = * parrTime;
450 refaTime[index2] = (Float_t)(mtime-etime);
451 number[index]++;
452 }
453
454 for (Int_t i=0;i<1;i++){
455 TH1F * hProf = Profile(i);
456 TF1* fGold = SetFitFunctions(hProf);
457 Int_t nfpar = fGold->GetNpar();
458 Float_t par[6];
459 for(Int_t kk=0;kk<6;kk++){
460 par[kk]=0;
461 }
462 for (Int_t kk = 0; kk< nfpar; kk++){
463 par[kk]=fGold->GetParameter(kk);
464 }
465 AliTOFChannel * CalChannel = fTOFCal->GetChannel(i);
466 CalChannel->SetSlewPar(par);
467 }
468 delete[] number;
469}
470
471//___________________________________________________________________________
472
473TH1F* AliTOFcalib::Profile(Int_t ich){
474 const Int_t nbinToT = 650;
475 Int_t nbinTime = 400;
476 Float_t minTime = -10.5; //ns
477 Float_t maxTime = 10.5; //ns
478 Float_t minToT = 7.5; //ns
479 Float_t maxToT = 40.; //ns
480 Float_t DeltaToT = (maxToT-minToT)/nbinToT;
481 Double_t mTime[nbinToT+1],mToT[nbinToT+1],meanTime[nbinToT+1], meanTime2[nbinToT+1],ToT[nbinToT+1], ToT2[nbinToT+1],meanToT[nbinToT+1],meanToT2[nbinToT+1],Time[nbinToT+1],Time2[nbinToT+1],xlow[nbinToT+1],sigmaTime[nbinToT+1];
482 Int_t n[nbinToT+1], nentrx[nbinToT+1];
483 Double_t sigmaToT[nbinToT+1];
484 for (Int_t i = 0; i < nbinToT+1 ; i++){
485 mTime[i]=0;
486 mToT[i]=0;
487 n[i]=0;
488 meanTime[i]=0;
489 meanTime2[i]=0;
490 ToT[i]=0;
491 ToT2[i]=0;
492 meanToT[i]=0;
493 meanToT2[i]=0;
494 Time[i]=0;
495 Time2[i]=0;
496 xlow[i]=0;
497 sigmaTime[i]=0;
498 sigmaToT[i]=0;
499 n[i]=0;
500 nentrx[i]=0;
501 }
502 TH2F* hSlewing = new TH2F("hSlewing", "hSlewing", nbinToT, minToT, maxToT, nbinTime, minTime, maxTime);
503 TArrayF * parrToT = fArrayToT->GetArray(ich);
504 TArrayF & refaToT = * parrToT;
505 TArrayF * parrTime = fArrayTime->GetArray(ich);
506 TArrayF & refaTime = * parrTime;
507 for (Int_t j = 0; j < AliTOFcalib::fgkchannel; j++){
508 if (refaToT[j] == 0) continue;
509 Int_t nx = (Int_t)((refaToT[j]-minToT)/DeltaToT)+1;
510 if ((refaToT[j] != 0) && (refaTime[j] != 0)){
511 Time[nx]+=refaTime[j];
512 Time2[nx]+=(refaTime[j])*(refaTime[j]);
513 ToT[nx]+=refaToT[j];
514 ToT2[nx]+=refaToT[j]*refaToT[j];
515 nentrx[nx]++;
516 hSlewing->Fill(refaToT[j],refaTime[j]);
517 }
518 }
519 Int_t nbinsToT=hSlewing->GetNbinsX();
520 if (nbinsToT != nbinToT) {
521 AliError("Profile :: incompatible numbers of bins");
522 return 0x0;
523 }
524
525 Int_t usefulBins=0;
526 TH1F *histo = new TH1F("histo", "1D Time vs ToT", nbinsToT, minToT, maxToT);
527 for (Int_t i=1;i<=nbinsToT;i++){
528 if (nentrx[i]!=0){
529 n[usefulBins]+=nentrx[i];
530 if (n[usefulBins]==0 && i == nbinsToT) {
531 break;
532 }
533 meanTime[usefulBins]+=Time[i];
534 meanTime2[usefulBins]+=Time2[i];
535 meanToT[usefulBins]+=ToT[i];
536 meanToT2[usefulBins]+=ToT2[i];
537 if (n[usefulBins]<20 && i!=nbinsToT) continue;
538 mTime[usefulBins]=meanTime[usefulBins]/n[usefulBins];
539 mToT[usefulBins]=meanToT[usefulBins]/n[usefulBins];
540 sigmaTime[usefulBins]=TMath::Sqrt(1./n[usefulBins]/n[usefulBins]
541 *(meanTime2[usefulBins]-meanTime[usefulBins]
542 *meanTime[usefulBins]/n[usefulBins]));
543 if ((1./n[usefulBins]/n[usefulBins]
544 *(meanToT2[usefulBins]-meanToT[usefulBins]
545 *meanToT[usefulBins]/n[usefulBins]))< 0) {
546 AliError(" too small radical" );
547 sigmaToT[usefulBins]=0;
548 }
549 else{
550 sigmaToT[usefulBins]=TMath::Sqrt(1./n[usefulBins]/n[usefulBins]
551 *(meanToT2[usefulBins]-meanToT[usefulBins]
552 *meanToT[usefulBins]/n[usefulBins]));
553 }
554 usefulBins++;
555 }
556 }
557 for (Int_t i=0;i<usefulBins;i++){
558 Int_t binN = (Int_t)((mToT[i]-minToT)/DeltaToT)+1;
559 histo->Fill(mToT[i],mTime[i]);
560 histo->SetBinError(binN,sigmaTime[i]);
561 }
562 return histo;
563}
564//_____________________________________________________________________________
565
566void AliTOFcalib::CorrectESDTime(){
567 Int_t nelements = fESDsel->GetEntries();
568 for (Int_t i=0; i< nelements; i++) {
569 AliTOFcalibESD *element=(AliTOFcalibESD*)fESDsel->At(i);
d4ad0d6b 570 Int_t index = element->GetTOFCalChannel();
571 Float_t ToT = element->GetToT();
6dc9348d 572 //select the correspondent channel with its simulated ToT spectrum
573 //summing up everything, index = 0 for all channels:
6dc9348d 574 Int_t ipid = element->GetCombID();
575 Double_t etime = 0; //expected time
576 Double_t expTime[10];
577 element->GetIntegratedTimes(expTime);
578 if (ipid == 0) etime = expTime[2]*1E-3; //ns
579 else if (ipid == 1) etime = expTime[3]*1E-3; //ns
580 else if (ipid == 2) etime = expTime[4]*1E-3; //ns
581 Float_t par[6];
582 AliTOFChannel * CalChannel = fTOFCal->GetChannel(index);
583 for (Int_t j = 0; j<6; j++){
584 par[j]=CalChannel->GetSlewPar(j);
585 }
d4ad0d6b 586 Float_t TimeCorr=0;
587 TimeCorr= par[0]+par[1]*ToT+par[2]*ToT*ToT+par[3]*ToT*ToT*ToT+par[4]*ToT*ToT*ToT*ToT+par[5]*ToT*ToT*ToT*ToT*ToT;
6dc9348d 588 }
589}
590//_____________________________________________________________________________
591
d4ad0d6b 592void AliTOFcalib::CorrectESDTime(AliESD *event){
593
6dc9348d 594 Int_t ntrk =0;
595 ntrk=event->GetNumberOfTracks();
596 for (Int_t itrk=0; itrk<ntrk; itrk++) {
597 AliESDtrack *t=event->GetTrack(itrk);
598 if ((t->GetStatus()&AliESDtrack::kTOFout)==0) {
599 continue;
600 }
601 //IsStartedTimeIntegral
602 if ((t->GetStatus()&AliESDtrack::kTIME)==0) {
603 continue;
604 }
605 UInt_t AssignedTOFcluster=t->GetTOFcluster();//index of the assigned TOF cluster, >0 ?
606 if(AssignedTOFcluster==0){ // not matched
607 continue;
608 }
d4ad0d6b 609 Int_t index = t->GetTOFCalChannel();
6dc9348d 610 AliTOFChannel * CalChannel = fTOFCal->GetChannel(index);
611 Float_t par[6];
612 for (Int_t j = 0; j<6; j++){
613 par[j]=CalChannel->GetSlewPar(j);
614 }
d4ad0d6b 615 Float_t ToT = t->GetTOFsignalToT();
6dc9348d 616 Float_t TimeCorr =0;
617 TimeCorr=par[0]+par[1]*ToT+par[2]*ToT*ToT+par[3]*ToT*ToT*ToT+par[4]*ToT*ToT*ToT*ToT+par[5]*ToT*ToT*ToT*ToT*ToT;
618 }
619}
620//_____________________________________________________________________________
621
d4ad0d6b 622void AliTOFcalib::WriteParOnCDB(Char_t *sel, Int_t minrun, Int_t maxrun){
6dc9348d 623 AliCDBManager *man = AliCDBManager::Instance();
d4ad0d6b 624 if(!man->IsDefaultStorageSet())man->SetDefaultStorage("local://$ALICE_ROOT");
28dd10b6 625 Char_t *sel1 = "Par" ;
d4ad0d6b 626 Char_t out[100];
627 sprintf(out,"%s/%s",sel,sel1);
628 AliCDBId id(out,minrun,maxrun);
6dc9348d 629 AliCDBMetaData *md = new AliCDBMetaData();
d4ad0d6b 630 md->SetResponsible("Chiara Zampolli");
6dc9348d 631 man->Put(fTOFCal,id,md);
632}
633//_____________________________________________________________________________
634
d4ad0d6b 635void AliTOFcalib::WriteParOnCDB(Char_t *sel, Int_t minrun, Int_t maxrun, AliTOFCal *cal){
6dc9348d 636 AliCDBManager *man = AliCDBManager::Instance();
d4ad0d6b 637 if(!man->IsDefaultStorageSet())man->SetDefaultStorage("local://$ALICE_ROOT");
28dd10b6 638 Char_t *sel1 = "Par" ;
d4ad0d6b 639 Char_t out[100];
640 sprintf(out,"%s/%s",sel,sel1);
641 AliCDBId id(out,minrun,maxrun);
642 AliCDBMetaData *md = new AliCDBMetaData();
643 md->SetResponsible("Chiara Zampolli");
644 man->Put(cal,id,md);
645}
646//_____________________________________________________________________________
647
648void AliTOFcalib::ReadParFromCDB(Char_t *sel, Int_t nrun){
649 AliCDBManager *man = AliCDBManager::Instance();
650 if(!man->IsDefaultStorageSet())man->SetDefaultStorage("local://$ALICE_ROOT");
28dd10b6 651 Char_t *sel1 = "Par" ;
d4ad0d6b 652 Char_t out[100];
653 sprintf(out,"%s/%s",sel,sel1);
654 AliCDBEntry *entry = man->Get(out,nrun);
6dc9348d 655 AliTOFCal *cal =(AliTOFCal*)entry->GetObject();
656 fTOFCal = cal;
657}
658//_____________________________________________________________________________
d4ad0d6b 659void AliTOFcalib::WriteSimParOnCDB(Char_t *sel, Int_t minrun, Int_t maxrun){
660
661
662 //for the time being, only one spectrum is used
663 TFile *spFile = new TFile("$ALICE_ROOT/TOF/data/spectrum.root","read");
664 TH1F * hToT;
665 // Retrieve ToT Spectrum
666 spFile->GetObject("ToT",hToT);
667
668 fTOFSimToT=hToT;
669
670 // Retrieve Time over TOT dependence
671
672 TH1F * h = (TH1F*)spFile->Get("TimeToTFit");
673 TList * list = (TList*)h->GetListOfFunctions();
674 TF1* f = (TF1*)list->At(0);
675 Float_t par[6] = {0,0,0,0,0,0};
676 Int_t npar=f->GetNpar();
677 for (Int_t ipar=0;ipar<npar;ipar++){
678 par[ipar]=f->GetParameter(ipar);
679 }
680
681 for(Int_t iTOFch=0; iTOFch<fTOFSimCal->NPads();iTOFch++){
682 AliTOFChannel * CalChannel = fTOFSimCal->GetChannel(iTOFch);
683 CalChannel->SetSlewPar(par);
684 }
685
686 // Store them in the CDB
687
688 AliCDBManager *man = AliCDBManager::Instance();
689 if(!man->IsDefaultStorageSet())man->SetDefaultStorage("local://$ALICE_ROOT");
690 AliCDBMetaData *md = new AliCDBMetaData();
691 md->SetResponsible("Chiara Zampolli");
28dd10b6 692 Char_t *sel1 = "SimPar" ;
d4ad0d6b 693 Char_t out[100];
694 sprintf(out,"%s/%s",sel,sel1);
695 AliCDBId id1(out,minrun,maxrun);
696 man->Put(fTOFSimCal,id1,md);
28dd10b6 697 Char_t *sel2 = "SimHisto" ;
d4ad0d6b 698 sprintf(out,"%s/%s",sel,sel2);
699 AliCDBId id2(out,minrun,maxrun);
700 man->Put(fTOFSimToT,id2,md);
701}
702
703//_____________________________________________________________________________
704void AliTOFcalib::WriteSimParOnCDB(Char_t *sel, Int_t minrun, Int_t maxrun, AliTOFCal *cal, TH1F * histo){
705
706 // Retrieve ToT Spectrum
707 fTOFSimToT=histo;
708 fTOFSimCal=cal;
709 AliCDBManager *man = AliCDBManager::Instance();
710 if(!man->IsDefaultStorageSet())man->SetDefaultStorage("local://$ALICE_ROOT");
711 AliCDBMetaData *md = new AliCDBMetaData();
712 md->SetResponsible("Chiara Zampolli");
28dd10b6 713 Char_t *sel1 = "SimPar" ;
d4ad0d6b 714 Char_t out[100];
715 sprintf(out,"%s/%s",sel,sel1);
716 AliCDBId id1(out,minrun,maxrun);
717 man->Put(fTOFSimCal,id1,md);
28dd10b6 718 Char_t *sel2 = "SimHisto" ;
d4ad0d6b 719 sprintf(out,"%s/%s",sel,sel2);
720 AliCDBId id2(out,minrun,maxrun);
721 man->Put(fTOFSimToT,id2,md);
722}
723//_____________________________________________________________________________
724void AliTOFcalib::ReadSimParFromCDB(Char_t *sel, Int_t nrun){
725 AliCDBManager *man = AliCDBManager::Instance();
726 if(!man->IsDefaultStorageSet())man->SetDefaultStorage("local://$ALICE_ROOT");
727 AliCDBMetaData *md = new AliCDBMetaData();
728 md->SetResponsible("Chiara Zampolli");
28dd10b6 729 Char_t *sel1 = "SimPar" ;
d4ad0d6b 730 Char_t out[100];
731 sprintf(out,"%s/%s",sel,sel1);
732 AliCDBEntry *entry1 = man->Get(out,nrun);
733 AliTOFCal *cal =(AliTOFCal*)entry1->GetObject();
734 fTOFSimCal=cal;
28dd10b6 735 Char_t *sel2 = "SimHisto" ;
d4ad0d6b 736 sprintf(out,"%s/%s",sel,sel2);
737 AliCDBEntry *entry2 = man->Get(out,nrun);
738 TH1F *histo =(TH1F*)entry2->GetObject();
739 fTOFSimToT=histo;
740}
741//_____________________________________________________________________________
6dc9348d 742
743Int_t AliTOFcalib::GetIndex(Int_t *detId){
744 Int_t isector = detId[0];
745 if (isector >= fNSector)
746 AliError(Form("Wrong sector number in TOF (%d) !",isector));
747 Int_t iplate = detId[1];
748 if (iplate >= fNPlate)
749 AliError(Form("Wrong plate number in TOF (%d) !",iplate));
750 Int_t istrip = detId[2];
751 Int_t ipadz = detId[3];
752 Int_t ipadx = detId[4];
753 Int_t stripOffset = 0;
754 switch (iplate) {
755 case 0:
756 stripOffset = 0;
757 break;
758 case 1:
759 stripOffset = fNStripC;
760 break;
761 case 2:
762 stripOffset = fNStripC+fNStripB;
763 break;
764 case 3:
765 stripOffset = fNStripC+fNStripB+fNStripA;
766 break;
767 case 4:
768 stripOffset = fNStripC+fNStripB+fNStripA+fNStripB;
769 break;
770 default:
771 AliError(Form("Wrong plate number in TOF (%d) !",iplate));
772 break;
773 };
774
775 Int_t idet = ((2*(fNStripC+fNStripB)+fNStripA)*fNpadZ*fNpadX)*isector +
776 (stripOffset*fNpadZ*fNpadX)+
777 (fNpadZ*fNpadX)*istrip+
778 (fNpadX)*ipadz+
779 ipadx;
780 return idet;
781}
782