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

Last change on this file since 4964 was 4960, checked in by gaug, 20 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->fAverageHiGainSectors->Expand(nsehi);
345 // cam->fAverageLoGainSectors->Expand(nselo);
346
347 /*
348 for (int i=0; i<nhi; i++)
349 {
350 delete (*cam->fHiGainArray)[i];
351 (*cam->fHiGainArray)[i] = (*fHiGainArray)[i]->Clone();
352 }
353 for (int i=0; i<nlo; i++)
354 {
355 delete (*cam->fLoGainArray)[i];
356 (*cam->fLoGainArray)[i] = (*fLoGainArray)[i]->Clone();
357 }
358 */
359
360 for (int i=0; i<navhi; i++)
361 (*cam->fAverageHiGainAreas)[i] = (*fAverageHiGainAreas)[i]->Clone();
362
363 // for (int i=0; i<nsehi; i++)
364 // (*cam->fAverageHiGainSectors)[i] = (*fAverageHiGainSectors)[i]->Clone();
365
366 if (IsLoGain())
367 {
368
369 const Int_t navlo = fAverageLoGainAreas->GetEntries();
370
371 cam->fAverageLoGainAreas->Expand(navlo);
372 // cam->fAverageLoGainSectors->Expand(nselo);
373 // const Int_t nselo = fAverageLoGainSectors->GetEntries();
374
375 for (int i=0; i<navlo; i++)
376 (*cam->fAverageLoGainAreas)[i] = (*fAverageLoGainAreas)[i]->Clone();
377
378 // for (int i=0; i<nselo; i++)
379 // (*cam->fAverageLoGainSectors)[i] = (*fAverageLoGainSectors)[i]->Clone();
380 }
381
382 cam->fAverageAreaNum = fAverageAreaNum;
383 cam->fAverageAreaSat = fAverageAreaSat;
384 cam->fAverageAreaSigma = fAverageAreaSigma;
385 cam->fAverageAreaSigmaVar = fAverageAreaSigmaVar;
386 cam->fAverageAreaRelSigma = fAverageAreaRelSigma;
387 cam->fAverageAreaRelSigmaVar = fAverageAreaRelSigmaVar;
388 cam->fAverageSectorNum = fAverageSectorNum;
389 cam->fRunNumbers = fRunNumbers;
390
391 cam->fPulserFrequency = fPulserFrequency;
392
393 return cam;
394}
395
396// --------------------------------------------------------------------------
397//
398// Gets the pointers to:
399// - MGeomCam
400//
401// Calls SetupHists(const MParList *pList)
402//
403// Calls Delete-Function of:
404// - MHCalibrationCam::fHiGainArray, MHCalibrationCam::fLoGainArray
405// - MHCalibrationCam::fAverageHiGainAreas, MHCalibrationCam::fAverageLoGainAreas
406// - MHCalibrationCam::fAverageHiGainSectors, MHCalibrationCam::fAverageLoGainSectors
407//
408Bool_t MHCalibrationCam::SetupFill(const MParList *pList)
409{
410
411 fGeom = (MGeomCam*)pList->FindObject("MGeomCam");
412 if (!fGeom)
413 {
414 *fLog << err << GetDescriptor()
415 << ": MGeomCam not found... aborting." << endl;
416 return kFALSE;
417 }
418
419 fRunHeader = (MRawRunHeader*)pList->FindObject("MRawRunHeader");
420 if (!fRunHeader)
421 {
422 *fLog << warn << GetDescriptor()
423 << ": MRawRunHeader not found... will not store run numbers." << endl;
424 }
425
426 fHiGainArray->Delete();
427 fLoGainArray->Delete();
428
429 fAverageHiGainAreas->Delete();
430 fAverageLoGainAreas->Delete();
431
432 fAverageHiGainSectors->Delete();
433 fAverageLoGainSectors->Delete();
434
435 return SetupHists(pList);
436}
437
438
439// --------------------------------------------------------------------------
440//
441// Gets or creates the pointers to:
442// - MBadPixelsCam
443//
444// Searches pointer to:
445// - MArrivalTimeCam
446//
447// Initializes, if empty to MArrivalTimeCam::GetSize() for:
448// - MHCalibrationCam::fHiGainArray, MHCalibrationCam::fLoGainArray
449//
450// Initializes, if empty to MGeomCam::GetNumAreas() for:
451// - MHCalibrationCam::fAverageHiGainAreas, MHCalibrationCam::fAverageLoGainAreas
452//
453// Initializes, if empty to MGeomCam::GetNumSectors() for:
454// - MHCalibrationCam::fAverageHiGainSectors, MHCalibrationCam::fAverageLoGainSectors
455//
456// Initializes TArray's to MGeomCam::GetNumAreas and MGeomCam::GetNumSectors, respectively
457// Fills with number of valid pixels (if !MBadPixelsPix::IsBad()):
458// - MHCalibrationCam::fAverageAreaNum[area index]
459// - MHCalibrationCam::fAverageSectorNum[area index]
460//
461// Calls InitializeHists() for every entry in:
462// - MHCalibrationCam::fHiGainArray
463// - MHCalibrationCam::fAverageHiGainAreas
464// - MHCalibrationCam::fAverageHiGainSectors
465//
466// Sets Titles and Names for the Histograms
467// - MHCalibrationCam::fAverageHiGainAreas
468// - MHCalibrationCam::fAverageHiGainSectors
469//
470// Retrieves the run numbers from MRawRunHeader and stores them in fRunNumbers
471//
472Bool_t MHCalibrationCam::ReInit(MParList *pList)
473{
474
475 const Int_t npixels = fGeom->GetNumPixels();
476 const Int_t nsectors = fGeom->GetNumSectors();
477 const Int_t nareas = fGeom->GetNumAreas();
478
479 fBadPixels = (MBadPixelsCam*)pList->FindObject("MBadPixelsCam");
480 if (!fBadPixels)
481 {
482
483 fBadPixels = (MBadPixelsCam*)pList->FindCreateObj(AddSerialNumber("MBadPixelsCam"));
484 if (!fBadPixels)
485 {
486 *fLog << err << "Cannot find nor create MBadPixelsCam ... abort." << endl;
487 return kFALSE;
488 }
489 else
490 fBadPixels->InitSize(npixels);
491 }
492
493 //
494 // The function TArrayF::Set() already sets all entries to 0.
495 //
496 fAverageAreaNum. Set(nareas);
497 fAverageAreaSat. Set(nareas);
498 fAverageAreaSigma. Set(nareas);
499 fAverageAreaSigmaVar. Set(nareas);
500 fAverageAreaRelSigma. Set(nareas);
501 fAverageAreaRelSigmaVar.Set(nareas);
502 fAverageSectorNum. Set(nsectors);
503 fRunNumbers. Set(fRunNumbers.GetSize()+1);
504
505 for (Int_t aidx=0; aidx<nareas; aidx++)
506 fAverageAreaNum[aidx] = 0;
507
508 for (Int_t sector=0; sector<nsectors; sector++)
509 fAverageSectorNum[sector] = 0;
510
511 for (Int_t i=0; i<npixels; i++)
512 {
513
514 if ((*fBadPixels)[i].IsBad())
515 continue;
516
517 fAverageAreaNum [(*fGeom)[i].GetAidx() ]++;
518 fAverageSectorNum[(*fGeom)[i].GetSector()]++;
519 }
520
521 if (fRunHeader)
522 {
523 fRunNumbers[fRunNumbers.GetSize()-1] = fRunHeader->GetRunNumber();
524 if (IsLoGain())
525 SetLoGain(fRunHeader->GetNumSamplesLoGain());
526 }
527
528 if (!ReInitHists(pList))
529 return kFALSE;
530
531 if (!fRunHeader)
532 return kTRUE;
533
534 for (Int_t i=0; i<fHiGainArray->GetEntries(); i++)
535 {
536 TH1F *h = (*this)[i].GetHGausHist();
537 h->SetTitle( Form("%s%i%s", h->GetTitle(),fRunNumbers[fRunNumbers.GetSize()-1]," "));
538 }
539
540 if (IsLoGain())
541 for (Int_t i=0; i<fLoGainArray->GetEntries(); i++)
542 {
543 TH1F *h = (*this)(i).GetHGausHist();
544 h->SetTitle( Form("%s%i%s", h->GetTitle(),fRunNumbers[fRunNumbers.GetSize()-1]," "));
545 }
546
547 for (Int_t j=0; j<nareas; j++)
548 {
549 TH1F *h = GetAverageHiGainArea(j).GetHGausHist();
550 h->SetTitle( Form("%s%i%s", h->GetTitle(),fRunNumbers[fRunNumbers.GetSize()-1]," "));
551 }
552
553 if (IsLoGain())
554 for (Int_t j=0; j<nareas; j++)
555 {
556 TH1F *h = GetAverageLoGainArea(j).GetHGausHist();
557 h->SetTitle( Form("%s%i%s", h->GetTitle(),fRunNumbers[fRunNumbers.GetSize()-1]," "));
558 }
559
560 for (Int_t j=0; j<nsectors; j++)
561 {
562 TH1F *h = GetAverageHiGainSector(j).GetHGausHist();
563 h->SetTitle( Form("%s%i%s", h->GetTitle(),fRunNumbers[fRunNumbers.GetSize()-1]," "));
564 }
565
566 if (IsLoGain())
567 for (Int_t j=0; j<nsectors; j++)
568 {
569 TH1F *h = GetAverageLoGainSector(j).GetHGausHist();
570 h->SetTitle( Form("%s%i%s", h->GetTitle(),fRunNumbers[fRunNumbers.GetSize()-1]," "));
571 }
572
573 return kTRUE;
574}
575
576//--------------------------------------------------------------------------------------
577//
578// Initializes the High Gain Arrays:
579//
580// - Expand fHiGainArrays to npixels
581// - Expand fAverageHiGainAreas to nareas
582// - Expand fAverageHiGainSectors to nsectors
583//
584// - For every entry in the expanded arrays:
585// * Initialize an MHCalibrationPix
586// * Set Binning from fNbins, fFirst and fLast
587// * Set Histgram names and titles from fHistName and fHistTitle
588// * Set X-axis and Y-axis titles with fHistXTitle and fHistYTitle
589// * Call InitHists
590//
591void MHCalibrationCam::InitHiGainArrays(const Int_t npixels, const Int_t nareas, const Int_t nsectors)
592{
593
594 if (fHiGainArray->GetEntries()==0)
595 {
596 fHiGainArray->Expand(npixels);
597 for (Int_t i=0; i<npixels; i++)
598 {
599 (*fHiGainArray)[i] = new MHCalibrationPix(Form("%s%s",fHistName.Data(),"HiGainPix"),
600 Form("%s%s",fHistTitle.Data()," High Gain Pixel"));
601
602 MHCalibrationPix &pix = (*this)[i];
603 pix.SetNbins(fNbins);
604 pix.SetFirst(fFirst);
605 pix.SetLast (fLast);
606
607 TH1F *h = pix.GetHGausHist();
608
609 h->SetName (Form("%s%s%s","H",fHistName.Data(),"HiGainPix"));
610 h->SetTitle(Form("%s%s",fHistTitle.Data()," High Gain Pixel "));
611 h->SetXTitle(fHistXTitle.Data());
612 h->SetYTitle(fHistYTitle.Data());
613
614 InitHists(pix,(*fBadPixels)[i],i);
615 }
616 }
617
618 if (fAverageHiGainAreas->GetEntries()==0)
619 {
620 fAverageHiGainAreas->Expand(nareas);
621
622 for (Int_t j=0; j<nareas; j++)
623 {
624 (*fAverageHiGainAreas)[j] = new MHCalibrationPix(Form("%s%s",fHistName.Data(),"HiGainArea"),
625 Form("%s%s",fHistTitle.Data()," High Gain Area Idx "));
626
627 MHCalibrationPix &pix = GetAverageHiGainArea(j);
628
629 pix.SetNbins(fNbins*(Int_t)TMath::Sqrt((Float_t)npixels/nareas));
630 pix.SetFirst(fFirst);
631 pix.SetLast (fLast);
632
633 TH1F *h = pix.GetHGausHist();
634
635 h->SetName (Form("%s%s%s","H",fHistName.Data(),"HiGainArea"));
636 h->SetXTitle(fHistXTitle.Data());
637 h->SetYTitle(fHistYTitle.Data());
638
639 if (fGeom && fGeom->InheritsFrom("MGeomCamMagic"))
640 {
641 h->SetTitle(Form("%s%s%s%s",fHistTitle.Data()," averaged on event-by-event basis ",
642 j==0 ? "Inner Pixels " : "Outer Pixels ","High Gain Runs: "));
643 pix.InitBins();
644 pix.SetEventFrequency(fPulserFrequency);
645 }
646 else
647 {
648 h->SetTitle(Form("%s%s",fHistTitle.Data()," averaged on event-by-event basis High Gain Area Idx "));
649 InitHists(pix,fIntensCam ? fIntensCam->GetAverageBadArea(j) : fCam->GetAverageBadArea(j),j);
650 }
651 }
652 }
653
654 if (fAverageHiGainSectors->GetEntries()==0)
655 {
656 fAverageHiGainSectors->Expand(nsectors);
657
658 for (Int_t j=0; j<nsectors; j++)
659 {
660 (*fAverageHiGainSectors)[j] = new MHCalibrationPix(Form("%s%s",fHistName.Data(),"HiGainSector"),
661 Form("%s%s",fHistTitle.Data()," High Gain Sector "));
662 MHCalibrationPix &pix = GetAverageHiGainSector(j);
663
664 pix.SetNbins(fNbins*(Int_t)TMath::Sqrt((Float_t)npixels/nsectors));
665 pix.SetFirst(fFirst);
666 pix.SetLast (fLast);
667
668 TH1F *h = pix.GetHGausHist();
669
670 h->SetName (Form("%s%s%s","H",fHistName.Data(),"HiGainSector"));
671 h->SetTitle(Form("%s%s",fHistTitle.Data()," High Gain Sector "));
672 h->SetXTitle(fHistXTitle.Data());
673 h->SetYTitle(fHistYTitle.Data());
674
675 InitHists(pix,fIntensCam ? fIntensCam->GetAverageBadSector(j) : fCam->GetAverageBadSector(j),j);
676 }
677 }
678}
679
680//--------------------------------------------------------------------------------------
681//
682// Return, if IsLoGain() is kFALSE
683//
684// Initializes the Low Gain Arrays:
685//
686// - Expand fLoGainArrays to npixels
687// - Expand fAverageLoGainAreas to nareas
688// - Expand fAverageLoGainSectors to nsectors
689//
690// - For every entry in the expanded arrays:
691// * Initialize an MHCalibrationPix
692// * Set Binning from fNbins, fFirst and fLast
693// * Set Histgram names and titles from fHistName and fHistTitle
694// * Set X-axis and Y-axis titles with fHistXTitle and fHistYTitle
695// * Call InitHists
696//
697void MHCalibrationCam::InitLoGainArrays(const Int_t npixels, const Int_t nareas, const Int_t nsectors)
698{
699
700 if (!IsLoGain())
701 return;
702
703 if (fLoGainArray->GetEntries()==0)
704 {
705 fLoGainArray->Expand(npixels);
706 for (Int_t i=0; i<npixels; i++)
707 {
708 (*fLoGainArray)[i] = new MHCalibrationPix(Form("%s%s",fHistName.Data(),"LoGainPix"),
709 Form("%s%s",fHistTitle.Data()," Low Gain Pixel"));
710
711 MHCalibrationPix &pix = (*this)[i];
712 pix.SetNbins(fNbins);
713 pix.SetFirst(fFirst);
714 pix.SetLast (fLast);
715
716 TH1F *h = pix.GetHGausHist();
717
718 h->SetName (Form("%s%s%s","H",fHistName.Data(),"LoGainPix"));
719 h->SetTitle(Form("%s%s",fHistTitle.Data()," Low Gain Pixel "));
720 h->SetXTitle(fHistXTitle.Data());
721 h->SetYTitle(fHistYTitle.Data());
722
723 InitHists(pix,(*fBadPixels)[i],i);
724 }
725 }
726
727 if (fAverageLoGainAreas->GetEntries()==0)
728 {
729 fAverageLoGainAreas->Expand(nareas);
730
731 for (Int_t j=0; j<nareas; j++)
732 {
733 (*fAverageLoGainAreas)[j] = new MHCalibrationPix(Form("%s%s",fHistName.Data(),"LoGainArea"),
734 Form("%s%s",fHistTitle.Data()," Low Gain Area Idx "));
735
736 MHCalibrationPix &pix = GetAverageLoGainArea(j);
737
738 pix.SetNbins(fNbins*(Int_t)TMath::Sqrt((Float_t)npixels/nareas));
739 pix.SetFirst(fFirst);
740 pix.SetLast (fLast);
741
742 TH1F *h = pix.GetHGausHist();
743
744 h->SetName (Form("%s%s%s","H",fHistName.Data(),"LoGainArea"));
745 h->SetXTitle(fHistXTitle.Data());
746 h->SetYTitle(fHistYTitle.Data());
747
748 if (fGeom && fGeom->InheritsFrom("MGeomCamMagic"))
749 {
750 h->SetTitle(Form("%s%s%s%s",fHistTitle.Data()," averaged on event-by-event basis ",
751 j==0 ? "Inner Pixels " : "Outer Pixels ","Low Gain Runs: "));
752 pix.InitBins();
753 pix.SetEventFrequency(fPulserFrequency);
754 }
755 else
756 {
757 h->SetTitle(Form("%s%s",fHistTitle.Data()," averaged on event-by-event basis Low Gain Area Idx "));
758 InitHists(pix,fIntensCam ? fIntensCam->GetAverageBadArea(j) : fCam->GetAverageBadArea(j),j);
759 }
760 }
761 }
762
763 if (fAverageLoGainSectors->GetEntries()==0)
764 {
765 fAverageLoGainSectors->Expand(nsectors);
766
767 for (Int_t j=0; j<nsectors; j++)
768 {
769 (*fAverageLoGainSectors)[j] = new MHCalibrationPix(Form("%s%s",fHistName.Data(),"LoGainSector"),
770 Form("%s%s",fHistTitle.Data()," Low Gain Sector "));
771 MHCalibrationPix &pix = GetAverageLoGainSector(j);
772
773 pix.SetNbins(fNbins*(Int_t)TMath::Sqrt((Float_t)npixels/nsectors));
774 pix.SetFirst(fFirst);
775 pix.SetLast (fLast);
776
777 TH1F *h = pix.GetHGausHist();
778
779 h->SetName (Form("%s%s%s","H",fHistName.Data(),"LoGainSector"));
780 h->SetTitle(Form("%s%s",fHistTitle.Data()," Low Gain Sector "));
781 h->SetXTitle(fHistXTitle.Data());
782 h->SetYTitle(fHistYTitle.Data());
783
784 InitHists(pix,fIntensCam ? fIntensCam->GetAverageBadSector(j) : fCam->GetAverageBadSector(j),j);
785 }
786 }
787}
788
789//--------------------------------------------------------------------------------
790//
791// Retrieves from MGeomCam:
792// - number of pixels
793// - number of pixel areas
794// - number of sectors
795//
796// Return kFALSE, if sizes of the TObjArrays do not match npixels, nareas or nsectors
797//
798// Call FillHists()
799//
800Bool_t MHCalibrationCam::Fill(const MParContainer *par, const Stat_t w)
801{
802
803 const Int_t npixels = fGeom->GetNumPixels();
804 const Int_t nareas = fGeom->GetNumAreas();
805 const Int_t nsectors = fGeom->GetNumSectors();
806
807 //
808 // Hi-Gain ObjArrays
809 //
810 if (fHiGainArray->GetEntries() != npixels)
811 {
812 *fLog << err << "ERROR - Size mismatch in number of pixels... abort." << endl;
813 return kFALSE;
814 }
815
816 if (fAverageHiGainAreas->GetEntries() != nareas)
817 {
818 *fLog << err << "ERROR - Size mismatch in number of areas ... abort." << endl;
819 return kFALSE;
820 }
821
822 if (fAverageHiGainSectors->GetEntries() != nsectors)
823 {
824 *fLog << err << "ERROR - Size mismatch in number of sectors ... abort." << endl;
825 return kFALSE;
826 }
827
828 if (IsLoGain())
829 {
830 if (fLoGainArray->GetEntries() != npixels)
831 {
832 *fLog << err << "ERROR - Size mismatch in number of pixels... abort." << endl;
833 return kFALSE;
834 }
835
836 if (fAverageLoGainAreas->GetEntries() != nareas)
837 {
838 *fLog << err << "ERROR - Size mismatch in number of areas ... abort." << endl;
839 return kFALSE;
840 }
841
842 if (fAverageLoGainSectors->GetEntries() != nsectors)
843 {
844 *fLog << err << "ERROR - Size mismatch in number of sectors ... abort." << endl;
845 return kFALSE;
846 }
847 }
848
849 return FillHists(par, w);
850}
851
852// --------------------------------------------------------------------------
853//
854// 0) Ask if fHiGainArray and fLoGainArray have been initialized,
855// otherwise return kFALSE.
856// 1) FinalizeHists()
857// 2) FinalizeBadPixels()
858// 3) CalcAverageSigma()
859//
860Bool_t MHCalibrationCam::Finalize()
861{
862
863 if (fHiGainArray->GetEntries() == 0 && fLoGainArray->GetEntries() == 0)
864 {
865 *fLog << err << GetDescriptor()
866 << ": ERROR: Both (HiGain and LoGain) histogram arrays have not been initialized... abort." << endl;
867 return kFALSE;
868 }
869
870 if (!FinalizeHists())
871 return kFALSE;
872
873 FinalizeBadPixels();
874 CalcAverageSigma();
875
876 return kTRUE;
877}
878
879// -------------------------------------------------------------
880//
881// If MBadPixelsPix::IsBad():
882// - calls MHCalibrationPix::SetExcluded()
883//
884// Calls:
885// - MHGausEvents::InitBins()
886// - MHCalibrationPix::ChangeHistId(i)
887// - MHCalibrationPix::SetEventFrequency(fPulserFrequency)
888//
889void MHCalibrationCam::InitHists(MHCalibrationPix &hist, MBadPixelsPix &bad, const Int_t i)
890{
891
892 if (bad.IsUnsuitable(MBadPixelsPix::kUnsuitableRun))
893 hist.SetExcluded();
894
895 hist.InitBins();
896 hist.ChangeHistId(i);
897 hist.SetEventFrequency(fPulserFrequency);
898
899 TH1F *h = hist.GetHGausHist();
900 h->SetTitle( Form("%s%s", h->GetTitle()," Runs: "));
901}
902
903// --------------------------------------------------------------------------
904//
905// Calls FitHiGainHists for every entry in:
906// - fHiGainArray
907// - fAverageHiGainAreas
908// - fAverageHiGainSectors
909//
910void MHCalibrationCam::FitHiGainArrays(MCalibrationCam &calcam, MBadPixelsCam &badcam,
911 MBadPixelsPix::UncalibratedType_t fittyp,
912 MBadPixelsPix::UncalibratedType_t osctyp)
913{
914
915 for (Int_t i=0; i<fHiGainArray->GetSize(); i++)
916 {
917
918 MHCalibrationPix &hist = (*this)[i];
919
920 if (hist.IsExcluded())
921 continue;
922
923 MCalibrationPix &pix = calcam[i];
924 MBadPixelsPix &bad = badcam[i];
925
926 FitHiGainHists(hist,pix,bad,fittyp,osctyp);
927
928 }
929
930 for (Int_t j=0; j<fAverageHiGainAreas->GetSize(); j++)
931 {
932
933 MHCalibrationPix &hist = GetAverageHiGainArea(j);
934 MCalibrationPix &pix = calcam.GetAverageArea(j);
935 MBadPixelsPix &bad = calcam.GetAverageBadArea(j);
936
937 FitHiGainHists(hist,pix,bad,fittyp,osctyp);
938 }
939
940
941 for (Int_t j=0; j<fAverageHiGainSectors->GetSize(); j++)
942 {
943
944 MHCalibrationPix &hist = GetAverageHiGainSector(j);
945 MCalibrationPix &pix = calcam.GetAverageSector(j);
946 MBadPixelsPix &bad = calcam.GetAverageBadSector(j);
947
948 FitHiGainHists(hist,pix,bad,fittyp,osctyp);
949 }
950
951}
952
953// --------------------------------------------------------------------------
954//
955// Calls FitLoGainHists for every entry in:
956// - fLoGainArray
957// - fAverageLoGainAreas
958// - fAverageLoGainSectors
959//
960void MHCalibrationCam::FitLoGainArrays(MCalibrationCam &calcam, MBadPixelsCam &badcam,
961 MBadPixelsPix::UncalibratedType_t fittyp,
962 MBadPixelsPix::UncalibratedType_t osctyp)
963{
964
965 if (!IsLoGain())
966 return;
967
968 for (Int_t i=0; i<fLoGainArray->GetSize(); i++)
969 {
970
971 MHCalibrationPix &hist = (*this)(i);
972
973 if (hist.IsExcluded())
974 continue;
975
976 MCalibrationPix &pix = calcam[i];
977 MBadPixelsPix &bad = badcam[i];
978
979 FitLoGainHists(hist,pix,bad,fittyp,osctyp);
980
981 }
982
983 for (Int_t j=0; j<fAverageLoGainAreas->GetSize(); j++)
984 {
985
986 MHCalibrationPix &hist = GetAverageLoGainArea(j);
987 MCalibrationPix &pix = calcam.GetAverageArea(j);
988 MBadPixelsPix &bad = calcam.GetAverageBadArea(j);
989
990 FitLoGainHists(hist,pix,bad,fittyp,osctyp);
991 }
992
993
994 for (Int_t j=0; j<fAverageLoGainSectors->GetSize(); j++)
995 {
996
997 MHCalibrationPix &hist = GetAverageLoGainSector(j);
998 MCalibrationPix &pix = calcam.GetAverageSector(j);
999 MBadPixelsPix &bad = calcam.GetAverageBadSector(j);
1000
1001 FitLoGainHists(hist,pix,bad,fittyp,osctyp);
1002 }
1003}
1004
1005//------------------------------------------------------------
1006//
1007// For all averaged areas, the fitted sigma is multiplied with the square root of
1008// the number involved pixels
1009//
1010void MHCalibrationCam::CalcAverageSigma()
1011{
1012
1013 if (!fCam)
1014 return;
1015
1016 for (UInt_t j=0; j<fGeom->GetNumAreas(); j++)
1017 {
1018
1019 MCalibrationPix &pix = fCam->GetAverageArea(j);
1020
1021 const Float_t numsqr = TMath::Sqrt((Float_t)fAverageAreaNum[j]);
1022 fAverageAreaSigma[j] = pix.GetSigma () * numsqr;
1023 fAverageAreaSigmaVar[j] = pix.GetSigmaErr () * pix.GetSigmaErr() * numsqr;
1024
1025 pix.SetSigma (fAverageAreaSigma[j]);
1026 pix.SetSigmaVar(fAverageAreaSigmaVar[j]);
1027
1028 fAverageAreaRelSigma [j] = fAverageAreaSigma[j] / pix.GetMean();
1029 fAverageAreaRelSigmaVar[j] = fAverageAreaSigmaVar[j] / (fAverageAreaSigma[j]*fAverageAreaSigma[j]);
1030 fAverageAreaRelSigmaVar[j] += pix.GetMeanRelVar();
1031 fAverageAreaRelSigmaVar[j] *= fAverageAreaRelSigma[j];
1032 }
1033}
1034
1035// ---------------------------------------------------------------------------
1036//
1037// Returns if the histogram is empty and sets the following flag:
1038// - MBadPixelsPix::SetUnsuitable(MBadPixelsPix::kUnsuitableRun)
1039//
1040// Fits the histograms with a Gaussian, in case of failure
1041// calls MHCalibrationPix::RepeatFit(), in case of repeated failure
1042// calls MHCalibrationPix::BypassFit() and sets the following flags:
1043// - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::UncalibratedType_t fittyp )
1044// - MBadPixelsPix::SetUnsuitable( MBadPixelsPix::kUnreliableRun )
1045//
1046// Creates the fourier spectrum and tests MHGausEvents::IsFourierSpectrumOK().
1047// In case no, sets the following flags:
1048// - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::UncalibratedType_t osctyp )
1049// - MBadPixelsPix::SetUnsuitable( MBadPixelsPix::kUnreliableRun )
1050//
1051// Retrieves the results and store them in MCalibrationPix
1052//
1053void MHCalibrationCam::FitHiGainHists(MHCalibrationPix &hist,
1054 MCalibrationPix &pix,
1055 MBadPixelsPix &bad,
1056 MBadPixelsPix::UncalibratedType_t fittyp,
1057 MBadPixelsPix::UncalibratedType_t osctyp)
1058{
1059
1060
1061 if (hist.IsEmpty() || hist.IsOnlyOverflow() || hist.IsOnlyUnderflow())
1062 return;
1063
1064 //
1065 // 2) Fit the Hi Gain histograms with a Gaussian
1066 //
1067 if (!hist.FitGaus())
1068 //
1069 // 3) In case of failure set the bit Fitted to false and take histogram means and RMS
1070 //
1071 if (!hist.RepeatFit())
1072 {
1073 hist.BypassFit();
1074 bad.SetUncalibrated( fittyp );
1075 }
1076
1077 hist.Renorm();
1078 //
1079 // 4) Check for oscillations
1080 //
1081 if (IsOscillations())
1082 {
1083 hist.CreateFourierSpectrum();
1084
1085 if (!hist.IsFourierSpectrumOK())
1086 bad.SetUncalibrated( osctyp );
1087 }
1088
1089 //
1090 // 5) Retrieve the results and store them in this class
1091 //
1092 pix.SetHiGainMean ( hist.GetMean() );
1093 pix.SetHiGainMeanVar ( hist.GetMeanErr() * hist.GetMeanErr() );
1094 pix.SetHiGainSigma ( hist.GetSigma() );
1095 pix.SetHiGainSigmaVar ( hist.GetSigmaErr()* hist.GetSigmaErr() );
1096 pix.SetHiGainProb ( hist.GetProb() );
1097 pix.SetHiGainNumBlackout( hist.GetBlackout() );
1098 pix.SetHiGainNumPickup ( hist.GetPickup() );
1099
1100 if (IsDebug())
1101 {
1102 *fLog << dbginf << GetDescriptor() << ": ID " << hist.GetPixId()
1103 << " HiGainSaturation: " << pix.IsHiGainSaturation()
1104 << " HiGainMean: " << hist.GetMean ()
1105 << " HiGainMeanErr: " << hist.GetMeanErr ()
1106 << " HiGainMeanSigma: " << hist.GetSigma ()
1107 << " HiGainMeanSigmaErr: " << hist.GetSigmaErr()
1108 << " HiGainMeanProb: " << hist.GetProb ()
1109 << " HiGainNumBlackout: " << hist.GetBlackout()
1110 << " HiGainNumPickup : " << hist.GetPickup ()
1111 << endl;
1112 }
1113
1114}
1115
1116
1117// ---------------------------------------------------------------------------
1118//
1119// Returns if the histogram is empty and sets the following flag:
1120// - MBadPixelsPix::SetUnsuitable(MBadPixelsPix::kUnsuitableRun)
1121//
1122// Fits the histograms with a Gaussian, in case of failure
1123// calls MHCalibrationPix::RepeatFit(), in case of repeated failure
1124// calls MHCalibrationPix::BypassFit() and sets the following flags:
1125// - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::UncalibratedType_t fittyp )
1126// - MBadPixelsPix::SetUnsuitable( MBadPixelsPix::kUnreliableRun )
1127//
1128// Creates the fourier spectrum and tests MHGausEvents::IsFourierSpectrumOK().
1129// In case no, sets the following flags:
1130// - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::UncalibratedType_t osctyp )
1131// - MBadPixelsPix::SetUnsuitable( MBadPixelsPix::kUnreliableRun )
1132//
1133// Retrieves the results and store them in MCalibrationPix
1134//
1135void MHCalibrationCam::FitLoGainHists(MHCalibrationPix &hist,
1136 MCalibrationPix &pix,
1137 MBadPixelsPix &bad,
1138 MBadPixelsPix::UncalibratedType_t fittyp,
1139 MBadPixelsPix::UncalibratedType_t osctyp)
1140{
1141
1142 if (hist.IsEmpty() || hist.IsOnlyOverflow() || hist.IsOnlyUnderflow())
1143 return;
1144
1145
1146 //
1147 // 2) Fit the Hi Gain histograms with a Gaussian
1148 //
1149 if (!hist.FitGaus())
1150 //
1151 // 3) In case of failure set the bit Fitted to false and take histogram means and RMS
1152 //
1153 if (!hist.RepeatFit())
1154 {
1155 hist.BypassFit();
1156 bad.SetUncalibrated( fittyp );
1157 }
1158
1159 //
1160 // 4) Check for oscillations
1161 //
1162 if (IsOscillations())
1163 {
1164 hist.CreateFourierSpectrum();
1165
1166 if (!hist.IsFourierSpectrumOK())
1167 bad.SetUncalibrated( osctyp );
1168 }
1169
1170 //
1171 // 5) Retrieve the results and store them in this class
1172 //
1173 pix.SetLoGainMean ( hist.GetMean() );
1174 pix.SetLoGainMeanVar ( hist.GetMeanErr() * hist.GetMeanErr() );
1175 pix.SetLoGainSigma ( hist.GetSigma() );
1176 pix.SetLoGainSigmaVar ( hist.GetSigmaErr() * hist.GetSigmaErr() );
1177 pix.SetLoGainProb ( hist.GetProb() );
1178 pix.SetLoGainNumBlackout( hist.GetBlackout() );
1179 pix.SetLoGainNumPickup ( hist.GetPickup() );
1180
1181 if (IsDebug())
1182 {
1183 *fLog << dbginf << GetDescriptor() << "ID: " << hist.GetPixId()
1184 << " HiGainSaturation: " << pix.IsHiGainSaturation()
1185 << " LoGainMean: " << hist.GetMean ()
1186 << " LoGainMeanErr: " << hist.GetMeanErr ()
1187 << " LoGainMeanSigma: " << hist.GetSigma ()
1188 << " LoGainMeanSigmaErr: " << hist.GetSigmaErr()
1189 << " LoGainMeanProb: " << hist.GetProb ()
1190 << " LoGainNumBlackout: " << hist.GetBlackout()
1191 << " LoGainNumPickup : " << hist.GetPickup ()
1192 << endl;
1193 }
1194
1195}
1196
1197
1198
1199// --------------------------------------------------------------------------
1200//
1201// Dummy, needed by MCamEvent
1202//
1203Bool_t MHCalibrationCam::GetPixelContent(Double_t &val, Int_t idx, const MGeomCam &cam, Int_t type) const
1204{
1205 return kTRUE;
1206}
1207
1208// --------------------------------------------------------------------------
1209//
1210// What MHCamera needs in order to draw an individual pixel in the camera
1211//
1212void MHCalibrationCam::DrawPixelContent(Int_t idx) const
1213{
1214}
1215
1216// -----------------------------------------------------------------------------
1217//
1218// Default draw:
1219//
1220// Displays the averaged areas, both High Gain and Low Gain
1221//
1222// Calls the Draw of the fAverageHiGainAreas and fAverageLoGainAreas objects with options
1223//
1224void MHCalibrationCam::Draw(const Option_t *opt)
1225{
1226
1227 const Int_t nareas = fAverageHiGainAreas->GetEntries();
1228 if (nareas == 0)
1229 return;
1230
1231 TVirtualPad *pad = gPad ? gPad : MH::MakeDefCanvas(this);
1232 pad->SetBorderMode(0);
1233
1234 pad->Divide(IsLoGain() ? 2 : 1,nareas);
1235
1236 for (Int_t i=0; i<nareas;i++)
1237 {
1238
1239 pad->cd(IsLoGain() ? 2*(i+1)-1 : i+1);
1240 GetAverageHiGainArea(i).Draw(opt);
1241
1242 if (!fAverageAreaSat[i])
1243 DrawAverageSigma(fAverageAreaSat[i], i,
1244 fAverageAreaSigma[i], fAverageAreaSigmaVar[i],
1245 fAverageAreaRelSigma[i], fAverageAreaRelSigmaVar[i]);
1246
1247 if (IsLoGain())
1248 {
1249 pad->cd(2*(i+1));
1250 GetAverageLoGainArea(i).Draw(opt);
1251 }
1252
1253 if (fAverageAreaSat[i])
1254 DrawAverageSigma(fAverageAreaSat[i], i,
1255 fAverageAreaSigma[i], fAverageAreaSigmaVar[i],
1256 fAverageAreaRelSigma[i], fAverageAreaRelSigmaVar[i]);
1257 }
1258
1259}
1260
1261// -----------------------------------------------------------------------------
1262//
1263// Default draw:
1264//
1265// Displays a TPaveText with the re-normalized sigmas of the average area
1266//
1267void MHCalibrationCam::DrawAverageSigma(Bool_t sat, Bool_t inner,
1268 Float_t sigma, Float_t sigmavar,
1269 Float_t relsigma, Float_t relsigmavar) const
1270{
1271
1272 if (sigma != 0 && sigmavar >= 0 && relsigmavar >= 0.)
1273 {
1274
1275 TPad *newpad = new TPad("newpad","transparent",0,0,1,1);
1276 newpad->SetFillStyle(4000);
1277 newpad->Draw();
1278 newpad->cd();
1279
1280 TPaveText *text = new TPaveText(sat? 0.1 : 0.35,0.7,sat ? 0.4 : 0.7,1.0);
1281 text->SetTextSize(0.07);
1282 const TString line1 = Form("%s%s%s",inner ? "Outer" : "Inner",
1283 " Pixels ", sat ? "Low Gain" : "High Gain");
1284 TText *txt1 = text->AddText(line1.Data());
1285 const TString line2 = Form("#sigma per pix: %2.2f #pm %2.2f",sigma,TMath::Sqrt(sigmavar));
1286 TText *txt2 = text->AddText(line2.Data());
1287 const TString line3 = Form("Rel. #sigma per pix: %2.2f #pm %2.2f",relsigma,TMath::Sqrt(relsigmavar));
1288 TText *txt3 = text->AddText(line3.Data());
1289 text->Draw("");
1290
1291 text->SetBit(kCanDelete);
1292 txt1->SetBit(kCanDelete);
1293 txt2->SetBit(kCanDelete);
1294 txt3->SetBit(kCanDelete);
1295 newpad->SetBit(kCanDelete);
1296 }
1297}
1298
1299Int_t MHCalibrationCam::ReadEnv(const TEnv &env, TString prefix, Bool_t print)
1300{
1301
1302 Bool_t rc = kFALSE;
1303 if (IsEnvDefined(env, prefix, "Debug", print))
1304 {
1305 SetDebug(GetEnvValue(env, prefix, "Debug", IsDebug()));
1306 rc = kTRUE;
1307 }
1308 if (IsEnvDefined(env, prefix, "LoGain", print))
1309 {
1310 SetDebug(GetEnvValue(env, prefix, "LoGain", IsLoGain()));
1311 rc = kTRUE;
1312 }
1313 if (IsEnvDefined(env, prefix, "Oscillations", print))
1314 {
1315 SetOscillations(GetEnvValue(env, prefix, "Oscillations", IsOscillations()));
1316 rc = kTRUE;
1317 }
1318
1319 if (IsEnvDefined(env, prefix, "Nbins", print))
1320 {
1321 SetNbins(GetEnvValue(env, prefix, "Nbins", fNbins));
1322 rc = kTRUE;
1323 }
1324 if (IsEnvDefined(env, prefix, "First", print))
1325 {
1326 SetFirst(GetEnvValue(env, prefix, "First", fFirst));
1327 rc = kTRUE;
1328 }
1329 if (IsEnvDefined(env, prefix, "Last", print))
1330 {
1331 SetLast(GetEnvValue(env, prefix, "Last", fLast));
1332 rc = kTRUE;
1333 }
1334
1335 return rc;
1336}
Note: See TracBrowser for help on using the repository browser.