source: trunk/MagicSoft/Mars/mhcalib/MHCalibrationPulseTimeCam.cc@ 8408

Last change on this file since 8408 was 8361, checked in by tbretz, 18 years ago
*** empty log message ***
File size: 21.7 KB
Line 
1/* ======================================================================== *\
2! $Name: not supported by cvs2svn $:$Id: MHCalibrationPulseTimeCam.cc,v 1.34 2007-03-04 11:55:55 tbretz Exp $
3! --------------------------------------------------------------------------
4!
5! *
6! * This file is part of MARS, the MAGIC Analysis and Reconstruction
7! * Software. It is distributed to you in the hope that it can be a useful
8! * and timesaving tool in analysing Data of imaging Cerenkov telescopes.
9! * It is distributed WITHOUT ANY WARRANTY.
10! *
11! * Permission to use, copy, modify and distribute this software and its
12! * documentation for any purpose is hereby granted without fee,
13! * provided that the above copyright notice appear in all copies and
14! * that both that copyright notice and this permission notice appear
15! * in supporting documentation. It is provided "as is" without express
16! * or implied warranty.
17! *
18!
19!
20! Author(s): Markus Gaug 02/2004 <mailto:markus@ifae.es>
21! Author(s): Thomas Bretz <mailto:tbretz@astro.uni-wuerzburg.de>
22!
23! Copyright: MAGIC Software Development, 2000-2007
24!
25!
26\* ======================================================================== */
27/////////////////////////////////////////////////////////////////////////////
28//
29// MHCalibrationPulseTimeCam
30//
31// Fills the extracted signals of MExtractedSignalCam into the MHCalibrationPix-classes
32// MHCalibrationPulseTimeHiGainPix and MHCalibrationPulseTimeLoGainPix for every:
33//
34// - Pixel, stored in the TOrdCollection's MHCalibrationCam::fHiGainArray and
35// MHCalibrationCam::fLoGainArray
36//
37// - Average pixel per AREA index (e.g. inner and outer for the MAGIC camera),
38// stored in the TOrdCollection's MHCalibrationCam::fAverageHiGainAreas and
39// MHCalibrationCam::fAverageLoGainAreas
40//
41// - Average pixel per camera SECTOR (e.g. sectors 1-6 for the MAGIC camera),
42// stored in the TOrdCollection's MHCalibrationCam::fAverageHiGainSectors and
43// MHCalibrationCam::fAverageLoGainSectors
44//
45// Every signal is taken from MExtractedSignalCam and filled into a histogram and
46// an array, in order to perform a Fourier analysis (see MHGausEvents).
47// The signals are moreover averaged on an event-by-event basis and written into
48// the corresponding average pixels.
49//
50// Additionally, the (FADC slice) position of the maximum is stored in an Absolute
51// Arrival Time histogram. This histogram serves for a rough cross-check if the
52// signal does not lie at or outside the edges of the extraction window.
53//
54// The PulseTime histograms are fitted to a Gaussian, mean and sigma with its errors
55// and the fit probability are extracted. If none of these values are NaN's and
56// if the probability is bigger than MHGausEvents::fProbLimit (default: 0.5%),
57// the fit is declared valid.
58// Otherwise, the fit is repeated within ranges of the previous mean
59// +- MHCalibrationPix::fPickupLimit (default: 5) sigma (see MHCalibrationPix::RepeatFit())
60// In case this does not make the fit valid, the histogram means and RMS's are
61// taken directly (see MHCalibrationPix::BypassFit()) and the following flags are set:
62// - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::kHiGainNotFitted ) or
63// - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::kLoGainNotFitted ) and
64// - MBadPixelsPix::SetUnsuitable( MBadPixelsPix::kUnreliableRun )
65//
66// Outliers of more than MHCalibrationPix::fPickupLimit (default: 5) sigmas
67// from the mean are counted as Pickup events (stored in MHCalibrationPix::fPickup)
68//
69// Unless more than fNumHiGainSaturationLimit (default: 1%) of the overall FADC
70// slices show saturation, the following flag is set:
71// - MCalibrationPulseTimePix::SetHiGainSaturation();
72// In that case, the calibration constants are derived from the low-gain results.
73//
74// If more than fNumLoGainSaturationLimit (default: 1%) of the overall
75// low-gain FADC slices saturate, the following flags are set:
76// - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::kLoGainSaturation ) and
77// - MBadPixelsPix::SetUnsuitable( MBadPixelsPix::kUnsuitableRun )
78//
79// The class also fills arrays with the signal vs. event number, creates a fourier
80// spectrum and investigates if the projected fourier components follow an exponential
81// distribution. In case that the probability of the exponential fit is less than
82// MHGausEvents::fProbLimit (default: 0.5%), the following flags are set:
83// - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::kHiGainOscillating ) or
84// - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::kLoGainOscillating ) and
85// - MBadPixelsPix::SetUnsuitable( MBadPixelsPix::kUnreliableRun )
86//
87// This same procedure is performed for the average pixels.
88//
89// The following results are written into MCalibrationPulseTimeCam:
90//
91// - MCalibrationPix::SetHiGainSaturation()
92// - MCalibrationPix::SetHiGainMean()
93// - MCalibrationPix::SetHiGainMeanErr()
94// - MCalibrationPix::SetHiGainSigma()
95// - MCalibrationPix::SetHiGainSigmaErr()
96// - MCalibrationPix::SetHiGainProb()
97// - MCalibrationPix::SetHiGainNumPickup()
98//
99// For all averaged areas, the fitted sigma is multiplied with the square root of
100// the number involved pixels in order to be able to compare it to the average of
101// sigmas in the camera.
102//
103/////////////////////////////////////////////////////////////////////////////
104#include "MHCalibrationPulseTimeCam.h"
105
106#include <TEnv.h>
107#include <TLine.h>
108#include <TGraph.h>
109#include <TLegend.h>
110#include <TCanvas.h>
111#include <TOrdCollection.h>
112
113#include "MLog.h"
114#include "MLogManip.h"
115
116#include "MParList.h"
117
118#include "MArrayD.h"
119#include "MRawRunHeader.h"
120#include "MExtractedSignalCam.h"
121#include "MPedestalSubtractedEvt.h"
122
123#include "MCalibrationPix.h"
124#include "MHCalibrationPix.h"
125#include "MCalibrationIntensityCam.h"
126#include "MBadPixelsIntensityCam.h"
127
128#include "MGeomCam.h"
129#include "MGeomPix.h"
130
131#include "MBadPixelsCam.h"
132
133ClassImp(MHCalibrationPulseTimeCam);
134
135using namespace std;
136
137const Byte_t MHCalibrationPulseTimeCam::fgSaturationLimit = 245;
138const Byte_t MHCalibrationPulseTimeCam::fgLowerSignalLimit = 85;
139const Int_t MHCalibrationPulseTimeCam::fgNumPixelsRequired = 2;
140const Int_t MHCalibrationPulseTimeCam::fgHiGainNbins = 20;
141const Axis_t MHCalibrationPulseTimeCam::fgHiGainFirst = -0.5;
142const Axis_t MHCalibrationPulseTimeCam::fgHiGainLast = 19.5;
143const Float_t MHCalibrationPulseTimeCam::fgProbLimit = 0.001;
144const TString MHCalibrationPulseTimeCam::gsHistName = "PulseTime";
145const TString MHCalibrationPulseTimeCam::gsHistTitle = "Extracted Times";
146const TString MHCalibrationPulseTimeCam::gsHistXTitle = "Time [FADC slices]";
147const TString MHCalibrationPulseTimeCam::gsHistYTitle = "Nr. events";
148const TString MHCalibrationPulseTimeCam::fgReferenceFile = "mjobs/signalref.rc";
149
150// --------------------------------------------------------------------------
151//
152// Default Constructor.
153//
154// Sets:
155// - all pointers to NULL
156//
157// - fNbins to fgHiGainNbins
158// - fFirst to fgHiGainFirst
159// - fLast to fgHiGainLast
160//
161// - fHistName to gsHistName
162// - fHistTitle to gsHistTitle
163// - fHistXTitle to gsHistXTitle
164// - fHistYTitle to gsHistYTitle
165//
166// - fSaturationLimit to fgSaturationLimit
167// - fLowerSignalLimit to fgLowerSignalLimit
168// - fNumPixelsRequired to fgNumPixelsRequired
169//
170MHCalibrationPulseTimeCam::MHCalibrationPulseTimeCam(const char *name, const char *title)
171 : fBadPixels(NULL)
172{
173
174 fName = name ? name : "MHCalibrationPulseTimeCam";
175 fTitle = title ? title : "Class to fill the extracted pulse times for cosmics ";
176
177 SetBinning(fgHiGainNbins, fgHiGainFirst, fgHiGainLast);
178
179 SetProbLimit(fgProbLimit);
180
181 SetHistName (gsHistName .Data());
182 SetHistTitle (gsHistTitle .Data());
183 SetHistXTitle(gsHistXTitle.Data());
184 SetHistYTitle(gsHistYTitle.Data());
185
186 SetReferenceFile();
187 SetLoGain(kFALSE);
188 SetOscillations(kFALSE);
189
190 SetSaturationLimit();
191 SetLowerSignalLimit();
192 SetNumPixelsRequired();
193
194 fInnerRefTime = 5.;
195 fOuterRefTime = 5.;
196}
197// --------------------------------------------------------------------------
198//
199// Creates new MHCalibrationPulseTimeCam only with the averaged areas:
200// the rest has to be retrieved directly, e.g. via:
201// MHCalibrationPulseTimeCam *cam = MParList::FindObject("MHCalibrationPulseTimeCam");
202// - cam->GetAverageSector(5).DrawClone();
203// - (*cam)[100].DrawClone()
204//
205TObject *MHCalibrationPulseTimeCam::Clone(const char *) const
206{
207
208 MHCalibrationPulseTimeCam *cam = new MHCalibrationPulseTimeCam();
209
210 //
211 // Copy the data members
212 //
213 cam->fColor = fColor;
214 cam->fRunNumbers = fRunNumbers;
215 cam->fPulserFrequency = fPulserFrequency;
216 cam->fFlags = fFlags;
217 cam->fNbins = fNbins;
218 cam->fFirst = fFirst;
219 cam->fLast = fLast;
220 cam->fReferenceFile = fReferenceFile;
221 cam->fInnerRefTime = fInnerRefTime;
222 cam->fOuterRefTime = fOuterRefTime;
223
224 if (!IsAverageing())
225 return cam;
226
227 const Int_t navhi = fAverageHiGainAreas->GetSize();
228
229 for (int i=0; i<navhi; i++)
230 cam->fAverageHiGainAreas->AddAt(GetAverageHiGainArea(i).Clone(),i);
231
232 return cam;
233}
234
235// --------------------------------------------------------------------------
236//
237// Gets the pointers to:
238// - MRawEvtData
239//
240Bool_t MHCalibrationPulseTimeCam::SetupHists(const MParList *pList)
241{
242
243 fBadPixels = (MBadPixelsCam*)pList->FindObject("MBadPixelsCam");
244 if (!fBadPixels)
245 {
246 *fLog << warn << GetDescriptor() << "MBadPixelsCam not found... " << endl;
247 }
248
249 return kTRUE;
250}
251
252// --------------------------------------------------------------------------
253//
254// Gets or creates the pointers to:
255// - MExtractedSignalCam
256// - MCalibrationPulseTimeCam or MCalibrationIntensityPulseTimeCam
257// - MBadPixelsCam
258//
259// Initializes the number of used FADC slices from MExtractedSignalCam
260// into MCalibrationPulseTimeCam and test for changes in that variable
261//
262// Calls:
263// - InitHiGainArrays()
264//
265// Sets:
266// - fSumhiarea to nareas
267// - fSumloarea to nareas
268// - fSumhisector to nsectors
269// - fSumlosector to nsectors
270//
271Bool_t MHCalibrationPulseTimeCam::ReInitHists(MParList *pList)
272{
273
274 MExtractedSignalCam *signal =
275 (MExtractedSignalCam*)pList->FindObject(AddSerialNumber("MExtractedSignalCam"));
276 if (!signal)
277 {
278 *fLog << err << "MExtractedSignalCam not found... abort." << endl;
279 return kFALSE;
280 }
281
282 if (!InitCams(pList,"PulseTime"))
283 return kFALSE;
284
285 const Int_t npixels = fGeom->GetNumPixels();
286 const Int_t nsectors = fGeom->GetNumSectors();
287 const Int_t nareas = fGeom->GetNumAreas();
288
289 InitHiGainArrays(npixels,nareas,nsectors);
290
291 return kTRUE;
292}
293
294void MHCalibrationPulseTimeCam::InitHiGainArrays(const Int_t npixels, const Int_t nareas, const Int_t nsectors)
295{
296
297 const Int_t samples = fRunHeader->GetNumSamples();
298 SetBinning(samples, -0.5, samples-0.5);
299
300 if (fHiGainArray->GetSize()==0)
301 {
302 for (Int_t i=0; i<npixels; i++)
303 {
304 fHiGainArray->AddAt(new MHCalibrationPix(Form("%sHiGainPix%04d",fHistName.Data(),i),
305 Form("%s High Gain Pixel %4d",fHistTitle.Data(),i)),i);
306
307 MHCalibrationPix &pix = (*this)[i];
308 pix.SetBinning(fNbins, fFirst, fLast);
309
310 MBadPixelsPix &bad = fIntensBad ? (*fIntensBad)[i] : (*fBadPixels)[i];
311 InitHists(pix,bad,i);
312
313 if (fCam)
314 (*fCam)[i].SetPixId(i);
315 }
316 }
317
318 if (!IsAverageing())
319 return;
320
321 if (fAverageHiGainAreas->GetSize()==0)
322 {
323 for (Int_t j=0; j<nareas; j++)
324 {
325 fAverageHiGainAreas->AddAt(new MHCalibrationPix(Form("%sHiGainArea%d",fHistName.Data(),j),
326 Form("%s High Gain Area Idx %d",fHistTitle.Data(),j)),j);
327
328 MHCalibrationPix &pix = GetAverageHiGainArea(j);
329 pix.SetBinning(fNbins, fFirst, fLast);
330
331 InitHists(pix,fIntensCam ? fIntensCam->GetAverageBadArea(j) : fCam->GetAverageBadArea(j),j);
332 }
333 }
334
335 if (fAverageHiGainSectors->GetSize()==0)
336 {
337 for (Int_t j=0; j<nsectors; j++)
338 {
339 fAverageHiGainSectors->AddAt(new MHCalibrationPix(Form("%sHiGainSector%02d",fHistName.Data(),j),
340 Form("%s High Gain Sector %02d",fHistTitle.Data(),j)),j);
341 MHCalibrationPix &pix = GetAverageHiGainSector(j);
342 pix.SetBinning(fNbins, fFirst, fLast);
343
344 InitHists(pix,fIntensCam ? fIntensCam->GetAverageBadSector(j) : fCam->GetAverageBadSector(j),j);
345 }
346 }
347}
348
349// --------------------------------------------------------------------------
350//
351// Retrieves from MExtractedSignalCam:
352// - first used LoGain FADC slice
353//
354// Retrieves from MGeomCam:
355// - number of pixels
356// - number of pixel areas
357// - number of sectors
358//
359// For all TOrdCollection's (including the averaged ones), the following steps are performed:
360//
361// 1) Fill PulseTimes histograms (MHGausEvents::FillHistAndArray()) with:
362// - MExtractedSignalPix::GetExtractedSignalHiGain();
363// - MExtractedSignalPix::GetExtractedSignalLoGain();
364//
365Bool_t MHCalibrationPulseTimeCam::FillHists(const MParContainer *par, const Stat_t w)
366{
367 MPedestalSubtractedEvt *evt = (MPedestalSubtractedEvt*)par;
368 if (!evt)
369 {
370 *fLog << err << "No argument in MHCalibrationPulseTimeCam::Fill... abort." << endl;
371 return kFALSE;
372 }
373
374 const UInt_t nareas = fGeom->GetNumAreas();
375 const UInt_t nsectors = fGeom->GetNumSectors();
376
377 MArrayD sumarea(nareas);
378 MArrayD sumsector(nsectors);
379
380 fAverageAreaNum.Reset();
381 fAverageSectorNum.Reset();
382
383 const Int_t npix = evt->GetNumPixels();
384 for (int idx=0; idx<npix; idx++)
385 {
386 if (fBadPixels)
387 if ((*fBadPixels)[idx].IsUnsuitable())
388 continue;
389
390 // Check for saturation
391 if (evt->GetSaturation(idx, fSaturationLimit)>0)
392 continue;
393
394 // Get position of maximum
395 const Int_t pos = evt->GetMax(idx/*, first, last*/);
396 const Float_t max = evt->GetSamples(idx)[pos];
397
398 // check if maximum is high enough
399 if (max<fLowerSignalLimit)
400 continue;
401
402 const Int_t maxpos = pos-1;
403
404 (*this)[idx].FillHist(maxpos);
405
406 const Int_t aidx = (*fGeom)[idx].GetAidx();
407 const Int_t sector = (*fGeom)[idx].GetSector();
408
409 sumarea[aidx] += maxpos;
410 sumsector[sector] += maxpos;
411
412 fAverageAreaNum[aidx]++;
413 fAverageSectorNum[sector]++;
414 }
415
416 for (UInt_t j=0; j<nareas; j++)
417 {
418 const Int_t npix = fAverageAreaNum[j];
419
420 if (npix > fNumPixelsRequired)
421 {
422 if (IsOscillations())
423 GetAverageHiGainArea(j).FillHistAndArray(sumarea[j]/npix);
424 else
425 GetAverageHiGainArea(j).FillHist(sumarea[j]/npix);
426
427 }
428 }
429
430 for (UInt_t j=0; j<nsectors; j++)
431 {
432 const Int_t npix = fAverageSectorNum[j];
433
434 if (npix > 0)
435 {
436 if (IsOscillations())
437 GetAverageHiGainSector(j).FillHistAndArray(sumsector[j]/npix);
438 else
439 GetAverageHiGainSector(j).FillHist(sumsector[j]/npix);
440 }
441 }
442
443 return kTRUE;
444}
445
446// --------------------------------------------------------------------------
447//
448// For all TOrdCollection's (including the averaged ones), the following steps are performed:
449//
450// 1) Returns if the pixel is excluded.
451// 2) Tests saturation. In case yes, set the flag: MCalibrationPix::SetHiGainSaturation()
452// or the flag: MBadPixelsPix::SetUncalibrated( MBadPixelsPix::kLoGainSaturated )
453// 3) Store the absolute arrival times in the MCalibrationPulseTimePix's. If flag
454// MCalibrationPix::IsHiGainSaturation() is set, the Low-Gain arrival times are stored,
455// otherwise the Hi-Gain ones.
456// 4) Calls to MHCalibrationCam::FitHiGainArrays() and MCalibrationCam::FitLoGainArrays()
457// with the flags:
458// - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::kHiGainNotFitted )
459// - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::kLoGainNotFitted )
460// - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::kHiGainOscillating )
461// - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::kLoGainOscillating )
462//
463Bool_t MHCalibrationPulseTimeCam::FinalizeHists()
464{
465
466 *fLog << endl;
467
468 MCalibrationCam &calcam = *(fIntensCam ? fIntensCam->GetCam() : fCam);
469 //
470 // Perform the fitting for the High Gain (done in MHCalibrationCam)
471 //
472 for (Int_t i=0; i<fHiGainArray->GetSize(); i++)
473 {
474
475 MHCalibrationPix &hist = (*this)[i];
476
477 if (hist.IsExcluded())
478 continue;
479
480 MCalibrationPix &pix = calcam[i];
481 CalcHists(hist,pix);
482 }
483
484 if (!IsAverageing())
485 return kTRUE;
486
487 for (Int_t j=0; j<fAverageHiGainAreas->GetSize(); j++)
488 {
489 MHCalibrationPix &hist = GetAverageHiGainArea(j);
490 MCalibrationPix &pix = calcam.GetAverageArea(j);
491 CalcHists(hist,pix);
492 }
493
494 for (Int_t j=0; j<fAverageHiGainSectors->GetSize(); j++)
495 {
496 MHCalibrationPix &hist = GetAverageHiGainSector(j);
497 MCalibrationPix &pix = calcam.GetAverageSector(j);
498 CalcHists(hist,pix);
499 }
500
501 return kTRUE;
502}
503
504void MHCalibrationPulseTimeCam::CalcHists(MHCalibrationPix &hist, MCalibrationPix &pix) const
505{
506
507 if (hist.IsEmpty())
508 {
509 *fLog << warn << hist.GetName() << ": Histogram empty." << endl;
510 return;
511 }
512 if (hist.IsOnlyOverflow())
513 {
514 *fLog << warn << hist.GetName() << ": Histogram contains only overflows." << endl;
515 return;
516 }
517 if (hist.IsOnlyUnderflow())
518 {
519 *fLog << warn << hist.GetName() << ": Histogram contains only underflows." << endl;
520 return;
521 }
522
523 hist.BypassFit();
524
525 pix.SetHiGainMean ( hist.GetMean() );
526 pix.SetHiGainMeanVar ( hist.GetMeanErr() * hist.GetMeanErr() );
527 pix.SetHiGainRms ( hist.GetHistRms() );
528 pix.SetHiGainSigma ( hist.GetSigma() );
529 pix.SetHiGainSigmaVar ( hist.GetSigmaErr()* hist.GetSigmaErr() );
530
531 if (IsDebug())
532 {
533 *fLog << dbginf << GetDescriptor() << ": ID " << GetName()
534 << " "<<pix.GetPixId()
535 << " Mean: " << hist.GetMean ()
536 << " MeanErr: " << hist.GetMeanErr ()
537 << " MeanSigma: " << hist.GetSigma ()
538 << " MeanSigmaErr: " << hist.GetSigmaErr()
539 << " Prob: " << hist.GetProb ()
540 << endl;
541 }
542}
543
544
545// --------------------------------------------------------------------------
546//
547// Calls MHCalibrationPix::DrawClone() for pixel idx
548//
549void MHCalibrationPulseTimeCam::DrawPixelContent(Int_t idx) const
550{
551 (*this)[idx].DrawClone();
552}
553
554
555// -----------------------------------------------------------------------------
556//
557// Default draw:
558//
559// Displays the averaged areas, both High Gain and Low Gain
560//
561// Calls the Draw of the fAverageHiGainAreas and fAverageLoGainAreas objects with options
562//
563void MHCalibrationPulseTimeCam::Draw(const Option_t *opt)
564{
565
566 const Int_t nareas = fAverageHiGainAreas->GetSize();
567 if (nareas == 0)
568 return;
569
570 TString option(opt);
571 option.ToLower();
572
573 TVirtualPad *pad = gPad ? gPad : MH::MakeDefCanvas(this);
574 pad->SetBorderMode(0);
575 pad->Divide(1,nareas, 1e-10, 1e-10);
576
577 //
578 // Loop over inner and outer pixels
579 //
580 for (Int_t i=0; i<nareas;i++)
581 {
582 pad->cd(i+1);
583
584 MHCalibrationPix &hipix = GetAverageHiGainArea(i);
585 DrawDataCheckPixel(hipix,i ? fOuterRefTime : fInnerRefTime);
586 }
587}
588
589// -----------------------------------------------------------------------------
590//
591// Draw the average pixel for the datacheck:
592//
593// Displays the averaged areas, both High Gain and Low Gain
594//
595// Calls the Draw of the fAverageHiGainAreas and fAverageLoGainAreas objects with options
596//
597void MHCalibrationPulseTimeCam::DrawDataCheckPixel(MHCalibrationPix &pix, const Float_t refline)
598{
599 gPad->SetBorderMode(0);
600 gPad->SetTicks();
601
602 TH1F *hist = pix.GetHGausHist();
603
604 //
605 // set the labels bigger
606 //
607 TAxis *xaxe = hist->GetXaxis();
608 TAxis *yaxe = hist->GetYaxis();
609 xaxe->CenterTitle();
610 yaxe->CenterTitle();
611 xaxe->SetTitleSize(0.06);
612 yaxe->SetTitleSize(0.076);
613 xaxe->SetTitleOffset(0.6);
614 yaxe->SetTitleOffset(0.65);
615 xaxe->SetLabelSize(0.06);
616 yaxe->SetLabelSize(0.06);
617 xaxe->SetTitle(hist->GetXaxis()->GetTitle());
618 yaxe->SetTitle(hist->GetYaxis()->GetTitle());
619 xaxe->SetRange(hist->GetMaximumBin()-30, hist->GetMaximumBin()+30);
620
621 hist->Draw();
622
623 DisplayRefLines(hist, refline);
624}
625
626void MHCalibrationPulseTimeCam::DisplayRefLines(const TH1F *hist, const Float_t refline) const
627{
628 TLine *line = new TLine(refline, 0, refline, hist->GetMaximum());
629 line->SetLineColor(106);
630 line->SetLineStyle(2);
631 line->SetLineWidth(3);
632 line->SetBit(kCanDelete);
633 line->Draw();
634
635 TLegend *leg = new TLegend(0.8,0.35,0.99,0.65);
636 leg->AddEntry(line, "Reference", "l");
637 leg->SetBit(kCanDelete);
638 leg->Draw();
639}
640
641Int_t MHCalibrationPulseTimeCam::ReadEnv(const TEnv &env, TString prefix, Bool_t print)
642{
643
644 Int_t rc = MHCalibrationCam::ReadEnv(env,prefix,print);
645 if (rc==kERROR)
646 return kERROR;
647
648 if (IsEnvDefined(env, prefix, "SaturationLimit", print))
649 {
650 SetSaturationLimit(GetEnvValue(env, prefix, "SaturationLimit", fSaturationLimit));
651 rc = kTRUE;
652 }
653
654 if (IsEnvDefined(env, prefix, "LowerSignalLimit", print))
655 {
656 SetLowerSignalLimit(GetEnvValue(env,prefix,"LowerSignalLimit",fLowerSignalLimit));
657 rc = kTRUE;
658 }
659
660 if (IsEnvDefined(env, prefix, "NumPixelsRequired", print))
661 {
662 SetNumPixelsRequired(GetEnvValue(env,prefix,"NumPixelsRequired",fNumPixelsRequired));
663 rc = kTRUE;
664 }
665
666 if (IsEnvDefined(env, prefix, "ReferenceFile", print))
667 {
668 SetReferenceFile(GetEnvValue(env,prefix,"ReferenceFile",fReferenceFile.Data()));
669 rc = kTRUE;
670 }
671
672 TEnv refenv(fReferenceFile);
673
674 fInnerRefTime = refenv.GetValue("InnerRefTime",fInnerRefTime);
675 fOuterRefTime = refenv.GetValue("OuterRefTime",fOuterRefTime);
676
677 return rc;
678}
Note: See TracBrowser for help on using the repository browser.