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

Last change on this file since 6595 was 6528, checked in by tbretz, 21 years ago
*** empty log message ***
File size: 46.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
84ClassImp(MHCalibrationCam);
85
86using namespace std;
87
88const Int_t MHCalibrationCam::fgPulserFrequency = 500;
89const Float_t MHCalibrationCam::fgProbLimit = 0.0001;
90const Float_t MHCalibrationCam::fgOverflowLimit = 0.005;
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// - fOverflowLimit to fgOverflowLimit
111//
112// - SetAveregeing (kTRUE);
113// - SetDebug (kFALSE);
114// - SetLoGain (kTRUE);
115//- SetOscillations(kTRUE);
116//- SetSizeCheck (kTRUE);
117//- SetInterlaced (kFALSE);
118//
119MHCalibrationCam::MHCalibrationCam(const char *name, const char *title)
120 : fHistName(gsHistName),fHistTitle(gsHistTitle),
121 fHistXTitle(gsHistXTitle),fHistYTitle(gsHistYTitle),
122 fColor(MCalibrationCam::kNONE), fIntensBad(NULL),
123 fBadPixels(NULL), fIntensCam(NULL), fCam(NULL), fGeom(NULL),
124 fRunHeader(NULL)
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 SetOverflowLimit();
148
149 SetAverageing (kTRUE);
150 SetDebug (kFALSE);
151 SetLoGain (kTRUE);
152 SetOscillations(kTRUE);
153 SetSizeCheck (kTRUE);
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 return SetupHists(pList);
518}
519
520
521// --------------------------------------------------------------------------
522//
523// Searches MRawEvtHeader to find the correct pulser colour
524//
525// Gets or creates the pointers to:
526// - MBadPixelsIntensityCam
527// - MBadPixelsCam
528//
529// Searches pointer to:
530// - MArrivalTimeCam
531//
532// Initializes, if empty to MArrivalTimeCam::GetSize() for:
533// - MHCalibrationCam::fHiGainArray, MHCalibrationCam::fLoGainArray
534//
535// Initializes, if empty to MGeomCam::GetNumAreas() for:
536// - MHCalibrationCam::fAverageHiGainAreas, MHCalibrationCam::fAverageLoGainAreas
537//
538// Initializes, if empty to MGeomCam::GetNumSectors() for:
539// - MHCalibrationCam::fAverageHiGainSectors, MHCalibrationCam::fAverageLoGainSectors
540//
541// Initializes TArray's to MGeomCam::GetNumAreas and MGeomCam::GetNumSectors, respectively
542// Fills with number of valid pixels (if !MBadPixelsPix::IsBad()):
543// - MHCalibrationCam::fAverageAreaNum[area index]
544// - MHCalibrationCam::fAverageSectorNum[area index]
545//
546// Calls InitializeHists() for every entry in:
547// - MHCalibrationCam::fHiGainArray
548// - MHCalibrationCam::fAverageHiGainAreas
549// - MHCalibrationCam::fAverageHiGainSectors
550//
551// Sets Titles and Names for the Histograms
552// - MHCalibrationCam::fAverageHiGainAreas
553// - MHCalibrationCam::fAverageHiGainSectors
554//
555// Retrieves the run numbers from MRawRunHeader and stores them in fRunNumbers
556//
557Bool_t MHCalibrationCam::ReInit(MParList *pList)
558{
559
560 const Int_t npixels = fGeom->GetNumPixels();
561 const Int_t nsectors = fGeom->GetNumSectors();
562 const Int_t nareas = fGeom->GetNumAreas();
563
564 fIntensBad = (MBadPixelsIntensityCam*)pList->FindObject(AddSerialNumber("MBadPixelsIntensityCam"));
565 if (fIntensBad)
566 *fLog << inf << "Found MBadPixelsIntensityCam ... " << endl;
567 else
568 {
569 fBadPixels = (MBadPixelsCam*)pList->FindObject("MBadPixelsCam");
570 if (!fBadPixels)
571 {
572
573 fBadPixels = (MBadPixelsCam*)pList->FindCreateObj(AddSerialNumber("MBadPixelsCam"));
574 if (!fBadPixels)
575 return kFALSE;
576
577 fBadPixels->InitSize(npixels);
578 }
579 }
580
581 if (IsAverageing())
582 {
583 //
584 // The function TArrayF::Set() already sets all entries to 0.
585 //
586 fAverageAreaNum. Set(nareas);
587 fAverageAreaSat. Set(nareas);
588 fAverageAreaSigma. Set(nareas);
589 fAverageAreaSigmaVar. Set(nareas);
590 fAverageAreaRelSigma. Set(nareas);
591 fAverageAreaRelSigmaVar.Set(nareas);
592 fAverageSectorNum. Set(nsectors);
593
594 for (Int_t aidx=0; aidx<nareas; aidx++)
595 fAverageAreaNum[aidx] = 0;
596
597 for (Int_t sector=0; sector<nsectors; sector++)
598 fAverageSectorNum[sector] = 0;
599
600 for (Int_t i=0; i<npixels; i++)
601 {
602
603 MBadPixelsPix &bad = fIntensBad ? (*fIntensBad)[i] : (*fBadPixels)[i];
604 if (bad.IsBad())
605 continue;
606
607 fAverageAreaNum [(*fGeom)[i].GetAidx() ]++;
608 fAverageSectorNum[(*fGeom)[i].GetSector()]++;
609 }
610 }
611
612 //
613 // Because ReInit has been called, a new run number is added
614 //
615 fRunNumbers.Set(fRunNumbers.GetSize()+1);
616
617 if (fRunHeader)
618 {
619 fRunNumbers[fRunNumbers.GetSize()-1] = fRunHeader->GetRunNumber();
620 if (IsLoGain())
621 SetLoGain(fRunHeader->GetNumSamplesLoGain());
622 }
623
624 if (!ReInitHists(pList))
625 return kFALSE;
626
627 ResetHistTitles();
628
629 if (!fRunHeader)
630 return kTRUE;
631
632 for (Int_t i=0; i<fHiGainArray->GetSize(); i++)
633 {
634 TH1F *h = (*this)[i].GetHGausHist();
635 h->SetTitle( Form("%s%i%s", h->GetTitle(),fRunNumbers[fRunNumbers.GetSize()-1]," "));
636 }
637
638 if (IsLoGain())
639 for (Int_t i=0; i<fLoGainArray->GetSize(); i++)
640 {
641 TH1F *h = (*this)(i).GetHGausHist();
642 h->SetTitle( Form("%s%i%s", h->GetTitle(),fRunNumbers[fRunNumbers.GetSize()-1]," "));
643 }
644
645 if (!IsAverageing())
646 return kTRUE;
647
648 for (Int_t j=0; j<nareas; j++)
649 {
650 TH1F *h = GetAverageHiGainArea(j).GetHGausHist();
651 h->SetTitle( Form("%s%i%s", h->GetTitle(),fRunNumbers[fRunNumbers.GetSize()-1]," "));
652 }
653
654 if (IsLoGain())
655 for (Int_t j=0; j<nareas; j++)
656 {
657 TH1F *h = GetAverageLoGainArea(j).GetHGausHist();
658 h->SetTitle( Form("%s%i%s", h->GetTitle(),fRunNumbers[fRunNumbers.GetSize()-1]," "));
659 }
660
661 for (Int_t j=0; j<nsectors; j++)
662 {
663 TH1F *h = GetAverageHiGainSector(j).GetHGausHist();
664 h->SetTitle( Form("%s%i%s", h->GetTitle(),fRunNumbers[fRunNumbers.GetSize()-1]," "));
665 }
666
667 if (IsLoGain())
668 for (Int_t j=0; j<nsectors; j++)
669 {
670 TH1F *h = GetAverageLoGainSector(j).GetHGausHist();
671 h->SetTitle( Form("%s%i%s", h->GetTitle(),fRunNumbers[fRunNumbers.GetSize()-1]," "));
672 }
673
674 return kTRUE;
675}
676
677//--------------------------------------------------------------------------------------
678//
679// Initializes the High Gain Arrays:
680//
681// - For every entry in the expanded arrays:
682// * Initialize an MHCalibrationPix
683// * Set Binning from fNbins, fFirst and fLast
684// * Set Histgram names and titles from fHistName and fHistTitle
685// * Set X-axis and Y-axis titles with fHistXTitle and fHistYTitle
686// * Call InitHists
687//
688void MHCalibrationCam::InitHiGainArrays(const Int_t npixels, const Int_t nareas, const Int_t nsectors)
689{
690
691 if (fHiGainArray->GetSize()==0)
692 {
693 for (Int_t i=0; i<npixels; i++)
694 {
695 fHiGainArray->AddAt(new MHCalibrationPix(Form("%sHiGainPix%04d",fHistName.Data(),i),
696 Form("%s High Gain Pixel %4d",fHistTitle.Data(),i)),i);
697
698 MHCalibrationPix &pix = (*this)[i];
699 pix.SetNbins(fNbins);
700 pix.SetFirst(fFirst);
701 pix.SetLast (fLast);
702 pix.SetProbLimit(fProbLimit);
703
704 MBadPixelsPix &bad = fIntensBad ? (*fIntensBad)[i] : (*fBadPixels)[i];
705 InitHists(pix,bad,i);
706 }
707 }
708
709 if (!IsAverageing())
710 return;
711
712 if (fAverageHiGainAreas->GetSize()==0)
713 {
714 for (Int_t j=0; j<nareas; j++)
715 {
716 fAverageHiGainAreas->AddAt(new MHCalibrationPix(Form("%sHiGainArea%d",fHistName.Data(),j),
717 Form("%s High Gain Area Idx %d",fHistTitle.Data(),j)),j);
718
719 MHCalibrationPix &pix = GetAverageHiGainArea(j);
720
721 pix.SetNbins(fNbins*(Int_t)TMath::Sqrt((Float_t)npixels/nareas));
722 pix.SetFirst(fFirst);
723 pix.SetLast (fLast);
724
725 if (fGeom && fGeom->InheritsFrom("MGeomCamMagic"))
726 {
727 pix.InitBins();
728 pix.SetEventFrequency(fPulserFrequency);
729 }
730 else
731 InitHists(pix,fIntensCam ? fIntensCam->GetAverageBadArea(j) : fCam->GetAverageBadArea(j),j);
732 }
733 }
734
735 if (fAverageHiGainSectors->GetSize()==0)
736 {
737 for (Int_t j=0; j<nsectors; j++)
738 {
739 fAverageHiGainSectors->AddAt(new MHCalibrationPix(Form("%sHiGainSector%02d",fHistName.Data(),j),
740 Form("%s High Gain Sector %02d",fHistTitle.Data(),j)),j);
741 MHCalibrationPix &pix = GetAverageHiGainSector(j);
742
743 pix.SetNbins(fNbins*(Int_t)TMath::Sqrt((Float_t)npixels/nsectors));
744 pix.SetFirst(fFirst);
745 pix.SetLast (fLast);
746
747 InitHists(pix,fIntensCam ? fIntensCam->GetAverageBadSector(j) : fCam->GetAverageBadSector(j),j);
748 }
749 }
750}
751
752//--------------------------------------------------------------------------------------
753//
754// Return, if IsLoGain() is kFALSE
755//
756// Initializes the Low Gain Arrays:
757//
758// - For every entry in the expanded arrays:
759// * Initialize an MHCalibrationPix
760// * Set Binning from fNbins, fFirst and fLast
761// * Set Histgram names and titles from fHistName and fHistTitle
762// * Set X-axis and Y-axis titles with fHistXTitle and fHistYTitle
763// * Call InitHists
764//
765void MHCalibrationCam::InitLoGainArrays(const Int_t npixels, const Int_t nareas, const Int_t nsectors)
766{
767
768 if (!IsLoGain())
769 return;
770
771 if (fLoGainArray->GetSize()==0)
772 {
773 for (Int_t i=0; i<npixels; i++)
774 {
775 fLoGainArray->AddAt(new MHCalibrationPix(Form("%sLoGainPix%04d",fHistName.Data(),i),
776 Form("%s Low Gain Pixel%04d",fHistTitle.Data(),i)),i);
777
778 MHCalibrationPix &pix = (*this)(i);
779 pix.SetNbins(fNbins);
780 pix.SetFirst(fFirst);
781 pix.SetLast (fLast);
782 pix.SetProbLimit(fProbLimit);
783
784 MBadPixelsPix &bad = fIntensBad ? (*fIntensBad)[i] : (*fBadPixels)[i];
785 InitHists(pix,bad,i);
786 }
787 }
788
789 if (!IsAverageing())
790 return;
791
792 if (fAverageLoGainAreas->GetSize()==0)
793 {
794 for (Int_t j=0; j<nareas; j++)
795 {
796 fAverageLoGainAreas->AddAt(new MHCalibrationPix(Form("%sLoGainArea%d",fHistName.Data(),j),
797 Form("%s Low Gain Area Idx %d",fHistTitle.Data(),j)),j);
798
799 MHCalibrationPix &pix = GetAverageLoGainArea(j);
800
801 pix.SetNbins(fNbins*(Int_t)TMath::Sqrt((Float_t)npixels/nareas));
802 pix.SetFirst(fFirst);
803 pix.SetLast (fLast);
804
805 if (fGeom && fGeom->InheritsFrom("MGeomCamMagic"))
806 {
807 pix.InitBins();
808 pix.SetEventFrequency(fPulserFrequency);
809 }
810 else
811 InitHists(pix,fIntensCam ? fIntensCam->GetAverageBadArea(j) : fCam->GetAverageBadArea(j),j);
812 }
813 }
814
815 if (fAverageLoGainSectors->GetSize()==0)
816 {
817 for (Int_t j=0; j<nsectors; j++)
818 {
819 fAverageLoGainSectors->AddAt(new MHCalibrationPix(Form("%sLoGainSector%02d",fHistName.Data(),j),
820 Form("%s Low Gain Sector %02d",fHistTitle.Data(),j)),j);
821 MHCalibrationPix &pix = GetAverageLoGainSector(j);
822
823 pix.SetNbins(fNbins*(Int_t)TMath::Sqrt((Float_t)npixels/nsectors));
824 pix.SetFirst(fFirst);
825 pix.SetLast (fLast);
826
827 InitHists(pix,fIntensCam ? fIntensCam->GetAverageBadSector(j) : fCam->GetAverageBadSector(j),j);
828 }
829 }
830}
831
832//--------------------------------------------------------------------------------
833//
834// Retrieves from MGeomCam:
835// - number of pixels
836// - number of pixel areas
837// - number of sectors
838//
839// Return kFALSE, if sizes of the TOrdCollections do not match npixels, nareas or nsectors
840//
841// Call FillHists()
842//
843Bool_t MHCalibrationCam::Fill(const MParContainer *par, const Stat_t w)
844{
845
846 if (!IsSizeCheck())
847 return FillHists(par,w);
848
849 const Int_t npixels = fGeom->GetNumPixels();
850 const Int_t nareas = fGeom->GetNumAreas();
851 const Int_t nsectors = fGeom->GetNumSectors();
852
853 //
854 // Hi-Gain OrdCollections
855 //
856 if (fHiGainArray->GetSize() != npixels)
857 {
858 *fLog << err << "ERROR - Size mismatch in number of pixels... abort." << endl;
859 return kFALSE;
860 }
861
862 if (IsLoGain())
863 {
864 if (fLoGainArray->GetSize() != npixels)
865 {
866 *fLog << err << "ERROR - Size mismatch in number of pixels... abort." << endl;
867 return kFALSE;
868 }
869 }
870
871 if (!IsAverageing())
872 return FillHists(par,w);
873
874 if (fAverageHiGainAreas->GetSize() != nareas)
875 {
876 *fLog << err << "ERROR - Size mismatch in number of areas ... abort." << endl;
877 return kFALSE;
878 }
879
880 if (fAverageHiGainSectors->GetSize() != nsectors)
881 {
882 *fLog << err << "ERROR - Size mismatch in number of sectors ... abort." << endl;
883 return kFALSE;
884 }
885
886 if (IsLoGain())
887 {
888
889 if (fAverageLoGainAreas->GetSize() != nareas)
890 {
891 *fLog << err << "ERROR - Size mismatch in number of areas ... abort." << endl;
892 return kFALSE;
893 }
894
895 if (fAverageLoGainSectors->GetSize() != nsectors)
896 {
897 *fLog << err << "ERROR - Size mismatch in number of sectors ... abort." << endl;
898 return kFALSE;
899 }
900 }
901
902 return FillHists(par, w);
903}
904
905// --------------------------------------------------------------------------
906//
907// 0) Ask if fHiGainArray and fLoGainArray have been initialized,
908// otherwise return kFALSE.
909// 1) FinalizeHists()
910// 2) FinalizeBadPixels()
911// 3) CalcAverageSigma()
912//
913Bool_t MHCalibrationCam::Finalize()
914{
915
916 if (GetNumExecutions() < 2)
917 return kTRUE;
918
919 if (fHiGainArray->GetSize() == 0 && fLoGainArray->GetSize() == 0)
920 {
921 *fLog << err << GetDescriptor()
922 << ": ERROR: Both (HiGain and LoGain) histogram arrays have not been initialized... abort." << endl;
923 return kFALSE;
924 }
925
926 for (Int_t i=0; i<fAverageHiGainAreas->GetSize(); i++)
927 {
928 TH1F *h = GetAverageHiGainArea(i).GetHGausHist();
929 switch (fColor)
930 {
931 case MCalibrationCam::kNONE:
932 break;
933 case MCalibrationCam::kBLUE:
934 h->SetTitle( Form("%s%s", h->GetTitle(),"BLUE "));
935 break;
936 case MCalibrationCam::kGREEN:
937 h->SetTitle( Form("%s%s", h->GetTitle(),"GREEN "));
938 break;
939 case MCalibrationCam::kUV:
940 h->SetTitle( Form("%s%s", h->GetTitle(),"UV "));
941 break;
942 case MCalibrationCam::kCT1:
943 h->SetTitle( Form("%s%s", h->GetTitle(),"CT1-Pulser "));
944 break;
945 }
946 }
947
948 for (Int_t i=0; i<fAverageLoGainAreas->GetSize(); i++)
949 {
950 TH1F *h = GetAverageLoGainArea(i).GetHGausHist();
951 switch (fColor)
952 {
953 case MCalibrationCam::kNONE:
954 break;
955 case MCalibrationCam::kBLUE:
956 h->SetTitle( Form("%s%s", h->GetTitle(),"BLUE "));
957 break;
958 case MCalibrationCam::kGREEN:
959 h->SetTitle( Form("%s%s", h->GetTitle(),"GREEN "));
960 break;
961 case MCalibrationCam::kUV:
962 h->SetTitle( Form("%s%s", h->GetTitle(),"UV "));
963 break;
964 case MCalibrationCam::kCT1:
965 h->SetTitle( Form("%s%s", h->GetTitle(),"CT1-Pulser "));
966 break;
967 }
968 }
969
970 if (!FinalizeHists())
971 return kFALSE;
972
973
974 FinalizeBadPixels();
975 CalcAverageSigma();
976
977 return kTRUE;
978}
979
980// -------------------------------------------------------------
981//
982// If MBadPixelsPix::IsUnsuitable(MBadPixelsPix::kUnsuitableRun):
983// - calls MHCalibrationPix::SetExcluded()
984//
985// Calls:
986// - MHGausEvents::InitBins()
987// - MHCalibrationPix::ChangeHistId(i)
988// - MHCalibrationPix::SetEventFrequency(fPulserFrequency)
989//
990void MHCalibrationCam::InitHists(MHCalibrationPix &hist, MBadPixelsPix &bad, const Int_t i)
991{
992
993 if (bad.IsUnsuitable(MBadPixelsPix::kUnsuitableRun))
994 hist.SetExcluded();
995
996 hist.InitBins();
997 hist.SetEventFrequency(fPulserFrequency);
998}
999
1000// -------------------------------------------------------------
1001//
1002// - Searches for the CalibrationIntensity*Cam corresponding to 'name'.
1003// - In case, it does not exist in the parameter list, it searches
1004// for the corresponding MCalibration*Cam.
1005// - Initializes the MCalibration*Cam, if not yet done.
1006//
1007Bool_t MHCalibrationCam::InitCams( MParList *plist, const TString name )
1008{
1009
1010 TString intensname = "MCalibrationIntensity";
1011 intensname += name;
1012 intensname += "Cam";
1013
1014 TString ordname = "MCalibration";
1015 ordname += name;
1016 ordname += "Cam";
1017
1018 fIntensCam = (MCalibrationIntensityCam*)plist->FindObject(AddSerialNumber(intensname));
1019 if (fIntensCam)
1020 *fLog << inf << "Found " << intensname << "... " << endl;
1021 else
1022 {
1023 fCam = (MCalibrationCam*)plist->FindObject(AddSerialNumber(ordname));
1024 if (!fCam)
1025 {
1026 fCam = (MCalibrationCam*)plist->FindCreateObj(AddSerialNumber(ordname));
1027 if (!fCam)
1028 return kFALSE;
1029
1030 fCam->Init(*fGeom);
1031 }
1032 }
1033 return kTRUE;
1034}
1035
1036// --------------------------------------------------------------------------
1037//
1038// Calls FitHiGainHists for every entry in:
1039// - fHiGainArray
1040// - fAverageHiGainAreas
1041// - fAverageHiGainSectors
1042//
1043void MHCalibrationCam::FitHiGainArrays(MCalibrationCam &calcam, MBadPixelsCam &badcam,
1044 MBadPixelsPix::UncalibratedType_t fittyp,
1045 MBadPixelsPix::UncalibratedType_t osctyp)
1046{
1047
1048 for (Int_t i=0; i<fHiGainArray->GetSize(); i++)
1049 {
1050
1051 MHCalibrationPix &hist = (*this)[i];
1052
1053 if (hist.IsExcluded())
1054 continue;
1055
1056 MCalibrationPix &pix = calcam[i];
1057 MBadPixelsPix &bad = badcam[i];
1058
1059 FitHiGainHists(hist,pix,bad,fittyp,osctyp);
1060 }
1061
1062 if (!IsAverageing())
1063 return;
1064
1065 for (Int_t j=0; j<fAverageHiGainAreas->GetSize(); j++)
1066 {
1067
1068 MHCalibrationPix &hist = GetAverageHiGainArea(j);
1069 MCalibrationPix &pix = calcam.GetAverageArea(j);
1070 MBadPixelsPix &bad = calcam.GetAverageBadArea(j);
1071
1072 FitHiGainHists(hist,pix,bad,fittyp,osctyp);
1073 }
1074
1075 for (Int_t j=0; j<fAverageHiGainSectors->GetSize(); j++)
1076 {
1077 MHCalibrationPix &hist = GetAverageHiGainSector(j);
1078 MCalibrationPix &pix = calcam.GetAverageSector(j);
1079 MBadPixelsPix &bad = calcam.GetAverageBadSector(j);
1080
1081 FitHiGainHists(hist,pix,bad,fittyp,osctyp);
1082 }
1083}
1084
1085// --------------------------------------------------------------------------
1086//
1087// Calls FitLoGainHists for every entry in:
1088// - fLoGainArray
1089// - fAverageLoGainAreas
1090// - fAverageLoGainSectors
1091//
1092void MHCalibrationCam::FitLoGainArrays(MCalibrationCam &calcam, MBadPixelsCam &badcam,
1093 MBadPixelsPix::UncalibratedType_t fittyp,
1094 MBadPixelsPix::UncalibratedType_t osctyp)
1095{
1096
1097 if (!IsLoGain())
1098 return;
1099
1100 for (Int_t i=0; i<fLoGainArray->GetSize(); i++)
1101 {
1102
1103 MHCalibrationPix &hist = (*this)(i);
1104
1105 if (hist.IsExcluded())
1106 continue;
1107
1108 MCalibrationPix &pix = calcam[i];
1109 MBadPixelsPix &bad = badcam[i];
1110
1111 FitLoGainHists(hist,pix,bad,fittyp,osctyp);
1112
1113 }
1114
1115 if (!IsAverageing())
1116 return;
1117
1118 for (Int_t j=0; j<fAverageLoGainAreas->GetSize(); j++)
1119 {
1120
1121 MHCalibrationPix &hist = GetAverageLoGainArea(j);
1122 MCalibrationPix &pix = calcam.GetAverageArea(j);
1123 MBadPixelsPix &bad = calcam.GetAverageBadArea(j);
1124
1125 FitLoGainHists(hist,pix,bad,fittyp,osctyp);
1126 }
1127
1128 for (Int_t j=0; j<fAverageLoGainSectors->GetSize(); j++)
1129 {
1130
1131 MHCalibrationPix &hist = GetAverageLoGainSector(j);
1132 MCalibrationPix &pix = calcam.GetAverageSector(j);
1133 MBadPixelsPix &bad = calcam.GetAverageBadSector(j);
1134
1135 FitLoGainHists(hist,pix,bad,fittyp,osctyp);
1136 }
1137}
1138
1139//------------------------------------------------------------
1140//
1141// For all averaged areas, the fitted sigma is multiplied with the square root of
1142// the number involved pixels
1143//
1144void MHCalibrationCam::CalcAverageSigma()
1145{
1146
1147 if (!fCam)
1148 return;
1149
1150 if (!IsAverageing())
1151 return;
1152
1153 for (UInt_t j=0; j<fGeom->GetNumAreas(); j++)
1154 {
1155
1156 MCalibrationPix &pix = fCam->GetAverageArea(j);
1157
1158 const Float_t numsqr = TMath::Sqrt((Float_t)fAverageAreaNum[j]);
1159 fAverageAreaSigma[j] = pix.GetSigma () * numsqr;
1160 fAverageAreaSigmaVar[j] = pix.GetSigmaErr () * pix.GetSigmaErr() * numsqr;
1161
1162 pix.SetSigma (fAverageAreaSigma[j]);
1163 pix.SetSigmaVar(fAverageAreaSigmaVar[j]);
1164
1165 fAverageAreaRelSigma [j] = fAverageAreaSigma[j] / pix.GetMean();
1166 fAverageAreaRelSigmaVar[j] = fAverageAreaSigmaVar[j] / (fAverageAreaSigma[j]*fAverageAreaSigma[j]);
1167 fAverageAreaRelSigmaVar[j] += pix.GetMeanRelVar();
1168 fAverageAreaRelSigmaVar[j] *= fAverageAreaRelSigma[j];
1169 }
1170}
1171
1172// ---------------------------------------------------------------------------
1173//
1174// Returns if the histogram is empty and sets the following flag:
1175// - MBadPixelsPix::SetUnsuitable(MBadPixelsPix::kUnsuitableRun)
1176//
1177// Fits the histograms with a Gaussian, in case of failure
1178// calls MHCalibrationPix::RepeatFit(), in case of repeated failure
1179// calls MHCalibrationPix::BypassFit() and sets the following flags:
1180// - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::UncalibratedType_t fittyp )
1181// - MBadPixelsPix::SetUnsuitable( MBadPixelsPix::kUnreliableRun )
1182//
1183// Creates the fourier spectrum and tests MHGausEvents::IsFourierSpectrumOK().
1184// In case no, sets the following flags:
1185// - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::UncalibratedType_t osctyp )
1186// - MBadPixelsPix::SetUnsuitable( MBadPixelsPix::kUnreliableRun )
1187//
1188// Retrieves the results and store them in MCalibrationPix
1189//
1190void MHCalibrationCam::FitHiGainHists(MHCalibrationPix &hist,
1191 MCalibrationPix &pix,
1192 MBadPixelsPix &bad,
1193 MBadPixelsPix::UncalibratedType_t fittyp,
1194 MBadPixelsPix::UncalibratedType_t osctyp)
1195{
1196
1197 if (hist.IsEmpty() || hist.IsOnlyOverflow() || hist.IsOnlyUnderflow())
1198 {
1199 *fLog << warn << GetDescriptor()
1200 << ": Only overflow or underflow in hi-gain pixel: " << pix.GetPixId() << endl;
1201 return;
1202 }
1203 //
1204 // 2) Fit the Hi Gain histograms with a Gaussian
1205 //
1206 if (!hist.FitGaus())
1207 //
1208 // 3) In case of failure set the bit Fitted to false and take histogram means and RMS
1209 //
1210 if (!hist.RepeatFit())
1211 {
1212 hist.BypassFit();
1213 bad.SetUncalibrated( fittyp );
1214 }
1215
1216 //
1217 // 4) Check for oscillations
1218 //
1219 if (IsOscillations())
1220 {
1221 hist.CreateFourierSpectrum();
1222
1223 if (!hist.IsFourierSpectrumOK())
1224 bad.SetUncalibrated( osctyp );
1225 }
1226
1227 //
1228 // 5) Retrieve the results and store them in this class
1229 //
1230 pix.SetHiGainMean ( hist.GetMean() );
1231 pix.SetHiGainMeanVar ( hist.GetMeanErr() * hist.GetMeanErr() );
1232 pix.SetHiGainRms ( hist.GetHistRms() );
1233 pix.SetHiGainSigma ( hist.GetSigma() );
1234 pix.SetHiGainSigmaVar ( hist.GetSigmaErr()* hist.GetSigmaErr() );
1235 pix.SetHiGainProb ( hist.GetProb() );
1236 pix.SetHiGainNumBlackout( hist.GetBlackout() );
1237 pix.SetHiGainNumPickup ( hist.GetPickup() );
1238
1239 if (IsDebug())
1240 {
1241 *fLog << dbginf << GetDescriptor() << ": ID " << GetName()
1242 << " HiGainSaturation: " << pix.IsHiGainSaturation()
1243 << " HiGainMean: " << hist.GetMean ()
1244 << " HiGainMeanErr: " << hist.GetMeanErr ()
1245 << " HiGainMeanSigma: " << hist.GetSigma ()
1246 << " HiGainMeanSigmaErr: " << hist.GetSigmaErr()
1247 << " HiGainMeanProb: " << hist.GetProb ()
1248 << " HiGainNumBlackout: " << hist.GetBlackout()
1249 << " HiGainNumPickup : " << hist.GetPickup ()
1250 << endl;
1251 }
1252
1253}
1254
1255
1256// ---------------------------------------------------------------------------
1257//
1258// Returns if the histogram is empty and sets the following flag:
1259// - MBadPixelsPix::SetUnsuitable(MBadPixelsPix::kUnsuitableRun)
1260//
1261// Fits the histograms with a Gaussian, in case of failure
1262// calls MHCalibrationPix::RepeatFit(), in case of repeated failure
1263// calls MHCalibrationPix::BypassFit() and sets the following flags:
1264// - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::UncalibratedType_t fittyp )
1265// - MBadPixelsPix::SetUnsuitable( MBadPixelsPix::kUnreliableRun )
1266//
1267// Creates the fourier spectrum and tests MHGausEvents::IsFourierSpectrumOK().
1268// In case no, sets the following flags:
1269// - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::UncalibratedType_t osctyp )
1270// - MBadPixelsPix::SetUnsuitable( MBadPixelsPix::kUnreliableRun )
1271//
1272// Retrieves the results and store them in MCalibrationPix
1273//
1274void MHCalibrationCam::FitLoGainHists(MHCalibrationPix &hist,
1275 MCalibrationPix &pix,
1276 MBadPixelsPix &bad,
1277 MBadPixelsPix::UncalibratedType_t fittyp,
1278 MBadPixelsPix::UncalibratedType_t osctyp)
1279{
1280
1281 if (hist.IsEmpty() || hist.IsOnlyOverflow() || hist.IsOnlyUnderflow())
1282 {
1283 *fLog << warn << GetDescriptor()
1284 << ": Only overflow or underflow in lo-gain pixel: " << pix.GetPixId() << endl;
1285 return;
1286 }
1287 //
1288 // 2) Fit the Hi Gain histograms with a Gaussian
1289 //
1290 if (!hist.FitGaus())
1291 //
1292 // 3) In case of failure set the bit Fitted to false and take histogram means and RMS
1293 //
1294 if (!hist.RepeatFit())
1295 {
1296 hist.BypassFit();
1297 if (pix.IsHiGainSaturation())
1298 bad.SetUncalibrated( fittyp );
1299 }
1300
1301 //
1302 // 4) Check for oscillations
1303 //
1304 if (IsOscillations())
1305 {
1306 hist.CreateFourierSpectrum();
1307
1308 if (!hist.IsFourierSpectrumOK())
1309 bad.SetUncalibrated( osctyp );
1310 }
1311
1312 //
1313 // 5) Retrieve the results and store them in this class
1314 //
1315 pix.SetLoGainMean ( hist.GetMean() );
1316 pix.SetLoGainMeanVar ( hist.GetMeanErr() * hist.GetMeanErr() );
1317 pix.SetLoGainRms ( hist.GetHistRms() );
1318 pix.SetLoGainSigma ( hist.GetSigma() );
1319 pix.SetLoGainSigmaVar ( hist.GetSigmaErr() * hist.GetSigmaErr() );
1320 pix.SetLoGainProb ( hist.GetProb() );
1321 pix.SetLoGainNumBlackout( hist.GetBlackout() );
1322 pix.SetLoGainNumPickup ( hist.GetPickup() );
1323
1324 if (IsDebug())
1325 {
1326 *fLog << dbginf << GetDescriptor() << "ID: " << hist.GetName()
1327 << " HiGainSaturation: " << pix.IsHiGainSaturation()
1328 << " LoGainMean: " << hist.GetMean ()
1329 << " LoGainMeanErr: " << hist.GetMeanErr ()
1330 << " LoGainMeanSigma: " << hist.GetSigma ()
1331 << " LoGainMeanSigmaErr: " << hist.GetSigmaErr()
1332 << " LoGainMeanProb: " << hist.GetProb ()
1333 << " LoGainNumBlackout: " << hist.GetBlackout()
1334 << " LoGainNumPickup : " << hist.GetPickup ()
1335 << endl;
1336 }
1337
1338}
1339
1340
1341
1342// -----------------------------------------------------------------------------
1343//
1344// Default draw:
1345//
1346// Displays the averaged areas, both High Gain and Low Gain
1347//
1348// Calls the Draw of the fAverageHiGainAreas and fAverageLoGainAreas objects with options
1349//
1350void MHCalibrationCam::Draw(const Option_t *opt)
1351{
1352
1353 if (!IsAverageing())
1354 return;
1355
1356 const Int_t nareas = fAverageHiGainAreas->GetSize();
1357 if (nareas == 0)
1358 return;
1359
1360 TVirtualPad *pad = gPad ? gPad : MH::MakeDefCanvas(this);
1361 pad->SetBorderMode(0);
1362
1363 pad->Divide(IsLoGain() ? 2 : 1,nareas);
1364
1365 for (Int_t i=0; i<nareas;i++)
1366 {
1367
1368 pad->cd(IsLoGain() ? 2*(i+1)-1 : i+1);
1369 GetAverageHiGainArea(i).Draw(opt);
1370
1371 if (!fAverageAreaSat[i])
1372 DrawAverageSigma(fAverageAreaSat[i], i,
1373 fAverageAreaSigma[i], fAverageAreaSigmaVar[i],
1374 fAverageAreaRelSigma[i], fAverageAreaRelSigmaVar[i]);
1375
1376 if (IsLoGain())
1377 {
1378 pad->cd(2*(i+1));
1379 GetAverageLoGainArea(i).Draw(opt);
1380 }
1381
1382 if (fAverageAreaSat[i])
1383 DrawAverageSigma(fAverageAreaSat[i], i,
1384 fAverageAreaSigma[i], fAverageAreaSigmaVar[i],
1385 fAverageAreaRelSigma[i], fAverageAreaRelSigmaVar[i]);
1386 }
1387
1388}
1389
1390// -----------------------------------------------------------------------------
1391//
1392// Default draw:
1393//
1394// Displays a TPaveText with the re-normalized sigmas of the average area
1395//
1396void MHCalibrationCam::DrawAverageSigma(Bool_t sat, Bool_t inner,
1397 Float_t sigma, Float_t sigmavar,
1398 Float_t relsigma, Float_t relsigmavar) const
1399{
1400
1401 if (sigma != 0 && sigmavar >= 0 && relsigmavar >= 0.)
1402 {
1403
1404 TPad *newpad = new TPad("newpad","transparent",0,0,1,1);
1405 newpad->SetFillStyle(4000);
1406 newpad->Draw();
1407 newpad->cd();
1408
1409 TPaveText *text = new TPaveText(sat? 0.1 : 0.35,0.7,sat ? 0.4 : 0.7,1.0);
1410 text->SetTextSize(0.07);
1411 const TString line1 = Form("%s%s%s",inner ? "Outer" : "Inner",
1412 " Pixels ", sat ? "Low Gain" : "High Gain");
1413 TText *txt1 = text->AddText(line1.Data());
1414 const TString line2 = Form("#sigma per pix: %2.2f #pm %2.2f",sigma,TMath::Sqrt(sigmavar));
1415 TText *txt2 = text->AddText(line2.Data());
1416 const TString line3 = Form("Rel. #sigma per pix: %2.2f #pm %2.2f",relsigma,TMath::Sqrt(relsigmavar));
1417 TText *txt3 = text->AddText(line3.Data());
1418 text->Draw("");
1419
1420 text->SetBit(kCanDelete);
1421 txt1->SetBit(kCanDelete);
1422 txt2->SetBit(kCanDelete);
1423 txt3->SetBit(kCanDelete);
1424 newpad->SetBit(kCanDelete);
1425 }
1426}
1427
1428Int_t MHCalibrationCam::ReadEnv(const TEnv &env, TString prefix, Bool_t print)
1429{
1430
1431 Bool_t rc = kFALSE;
1432 if (IsEnvDefined(env, prefix, "Debug", print))
1433 {
1434 SetDebug(GetEnvValue(env, prefix, "Debug", IsDebug()));
1435 rc = kTRUE;
1436 }
1437 if (IsEnvDefined(env, prefix, "LoGain", print))
1438 {
1439 SetDebug(GetEnvValue(env, prefix, "LoGain", IsLoGain()));
1440 rc = kTRUE;
1441 }
1442 if (IsEnvDefined(env, prefix, "Oscillations", print))
1443 {
1444 SetOscillations(GetEnvValue(env, prefix, "Oscillations", IsOscillations()));
1445 rc = kTRUE;
1446 }
1447 if (IsEnvDefined(env, prefix, "SizeCheck", print))
1448 {
1449 SetSizeCheck(GetEnvValue(env, prefix, "SizeCheck", IsSizeCheck()));
1450 rc = kTRUE;
1451 }
1452 if (IsEnvDefined(env, prefix, "Averageing", print))
1453 {
1454 SetAverageing(GetEnvValue(env, prefix, "Averageing", IsAverageing()));
1455 rc = kTRUE;
1456 }
1457
1458 if (IsEnvDefined(env, prefix, "Nbins", print))
1459 {
1460 SetNbins(GetEnvValue(env, prefix, "Nbins", fNbins));
1461 rc = kTRUE;
1462 }
1463 if (IsEnvDefined(env, prefix, "First", print))
1464 {
1465 SetFirst(GetEnvValue(env, prefix, "First", fFirst));
1466 rc = kTRUE;
1467 }
1468
1469 if (IsEnvDefined(env, prefix, "Last", print))
1470 {
1471 SetLast(GetEnvValue(env, prefix, "Last", fLast));
1472 rc = kTRUE;
1473 }
1474
1475 if (IsEnvDefined(env, prefix, "ProbLimit", print))
1476 {
1477 SetProbLimit(GetEnvValue(env, prefix, "ProbLimit", fProbLimit));
1478 rc = kTRUE;
1479 }
1480
1481 if (IsEnvDefined(env, prefix, "OverflowLimit", print))
1482 {
1483 SetOverflowLimit(GetEnvValue(env, prefix, "OverflowLimit", fOverflowLimit));
1484 rc = kTRUE;
1485 }
1486
1487 if (IsEnvDefined(env, prefix, "PulserFrequency", print))
1488 {
1489 SetPulserFrequency(GetEnvValue(env, prefix, "PulserFrequency", fPulserFrequency));
1490 rc = kTRUE;
1491 }
1492
1493 return rc;
1494}
1495
Note: See TracBrowser for help on using the repository browser.