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

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