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

Last change on this file since 4941 was 4929, checked in by gaug, 20 years ago
*** empty log message ***
File size: 33.2 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
53#include <TVirtualPad.h>
54#include <TCanvas.h>
55#include <TPad.h>
56#include <TText.h>
57#include <TPaveText.h>
58
59#include "MLog.h"
60#include "MLogManip.h"
61
62#include "MCalibrationPix.h"
63#include "MCalibrationCam.h"
64
65#include "MHCalibrationPix.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::fgAverageNbins = 2000;
82const Int_t MHCalibrationCam::fgPulserFrequency = 500;
83
84// --------------------------------------------------------------------------
85//
86// Default Constructor.
87//
88// Sets:
89// - all pointers to NULL
90//
91// Initializes and sets owner of:
92// - fHiGainArray, fLoGainArray
93// - fAverageHiGainAreas, fAverageLoGainAreas
94// - fAverageHiGainSectors, fAverageLoGainSectors
95//
96// Initializes:
97// - fPulserFrequency to fgPulserFrequency
98//
99MHCalibrationCam::MHCalibrationCam(const char *name, const char *title)
100 : fColor(MCalibrationCam::kNONE),
101 fBadPixels(NULL), fIntensCam(NULL), fCam(NULL), fGeom(NULL), fRunHeader(NULL)
102{
103
104 fHiGainArray = new TObjArray;
105 fHiGainArray->SetOwner();
106
107 fLoGainArray = new TObjArray;
108 fLoGainArray->SetOwner();
109
110 fAverageHiGainAreas = new TObjArray;
111 fAverageHiGainAreas->SetOwner();
112
113 fAverageLoGainAreas = new TObjArray;
114 fAverageLoGainAreas->SetOwner();
115
116 fAverageHiGainSectors = new TObjArray;
117 fAverageHiGainSectors->SetOwner();
118
119 fAverageLoGainSectors = new TObjArray;
120 fAverageLoGainSectors->SetOwner();
121
122 SetAverageNbins();
123 SetPulserFrequency();
124
125 SetDebug (kFALSE);
126 SetLoGain (kTRUE);
127 SetOscillations(kTRUE);
128}
129
130const Int_t MHCalibrationCam::GetSize() const
131{
132 return fHiGainArray->GetSize();
133}
134
135// --------------------------------------------------------------------------
136//
137// Deletes the TClonesArray of:
138// - fHiGainArray, fLoGainArray
139// - fAverageHiGainAreas, fAverageLoGainAreas
140// - fAverageHiGainSectors, fAverageLoGainSectors
141//
142MHCalibrationCam::~MHCalibrationCam()
143{
144
145 delete fHiGainArray;
146 delete fLoGainArray;
147
148 delete fAverageHiGainAreas;
149 delete fAverageLoGainAreas;
150
151 delete fAverageHiGainSectors;
152 delete fAverageLoGainSectors;
153}
154
155void MHCalibrationCam::ResetHists()
156{
157
158 if (fHiGainArray)
159 { fHiGainArray->ForEach(MHCalibrationPix,Reset)(); }
160 if (fLoGainArray)
161 { fLoGainArray->ForEach(MHCalibrationPix,Reset)(); }
162
163 if (fAverageHiGainAreas)
164 { fAverageHiGainAreas->ForEach(MHCalibrationPix,Reset)(); }
165 if (fAverageLoGainAreas)
166 { fAverageLoGainAreas->ForEach(MHCalibrationPix,Reset)(); }
167 if (fAverageHiGainSectors)
168 { fAverageHiGainSectors->ForEach(MHCalibrationPix,Reset)(); }
169 if (fAverageLoGainSectors)
170 { fAverageLoGainSectors->ForEach(MHCalibrationPix,Reset)(); }
171
172}
173
174
175
176// --------------------------------------------------------------------------
177//
178// Get i-th High Gain pixel (pixel number)
179//
180MHCalibrationPix &MHCalibrationCam::operator[](UInt_t i)
181{
182 return *static_cast<MHCalibrationPix*>(fHiGainArray->UncheckedAt(i));
183}
184
185// --------------------------------------------------------------------------
186//
187// Get i-th High Gain pixel (pixel number)
188//
189const MHCalibrationPix &MHCalibrationCam::operator[](UInt_t i) const
190{
191 return *static_cast<MHCalibrationPix*>(fHiGainArray->UncheckedAt(i));
192}
193
194// --------------------------------------------------------------------------
195//
196// Get i-th Low Gain pixel (pixel number)
197//
198MHCalibrationPix &MHCalibrationCam::operator()(UInt_t i)
199{
200 return *static_cast<MHCalibrationPix*>(fLoGainArray->UncheckedAt(i));
201}
202
203// --------------------------------------------------------------------------
204//
205// Get i-th Low Gain pixel (pixel number)
206//
207const MHCalibrationPix &MHCalibrationCam::operator()(UInt_t i) const
208{
209 return *static_cast<MHCalibrationPix*>(fLoGainArray->UncheckedAt(i));
210}
211
212// --------------------------------------------------------------------------
213//
214// Returns the current size of the TObjArray fAverageHiGainAreas
215// independently if the MHCalibrationPix is filled with values or not.
216//
217const Int_t MHCalibrationCam::GetAverageAreas() const
218{
219 return fAverageHiGainAreas->GetEntries();
220}
221
222// --------------------------------------------------------------------------
223//
224// Get i-th High Gain pixel Area (area number)
225//
226MHCalibrationPix &MHCalibrationCam::GetAverageHiGainArea(UInt_t i)
227{
228 return *static_cast<MHCalibrationPix*>(fAverageHiGainAreas->UncheckedAt(i));
229}
230
231// --------------------------------------------------------------------------
232//
233// Get i-th High Gain pixel Area (area number)
234//
235const MHCalibrationPix &MHCalibrationCam::GetAverageHiGainArea(UInt_t i) const
236{
237 return *static_cast<MHCalibrationPix *>(fAverageHiGainAreas->UncheckedAt(i));
238}
239
240// --------------------------------------------------------------------------
241//
242// Get i-th Low Gain pixel Area (area number)
243//
244MHCalibrationPix &MHCalibrationCam::GetAverageLoGainArea(UInt_t i)
245{
246 return *static_cast<MHCalibrationPix*>(fAverageLoGainAreas->UncheckedAt(i));
247}
248
249// --------------------------------------------------------------------------
250//
251// Get i-th Low Gain pixel Area (area number)
252//
253const MHCalibrationPix &MHCalibrationCam::GetAverageLoGainArea(UInt_t i) const
254{
255 return *static_cast<MHCalibrationPix*>(fAverageLoGainAreas->UncheckedAt(i));
256}
257
258// --------------------------------------------------------------------------
259//
260// Returns the current size of the TObjArray fAverageHiGainSectors
261// independently if the MHCalibrationPix is filled with values or not.
262//
263const Int_t MHCalibrationCam::GetAverageSectors() const
264{
265 return fAverageHiGainSectors->GetEntries();
266}
267
268// --------------------------------------------------------------------------
269//
270// Get i-th High Gain Sector (sector number)
271//
272MHCalibrationPix &MHCalibrationCam::GetAverageHiGainSector(UInt_t i)
273{
274 return *static_cast<MHCalibrationPix*>(fAverageHiGainSectors->UncheckedAt(i));
275}
276
277// --------------------------------------------------------------------------
278//
279// Get i-th High Gain Sector (sector number)
280//
281const MHCalibrationPix &MHCalibrationCam::GetAverageHiGainSector(UInt_t i) const
282{
283 return *static_cast<MHCalibrationPix*>(fAverageHiGainSectors->UncheckedAt(i));
284}
285
286// --------------------------------------------------------------------------
287//
288// Get i-th Low Gain Sector (sector number)
289//
290MHCalibrationPix &MHCalibrationCam::GetAverageLoGainSector(UInt_t i)
291{
292 return *static_cast<MHCalibrationPix*>(fAverageLoGainSectors->UncheckedAt(i));
293}
294
295// --------------------------------------------------------------------------
296//
297// Get i-th Low Gain Sector (sector number)
298//
299const MHCalibrationPix &MHCalibrationCam::GetAverageLoGainSector(UInt_t i) const
300{
301 return *static_cast<MHCalibrationPix*>(fAverageLoGainSectors->UncheckedAt(i));
302}
303
304
305const TArrayI &MHCalibrationCam::GetRunNumbers() const
306{
307 return fRunNumbers;
308}
309
310// --------------------------------------------------------------------------
311//
312// Our own clone function is necessary since root 3.01/06 or Mars 0.4
313// I don't know the reason.
314//
315// Creates new MHCalibrationCam
316// Deletes the TObjArray's and Clones them individually
317// Copies the TArray's
318// Copies the fPulserFrequency
319//
320TObject *MHCalibrationCam::Clone(const char *) const
321{
322
323 // const Int_t nhi = fHiGainArray->GetEntries();
324 // const Int_t nlo = fLoGainArray->GetEntries();
325 const Int_t navhi = fAverageHiGainAreas->GetEntries();
326 const Int_t navlo = fAverageLoGainAreas->GetEntries();
327 const Int_t nsehi = fAverageHiGainSectors->GetEntries();
328 const Int_t nselo = fAverageLoGainSectors->GetEntries();
329
330 //
331 // FIXME, this might be done faster and more elegant, by direct copy.
332 //
333 MHCalibrationCam *cam = new MHCalibrationCam();
334
335 // cam->fHiGainArray->Expand(nhi);
336 // cam->fLoGainArray->Expand(nlo);
337 cam->fAverageHiGainAreas->Expand(navhi);
338 cam->fAverageLoGainAreas->Expand(navlo);
339 cam->fAverageHiGainSectors->Expand(nsehi);
340 cam->fAverageLoGainSectors->Expand(nselo);
341
342 /*
343 for (int i=0; i<nhi; i++)
344 {
345 delete (*cam->fHiGainArray)[i];
346 (*cam->fHiGainArray)[i] = (*fHiGainArray)[i]->Clone();
347 }
348 for (int i=0; i<nlo; i++)
349 {
350 delete (*cam->fLoGainArray)[i];
351 (*cam->fLoGainArray)[i] = (*fLoGainArray)[i]->Clone();
352 }
353 */
354
355 for (int i=0; i<navhi; i++)
356 {
357 // delete (*cam->fAverageHiGainAreas)[i];
358 (*cam->fAverageHiGainAreas)[i] = (*fAverageHiGainAreas)[i]->Clone();
359 }
360 for (int i=0; i<navlo; i++)
361 {
362 // delete (*cam->fAverageLoGainAreas)[i];
363 (*cam->fAverageLoGainAreas)[i] = (*fAverageLoGainAreas)[i]->Clone();
364 }
365 for (int i=0; i<nsehi; i++)
366 {
367 // delete (*cam->fAverageHiGainSectors)[i];
368 (*cam->fAverageHiGainSectors)[i] = (*fAverageHiGainSectors)[i]->Clone();
369 }
370 for (int i=0; i<nselo; i++)
371 {
372 // delete (*cam->fAverageLoGainSectors)[i];
373 (*cam->fAverageLoGainSectors)[i] = (*fAverageLoGainSectors)[i]->Clone();
374 }
375
376 cam->fAverageAreaNum = fAverageAreaNum;
377 cam->fAverageAreaSat = fAverageAreaSat;
378 cam->fAverageAreaSigma = fAverageAreaSigma;
379 cam->fAverageAreaSigmaVar = fAverageAreaSigmaVar;
380 cam->fAverageAreaRelSigma = fAverageAreaRelSigma;
381 cam->fAverageAreaRelSigmaVar = fAverageAreaRelSigmaVar;
382 cam->fAverageSectorNum = fAverageSectorNum;
383 cam->fRunNumbers = fRunNumbers;
384
385 cam->fPulserFrequency = fPulserFrequency;
386 cam->fAverageNbins = fAverageNbins;
387
388 return cam;
389}
390
391// --------------------------------------------------------------------------
392//
393// Gets the pointers to:
394// - MGeomCam
395//
396// Calls SetupHists(const MParList *pList)
397//
398// Calls Delete-Function of:
399// - MHCalibrationCam::fHiGainArray, MHCalibrationCam::fLoGainArray
400// - MHCalibrationCam::fAverageHiGainAreas, MHCalibrationCam::fAverageLoGainAreas
401// - MHCalibrationCam::fAverageHiGainSectors, MHCalibrationCam::fAverageLoGainSectors
402//
403Bool_t MHCalibrationCam::SetupFill(const MParList *pList)
404{
405
406 fGeom = (MGeomCam*)pList->FindObject("MGeomCam");
407 if (!fGeom)
408 {
409 *fLog << err << GetDescriptor()
410 << ": MGeomCam not found... aborting." << endl;
411 return kFALSE;
412 }
413
414 fRunHeader = (MRawRunHeader*)pList->FindObject("MRawRunHeader");
415 if (!fRunHeader)
416 {
417 *fLog << warn << GetDescriptor()
418 << ": MRawRunHeader not found... will not store run numbers." << endl;
419 }
420
421 fHiGainArray->Delete();
422 fLoGainArray->Delete();
423
424 fAverageHiGainAreas->Delete();
425 fAverageLoGainAreas->Delete();
426
427 fAverageHiGainSectors->Delete();
428 fAverageLoGainSectors->Delete();
429
430 return SetupHists(pList);
431}
432
433
434Bool_t MHCalibrationCam::SetupHists(const MParList *pList)
435{
436 return kTRUE;
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
578Bool_t MHCalibrationCam::ReInitHists(MParList *pList)
579{
580 return kTRUE;
581}
582
583//--------------------------------------------------------------------------------
584//
585// Retrieves from MGeomCam:
586// - number of pixels
587// - number of pixel areas
588// - number of sectors
589//
590// For all TObjArray's (including the averaged ones), the following steps are performed:
591//
592// 1) Test size and return kFALSE if not matching
593// 2)
594//
595Bool_t MHCalibrationCam::Fill(const MParContainer *par, const Stat_t w)
596{
597
598 const Int_t npixels = fGeom->GetNumPixels();
599 const Int_t nareas = fGeom->GetNumAreas();
600 const Int_t nsectors = fGeom->GetNumSectors();
601
602 if (fHiGainArray->GetEntries() != npixels)
603 {
604 *fLog << err << "ERROR - Size mismatch... abort." << endl;
605 return kFALSE;
606 }
607
608 if (IsLoGain())
609 if (fLoGainArray->GetEntries() != npixels)
610 {
611 *fLog << err << "ERROR - Size mismatch... abort." << endl;
612 return kFALSE;
613 }
614
615 if (fAverageHiGainAreas->GetEntries() != nareas)
616 {
617 *fLog << err << "ERROR - Size mismatch in number of areas ... abort." << endl;
618 return kFALSE;
619 }
620
621 if (IsLoGain())
622 if (fAverageLoGainAreas->GetEntries() != nareas)
623 {
624 *fLog << err << "ERROR - Size mismatch in number of areas ... abort." << endl;
625 return kFALSE;
626 }
627
628 if (fAverageHiGainSectors->GetEntries() != nsectors)
629 {
630 *fLog << err << "ERROR - Size mismatch in number of sectors ... abort." << endl;
631 return kFALSE;
632 }
633
634 if (IsLoGain())
635 if (fAverageLoGainSectors->GetEntries() != nsectors)
636 {
637 *fLog << err << "ERROR - Size mismatch in number of sectors ... abort." << endl;
638 return kFALSE;
639 }
640
641 return FillHists(par, w);
642}
643
644Bool_t MHCalibrationCam::FillHists(const MParContainer *par, const Stat_t w)
645{
646 *fLog << warn << GetDescriptor() << "FillHists not overloaded! Can't be used!" << endl;
647 return kFALSE;
648}
649
650// --------------------------------------------------------------------------
651//
652// 0) Ask if fHiGainArray and fLoGainArray have been initialized,
653// otherwise return kFALSE.
654// 1) FinalizeHists()
655// 2) FinalizeBadPixels()
656// 3) CalcAverageSigma()
657//
658Bool_t MHCalibrationCam::Finalize()
659{
660
661 if (fHiGainArray->GetEntries() == 0 && fLoGainArray->GetEntries() == 0)
662 {
663 *fLog << err << GetDescriptor()
664 << ": ERROR: Both (HiGain and LoGain) histogram arrays have not been initialized... abort." << endl;
665 return kFALSE;
666 }
667
668 if (!FinalizeHists())
669 return kFALSE;
670
671 FinalizeBadPixels();
672 CalcAverageSigma();
673
674 return kTRUE;
675}
676
677Bool_t MHCalibrationCam::FinalizeHists()
678{
679 return kTRUE;
680}
681
682void MHCalibrationCam::FinalizeBadPixels()
683{
684}
685
686
687// -------------------------------------------------------------
688//
689// If MBadPixelsPix::IsBad():
690// - calls MHCalibrationPix::SetExcluded()
691//
692// Calls:
693// - MHGausEvents::InitBins()
694// - MHCalibrationPix::ChangeHistId(i)
695// - MHCalibrationPix::SetEventFrequency(fPulserFrequency)
696//
697void MHCalibrationCam::InitHists(MHCalibrationPix &hist, MBadPixelsPix &bad, const Int_t i)
698{
699
700 if (bad.IsUnsuitable(MBadPixelsPix::kUnsuitableRun))
701 hist.SetExcluded();
702
703 hist.InitBins();
704 hist.ChangeHistId(i);
705 hist.SetEventFrequency(fPulserFrequency);
706
707 TH1F *h = hist.GetHGausHist();
708 h->SetTitle( Form("%s%s", h->GetTitle()," Runs: "));
709}
710
711void MHCalibrationCam::FitHiGainArrays(MCalibrationCam &calcam, MBadPixelsCam &badcam,
712 MBadPixelsPix::UncalibratedType_t fittyp,
713 MBadPixelsPix::UncalibratedType_t osctyp)
714{
715
716 for (Int_t i=0; i<fHiGainArray->GetSize(); i++)
717 {
718
719 MHCalibrationPix &hist = (*this)[i];
720
721 if (hist.IsExcluded())
722 continue;
723
724 MCalibrationPix &pix = calcam[i];
725 MBadPixelsPix &bad = badcam[i];
726
727 FitHiGainHists(hist,pix,bad,fittyp,osctyp);
728
729 }
730
731 for (Int_t j=0; j<fAverageHiGainAreas->GetSize(); j++)
732 {
733
734 MHCalibrationPix &hist = GetAverageHiGainArea(j);
735 MCalibrationPix &pix = calcam.GetAverageArea(j);
736 MBadPixelsPix &bad = calcam.GetAverageBadArea(j);
737
738 FitHiGainHists(hist,pix,bad,fittyp,osctyp);
739 }
740
741
742 for (Int_t j=0; j<fAverageHiGainSectors->GetSize(); j++)
743 {
744
745 MHCalibrationPix &hist = GetAverageHiGainSector(j);
746 MCalibrationPix &pix = calcam.GetAverageSector(j);
747 MBadPixelsPix &bad = calcam.GetAverageBadSector(j);
748
749 FitHiGainHists(hist,pix,bad,fittyp,osctyp);
750 }
751
752}
753
754void MHCalibrationCam::FitLoGainArrays(MCalibrationCam &calcam, MBadPixelsCam &badcam,
755 MBadPixelsPix::UncalibratedType_t fittyp,
756 MBadPixelsPix::UncalibratedType_t osctyp)
757{
758
759 for (Int_t i=0; i<fLoGainArray->GetSize(); i++)
760 {
761
762 MHCalibrationPix &hist = (*this)(i);
763
764 if (hist.IsExcluded())
765 continue;
766
767 MCalibrationPix &pix = calcam[i];
768 MBadPixelsPix &bad = badcam[i];
769
770 FitLoGainHists(hist,pix,bad,fittyp,osctyp);
771
772 }
773
774 for (Int_t j=0; j<fAverageLoGainAreas->GetSize(); j++)
775 {
776
777 MHCalibrationPix &hist = GetAverageLoGainArea(j);
778 MCalibrationPix &pix = calcam.GetAverageArea(j);
779 MBadPixelsPix &bad = calcam.GetAverageBadArea(j);
780
781 FitLoGainHists(hist,pix,bad,fittyp,osctyp);
782 }
783
784
785 for (Int_t j=0; j<fAverageLoGainSectors->GetSize(); j++)
786 {
787
788 MHCalibrationPix &hist = GetAverageLoGainSector(j);
789 MCalibrationPix &pix = calcam.GetAverageSector(j);
790 MBadPixelsPix &bad = calcam.GetAverageBadSector(j);
791
792 FitLoGainHists(hist,pix,bad,fittyp,osctyp);
793 }
794}
795
796//------------------------------------------------------------
797//
798// For all averaged areas, the fitted sigma is multiplied with the square root of
799// the number involved pixels
800//
801void MHCalibrationCam::CalcAverageSigma()
802{
803
804 if (!fCam)
805 return;
806
807 for (UInt_t j=0; j<fGeom->GetNumAreas(); j++)
808 {
809
810 MCalibrationPix &pix = fCam->GetAverageArea(j);
811
812 const Float_t numsqr = TMath::Sqrt((Float_t)fAverageAreaNum[j]);
813 fAverageAreaSigma[j] = pix.GetSigma () * numsqr;
814 fAverageAreaSigmaVar[j] = pix.GetSigmaErr () * pix.GetSigmaErr() * numsqr;
815
816 pix.SetSigma (fAverageAreaSigma[j]);
817 pix.SetSigmaVar(fAverageAreaSigmaVar[j]);
818
819 fAverageAreaRelSigma [j] = fAverageAreaSigma[j] / pix.GetMean();
820 fAverageAreaRelSigmaVar[j] = fAverageAreaSigmaVar[j] / (fAverageAreaSigma[j]*fAverageAreaSigma[j]);
821 fAverageAreaRelSigmaVar[j] += pix.GetMeanRelVar();
822 fAverageAreaRelSigmaVar[j] *= fAverageAreaRelSigma[j];
823 }
824}
825
826// ---------------------------------------------------------------------------
827//
828// Returns if the histogram is empty and sets the following flag:
829// - MBadPixelsPix::SetUnsuitable(MBadPixelsPix::kUnsuitableRun)
830//
831// Fits the histograms with a Gaussian, in case of failure
832// calls MHCalibrationPix::RepeatFit(), in case of repeated failure
833// calls MHCalibrationPix::BypassFit() and sets the following flags:
834// - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::UncalibratedType_t fittyp )
835// - MBadPixelsPix::SetUnsuitable( MBadPixelsPix::kUnreliableRun )
836//
837// Creates the fourier spectrum and tests MHGausEvents::IsFourierSpectrumOK().
838// In case no, sets the following flags:
839// - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::UncalibratedType_t osctyp )
840// - MBadPixelsPix::SetUnsuitable( MBadPixelsPix::kUnreliableRun )
841//
842// Retrieves the results and store them in MCalibrationPix
843//
844void MHCalibrationCam::FitHiGainHists(MHCalibrationPix &hist,
845 MCalibrationPix &pix,
846 MBadPixelsPix &bad,
847 MBadPixelsPix::UncalibratedType_t fittyp,
848 MBadPixelsPix::UncalibratedType_t osctyp)
849{
850
851
852 if (hist.IsEmpty() || hist.IsOnlyOverflow() || hist.IsOnlyUnderflow())
853 return;
854
855 //
856 // 2) Fit the Hi Gain histograms with a Gaussian
857 //
858 if (!hist.FitGaus())
859 //
860 // 3) In case of failure set the bit Fitted to false and take histogram means and RMS
861 //
862 if (!hist.RepeatFit())
863 {
864 hist.BypassFit();
865 bad.SetUncalibrated( fittyp );
866 }
867
868 hist.Renorm();
869 //
870 // 4) Check for oscillations
871 //
872 if (IsOscillations())
873 {
874 hist.CreateFourierSpectrum();
875
876 if (!hist.IsFourierSpectrumOK())
877 bad.SetUncalibrated( osctyp );
878 }
879
880 //
881 // 5) Retrieve the results and store them in this class
882 //
883 pix.SetHiGainMean ( hist.GetMean() );
884 pix.SetHiGainMeanVar ( hist.GetMeanErr() * hist.GetMeanErr() );
885 pix.SetHiGainSigma ( hist.GetSigma() );
886 pix.SetHiGainSigmaVar ( hist.GetSigmaErr()* hist.GetSigmaErr() );
887 pix.SetHiGainProb ( hist.GetProb() );
888 pix.SetHiGainNumBlackout( hist.GetBlackout() );
889 pix.SetHiGainNumPickup ( hist.GetPickup() );
890
891 if (IsDebug())
892 {
893 *fLog << dbginf << GetDescriptor() << ": ID " << hist.GetPixId()
894 << " HiGainSaturation: " << pix.IsHiGainSaturation()
895 << " HiGainMean: " << hist.GetMean ()
896 << " HiGainMeanErr: " << hist.GetMeanErr ()
897 << " HiGainMeanSigma: " << hist.GetSigma ()
898 << " HiGainMeanSigmaErr: " << hist.GetSigmaErr()
899 << " HiGainMeanProb: " << hist.GetProb ()
900 << " HiGainNumBlackout: " << hist.GetBlackout()
901 << " HiGainNumPickup : " << hist.GetPickup ()
902 << endl;
903 }
904
905}
906
907
908// ---------------------------------------------------------------------------
909//
910// Returns if the histogram is empty and sets the following flag:
911// - MBadPixelsPix::SetUnsuitable(MBadPixelsPix::kUnsuitableRun)
912//
913// Fits the histograms with a Gaussian, in case of failure
914// calls MHCalibrationPix::RepeatFit(), in case of repeated failure
915// calls MHCalibrationPix::BypassFit() and sets the following flags:
916// - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::UncalibratedType_t fittyp )
917// - MBadPixelsPix::SetUnsuitable( MBadPixelsPix::kUnreliableRun )
918//
919// Creates the fourier spectrum and tests MHGausEvents::IsFourierSpectrumOK().
920// In case no, sets the following flags:
921// - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::UncalibratedType_t osctyp )
922// - MBadPixelsPix::SetUnsuitable( MBadPixelsPix::kUnreliableRun )
923//
924// Retrieves the results and store them in MCalibrationPix
925//
926void MHCalibrationCam::FitLoGainHists(MHCalibrationPix &hist,
927 MCalibrationPix &pix,
928 MBadPixelsPix &bad,
929 MBadPixelsPix::UncalibratedType_t fittyp,
930 MBadPixelsPix::UncalibratedType_t osctyp)
931{
932
933 if (hist.IsEmpty() || hist.IsOnlyOverflow() || hist.IsOnlyUnderflow())
934 return;
935
936
937 //
938 // 2) Fit the Hi Gain histograms with a Gaussian
939 //
940 if (!hist.FitGaus())
941 //
942 // 3) In case of failure set the bit Fitted to false and take histogram means and RMS
943 //
944 if (!hist.RepeatFit())
945 {
946 hist.BypassFit();
947 bad.SetUncalibrated( fittyp );
948 }
949
950 //
951 // 4) Check for oscillations
952 //
953 if (IsOscillations())
954 {
955 hist.CreateFourierSpectrum();
956
957 if (!hist.IsFourierSpectrumOK())
958 bad.SetUncalibrated( osctyp );
959 }
960
961 //
962 // 5) Retrieve the results and store them in this class
963 //
964 pix.SetLoGainMean ( hist.GetMean() );
965 pix.SetLoGainMeanVar ( hist.GetMeanErr() * hist.GetMeanErr() );
966 pix.SetLoGainSigma ( hist.GetSigma() );
967 pix.SetLoGainSigmaVar ( hist.GetSigmaErr() * hist.GetSigmaErr() );
968 pix.SetLoGainProb ( hist.GetProb() );
969 pix.SetLoGainNumBlackout( hist.GetBlackout() );
970 pix.SetLoGainNumPickup ( hist.GetPickup() );
971
972 if (IsDebug())
973 {
974 *fLog << dbginf << GetDescriptor() << "ID: " << hist.GetPixId()
975 << " HiGainSaturation: " << pix.IsHiGainSaturation()
976 << " LoGainMean: " << hist.GetMean ()
977 << " LoGainMeanErr: " << hist.GetMeanErr ()
978 << " LoGainMeanSigma: " << hist.GetSigma ()
979 << " LoGainMeanSigmaErr: " << hist.GetSigmaErr()
980 << " LoGainMeanProb: " << hist.GetProb ()
981 << " LoGainNumBlackout: " << hist.GetBlackout()
982 << " LoGainNumPickup : " << hist.GetPickup ()
983 << endl;
984 }
985
986}
987
988
989
990// --------------------------------------------------------------------------
991//
992// Dummy, needed by MCamEvent
993//
994Bool_t MHCalibrationCam::GetPixelContent(Double_t &val, Int_t idx, const MGeomCam &cam, Int_t type) const
995{
996 return kTRUE;
997}
998
999// --------------------------------------------------------------------------
1000//
1001// What MHCamera needs in order to draw an individual pixel in the camera
1002//
1003void MHCalibrationCam::DrawPixelContent(Int_t idx) const
1004{
1005}
1006
1007// -----------------------------------------------------------------------------
1008//
1009// Default draw:
1010//
1011// Displays the averaged areas, both High Gain and Low Gain
1012//
1013// Calls the Draw of the fAverageHiGainAreas and fAverageLoGainAreas objects with options
1014//
1015void MHCalibrationCam::Draw(const Option_t *opt)
1016{
1017
1018 const Int_t nareas = fAverageHiGainAreas->GetEntries();
1019 if (nareas == 0)
1020 return;
1021
1022 TVirtualPad *pad = gPad ? gPad : MH::MakeDefCanvas(this);
1023 pad->SetBorderMode(0);
1024
1025 pad->Divide(IsLoGain() ? 2 : 1,nareas);
1026
1027 for (Int_t i=0; i<nareas;i++)
1028 {
1029
1030 pad->cd(IsLoGain() ? 2*(i+1)-1 : i+1);
1031 GetAverageHiGainArea(i).Draw(opt);
1032
1033 if (!fAverageAreaSat[i])
1034 DrawAverageSigma(fAverageAreaSat[i], i,
1035 fAverageAreaSigma[i], fAverageAreaSigmaVar[i],
1036 fAverageAreaRelSigma[i], fAverageAreaRelSigmaVar[i]);
1037
1038 if (IsLoGain())
1039 {
1040 pad->cd(2*(i+1));
1041 GetAverageLoGainArea(i).Draw(opt);
1042 }
1043
1044 if (fAverageAreaSat[i])
1045 DrawAverageSigma(fAverageAreaSat[i], i,
1046 fAverageAreaSigma[i], fAverageAreaSigmaVar[i],
1047 fAverageAreaRelSigma[i], fAverageAreaRelSigmaVar[i]);
1048 }
1049
1050}
1051
1052// -----------------------------------------------------------------------------
1053//
1054// Default draw:
1055//
1056// Displays a TPaveText with the re-normalized sigmas of the average area
1057//
1058void MHCalibrationCam::DrawAverageSigma(Bool_t sat, Bool_t inner,
1059 Float_t sigma, Float_t sigmavar,
1060 Float_t relsigma, Float_t relsigmavar) const
1061{
1062
1063 if (sigma != 0 && sigmavar >= 0 && relsigmavar >= 0.)
1064 {
1065
1066 TPad *newpad = new TPad("newpad","transparent",0,0,1,1);
1067 newpad->SetFillStyle(4000);
1068 newpad->Draw();
1069 newpad->cd();
1070
1071 TPaveText *text = new TPaveText(sat? 0.1 : 0.35,0.7,sat ? 0.4 : 0.7,1.0);
1072 text->SetTextSize(0.07);
1073 const TString line1 = Form("%s%s%s",inner ? "Outer" : "Inner",
1074 " Pixels ", sat ? "Low Gain" : "High Gain");
1075 TText *txt1 = text->AddText(line1.Data());
1076 const TString line2 = Form("#sigma per pix: %2.2f #pm %2.2f",sigma,TMath::Sqrt(sigmavar));
1077 TText *txt2 = text->AddText(line2.Data());
1078 const TString line3 = Form("Rel. #sigma per pix: %2.2f #pm %2.2f",relsigma,TMath::Sqrt(relsigmavar));
1079 TText *txt3 = text->AddText(line3.Data());
1080 text->Draw("");
1081
1082 text->SetBit(kCanDelete);
1083 txt1->SetBit(kCanDelete);
1084 txt2->SetBit(kCanDelete);
1085 txt3->SetBit(kCanDelete);
1086 newpad->SetBit(kCanDelete);
1087 }
1088}
1089
1090Int_t MHCalibrationCam::ReadEnv(const TEnv &env, TString prefix, Bool_t print)
1091{
1092 Bool_t rc = kFALSE;
1093 if (IsEnvDefined(env, prefix, "Debug", print))
1094 {
1095 SetDebug(GetEnvValue(env, prefix, "Debug", IsDebug()));
1096 rc = kTRUE;
1097 }
1098 if (IsEnvDefined(env, prefix, "LoGain", print))
1099 {
1100 SetDebug(GetEnvValue(env, prefix, "LoGain", IsLoGain()));
1101 rc = kTRUE;
1102 }
1103 if (IsEnvDefined(env, prefix, "Oscillations", print))
1104 {
1105 SetDebug(GetEnvValue(env, prefix, "Oscillations", IsOscillations()));
1106 rc = kTRUE;
1107 }
1108
1109
1110 return rc;
1111}
Note: See TracBrowser for help on using the repository browser.