source: trunk/MagicSoft/Mars/mhcalib/MHCalibrationCam.cc@ 5851

Last change on this file since 5851 was 5851, checked in by gaug, 20 years ago
*** empty log message ***
File size: 45.0 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!
20! Copyright: MAGIC Software Development, 2000-2004
21!
22!
23\* ======================================================================== */
24/////////////////////////////////////////////////////////////////////////////
25//
26// MHCalibrationCam
27//
28// Base class for camera calibration classes. Incorporates the TOrdCollection's:
29// - fHiGainArray (for calibrated High Gains per pixel)
30// - fLoGainArray (for calibrated Low Gains per pixel)
31// - fAverageHiGainAreas (for averaged High Gains events per camera area index)
32// - fAverageLoGainAreas (for averaged High Gains events per camera area index)
33// - fAverageHiGainSectors (for averaged High Gains events per camera sector )
34// - fAverageLoGainSectors (for averaged High Gains events per camera sector )
35// These TOrdCollection's are called by their default constructors, thus no objects
36// are created, until the derived class does so.
37//
38// The corresponding operators: [],() and the operators GetAverageHiGainArea(),
39// GetAverageLoGainArea(), GetAverageHiGainSector() and GetAverageLoGainSector()
40// have to be cast to the corresponding class. It is assumed that all classes
41// dealing with calibration pixels derive from MHCalibrationPix.
42//
43// The following flag can be set:
44// - SetAverageing() for calculating the event-by-event averages.
45// - SetDebug() for debug output.
46// - SetLoGain() for the case that low-gain slices are available, but
47// MRawRunHeader::GetNumLoGainSlices() gives still 0.
48// - SetCheckSize() for testing the sizes of the f*Arrays in ReInit() if they
49// are coinciding with the equiv. sizes in MGeomCam
50//
51// The flag kLoGain steers if the low-gain signal is treated at all or not.
52// The flag kAverageing steers if the event-by-event averages are treated at all.
53//
54/////////////////////////////////////////////////////////////////////////////
55#include "MHCalibrationCam.h"
56#include "MHCalibrationPix.h"
57
58#include <TVirtualPad.h>
59#include <TCanvas.h>
60#include <TPad.h>
61#include <TText.h>
62#include <TPaveText.h>
63#include <TOrdCollection.h>
64#include <TROOT.h>
65
66#include "MLog.h"
67#include "MLogManip.h"
68
69#include "MCalibrationIntensityCam.h"
70#include "MCalibrationCam.h"
71#include "MCalibrationPix.h"
72
73#include "MBadPixelsIntensityCam.h"
74#include "MBadPixelsCam.h"
75#include "MBadPixelsPix.h"
76
77#include "MGeomCam.h"
78#include "MGeomPix.h"
79
80#include "MParList.h"
81
82#include "MRawRunHeader.h"
83#include "MTriggerPattern.h"
84
85ClassImp(MHCalibrationCam);
86
87using namespace std;
88
89const Int_t MHCalibrationCam::fgPulserFrequency = 500;
90const Float_t MHCalibrationCam::fgProbLimit = 0.0001;
91const TString MHCalibrationCam::gsHistName = "Hist";
92const TString MHCalibrationCam::gsHistTitle = "";
93const TString MHCalibrationCam::gsHistXTitle = "";
94const TString MHCalibrationCam::gsHistYTitle = "Nr. events";
95// --------------------------------------------------------------------------
96//
97// Default Constructor.
98//
99// Sets:
100// - all pointers to NULL
101//
102// Initializes and sets owner of:
103// - fHiGainArray, fLoGainArray
104// - fAverageHiGainAreas, fAverageLoGainAreas
105// - fAverageHiGainSectors, fAverageLoGainSectors
106//
107// Initializes:
108// - fPulserFrequency to fgPulserFrequency
109// - fProbLimit to fgProbLimit
110//
111// - SetAveregeing (kTRUE);
112// - SetDebug (kFALSE);
113// - SetLoGain (kTRUE);
114//- SetOscillations(kTRUE);
115//- SetSizeCheck (kTRUE);
116//- SetInterlaced (kFALSE);
117//
118MHCalibrationCam::MHCalibrationCam(const char *name, const char *title)
119 : fHistName(gsHistName),fHistTitle(gsHistTitle),
120 fHistXTitle(gsHistXTitle),fHistYTitle(gsHistYTitle),
121 fColor(MCalibrationCam::kNONE), fIntensBad(NULL),
122 fBadPixels(NULL), fIntensCam(NULL), fCam(NULL), fGeom(NULL),
123 fRunHeader(NULL), fTrigPattern(NULL), fInterlacedFlags(0),
124 fMaskInterlaced(0)
125{
126
127 fHiGainArray = new TOrdCollection;
128 fHiGainArray->SetOwner();
129
130 fLoGainArray = new TOrdCollection;
131 fLoGainArray->SetOwner();
132
133 fAverageHiGainAreas = new TOrdCollection;
134 fAverageHiGainAreas->SetOwner();
135
136 fAverageLoGainAreas = new TOrdCollection;
137 fAverageLoGainAreas->SetOwner();
138
139 fAverageHiGainSectors = new TOrdCollection;
140 fAverageHiGainSectors->SetOwner();
141
142 fAverageLoGainSectors = new TOrdCollection;
143 fAverageLoGainSectors->SetOwner();
144
145 SetPulserFrequency();
146 SetProbLimit();
147
148 SetAverageing (kTRUE);
149 SetDebug (kFALSE);
150 SetLoGain (kTRUE);
151 SetOscillations(kTRUE);
152 SetSizeCheck (kTRUE);
153 SetInterlaced (kFALSE);
154}
155
156// --------------------------------------------------------------------------
157//
158// Deletes the TOrdCollection of:
159// - fHiGainArray, fLoGainArray
160// - fAverageHiGainAreas, fAverageLoGainAreas
161// - fAverageHiGainSectors, fAverageLoGainSectors
162//
163MHCalibrationCam::~MHCalibrationCam()
164{
165
166 delete fHiGainArray;
167 delete fLoGainArray;
168
169 delete fAverageHiGainAreas;
170 delete fAverageLoGainAreas;
171
172 delete fAverageHiGainSectors;
173 delete fAverageLoGainSectors;
174 /*
175
176 Remove ( fHiGainArray );
177 Remove ( fLoGainArray );
178
179 Remove ( fAverageHiGainAreas );
180 Remove ( fAverageLoGainAreas );
181
182 Remove ( fAverageHiGainSectors );
183 Remove ( fAverageLoGainSectors );
184 */
185
186}
187
188void MHCalibrationCam::Remove(TOrdCollection *col)
189{
190
191 if (!col)
192 return;
193
194 TOrdCollectionIter Next(col);
195
196 Int_t count = 0;
197
198 while (MHCalibrationPix *obj = (MHCalibrationPix*)Next())
199 {
200 *fLog << ++count << " " << obj << flush;
201 if (obj && obj->IsOnHeap())
202 {
203 obj->Draw();
204 delete obj;
205 }
206 }
207
208 delete col;
209}
210
211// --------------------------------------------------------------------------
212//
213// Returns size of fHiGainArray
214//
215const Int_t MHCalibrationCam::GetSize() const
216{
217 return fHiGainArray->GetSize();
218}
219
220// --------------------------------------------------------------------------
221//
222// Get i-th High Gain pixel (pixel number)
223//
224MHCalibrationPix &MHCalibrationCam::operator[](UInt_t i)
225{
226 return *static_cast<MHCalibrationPix*>(fHiGainArray->At(i));
227}
228
229// --------------------------------------------------------------------------
230//
231// Get i-th High Gain pixel (pixel number)
232//
233const MHCalibrationPix &MHCalibrationCam::operator[](UInt_t i) const
234{
235 return *static_cast<MHCalibrationPix*>(fHiGainArray->At(i));
236}
237
238// --------------------------------------------------------------------------
239//
240// Get i-th Low Gain pixel (pixel number)
241//
242MHCalibrationPix &MHCalibrationCam::operator()(UInt_t i)
243{
244 return *static_cast<MHCalibrationPix*>(fLoGainArray->At(i));
245}
246
247// --------------------------------------------------------------------------
248//
249// Get i-th Low Gain pixel (pixel number)
250//
251const MHCalibrationPix &MHCalibrationCam::operator()(UInt_t i) const
252{
253 return *static_cast<MHCalibrationPix*>(fLoGainArray->At(i));
254}
255
256// --------------------------------------------------------------------------
257//
258// Returns the current size of the TOrdCollection fAverageHiGainAreas
259// independently if the MHCalibrationPix is filled with values or not.
260//
261const Int_t MHCalibrationCam::GetAverageAreas() const
262{
263 return fAverageHiGainAreas->GetSize();
264}
265
266// --------------------------------------------------------------------------
267//
268// Get i-th High Gain pixel Area (area number)
269//
270MHCalibrationPix &MHCalibrationCam::GetAverageHiGainArea(UInt_t i)
271{
272 return *static_cast<MHCalibrationPix*>(fAverageHiGainAreas->At(i));
273}
274
275// --------------------------------------------------------------------------
276//
277// Get i-th High Gain pixel Area (area number)
278//
279const MHCalibrationPix &MHCalibrationCam::GetAverageHiGainArea(UInt_t i) const
280{
281 return *static_cast<MHCalibrationPix *>(fAverageHiGainAreas->At(i));
282}
283
284// --------------------------------------------------------------------------
285//
286// Get i-th Low Gain pixel Area (area number)
287//
288MHCalibrationPix &MHCalibrationCam::GetAverageLoGainArea(UInt_t i)
289{
290 return *static_cast<MHCalibrationPix*>(fAverageLoGainAreas->At(i));
291}
292
293// --------------------------------------------------------------------------
294//
295// Get i-th Low Gain pixel Area (area number)
296//
297const MHCalibrationPix &MHCalibrationCam::GetAverageLoGainArea(UInt_t i) const
298{
299 return *static_cast<MHCalibrationPix*>(fAverageLoGainAreas->At(i));
300}
301
302// --------------------------------------------------------------------------
303//
304// Returns the current size of the TOrdCollection fAverageHiGainSectors
305// independently if the MHCalibrationPix is filled with values or not.
306//
307const Int_t MHCalibrationCam::GetAverageSectors() const
308{
309 return fAverageHiGainSectors->GetSize();
310}
311
312// --------------------------------------------------------------------------
313//
314// Get i-th High Gain Sector (sector number)
315//
316MHCalibrationPix &MHCalibrationCam::GetAverageHiGainSector(UInt_t i)
317{
318 return *static_cast<MHCalibrationPix*>(fAverageHiGainSectors->At(i));
319}
320
321// --------------------------------------------------------------------------
322//
323// Get i-th High Gain Sector (sector number)
324//
325const MHCalibrationPix &MHCalibrationCam::GetAverageHiGainSector(UInt_t i) const
326{
327 return *static_cast<MHCalibrationPix*>(fAverageHiGainSectors->At(i));
328}
329
330// --------------------------------------------------------------------------
331//
332// Get i-th Low Gain Sector (sector number)
333//
334MHCalibrationPix &MHCalibrationCam::GetAverageLoGainSector(UInt_t i)
335{
336 return *static_cast<MHCalibrationPix*>(fAverageLoGainSectors->At(i));
337}
338
339// --------------------------------------------------------------------------
340//
341// Get i-th Low Gain Sector (sector number)
342//
343const MHCalibrationPix &MHCalibrationCam::GetAverageLoGainSector(UInt_t i) const
344{
345 return *static_cast<MHCalibrationPix*>(fAverageLoGainSectors->At(i));
346}
347
348// --------------------------------------------------------------------------
349//
350// Calls ResetHistTitles()
351//
352// Calls Reset() for each entry in:
353// - fHiGainArray, fLoGainArray
354// - fAverageHiGainAreas, fAverageLoGainAreas
355// - fAverageHiGainSectors, fAverageLoGainSectors
356//
357void MHCalibrationCam::ResetHists()
358{
359
360 ResetHistTitles();
361
362 if (fHiGainArray)
363 { fHiGainArray->ForEach(MHCalibrationPix,Reset)(); }
364
365 if (IsAverageing())
366 {
367 if (fAverageHiGainAreas)
368 { fAverageHiGainAreas->ForEach(MHCalibrationPix,Reset)(); }
369 if (fAverageHiGainSectors)
370 { fAverageHiGainSectors->ForEach(MHCalibrationPix,Reset)(); }
371 }
372
373 if (!IsLoGain())
374 return;
375
376 if (fLoGainArray)
377 { fLoGainArray->ForEach(MHCalibrationPix,Reset)(); }
378 if (IsAverageing())
379 {
380 if (fAverageLoGainAreas)
381 { fAverageLoGainAreas->ForEach(MHCalibrationPix,Reset)(); }
382 if (fAverageLoGainSectors)
383 { fAverageLoGainSectors->ForEach(MHCalibrationPix,Reset)(); }
384 }
385}
386
387// --------------------------------------------------------------------------
388//
389// Resets the histogram titles for each entry in:
390// - fHiGainArray, fLoGainArray
391// - fAverageHiGainAreas, fAverageLoGainAreas
392// - fAverageHiGainSectors, fAverageLoGainSectors
393//
394void MHCalibrationCam::ResetHistTitles()
395{
396
397 TH1F *h;
398
399 if (fHiGainArray)
400 for (Int_t i=0;i<fHiGainArray->GetSize(); i++)
401 {
402 MHCalibrationPix &pix = (*this)[i];
403 h = pix.GetHGausHist();
404 h->SetName (Form("%s%s%s%4i","H",fHistName.Data(),"HiGainPix",i));
405 h->SetTitle(Form("%s%s%4i%s",fHistTitle.Data()," High Gain Pixel ",i," Runs: "));
406 h->SetXTitle(fHistXTitle.Data());
407 h->SetYTitle(fHistYTitle.Data());
408 }
409
410 if (IsAverageing())
411 {
412 if (fAverageHiGainAreas)
413 for (Int_t j=0; j<fAverageHiGainAreas->GetSize(); j++)
414 {
415 MHCalibrationPix &pix = GetAverageHiGainArea(j);
416 h = pix.GetHGausHist();
417 h->SetName (Form("%s%s%s%d","H",fHistName.Data(),"HiGainArea",j));
418 h->SetXTitle(fHistXTitle.Data());
419 h->SetYTitle(fHistYTitle.Data());
420 if (fGeom->InheritsFrom("MGeomCamMagic"))
421 h->SetTitle(Form("%s%s%s%s",fHistTitle.Data()," averaged on event-by-event basis ",
422 j==0 ? "Inner Pixels " : "Outer Pixels ","High Gain Runs: "));
423 else
424 h->SetTitle(Form("%s%s%d%s",fHistTitle.Data(),
425 " averaged on event-by-event basis High Gain Area Idx ",j," Runs: "));
426 }
427
428 if (fAverageHiGainSectors)
429 for (Int_t j=0; j<fAverageHiGainSectors->GetSize(); j++)
430 {
431 MHCalibrationPix &pix = GetAverageHiGainSector(j);
432 h = pix.GetHGausHist();
433 h->SetName (Form("%s%s%s%2i","H",fHistName.Data(),"HiGainSector",j));
434 h->SetTitle(Form("%s%s%2i%s",fHistTitle.Data(),
435 " averaged on event-by-event basis High Gain Sector ",j," Runs: "));
436 h->SetXTitle(fHistXTitle.Data());
437 h->SetYTitle(fHistYTitle.Data());
438 }
439 }
440
441 if (!IsLoGain())
442 return;
443
444 if (fLoGainArray)
445 for (Int_t i=0;i<fLoGainArray->GetSize(); i++)
446 {
447 MHCalibrationPix &pix = (*this)(i);
448 h = pix.GetHGausHist();
449 h->SetName (Form("%s%s%s%4i","H",fHistName.Data(),"LoGainPix",i));
450 h->SetTitle(Form("%s%s%4i%s",fHistTitle.Data()," Low Gain Pixel ",i," Runs: "));
451 h->SetXTitle(fHistXTitle.Data());
452 h->SetYTitle(fHistYTitle.Data());
453 }
454
455 if (IsAverageing())
456 {
457 if (fAverageLoGainAreas)
458 for (Int_t j=0; j<fAverageLoGainAreas->GetSize(); j++)
459 {
460 MHCalibrationPix &pix = GetAverageLoGainArea(j);
461 h = pix.GetHGausHist();
462 h->SetName (Form("%s%s%s%d","H",fHistName.Data(),"LoGainArea",j));
463 h->SetXTitle(fHistXTitle.Data());
464 h->SetYTitle(fHistYTitle.Data());
465 if (fGeom->InheritsFrom("MGeomCamMagic"))
466 h->SetTitle(Form("%s%s%s%s",fHistTitle.Data()," averaged on event-by-event basis ",
467 j==0 ? "Inner Pixels " : "Outer Pixels ","Low Gain Runs: "));
468 else
469 h->SetTitle(Form("%s%s%d%s",fHistTitle.Data(),
470 " averaged on event-by-event basis Low Gain Area Idx ",j," Runs: "));
471 }
472
473 if (fAverageLoGainSectors)
474 for (Int_t j=0; j<fAverageLoGainSectors->GetSize(); j++)
475 {
476 MHCalibrationPix &pix = GetAverageLoGainSector(j);
477 h = pix.GetHGausHist();
478 h->SetName (Form("%s%s%s%2i","H",fHistName.Data(),"LoGainSector",j));
479 h->SetTitle(Form("%s%s%2i%s",fHistTitle.Data(),
480 " averaged on event-by-event basis Low Gain Sector ",j," Runs: "));
481 h->SetXTitle(fHistXTitle.Data());
482 h->SetYTitle(fHistYTitle.Data());
483 }
484 }
485}
486
487// --------------------------------------------------------------------------
488//
489// Gets the pointers to:
490// - MGeomCam
491//
492// Calls SetupHists(const MParList *pList)
493//
494// Calls Delete-Function of:
495// - MHCalibrationCam::fHiGainArray, MHCalibrationCam::fLoGainArray
496// - MHCalibrationCam::fAverageHiGainAreas, MHCalibrationCam::fAverageLoGainAreas
497// - MHCalibrationCam::fAverageHiGainSectors, MHCalibrationCam::fAverageLoGainSectors
498//
499Bool_t MHCalibrationCam::SetupFill(const MParList *pList)
500{
501
502 fGeom = (MGeomCam*)pList->FindObject("MGeomCam");
503 if (!fGeom)
504 {
505 *fLog << err << GetDescriptor()
506 << ": MGeomCam not found... aborting." << endl;
507 return kFALSE;
508 }
509
510 fRunHeader = (MRawRunHeader*)pList->FindObject("MRawRunHeader");
511 if (!fRunHeader)
512 {
513 *fLog << warn << GetDescriptor()
514 << ": MRawRunHeader not found... will not store run numbers." << endl;
515 }
516
517 if (IsInterlaced())
518 {
519 fTrigPattern = (MTriggerPattern*)pList->FindObject("MTriggerPattern");
520 if (!fTrigPattern)
521 {
522 *fLog << err << "MTriggerPattern not found... abort." << endl;
523 return kFALSE;
524 }
525 if (IsInterlacedTypePed())
526 fMaskInterlaced |= MTriggerPattern::kPedestal;
527 if (IsInterlacedTypeCal())
528 fMaskInterlaced |= MTriggerPattern::kCalibration;
529 if (IsInterlacedTypePin())
530 fMaskInterlaced |= MTriggerPattern::kPinDiode;
531 }
532
533 return SetupHists(pList);
534}
535
536
537// --------------------------------------------------------------------------
538//
539// Searches MRawEvtHeader to find the correct pulser colour
540//
541// Gets or creates the pointers to:
542// - MBadPixelsIntensityCam
543// - MBadPixelsCam
544//
545// Searches pointer to:
546// - MArrivalTimeCam
547//
548// Initializes, if empty to MArrivalTimeCam::GetSize() for:
549// - MHCalibrationCam::fHiGainArray, MHCalibrationCam::fLoGainArray
550//
551// Initializes, if empty to MGeomCam::GetNumAreas() for:
552// - MHCalibrationCam::fAverageHiGainAreas, MHCalibrationCam::fAverageLoGainAreas
553//
554// Initializes, if empty to MGeomCam::GetNumSectors() for:
555// - MHCalibrationCam::fAverageHiGainSectors, MHCalibrationCam::fAverageLoGainSectors
556//
557// Initializes TArray's to MGeomCam::GetNumAreas and MGeomCam::GetNumSectors, respectively
558// Fills with number of valid pixels (if !MBadPixelsPix::IsBad()):
559// - MHCalibrationCam::fAverageAreaNum[area index]
560// - MHCalibrationCam::fAverageSectorNum[area index]
561//
562// Calls InitializeHists() for every entry in:
563// - MHCalibrationCam::fHiGainArray
564// - MHCalibrationCam::fAverageHiGainAreas
565// - MHCalibrationCam::fAverageHiGainSectors
566//
567// Sets Titles and Names for the Histograms
568// - MHCalibrationCam::fAverageHiGainAreas
569// - MHCalibrationCam::fAverageHiGainSectors
570//
571// Retrieves the run numbers from MRawRunHeader and stores them in fRunNumbers
572//
573Bool_t MHCalibrationCam::ReInit(MParList *pList)
574{
575
576 const Int_t npixels = fGeom->GetNumPixels();
577 const Int_t nsectors = fGeom->GetNumSectors();
578 const Int_t nareas = fGeom->GetNumAreas();
579
580 fIntensBad = (MBadPixelsIntensityCam*)pList->FindObject(AddSerialNumber("MBadPixelsIntensityCam"));
581 if (fIntensBad)
582 *fLog << inf << "Found MBadPixelsIntensityCam ... " << endl;
583 else
584 {
585 fBadPixels = (MBadPixelsCam*)pList->FindObject("MBadPixelsCam");
586 if (!fBadPixels)
587 {
588
589 fBadPixels = (MBadPixelsCam*)pList->FindCreateObj(AddSerialNumber("MBadPixelsCam"));
590 if (!fBadPixels)
591 {
592 *fLog << err << "Cannot find nor create MBadPixelsCam ... abort." << endl;
593 return kFALSE;
594 }
595 else
596 fBadPixels->InitSize(npixels);
597 }
598 }
599
600 if (IsAverageing())
601 {
602 //
603 // The function TArrayF::Set() already sets all entries to 0.
604 //
605 fAverageAreaNum. Set(nareas);
606 fAverageAreaSat. Set(nareas);
607 fAverageAreaSigma. Set(nareas);
608 fAverageAreaSigmaVar. Set(nareas);
609 fAverageAreaRelSigma. Set(nareas);
610 fAverageAreaRelSigmaVar.Set(nareas);
611 fAverageSectorNum. Set(nsectors);
612
613 for (Int_t aidx=0; aidx<nareas; aidx++)
614 fAverageAreaNum[aidx] = 0;
615
616 for (Int_t sector=0; sector<nsectors; sector++)
617 fAverageSectorNum[sector] = 0;
618
619 for (Int_t i=0; i<npixels; i++)
620 {
621
622 MBadPixelsPix &bad = fIntensBad ? (*fIntensBad)[i] : (*fBadPixels)[i];
623 if (bad.IsBad())
624 continue;
625
626 fAverageAreaNum [(*fGeom)[i].GetAidx() ]++;
627 fAverageSectorNum[(*fGeom)[i].GetSector()]++;
628 }
629 }
630
631 //
632 // Because ReInit has been called, a new run number is added
633 //
634 fRunNumbers.Set(fRunNumbers.GetSize()+1);
635
636 if (fRunHeader)
637 {
638 fRunNumbers[fRunNumbers.GetSize()-1] = fRunHeader->GetRunNumber();
639 if (IsLoGain())
640 SetLoGain(fRunHeader->GetNumSamplesLoGain());
641 }
642
643 if (!ReInitHists(pList))
644 return kFALSE;
645
646 ResetHistTitles();
647
648 if (!fRunHeader)
649 return kTRUE;
650
651 for (Int_t i=0; i<fHiGainArray->GetSize(); i++)
652 {
653 TH1F *h = (*this)[i].GetHGausHist();
654 h->SetTitle( Form("%s%i%s", h->GetTitle(),fRunNumbers[fRunNumbers.GetSize()-1]," "));
655 }
656
657 if (IsLoGain())
658 for (Int_t i=0; i<fLoGainArray->GetSize(); i++)
659 {
660 TH1F *h = (*this)(i).GetHGausHist();
661 h->SetTitle( Form("%s%i%s", h->GetTitle(),fRunNumbers[fRunNumbers.GetSize()-1]," "));
662 }
663
664 if (!IsAverageing())
665 return kTRUE;
666
667 for (Int_t j=0; j<nareas; j++)
668 {
669 TH1F *h = GetAverageHiGainArea(j).GetHGausHist();
670 h->SetTitle( Form("%s%i%s", h->GetTitle(),fRunNumbers[fRunNumbers.GetSize()-1]," "));
671 }
672
673 if (IsLoGain())
674 for (Int_t j=0; j<nareas; j++)
675 {
676 TH1F *h = GetAverageLoGainArea(j).GetHGausHist();
677 h->SetTitle( Form("%s%i%s", h->GetTitle(),fRunNumbers[fRunNumbers.GetSize()-1]," "));
678 }
679
680 for (Int_t j=0; j<nsectors; j++)
681 {
682 TH1F *h = GetAverageHiGainSector(j).GetHGausHist();
683 h->SetTitle( Form("%s%i%s", h->GetTitle(),fRunNumbers[fRunNumbers.GetSize()-1]," "));
684 }
685
686 if (IsLoGain())
687 for (Int_t j=0; j<nsectors; j++)
688 {
689 TH1F *h = GetAverageLoGainSector(j).GetHGausHist();
690 h->SetTitle( Form("%s%i%s", h->GetTitle(),fRunNumbers[fRunNumbers.GetSize()-1]," "));
691 }
692
693 return kTRUE;
694}
695
696//--------------------------------------------------------------------------------------
697//
698// Initializes the High Gain Arrays:
699//
700// - For every entry in the expanded arrays:
701// * Initialize an MHCalibrationPix
702// * Set Binning from fNbins, fFirst and fLast
703// * Set Histgram names and titles from fHistName and fHistTitle
704// * Set X-axis and Y-axis titles with fHistXTitle and fHistYTitle
705// * Call InitHists
706//
707void MHCalibrationCam::InitHiGainArrays(const Int_t npixels, const Int_t nareas, const Int_t nsectors)
708{
709
710 if (fHiGainArray->GetSize()==0)
711 {
712 for (Int_t i=0; i<npixels; i++)
713 {
714 fHiGainArray->AddAt(new MHCalibrationPix(Form("%sHiGainPix%04d",fHistName.Data(),i),
715 Form("%s High Gain Pixel %4d",fHistTitle.Data(),i)),i);
716
717 MHCalibrationPix &pix = (*this)[i];
718 pix.SetNbins(fNbins);
719 pix.SetFirst(fFirst);
720 pix.SetLast (fLast);
721 pix.SetProbLimit(fProbLimit);
722
723 MBadPixelsPix &bad = fIntensBad ? (*fIntensBad)[i] : (*fBadPixels)[i];
724 InitHists(pix,bad,i);
725 }
726 }
727
728 if (!IsAverageing())
729 return;
730
731 if (fAverageHiGainAreas->GetSize()==0)
732 {
733 for (Int_t j=0; j<nareas; j++)
734 {
735 fAverageHiGainAreas->AddAt(new MHCalibrationPix(Form("%sHiGainArea%d",fHistName.Data(),j),
736 Form("%s High Gain Area Idx %d",fHistTitle.Data(),j)),j);
737
738 MHCalibrationPix &pix = GetAverageHiGainArea(j);
739
740 pix.SetNbins(fNbins*(Int_t)TMath::Sqrt((Float_t)npixels/nareas));
741 pix.SetFirst(fFirst);
742 pix.SetLast (fLast);
743
744 if (fGeom && fGeom->InheritsFrom("MGeomCamMagic"))
745 {
746 pix.InitBins();
747 pix.SetEventFrequency(fPulserFrequency);
748 }
749 else
750 InitHists(pix,fIntensCam ? fIntensCam->GetAverageBadArea(j) : fCam->GetAverageBadArea(j),j);
751 }
752 }
753
754 if (fAverageHiGainSectors->GetSize()==0)
755 {
756 for (Int_t j=0; j<nsectors; j++)
757 {
758 fAverageHiGainSectors->AddAt(new MHCalibrationPix(Form("%sHiGainSector%02d",fHistName.Data(),j),
759 Form("%s High Gain Sector %02d",fHistTitle.Data(),j)),j);
760 MHCalibrationPix &pix = GetAverageHiGainSector(j);
761
762 pix.SetNbins(fNbins*(Int_t)TMath::Sqrt((Float_t)npixels/nsectors));
763 pix.SetFirst(fFirst);
764 pix.SetLast (fLast);
765
766 InitHists(pix,fIntensCam ? fIntensCam->GetAverageBadSector(j) : fCam->GetAverageBadSector(j),j);
767 }
768 }
769}
770
771//--------------------------------------------------------------------------------------
772//
773// Return, if IsLoGain() is kFALSE
774//
775// Initializes the Low Gain Arrays:
776//
777// - For every entry in the expanded arrays:
778// * Initialize an MHCalibrationPix
779// * Set Binning from fNbins, fFirst and fLast
780// * Set Histgram names and titles from fHistName and fHistTitle
781// * Set X-axis and Y-axis titles with fHistXTitle and fHistYTitle
782// * Call InitHists
783//
784void MHCalibrationCam::InitLoGainArrays(const Int_t npixels, const Int_t nareas, const Int_t nsectors)
785{
786
787 if (!IsLoGain())
788 return;
789
790 if (fLoGainArray->GetSize()==0)
791 {
792 for (Int_t i=0; i<npixels; i++)
793 {
794 fLoGainArray->AddAt(new MHCalibrationPix(Form("%sLoGainPix%04d",fHistName.Data(),i),
795 Form("%s Low Gain Pixel%04d",fHistTitle.Data(),i)),i);
796
797 MHCalibrationPix &pix = (*this)(i);
798 pix.SetNbins(fNbins);
799 pix.SetFirst(fFirst);
800 pix.SetLast (fLast);
801 pix.SetProbLimit(fProbLimit);
802
803 MBadPixelsPix &bad = fIntensBad ? (*fIntensBad)[i] : (*fBadPixels)[i];
804 InitHists(pix,bad,i);
805 }
806 }
807
808 if (!IsAverageing())
809 return;
810
811 if (fAverageLoGainAreas->GetSize()==0)
812 {
813 for (Int_t j=0; j<nareas; j++)
814 {
815 fAverageLoGainAreas->AddAt(new MHCalibrationPix(Form("%sLoGainArea%d",fHistName.Data(),j),
816 Form("%s Low Gain Area Idx %d",fHistTitle.Data(),j)),j);
817
818 MHCalibrationPix &pix = GetAverageLoGainArea(j);
819
820 pix.SetNbins(fNbins*(Int_t)TMath::Sqrt((Float_t)npixels/nareas));
821 pix.SetFirst(fFirst);
822 pix.SetLast (fLast);
823
824 if (fGeom && fGeom->InheritsFrom("MGeomCamMagic"))
825 {
826 pix.InitBins();
827 pix.SetEventFrequency(fPulserFrequency);
828 }
829 else
830 InitHists(pix,fIntensCam ? fIntensCam->GetAverageBadArea(j) : fCam->GetAverageBadArea(j),j);
831 }
832 }
833
834 if (fAverageLoGainSectors->GetSize()==0)
835 {
836 for (Int_t j=0; j<nsectors; j++)
837 {
838 fAverageLoGainSectors->AddAt(new MHCalibrationPix(Form("%sLoGainSector%02d",fHistName.Data(),j),
839 Form("%s Low Gain Sector %02d",fHistTitle.Data(),j)),j);
840 MHCalibrationPix &pix = GetAverageLoGainSector(j);
841
842 pix.SetNbins(fNbins*(Int_t)TMath::Sqrt((Float_t)npixels/nsectors));
843 pix.SetFirst(fFirst);
844 pix.SetLast (fLast);
845
846 InitHists(pix,fIntensCam ? fIntensCam->GetAverageBadSector(j) : fCam->GetAverageBadSector(j),j);
847 }
848 }
849}
850
851//--------------------------------------------------------------------------------
852//
853// Retrieves from MGeomCam:
854// - number of pixels
855// - number of pixel areas
856// - number of sectors
857//
858// Return kFALSE, if sizes of the TOrdCollections do not match npixels, nareas or nsectors
859//
860// Call FillHists()
861//
862Bool_t MHCalibrationCam::Fill(const MParContainer *par, const Stat_t w)
863{
864
865 if (IsInterlaced())
866 // check whether one of the bits in fMaskInterlaced is set in the trigger pattern
867 if (!(fTrigPattern->GetPrescaled() & fMaskInterlaced))
868 return kTRUE;
869
870 if (!IsSizeCheck())
871 return FillHists(par,w);
872
873 const Int_t npixels = fGeom->GetNumPixels();
874 const Int_t nareas = fGeom->GetNumAreas();
875 const Int_t nsectors = fGeom->GetNumSectors();
876
877 //
878 // Hi-Gain OrdCollections
879 //
880 if (fHiGainArray->GetSize() != npixels)
881 {
882 *fLog << err << "ERROR - Size mismatch in number of pixels... abort." << endl;
883 return kFALSE;
884 }
885
886 if (IsLoGain())
887 {
888 if (fLoGainArray->GetSize() != npixels)
889 {
890 *fLog << err << "ERROR - Size mismatch in number of pixels... abort." << endl;
891 return kFALSE;
892 }
893 }
894
895 if (!IsAverageing())
896 return FillHists(par,w);
897
898 if (fAverageHiGainAreas->GetSize() != nareas)
899 {
900 *fLog << err << "ERROR - Size mismatch in number of areas ... abort." << endl;
901 return kFALSE;
902 }
903
904 if (fAverageHiGainSectors->GetSize() != nsectors)
905 {
906 *fLog << err << "ERROR - Size mismatch in number of sectors ... abort." << endl;
907 return kFALSE;
908 }
909
910 if (IsLoGain())
911 {
912
913 if (fAverageLoGainAreas->GetSize() != nareas)
914 {
915 *fLog << err << "ERROR - Size mismatch in number of areas ... abort." << endl;
916 return kFALSE;
917 }
918
919 if (fAverageLoGainSectors->GetSize() != nsectors)
920 {
921 *fLog << err << "ERROR - Size mismatch in number of sectors ... abort." << endl;
922 return kFALSE;
923 }
924 }
925
926 return FillHists(par, w);
927}
928
929// --------------------------------------------------------------------------
930//
931// 0) Ask if fHiGainArray and fLoGainArray have been initialized,
932// otherwise return kFALSE.
933// 1) FinalizeHists()
934// 2) FinalizeBadPixels()
935// 3) CalcAverageSigma()
936//
937Bool_t MHCalibrationCam::Finalize()
938{
939
940 if (fHiGainArray->GetSize() == 0 && fLoGainArray->GetSize() == 0)
941 {
942 *fLog << err << GetDescriptor()
943 << ": ERROR: Both (HiGain and LoGain) histogram arrays have not been initialized... abort." << endl;
944 return kFALSE;
945 }
946
947 for (Int_t i=0; i<fAverageHiGainAreas->GetSize(); i++)
948 {
949 TH1F *h = GetAverageHiGainArea(i).GetHGausHist();
950 switch (fColor)
951 {
952 case MCalibrationCam::kNONE:
953 break;
954 case MCalibrationCam::kBLUE:
955 h->SetTitle( Form("%s%s", h->GetTitle(),"BLUE "));
956 break;
957 case MCalibrationCam::kGREEN:
958 h->SetTitle( Form("%s%s", h->GetTitle(),"GREEN "));
959 break;
960 case MCalibrationCam::kUV:
961 h->SetTitle( Form("%s%s", h->GetTitle(),"UV "));
962 break;
963 case MCalibrationCam::kCT1:
964 h->SetTitle( Form("%s%s", h->GetTitle(),"CT1-Pulser "));
965 break;
966 }
967 }
968
969 for (Int_t i=0; i<fAverageLoGainAreas->GetSize(); i++)
970 {
971 TH1F *h = GetAverageLoGainArea(i).GetHGausHist();
972 switch (fColor)
973 {
974 case MCalibrationCam::kNONE:
975 break;
976 case MCalibrationCam::kBLUE:
977 h->SetTitle( Form("%s%s", h->GetTitle(),"BLUE "));
978 break;
979 case MCalibrationCam::kGREEN:
980 h->SetTitle( Form("%s%s", h->GetTitle(),"GREEN "));
981 break;
982 case MCalibrationCam::kUV:
983 h->SetTitle( Form("%s%s", h->GetTitle(),"UV "));
984 break;
985 case MCalibrationCam::kCT1:
986 h->SetTitle( Form("%s%s", h->GetTitle(),"CT1-Pulser "));
987 break;
988 }
989 }
990
991 if (!FinalizeHists())
992 return kFALSE;
993
994
995 FinalizeBadPixels();
996 CalcAverageSigma();
997
998 return kTRUE;
999}
1000
1001// -------------------------------------------------------------
1002//
1003// If MBadPixelsPix::IsBad():
1004// - calls MHCalibrationPix::SetExcluded()
1005//
1006// Calls:
1007// - MHGausEvents::InitBins()
1008// - MHCalibrationPix::ChangeHistId(i)
1009// - MHCalibrationPix::SetEventFrequency(fPulserFrequency)
1010//
1011void MHCalibrationCam::InitHists(MHCalibrationPix &hist, MBadPixelsPix &bad, const Int_t i)
1012{
1013
1014 if (bad.IsUnsuitable(MBadPixelsPix::kUnsuitableRun))
1015 hist.SetExcluded();
1016
1017 hist.InitBins();
1018 hist.SetEventFrequency(fPulserFrequency);
1019}
1020
1021// --------------------------------------------------------------------------
1022//
1023// Calls FitHiGainHists for every entry in:
1024// - fHiGainArray
1025// - fAverageHiGainAreas
1026// - fAverageHiGainSectors
1027//
1028void MHCalibrationCam::FitHiGainArrays(MCalibrationCam &calcam, MBadPixelsCam &badcam,
1029 MBadPixelsPix::UncalibratedType_t fittyp,
1030 MBadPixelsPix::UncalibratedType_t osctyp)
1031{
1032
1033 for (Int_t i=0; i<fHiGainArray->GetSize(); i++)
1034 {
1035
1036 MHCalibrationPix &hist = (*this)[i];
1037
1038 if (hist.IsExcluded())
1039 continue;
1040
1041 MCalibrationPix &pix = calcam[i];
1042 MBadPixelsPix &bad = badcam[i];
1043
1044 FitHiGainHists(hist,pix,bad,fittyp,osctyp);
1045 }
1046
1047 if (!IsAverageing())
1048 return;
1049
1050 for (Int_t j=0; j<fAverageHiGainAreas->GetSize(); j++)
1051 {
1052
1053 MHCalibrationPix &hist = GetAverageHiGainArea(j);
1054 MCalibrationPix &pix = calcam.GetAverageArea(j);
1055 MBadPixelsPix &bad = calcam.GetAverageBadArea(j);
1056
1057 FitHiGainHists(hist,pix,bad,fittyp,osctyp);
1058 }
1059
1060 for (Int_t j=0; j<fAverageHiGainSectors->GetSize(); j++)
1061 {
1062 MHCalibrationPix &hist = GetAverageHiGainSector(j);
1063 MCalibrationPix &pix = calcam.GetAverageSector(j);
1064 MBadPixelsPix &bad = calcam.GetAverageBadSector(j);
1065
1066 FitHiGainHists(hist,pix,bad,fittyp,osctyp);
1067 }
1068}
1069
1070// --------------------------------------------------------------------------
1071//
1072// Calls FitLoGainHists for every entry in:
1073// - fLoGainArray
1074// - fAverageLoGainAreas
1075// - fAverageLoGainSectors
1076//
1077void MHCalibrationCam::FitLoGainArrays(MCalibrationCam &calcam, MBadPixelsCam &badcam,
1078 MBadPixelsPix::UncalibratedType_t fittyp,
1079 MBadPixelsPix::UncalibratedType_t osctyp)
1080{
1081
1082 if (!IsLoGain())
1083 return;
1084
1085 for (Int_t i=0; i<fLoGainArray->GetSize(); i++)
1086 {
1087
1088 MHCalibrationPix &hist = (*this)(i);
1089
1090 if (hist.IsExcluded())
1091 continue;
1092
1093 MCalibrationPix &pix = calcam[i];
1094 MBadPixelsPix &bad = badcam[i];
1095
1096 FitLoGainHists(hist,pix,bad,fittyp,osctyp);
1097
1098 }
1099
1100 if (!IsAverageing())
1101 return;
1102
1103 for (Int_t j=0; j<fAverageLoGainAreas->GetSize(); j++)
1104 {
1105
1106 MHCalibrationPix &hist = GetAverageLoGainArea(j);
1107 MCalibrationPix &pix = calcam.GetAverageArea(j);
1108 MBadPixelsPix &bad = calcam.GetAverageBadArea(j);
1109
1110 FitLoGainHists(hist,pix,bad,fittyp,osctyp);
1111 }
1112
1113 for (Int_t j=0; j<fAverageLoGainSectors->GetSize(); j++)
1114 {
1115
1116 MHCalibrationPix &hist = GetAverageLoGainSector(j);
1117 MCalibrationPix &pix = calcam.GetAverageSector(j);
1118 MBadPixelsPix &bad = calcam.GetAverageBadSector(j);
1119
1120 FitLoGainHists(hist,pix,bad,fittyp,osctyp);
1121 }
1122}
1123
1124//------------------------------------------------------------
1125//
1126// For all averaged areas, the fitted sigma is multiplied with the square root of
1127// the number involved pixels
1128//
1129void MHCalibrationCam::CalcAverageSigma()
1130{
1131
1132 if (!fCam)
1133 return;
1134
1135 if (!IsAverageing())
1136 return;
1137
1138 for (UInt_t j=0; j<fGeom->GetNumAreas(); j++)
1139 {
1140
1141 MCalibrationPix &pix = fCam->GetAverageArea(j);
1142
1143 const Float_t numsqr = TMath::Sqrt((Float_t)fAverageAreaNum[j]);
1144 fAverageAreaSigma[j] = pix.GetSigma () * numsqr;
1145 fAverageAreaSigmaVar[j] = pix.GetSigmaErr () * pix.GetSigmaErr() * numsqr;
1146
1147 pix.SetSigma (fAverageAreaSigma[j]);
1148 pix.SetSigmaVar(fAverageAreaSigmaVar[j]);
1149
1150 fAverageAreaRelSigma [j] = fAverageAreaSigma[j] / pix.GetMean();
1151 fAverageAreaRelSigmaVar[j] = fAverageAreaSigmaVar[j] / (fAverageAreaSigma[j]*fAverageAreaSigma[j]);
1152 fAverageAreaRelSigmaVar[j] += pix.GetMeanRelVar();
1153 fAverageAreaRelSigmaVar[j] *= fAverageAreaRelSigma[j];
1154 }
1155}
1156
1157// ---------------------------------------------------------------------------
1158//
1159// Returns if the histogram is empty and sets the following flag:
1160// - MBadPixelsPix::SetUnsuitable(MBadPixelsPix::kUnsuitableRun)
1161//
1162// Fits the histograms with a Gaussian, in case of failure
1163// calls MHCalibrationPix::RepeatFit(), in case of repeated failure
1164// calls MHCalibrationPix::BypassFit() and sets the following flags:
1165// - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::UncalibratedType_t fittyp )
1166// - MBadPixelsPix::SetUnsuitable( MBadPixelsPix::kUnreliableRun )
1167//
1168// Creates the fourier spectrum and tests MHGausEvents::IsFourierSpectrumOK().
1169// In case no, sets the following flags:
1170// - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::UncalibratedType_t osctyp )
1171// - MBadPixelsPix::SetUnsuitable( MBadPixelsPix::kUnreliableRun )
1172//
1173// Retrieves the results and store them in MCalibrationPix
1174//
1175void MHCalibrationCam::FitHiGainHists(MHCalibrationPix &hist,
1176 MCalibrationPix &pix,
1177 MBadPixelsPix &bad,
1178 MBadPixelsPix::UncalibratedType_t fittyp,
1179 MBadPixelsPix::UncalibratedType_t osctyp)
1180{
1181
1182 if (hist.IsEmpty() || hist.IsOnlyOverflow() || hist.IsOnlyUnderflow())
1183 return;
1184
1185 //
1186 // 2) Fit the Hi Gain histograms with a Gaussian
1187 //
1188 if (!hist.FitGaus())
1189 //
1190 // 3) In case of failure set the bit Fitted to false and take histogram means and RMS
1191 //
1192 if (!hist.RepeatFit())
1193 {
1194 hist.BypassFit();
1195 bad.SetUncalibrated( fittyp );
1196 }
1197
1198 //
1199 // 4) Check for oscillations
1200 //
1201 if (IsOscillations())
1202 {
1203 hist.CreateFourierSpectrum();
1204
1205 if (!hist.IsFourierSpectrumOK())
1206 bad.SetUncalibrated( osctyp );
1207 }
1208
1209 //
1210 // 5) Retrieve the results and store them in this class
1211 //
1212 pix.SetHiGainMean ( hist.GetMean() );
1213 pix.SetHiGainMeanVar ( hist.GetMeanErr() * hist.GetMeanErr() );
1214 pix.SetHiGainRms ( hist.GetHistRms() );
1215 pix.SetHiGainSigma ( hist.GetSigma() );
1216 pix.SetHiGainSigmaVar ( hist.GetSigmaErr()* hist.GetSigmaErr() );
1217 pix.SetHiGainProb ( hist.GetProb() );
1218 pix.SetHiGainNumBlackout( hist.GetBlackout() );
1219 pix.SetHiGainNumPickup ( hist.GetPickup() );
1220
1221 if (IsDebug())
1222 {
1223 *fLog << dbginf << GetDescriptor() << ": ID " << GetName()
1224 << " HiGainSaturation: " << pix.IsHiGainSaturation()
1225 << " HiGainMean: " << hist.GetMean ()
1226 << " HiGainMeanErr: " << hist.GetMeanErr ()
1227 << " HiGainMeanSigma: " << hist.GetSigma ()
1228 << " HiGainMeanSigmaErr: " << hist.GetSigmaErr()
1229 << " HiGainMeanProb: " << hist.GetProb ()
1230 << " HiGainNumBlackout: " << hist.GetBlackout()
1231 << " HiGainNumPickup : " << hist.GetPickup ()
1232 << endl;
1233 }
1234
1235}
1236
1237
1238// ---------------------------------------------------------------------------
1239//
1240// Returns if the histogram is empty and sets the following flag:
1241// - MBadPixelsPix::SetUnsuitable(MBadPixelsPix::kUnsuitableRun)
1242//
1243// Fits the histograms with a Gaussian, in case of failure
1244// calls MHCalibrationPix::RepeatFit(), in case of repeated failure
1245// calls MHCalibrationPix::BypassFit() and sets the following flags:
1246// - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::UncalibratedType_t fittyp )
1247// - MBadPixelsPix::SetUnsuitable( MBadPixelsPix::kUnreliableRun )
1248//
1249// Creates the fourier spectrum and tests MHGausEvents::IsFourierSpectrumOK().
1250// In case no, sets the following flags:
1251// - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::UncalibratedType_t osctyp )
1252// - MBadPixelsPix::SetUnsuitable( MBadPixelsPix::kUnreliableRun )
1253//
1254// Retrieves the results and store them in MCalibrationPix
1255//
1256void MHCalibrationCam::FitLoGainHists(MHCalibrationPix &hist,
1257 MCalibrationPix &pix,
1258 MBadPixelsPix &bad,
1259 MBadPixelsPix::UncalibratedType_t fittyp,
1260 MBadPixelsPix::UncalibratedType_t osctyp)
1261{
1262
1263 if (hist.IsEmpty() || hist.IsOnlyOverflow() || hist.IsOnlyUnderflow())
1264 return;
1265
1266
1267 //
1268 // 2) Fit the Hi Gain histograms with a Gaussian
1269 //
1270 if (!hist.FitGaus())
1271 //
1272 // 3) In case of failure set the bit Fitted to false and take histogram means and RMS
1273 //
1274 if (!hist.RepeatFit())
1275 {
1276 hist.BypassFit();
1277 if (pix.IsHiGainSaturation())
1278 bad.SetUncalibrated( fittyp );
1279 }
1280
1281 //
1282 // 4) Check for oscillations
1283 //
1284 if (IsOscillations())
1285 {
1286 hist.CreateFourierSpectrum();
1287
1288 if (!hist.IsFourierSpectrumOK())
1289 bad.SetUncalibrated( osctyp );
1290 }
1291
1292 //
1293 // 5) Retrieve the results and store them in this class
1294 //
1295 pix.SetLoGainMean ( hist.GetMean() );
1296 pix.SetLoGainMeanVar ( hist.GetMeanErr() * hist.GetMeanErr() );
1297 pix.SetLoGainRms ( hist.GetHistRms() );
1298 pix.SetLoGainSigma ( hist.GetSigma() );
1299 pix.SetLoGainSigmaVar ( hist.GetSigmaErr() * hist.GetSigmaErr() );
1300 pix.SetLoGainProb ( hist.GetProb() );
1301 pix.SetLoGainNumBlackout( hist.GetBlackout() );
1302 pix.SetLoGainNumPickup ( hist.GetPickup() );
1303
1304 if (IsDebug())
1305 {
1306 *fLog << dbginf << GetDescriptor() << "ID: " << hist.GetName()
1307 << " HiGainSaturation: " << pix.IsHiGainSaturation()
1308 << " LoGainMean: " << hist.GetMean ()
1309 << " LoGainMeanErr: " << hist.GetMeanErr ()
1310 << " LoGainMeanSigma: " << hist.GetSigma ()
1311 << " LoGainMeanSigmaErr: " << hist.GetSigmaErr()
1312 << " LoGainMeanProb: " << hist.GetProb ()
1313 << " LoGainNumBlackout: " << hist.GetBlackout()
1314 << " LoGainNumPickup : " << hist.GetPickup ()
1315 << endl;
1316 }
1317
1318}
1319
1320
1321
1322// -----------------------------------------------------------------------------
1323//
1324// Default draw:
1325//
1326// Displays the averaged areas, both High Gain and Low Gain
1327//
1328// Calls the Draw of the fAverageHiGainAreas and fAverageLoGainAreas objects with options
1329//
1330void MHCalibrationCam::Draw(const Option_t *opt)
1331{
1332
1333 if (!IsAverageing())
1334 return;
1335
1336 const Int_t nareas = fAverageHiGainAreas->GetSize();
1337 if (nareas == 0)
1338 return;
1339
1340 TVirtualPad *pad = gPad ? gPad : MH::MakeDefCanvas(this);
1341 pad->SetBorderMode(0);
1342
1343 pad->Divide(IsLoGain() ? 2 : 1,nareas);
1344
1345 for (Int_t i=0; i<nareas;i++)
1346 {
1347
1348 pad->cd(IsLoGain() ? 2*(i+1)-1 : i+1);
1349 GetAverageHiGainArea(i).Draw(opt);
1350
1351 if (!fAverageAreaSat[i])
1352 DrawAverageSigma(fAverageAreaSat[i], i,
1353 fAverageAreaSigma[i], fAverageAreaSigmaVar[i],
1354 fAverageAreaRelSigma[i], fAverageAreaRelSigmaVar[i]);
1355
1356 if (IsLoGain())
1357 {
1358 pad->cd(2*(i+1));
1359 GetAverageLoGainArea(i).Draw(opt);
1360 }
1361
1362 if (fAverageAreaSat[i])
1363 DrawAverageSigma(fAverageAreaSat[i], i,
1364 fAverageAreaSigma[i], fAverageAreaSigmaVar[i],
1365 fAverageAreaRelSigma[i], fAverageAreaRelSigmaVar[i]);
1366 }
1367
1368}
1369
1370// -----------------------------------------------------------------------------
1371//
1372// Default draw:
1373//
1374// Displays a TPaveText with the re-normalized sigmas of the average area
1375//
1376void MHCalibrationCam::DrawAverageSigma(Bool_t sat, Bool_t inner,
1377 Float_t sigma, Float_t sigmavar,
1378 Float_t relsigma, Float_t relsigmavar) const
1379{
1380
1381 if (sigma != 0 && sigmavar >= 0 && relsigmavar >= 0.)
1382 {
1383
1384 TPad *newpad = new TPad("newpad","transparent",0,0,1,1);
1385 newpad->SetFillStyle(4000);
1386 newpad->Draw();
1387 newpad->cd();
1388
1389 TPaveText *text = new TPaveText(sat? 0.1 : 0.35,0.7,sat ? 0.4 : 0.7,1.0);
1390 text->SetTextSize(0.07);
1391 const TString line1 = Form("%s%s%s",inner ? "Outer" : "Inner",
1392 " Pixels ", sat ? "Low Gain" : "High Gain");
1393 TText *txt1 = text->AddText(line1.Data());
1394 const TString line2 = Form("#sigma per pix: %2.2f #pm %2.2f",sigma,TMath::Sqrt(sigmavar));
1395 TText *txt2 = text->AddText(line2.Data());
1396 const TString line3 = Form("Rel. #sigma per pix: %2.2f #pm %2.2f",relsigma,TMath::Sqrt(relsigmavar));
1397 TText *txt3 = text->AddText(line3.Data());
1398 text->Draw("");
1399
1400 text->SetBit(kCanDelete);
1401 txt1->SetBit(kCanDelete);
1402 txt2->SetBit(kCanDelete);
1403 txt3->SetBit(kCanDelete);
1404 newpad->SetBit(kCanDelete);
1405 }
1406}
1407
1408Int_t MHCalibrationCam::ReadEnv(const TEnv &env, TString prefix, Bool_t print)
1409{
1410
1411 Bool_t rc = kFALSE;
1412 if (IsEnvDefined(env, prefix, "Debug", print))
1413 {
1414 SetDebug(GetEnvValue(env, prefix, "Debug", IsDebug()));
1415 rc = kTRUE;
1416 }
1417 if (IsEnvDefined(env, prefix, "LoGain", print))
1418 {
1419 SetDebug(GetEnvValue(env, prefix, "LoGain", IsLoGain()));
1420 rc = kTRUE;
1421 }
1422 if (IsEnvDefined(env, prefix, "Oscillations", print))
1423 {
1424 SetOscillations(GetEnvValue(env, prefix, "Oscillations", IsOscillations()));
1425 rc = kTRUE;
1426 }
1427 if (IsEnvDefined(env, prefix, "SizeCheck", print))
1428 {
1429 SetSizeCheck(GetEnvValue(env, prefix, "SizeCheck", IsSizeCheck()));
1430 rc = kTRUE;
1431 }
1432 if (IsEnvDefined(env, prefix, "Averageing", print))
1433 {
1434 SetAverageing(GetEnvValue(env, prefix, "Averageing", IsAverageing()));
1435 rc = kTRUE;
1436 }
1437
1438 if (IsEnvDefined(env, prefix, "Nbins", print))
1439 {
1440 SetNbins(GetEnvValue(env, prefix, "Nbins", fNbins));
1441 rc = kTRUE;
1442 }
1443 if (IsEnvDefined(env, prefix, "First", print))
1444 {
1445 SetFirst(GetEnvValue(env, prefix, "First", fFirst));
1446 rc = kTRUE;
1447 }
1448 if (IsEnvDefined(env, prefix, "Last", print))
1449 {
1450 SetLast(GetEnvValue(env, prefix, "Last", fLast));
1451 rc = kTRUE;
1452 }
1453
1454 return rc;
1455}
1456
Note: See TracBrowser for help on using the repository browser.