Line data Source code
1 :
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 : // //
21 : // class for T0 calibration TM-AC-AM_6-02-2006
22 : // equalize time shift for each time CFD channel
23 : // //
24 : ///////////////////////////////////////////////////////////////////////////////
25 :
26 : #include "AliT0CalibTimeEq.h"
27 : #include "AliLog.h"
28 : #include <TFile.h>
29 : #include <TMath.h>
30 : #include <TF1.h>
31 : #include <TSpectrum.h>
32 : #include <TProfile.h>
33 : #include <iostream>
34 :
35 20 : ClassImp(AliT0CalibTimeEq)
36 :
37 : //________________________________________________________________
38 3 : AliT0CalibTimeEq::AliT0CalibTimeEq():TNamed(),
39 3 : fMeanVertex(0),
40 3 : fRmsVertex(0)
41 15 : {
42 : //
43 150 : for(Int_t i=0; i<24; i++) {
44 72 : fTimeEq[i] = 0; // Time Equalized for OCDB
45 72 : fTimeEqRms[i] = -1; // RMS of Time Equalized for OCDB
46 864 : for (Int_t ih=0; ih<5; ih++) fCFDvalue[i][ih] = 0;
47 : }
48 6 : }
49 :
50 : //________________________________________________________________
51 0 : AliT0CalibTimeEq::AliT0CalibTimeEq(const char* name):TNamed(),
52 0 : fMeanVertex(0),
53 0 : fRmsVertex(0)
54 0 : {
55 : //constructor
56 :
57 0 : TString namst = "Calib_";
58 0 : namst += name;
59 0 : SetName(namst.Data());
60 0 : SetTitle(namst.Data());
61 0 : for(Int_t i=0; i<24; i++) {
62 0 : fTimeEq[i] = 0; // Time Equalized for OCDB
63 0 : fTimeEqRms[i] = -1; // RMS of Time Equalized for OCDB
64 0 : for (Int_t ih=0; ih<5; ih++) fCFDvalue[i][ih] = 0;
65 : }
66 0 : }
67 :
68 : //________________________________________________________________
69 0 : AliT0CalibTimeEq::AliT0CalibTimeEq(const AliT0CalibTimeEq& calibda):TNamed(calibda),
70 0 : fMeanVertex(0),
71 0 : fRmsVertex(0)
72 0 : {
73 : // copy constructor
74 0 : SetName(calibda.GetName());
75 0 : SetTitle(calibda.GetName());
76 0 : ((AliT0CalibTimeEq &) calibda).Copy(*this);
77 :
78 0 : }
79 :
80 : //________________________________________________________________
81 : AliT0CalibTimeEq &AliT0CalibTimeEq::operator =(const AliT0CalibTimeEq& calibda)
82 : {
83 : // assignment operator
84 0 : SetName(calibda.GetName());
85 0 : SetTitle(calibda.GetName());
86 :
87 0 : if (this != &calibda) (( AliT0CalibTimeEq &) calibda).Copy(*this);
88 0 : return *this;
89 : }
90 :
91 : //________________________________________________________________
92 : AliT0CalibTimeEq::~AliT0CalibTimeEq()
93 0 : {
94 : //
95 : // destrictor
96 0 : }
97 : //________________________________________________________________
98 : void AliT0CalibTimeEq::Reset()
99 : {
100 : //reset values
101 :
102 0 : memset(fCFDvalue,0,120*sizeof(Float_t));
103 0 : memset(fTimeEq,1,24*sizeof(Float_t));
104 0 : }
105 :
106 :
107 : //________________________________________________________________
108 : void AliT0CalibTimeEq::Print(Option_t*) const
109 : {
110 : // print time values
111 :
112 0 : printf("\n ---- PM Arrays ----\n\n");
113 0 : printf(" Time delay CFD \n");
114 0 : for (Int_t i=0; i<24; i++)
115 0 : printf(" CFD %f diff %f qt1 %f pedestal %f \n",fCFDvalue[i][0],fTimeEq[i],fCFDvalue[i][1], fCFDvalue[i][3]);
116 0 : printf(" TVDC %f OrA %f OrC %f \n",fMeanVertex,fCFDvalue[0][2], fCFDvalue[1][2]);
117 :
118 :
119 0 : }
120 :
121 :
122 : //________________________________________________________________
123 : Bool_t AliT0CalibTimeEq::ComputeOnlineParams(const char* filePhys)
124 : {
125 : // compute online equalized time
126 0 : Float_t meandiff, sigmadiff, meanver, meancfdtime, sigmacfdtime;
127 0 : Float_t ora,orc, meanqt, sigmaor,sigmaver, meanEstimate;
128 0 : Float_t meanpedold=0, sigmaped=0;
129 :
130 0 : meandiff = sigmadiff = meanver = meancfdtime = sigmacfdtime = ora = orc = meanqt =0;
131 0 : meanEstimate=sigmaver=sigmaor = 0;
132 : Int_t maxBin=0;
133 : Int_t nent=0;
134 : Int_t okdiff=0;
135 : Int_t oktime=0;
136 : Bool_t ok=false;
137 : //
138 0 : gFile = TFile::Open(filePhys);
139 0 : if(!gFile) {
140 0 : AliError("No input PHYS data found ");
141 0 : }
142 : else
143 : {
144 : ok=true;
145 0 : for (Int_t i=0; i<24; i++)
146 : {
147 0 : meandiff = sigmadiff = meanver = meancfdtime = sigmacfdtime =0;
148 0 : TH1F *hcfd = (TH1F*) gFile->Get(Form("CFD1minCFD%d",i+1));
149 0 : TH1F *hcfdtime = (TH1F*) gFile->Get(Form("CFD%d",i+1));
150 0 : TH1F *hqt1 = (TH1F*) gFile->Get(Form("QT1%d",i+1));
151 0 : TH1F *hPedOld = (TH1F*) gFile->Get(Form("hPed%i",i+1));
152 :
153 0 : if(!hcfd) {
154 0 : AliWarning(Form("no Diff histograms collected by PHYS DA for channel %i", i));
155 0 : okdiff++;
156 0 : if(okdiff<4) {
157 0 : meandiff = 0;
158 0 : sigmadiff = 0;
159 0 : }
160 : else
161 : ok = false;
162 : }
163 0 : if(!hcfdtime) {
164 0 : AliWarning(Form("no CFD histograms collected by PHYS DA for channel %i", i));
165 0 : oktime++;
166 0 : if(oktime<4) {
167 0 : meancfdtime = 0;
168 0 : }
169 : else
170 : ok = false;
171 : }
172 0 : if(hcfd) {
173 0 : nent = Int_t(hcfd->GetEntries());
174 0 : if( nent<=50) {
175 0 : okdiff++;
176 : // printf(" pmt %i nent %i cfd->GetRMS() %f cfd->GetMean() %f \n",
177 : // i, nent, hcfd->GetRMS(), hcfd->GetMean() );
178 0 : if(okdiff<4) {
179 0 : meandiff = 0;
180 0 : sigmadiff = 0;
181 0 : }
182 : else
183 : {
184 : // printf(" OK fsle:: pmt %i nent %i cfd->GetRMS() %f cfd->GetMean() %f \n",
185 : // i, nent, cfd->GetRMS(), cfd->GetMean() );
186 0 : AliWarning(Form(" Not enouph data in PMT %i- PMT1: %i ", i, nent));
187 : ok = false;
188 : }
189 : }
190 0 : if(nent>50) {
191 0 : if(hcfd->GetRMS()>1.5 )
192 0 : GetMeanAndSigma(hcfd, meandiff, sigmadiff);
193 0 : if(hcfd->GetRMS()<=1.5)
194 : {
195 0 : meandiff = hcfd->GetMean();
196 0 : sigmadiff=hcfd->GetRMS();
197 0 : }
198 0 : maxBin = hcfd->GetMaximumBin();
199 0 : meanEstimate = hcfd->GetBinCenter( maxBin);
200 0 : if(TMath::Abs(meanEstimate - meandiff) > 20 ) meandiff = meanEstimate;
201 : }
202 : }
203 :
204 0 : if(hcfdtime) {
205 0 : nent = Int_t(hcfdtime->GetEntries());
206 0 : if( nent<=50 ) {
207 0 : oktime++;
208 0 : if(oktime<4) {
209 0 : meancfdtime = 0;
210 0 : }
211 : else
212 : {
213 0 : AliWarning(Form(" Not enouph data in PMT %i: %i ", i, nent));
214 : ok = false;
215 : }
216 : }
217 0 : if(nent > 50 ) { //!!!!!!!!!!
218 0 : if(hcfdtime->GetRMS()>1.5 )
219 0 : GetMeanAndSigma(hcfdtime,meancfdtime, sigmacfdtime);
220 0 : if(hcfdtime->GetRMS()<=1.5)
221 : {
222 0 : if(hcfdtime->GetRMS()==0 ||hcfdtime->GetMean()==0 ) {
223 : ok = false;
224 0 : }
225 0 : meancfdtime = hcfdtime->GetMean();
226 0 : sigmacfdtime = hcfdtime->GetRMS();
227 0 : }
228 : }
229 0 : if(hqt1) GetMeanAndSigma(hqt1,meanqt, sigmaor);
230 : //Pedestals
231 0 : if(hPedOld ) GetMeanAndSigma(hPedOld ,meanpedold, sigmaped);
232 : } //cycle 24 PMT
233 0 : SetTimeEq(i,meandiff);
234 0 : SetTimeEqRms(i,sigmadiff);
235 0 : SetCFDvalue(i,0,meancfdtime);
236 0 : SetPedestalOld(i,meanpedold);
237 0 : SetCFDvalue(i,1,meanqt);
238 0 : if (hcfd) delete hcfd;
239 0 : if (hcfdtime) delete hcfdtime;
240 0 : if(hqt1) delete hqt1;
241 0 : if(hPedOld) delete hPedOld;
242 : }
243 0 : TH1F *hver = (TH1F*) gFile->Get("hVertex") ;
244 0 : TH1F *hora = (TH1F*) gFile->Get("hOrA") ;
245 0 : TH1F *horc = (TH1F*) gFile->Get("hOrC") ;
246 0 : if(hver) GetMeanAndSigma(hver,meanver, sigmaver);
247 0 : if(hora) GetMeanAndSigma(hora,ora, sigmaor);
248 0 : if(horc) GetMeanAndSigma(horc,orc, sigmaor);
249 0 : SetMeanVertex(meanver);
250 0 : SetRmsVertex(sigmaver);
251 0 : SetOrA(ora);
252 0 : SetOrC(orc);
253 :
254 0 : gFile->Close();
255 0 : delete gFile;
256 : }
257 0 : return ok;
258 0 : }
259 :
260 : //________________________________________________________________
261 : Int_t AliT0CalibTimeEq::ComputeOfflineParams(const char* filePhys, Float_t *timecdb, Float_t *cfdvalue, Int_t badpmt)
262 : {
263 : // fStatus implementation:
264 : //ok means
265 : // 1000 - can not open file
266 : // for timediff
267 : // 20 >3 histos are empty
268 : // -11 low statistics oe empty histos in few channels; for these channels OCDB value will be written back WARNING
269 : // for cfd
270 : // 20 >2 histos are empty or with low statistic
271 : // -11 if less 3 channels are empty or low statistics in few channels ;for these channels OCDB value will be written back WARNING
272 : //
273 : // compute offline equalized time
274 0 : Float_t meandiff, sigmadiff, meanver, meancfdtime, sigmacfdtime;
275 0 : meandiff = sigmadiff = meanver = meancfdtime = sigmacfdtime =0;
276 : Int_t nent=0;
277 : Int_t ok = 0;
278 : Int_t okcfd=0;
279 : TH1F *cfddiff = NULL;
280 : TH1F *cfdtime = NULL;
281 : TObjArray * tzeroObj = NULL;
282 0 : Float_t qt1[24], ped[24], orA, orC, tvdc;
283 :
284 0 : gFile = TFile::Open(filePhys);
285 0 : if(!gFile) {
286 0 : AliError("No input PHYS data found ");
287 : ok = 1000;
288 0 : return ok;
289 : }
290 : else {
291 0 : meandiff = sigmadiff = meanver = meancfdtime = sigmacfdtime =0;
292 : // TDirectory *dr = (TDirectory*) gFile->Get("T0Calib");
293 0 : tzeroObj = dynamic_cast<TObjArray*>(gFile->Get("T0Calib"));
294 0 : for (Int_t i=0; i<24; i++)
295 : {
296 0 : if (i != badpmt) {
297 0 : if(tzeroObj) {
298 0 : cfddiff = (TH1F*) tzeroObj->FindObject(Form("CFD1minCFD%d",i+1));
299 0 : cfdtime = (TH1F*)tzeroObj->FindObject(Form("CFD%d",i+1));
300 0 : }
301 : else
302 : {
303 0 : cfddiff = (TH1F*)gFile->Get(Form("CFD1minCFD%d",i+1));
304 0 : cfdtime = (TH1F*)gFile->Get(Form("CFD%d",i+1));
305 :
306 : }
307 0 : if(!cfddiff ) {
308 0 : AliWarning(Form("no histograms collected by pass0 for diff channel %i\n", i));
309 0 : meandiff = timecdb[i];
310 0 : sigmadiff = 0;
311 0 : }
312 0 : if(cfddiff) {
313 0 : nent = Int_t(cfddiff->GetEntries());
314 0 : if ( nent == 0 ) {
315 0 : AliWarning(Form("no entries in histogram for diff channel %i\n", i));
316 0 : meandiff = timecdb[i];
317 0 : sigmadiff = 0;
318 0 : }
319 0 : if(nent<=100 && nent>0) { //!!!!!
320 0 : AliWarning(Form(" Not enouph data in PMT %i- PMT1: %i \n", i, nent));
321 0 : meandiff = timecdb[i];
322 0 : sigmadiff = 0;
323 0 : }
324 0 : if(nent>=100 ) { //!!!!!
325 : {
326 0 : if(cfddiff->GetRMS()>1.5 )
327 0 : GetMeanAndSigma(cfddiff, meandiff, sigmadiff);
328 0 : if(cfddiff->GetRMS()<=1.5)
329 : {
330 0 : meandiff = cfddiff->GetMean();
331 0 : sigmadiff = cfddiff->GetRMS();
332 0 : }
333 0 : Int_t maxBin = cfddiff->GetMaximumBin();
334 0 : Double_t meanEstimate = cfddiff->GetBinCenter( maxBin);
335 0 : if(TMath::Abs(meanEstimate - meandiff) > 20 ) meandiff = meanEstimate;
336 : }
337 0 : }
338 : }
339 0 : if(!cfdtime ) {
340 0 : AliWarning(Form("no histograms collected by pass0 for time channel %i", i));
341 0 : meancfdtime = cfdvalue[i];
342 0 : okcfd++;
343 0 : if(okcfd<2) {
344 0 : meancfdtime = cfdvalue[i];
345 : // ok = -11;
346 : }
347 : else {
348 0 : AliError(Form("no histograms collected by pass0 for time %i channels ", okcfd));
349 0 : if (tzeroObj) delete tzeroObj;
350 0 : return 20;
351 : }
352 0 : }
353 0 : if(cfdtime) {
354 0 : nent = Int_t(cfdtime->GetEntries());
355 0 : if (nent == 0 ||
356 0 : (cfdtime->GetRMS() == 0 || cfdtime->GetMean() ==0 ) )
357 : {
358 0 : okcfd++;
359 0 : if(okcfd<2) {
360 0 : meancfdtime = cfdvalue[i];
361 : // ok = -11;
362 0 : printf("!!!!bad data:: pmt %i nent%i RMS %f mean %f cdbtime %f \n",
363 0 : i, nent, cfdtime->GetRMS(), cfdtime->GetMean(), cfdvalue[i] );
364 : }
365 : else
366 : {
367 0 : printf("!!!!fatal data:: pmt %i nent%i RMS %f mean %f cdbtime %f \n",
368 0 : i, nent, cfdtime->GetRMS(), cfdtime->GetMean(), cfdvalue[i]);
369 0 : AliError(Form(" histograms collected by pass0 for time %i channels are empty", okcfd));
370 0 : if (tzeroObj) delete tzeroObj;
371 0 : return 20;
372 : }
373 0 : }
374 :
375 :
376 0 : if(nent<=100 && nent>0 )
377 : {
378 0 : okcfd++;
379 0 : AliWarning(Form(" Not enouph data in PMT in CFD peak %i - %i ", i, nent));
380 0 : meancfdtime = cfdvalue[i];
381 : // ok = -11;
382 0 : printf("!!!!low statstics:: pmt %i nent%i RMS %f mean %f cdbtime %f \n",
383 0 : i, nent, cfdtime->GetRMS(), cfdtime->GetMean(), cfdvalue[i]);
384 0 : if (okcfd>2) {
385 : ok = -11;
386 0 : if (tzeroObj) delete tzeroObj;
387 0 : return ok;
388 : }
389 : }
390 :
391 0 : if( nent>100 ) { //!!!!!
392 0 : if(cfdtime->GetRMS()>1.5 )
393 0 : GetMeanAndSigma(cfdtime,meancfdtime, sigmacfdtime);
394 0 : if(cfdtime->GetRMS()<=1.5)
395 : {
396 0 : meancfdtime = cfdtime->GetMean();
397 0 : sigmacfdtime=cfdtime->GetRMS();
398 0 : }
399 0 : Int_t maxBin = cfdtime->GetMaximumBin();
400 0 : Double_t meanEstimate = cfdtime->GetBinCenter( maxBin);
401 0 : if(TMath::Abs(meanEstimate - meancfdtime) > 20 ) meancfdtime = meanEstimate;
402 0 : }
403 : }
404 :
405 0 : SetTimeEq(i,meandiff);
406 0 : SetTimeEqRms(i,sigmadiff);
407 0 : SetCFDvalue(i,0, meancfdtime );
408 0 : qt1[i]=cfdvalue[24+i];
409 0 : SetCFDvalue(i,1,qt1[i]);
410 0 : ped[i]=cfdvalue[52+i];
411 0 : SetCFDvalue(i,3,ped[i]);
412 0 : AliInfo(Form(" !!! AliT0CalibTimeEq pmt %i pedestal %f \n ", i, ped[i]) );
413 0 : if (cfddiff) cfddiff->Reset();
414 0 : if (cfdtime) cfdtime->Reset();
415 : } //bad pmt
416 : }
417 0 : SetMeanVertex(cfdvalue[48]);
418 0 : SetOrA(cfdvalue[49]);
419 0 : SetOrC(cfdvalue[50]);
420 :
421 0 : gFile->Close();
422 0 : delete gFile;
423 : }
424 0 : if (tzeroObj) delete tzeroObj;
425 0 : return ok;
426 0 : }
427 :
428 : //________________________________________________________________________
429 : void AliT0CalibTimeEq::GetMeanAndSigma(TH1F* hist, Float_t &mean, Float_t &sigma) {
430 :
431 : const double window = 3.; //fit window
432 : double meanEstimate, sigmaEstimate;
433 : int maxBin;
434 0 : maxBin = hist->GetMaximumBin(); //position of maximum
435 0 : meanEstimate = hist->GetBinCenter( maxBin); // mean of gaussian sitting in maximum
436 0 : sigmaEstimate = hist->GetRMS();
437 : // sigmaEstimate = 10;
438 0 : TF1* fit= new TF1("fit","gaus", meanEstimate - window*sigmaEstimate, meanEstimate + window*sigmaEstimate);
439 0 : fit->SetParameters(hist->GetBinContent(maxBin), meanEstimate, sigmaEstimate);
440 0 : hist->Fit("fit","RQ","");
441 :
442 0 : mean = (Float_t) fit->GetParameter(1);
443 0 : sigma = (Float_t) fit->GetParameter(2);
444 :
445 0 : if(TMath::Abs(meanEstimate - mean) > 20 ) mean = meanEstimate;
446 :
447 :
448 0 : delete fit;
449 0 : }
450 :
|
