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

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