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

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