source: trunk/MagicSoft/Mars/mhcalib/MHCalibrationRelTimeCam.cc@ 8433

Last change on this file since 8433 was 8428, checked in by tbretz, 18 years ago
*** empty log message ***
File size: 25.3 KB
Line 
1/* ======================================================================== *\
2!
3! *
4! * This file is part of MARS, the MAGIC Analysis and Reconstruction
5! * Software. It is distributed to you in the hope that it can be a useful
6! * and timesaving tool in analysing Data of imaging Cerenkov telescopes.
7! * It is distributed WITHOUT ANY WARRANTY.
8! *
9! * Permission to use, copy, modify and distribute this software and its
10! * documentation for any purpose is hereby granted without fee,
11! * provided that the above copyright notice appear in all copies and
12! * that both that copyright notice and this permission notice appear
13! * in supporting documentation. It is provided "as is" without express
14! * or implied warranty.
15! *
16!
17!
18! Author(s): Markus Gaug 02/2004 <mailto:markus@ifae.es>
19! Author(s): Thomas Bretz <mailto:tbretz@astro.uni-wuerzburg.de>
20!
21! Copyright: MAGIC Software Development, 2000-2007
22!
23!
24\* ======================================================================== */
25
26/////////////////////////////////////////////////////////////////////////////
27//
28// MHCalibrationRelTimeCam
29//
30// Fills the extracted relative arrival times of MArrivalTimeCam into
31// the MHCalibrationPix-classes MHCalibrationPix for every:
32//
33// - Pixel, stored in the TObjArray's MHCalibrationCam::fHiGainArray
34// or MHCalibrationCam::fHiGainArray, respectively, depending if
35// MArrivalTimePix::IsLoGainUsed() is set.
36//
37// - Average pixel per AREA index (e.g. inner and outer for the MAGIC camera),
38// stored in the TObjArray's MHCalibrationCam::fAverageHiGainAreas and
39// MHCalibrationCam::fAverageHiGainAreas
40//
41// - Average pixel per camera SECTOR (e.g. sectors 1-6 for the MAGIC camera),
42// stored in the TObjArray's MHCalibrationCam::fAverageHiGainSectors
43// and MHCalibrationCam::fAverageHiGainSectors
44//
45// Every relative time is calculated as the difference between the individual
46// pixel arrival time and the one of pixel 1 (hardware number: 2).
47// The relative times are filled into a histogram and an array, in order to perform
48// a Fourier analysis (see MHGausEvents). The signals are moreover averaged on an
49// event-by-event basis and written into the corresponding average pixels.
50//
51// The histograms are fitted to a Gaussian, mean and sigma with its errors
52// and the fit probability are extracted. If none of these values are NaN's and
53// if the probability is bigger than MHGausEvents::fProbLimit (default: 0.5%),
54// the fit is declared valid.
55// Otherwise, the fit is repeated within ranges of the previous mean
56// - MHCalibrationPix::fPickupLimit (default: 5) sigma (see MHCalibrationPix::RepeatFit())
57// In case this does not make the fit valid, the histogram means and RMS's are
58// taken directly (see MHCalibrationPix::BypassFit()) and the following flags are set:
59// - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::kRelTimeNotFitted ) and
60// - MBadPixelsPix::SetUnsuitable( MBadPixelsPix::kUnreliableRun )
61//
62// Outliers of more than MHCalibrationPix::fPickupLimit (default: 5) sigmas
63// from the mean are counted as Pickup events (stored in MHCalibrationPix::fPickup)
64//
65// The class also fills arrays with the signal vs. event number, creates a fourier
66// spectrum (see MHGausEvents::CreateFourierSpectrum()) and investigates if the
67// projected fourier components follow an exponential distribution.
68// In case that the probability of the exponential fit is less than
69// MHGausEvents::fProbLimit (default: 0.5%), the following flags are set:
70// - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::kRelTimeOscillating ) and
71// - MBadPixelsPix::SetUnsuitable( MBadPixelsPix::kUnreliableRun )
72//
73// This same procedure is performed for the average pixels.
74//
75// The following results are written into MCalibrationRelTimeCam:
76//
77// - MCalibrationPix::SetMean()
78// - MCalibrationPix::SetMeanErr()
79// - MCalibrationPix::SetSigma()
80// - MCalibrationPix::SetSigmaErr()
81// - MCalibrationPix::SetProb()
82// - MCalibrationPix::SetNumPickup()
83//
84// For all averaged areas, the fitted sigma is multiplied with the square root of
85// the number involved pixels in order to be able to compare it to the average of
86// sigmas in the camera.
87//
88/////////////////////////////////////////////////////////////////////////////
89#include "MHCalibrationRelTimeCam.h"
90#include "MHCalibrationPix.h"
91
92#include "MLog.h"
93#include "MLogManip.h"
94
95#include "MParList.h"
96
97#include "MCalibrationRelTimeCam.h"
98#include "MCalibrationRelTimePix.h"
99#include "MCalibrationPix.h"
100
101#include "MArrivalTimeCam.h"
102#include "MArrivalTimePix.h"
103
104#include "MGeomCam.h"
105#include "MGeomPix.h"
106
107#include "MBadPixelsCam.h"
108#include "MBadPixelsPix.h"
109
110#include <TOrdCollection.h>
111#include <TPad.h>
112#include <TVirtualPad.h>
113#include <TCanvas.h>
114#include <TStyle.h>
115#include <TF1.h>
116#include <TLine.h>
117#include <TLatex.h>
118#include <TLegend.h>
119#include <TGraph.h>
120#include <TEnv.h>
121
122ClassImp(MHCalibrationRelTimeCam);
123
124using namespace std;
125
126const Float_t MHCalibrationRelTimeCam::fgNumHiGainSaturationLimit = 0.25;
127const UInt_t MHCalibrationRelTimeCam::fgReferencePixel = 1;
128const Int_t MHCalibrationRelTimeCam::fgNbins = 400;
129const Axis_t MHCalibrationRelTimeCam::fgFirst = -9.975;
130const Axis_t MHCalibrationRelTimeCam::fgLast = 10.025;
131const Float_t MHCalibrationRelTimeCam::fgProbLimit = 0.0;
132const TString MHCalibrationRelTimeCam::gsHistName = "RelTime";
133const TString MHCalibrationRelTimeCam::gsHistTitle = "Arr. Times";
134const TString MHCalibrationRelTimeCam::gsHistXTitle = "Arr. Time [FADC slices]";
135const TString MHCalibrationRelTimeCam::gsHistYTitle = "Nr. events";
136const TString MHCalibrationRelTimeCam::fgReferenceFile = "mjobs/calibrationref.rc";
137
138// --------------------------------------------------------------------------
139//
140// Default Constructor.
141//
142// Sets:
143// - fReferencePixel to fgReferencePixel
144// - fNbins to fgNbins
145// - fFirst to fgFirst
146// - fLast to fgLast
147//
148// - fHistName to gsHistName
149// - fHistTitle to gsHistTitle
150// - fHistXTitle to gsHistXTitle
151// - fHistYTitle to gsHistYTitle
152//
153MHCalibrationRelTimeCam::MHCalibrationRelTimeCam(const char *name, const char *title)
154{
155
156 fName = name ? name : "MHCalibrationRelTimeCam";
157 fTitle = title ? title : "Histogram class for the relative time calibration of the camera";
158
159 SetNumHiGainSaturationLimit(fgNumHiGainSaturationLimit);
160
161 SetReferencePixel();
162
163 SetBinning(fgNbins, fgFirst, fgLast);
164
165 SetProbLimit(fgProbLimit);
166
167 SetHistName (gsHistName .Data());
168 SetHistTitle (gsHistTitle .Data());
169 SetHistXTitle(gsHistXTitle.Data());
170 SetHistYTitle(gsHistYTitle.Data());
171
172 SetReferenceFile();
173
174 fInnerRefTime = 2.95;
175 fOuterRefTime = 3.6;
176}
177
178// --------------------------------------------------------------------------
179//
180// Creates new MHCalibrationRelTimeCam only with the averaged areas:
181// the rest has to be retrieved directly, e.g. via:
182// MHCalibrationRelTimeCam *cam = MParList::FindObject("MHCalibrationRelTimeCam");
183// - cam->GetAverageSector(5).DrawClone();
184// - (*cam)[100].DrawClone()
185//
186TObject *MHCalibrationRelTimeCam::Clone(const char *) const
187{
188
189 MHCalibrationRelTimeCam *cam = new MHCalibrationRelTimeCam();
190
191 //
192 // Copy the data members
193 //
194 cam->fColor = fColor;
195 cam->fRunNumbers = fRunNumbers;
196 cam->fPulserFrequency = fPulserFrequency;
197 cam->fFlags = fFlags;
198 cam->fNbins = fNbins;
199 cam->fFirst = fFirst;
200 cam->fLast = fLast;
201
202 cam->fReferenceFile = fReferenceFile;
203 cam->fInnerRefTime = fInnerRefTime;
204 cam->fOuterRefTime = fOuterRefTime;
205
206 //
207 // Copy the MArrays
208 //
209 cam->fAverageAreaRelSigma = fAverageAreaRelSigma;
210 cam->fAverageAreaRelSigmaVar = fAverageAreaRelSigmaVar;
211 cam->fAverageAreaSat = fAverageAreaSat;
212 cam->fAverageAreaSigma = fAverageAreaSigma;
213 cam->fAverageAreaSigmaVar = fAverageAreaSigmaVar;
214 cam->fAverageAreaNum = fAverageAreaNum;
215 cam->fAverageSectorNum = fAverageSectorNum;
216
217 if (!IsAverageing())
218 return cam;
219
220 const Int_t navhi = fAverageHiGainAreas->GetSize();
221
222 for (int i=0; i<navhi; i++)
223 cam->fAverageHiGainAreas->AddAt(GetAverageHiGainArea(i).Clone(),i);
224
225 if (IsLoGain())
226 {
227
228 const Int_t navlo = fAverageLoGainAreas->GetSize();
229 for (int i=0; i<navlo; i++)
230 cam->fAverageLoGainAreas->AddAt(GetAverageLoGainArea(i).Clone(),i);
231
232 }
233
234 return cam;
235}
236
237// --------------------------------------------------------------------------
238//
239// Gets or creates the pointers to:
240// - MCalibrationRelTimeCam
241//
242// Searches pointer to:
243// - MArrivalTimeCam
244//
245// Calls:
246// - MHCalibrationCam::InitHiGainArrays()
247// - MHCalibrationCam::InitLoGainArrays()
248//
249// Sets:
250// - fSumareahi to nareas
251// - fSumarealo to nareas
252// - fSumsectorhi to nareas
253// - fSumsectorlo to nareas
254// - fNumareahi to nareas
255// - fNumarealo to nareas
256// - fNumsectorhi to nareas
257// - fNumsectorlo to nareas
258//
259Bool_t MHCalibrationRelTimeCam::ReInitHists(MParList *pList)
260{
261
262 if (!InitCams(pList,"RelTime"))
263 return kFALSE;
264
265 MArrivalTimeCam *signal = (MArrivalTimeCam*)pList->FindObject("MArrivalTimeCam");
266 if (!signal)
267 {
268 *fLog << err << "MArrivalTimeCam not found... abort." << endl;
269 return kFALSE;
270 }
271
272 const Int_t npixels = fGeom->GetNumPixels();
273 const Int_t nsectors = fGeom->GetNumSectors();
274 const Int_t nareas = fGeom->GetNumAreas();
275
276 InitHiGainArrays(npixels,nareas,nsectors);
277 InitLoGainArrays(npixels,nareas,nsectors);
278
279 fSumareahi .Set(nareas);
280 fSumarealo .Set(nareas);
281 fSumsectorhi.Set(nsectors);
282 fSumsectorlo.Set(nsectors);
283 fNumareahi .Set(nareas);
284 fNumarealo .Set(nareas);
285 fNumsectorhi.Set(nsectors);
286 fNumsectorlo.Set(nsectors);
287
288 return kTRUE;
289}
290
291
292// -------------------------------------------------------------------------------
293//
294// Retrieves pointer to MArrivalTimeCam:
295//
296// Retrieves from MGeomCam:
297// - number of pixels
298// - number of pixel areas
299// - number of sectors
300//
301// Fills HiGain or LoGain histograms (MHGausEvents::FillHistAndArray()), respectively
302// depending on MArrivalTimePix::IsLoGainUsed(), with:
303// - MArrivalTimePix::GetArrivalTime(pixid) - MArrivalTimePix::GetArrivalTime(1);
304// (i.e. the time difference between pixel i and pixel 1 (hardware number: 2) )
305//
306Bool_t MHCalibrationRelTimeCam::FillHists(const MParContainer *par, const Stat_t w)
307{
308
309 MArrivalTimeCam *arrtime = (MArrivalTimeCam*)par;
310 if (!arrtime)
311 {
312 gLog << err << "No argument in MArrivalTime::Fill... abort." << endl;
313 return kFALSE;
314 }
315
316 const Int_t npixels = fGeom->GetNumPixels();
317 const Int_t nareas = fGeom->GetNumAreas();
318 const Int_t nsectors = fGeom->GetNumSectors();
319
320 fSumareahi .Reset();
321 fSumarealo .Reset();
322 fSumsectorhi.Reset();
323 fSumsectorlo.Reset();
324 fNumareahi .Reset();
325 fNumarealo .Reset();
326 fNumsectorhi.Reset();
327 fNumsectorlo.Reset();
328
329 const MArrivalTimePix &refpix = (*arrtime)[fReferencePixel];
330 // FIXME: What do we do if pixel is invalid?
331 if (!refpix.IsArrivalTimeValid())
332 {
333 gLog << warn << "WARNING - Arrival time in refrence pixel " << fReferencePixel << " invalid." << endl;
334 return kTRUE;
335 }
336
337 const Float_t reftime = refpix.GetArrivalTime();
338// refpix.IsHiGainSaturated()
339// ? refpix.GetArrivalTimeLoGain() : refpix.GetArrivalTimeHiGain();
340
341 for (Int_t i=0; i<npixels; i++)
342 {
343 MHCalibrationPix &histhi = (*this)[i];
344 if (histhi.IsExcluded())
345 continue;
346
347 const Int_t aidx = (*fGeom)[i].GetAidx();
348 const Int_t sector = (*fGeom)[i].GetSector();
349
350 const MArrivalTimePix &pix = (*arrtime)[i];
351
352 // If hi-gain arrival time has been extracted successfully
353 // fill hi-gain histograms and arrays
354 if (pix.IsHiGainValid() && !pix.IsHiGainSaturated())
355 {
356 const Float_t time = pix.GetArrivalTimeHiGain();
357
358 if (IsOscillations())
359 histhi.FillHistAndArray(time-reftime);
360 else
361 histhi.FillHist(time-reftime);
362
363 fSumareahi [aidx] += time;
364 fNumareahi [aidx] ++;
365 fSumsectorhi[sector] += time;
366 fNumsectorhi[sector] ++;
367 }
368
369 if (!pix.IsHiGainSaturated())
370 continue;
371
372 histhi.AddSaturated(1);
373/*
374 // If lo-gain arrival time has been extracted successfully,
375 // the hi-gain has saturateed and the lo-gain is switched on
376 // fill hi-gain histograms and arrays
377 if (pix.IsLoGainValid() && IsLoGain())
378 {
379 const Float_t time = pix.GetArrivalTimeLoGain();
380
381 MHCalibrationPix &histlo = (*this)(i);
382 if (IsOscillations())
383 histlo.FillHistAndArray(time-reftime);
384 else
385 histlo.FillHist(time-reftime);
386
387 fSumarealo [aidx] += time;
388 fNumarealo [aidx] ++;
389 fSumsectorlo[sector] += time;
390 fNumsectorlo[sector] ++;
391 }*/
392 }
393
394 for (Int_t j=0; j<nareas; j++)
395 {
396 MHCalibrationPix &histhi = GetAverageHiGainArea(j);
397 if (IsOscillations())
398 histhi.FillHistAndArray(fNumareahi[j] == 0 ? 0. : fSumareahi[j]/fNumareahi[j]);
399 else
400 histhi.FillHist(fNumareahi[j] == 0 ? 0. : fSumareahi[j]/fNumareahi[j]);
401
402 if (IsLoGain())
403 {
404 MHCalibrationPix &histlo = GetAverageLoGainArea(j);
405 if (IsOscillations())
406 histlo.FillHistAndArray(fNumarealo[j] == 0 ? 0. : fSumarealo[j]/fNumarealo[j]);
407 else
408 histlo.FillHist(fNumarealo[j] == 0 ? 0. : fSumarealo[j]/fNumarealo[j]);
409 }
410 }
411
412 for (Int_t j=0; j<nsectors; j++)
413 {
414 MHCalibrationPix &histhi = GetAverageHiGainSector(j);
415 if (IsOscillations())
416 histhi.FillHistAndArray(fNumsectorhi[j] == 0 ? 0. : fSumsectorhi[j]/fNumsectorhi[j]);
417 else
418 histhi.FillHist(fNumsectorhi[j] == 0 ? 0. : fSumsectorhi[j]/fNumsectorhi[j]);
419
420 if (IsLoGain())
421 {
422 MHCalibrationPix &histlo = GetAverageLoGainSector(j);
423 if (IsOscillations())
424 histlo.FillHistAndArray(fNumsectorlo[j] == 0 ? 0. : fSumsectorlo[j]/fNumsectorlo[j]);
425 else
426 histlo.FillHist(fNumsectorlo[j] == 0 ? 0. : fSumsectorlo[j]/fNumsectorlo[j]);
427 }
428 }
429
430 return kTRUE;
431}
432
433// --------------------------------------------------------------------------
434//
435// Calls:
436// - MHCalibrationCam::FitHiGainArrays() with flags:
437// MBadPixelsPix::kRelTimeNotFitted and MBadPixelsPix::kRelTimeOscillating
438// - MHCalibrationCam::FitLoGainArrays() with flags:
439// MBadPixelsPix::kRelTimeNotFitted and MBadPixelsPix::kRelTimeOscillating
440//
441Bool_t MHCalibrationRelTimeCam::FinalizeHists()
442{
443
444 *fLog << endl;
445
446 const Int_t nareas = fAverageHiGainAreas->GetSize();
447 const Int_t nsectors = fAverageHiGainSectors->GetSize();
448
449 TArrayI satarea(nareas);
450 TArrayI satsect(nsectors);
451 fNumareahi .Reset();
452 fNumsectorhi.Reset();
453
454 for (Int_t i=0; i<fHiGainArray->GetSize(); i++)
455 {
456
457 MHCalibrationPix &histhi = (*this)[i];
458
459 if (histhi.IsExcluded())
460 continue;
461
462 const Int_t aidx = (*fGeom)[i].GetAidx();
463 const Int_t sector = (*fGeom)[i].GetSector();
464
465 fNumareahi[aidx]++;
466 fNumsectorhi[sector]++;
467
468 //
469 // Check saturation
470 //
471 MCalibrationRelTimePix &pix = (MCalibrationRelTimePix&)(*fCam)[i] ;
472 if (histhi.GetSaturated() > fNumHiGainSaturationLimit*histhi.GetHGausHist()->GetEntries())
473 {
474 pix.SetHiGainSaturation();
475 histhi.SetExcluded();
476 satarea[aidx]++;
477 satsect[sector]++;
478 }
479 else
480 if (IsLoGain())
481 (*this)(i).SetExcluded();
482
483 //
484 // Check histogram overflow
485 //
486 CheckOverflow(histhi);
487 if (IsLoGain())
488 CheckOverflow((*this)(i));
489
490 }
491
492 for (Int_t j=0; j<nareas; j++)
493 {
494
495 MHCalibrationPix &histhi = GetAverageHiGainArea(j);
496 MCalibrationRelTimePix &pix = (MCalibrationRelTimePix&)fCam->GetAverageArea(j);
497
498 if (satarea[j] > 0.5*fNumareahi[j])
499 {
500 pix.SetHiGainSaturation();
501 histhi.SetExcluded();
502 }
503 else
504 if (IsLoGain())
505 GetAverageLoGainArea(j).SetExcluded();
506
507 //
508 // Check histogram overflow
509 //
510 CheckOverflow(histhi);
511 if (IsLoGain())
512 CheckOverflow(GetAverageLoGainArea(j));
513 }
514
515 for (Int_t j=0; j<fAverageHiGainSectors->GetSize(); j++)
516 {
517
518 MHCalibrationPix &histhi = GetAverageHiGainSector(j);
519
520 if (satsect[j] > 0.5*fNumsectorhi[j])
521 {
522 MCalibrationRelTimePix &pix = (MCalibrationRelTimePix&)fCam->GetAverageSector(j) ;
523 pix.SetHiGainSaturation();
524 histhi.SetExcluded();
525 }
526 else
527 if (IsLoGain())
528 GetAverageLoGainSector(j).SetExcluded();
529
530 //
531 // Check histogram overflow
532 //
533 CheckOverflow(histhi);
534 if (IsLoGain())
535 CheckOverflow(GetAverageLoGainSector(j));
536 }
537
538 FitHiGainArrays(*fCam, *fBadPixels,
539 MBadPixelsPix::kRelTimeNotFitted,
540 MBadPixelsPix::kRelTimeOscillating);
541
542 if (IsLoGain())
543 FitLoGainArrays(*fCam, *fBadPixels,
544 MBadPixelsPix::kRelTimeNotFitted,
545 MBadPixelsPix::kRelTimeOscillating);
546
547 return kTRUE;
548}
549
550// --------------------------------------------------------------------------
551//
552// Sets all pixels to MBadPixelsPix::kUnreliableRun, if following flags are set:
553// - MBadPixelsPix::kRelTimeNotFitted
554// - MBadPixelsPix::kRelTimeOscillating
555//
556void MHCalibrationRelTimeCam::FinalizeBadPixels()
557{
558
559 for (Int_t i=0; i<fBadPixels->GetSize(); i++)
560 {
561 MBadPixelsPix &bad = (*fBadPixels)[i];
562
563 if (bad.IsUncalibrated( MBadPixelsPix::kRelTimeNotFitted ))
564 bad.SetUnsuitable( MBadPixelsPix::kUnreliableRun );
565
566 if (bad.IsUncalibrated( MBadPixelsPix::kRelTimeOscillating))
567 bad.SetUnsuitable( MBadPixelsPix::kUnreliableRun );
568
569 }
570}
571
572// --------------------------------------------------------------------------
573//
574// The types are as follows:
575//
576// Fitted values:
577// ==============
578//
579// 0: Fitted Mean Relative Arrival Time in FADC slices (MHGausEvents::GetMean()
580// 1: Error Mean Relative Arrival Time in FADC slices (MHGausEvents::GetMeanErr()
581// 2: Sigma fitted Relative Arrival Time in FADC slices (MHGausEvents::GetSigma()
582// 3: Error Sigma Relative Arrival Time in FADC slices (MHGausEvents::GetSigmaErr()
583//
584// Useful variables derived from the fit results:
585// =============================================
586//
587// 4: Returned probability of Gauss fit (calls: MHGausEvents::GetProb())
588//
589// Localized defects:
590// ==================
591//
592// 5: Gaus fit not OK (calls: MHGausEvents::IsGausFitOK())
593// 6: Fourier spectrum not OK (calls: MHGausEvents::IsFourierSpectrumOK())
594//
595Bool_t MHCalibrationRelTimeCam::GetPixelContent(Double_t &val, Int_t idx, const MGeomCam &cam, Int_t type) const
596{
597
598 if (fHiGainArray->GetSize() <= idx)
599 return kFALSE;
600
601 const MHCalibrationPix &pix = (*this)[idx];
602
603 switch (type)
604 {
605 case 0:
606 val = pix.GetMean();
607 break;
608 case 1:
609 val = pix.GetMeanErr();
610 break;
611 case 2:
612 val = pix.GetSigma();
613 break;
614 case 3:
615 val = pix.GetSigmaErr();
616 break;
617 case 4:
618 val = pix.GetProb();
619 break;
620 case 5:
621 if (!pix.IsGausFitOK())
622 val = 1.;
623 break;
624 case 6:
625 if (!pix.IsFourierSpectrumOK())
626 val = 1.;
627 break;
628 default:
629 return kFALSE;
630 }
631 return kTRUE;
632}
633
634// --------------------------------------------------------------------------
635//
636// Calls MHCalibrationPix::DrawClone() for pixel idx
637//
638void MHCalibrationRelTimeCam::DrawPixelContent(Int_t idx) const
639{
640 (*this)[idx].DrawClone();
641}
642
643// -----------------------------------------------------------------------------
644//
645// Default draw:
646//
647// Displays the averaged areas, both High Gain and Low Gain
648//
649// Calls the Draw of the fAverageHiGainAreas and fAverageLoGainAreas objects with options
650//
651void MHCalibrationRelTimeCam::Draw(const Option_t *opt)
652{
653
654 const Int_t nareas = fAverageHiGainAreas->GetSize();
655 if (nareas == 0)
656 return;
657
658 TString option(opt);
659 option.ToLower();
660
661 if (!option.Contains("datacheck"))
662 {
663 MHCalibrationCam::Draw(opt);
664 return;
665 }
666
667 //
668 // From here on , the datacheck - Draw
669 //
670 TVirtualPad *pad = gPad ? gPad : MH::MakeDefCanvas(this);
671 pad->SetBorderMode(0);
672 pad->Divide(1,nareas);
673
674 //
675 // Loop over inner and outer pixels
676 //
677 for (Int_t i=0; i<nareas;i++)
678 {
679
680 pad->cd(i+1);
681
682 MHCalibrationPix &hipix = GetAverageHiGainArea(i);
683 //
684 // Ask for Hi-Gain saturation
685 //
686 if (hipix.IsExcluded() && IsLoGain())
687 {
688 MHCalibrationPix &lopix = GetAverageLoGainArea(i);
689 DrawDataCheckPixel(lopix,i ? fOuterRefTime+1.5 : fInnerRefTime+1.5);
690 }
691 else
692 DrawDataCheckPixel(hipix,i ? fOuterRefTime : fInnerRefTime);
693 }
694}
695
696void MHCalibrationRelTimeCam::CheckOverflow( MHCalibrationPix &pix ) const
697{
698 if (pix.IsExcluded())
699 return;
700
701 const TH1F &hist = *pix.GetHGausHist();
702
703 const Int_t n = hist.GetNbinsX();
704 const Float_t max = fOverflowLimit*hist.GetEntries();
705
706 const Stat_t overflow = hist.GetBinContent(n+1);
707 if (overflow > max)
708 {
709 *fLog << warn << overflow << " overflows above " << hist.GetBinLowEdge(n);
710 *fLog << " in " << pix.GetName() << " (w/o saturation!) " << endl;
711 }
712
713 const Stat_t underflow = hist.GetBinContent(0);
714 if (underflow > max)
715 {
716 *fLog << warn << underflow << " underflows below " << hist.GetBinLowEdge(1);
717 *fLog << " in " << pix.GetName() << " (w/o saturation!) " << endl;
718 }
719}
720
721
722// -----------------------------------------------------------------------------
723//
724// Draw the average pixel for the datacheck:
725//
726// Displays the averaged areas, both High Gain and Low Gain
727//
728// Calls the Draw of the fAverageHiGainAreas and fAverageLoGainAreas objects with options
729//
730void MHCalibrationRelTimeCam::DrawDataCheckPixel(MHCalibrationPix &pix, const Float_t refline)
731{
732 if (pix.IsEmpty())
733 return;
734
735 TVirtualPad *pad = gPad;
736 pad->Divide(1,2, 1e-10, 1e-10);
737 pad->cd(1);
738
739 gPad->SetBorderMode(0);
740
741 gPad->SetTicks();
742 if (!pix.IsEmpty() && !pix.IsOnlyOverflow() && !pix.IsOnlyUnderflow())
743 gPad->SetLogy();
744
745 TH1F *hist = pix.GetHGausHist();
746
747 //
748 // set the labels bigger
749 //
750 TAxis *xaxe = hist->GetXaxis();
751 TAxis *yaxe = hist->GetYaxis();
752 xaxe->CenterTitle();
753 yaxe->CenterTitle();
754 xaxe->SetTitleSize(0.07);
755 yaxe->SetTitleSize(0.07);
756 xaxe->SetTitleOffset(0.65);
757 yaxe->SetTitleOffset(0.55);
758 xaxe->SetLabelSize(0.06);
759 yaxe->SetLabelSize(0.06);
760 xaxe->SetTitle(hist->GetXaxis()->GetTitle());
761 yaxe->SetTitle(hist->GetYaxis()->GetTitle());
762 xaxe->SetRange(hist->GetMaximumBin()-1500, hist->GetMaximumBin()+1500);
763
764 hist->Draw();
765
766 gStyle->SetOptFit();
767
768 TF1 *fit = pix.GetFGausFit();
769
770 if (fit)
771 {
772 switch (fColor)
773 {
774 case MCalibrationCam::kGREEN:
775 fit->SetLineColor(kGreen);
776 break;
777 case MCalibrationCam::kBLUE:
778 fit->SetLineColor(kBlue);
779 break;
780 case MCalibrationCam::kUV:
781 fit->SetLineColor(106);
782 break;
783 case MCalibrationCam::kCT1:
784 fit->SetLineColor(006);
785 break;
786 default:
787 fit->SetLineColor(kRed);
788 }
789 fit->Draw("same");
790 }
791
792 DisplayRefLines(hist, refline);
793
794 pad->cd(2);
795 gPad->SetBorderMode(0);
796 gPad->SetTicks();
797
798 pix.CreateGraphEvents();
799 TGraph *gr = pix.GetGraphEvents();
800 if (gr)
801 {
802 TH1F *null2 = gr->GetHistogram();
803
804 null2->SetMinimum(pix.GetMean()-10.*pix.GetSigma());
805 null2->SetMaximum(pix.GetMean()+10.*pix.GetSigma());
806
807 //
808 // set the labels bigger
809 //
810 TAxis *xaxe2 = null2->GetXaxis();
811 TAxis *yaxe2 = null2->GetYaxis();
812 xaxe2->CenterTitle();
813 yaxe2->CenterTitle();
814 xaxe2->SetTitleSize(0.07);
815 yaxe2->SetTitleSize(0.07);
816 xaxe2->SetTitleOffset(0.65);
817 yaxe2->SetTitleOffset(0.55);
818 xaxe2->SetLabelSize(0.06);
819 yaxe2->SetLabelSize(0.06);
820 xaxe2->SetRangeUser(0.,pix.GetEvents()->GetSize()/pix.GetEventFrequency());
821 }
822
823 pix.DrawEvents();
824}
825
826void MHCalibrationRelTimeCam::DisplayRefLines(const TH1F *hist, const Float_t refline) const
827{
828 TLine *line = new TLine(refline, 0, refline, hist->GetMaximum());
829 line->SetLineColor(kGreen);
830 line->SetLineStyle(2);
831 line->SetLineWidth(3);
832 line->SetBit(kCanDelete);
833 line->Draw();
834
835 TLegend *leg = new TLegend(0.75,0.01,0.99,0.3);
836 leg->SetBit(kCanDelete);
837 leg->AddEntry(line, "Trigger Calibration", "l");
838 leg->Draw();
839}
840
841Int_t MHCalibrationRelTimeCam::ReadEnv(const TEnv &env, TString prefix, Bool_t print)
842{
843 Bool_t rc = kFALSE;
844
845 if (IsEnvDefined(env, prefix, "ReferenceFile", print))
846 {
847 SetReferenceFile(GetEnvValue(env, prefix, "ReferenceFile", fReferenceFile.Data()));
848 rc = kTRUE;
849 }
850
851 if (IsEnvDefined(env, prefix, "ReferencePixel", print))
852 {
853 SetReferencePixel(GetEnvValue(env, prefix, "ReferencePixel", (Int_t)fReferencePixel));
854 rc = kTRUE;
855 }
856
857 TEnv refenv(fReferenceFile);
858
859 fInnerRefTime = refenv.GetValue("InnerRefTime", fInnerRefTime);
860 fOuterRefTime = refenv.GetValue("OuterRefTime", fOuterRefTime);
861
862 return MHCalibrationCam::ReadEnv(env,prefix,print) ? kTRUE : rc;
863}
Note: See TracBrowser for help on using the repository browser.