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