Allow writing in OCDB TimeDelay if less then 2 PMTs are empty
[u/mrichter/AliRoot.git] / T0 / AliT0CalibTimeEq.cxx
CommitLineData
18086258 1
a2ad8166 2/**************************************************************************
3 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * *
5 * Author: The ALICE Off-line Project. *
6 * Contributors are mentioned in the code where appropriate. *
7 * *
8 * Permission to use, copy, modify and distribute this software and its *
9 * documentation strictly for non-commercial purposes is hereby granted *
10 * without fee, provided that the above copyright notice appears in all *
11 * copies and that both the copyright notice and this permission notice *
12 * appear in the supporting documentation. The authors make no claims *
13 * about the suitability of this software for any purpose. It is *
14 * provided "as is" without express or implied warranty. *
15 **************************************************************************/
16
17/* $Id$ */
18
19///////////////////////////////////////////////////////////////////////////////
20// //
59b8a5fc 21// class for T0 calibration TM-AC-AM_6-02-2006
22// equalize time shift for each time CFD channel
a2ad8166 23// //
24///////////////////////////////////////////////////////////////////////////////
25
26#include "AliT0CalibTimeEq.h"
e47b556f 27#include "AliLog.h"
a2ad8166 28#include <TFile.h>
29#include <TMath.h>
30#include <TF1.h>
a2ad8166 31#include <TSpectrum.h>
a2ad8166 32#include <TProfile.h>
c61a7285 33#include <iostream>
a2ad8166 34
35ClassImp(AliT0CalibTimeEq)
36
37//________________________________________________________________
eef60a6a 38 AliT0CalibTimeEq::AliT0CalibTimeEq():TNamed(),
39 fMeanVertex(0),
40 fRmsVertex(0)
a2ad8166 41{
42 //
57b90263 43 for(Int_t i=0; i<24; i++) {
44 fTimeEq[i] = 0; // Time Equalized for OCDB
63c69eb3 45 fTimeEqRms[i] = -1; // RMS of Time Equalized for OCDB
57b90263 46 for (Int_t ih=0; ih<5; ih++) fCFDvalue[i][ih] = 0;
47 }
a2ad8166 48}
49
50//________________________________________________________________
eef60a6a 51AliT0CalibTimeEq::AliT0CalibTimeEq(const char* name):TNamed(),
52 fMeanVertex(0),
53 fRmsVertex(0)
a2ad8166 54{
59b8a5fc 55 //constructor
56
a2ad8166 57 TString namst = "Calib_";
58 namst += name;
59 SetName(namst.Data());
60 SetTitle(namst.Data());
57b90263 61 for(Int_t i=0; i<24; i++) {
62 fTimeEq[i] = 0; // Time Equalized for OCDB
63c69eb3 63 fTimeEqRms[i] = -1; // RMS of Time Equalized for OCDB
57b90263 64 for (Int_t ih=0; ih<5; ih++) fCFDvalue[i][ih] = 0;
65 }
a2ad8166 66}
67
68//________________________________________________________________
eef60a6a 69AliT0CalibTimeEq::AliT0CalibTimeEq(const AliT0CalibTimeEq& calibda):TNamed(calibda),
70 fMeanVertex(0),
71 fRmsVertex(0)
a2ad8166 72{
73// copy constructor
74 SetName(calibda.GetName());
75 SetTitle(calibda.GetName());
57b90263 76 ((AliT0CalibTimeEq &) calibda).Copy(*this);
a2ad8166 77
78}
79
80//________________________________________________________________
81AliT0CalibTimeEq &AliT0CalibTimeEq::operator =(const AliT0CalibTimeEq& calibda)
82{
83// assignment operator
84 SetName(calibda.GetName());
85 SetTitle(calibda.GetName());
86
57b90263 87 if (this != &calibda) (( AliT0CalibTimeEq &) calibda).Copy(*this);
a2ad8166 88 return *this;
89}
90
91//________________________________________________________________
92AliT0CalibTimeEq::~AliT0CalibTimeEq()
93{
94 //
59b8a5fc 95 // destrictor
a2ad8166 96}
97//________________________________________________________________
98void AliT0CalibTimeEq::Reset()
99{
59b8a5fc 100 //reset values
101
a2ad8166 102 memset(fCFDvalue,0,120*sizeof(Float_t));
103 memset(fTimeEq,1,24*sizeof(Float_t));
104}
105
106
107//________________________________________________________________
108void AliT0CalibTimeEq::Print(Option_t*) const
109{
59b8a5fc 110 // print time values
a2ad8166 111
112 printf("\n ---- PM Arrays ----\n\n");
113 printf(" Time delay CFD \n");
2a1ccb90 114 for (Int_t i=0; i<24; i++)
115 printf(" CFD %f diff %f \n",fCFDvalue[i][0],fTimeEq[i]);
a2ad8166 116}
117
118
119//________________________________________________________________
02888358 120Bool_t AliT0CalibTimeEq::ComputeOnlineParams(const char* filePhys)
a2ad8166 121{
59b8a5fc 122 // compute online equalized time
1e7108ce 123 Float_t meandiff, sigmadiff, meanver, meancfdtime, sigmacfdtime;
c4ba15c3 124 meandiff = sigmadiff = meanver = meancfdtime = sigmacfdtime =0;
c4ba15c3 125 Double_t rmsver=0;
b5e0299c 126 Int_t nent=0;
32d3e409 127 Int_t okdiff=0;
128 Int_t oktime=0;
02888358 129 Bool_t ok=false;
61bd03f6 130 //
c2337900 131 gFile = TFile::Open(filePhys);
370867d2 132 if(!gFile) {
6318dd46 133 AliError("No input PHYS data found ");
134 }
135 else
136 {
02888358 137 ok=true;
6318dd46 138 for (Int_t i=0; i<24; i++)
139 {
e691fec9 140 meandiff = sigmadiff = meanver = meancfdtime = sigmacfdtime =0;
61bd03f6 141 TH1F *cfd = (TH1F*) gFile->Get(Form("CFD1minCFD%d",i+1));
1e7108ce 142 TH1F *cfdtime = (TH1F*) gFile->Get(Form("CFD%d",i+1));
dd92f5b8 143 if(!cfd) {
144 AliWarning(Form("no Diff histograms collected by PHYS DA for channel %i", i));
32d3e409 145 okdiff++;
146 if(okdiff<4) {
147 meandiff = 0;
148 sigmadiff = 0;
149 }
150 else
151 ok = false;
dd92f5b8 152 }
153 if(!cfdtime) {
154 AliWarning(Form("no CFD histograms collected by PHYS DA for channel %i", i));
32d3e409 155 oktime++;
156 if(oktime<4) {
157 meancfdtime = 0;
158 }
159 else
160 ok = false;
dd92f5b8 161 }
eef60a6a 162 if(cfd) {
e691fec9 163 nent = Int_t(cfd->GetEntries());
11c0962a 164 if( nent<=50) {
32d3e409 165 okdiff++;
2a1ccb90 166 // printf(" pmt %i nent %i cfd->GetRMS() %f cfd->GetMean() %f \n",
167 // i, nent, cfd->GetRMS(), cfd->GetMean() );
32d3e409 168 if(okdiff<4) {
169 meandiff = 0;
170 sigmadiff = 0;
171 }
172 else
173 {
2a1ccb90 174 // printf(" OK fsle:: pmt %i nent %i cfd->GetRMS() %f cfd->GetMean() %f \n",
175 // i, nent, cfd->GetRMS(), cfd->GetMean() );
32d3e409 176 AliWarning(Form(" Not enouph data in PMT %i- PMT1: %i ", i, nent));
177 ok = false;
178 }
179 }
18086258 180 if(nent>50) {
612737bb 181 if(cfd->GetRMS()>1.5 )
e691fec9 182 GetMeanAndSigma(cfd, meandiff, sigmadiff);
183 if(cfd->GetRMS()<=1.5)
184 {
185 meandiff = cfd->GetMean();
186 sigmadiff=cfd->GetRMS();
187 }
dd92f5b8 188 Int_t maxBin = cfd->GetMaximumBin();
189 Double_t meanEstimate = cfd->GetBinCenter( maxBin);
190 if(TMath::Abs(meanEstimate - meandiff) > 20 ) meandiff = meanEstimate;
1e7108ce 191 }
1e7108ce 192 }
32d3e409 193
1e7108ce 194 if(cfdtime) {
e691fec9 195 nent = Int_t(cfdtime->GetEntries());
11c0962a 196 if( nent<=50 ) {
32d3e409 197 oktime++;
198 if(oktime<4) {
199 meancfdtime = 0;
200 }
201 else
202 {
203 AliWarning(Form(" Not enouph data in PMT %i: %i ", i, nent));
204 ok = false;
205 }
206 }
11c0962a 207 if(nent > 50 ) { //!!!!!!!!!!
b5e0299c 208 if(cfdtime->GetRMS()>1.5 )
e691fec9 209 GetMeanAndSigma(cfdtime,meancfdtime, sigmacfdtime);
210 if(cfdtime->GetRMS()<=1.5)
211 {
2a1ccb90 212 if(cfdtime->GetRMS()==0 ||cfdtime->GetMean()==0 ) {
213 ok = false;
214 }
e691fec9 215 meancfdtime = cfdtime->GetMean();
216 sigmacfdtime = cfdtime->GetRMS();
217 }
dd92f5b8 218 Int_t maxBin = cfdtime->GetMaximumBin();
219 Double_t meanEstimate = cfdtime->GetBinCenter( maxBin);
220 if(TMath::Abs(meanEstimate - meancfdtime) > 20 ) meancfdtime = meanEstimate;
1e7108ce 221 }
222 }
1e7108ce 223 SetTimeEq(i,meandiff);
224 SetTimeEqRms(i,sigmadiff);
225 SetCFDvalue(i,0,meancfdtime);
1e7108ce 226 if (cfd) delete cfd;
227 if (cfdtime) delete cfdtime;
1e7108ce 228 }
61bd03f6 229 TH1F *ver = (TH1F*) gFile->Get("hVertex") ;
1e7108ce 230 if(ver) {
231 meanver = ver->GetMean();
232 rmsver = ver->GetRMS();
233 }
234 SetMeanVertex(meanver);
235 SetRmsVertex(rmsver);
1e7108ce 236 gFile->Close();
237 delete gFile;
1e7108ce 238 }
239 return ok;
240}
241
242//________________________________________________________________
61bd03f6 243Int_t AliT0CalibTimeEq::ComputeOfflineParams(const char* filePhys, Float_t *timecdb, Float_t *cfdvalue, Int_t badpmt)
1e7108ce 244{
32d3e409 245 // fStatus implementation:
246 //ok means
247 // 1000 - can not open file
248 // for timediff
11c0962a 249 // 20 >3 histos are empty
250 // -11 low statistics oe empty histos in few channels; for these channels OCDB value will be written back WARNING
32d3e409 251 // for cfd
11c0962a 252 // 20 >2 histos are empty or with low statistic
253 // -11 if less 3 channels are empty or low statistics in few channels ;for these channels OCDB value will be written back WARNING
32d3e409 254 //
dd92f5b8 255 // compute offline equalized time
1e7108ce 256 Float_t meandiff, sigmadiff, meanver, meancfdtime, sigmacfdtime;
c4ba15c3 257 meandiff = sigmadiff = meanver = meancfdtime = sigmacfdtime =0;
1e7108ce 258 Int_t nent=0;
32d3e409 259 Int_t ok = 0;
32d3e409 260 Int_t okcfd=0;
bb54f817 261 TH1F *cfddiff = NULL;
262 TH1F *cfdtime = NULL;
61bd03f6 263 TObjArray * tzeroObj = NULL;
bb54f817 264
c4ba15c3 265 gFile = TFile::Open(filePhys);
1e7108ce 266 if(!gFile) {
267 AliError("No input PHYS data found ");
32d3e409 268 ok = 1000;
612737bb 269 return ok;
1e7108ce 270 }
271 else
272 {
e691fec9 273 meandiff = sigmadiff = meanver = meancfdtime = sigmacfdtime =0;
61bd03f6 274 // TDirectory *dr = (TDirectory*) gFile->Get("T0Calib");
275 tzeroObj = dynamic_cast<TObjArray*>(gFile->Get("T0Calib"));
612737bb 276 for (Int_t i=0; i<24; i++)
1e7108ce 277 {
61bd03f6 278 if (i != badpmt) {
279 if(tzeroObj) {
280 cfddiff = (TH1F*) tzeroObj->FindObject(Form("CFD1minCFD%d",i+1));
281 cfdtime = (TH1F*)tzeroObj->FindObject(Form("CFD%d",i+1));
bb54f817 282 }
283 else
284 {
dd92f5b8 285 cfddiff = (TH1F*)gFile->Get(Form("CFD1minCFD%d",i+1));
18086258 286 cfdtime = (TH1F*)gFile->Get(Form("CFD%d",i+1));
287
bb54f817 288 }
18086258 289 if(!cfddiff ) {
11c0962a 290 AliWarning(Form("no histograms collected by pass0 for diff channel %i\n", i));
32d3e409 291 meandiff = timecdb[i];
292 sigmadiff = 0;
32d3e409 293 }
dd92f5b8 294 if(cfddiff) {
18086258 295 nent = Int_t(cfddiff->GetEntries());
11c0962a 296 if ( nent == 0 ) {
297 AliWarning(Form("no entries in histogram for diff channel %i\n", i));
298 meandiff = timecdb[i];
299 sigmadiff = 0;
300 }
301 if(nent<=100 && nent>0) { //!!!!!
302 AliWarning(Form(" Not enouph data in PMT %i- PMT1: %i \n", i, nent));
32d3e409 303 meandiff = timecdb[i];
304 sigmadiff = 0;
f824d160 305 }
11c0962a 306 if(nent>=100 ) { //!!!!!
e691fec9 307 {
32d3e409 308 if(cfddiff->GetRMS()>1.5 )
309 GetMeanAndSigma(cfddiff, meandiff, sigmadiff);
310 if(cfddiff->GetRMS()<=1.5)
311 {
312 meandiff = cfddiff->GetMean();
313 sigmadiff = cfddiff->GetRMS();
314 }
315 Int_t maxBin = cfddiff->GetMaximumBin();
316 Double_t meanEstimate = cfddiff->GetBinCenter( maxBin);
317 if(TMath::Abs(meanEstimate - meandiff) > 20 ) meandiff = meanEstimate;
18086258 318 }
32d3e409 319 }
320 }
18086258 321 if(!cfdtime ) {
322 AliWarning(Form("no histograms collected by pass0 for time channel %i", i));
db082d76 323 meancfdtime = cfdvalue[i];
32d3e409 324 okcfd++;
11c0962a 325 if(okcfd<2) {
32d3e409 326 meancfdtime = cfdvalue[i];
11c0962a 327 // ok = -11;
32d3e409 328 }
329 else {
330 AliError(Form("no histograms collected by pass0 for time %i channels ", okcfd));
11c0962a 331 return 20;
32d3e409 332 }
18086258 333 }
334 if(cfdtime) {
335 nent = Int_t(cfdtime->GetEntries());
32d3e409 336 if (nent == 0 ||
337 (cfdtime->GetRMS() == 0 || cfdtime->GetMean() ==0 ) )
338 {
339 okcfd++;
11c0962a 340 if(okcfd<2) {
32d3e409 341 meancfdtime = cfdvalue[i];
11c0962a 342 // ok = -11;
32d3e409 343 printf("!!!!bad data:: pmt %i nent%i RMS %f mean %f cdbtime %f \n",
11c0962a 344 i, nent, cfdtime->GetRMS(), cfdtime->GetMean(), cfdvalue[i] );
32d3e409 345 }
346 else
347 {
348 printf("!!!!fatal data:: pmt %i nent%i RMS %f mean %f cdbtime %f \n",
11c0962a 349 i, nent, cfdtime->GetRMS(), cfdtime->GetMean(), cfdvalue[i]);
32d3e409 350 AliError(Form(" histograms collected by pass0 for time %i channels are empty", okcfd));
11c0962a 351 return 20;
32d3e409 352 }
353 }
354
355
11c0962a 356 if(nent<=100 && nent>0 )
32d3e409 357 {
11c0962a 358 okcfd++;
32d3e409 359 AliWarning(Form(" Not enouph data in PMT in CFD peak %i - %i ", i, nent));
32d3e409 360 meancfdtime = cfdvalue[i];
d09f51a3 361 // ok = -11;
32d3e409 362 printf("!!!!low statstics:: pmt %i nent%i RMS %f mean %f cdbtime %f \n",
11c0962a 363 i, nent, cfdtime->GetRMS(), cfdtime->GetMean(), cfdvalue[i]);
d09f51a3 364 if (okcfd>2) {
365 ok = -11;
366 return ok;
367 }
32d3e409 368 }
d09f51a3 369
11c0962a 370 if( nent>100 ) { //!!!!!
18086258 371 if(cfdtime->GetRMS()>1.5 )
372 GetMeanAndSigma(cfdtime,meancfdtime, sigmacfdtime);
373 if(cfdtime->GetRMS()<=1.5)
374 {
375 meancfdtime = cfdtime->GetMean();
376 sigmacfdtime=cfdtime->GetRMS();
32d3e409 377 }
657abde7 378 Int_t maxBin = cfdtime->GetMaximumBin();
379 Double_t meanEstimate = cfdtime->GetBinCenter( maxBin);
380 if(TMath::Abs(meanEstimate - meancfdtime) > 20 ) meancfdtime = meanEstimate;
1e7108ce 381 }
eef60a6a 382 }
bb54f817 383
b9dd81a7 384 SetTimeEq(i,meandiff);
1e7108ce 385 SetTimeEqRms(i,sigmadiff);
18086258 386 SetCFDvalue(i,0, meancfdtime );
bb54f817 387 if (cfddiff) cfddiff->Reset();
388 if (cfdtime) cfdtime->Reset();
b5e0299c 389 } //bad pmt
61bd03f6 390 }
6318dd46 391 gFile->Close();
392 delete gFile;
6318dd46 393 }
02888358 394 return ok;
e691fec9 395 }
a2ad8166 396
1e7108ce 397//________________________________________________________________________
398void AliT0CalibTimeEq::GetMeanAndSigma(TH1F* hist, Float_t &mean, Float_t &sigma) {
612737bb 399
6c18f84e 400 const double window = 5.; //fit window
1e7108ce 401 double meanEstimate, sigmaEstimate;
402 int maxBin;
403 maxBin = hist->GetMaximumBin(); //position of maximum
404 meanEstimate = hist->GetBinCenter( maxBin); // mean of gaussian sitting in maximum
6c18f84e 405 sigmaEstimate = hist->GetRMS();
406 // sigmaEstimate = 10;
1e7108ce 407 TF1* fit= new TF1("fit","gaus", meanEstimate - window*sigmaEstimate, meanEstimate + window*sigmaEstimate);
408 fit->SetParameters(hist->GetBinContent(maxBin), meanEstimate, sigmaEstimate);
11c0962a 409 hist->Fit("fit","QR","");
1e7108ce 410
411 mean = (Float_t) fit->GetParameter(1);
412 sigma = (Float_t) fit->GetParameter(2);
413
414 delete fit;
415}
416
8bfd9a3e 417