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

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