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