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