fb754b6f98ff8da920836f044ae379c4dc1efd23
[u/mrichter/AliRoot.git] / TOF / AliTOFAnalysisTaskCalibPass0.cxx
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
16 ///////////////////////////////////////////////////////////////
17 //                                                           //
18 //   This class provides TOF pass0/passX calibration tools   //
19 //                                                           //
20 ///////////////////////////////////////////////////////////////
21
22 #include "AliTOFAnalysisTaskCalibPass0.h"
23 #include "AliESDtrack.h"
24 #include "AliESDEvent.h"
25 #include "AliPhysicsSelection.h"
26 #include "AliESDtrackCuts.h"
27 #include "AliCDBManager.h"
28 #include "AliGRPManager.h"
29 #include "AliGRPObject.h"
30 #include "AliTOFcalib.h"
31 #include "TList.h"
32 #include "TH2F.h"
33 #include "TFile.h"
34 #include "TH1D.h"
35 #include "TObjArray.h"
36 #include "TString.h"
37 #include "TObjString.h"
38 #include "TROOT.h"
39 #include "TF1.h"
40 #include "TTime.h"
41 #include "TGrid.h"
42 #include "TTimeStamp.h"
43 #include "AliTOFRunParams.h"
44 #include "AliCDBStorage.h"
45 #include "AliCDBId.h"
46 #include "AliCDBMetaData.h"
47 #include "TSystem.h"
48 #include "AliLog.h"
49
50 ClassImp(AliTOFAnalysisTaskCalibPass0)
51   
52 //_______________________________________________________
53
54 const Int_t AliTOFAnalysisTaskCalibPass0::fgkMaxNumberOfPoints = 10000; // max number of points
55 const Double_t AliTOFAnalysisTaskCalibPass0::fgkMinVertexIntegral = 1000.;
56 const Double_t AliTOFAnalysisTaskCalibPass0::fgkMinDeltatIntegral = 20000.;
57
58 //_______________________________________________________
59   
60 AliTOFAnalysisTaskCalibPass0::AliTOFAnalysisTaskCalibPass0() :
61   AliAnalysisTaskSE("TOFCalib-Pass0"),
62   fInitFlag(kFALSE),
63   fEventSelectionFlag(kFALSE),
64   fVertexSelectionFlag(kFALSE),
65   fVertexCut(10.),
66   fRunNumber(0),
67   fESDEvent(NULL),
68   fEventCuts(new AliPhysicsSelection()),
69   fTrackCuts(new AliESDtrackCuts()),
70   fStartTime(0),
71   fEndTime(0),
72   fEventTime(0),
73   fElapsedTime(0),
74   fkVertex(NULL),
75   fGRPManager(new AliGRPManager()),
76   fkGRPObject(NULL),
77   fTOFcalib(new AliTOFcalib()),
78   fHistoList(new TList()),
79   fHistoVertexTimestamp(NULL),
80   fHistoDeltatTimestamp(NULL),
81   fHistoDeltazEta(NULL),
82   fHistoDeltazCosTheta(NULL),
83   fHistoAcceptedTracksEtaPt(NULL),
84   fHistoMatchedTracksEtaPt(NULL)
85 {
86   /* 
87    * default constructor 
88    */
89
90   /* define output */
91   fHistoList->SetOwner(kTRUE);
92   DefineOutput(1, TList::Class());
93 }
94
95 //_______________________________________________________
96
97 AliTOFAnalysisTaskCalibPass0::~AliTOFAnalysisTaskCalibPass0()
98 {
99   /*
100    * default destructor
101    */
102
103   if (fHistoList) delete fHistoList;
104 }
105
106 //_______________________________________________________
107
108 Bool_t
109 AliTOFAnalysisTaskCalibPass0::InitRun()
110 {
111   /*
112    * init run
113    */
114
115   /* get ESD event */
116   fESDEvent = dynamic_cast<AliESDEvent *>(InputEvent());
117   if (!fESDEvent) {
118     AliError("cannot get ESD event");
119     return kFALSE;
120   }
121   /* get run number */
122   Int_t runNb = fESDEvent->GetRunNumber();
123   /* check run already initialized */
124   if (fInitFlag && fRunNumber == runNb) return kTRUE;
125   if (fInitFlag && fRunNumber != runNb) {
126     AliFatal(Form("run number has changed within same job and this is not allowed: %d -> %d", fRunNumber, runNb));
127     return kTRUE;
128   }
129   fInitFlag = kFALSE;
130   /* init cdb if not yet done */
131   AliCDBManager *cdb = AliCDBManager::Instance();
132   if (!cdb->IsDefaultStorageSet())
133     cdb->SetDefaultStorage("raw://");
134   if (cdb->GetRun() < 0)
135     cdb->SetRun(runNb);
136     /* init TOF calib */
137   if (!fTOFcalib->Init()) {
138     AliError("cannot init TOF calib");
139     return kFALSE;
140   }
141   /* get GRP data */
142   if (!fGRPManager->ReadGRPEntry()) {
143     AliError("error while reading \"GRPEntry\" from OCDB");
144     return kFALSE;
145   }
146   fkGRPObject = fGRPManager->GetGRPData();
147   if (!fkGRPObject) {
148     AliError("cannot get \"GRPData\" from GRP manager");
149     return kFALSE;
150   }
151   fStartTime = fkGRPObject->GetTimeStart();
152   fEndTime = fkGRPObject->GetTimeEnd();
153   AliInfo(Form("got \"GRPData\": startTime=%d, endTime=%d", fStartTime, fEndTime));
154   
155   AliInfo(Form("initialized for run %d", runNb));
156   fInitFlag = kTRUE;
157   fRunNumber = runNb;
158   
159   /* set histo title with run-number and start-time */
160   fHistoVertexTimestamp->SetTitle(Form("run: %d, startTimestamp: %u", fRunNumber, fStartTime));
161   fHistoDeltatTimestamp->SetTitle(Form("run: %d, startTimestamp: %u", fRunNumber, fStartTime));
162   fHistoDeltazEta->SetTitle(Form("run: %d, startTimestamp: %u", fRunNumber, fStartTime));
163   fHistoDeltazCosTheta->SetTitle(Form("run: %d, startTimestamp: %u", fRunNumber, fStartTime));
164   fHistoAcceptedTracksEtaPt->SetTitle(Form("run: %d, startTimestamp: %u", fRunNumber, fStartTime));
165   fHistoMatchedTracksEtaPt->SetTitle(Form("run: %d, startTimestamp: %u", fRunNumber, fStartTime));
166   
167   return kTRUE;
168 }
169
170 //_______________________________________________________
171
172 Bool_t
173 AliTOFAnalysisTaskCalibPass0::InitEvent()
174 {
175   /*
176    * init event
177    */
178
179   /* get ESD event */
180   fESDEvent = dynamic_cast<AliESDEvent *>(InputEvent());
181   if (!fESDEvent) return kFALSE;
182   /* set event time and elapsed time */
183   fEventTime = fESDEvent->GetTimeStamp();
184   fElapsedTime = fESDEvent->GetTimeStamp() - fStartTime;
185   /* event selection */
186   if (fEventSelectionFlag && !fEventCuts->IsCollisionCandidate(fESDEvent)) return kFALSE;
187   /* vertex selection */
188   fkVertex = fESDEvent->GetPrimaryVertexTracks();
189   if (fVertexSelectionFlag && fkVertex->GetNContributors() < 1) {
190     fkVertex = fESDEvent->GetPrimaryVertexSPD();
191     if (fkVertex->GetNContributors() < 1) return kFALSE;
192   }
193   if (fVertexSelectionFlag && TMath::Abs(fkVertex->GetZ()) > fVertexCut) return kFALSE;
194   /* calibrate ESD if requested */
195   fTOFcalib->CalibrateESD(fESDEvent);
196
197   return kTRUE;
198 }
199
200 //_______________________________________________________
201
202 Bool_t
203 AliTOFAnalysisTaskCalibPass0::HasTOFMeasurement(const AliESDtrack *track) const
204 {
205   /*
206    * has TOF measurement
207    */
208
209   /* check TOF status flags */
210   if (!(track->GetStatus() & AliESDtrack::kTOFout) ||
211       !(track->GetStatus() & AliESDtrack::kTIME)) return kFALSE;
212   /* check integrated length */
213   if (track->GetIntegratedLength() < 350.) return kFALSE;
214
215   /* TOF measurement ok */
216   return kTRUE;
217 }
218
219 //_______________________________________________________
220
221 void
222 AliTOFAnalysisTaskCalibPass0::UserCreateOutputObjects()
223 {
224   /*
225    * user create output objects
226    */
227
228   /* time binning */
229   Int_t timeBins = 288;
230   Float_t timeMin = 0.;
231   Float_t timeMax = 24. * 3600.;
232   /* vertex binning */
233   Int_t vertexBins = 200;
234   Float_t vertexMin = -50.;
235   Float_t vertexMax = 50.;
236   /* deltat binning */
237   Int_t deltatBins = 2000;
238   Float_t deltatMin = -24400.;
239   Float_t deltatMax = 24400.;
240   /* deltaz binning */
241   Int_t deltazBins = 200;
242   Float_t deltazMin = -10.;
243   Float_t deltazMax = 10.;
244   /* eta binning */
245   Int_t etaBins = 200;
246   Float_t etaMin = -1.;
247   Float_t etaMax = 1.;
248   /* p binning */
249   Int_t pBins = 500;
250   Float_t pMin = 0.;
251   Float_t pMax = 5.;
252   
253   fHistoVertexTimestamp = new TH2F("hHistoVertexTimestamp", "Vertex position;elapsed time (s);z (cm);", timeBins, timeMin, timeMax, vertexBins, vertexMin, vertexMax);
254   fHistoList->Add(fHistoVertexTimestamp);
255
256   fHistoDeltatTimestamp = new TH2F("hHistoDeltatTimestamp", "Global time shift (T0-fill);elapsed time (s);t - t_{exp}^{(#pi)} (ps);", timeBins, timeMin, timeMax, deltatBins, deltatMin, deltatMax);
257   fHistoList->Add(fHistoDeltatTimestamp);
258
259   fHistoDeltazEta = new TH2F("hHistoDeltazEta", "Matching residuals (longitudinal);#eta;#Deltaz (cm);", etaBins, etaMin, etaMax, deltazBins, deltazMin, deltazMax);
260   fHistoList->Add(fHistoDeltazEta);
261
262   fHistoDeltazCosTheta = new TH2F("hHistoDeltazCosTheta", "Matching residuals (longitudinal);cos #theta;#Deltaz (cm);", etaBins, etaMin, etaMax, deltazBins, deltazMin, deltazMax);
263   fHistoList->Add(fHistoDeltazCosTheta);
264
265   fHistoAcceptedTracksEtaPt = new TH2F("hHistoAcceptedTracksEtaPt", "Accepted tracks;#eta;#p_{T} (GeV/c);", etaBins, etaMin, etaMax, pBins, pMin, pMax);
266   fHistoList->Add(fHistoAcceptedTracksEtaPt);
267
268   fHistoMatchedTracksEtaPt = new TH2F("hHistoMatchedTracksEtaPt", "Matched tracks;#eta;p_{T} (GeV/c);", etaBins, etaMin, etaMax, pBins, pMin, pMax);
269   fHistoList->Add(fHistoMatchedTracksEtaPt);
270
271   /* post data */
272   PostData(1, fHistoList);
273 }
274
275 //_______________________________________________________
276
277 void
278 AliTOFAnalysisTaskCalibPass0::UserExec(Option_t *)
279 {
280   /*
281    * user exec
282    */
283
284   /* init run */
285   if (!InitRun()) return;
286   /* init event */
287   if (!InitEvent()) return;
288
289   /*** ACCEPTED EVENT ***/
290   
291   /* fill vertex histo */
292   fHistoVertexTimestamp->Fill(fElapsedTime, fkVertex->GetZ());
293
294   /* loop over ESD tracks */
295   Int_t nTracks = fESDEvent->GetNumberOfTracks();
296   AliESDtrack *track;
297   Double_t eta, costheta, pt, time, timei[AliPID::kSPECIES], deltat, deltaz;
298   for (Int_t itrk = 0; itrk < nTracks; itrk++) {
299     /* get track */
300     track = fESDEvent->GetTrack(itrk);
301     if (!track) continue;
302     /* check accept track */
303     if (!fTrackCuts->AcceptTrack(track)) continue;
304     /* get track info */
305     eta = track->Eta();
306     costheta = TMath::Cos(track->Theta());
307     pt = track->Pt();
308     /* fill accepted tracks histo */
309     fHistoAcceptedTracksEtaPt->Fill(eta, pt);
310     /* check TOF measurement */
311     if (!HasTOFMeasurement(track)) continue;
312
313     /*** ACCEPTED TRACK WITH TOF MEASUREMENT ***/
314
315     /* fill matched tracks histo */
316     fHistoMatchedTracksEtaPt->Fill(eta, pt);
317     /* get TOF info */
318     time = track->GetTOFsignal();
319     track->GetIntegratedTimes(timei);
320     deltat = time - timei[AliPID::kPion];
321     deltaz = track->GetTOFsignalDz();
322     
323     /* fill histos */
324     fHistoDeltatTimestamp->Fill(fElapsedTime, deltat);
325     fHistoDeltazEta->Fill(eta, deltaz);
326     fHistoDeltazCosTheta->Fill(costheta, deltaz);
327     
328   } /* end of loop over ESD tracks */
329
330   /* post data */
331   PostData(1, fHistoList);
332   
333 }
334
335 //_______________________________________________________
336
337 Bool_t
338 AliTOFAnalysisTaskCalibPass0::ProcessOutput(const Char_t *filename, const Char_t *dbString)
339 {
340   /*
341    * process output
342    */
343
344   /* open file */
345   TFile *file = TFile::Open(filename);
346   if (!file || !file->IsOpen()) {
347     AliError(Form("cannot open output file %s", filename));
348     return kFALSE;
349   }
350   /* get histo list */
351   TList *list = (TList *)file->Get("Histos");
352   if (!list) {
353     AliError(Form("cannot get \"Histos\" list from file %s", filename));
354     return kFALSE;
355   }
356   /* get histos */
357   TH2F *histoVertexTimestamp = (TH2F *)list->FindObject("hHistoVertexTimestamp");
358   if (!histoVertexTimestamp) {
359     AliError(Form("cannot get \"hHistoVertexTimestamp\" object from file %s", filename));
360     return kFALSE;
361   }
362   TH2F *histoDeltatTimestamp = (TH2F *)list->FindObject("hHistoDeltatTimestamp");
363   if (!histoDeltatTimestamp) {
364     AliError(Form("cannot get \"hHistoDeltatTimestamp\" object from file %s", filename));
365     return kFALSE;
366   }
367   TH2F *histoDeltazEta = (TH2F *)list->FindObject("hHistoDeltazEta");
368   if (!histoDeltazEta) {
369     AliError(Form("cannot get \"hHistoDeltazEta\" object from file %s", filename));
370     return kFALSE;
371   }
372   TH2F *histoDeltazCosTheta = (TH2F *)list->FindObject("hHistoDeltazCosTheta");
373   if (!histoDeltazCosTheta) {
374     AliError(Form("cannot get \"hHistoDeltazCosTheta\" object from file %s", filename));
375     return kFALSE;
376   }
377   TH2F *histoAcceptedTracksEtaPt = (TH2F *)list->FindObject("hHistoAcceptedTracksEtaPt");
378   if (!histoAcceptedTracksEtaPt) {
379     AliError(Form("cannot get \"hHistoAccptedTracksEtaPt\" object from file %s", filename));
380     return kFALSE;
381   }
382   TH2F *histoMatchedTracksEtaPt = (TH2F *)list->FindObject("hHistoMatchedTracksEtaPt");
383   if (!histoMatchedTracksEtaPt) {
384     AliError(Form("cannot get \"hHistoMatchedTracksEtaPt\" object from file %s", filename));
385     return kFALSE;
386   }
387   /* check matching performance */
388   if (!CheckMatchingPerformance(histoDeltazEta, histoAcceptedTracksEtaPt, histoMatchedTracksEtaPt)) {
389     AliError("error while checking matching efficiency");
390     return kFALSE;
391   }
392   /* calibrate and store */
393   if (!CalibrateAndStore(histoVertexTimestamp, histoDeltatTimestamp, dbString)) {
394     AliError("error while calibrating and storing");
395     return kFALSE;
396   }
397
398   /* success */
399   return kTRUE;
400 }
401
402 //_______________________________________________________
403
404 Bool_t
405 AliTOFAnalysisTaskCalibPass0::CheckMatchingPerformance(const TH2F *histoDeltazCosTheta, const TH2F *histoAcceptedTracksEtaPt, const TH2F *histoMatchedTracksEtaPt) const
406 {
407   /*
408    * check matching performance
409    */
410
411   /* check pointers */
412   if (!histoDeltazCosTheta || !histoAcceptedTracksEtaPt || !histoMatchedTracksEtaPt)
413     return kFALSE;
414   /* dummy for the time being */
415   return kTRUE;
416 }
417
418 //_______________________________________________________
419
420 Bool_t
421 AliTOFAnalysisTaskCalibPass0::CalibrateAndStore(TH2F *histoVertexTimestamp, TH2F *histoDeltatTimestamp, const Char_t *dbString)
422 {
423   /*
424    * calibrate and store
425    */
426
427   /* check pointers */
428   if (!histoVertexTimestamp || !histoDeltatTimestamp)
429     return kFALSE;
430
431   /*** GET RUN-NUMBER AND START-TIMESTAMP ***/
432
433   TString str;
434   TObjArray *strarr = NULL;
435   TObjString *ostr = NULL;
436   str = histoVertexTimestamp->GetTitle();
437   strarr = str.Tokenize(",");
438   if (!strarr) {
439     AliError("problems whith tokenize histogram title");
440     return kFALSE;
441   }
442   
443   /* run number */
444   ostr = (TObjString *)strarr->At(0);
445   if (!ostr) {
446     AliError("problems while getting run number from histogram title");
447     return kFALSE;
448   }
449   str = ostr->GetString();
450   if (!str.BeginsWith("run:")) {
451     AliError("problems while getting run number from histogram title");
452     return kFALSE;
453   }
454   str.Remove(0, 5);
455   Int_t runNb = atoi(str.Data());
456   if (runNb <= 0) {
457     AliError(Form("bad run number: %d", runNb));
458     return kFALSE;
459   }
460   AliInfo(Form("got run number: %d", runNb));
461
462   /* start timestamp */
463   ostr = (TObjString *)strarr->At(1);
464   if (!ostr) {
465     AliError("problems while getting start timestamp from histogram title");
466     return kFALSE;
467   }
468   str = ostr->GetString();
469   str.Remove(0, 1); /* remove empty space at the beginning */
470   if (!str.BeginsWith("startTimestamp:")) {
471     AliError("problems while getting start timestamp from histogram title");
472     return kFALSE;
473   }
474   str.Remove(0, 16);
475   UInt_t startTimestamp = atoi(str.Data());
476   if (startTimestamp <= 0) {
477     AliError(Form("bad start timestamp: %d", startTimestamp));
478     return kFALSE;
479   }
480   TTimeStamp ts = startTimestamp;
481   AliInfo(Form("got start timestamp: %d (%s)", startTimestamp, ts.AsString()));
482
483   /*** CALIBRATION STAGE ***/
484
485   /* get fit function */
486   TF1 *fitFunc = (TF1 *)gROOT->GetFunction("gaus");
487
488   /* projection-x */
489   TH1D *histoVertexTimestamppx = histoVertexTimestamp->ProjectionX("histoVertexTimestamppx");
490   TH1D *histoDeltatTimestamppx = histoDeltatTimestamp->ProjectionX("histoDeltatTimestamppx");
491
492   /* define mix and max time bin */
493   Int_t minBin = histoVertexTimestamppx->FindFirstBinAbove(0);
494   Int_t maxBin = histoVertexTimestamppx->FindLastBinAbove(0);
495   Float_t minTime = histoVertexTimestamppx->GetBinLowEdge(minBin);
496   Float_t maxTime = histoVertexTimestamppx->GetBinLowEdge(maxBin + 1);
497   AliInfo(Form("min/max time defined: %d < t < %d s [%d, %d]", (Int_t)minTime, (Int_t)maxTime, minBin, maxBin));
498
499   /* loop over time bins */
500   Int_t nPoints = 0;
501   Float_t time[fgkMaxNumberOfPoints], timeerr[fgkMaxNumberOfPoints];
502   Float_t vertexMean[fgkMaxNumberOfPoints], vertexMeanerr[fgkMaxNumberOfPoints];
503   Float_t vertexSigma[fgkMaxNumberOfPoints], vertexSigmaerr[fgkMaxNumberOfPoints];
504   Float_t timeZeroMean[fgkMaxNumberOfPoints], timeZeroMeanerr[fgkMaxNumberOfPoints];
505   Float_t timeZeroSigma[fgkMaxNumberOfPoints], timeZeroSigmaerr[fgkMaxNumberOfPoints];
506   Float_t averageTimeZero = 0.;
507   Float_t averageTimeSigma = 0.;
508   Float_t averageVertexSpread = 0.;
509   for (Int_t ibin = minBin; ibin <= maxBin; ibin++) {
510
511     /* define time window */
512     Int_t startBin = ibin;
513     Int_t endBin = ibin;
514     while(histoVertexTimestamppx->Integral(startBin, endBin) < fgkMinVertexIntegral ||
515           histoDeltatTimestamppx->Integral(startBin, endBin) < fgkMinDeltatIntegral) {
516       if (endBin < maxBin) endBin++;
517       else if (startBin > minBin) startBin--;
518       else break;
519     }
520     if (histoVertexTimestamppx->Integral(startBin, endBin) <= 0 ||
521         histoDeltatTimestamppx->Integral(startBin, endBin) <= 0) continue;
522     Float_t startTime = histoVertexTimestamppx->GetBinLowEdge(startBin);
523     Float_t endTime = histoVertexTimestamppx->GetBinLowEdge(endBin + 1);
524     Float_t vertexIntegral = histoVertexTimestamppx->Integral(startBin, endBin);
525     Float_t deltatIntegral = histoDeltatTimestamppx->Integral(startBin, endBin);
526     AliInfo(Form("time window defined: %d < t < %d s [%d, %d]: %d vertices, %d tracks", (Int_t)startTime, (Int_t)endTime, startBin, endBin, (Int_t)vertexIntegral, (Int_t)deltatIntegral));
527
528     /* projection-y */
529     TH1D *histoVertexTimestamppy = histoVertexTimestamp->ProjectionY("histoVertexTimestamppy", startBin, endBin);
530     TH1D *histoDeltatTimestamppy = histoDeltatTimestamp->ProjectionY("histoDeltatTimestamppy", startBin, endBin);
531
532     /* average time */
533     histoVertexTimestamppx->GetXaxis()->SetRange(startBin, endBin);
534     time[nPoints] = histoVertexTimestamppx->GetMean();
535     timeerr[nPoints] = histoVertexTimestamppx->GetMeanError();
536
537     /* fit vertex */
538     if (FitPeak(fitFunc, histoVertexTimestamppy, 10., 3., 3.) != 0) {
539       AliError("troubles fitting vertex, skip");
540       delete histoVertexTimestamppy;
541       delete histoDeltatTimestamppy;
542       continue;
543     }
544     vertexMean[nPoints] = fitFunc->GetParameter(1);
545     vertexMeanerr[nPoints] = fitFunc->GetParError(1);
546     vertexSigma[nPoints] = fitFunc->GetParameter(2);
547     vertexSigmaerr[nPoints] = fitFunc->GetParError(2);
548     averageVertexSpread += fitFunc->GetParameter(2);
549
550     /* fit time-zero */
551     if (FitPeak(fitFunc, histoDeltatTimestamppy, 500., 2., 1.) != 0) {
552       AliError("troubles fitting time-zero TRACKS, skip");
553       delete histoVertexTimestamppy;
554       delete histoDeltatTimestamppy;
555       continue;
556     }
557     timeZeroMean[nPoints] = fitFunc->GetParameter(1);
558     timeZeroMeanerr[nPoints] = fitFunc->GetParError(1);
559     timeZeroSigma[nPoints] = fitFunc->GetParameter(2);
560     timeZeroSigmaerr[nPoints] = fitFunc->GetParError(2);
561     averageTimeZero += fitFunc->GetParameter(1);
562     averageTimeSigma += fitFunc->GetParameter(2);
563
564     /* delete projection-y */
565     delete histoVertexTimestamppy;
566     delete histoDeltatTimestamppy;
567
568     /* increment n points */
569     nPoints++;
570
571     /* set current bin */
572     ibin = endBin;
573
574   } /* end of loop over time bins */
575
576   /* delete projection-x */
577   delete histoVertexTimestamppx;
578   delete histoDeltatTimestamppx;
579
580   /* check points */
581   if (nPoints <= 0) {
582     AliError("no measurement available, quit");
583     return kFALSE;
584   }
585   AliInfo(Form("average time-zero  = %f", averageTimeZero / nPoints));
586   AliInfo(Form("average time-sigma = %f", averageTimeSigma / nPoints));
587   AliInfo(Form("average v-spread   = %f", averageVertexSpread / nPoints));
588
589   /*** CREATE RUN PARAMS OBJECT ***/
590
591 #if 0
592   /* get start time from GRP */
593   TGrid::Connect("alien://", 0, 0, "t");
594   AliCDBManager *cdb = AliCDBManager::Instance();
595   cdb->SetDefaultStorage("raw://");
596   cdb->SetRun(runNb);
597   AliGRPManager grp;
598   if (!grp.ReadGRPEntry()) {
599     AliError("error while reading GRP entry");
600     return kFALSE;
601   }
602   UInt_t startTimestamp = grp.GetGRPData()->GetTimeStart();
603   TTimeStamp ts;
604   ts = startTimestamp;
605   AliInfo(Form("got start time from GRP: %s", ts.AsString()));
606 #endif
607     
608   /* create arrays */
609   UInt_t timestamp[fgkMaxNumberOfPoints];
610   Float_t t0[fgkMaxNumberOfPoints];
611   Float_t tofReso[fgkMaxNumberOfPoints];
612   Float_t t0Spread[fgkMaxNumberOfPoints];
613   for (Int_t ipoint = 0; ipoint < nPoints; ipoint++) {
614     timestamp[ipoint] = time[ipoint] + (Float_t)startTimestamp;
615     t0[ipoint] = timeZeroMean[ipoint];
616     tofReso[ipoint] = timeZeroSigma[ipoint];
617     t0Spread[ipoint] = vertexSigma[ipoint] / 2.99792457999999984e-02;
618   }
619   UInt_t run[1] = {runNb};
620   UInt_t runFirstPoint[1] = {0};
621   UInt_t runLastPoint[1] = {nPoints - 1};
622   
623   /* create run params object */
624   AliTOFRunParams obj(nPoints, 1);
625   AliInfo(Form("create run params object for run %d with %d points", runNb, nPoints));
626   obj.SetTimestamp(timestamp);
627   obj.SetT0(t0);
628   obj.SetTOFResolution(tofReso);
629   obj.SetT0Spread(t0Spread);
630   obj.SetRunNb(run);
631   obj.SetRunFirstPoint(runFirstPoint);
632   obj.SetRunLastPoint(runLastPoint);
633   
634   /*** CREATE OCDB ENTRY ***/
635
636   if (!dbString) {
637     AliError("cannot store object because of NULL string");
638     return kFALSE;
639   }
640
641   /* install run params object in OCDB */
642   AliCDBManager *cdb = AliCDBManager::Instance();
643   AliCDBStorage *sto = cdb->GetStorage(dbString);
644   if (!sto) {
645     AliError(Form("cannot get storage %s", dbString));
646     return kFALSE;
647   }
648   AliCDBId id("TOF/Calib/RunParams", runNb, runNb);
649   AliCDBMetaData md;
650   md.SetResponsible("Roberto Preghenella");
651   md.SetComment("offline TOF run parameters");
652   md.SetAliRootVersion(gSystem->Getenv("ARVERSION"));
653   md.SetBeamPeriod(0);
654   if (!sto->Put(&obj, id, &md)) {
655     AliError(Form("error while putting object in storage %s", dbString));
656   }
657
658   /* success */
659   return kTRUE;
660 }
661
662 //_____________________________________________________________
663
664 Int_t
665 AliTOFAnalysisTaskCalibPass0::FitPeak(TF1 *fitFunc, TH1D *h, Float_t startSigma, Float_t nSigmaMin, Float_t nSigmaMax)
666 {
667   /*
668    * fit peak
669    */
670
671   Double_t fitCent = h->GetBinCenter(h->GetMaximumBin());
672   Double_t fitMin = fitCent - nSigmaMin * startSigma;
673   Double_t fitMax = fitCent + nSigmaMax * startSigma;
674   if (fitMin < h->GetXaxis()->GetXmin()) fitMin = h->GetXaxis()->GetXmin();
675   if (fitMax > h->GetXaxis()->GetXmax()) fitMax = h->GetXaxis()->GetXmax();
676   fitFunc->SetParameter(1, fitCent);
677   fitFunc->SetParameter(2, startSigma);
678   Int_t fitres = h->Fit(fitFunc, "WWq0", "", fitMin, fitMax);
679   if (fitres != 0) return fitres;
680   /* refit with better range */
681   for (Int_t i = 0; i < 3; i++) {
682     fitCent = fitFunc->GetParameter(1);
683     fitMin = fitCent - nSigmaMin * fitFunc->GetParameter(2);
684     fitMax = fitCent + nSigmaMax * fitFunc->GetParameter(2);
685     if (fitMin < h->GetXaxis()->GetXmin()) fitMin = h->GetXaxis()->GetXmin();
686     if (fitMax > h->GetXaxis()->GetXmax()) fitMax = h->GetXaxis()->GetXmax();
687     fitres = h->Fit(fitFunc, "q0", "", fitMin, fitMax);
688     if (fitres != 0) return fitres;
689   }
690   return fitres;
691 }