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

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