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