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

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