]> git.uio.no Git - u/mrichter/AliRoot.git/blame - TOF/AliTOFcalib.cxx
changes for proper protection against failed retrieval of CDB Reco object (moved...
[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$
7037bd93 18Revision 1.11 2007/02/28 18:08:26 arcelli
19Add protection against failed retrieval of the CDB cal object
20
1b20c168 21Revision 1.10 2006/08/22 13:30:49 arcelli
22removal of effective c++ warnings (C.Zampolli)
23
655e379f 24Revision 1.9 2006/04/20 22:30:50 hristov
25Coding conventions (Annalisa)
26
0e46b9ae 27Revision 1.8 2006/04/16 22:29:05 hristov
28Coding conventions (Annalisa)
29
7aeeaf38 30Revision 1.7 2006/04/16 20:12:46 hristov
31Removing memory leak in case of cached CDB entries
32
85fc78e3 33Revision 1.6 2006/04/11 15:28:32 hristov
34Checks on cache status before deleting calibration objects (A.Colla)
35
5eaae242 36Revision 1.5 2006/04/05 08:35:38 hristov
37Coding conventions (S.Arcelli, C.Zampolli)
38
340693af 39Revision 1.4 2006/03/31 11:26:46 arcelli
40 changing CDB Ids according to standard convention
41
28dd10b6 42Revision 1.3 2006/03/28 14:57:02 arcelli
43updates to handle new V5 geometry & some re-arrangements
44
d4ad0d6b 45Revision 1.2 2006/02/13 17:22:26 arcelli
46just Fixing Log info
47
762446e0 48Revision 1.1 2006/02/13 16:10:48 arcelli
49Add classes for TOF Calibration (C.Zampolli)
50
6dc9348d 51author: Chiara Zampolli, zampolli@bo.infn.it
762446e0 52*/
6dc9348d 53
54///////////////////////////////////////////////////////////////////////////////
55// //
56// class for TOF calibration //
57// //
58///////////////////////////////////////////////////////////////////////////////
59
6dc9348d 60#include "TF1.h"
0e46b9ae 61#include "TFile.h"
6dc9348d 62#include "TH1F.h"
63#include "TH2F.h"
6dc9348d 64#include "TList.h"
0e46b9ae 65#include "TROOT.h"
66#include "TStyle.h"
67
68#include "AliCDBEntry.h"
69#include "AliCDBId.h"
6dc9348d 70#include "AliCDBManager.h"
71#include "AliCDBMetaData.h"
0e46b9ae 72#include "AliESDtrack.h"
73#include "AliESD.h"
74#include "AliLog.h"
75
76#include "AliTOFCal.h"
77#include "AliTOFcalibESD.h"
78#include "AliTOFcalib.h"
79#include "AliTOFChannel.h"
80#include "AliTOFGeometryV5.h"
81#include "AliTOFGeometry.h"
7037bd93 82#include "AliTOFRecoParam.h"
6dc9348d 83
84extern TROOT *gROOT;
85extern TStyle *gStyle;
86
87ClassImp(AliTOFcalib)
88
89const Int_t AliTOFcalib::fgkchannel = 5000;
6dc9348d 90//_______________________________________________________________________
655e379f 91AliTOFcalib::AliTOFcalib():
92 TTask("AliTOFcalib",""),
93 fNChannels(-1),
94 fNSector(-1),
95 fNPlate(-1),
96 fNStripA(-1),
97 fNStripB(-1),
98 fNStripC(-1),
99 fNpadZ(-1),
100 fNpadX(-1),
101 fNevents(0),
102 fESDsel(0x0),
103 fArrayToT(0x0),
104 fArrayTime(0x0),
105 fTOFCal(0x0),
106 fTOFSimCal(0x0),
107 fTOFSimToT(0x0)
108{
340693af 109 //TOF Calibration Class ctor
6dc9348d 110 fArrayToT = 0x0;
111 fArrayTime = 0x0;
6dc9348d 112 fESDsel = 0x0;
d4ad0d6b 113 AliTOFGeometry *geom=new AliTOFGeometryV5();
114 AliInfo("V5 TOF Geometry is taken as the default");
115 fNSector = geom->NSectors();
116 fNPlate = geom->NPlates();
117 fNStripA = geom->NStripA();
118 fNStripB = geom->NStripB();
119 fNStripC = geom->NStripC();
120 fNpadZ = geom->NpadZ();
121 fNpadX = geom->NpadX();
122 fNChannels = fNSector*(2*(fNStripC+fNStripB)+fNStripA)*fNpadZ*fNpadX; //generalized version
85fc78e3 123 fTOFCal = 0x0;
124 fTOFSimCal = 0x0;
125 fTOFSimToT = 0x0;
d4ad0d6b 126 delete geom;
6dc9348d 127}
d4ad0d6b 128//_______________________________________________________________________
655e379f 129AliTOFcalib::AliTOFcalib(AliTOFGeometry *geom):
130 TTask("AliTOFcalib",""),
131 fNChannels(-1),
132 fNSector(-1),
133 fNPlate(-1),
134 fNStripA(-1),
135 fNStripB(-1),
136 fNStripC(-1),
137 fNpadZ(-1),
138 fNpadX(-1),
139 fNevents(0),
140 fESDsel(0x0),
141 fArrayToT(0x0),
142 fArrayTime(0x0),
143 fTOFCal(0x0),
144 fTOFSimCal(0x0),
145 fTOFSimToT(0x0)
146{
340693af 147 //TOF Calibration Class ctor, taking the TOF geometry as input
6dc9348d 148 fArrayToT = 0x0;
149 fArrayTime = 0x0;
6dc9348d 150 fESDsel = 0x0;
d4ad0d6b 151 fNSector = geom->NSectors();
152 fNPlate = geom->NPlates();
153 fNStripA = geom->NStripA();
154 fNStripB = geom->NStripB();
155 fNStripC = geom->NStripC();
156 fNpadZ = geom->NpadZ();
157 fNpadX = geom->NpadX();
158 fNChannels = fNSector*(2*(fNStripC+fNStripB)+fNStripA)*fNpadZ*fNpadX; //generalized version
85fc78e3 159 fTOFCal = 0x0;
160 fTOFSimCal = 0x0;
161 fTOFSimToT = 0x0;
6dc9348d 162}
163//____________________________________________________________________________
164
655e379f 165AliTOFcalib::AliTOFcalib(const AliTOFcalib & calib):
166 TTask("AliTOFcalib",""),
167 fNChannels(-1),
168 fNSector(-1),
169 fNPlate(-1),
170 fNStripA(-1),
171 fNStripB(-1),
172 fNStripC(-1),
173 fNpadZ(-1),
174 fNpadX(-1),
175 fNevents(0),
176 fESDsel(0x0),
177 fArrayToT(0x0),
178 fArrayTime(0x0),
179 fTOFCal(0x0),
180 fTOFSimCal(0x0),
181 fTOFSimToT(0x0)
6dc9348d 182{
340693af 183 //TOF Calibration Class copy ctor
6dc9348d 184 fNSector = calib.fNSector;
185 fNPlate = calib.fNPlate;
186 fNStripA = calib.fNStripA;
187 fNStripB = calib.fNStripB;
188 fNStripC = calib.fNStripC;
189 fNpadZ = calib.fNpadZ;
190 fNpadX = calib.fNpadX;
d4ad0d6b 191 fNChannels = calib.fNChannels;
6dc9348d 192 fArrayToT = calib.fArrayToT;
193 fArrayTime = calib.fArrayTime;
6dc9348d 194 fTOFCal=calib.fTOFCal;
d4ad0d6b 195 fTOFSimCal = calib.fTOFSimCal;
340693af 196 fTOFSimToT=calib.fTOFSimToT;
6dc9348d 197}
198
199//____________________________________________________________________________
7aeeaf38 200
201AliTOFcalib& AliTOFcalib::operator=(const AliTOFcalib &calib)
202{
203 //TOF Calibration Class assignment operator
204 this->fNSector = calib.fNSector;
205 this->fNPlate = calib.fNPlate;
206 this->fNStripA = calib.fNStripA;
207 this->fNStripB = calib.fNStripB;
208 this->fNStripC = calib.fNStripC;
209 this->fNpadZ = calib.fNpadZ;
210 this->fNpadX = calib.fNpadX;
211 this->fNChannels = calib.fNChannels;
212 this->fArrayToT = calib.fArrayToT;
213 this->fArrayTime = calib.fArrayTime;
214 this->fTOFCal=calib.fTOFCal;
215 this->fTOFSimCal = calib.fTOFSimCal;
216 this->fTOFSimToT=calib.fTOFSimToT;
217 return *this;
218}
219
220//____________________________________________________________________________
6dc9348d 221
222AliTOFcalib::~AliTOFcalib()
223{
340693af 224 //TOF Calibration Class dtor
6dc9348d 225 delete fArrayToT;
226 delete fArrayTime;
5eaae242 227
228 if(!(AliCDBManager::Instance()->GetCacheFlag())){ // CDB objects must NOT be deleted if cache is active!
229 delete fTOFCal;
230 delete fTOFSimCal;
231 }
232
6dc9348d 233 delete fESDsel;
234}
6dc9348d 235//__________________________________________________________________________
236
340693af 237TF1* AliTOFcalib::SetFitFunctions(TH1F *histo)
238{
239 //Define Fit Functions for Slewing Correction
6dc9348d 240 TF1 * fpol[3];
340693af 241 const Int_t knbins = histo->GetNbinsX();
242 Float_t delta = histo->GetBinWidth(1); //all the bins have the same width
243 Double_t max = histo->GetBinLowEdge(knbins)+delta;
6dc9348d 244 max = 15;
245 fpol[0]=new TF1("poly3","pol3",5,max);
246 fpol[1]=new TF1("poly4","pol4",5,max);
247 fpol[2]=new TF1("poly5","pol5",5,max);
248 char npoly[10];
249 Double_t chi[3]={1E6,1E6,1E6};
250 Int_t ndf[3]={-1,-1,-1};
340693af 251 Double_t nchi[3]={1E6,1E6,1E6};
6dc9348d 252 Double_t bestchi=1E6;
253 TF1 * fGold=0x0;
254 Int_t nonzero =0;
255 Int_t numberOfpar =0;
340693af 256 for (Int_t j=0; j<knbins; j++){
6dc9348d 257 if (histo->GetBinContent(j)!=0) {
258 nonzero++;
259 }
260 }
261 Int_t norderfit = 0;
262 if (nonzero<=4) {
263 AliError(" Too few points in the histo. No fit performed.");
264 return 0x0;
265 }
266 else if (nonzero<=5) {
267 norderfit = 3;
268 AliInfo(" Only 3rd order polynomial fit possible.");
269 }
270 else if (nonzero<=6) {
271 norderfit = 4;
272 AliInfo(" Only 3rd and 4th order polynomial fit possible.");
273 }
274 else {
275 norderfit = 5;
276 AliInfo(" All 3rd, 4th and 5th order polynomial fit possible.");
277 }
278 for (Int_t ifun=norderfit-3;ifun<norderfit-2;ifun++){
279 sprintf(npoly,"poly%i",ifun+3);
280 histo->Fit(npoly, "ERN", " ", 5.,14.);
281 chi[ifun] = fpol[ifun]->GetChisquare();
282 ndf[ifun] = fpol[ifun]->GetNDF();
340693af 283 nchi[ifun] = (Double_t)chi[ifun]/ndf[ifun];
284 if (nchi[ifun]<bestchi) {
285 bestchi=nchi[ifun];
6dc9348d 286 fGold = fpol[ifun];
287 numberOfpar = fGold->GetNpar();
288 }
289 }
290 fGold=fpol[2]; //Gold fit function set to pol5 in any case
291 histo->Fit(fGold,"ER " ," ",5.,15.);
292 return fGold;
293}
294//____________________________________________________________________________
295
d4ad0d6b 296void AliTOFcalib::SelectESD(AliESD *event)
6dc9348d 297{
340693af 298 //track selection for Calibration
299 Float_t lowerMomBound=0.8; // [GeV/c] default value Pb-Pb
300 Float_t upperMomBound=1.8 ; // [GeV/c] default value Pb-Pb
6dc9348d 301 Int_t ntrk =0;
302 Int_t ngoodtrkfinalToT = 0;
303 ntrk=event->GetNumberOfTracks();
304 fESDsel = new TObjArray(ntrk);
305 fESDsel->SetOwner();
340693af 306 TObjArray uCdatatemp(ntrk);
6dc9348d 307 Int_t ngoodtrk = 0;
308 Int_t ngoodtrkfinal = 0;
309 Float_t mintime =1E6;
310 for (Int_t itrk=0; itrk<ntrk; itrk++) {
311 AliESDtrack *t=event->GetTrack(itrk);
312 //track selection: reconstrution to TOF:
313 if ((t->GetStatus()&AliESDtrack::kTOFout)==0) {
314 continue;
315 }
316 //IsStartedTimeIntegral
317 if ((t->GetStatus()&AliESDtrack::kTIME)==0) {
318 continue;
319 }
320 Double_t time=t->GetTOFsignal();
321 time*=1.E-3; // tof given in nanoseconds
322 if(time <= mintime)mintime=time;
323 Double_t mom=t->GetP();
340693af 324 if (!(mom<=upperMomBound && mom>=lowerMomBound))continue;
325 UInt_t assignedTOFcluster=t->GetTOFcluster();//index of the assigned TOF cluster, >0 ?
326 if(assignedTOFcluster==0){ // not matched
6dc9348d 327 continue;
328 }
6dc9348d 329 AliTOFcalibESD *unc = new AliTOFcalibESD;
330 unc->CopyFromAliESD(t);
6dc9348d 331 Double_t c1[15];
332 unc->GetExternalCovariance(c1);
340693af 333 uCdatatemp.Add(unc);
6dc9348d 334 ngoodtrk++;
335 }
336 for (Int_t i = 0; i < ngoodtrk ; i ++){
340693af 337 AliTOFcalibESD *unc = (AliTOFcalibESD*)uCdatatemp.At(i);
6dc9348d 338 if((unc->GetTOFsignal()-mintime*1.E3)<5.E3){
339 fESDsel->Add(unc);
340 ngoodtrkfinal++;
341 ngoodtrkfinalToT++;
342 }
343 }
344 fESDsel->Sort();
345}
346//_____________________________________________________________________________
347
340693af 348void AliTOFcalib::CombESDId()
349{
350 //track PID for calibration
6dc9348d 351 Float_t t0offset=0;
352 Float_t loffset=0;
353 Int_t ntracksinset=6;
354 Float_t exptof[6][3];
355 Float_t momentum[6]={0.,0.,0.,0.,0.,0.};
356 Int_t assparticle[6]={3,3,3,3,3,3};
357 Float_t massarray[3]={0.13957,0.493677,0.9382723};
358 Float_t timeofflight[6]={0.,0.,0.,0.,0.,0.};
359 Float_t beta[6]={0.,0.,0.,0.,0.,0.};
360 Float_t texp[6]={0.,0.,0.,0.,0.,0.};
361 Float_t sqMomError[6]={0.,0.,0.,0.,0.,0.};
362 Float_t sqTrackError[6]={0.,0.,0.,0.,0.,0.};
363 Float_t tracktoflen[6]={0.,0.,0.,0.,0.,0.};
340693af 364 Float_t timeResolution = 0.90e-10; // 90 ps by default
365 Float_t timeresolutioninns=timeResolution*(1.e+9); // convert in [ns]
6dc9348d 366 Float_t timezero[6]={0.,0.,0.,0.,0.,0.};
367 Float_t weightedtimezero[6]={0.,0.,0.,0.,0.,0.};
368 Float_t besttimezero[6]={0.,0.,0.,0.,0.,0.};
369 Float_t bestchisquare[6]={0.,0.,0.,0.,0.,0.};
370 Float_t bestweightedtimezero[6]={0.,0.,0.,0.,0.,0.};
371 Float_t bestsqTrackError[6]={0.,0.,0.,0.,0.,0.};
372
373 Int_t nelements = fESDsel->GetEntries();
374 Int_t nset= (Int_t)(nelements/ntracksinset);
375 for (Int_t i=0; i< nset; i++) {
376
377 AliTOFcalibESD **unc=new AliTOFcalibESD*[ntracksinset];
378 for (Int_t itrk=0; itrk<ntracksinset; itrk++) {
379 Int_t index = itrk+i*ntracksinset;
380 AliTOFcalibESD *element=(AliTOFcalibESD*)fESDsel->At(index);
381 unc[itrk]=element;
382 }
383
384 for (Int_t j=0; j<ntracksinset; j++) {
385 AliTOFcalibESD *element=unc[j];
386 Double_t mom=element->GetP();
387 Double_t time=element->GetTOFsignal()*1.E-3; // in ns
388 Double_t exptime[10];
389 element->GetIntegratedTimes(exptime);
390 Double_t toflen=element->GetIntegratedLength()/100.; // in m
391 timeofflight[j]=time+t0offset;
392 tracktoflen[j]=toflen+loffset;
393 exptof[j][0]=exptime[2]*1.E-3+0.005;
394 exptof[j][1]=exptime[3]*1.E-3+0.005;
395 exptof[j][2]=exptime[4]*1.E-3+0.005;
396 momentum[j]=mom;
397 }
398 Float_t t0best=999.;
340693af 399 Float_t et0best=999.;
6dc9348d 400 Float_t chisquarebest=999.;
401 for (Int_t i1=0; i1<3;i1++) {
402 beta[0]=momentum[0]/sqrt(massarray[i1]*massarray[i1]+momentum[0]*momentum[0]);
403 texp[0]=exptof[0][i1];
404 for (Int_t i2=0; i2<3;i2++) {
405 beta[1]=momentum[1]/sqrt(massarray[i2]*massarray[i2]+momentum[1]*momentum[1]);
406 texp[1]=exptof[1][i2];
407 for (Int_t i3=0; i3<3;i3++) {
408 beta[2]=momentum[2]/sqrt(massarray[i3]*massarray[i3]+momentum[2]*momentum[2]);
409 texp[2]=exptof[2][i3];
410 for (Int_t i4=0; i4<3;i4++) {
411 beta[3]=momentum[3]/sqrt(massarray[i4]*massarray[i4]+momentum[3]*momentum[3]);
412 texp[3]=exptof[3][i4];
413
414 for (Int_t i5=0; i5<3;i5++) {
415 beta[4]=momentum[4]/sqrt(massarray[i5]*massarray[i5]+momentum[4]*momentum[4]);
416 texp[4]=exptof[4][i5];
417 for (Int_t i6=0; i6<3;i6++) {
418 beta[5]=momentum[5]/sqrt(massarray[i6]*massarray[i6]+momentum[5]*momentum[5]);
419 texp[5]=exptof[5][i6];
420
421 Float_t sumAllweights=0.;
422 Float_t meantzero=0.;
340693af 423 Float_t emeantzero=0.;
6dc9348d 424
425 for (Int_t itz=0; itz<ntracksinset;itz++) {
426 sqMomError[itz]=
427 ((1.-beta[itz]*beta[itz])*0.025)*
428 ((1.-beta[itz]*beta[itz])*0.025)*
429 (tracktoflen[itz]/
430 (0.299792*beta[itz]))*
431 (tracktoflen[itz]/
432 (0.299792*beta[itz]));
433 sqTrackError[itz]=
434 (timeresolutioninns*
435 timeresolutioninns
436 +sqMomError[itz]);
437
438 timezero[itz]=texp[itz]-timeofflight[itz];
439 weightedtimezero[itz]=timezero[itz]/sqTrackError[itz];
440 sumAllweights+=1./sqTrackError[itz];
441 meantzero+=weightedtimezero[itz];
442
443 } // end loop for (Int_t itz=0; itz<15;itz++)
444
445 meantzero=meantzero/sumAllweights; // it is given in [ns]
340693af 446 emeantzero=sqrt(1./sumAllweights); // it is given in [ns]
6dc9348d 447
448 // calculate chisquare
449
450 Float_t chisquare=0.;
451 for (Int_t icsq=0; icsq<ntracksinset;icsq++) {
452 chisquare+=(timezero[icsq]-meantzero)*(timezero[icsq]-meantzero)/sqTrackError[icsq];
453 } // end loop for (Int_t icsq=0; icsq<15;icsq++)
454 // cout << " chisquare " << chisquare << endl;
455
456 Int_t npion=0;
457 if(i1==0)npion++;
458 if(i2==0)npion++;
459 if(i3==0)npion++;
460 if(i4==0)npion++;
461 if(i5==0)npion++;
462 if(i6==0)npion++;
463
464 if(chisquare<=chisquarebest && ((Float_t) npion/ ((Float_t) ntracksinset)>0.3)){
465 // if(chisquare<=chisquarebest){
466
467 for(Int_t iqsq = 0; iqsq<ntracksinset; iqsq++) {
468 bestsqTrackError[iqsq]=sqTrackError[iqsq];
469 besttimezero[iqsq]=timezero[iqsq];
470 bestweightedtimezero[iqsq]=weightedtimezero[iqsq];
471 bestchisquare[iqsq]=(timezero[iqsq]-meantzero)*(timezero[iqsq]-meantzero)/sqTrackError[iqsq];
472 }
473
474 assparticle[0]=i1;
475 assparticle[1]=i2;
476 assparticle[2]=i3;
477 assparticle[3]=i4;
478 assparticle[4]=i5;
479 assparticle[5]=i6;
480
481 chisquarebest=chisquare;
482 t0best=meantzero;
340693af 483 et0best=emeantzero;
6dc9348d 484 } // close if(dummychisquare<=chisquare)
485 } // end loop on i6
486 } // end loop on i5
487 } // end loop on i4
488 } // end loop on i3
489 } // end loop on i2
490 } // end loop on i1
491
492
493 Float_t confLevel=999;
494 if(chisquarebest<999.){
495 Double_t dblechisquare=(Double_t)chisquarebest;
496 confLevel=(Float_t)TMath::Prob(dblechisquare,ntracksinset-1);
497 }
498 // assume they are all pions for fake sets
499 if(confLevel<0.01 || confLevel==999. ){
500 for (Int_t itrk=0; itrk<ntracksinset; itrk++)assparticle[itrk]=0;
501 }
502 for (Int_t itrk=0; itrk<ntracksinset; itrk++) {
503 Int_t index = itrk+i*ntracksinset;
504 AliTOFcalibESD *element=(AliTOFcalibESD*)fESDsel->At(index);
505 element->SetCombID(assparticle[itrk]);
506 }
507 }
508}
509
510//_____________________________________________________________________________
511
512void AliTOFcalib::CalibrateESD(){
340693af 513 //Calibrate selected ESD times
6dc9348d 514 Int_t nelements = fESDsel->GetEntries();
d4ad0d6b 515 Int_t *number=new Int_t[fNChannels];
516 fArrayToT = new AliTOFArray(fNChannels);
517 fArrayTime = new AliTOFArray(fNChannels);
518 for (Int_t i=0; i<fNChannels; i++){
6dc9348d 519 number[i]=0;
520 fArrayToT->AddArray(i, new TArrayF(fgkchannel));
521 TArrayF * parrToT = fArrayToT->GetArray(i);
522 TArrayF & refaToT = * parrToT;
523 fArrayTime->AddArray(i, new TArrayF(fgkchannel));
524 TArrayF * parrTime = fArrayToT->GetArray(i);
525 TArrayF & refaTime = * parrTime;
526 for (Int_t j = 0;j<AliTOFcalib::fgkchannel;j++){
527 refaToT[j]=0.; //ToT[i][j]=j;
528 refaTime[j]=0.; //Time[i][j]=j;
529 }
530 }
531
532 for (Int_t i=0; i< nelements; i++) {
533 AliTOFcalibESD *element=(AliTOFcalibESD*)fESDsel->At(i);
534 Int_t ipid = element->GetCombID();
535 Double_t etime = 0; //expected time
536 Double_t expTime[10];
537 element->GetIntegratedTimes(expTime);
538 if (ipid == 0) etime = expTime[2]*1E-3; //ns
539 else if (ipid == 1) etime = expTime[3]*1E-3; //ns
540 else if (ipid == 2) etime = expTime[4]*1E-3; //ns
541 else AliError("No pid from combinatorial algo for this track");
542 Double_t mtime = (Double_t)element->GetTOFsignal()*1E-3; //measured time
543 Double_t mToT = (Double_t) element->GetToT(); //measured ToT, ns
6dc9348d 544 //select the correspondent channel with its simulated ToT spectrum
545 //summing up everything, index = 0 for all channels:
d4ad0d6b 546 Int_t index = element->GetTOFCalChannel();
6dc9348d 547 Int_t index2 = number[index];
548 TArrayF * parrToT = fArrayToT->GetArray(index);
549 TArrayF & refaToT = * parrToT;
550 refaToT[index2] = (Float_t)mToT;
551 TArrayF * parrTime = fArrayTime->GetArray(index);
552 TArrayF & refaTime = * parrTime;
553 refaTime[index2] = (Float_t)(mtime-etime);
554 number[index]++;
555 }
556
557 for (Int_t i=0;i<1;i++){
558 TH1F * hProf = Profile(i);
559 TF1* fGold = SetFitFunctions(hProf);
560 Int_t nfpar = fGold->GetNpar();
561 Float_t par[6];
562 for(Int_t kk=0;kk<6;kk++){
563 par[kk]=0;
564 }
565 for (Int_t kk = 0; kk< nfpar; kk++){
566 par[kk]=fGold->GetParameter(kk);
567 }
85fc78e3 568 if (!fTOFCal) {
569 AliTOFGeometry *geom=new AliTOFGeometryV5();
570 fTOFCal = new AliTOFCal(geom);
571 fTOFCal->CreateArray();
572 delete geom;
573 }
340693af 574 AliTOFChannel * calChannel = fTOFCal->GetChannel(i);
575 calChannel->SetSlewPar(par);
6dc9348d 576 }
577 delete[] number;
578}
579
580//___________________________________________________________________________
581
340693af 582TH1F* AliTOFcalib::Profile(Int_t ich)
583{
584 //Prepare histograms for Slewing Correction
585 const Int_t knbinToT = 650;
6dc9348d 586 Int_t nbinTime = 400;
587 Float_t minTime = -10.5; //ns
588 Float_t maxTime = 10.5; //ns
589 Float_t minToT = 7.5; //ns
590 Float_t maxToT = 40.; //ns
340693af 591 Float_t deltaToT = (maxToT-minToT)/knbinToT;
592 Double_t mTime[knbinToT+1],mToT[knbinToT+1],meanTime[knbinToT+1], meanTime2[knbinToT+1],vToT[knbinToT+1], vToT2[knbinToT+1],meanToT[knbinToT+1],meanToT2[knbinToT+1],vTime[knbinToT+1],vTime2[knbinToT+1],xlow[knbinToT+1],sigmaTime[knbinToT+1];
593 Int_t n[knbinToT+1], nentrx[knbinToT+1];
594 Double_t sigmaToT[knbinToT+1];
595 for (Int_t i = 0; i < knbinToT+1 ; i++){
6dc9348d 596 mTime[i]=0;
597 mToT[i]=0;
598 n[i]=0;
599 meanTime[i]=0;
600 meanTime2[i]=0;
340693af 601 vToT[i]=0;
602 vToT2[i]=0;
6dc9348d 603 meanToT[i]=0;
604 meanToT2[i]=0;
340693af 605 vTime[i]=0;
606 vTime2[i]=0;
6dc9348d 607 xlow[i]=0;
608 sigmaTime[i]=0;
609 sigmaToT[i]=0;
610 n[i]=0;
611 nentrx[i]=0;
612 }
340693af 613 TH2F* hSlewing = new TH2F("hSlewing", "hSlewing", knbinToT, minToT, maxToT, nbinTime, minTime, maxTime);
6dc9348d 614 TArrayF * parrToT = fArrayToT->GetArray(ich);
615 TArrayF & refaToT = * parrToT;
616 TArrayF * parrTime = fArrayTime->GetArray(ich);
617 TArrayF & refaTime = * parrTime;
618 for (Int_t j = 0; j < AliTOFcalib::fgkchannel; j++){
619 if (refaToT[j] == 0) continue;
340693af 620 Int_t nx = (Int_t)((refaToT[j]-minToT)/deltaToT)+1;
6dc9348d 621 if ((refaToT[j] != 0) && (refaTime[j] != 0)){
340693af 622 vTime[nx]+=refaTime[j];
623 vTime2[nx]+=(refaTime[j])*(refaTime[j]);
624 vToT[nx]+=refaToT[j];
625 vToT2[nx]+=refaToT[j]*refaToT[j];
6dc9348d 626 nentrx[nx]++;
627 hSlewing->Fill(refaToT[j],refaTime[j]);
628 }
629 }
630 Int_t nbinsToT=hSlewing->GetNbinsX();
340693af 631 if (nbinsToT != knbinToT) {
6dc9348d 632 AliError("Profile :: incompatible numbers of bins");
633 return 0x0;
634 }
635
636 Int_t usefulBins=0;
637 TH1F *histo = new TH1F("histo", "1D Time vs ToT", nbinsToT, minToT, maxToT);
638 for (Int_t i=1;i<=nbinsToT;i++){
639 if (nentrx[i]!=0){
640 n[usefulBins]+=nentrx[i];
641 if (n[usefulBins]==0 && i == nbinsToT) {
642 break;
643 }
340693af 644 meanTime[usefulBins]+=vTime[i];
645 meanTime2[usefulBins]+=vTime2[i];
646 meanToT[usefulBins]+=vToT[i];
647 meanToT2[usefulBins]+=vToT2[i];
6dc9348d 648 if (n[usefulBins]<20 && i!=nbinsToT) continue;
649 mTime[usefulBins]=meanTime[usefulBins]/n[usefulBins];
650 mToT[usefulBins]=meanToT[usefulBins]/n[usefulBins];
651 sigmaTime[usefulBins]=TMath::Sqrt(1./n[usefulBins]/n[usefulBins]
652 *(meanTime2[usefulBins]-meanTime[usefulBins]
653 *meanTime[usefulBins]/n[usefulBins]));
654 if ((1./n[usefulBins]/n[usefulBins]
655 *(meanToT2[usefulBins]-meanToT[usefulBins]
656 *meanToT[usefulBins]/n[usefulBins]))< 0) {
657 AliError(" too small radical" );
658 sigmaToT[usefulBins]=0;
659 }
660 else{
661 sigmaToT[usefulBins]=TMath::Sqrt(1./n[usefulBins]/n[usefulBins]
662 *(meanToT2[usefulBins]-meanToT[usefulBins]
663 *meanToT[usefulBins]/n[usefulBins]));
664 }
665 usefulBins++;
666 }
667 }
668 for (Int_t i=0;i<usefulBins;i++){
340693af 669 Int_t binN = (Int_t)((mToT[i]-minToT)/deltaToT)+1;
6dc9348d 670 histo->Fill(mToT[i],mTime[i]);
671 histo->SetBinError(binN,sigmaTime[i]);
672 }
673 return histo;
674}
675//_____________________________________________________________________________
676
340693af 677void AliTOFcalib::CorrectESDTime()
678{
679 //Calculate the corrected TOF time
6dc9348d 680 Int_t nelements = fESDsel->GetEntries();
681 for (Int_t i=0; i< nelements; i++) {
682 AliTOFcalibESD *element=(AliTOFcalibESD*)fESDsel->At(i);
d4ad0d6b 683 Int_t index = element->GetTOFCalChannel();
340693af 684 Float_t tToT = element->GetToT();
6dc9348d 685 //select the correspondent channel with its simulated ToT spectrum
686 //summing up everything, index = 0 for all channels:
6dc9348d 687 Int_t ipid = element->GetCombID();
688 Double_t etime = 0; //expected time
689 Double_t expTime[10];
690 element->GetIntegratedTimes(expTime);
691 if (ipid == 0) etime = expTime[2]*1E-3; //ns
692 else if (ipid == 1) etime = expTime[3]*1E-3; //ns
693 else if (ipid == 2) etime = expTime[4]*1E-3; //ns
694 Float_t par[6];
85fc78e3 695 if (!fTOFCal) {
696 AliTOFGeometry *geom=new AliTOFGeometryV5();
697 fTOFCal = new AliTOFCal(geom);
698 fTOFCal->CreateArray();
699 delete geom;
700 }
340693af 701 AliTOFChannel * calChannel = fTOFCal->GetChannel(index);
6dc9348d 702 for (Int_t j = 0; j<6; j++){
340693af 703 par[j]=calChannel->GetSlewPar(j);
6dc9348d 704 }
340693af 705 Float_t timeCorr=0;
706 timeCorr= par[0]+par[1]*tToT+par[2]*tToT*tToT+par[3]*tToT*tToT*tToT+par[4]*tToT*tToT*tToT*tToT+par[5]*tToT*tToT*tToT*tToT*tToT;
6dc9348d 707 }
708}
709//_____________________________________________________________________________
710
d4ad0d6b 711void AliTOFcalib::CorrectESDTime(AliESD *event){
340693af 712 //Calculate the corrected TOF time
d4ad0d6b 713
6dc9348d 714 Int_t ntrk =0;
715 ntrk=event->GetNumberOfTracks();
716 for (Int_t itrk=0; itrk<ntrk; itrk++) {
717 AliESDtrack *t=event->GetTrack(itrk);
718 if ((t->GetStatus()&AliESDtrack::kTOFout)==0) {
719 continue;
720 }
721 //IsStartedTimeIntegral
722 if ((t->GetStatus()&AliESDtrack::kTIME)==0) {
723 continue;
724 }
340693af 725 UInt_t assignedTOFcluster=t->GetTOFcluster();//index of the assigned TOF cluster, >0 ?
726 if(assignedTOFcluster==0){ // not matched
6dc9348d 727 continue;
728 }
d4ad0d6b 729 Int_t index = t->GetTOFCalChannel();
85fc78e3 730 if (!fTOFCal) {
731 AliTOFGeometry *geom=new AliTOFGeometryV5();
732 fTOFCal = new AliTOFCal(geom);
733 fTOFCal->CreateArray();
734 delete geom;
735 }
340693af 736 AliTOFChannel * calChannel = fTOFCal->GetChannel(index);
6dc9348d 737 Float_t par[6];
738 for (Int_t j = 0; j<6; j++){
340693af 739 par[j]=calChannel->GetSlewPar(j);
6dc9348d 740 }
340693af 741 Float_t tToT = t->GetTOFsignalToT();
742 Float_t timeCorr =0;
743 timeCorr=par[0]+par[1]*tToT+par[2]*tToT*tToT+par[3]*tToT*tToT*tToT+par[4]*tToT*tToT*tToT*tToT+par[5]*tToT*tToT*tToT*tToT*tToT;
6dc9348d 744 }
745}
746//_____________________________________________________________________________
747
340693af 748void AliTOFcalib::WriteParOnCDB(Char_t *sel, Int_t minrun, Int_t maxrun)
749{
750 //Write calibration parameters to the CDB
6dc9348d 751 AliCDBManager *man = AliCDBManager::Instance();
d4ad0d6b 752 if(!man->IsDefaultStorageSet())man->SetDefaultStorage("local://$ALICE_ROOT");
28dd10b6 753 Char_t *sel1 = "Par" ;
d4ad0d6b 754 Char_t out[100];
755 sprintf(out,"%s/%s",sel,sel1);
756 AliCDBId id(out,minrun,maxrun);
6dc9348d 757 AliCDBMetaData *md = new AliCDBMetaData();
d4ad0d6b 758 md->SetResponsible("Chiara Zampolli");
85fc78e3 759 if (!fTOFCal) {
760 AliTOFGeometry *geom=new AliTOFGeometryV5();
761 fTOFCal = new AliTOFCal(geom);
762 fTOFCal->CreateArray();
763 delete geom;
764 }
6dc9348d 765 man->Put(fTOFCal,id,md);
5eaae242 766 delete md;
6dc9348d 767}
768//_____________________________________________________________________________
769
d4ad0d6b 770void AliTOFcalib::WriteParOnCDB(Char_t *sel, Int_t minrun, Int_t maxrun, AliTOFCal *cal){
340693af 771 //Write calibration parameters to the CDB
6dc9348d 772 AliCDBManager *man = AliCDBManager::Instance();
d4ad0d6b 773 if(!man->IsDefaultStorageSet())man->SetDefaultStorage("local://$ALICE_ROOT");
28dd10b6 774 Char_t *sel1 = "Par" ;
d4ad0d6b 775 Char_t out[100];
776 sprintf(out,"%s/%s",sel,sel1);
777 AliCDBId id(out,minrun,maxrun);
778 AliCDBMetaData *md = new AliCDBMetaData();
779 md->SetResponsible("Chiara Zampolli");
780 man->Put(cal,id,md);
5eaae242 781 delete md;
d4ad0d6b 782}
783//_____________________________________________________________________________
784
1b20c168 785Bool_t AliTOFcalib::ReadParFromCDB(Char_t *sel, Int_t nrun)
340693af 786{
787 //Read calibration parameters from the CDB
d4ad0d6b 788 AliCDBManager *man = AliCDBManager::Instance();
789 if(!man->IsDefaultStorageSet())man->SetDefaultStorage("local://$ALICE_ROOT");
28dd10b6 790 Char_t *sel1 = "Par" ;
d4ad0d6b 791 Char_t out[100];
792 sprintf(out,"%s/%s",sel,sel1);
793 AliCDBEntry *entry = man->Get(out,nrun);
1b20c168 794 if (!entry) {
795 AliError(Form("Failed to get entry: %s",out));
796 return kFALSE;
797 }
6dc9348d 798 AliTOFCal *cal =(AliTOFCal*)entry->GetObject();
799 fTOFCal = cal;
1b20c168 800 return kTRUE;
801
6dc9348d 802}
803//_____________________________________________________________________________
340693af 804void AliTOFcalib::WriteSimParOnCDB(Char_t *sel, Int_t minrun, Int_t maxrun)
805{
806 //Write Sim miscalibration parameters to the CDB
d4ad0d6b 807
808
809 //for the time being, only one spectrum is used
810 TFile *spFile = new TFile("$ALICE_ROOT/TOF/data/spectrum.root","read");
811 TH1F * hToT;
812 // Retrieve ToT Spectrum
813 spFile->GetObject("ToT",hToT);
814
815 fTOFSimToT=hToT;
816
817 // Retrieve Time over TOT dependence
818
819 TH1F * h = (TH1F*)spFile->Get("TimeToTFit");
820 TList * list = (TList*)h->GetListOfFunctions();
821 TF1* f = (TF1*)list->At(0);
822 Float_t par[6] = {0,0,0,0,0,0};
823 Int_t npar=f->GetNpar();
824 for (Int_t ipar=0;ipar<npar;ipar++){
825 par[ipar]=f->GetParameter(ipar);
826 }
85fc78e3 827 if (!fTOFSimCal) {
828 AliTOFGeometry *geom=new AliTOFGeometryV5();
829 fTOFSimCal = new AliTOFCal(geom);
830 fTOFSimCal->CreateArray();
831 delete geom;
832 }
d4ad0d6b 833 for(Int_t iTOFch=0; iTOFch<fTOFSimCal->NPads();iTOFch++){
340693af 834 AliTOFChannel * calChannel = fTOFSimCal->GetChannel(iTOFch);
835 calChannel->SetSlewPar(par);
d4ad0d6b 836 }
837
838 // Store them in the CDB
839
840 AliCDBManager *man = AliCDBManager::Instance();
841 if(!man->IsDefaultStorageSet())man->SetDefaultStorage("local://$ALICE_ROOT");
842 AliCDBMetaData *md = new AliCDBMetaData();
843 md->SetResponsible("Chiara Zampolli");
28dd10b6 844 Char_t *sel1 = "SimPar" ;
d4ad0d6b 845 Char_t out[100];
846 sprintf(out,"%s/%s",sel,sel1);
847 AliCDBId id1(out,minrun,maxrun);
848 man->Put(fTOFSimCal,id1,md);
28dd10b6 849 Char_t *sel2 = "SimHisto" ;
d4ad0d6b 850 sprintf(out,"%s/%s",sel,sel2);
851 AliCDBId id2(out,minrun,maxrun);
852 man->Put(fTOFSimToT,id2,md);
5eaae242 853 delete md;
d4ad0d6b 854}
855
856//_____________________________________________________________________________
857void AliTOFcalib::WriteSimParOnCDB(Char_t *sel, Int_t minrun, Int_t maxrun, AliTOFCal *cal, TH1F * histo){
340693af 858 //Write Sim miscalibration parameters to the CDB
d4ad0d6b 859
d4ad0d6b 860 fTOFSimToT=histo;
861 fTOFSimCal=cal;
862 AliCDBManager *man = AliCDBManager::Instance();
863 if(!man->IsDefaultStorageSet())man->SetDefaultStorage("local://$ALICE_ROOT");
864 AliCDBMetaData *md = new AliCDBMetaData();
865 md->SetResponsible("Chiara Zampolli");
28dd10b6 866 Char_t *sel1 = "SimPar" ;
d4ad0d6b 867 Char_t out[100];
868 sprintf(out,"%s/%s",sel,sel1);
869 AliCDBId id1(out,minrun,maxrun);
870 man->Put(fTOFSimCal,id1,md);
28dd10b6 871 Char_t *sel2 = "SimHisto" ;
d4ad0d6b 872 sprintf(out,"%s/%s",sel,sel2);
873 AliCDBId id2(out,minrun,maxrun);
874 man->Put(fTOFSimToT,id2,md);
5eaae242 875 delete md;
d4ad0d6b 876}
877//_____________________________________________________________________________
340693af 878void AliTOFcalib::ReadSimParFromCDB(Char_t *sel, Int_t nrun)
879{
880 //Read miscalibration parameters from the CDB
d4ad0d6b 881 AliCDBManager *man = AliCDBManager::Instance();
882 if(!man->IsDefaultStorageSet())man->SetDefaultStorage("local://$ALICE_ROOT");
28dd10b6 883 Char_t *sel1 = "SimPar" ;
d4ad0d6b 884 Char_t out[100];
885 sprintf(out,"%s/%s",sel,sel1);
886 AliCDBEntry *entry1 = man->Get(out,nrun);
887 AliTOFCal *cal =(AliTOFCal*)entry1->GetObject();
888 fTOFSimCal=cal;
28dd10b6 889 Char_t *sel2 = "SimHisto" ;
d4ad0d6b 890 sprintf(out,"%s/%s",sel,sel2);
891 AliCDBEntry *entry2 = man->Get(out,nrun);
892 TH1F *histo =(TH1F*)entry2->GetObject();
893 fTOFSimToT=histo;
894}
895//_____________________________________________________________________________
6dc9348d 896
7037bd93 897//_____________________________________________________________________________
898void AliTOFcalib::WriteRecParOnCDB(Char_t *sel, Int_t minrun, Int_t maxrun, AliTOFRecoParam *param){
899 //Write reconstruction parameters to the CDB
900
901 AliCDBManager *man = AliCDBManager::Instance();
902 if(!man->IsDefaultStorageSet())man->SetDefaultStorage("local://$ALICE_ROOT");
903 AliCDBMetaData *md = new AliCDBMetaData();
904 md->SetResponsible("Silvia Arcelli");
905 Char_t *sel1 = "RecPar" ;
906 Char_t out[100];
907 sprintf(out,"%s/%s",sel,sel1);
908 AliCDBId id(out,minrun,maxrun);
909 man->Put(param,id,md);
910 delete md;
911}
912//_____________________________________________________________________________
913AliTOFRecoParam * AliTOFcalib::ReadRecParFromCDB(Char_t *sel, Int_t nrun)
914{
915 //Read reconstruction parameters from the CDB
916 AliCDBManager *man = AliCDBManager::Instance();
917 if(!man->IsDefaultStorageSet())man->SetDefaultStorage("local://$ALICE_ROOT");
918 Char_t *sel1 = "RecPar" ;
919 Char_t out[100];
920 sprintf(out,"%s/%s",sel,sel1);
921 AliCDBEntry *entry = man->Get(out,nrun);
922 AliTOFRecoParam *param=(AliTOFRecoParam*)entry->GetObject();
923 return param;
924}
925//_____________________________________________________________________________
926
340693af 927Int_t AliTOFcalib::GetIndex(Int_t *detId)
928{
929 //Retrieve calibration channel index
6dc9348d 930 Int_t isector = detId[0];
931 if (isector >= fNSector)
932 AliError(Form("Wrong sector number in TOF (%d) !",isector));
933 Int_t iplate = detId[1];
934 if (iplate >= fNPlate)
935 AliError(Form("Wrong plate number in TOF (%d) !",iplate));
936 Int_t istrip = detId[2];
937 Int_t ipadz = detId[3];
938 Int_t ipadx = detId[4];
939 Int_t stripOffset = 0;
940 switch (iplate) {
941 case 0:
942 stripOffset = 0;
943 break;
944 case 1:
945 stripOffset = fNStripC;
946 break;
947 case 2:
948 stripOffset = fNStripC+fNStripB;
949 break;
950 case 3:
951 stripOffset = fNStripC+fNStripB+fNStripA;
952 break;
953 case 4:
954 stripOffset = fNStripC+fNStripB+fNStripA+fNStripB;
955 break;
956 default:
957 AliError(Form("Wrong plate number in TOF (%d) !",iplate));
958 break;
959 };
960
961 Int_t idet = ((2*(fNStripC+fNStripB)+fNStripA)*fNpadZ*fNpadX)*isector +
962 (stripOffset*fNpadZ*fNpadX)+
963 (fNpadZ*fNpadX)*istrip+
964 (fNpadX)*ipadz+
965 ipadx;
966 return idet;
967}
968