ProcessOutput update: deal with both TFileMerged merged data (histos inside TList...
[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
55const Double_t AliTOFAnalysisTaskCalibPass0::fgkMinVertexIntegral = 1000.;
56const Double_t AliTOFAnalysisTaskCalibPass0::fgkMinDeltatIntegral = 20000.;
57
58//_______________________________________________________
59
60AliTOFAnalysisTaskCalibPass0::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
97AliTOFAnalysisTaskCalibPass0::~AliTOFAnalysisTaskCalibPass0()
98{
99 /*
100 * default destructor
101 */
102
103 if (fHistoList) delete fHistoList;
104}
105
106//_______________________________________________________
107
108Bool_t
109AliTOFAnalysisTaskCalibPass0::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
172Bool_t
173AliTOFAnalysisTaskCalibPass0::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
202Bool_t
203AliTOFAnalysisTaskCalibPass0::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
221void
222AliTOFAnalysisTaskCalibPass0::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
277void
278AliTOFAnalysisTaskCalibPass0::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
337Bool_t
338AliTOFAnalysisTaskCalibPass0::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 }
36c9ca5c 350 /* get histograms */
3c8efc07 351 TList *list = (TList *)file->Get("Histos");
36c9ca5c 352 TH2F *histoVertexTimestamp = NULL;
353 TH2F *histoDeltatTimestamp = NULL;
354 TH2F *histoDeltazEta = NULL;
355 TH2F *histoDeltazCosTheta = NULL;
356 TH2F *histoAcceptedTracksEtaPt = NULL;
357 TH2F *histoMatchedTracksEtaPt = NULL;
358 if (list) {
359 AliInfo(Form("getting histograms from \"Histos\" list from file %s", filename));
360 histoVertexTimestamp = (TH2F *)list->FindObject("hHistoVertexTimestamp");
361 histoDeltatTimestamp = (TH2F *)list->FindObject("hHistoDeltatTimestamp");
362 histoDeltazEta = (TH2F *)list->FindObject("hHistoDeltazEta");
363 histoDeltazCosTheta = (TH2F *)list->FindObject("hHistoDeltazCosTheta");
364 histoAcceptedTracksEtaPt = (TH2F *)list->FindObject("hHistoAcceptedTracksEtaPt");
365 histoMatchedTracksEtaPt = (TH2F *)list->FindObject("hHistoMatchedTracksEtaPt");
366 }
367 else {
368 AliInfo(Form("getting histograms directly from file %s", filename));
369 histoVertexTimestamp = (TH2F *)file->Get("hHistoVertexTimestamp");
370 histoDeltatTimestamp = (TH2F *)file->Get("hHistoDeltatTimestamp");
371 histoDeltazEta = (TH2F *)file->Get("hHistoDeltazEta");
372 histoDeltazCosTheta = (TH2F *)file->Get("hHistoDeltazCosTheta");
373 histoAcceptedTracksEtaPt = (TH2F *)file->Get("hHistoAcceptedTracksEtaPt");
374 histoMatchedTracksEtaPt = (TH2F *)file->Get("hHistoMatchedTracksEtaPt");
3c8efc07 375 }
36c9ca5c 376 /* check histos */
3c8efc07 377 if (!histoVertexTimestamp) {
378 AliError(Form("cannot get \"hHistoVertexTimestamp\" object from file %s", filename));
379 return kFALSE;
380 }
3c8efc07 381 if (!histoDeltatTimestamp) {
382 AliError(Form("cannot get \"hHistoDeltatTimestamp\" object from file %s", filename));
383 return kFALSE;
384 }
3c8efc07 385 if (!histoDeltazEta) {
386 AliError(Form("cannot get \"hHistoDeltazEta\" object from file %s", filename));
387 return kFALSE;
388 }
3c8efc07 389 if (!histoDeltazCosTheta) {
390 AliError(Form("cannot get \"hHistoDeltazCosTheta\" object from file %s", filename));
391 return kFALSE;
392 }
3c8efc07 393 if (!histoAcceptedTracksEtaPt) {
394 AliError(Form("cannot get \"hHistoAccptedTracksEtaPt\" object from file %s", filename));
395 return kFALSE;
396 }
3c8efc07 397 if (!histoMatchedTracksEtaPt) {
398 AliError(Form("cannot get \"hHistoMatchedTracksEtaPt\" object from file %s", filename));
399 return kFALSE;
400 }
401 /* check matching performance */
402 if (!CheckMatchingPerformance(histoDeltazEta, histoAcceptedTracksEtaPt, histoMatchedTracksEtaPt)) {
403 AliError("error while checking matching efficiency");
404 return kFALSE;
405 }
406 /* calibrate and store */
407 if (!CalibrateAndStore(histoVertexTimestamp, histoDeltatTimestamp, dbString)) {
408 AliError("error while calibrating and storing");
409 return kFALSE;
410 }
411
412 /* success */
413 return kTRUE;
414}
415
416//_______________________________________________________
417
418Bool_t
419AliTOFAnalysisTaskCalibPass0::CheckMatchingPerformance(const TH2F *histoDeltazCosTheta, const TH2F *histoAcceptedTracksEtaPt, const TH2F *histoMatchedTracksEtaPt) const
420{
421 /*
422 * check matching performance
423 */
424
425 /* check pointers */
426 if (!histoDeltazCosTheta || !histoAcceptedTracksEtaPt || !histoMatchedTracksEtaPt)
427 return kFALSE;
428 /* dummy for the time being */
429 return kTRUE;
430}
431
432//_______________________________________________________
433
434Bool_t
435AliTOFAnalysisTaskCalibPass0::CalibrateAndStore(TH2F *histoVertexTimestamp, TH2F *histoDeltatTimestamp, const Char_t *dbString)
436{
437 /*
438 * calibrate and store
439 */
440
441 /* check pointers */
442 if (!histoVertexTimestamp || !histoDeltatTimestamp)
443 return kFALSE;
444
445 /*** GET RUN-NUMBER AND START-TIMESTAMP ***/
446
447 TString str;
448 TObjArray *strarr = NULL;
449 TObjString *ostr = NULL;
450 str = histoVertexTimestamp->GetTitle();
451 strarr = str.Tokenize(",");
452 if (!strarr) {
453 AliError("problems whith tokenize histogram title");
454 return kFALSE;
455 }
456
457 /* run number */
458 ostr = (TObjString *)strarr->At(0);
459 if (!ostr) {
460 AliError("problems while getting run number from histogram title");
461 return kFALSE;
462 }
463 str = ostr->GetString();
464 if (!str.BeginsWith("run:")) {
465 AliError("problems while getting run number from histogram title");
466 return kFALSE;
467 }
468 str.Remove(0, 5);
469 Int_t runNb = atoi(str.Data());
470 if (runNb <= 0) {
471 AliError(Form("bad run number: %d", runNb));
472 return kFALSE;
473 }
474 AliInfo(Form("got run number: %d", runNb));
475
476 /* start timestamp */
477 ostr = (TObjString *)strarr->At(1);
478 if (!ostr) {
479 AliError("problems while getting start timestamp from histogram title");
480 return kFALSE;
481 }
482 str = ostr->GetString();
483 str.Remove(0, 1); /* remove empty space at the beginning */
484 if (!str.BeginsWith("startTimestamp:")) {
485 AliError("problems while getting start timestamp from histogram title");
486 return kFALSE;
487 }
488 str.Remove(0, 16);
489 UInt_t startTimestamp = atoi(str.Data());
490 if (startTimestamp <= 0) {
491 AliError(Form("bad start timestamp: %d", startTimestamp));
492 return kFALSE;
493 }
494 TTimeStamp ts = startTimestamp;
495 AliInfo(Form("got start timestamp: %d (%s)", startTimestamp, ts.AsString()));
496
497 /*** CALIBRATION STAGE ***/
498
499 /* get fit function */
500 TF1 *fitFunc = (TF1 *)gROOT->GetFunction("gaus");
501
502 /* projection-x */
503 TH1D *histoVertexTimestamppx = histoVertexTimestamp->ProjectionX("histoVertexTimestamppx");
504 TH1D *histoDeltatTimestamppx = histoDeltatTimestamp->ProjectionX("histoDeltatTimestamppx");
505
506 /* define mix and max time bin */
507 Int_t minBin = histoVertexTimestamppx->FindFirstBinAbove(0);
508 Int_t maxBin = histoVertexTimestamppx->FindLastBinAbove(0);
509 Float_t minTime = histoVertexTimestamppx->GetBinLowEdge(minBin);
510 Float_t maxTime = histoVertexTimestamppx->GetBinLowEdge(maxBin + 1);
511 AliInfo(Form("min/max time defined: %d < t < %d s [%d, %d]", (Int_t)minTime, (Int_t)maxTime, minBin, maxBin));
512
513 /* loop over time bins */
514 Int_t nPoints = 0;
515 Float_t time[fgkMaxNumberOfPoints], timeerr[fgkMaxNumberOfPoints];
516 Float_t vertexMean[fgkMaxNumberOfPoints], vertexMeanerr[fgkMaxNumberOfPoints];
517 Float_t vertexSigma[fgkMaxNumberOfPoints], vertexSigmaerr[fgkMaxNumberOfPoints];
518 Float_t timeZeroMean[fgkMaxNumberOfPoints], timeZeroMeanerr[fgkMaxNumberOfPoints];
519 Float_t timeZeroSigma[fgkMaxNumberOfPoints], timeZeroSigmaerr[fgkMaxNumberOfPoints];
520 Float_t averageTimeZero = 0.;
521 Float_t averageTimeSigma = 0.;
522 Float_t averageVertexSpread = 0.;
523 for (Int_t ibin = minBin; ibin <= maxBin; ibin++) {
524
525 /* define time window */
526 Int_t startBin = ibin;
527 Int_t endBin = ibin;
528 while(histoVertexTimestamppx->Integral(startBin, endBin) < fgkMinVertexIntegral ||
529 histoDeltatTimestamppx->Integral(startBin, endBin) < fgkMinDeltatIntegral) {
530 if (endBin < maxBin) endBin++;
531 else if (startBin > minBin) startBin--;
532 else break;
533 }
534 if (histoVertexTimestamppx->Integral(startBin, endBin) <= 0 ||
535 histoDeltatTimestamppx->Integral(startBin, endBin) <= 0) continue;
536 Float_t startTime = histoVertexTimestamppx->GetBinLowEdge(startBin);
537 Float_t endTime = histoVertexTimestamppx->GetBinLowEdge(endBin + 1);
538 Float_t vertexIntegral = histoVertexTimestamppx->Integral(startBin, endBin);
539 Float_t deltatIntegral = histoDeltatTimestamppx->Integral(startBin, endBin);
540 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));
541
542 /* projection-y */
543 TH1D *histoVertexTimestamppy = histoVertexTimestamp->ProjectionY("histoVertexTimestamppy", startBin, endBin);
544 TH1D *histoDeltatTimestamppy = histoDeltatTimestamp->ProjectionY("histoDeltatTimestamppy", startBin, endBin);
545
546 /* average time */
547 histoVertexTimestamppx->GetXaxis()->SetRange(startBin, endBin);
548 time[nPoints] = histoVertexTimestamppx->GetMean();
549 timeerr[nPoints] = histoVertexTimestamppx->GetMeanError();
550
551 /* fit vertex */
552 if (FitPeak(fitFunc, histoVertexTimestamppy, 10., 3., 3.) != 0) {
553 AliError("troubles fitting vertex, skip");
554 delete histoVertexTimestamppy;
555 delete histoDeltatTimestamppy;
556 continue;
557 }
558 vertexMean[nPoints] = fitFunc->GetParameter(1);
559 vertexMeanerr[nPoints] = fitFunc->GetParError(1);
560 vertexSigma[nPoints] = fitFunc->GetParameter(2);
561 vertexSigmaerr[nPoints] = fitFunc->GetParError(2);
562 averageVertexSpread += fitFunc->GetParameter(2);
563
564 /* fit time-zero */
565 if (FitPeak(fitFunc, histoDeltatTimestamppy, 500., 2., 1.) != 0) {
566 AliError("troubles fitting time-zero TRACKS, skip");
567 delete histoVertexTimestamppy;
568 delete histoDeltatTimestamppy;
569 continue;
570 }
571 timeZeroMean[nPoints] = fitFunc->GetParameter(1);
572 timeZeroMeanerr[nPoints] = fitFunc->GetParError(1);
573 timeZeroSigma[nPoints] = fitFunc->GetParameter(2);
574 timeZeroSigmaerr[nPoints] = fitFunc->GetParError(2);
575 averageTimeZero += fitFunc->GetParameter(1);
576 averageTimeSigma += fitFunc->GetParameter(2);
577
578 /* delete projection-y */
579 delete histoVertexTimestamppy;
580 delete histoDeltatTimestamppy;
581
582 /* increment n points */
583 nPoints++;
584
585 /* set current bin */
586 ibin = endBin;
587
588 } /* end of loop over time bins */
589
590 /* delete projection-x */
591 delete histoVertexTimestamppx;
592 delete histoDeltatTimestamppx;
593
594 /* check points */
595 if (nPoints <= 0) {
596 AliError("no measurement available, quit");
597 return kFALSE;
598 }
599 AliInfo(Form("average time-zero = %f", averageTimeZero / nPoints));
600 AliInfo(Form("average time-sigma = %f", averageTimeSigma / nPoints));
601 AliInfo(Form("average v-spread = %f", averageVertexSpread / nPoints));
602
603 /*** CREATE RUN PARAMS OBJECT ***/
604
605#if 0
606 /* get start time from GRP */
607 TGrid::Connect("alien://", 0, 0, "t");
608 AliCDBManager *cdb = AliCDBManager::Instance();
609 cdb->SetDefaultStorage("raw://");
610 cdb->SetRun(runNb);
611 AliGRPManager grp;
612 if (!grp.ReadGRPEntry()) {
613 AliError("error while reading GRP entry");
614 return kFALSE;
615 }
616 UInt_t startTimestamp = grp.GetGRPData()->GetTimeStart();
617 TTimeStamp ts;
618 ts = startTimestamp;
619 AliInfo(Form("got start time from GRP: %s", ts.AsString()));
620#endif
621
622 /* create arrays */
623 UInt_t timestamp[fgkMaxNumberOfPoints];
624 Float_t t0[fgkMaxNumberOfPoints];
625 Float_t tofReso[fgkMaxNumberOfPoints];
626 Float_t t0Spread[fgkMaxNumberOfPoints];
627 for (Int_t ipoint = 0; ipoint < nPoints; ipoint++) {
7ca3d6d6 628 timestamp[ipoint] = (UInt_t)time[ipoint] + startTimestamp;
3c8efc07 629 t0[ipoint] = timeZeroMean[ipoint];
630 tofReso[ipoint] = timeZeroSigma[ipoint];
631 t0Spread[ipoint] = vertexSigma[ipoint] / 2.99792457999999984e-02;
632 }
633 UInt_t run[1] = {runNb};
634 UInt_t runFirstPoint[1] = {0};
635 UInt_t runLastPoint[1] = {nPoints - 1};
636
637 /* create run params object */
638 AliTOFRunParams obj(nPoints, 1);
639 AliInfo(Form("create run params object for run %d with %d points", runNb, nPoints));
640 obj.SetTimestamp(timestamp);
641 obj.SetT0(t0);
642 obj.SetTOFResolution(tofReso);
643 obj.SetT0Spread(t0Spread);
644 obj.SetRunNb(run);
645 obj.SetRunFirstPoint(runFirstPoint);
646 obj.SetRunLastPoint(runLastPoint);
647
648 /*** CREATE OCDB ENTRY ***/
649
650 if (!dbString) {
651 AliError("cannot store object because of NULL string");
652 return kFALSE;
653 }
654
655 /* install run params object in OCDB */
656 AliCDBManager *cdb = AliCDBManager::Instance();
657 AliCDBStorage *sto = cdb->GetStorage(dbString);
658 if (!sto) {
659 AliError(Form("cannot get storage %s", dbString));
660 return kFALSE;
661 }
662 AliCDBId id("TOF/Calib/RunParams", runNb, runNb);
663 AliCDBMetaData md;
664 md.SetResponsible("Roberto Preghenella");
665 md.SetComment("offline TOF run parameters");
666 md.SetAliRootVersion(gSystem->Getenv("ARVERSION"));
667 md.SetBeamPeriod(0);
668 if (!sto->Put(&obj, id, &md)) {
669 AliError(Form("error while putting object in storage %s", dbString));
670 }
671
672 /* success */
673 return kTRUE;
674}
675
676//_____________________________________________________________
677
678Int_t
679AliTOFAnalysisTaskCalibPass0::FitPeak(TF1 *fitFunc, TH1D *h, Float_t startSigma, Float_t nSigmaMin, Float_t nSigmaMax)
680{
681 /*
682 * fit peak
683 */
684
685 Double_t fitCent = h->GetBinCenter(h->GetMaximumBin());
686 Double_t fitMin = fitCent - nSigmaMin * startSigma;
687 Double_t fitMax = fitCent + nSigmaMax * startSigma;
688 if (fitMin < h->GetXaxis()->GetXmin()) fitMin = h->GetXaxis()->GetXmin();
689 if (fitMax > h->GetXaxis()->GetXmax()) fitMax = h->GetXaxis()->GetXmax();
690 fitFunc->SetParameter(1, fitCent);
691 fitFunc->SetParameter(2, startSigma);
692 Int_t fitres = h->Fit(fitFunc, "WWq0", "", fitMin, fitMax);
693 if (fitres != 0) return fitres;
694 /* refit with better range */
695 for (Int_t i = 0; i < 3; i++) {
696 fitCent = fitFunc->GetParameter(1);
697 fitMin = fitCent - nSigmaMin * fitFunc->GetParameter(2);
698 fitMax = fitCent + nSigmaMax * fitFunc->GetParameter(2);
699 if (fitMin < h->GetXaxis()->GetXmin()) fitMin = h->GetXaxis()->GetXmin();
700 if (fitMax > h->GetXaxis()->GetXmax()) fitMax = h->GetXaxis()->GetXmax();
701 fitres = h->Fit(fitFunc, "q0", "", fitMin, fitMax);
702 if (fitres != 0) return fitres;
703 }
704 return fitres;
705}