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