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

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