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

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