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