source: branches/MarsMoreSimulationTruth/mcalib/MCalibrationChargeCam.cc

Last change on this file was 10166, checked in by tbretz, 14 years ago
Removed the old obsolete cvs header line.
File size: 37.1 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 11/2003 <mailto:markus@ifae.es>
19!
20! Copyright: MAGIC Software Development, 2000-2008
21!
22!
23\* ======================================================================== */
24/////////////////////////////////////////////////////////////////////////////
25//
26// MCalibrationChargeCam
27//
28// Storage container for charge calibration results from the signal distribution
29// fits (see MHCalibrationChargeCam and MHCalibrationChargePix), the calculation
30// of reduced sigmas and number of photo-electrons (this class) and conversion
31// factors sum FADC slices to photo-electrons (see MCalibrationChargeCalc)
32//
33// Individual pixels have to be cast when retrieved e.g.:
34// MCalibrationChargePix &avpix = (MCalibrationChargePix&)(*fChargeCam)[i]
35//
36// Averaged values over one whole area index (e.g. inner or outer pixels for
37// the MAGIC camera), can be retrieved via:
38// MCalibrationChargePix &avpix = (MCalibrationChargePix&)fChargeCam->GetAverageArea(i)
39//
40// Averaged values over one whole camera sector can be retrieved via:
41// MCalibrationChargePix &avpix = (MCalibrationChargePix&)fChargeCam->GetAverageSector(i)
42//
43// Note the averageing has been done on an event-by-event basis. Resulting
44// Sigma's of the Gauss fit have been multiplied with the square root of the number
45// of involved pixels in order to make a direct comparison possible with the mean of
46// sigmas.
47//
48// Final numbers of uncalibrated or unreliable pixels can be retrieved via the commands:
49// GetNumUncalibrated(aidx) and GetNumUnreliable(aidx) where aidx is the area index (0 for
50// inner and 1 for outer pixels in the MAGIC camera).
51//
52// The following "calibration" constants are used for the calibration of each pixel
53// (see MCalibrate):
54//
55// - MCalibrationQEPix::GetMeanConvFADC2Phe(): The mean conversion factor from the
56// summed FADC slices to the number of photo-electrons (in first order independent
57// on colour and intensity)
58// - MCalibrationQEPix::GetMeanFFactorFADC2Phot(): The mean F-Factor of the total
59// readout chain dividing the signal to noise of the incoming number of photons
60// (= sqrt(number photons)) by the signal to noise of the outgoing summed FADC slices
61// signal (= MCalibrationPix::GetMean() / MCalibrationChargePix::GetRSigma() )
62//
63// The following calibration constants can be retrieved directly from this class:
64//
65// - GetConversionFactorFFactor ( Int_t idx, Float_t &mean, Float_t &err, Float_t &sigma );
66//
67// where:
68// - idx is the pixel software ID
69// - "mean" is the mean conversion constant, to be multiplied with the retrieved signal
70// in order to get a calibrated number of PHOTONS.
71// - "err" is the pure statistical uncertainty about the MEAN
72// - "sigma", if mulitplied with the square root of signal, gives the approximate sigma of the
73// retrieved mean number of incident Cherenkov photons.
74//
75// Note, Conversion is ALWAYS (included the F-Factor method!) from summed FADC slices to PHOTONS.
76//
77// See also: MCalibrationChargePix, MCalibrationChargeCalc, MCalibrationQECam
78// MHCalibrationChargePix, MHCalibrationChargeCam
79// MCalibrationBlindPix, MCalibrationChargePINDiode
80//
81/////////////////////////////////////////////////////////////////////////////
82#include "MCalibrationChargeCam.h"
83#include "MCalibrationChargePix.h"
84
85#include <TMath.h>
86
87#include <TH1.h>
88#include <TF1.h>
89#include <TOrdCollection.h>
90
91#include "MLog.h"
92#include "MLogManip.h"
93
94#include "MGeomCam.h"
95#include "MGeom.h"
96
97#include "MBadPixelsCam.h"
98#include "MBadPixelsPix.h"
99
100#include "MCalibrationQECam.h"
101#include "MCalibrationQEPix.h"
102
103#include "MHCamera.h"
104
105ClassImp(MCalibrationChargeCam);
106
107using namespace std;
108// --------------------------------------------------------------------------
109//
110// Default constructor.
111//
112// Calls:
113// - Clear()
114//
115MCalibrationChargeCam::MCalibrationChargeCam(const char *name, const char *title)
116{
117 fName = name ? name : "MCalibrationChargeCam";
118 fTitle = title ? title : "Storage container for the Calibration Information in the camera";
119
120 Clear();
121}
122
123/*
124// --------------------------------------------------------------------------
125//
126// Creates new MCalibrationCam only for the averaged areas:
127// the rest has to be retrieved directly, e.g. via:
128//
129TObject *MCalibrationChargeCam::Clone(const char *) const
130{
131
132 //
133 // FIXME, this might be done faster and more elegant, by direct copy.
134 //
135 MCalibrationChargeCam *cam = new MCalibrationChargeCam(fName,fTitle);
136
137 cam->fNumUnsuitable = fNumUnsuitable;
138 cam->fNumUnreliable = fNumUnreliable;
139 cam->fNumHiGainFADCSlices = fNumHiGainFADCSlices;
140 cam->fNumLoGainFADCSlices = fNumLoGainFADCSlices;
141 cam->fPulserColor = fPulserColor;
142
143 cam->fFlags = fFlags;
144
145 cam->fNumPhotonsBlindPixelMethod = fNumPhotonsBlindPixelMethod;
146 cam->fNumPhotonsFFactorMethod = fNumPhotonsFFactorMethod;
147 cam->fNumPhotonsPINDiodeMethod = fNumPhotonsPINDiodeMethod;
148 cam->fNumPhotonsBlindPixelMethodErr = fNumPhotonsBlindPixelMethodErr;
149 cam->fNumPhotonsFFactorMethodErr = fNumPhotonsFFactorMethodErr;
150 cam->fNumPhotonsPINDiodeMethodErr = fNumPhotonsPINDiodeMethodErr;
151
152 for (Int_t i=0; i<GetSize(); i++)
153 cam->fPixels->AddAt((*this)[i].Clone(),i);
154
155 for (Int_t i=0; i<GetAverageAreas(); i++)
156 {
157 cam->fAverageAreas->AddAt(GetAverageArea(i).Clone(),i);
158 cam->fAverageBadAreas->AddAt(GetAverageBadArea(i).Clone(),i);
159 }
160 for (Int_t i=0; i<GetAverageSectors(); i++)
161 {
162 cam->fAverageSectors->AddAt(GetAverageSector(i).Clone(),i);
163 cam->fAverageBadSectors->AddAt(GetAverageBadSector(i).Clone(),i);
164 }
165
166 return cam;
167}
168*/
169
170// -------------------------------------------------------------------
171//
172// Add MCalibrationChargePix's in the ranges from - to to fPixels
173//
174void MCalibrationChargeCam::Add(const UInt_t a, const UInt_t b)
175{
176 for (UInt_t i=a; i<b; i++)
177 {
178 fPixels->AddAt(new MCalibrationChargePix,i);
179 (*this)[i].SetPixId(i);
180 }
181}
182
183// -------------------------------------------------------------------
184//
185// Add MCalibrationChargePix's in the ranges from - to to fAverageAreas
186//
187void MCalibrationChargeCam::AddArea(const UInt_t a, const UInt_t b)
188{
189 for (UInt_t i=a; i<b; i++)
190 {
191 fAverageAreas->AddAt(new MCalibrationChargePix,i);
192 GetAverageArea(i).SetPixId(i);
193 }
194}
195
196// -------------------------------------------------------------------
197//
198// Add MCalibrationChargePix's in the ranges from - to to fAverageSectors
199//
200void MCalibrationChargeCam::AddSector(const UInt_t a, const UInt_t b)
201{
202 for (UInt_t i=a; i<b; i++)
203 {
204 fAverageSectors->AddAt(new MCalibrationChargePix,i);
205 GetAverageSector(i).SetPixId(i);
206 }
207}
208
209
210// --------------------------------------
211//
212// Sets all variable to 0.
213// Sets all flags to kFALSE
214// Calls MCalibrationCam::Clear()
215//
216void MCalibrationChargeCam::Clear(Option_t *o)
217{
218 // DO NOT RESET THE HI-/LO-GAIN CONSTANTS
219 SetFFactorMethodValid ( kFALSE );
220
221 fNumPhotonsBlindPixelMethod = 0.;
222 fNumPhotonsFFactorMethod = 0.;
223 fNumPhotonsPINDiodeMethod = 0.;
224 fNumPhotonsBlindPixelMethodErr = 0.;
225 fNumPhotonsFFactorMethodErr = 0.;
226 fNumPhotonsPINDiodeMethodErr = 0.;
227
228 MCalibrationCam::Clear();
229}
230
231// -----------------------------------------------
232//
233// Sets the kFFactorMethodValid bit from outside
234//
235void MCalibrationChargeCam::SetFFactorMethodValid(const Bool_t b)
236{
237 b ? SETBIT(fFlags, kFFactorMethodValid) : CLRBIT(fFlags, kFFactorMethodValid);
238}
239
240
241// --------------------------------------------------------------------------
242//
243// Test bit kFFactorMethodValid
244//
245Bool_t MCalibrationChargeCam::IsFFactorMethodValid() const
246{
247 return TESTBIT(fFlags,kFFactorMethodValid);
248}
249
250// --------------------------------------------------------------------------
251//
252// Copy High-gain vs. low-gain conversion factors from cam to this.
253//
254Bool_t MCalibrationChargeCam::MergeHiLoConversionFactors(const MCalibrationChargeCam &cam) const
255{
256
257 if (GetSize() != cam.GetSize())
258 {
259 *fLog << warn << "Sizes mismatch! Cannot merge high-gain vs. low-gain convertion factors" << endl;
260 return kFALSE;
261 }
262
263 for (Int_t i=0; i<GetSize(); i++)
264 {
265 ((MCalibrationChargePix&)(*this)[i]).SetConversionHiLo (((MCalibrationChargePix&)cam[i]).GetConversionHiLo());
266 ((MCalibrationChargePix&)(*this)[i]).SetConversionHiLoErr(((MCalibrationChargePix&)cam[i]).GetConversionHiLoErr());
267 }
268
269 return kTRUE;
270}
271
272// --------------------------------------------------------------------------
273//
274// Print first the well fitted pixels
275// and then the ones which are not FitValid
276//
277void MCalibrationChargeCam::Print(Option_t *o) const
278{
279
280 *fLog << all << GetDescriptor() << ":" << endl;
281 int id = 0;
282
283 *fLog << all << "Calibrated pixels:" << endl;
284 *fLog << all << endl;
285
286 TIter Next(fPixels);
287 MCalibrationChargePix *pix;
288 while ((pix=(MCalibrationChargePix*)Next()))
289 {
290
291 if (pix->IsFFactorMethodValid())
292 {
293
294 *fLog << all
295 << Form("%s%3i","Pixel: ",pix->GetPixId())
296 << Form("%s%4.2f%s%4.2f"," Ped.Rms: ",pix->GetPedRms(),"+-",pix->GetPedRmsErr())
297 << Form("%s%4.2f%s%4.2f"," Charge: " ,pix->GetConvertedMean(),"+-",pix->GetConvertedSigma())
298 << Form("%s%4.2f%s%4.2f"," Red.Sigma: ",pix->GetConvertedRSigma(),"+-",pix->GetConvertedRSigmaErr())
299 << Form("%s%4.2f%s%4.2f"," Num.Phes: ",pix->GetPheFFactorMethod(),"+-",pix->GetPheFFactorMethodErr())
300 << Form("%s%4.2f%s%4.2f"," Conv.FADC2Phe: ",pix->GetMeanConvFADC2Phe(),"+-",pix->GetMeanConvFADC2PheErr())
301 << " Saturated? :" << pix->IsHiGainSaturation()
302 << Form("%s%4.2f%s%4.2f"," Conv.HiLo: ",pix->GetConversionHiLo(),"+-",pix->GetConversionHiLoErr())
303 << endl;
304 id++;
305 }
306 }
307
308 *fLog << all << id << " pixels" << endl;
309 id = 0;
310
311
312 *fLog << all << endl;
313 *fLog << all << "Previously Excluded pixels:" << endl;
314 *fLog << all << endl;
315
316 id = 0;
317
318 TIter Next4(fPixels);
319 while ((pix=(MCalibrationChargePix*)Next4()))
320 {
321 if (pix->IsExcluded())
322 {
323 *fLog << all << pix->GetPixId() << " ";
324 id++;
325
326 if (!(id % 25))
327 *fLog << endl;
328 }
329 }
330
331 *fLog << all << endl;
332 *fLog << all << "New Excluded pixels:" << endl;
333 *fLog << all << endl;
334
335 TIter Next5(fPixels);
336 while ((pix=(MCalibrationChargePix*)Next5()))
337 {
338 if (!pix->IsFFactorMethodValid() && !pix->IsExcluded())
339 {
340 *fLog << all << pix->GetPixId() << " ";
341 id++;
342
343 if (!(id % 25))
344 *fLog << endl;
345 }
346 }
347
348 *fLog << endl;
349 *fLog << all << id << " Excluded pixels " << endl;
350 *fLog << endl;
351
352 *fLog << all << endl;
353 *fLog << all << "Averaged Areas:" << endl;
354 *fLog << all << endl;
355
356 TIter Next6(fAverageAreas);
357 while ((pix=(MCalibrationChargePix*)Next6()))
358 {
359 *fLog << all << Form("%s%3i","Area Idx: ",pix->GetPixId())
360 << " Ped. Rms: " << pix->GetPedRms() << " +- " << pix->GetPedRmsErr()
361 << " Mean signal: " << pix->GetMean() << " +- " << pix->GetMeanErr()
362 << " Sigma signal: " << pix->GetSigma() << " +- "<< pix->GetSigmaErr()
363 << " Reduced Sigma: " << pix->GetRSigma()
364 << " Nr Phe's: " << pix->GetPheFFactorMethod()
365 << endl;
366 }
367
368 *fLog << all << endl;
369 *fLog << all << "Averaged Sectors:" << endl;
370 *fLog << all << endl;
371
372 TIter Next7(fAverageSectors);
373 while ((pix=(MCalibrationChargePix*)Next7()))
374 {
375 *fLog << all << Form("%s%3i","Sector: ",pix->GetPixId())
376 << " Ped. Rms: " << pix->GetPedRms() << " +- " << pix->GetPedRmsErr()
377 << " Mean signal: " << pix->GetMean() << " +- " << pix->GetMeanErr()
378 << " Sigma signal: " << pix->GetSigma() << " +- "<< pix->GetSigmaErr()
379 << " Reduced Sigma: " << pix->GetRSigma()
380 << " Nr Phe's: " << pix->GetPheFFactorMethod()
381 << endl;
382 }
383 *fLog << all << endl;
384}
385
386
387// --------------------------------------------------------------------------
388//
389// The types are as follows:
390//
391// Fitted values:
392// ==============
393//
394// 0: Fitted Charge (see MCalibrationPix::GetMean())
395// 1: Error of fitted Charge (see MCalibrationPix::GetMeanErr())
396// 2: Sigma of fitted Charge (see MCalibrationPix::GetSigma())
397// 3: Error of Sigma of fitted Charge (see MCalibrationPix::GetSigmaErr())
398//
399// Useful variables derived from the fit results:
400// =============================================
401//
402// 4: Probability Gauss fit Charge distribution (see MCalibrationPix::GetProb())
403// 5: Reduced Sigma of fitted Charge (see MCalibrationChargePix::GetRSigma())
404// 6: Error Reduced Sigma of fitted Charge (see MCalibrationChargePix::GetRSigmaErr())
405// 7: Reduced Sigma per Charge (see MCalibrationChargePix::GetRSigmaPerCharge())
406// 8: Error of Reduced Sigma per Charge (see MCalibrationChargePix::GetRSigmaPerChargeErr())
407//
408// Results of the F-Factor calibration Method:
409// ===========================================
410//
411// 9: Nr. Photo-electrons from F-Factor Method (see MCalibrationChargePix::GetPheFFactorMethod())
412// 10: Error Nr. Photo-el. from F-Factor Method (see MCalibrationChargePix::GetPheFFactorMethodErr())
413// 11: Conversion factor from F-Factor Method (see MCalibrationChargePix::GetMeanConvFADC2Phe()
414// 12: Error conv. factor from F-Factor Method (see MCalibrationChargePix::GetMeanConvFADC2PheErr()
415// 13: Overall F-Factor from F-Factor Method (see MCalibrationChargePix::GetMeanFFactorFADC2Phot()
416// 14: Error F-Factor from F-Factor Method (see MCalibrationChargePix::GetMeanFFactorFADC2PhotErr()
417// 15: Pixels valid calibration F-Factor-Method (see MCalibrationChargePix::IsFFactorMethodValid())
418//
419// Results of the Low-Gain vs. High-Gain Conversion:
420// =================================================
421//
422// 16: Mean Signal Hi-Gain / Mean Signal Lo-Gain (see MCalibrationPix::GetHiLoMeansDivided())
423// 17: Error Signal High-Gain / Signal Low Gain (see MCalibrationPix::GetHiLoMeansDividedErr())
424// 18: Sigma High-Gain / Sigma Low Gain (see MCalibrationPix::GetHiLoSigmasDivided())
425// 19: Error Sigma High-Gain / Sigma Low Gain (see MCalibrationPix::GetHiLoSigmasDividedErr())
426//
427// Localized defects:
428// ==================
429//
430// 20: Excluded Pixels
431// 21: Number of pickup events in the Hi Gain (see MCalibrationPix::GetHiGainNumPickup())
432// 22: Number of pickup events in the Lo Gain (see MCalibrationPix::GetLoGainNumPickup())
433// 23: Number of blackout events in the Hi Gain (see MCalibrationPix::GetHiGainNumBlackout())
434// 24: Number of blackout events in the Lo Gain (see MCalibrationPix::GetLoGainNumBlackout())
435//
436// Other classifications of pixels:
437// ================================
438//
439// 25: Pixels with saturated High-Gain (see MCalibrationPix::IsHiGainSaturation())
440//
441// Calculated absolute arrival times (very low precision!):
442// ========================================================
443//
444// 26: Absolute Arrival time of the signal (see MCalibrationChargePix::GetAbsTimeMean())
445// 27: RMS Ab. Arrival time of the signal (see MCalibrationChargePix::GetAbsTimeRms())
446//
447// Used Pedestals:
448// ===============
449//
450// 28: Pedestal for entire signal extr. range (see MCalibrationChargePix::Ped())
451// 29: Error Pedestal entire signal extr. range (see MCalibrationChargePix::PedErr())
452// 30: Ped. RMS entire signal extraction range (see MCalibrationChargePix::PedRms())
453// 31: Error Ped. RMS entire signal extr. range (see MCalibrationChargePix::PedRmsErr())
454//
455// Special variables (for data check):
456// ====================================
457//
458// 32: HiGain RMS divided by Mean for every pixel (with inclusion of the excluded pixels)
459//
460Bool_t MCalibrationChargeCam::GetPixelContent(Double_t &val, Int_t idx, const MGeomCam &cam, Int_t type) const
461{
462 if (idx > GetSize())
463 return kFALSE;
464
465 const Float_t area = cam[idx].GetA();
466 if (area == 0)
467 return kFALSE;
468
469 const MCalibrationChargePix &pix = (MCalibrationChargePix&)(*this)[idx];
470 switch (type)
471 {
472 case 0:
473 if (pix.IsExcluded() || pix.GetConvertedMeanErr()<0)
474 return kFALSE;
475 val = pix.GetConvertedMean();
476 return kTRUE;
477 case 1:
478 if (pix.IsExcluded() || pix.GetConvertedMeanErr()<0)
479 return kFALSE;
480 val = pix.GetConvertedMeanErr();
481 return kTRUE;
482 case 2:
483 if (pix.IsExcluded() || pix.GetConvertedSigmaErr()<0)
484 return kFALSE;
485 val = pix.GetConvertedSigma();
486 return kTRUE;
487 case 3:
488 if (pix.IsExcluded() || pix.GetConvertedSigmaErr()<0)
489 return kFALSE;
490 val = pix.GetConvertedSigmaErr();
491 return kTRUE;
492 case 4:
493 if (pix.IsExcluded())
494 return kFALSE;
495 val = pix.GetProb();
496 return val>=0;
497 case 5:
498 if (!pix.IsFFactorMethodValid())
499 return kFALSE;
500 if (pix.GetRSigma() == -1.)
501 return kFALSE;
502 val = pix.GetConvertedRSigma();
503 break;
504 case 6:
505 if (!pix.IsFFactorMethodValid())
506 return kFALSE;
507 if (pix.GetRSigma() == -1.)
508 return kFALSE;
509 val = pix.GetConvertedRSigmaErr();
510 break;
511 case 7:
512 if (!pix.IsFFactorMethodValid())
513 return kFALSE;
514 val = pix.GetRSigmaPerCharge();
515 break;
516 case 8:
517 if (!pix.IsFFactorMethodValid())
518 return kFALSE;
519 val = pix.GetRSigmaPerChargeErr();
520 break;
521 case 9:
522 // if (!pix.IsFFactorMethodValid())
523 // return kFALSE;
524 if (pix.IsExcluded() || pix.GetPheFFactorMethodErr()<0)
525 return kFALSE;
526 val = pix.GetPheFFactorMethod();
527 return kTRUE;
528 case 10:
529 // if (!pix.IsFFactorMethodValid())
530 // return kFALSE;
531 if (pix.IsExcluded() || pix.GetPheFFactorMethodErr()<=0)
532 return kFALSE;
533 val = pix.GetPheFFactorMethodErr();
534 return kTRUE;
535 case 11:
536 // if (!pix.IsFFactorMethodValid())
537 // return kFALSE;
538 if (pix.IsExcluded() || pix.GetMeanConvFADC2PheErr()<0)
539 return kFALSE;
540 val = pix.GetMeanConvFADC2Phe();
541 return kTRUE;
542 case 12:
543 // if (!pix.IsFFactorMethodValid())
544 // return kFALSE;
545 if (pix.IsExcluded() || pix.GetMeanConvFADC2PheErr()<=0)
546 return kFALSE;
547 val = pix.GetMeanConvFADC2PheErr();
548 return kTRUE;
549 case 13:
550 // if (!pix.IsFFactorMethodValid())
551 // return kFALSE;
552 if (pix.IsExcluded())
553 return kFALSE;
554 val = pix.GetMeanFFactorFADC2Phot();
555 break;
556 case 14:
557 // if (!pix.IsFFactorMethodValid())
558 // return kFALSE;
559 if (pix.IsExcluded())
560 return kFALSE;
561 val = pix.GetMeanFFactorFADC2PhotErr();
562 if (val <= 0.)
563 val = 0.00001;
564 break;
565 case 15:
566 if (pix.IsExcluded())
567 return kFALSE;
568 if (!pix.IsFFactorMethodValid())
569 return kFALSE;
570 val = 1;
571 return kTRUE;
572 case 16:
573 if (pix.IsExcluded())
574 return kFALSE;
575 val = pix.GetHiLoMeansDivided();
576 break;
577 case 17:
578 if (pix.IsExcluded())
579 return kFALSE;
580 val = pix.GetHiLoMeansDividedErr();
581 break;
582 case 18:
583 if (pix.IsExcluded())
584 return kFALSE;
585 val = pix.GetHiLoSigmasDivided();
586 break;
587 case 19:
588 if (pix.IsExcluded())
589 return kFALSE;
590 val = pix.GetHiLoSigmasDividedErr();
591 break;
592 case 20:
593 if (!pix.IsExcluded())
594 return kFALSE;
595 val = 1.;
596 return kTRUE;
597 case 21:
598 if (pix.IsExcluded())
599 return kFALSE;
600 val = pix.GetHiGainNumPickup();
601 break;
602 case 22:
603 if (pix.IsExcluded())
604 return kFALSE;
605 val = pix.GetLoGainNumPickup();
606 break;
607 case 23:
608 if (pix.IsExcluded())
609 return kFALSE;
610 val = pix.GetHiGainNumBlackout();
611 break;
612 case 24:
613 if (pix.IsExcluded())
614 return kFALSE;
615 val = pix.GetLoGainNumBlackout();
616 break;
617 case 25:
618 if (pix.IsExcluded())
619 return kFALSE;
620 val = pix.IsHiGainSaturation();
621 break;
622 case 26:
623 if (pix.IsExcluded())
624 return kFALSE;
625 val = pix.GetAbsTimeMean();
626 break;
627 case 27:
628 if (pix.IsExcluded())
629 return kFALSE;
630 val = pix.GetAbsTimeRms();
631 break;
632 case 28:
633 if (pix.IsExcluded())
634 return kFALSE;
635 val = pix.GetPed();
636 break;
637 case 29:
638 if (pix.IsExcluded())
639 return kFALSE;
640 val = pix.GetPedErr();
641 break;
642 case 30:
643 if (pix.IsExcluded())
644 return kFALSE;
645 val = pix.GetPedRms();
646 break;
647 case 31:
648 if (pix.IsExcluded())
649 return kFALSE;
650 val = pix.GetPedErr()/2.;
651 break;
652 case 32:
653 if (pix.IsExcluded() || pix.GetRms()<0)
654 return kFALSE;
655 val = pix.GetMean() == 0. ? 0. : pix.GetRms()/pix.GetMean();
656 return kTRUE;
657 case 33:
658 if (pix.IsExcluded() || pix.GetRms()<0)
659 return kFALSE;
660 if (pix.GetMean() == 0.)
661 val = 0.;
662 else
663 val = pix.GetSigmaErr()/pix.GetMean() + pix.GetRms()*pix.GetMeanErr()/pix.GetMean()/pix.GetMean();
664 return kTRUE;
665 default:
666 return kFALSE;
667 }
668
669 return val!=-1.;
670}
671
672Bool_t MCalibrationChargeCam::GetConversionFactorFFactor(Int_t ipx, Float_t &mean, Float_t &ferr, Float_t &ffactor)
673{
674
675 MCalibrationChargePix &pix = (MCalibrationChargePix&)(*this)[ipx];
676
677 Float_t conv = pix.GetMeanConvFADC2Phe();
678
679 if (conv < 0.)
680 return kFALSE;
681
682 mean = conv;
683 ferr = pix.GetMeanConvFADC2PheErr();
684 ffactor = pix.GetMeanFFactorFADC2Phot();
685
686 return kTRUE;
687}
688
689
690// --------------------------------------------------------------------------
691//
692// Calculates the average conversion factor FADC counts to photons for pixel sizes.
693// The geometry container is used to get the necessary
694// geometry information (area index). The MCalibrationQECam container is necessary for
695// the quantum efficiency information.
696// If the bad pixel container is given all pixels which have the flag 'kUnsuitableRun' are ignored
697// in the calculation of the size average.
698//
699// Returns a TArrayF of dimension 2:
700// arr[0]: averaged conversion factors (default: -1.)
701// arr[1]: Error (rms) of averaged conversion factors (default: 0.)
702//
703TArrayF MCalibrationChargeCam::GetAveragedConvFADC2PhotPerArea ( const MGeomCam &geom, const MCalibrationQECam &qecam,
704 const UInt_t ai, MBadPixelsCam *bad)
705{
706
707 const Int_t np = GetSize();
708
709 Double_t mean = 0.;
710 Double_t mean2 = 0.;
711 Int_t nr = 0;
712
713 MHCamera convcam(geom,"ConvFactors","Conversion Factors;Conv Factor [phot/FADC cnts];channels");
714
715 for (int i=0; i<np; i++)
716 {
717 if (bad && (*bad)[i].IsUnsuitable(MBadPixelsPix::kUnsuitableRun))
718 continue;
719
720 if (bad && (*bad)[i].IsUncalibrated(MBadPixelsPix::kDeviatingNumPhes))
721 continue;
722
723 const UInt_t aidx = geom[i].GetAidx();
724
725 if (ai != aidx)
726 continue;
727
728 const MCalibrationQEPix &qepix = (MCalibrationQEPix&)qecam[i];
729 if (!qepix.IsAverageQEFFactorAvailable())
730 continue;
731
732 const MCalibrationChargePix &pix = (MCalibrationChargePix&)(*this)[i];
733 const Float_t conv = pix.GetMeanConvFADC2Phe()/qepix.GetQECascadesFFactor();
734
735 mean += conv;
736 mean2 += conv*conv;
737 nr ++;
738
739 convcam.Fill(i,conv);
740 convcam.SetUsed(i);
741 }
742
743 Float_t mn = nr ? mean/nr : -1.;
744 Float_t sg = nr>1 ? TMath::Sqrt((mean2 - mean*mean/nr)/(nr-1)) : 0.;
745
746 const Int_t aidx = (Int_t)ai;
747
748 TH1D *h = convcam.ProjectionS(TArrayI(),TArrayI(1,&aidx),"_py",750);
749 h->SetDirectory(NULL);
750
751 TF1 *fit = NULL;
752
753 if (geom.InheritsFrom("MGeomCamMagic"))
754 {
755
756 fit = new TF1("fit","gaus",0.4,5.);
757
758 // Fix the ranges, as found by Nadia
759 if(aidx == 0)
760 h->Fit(fit, "REQ0", "",0.4,1.5);
761 else
762 h->Fit(fit ,"REQ0", "",1.,5.);
763 }
764 else
765 {
766 h->Fit("gaus","Q0");
767 fit = h->GetFunction("gaus");
768 }
769
770 Float_t ci2 = fit->GetChisquare();
771 Float_t sigma = fit->GetParameter(2);
772
773 if (ci2 > 500. || sigma > sg)
774 {
775 if (geom.InheritsFrom("MGeomCamMagic"))
776 {
777 // Fix the ranges, as found by Nadia
778 if(aidx == 0)
779 h->Fit(fit, "REQ0", "",0.4,1.5);
780 else
781 h->Fit(fit, "REQ0", "",1.,5.);
782 }
783 else
784 {
785 h->Fit("gaus","MREQ0");
786 fit = h->GetFunction("gaus");
787 }
788
789 ci2 = fit->GetChisquare();
790 sigma = fit->GetParameter(2);
791 }
792
793 const Int_t ndf = fit->GetNDF();
794
795 if (ci2 < 500. && sigma < sg && ndf > 2)
796 {
797 mn = fit->GetParameter(1);
798 sg = sigma;
799 }
800
801 *fLog << inf << "Conversion Factors to photons area idx: " << aidx << ":" << endl;
802 *fLog << inf << "Mean: " << Form("%4.3f",mn)
803 << "+-" << Form("%4.3f",fit->GetParError(1))
804 << " Sigma: " << Form("%4.3f",sg) << "+-" << Form("%4.3f",fit->GetParError(2))
805 << " Chisquare: " << Form("%4.3f",fit->GetChisquare()) << " NDF : " << ndf << endl;
806
807 delete fit;
808 delete h;
809 gROOT->GetListOfFunctions()->Remove(fit);
810
811 TArrayF arr(2);
812 arr[0] = mn;
813 arr[1] = sg;
814
815 return arr;
816}
817
818// --------------------------------------------------------------------------
819//
820// Calculates the average conversion factor FADC counts to equiv. photo-electrons for pixel sizes.
821// The geometry container is used to get the necessary
822// geometry information (area index). The MCalibrationQECam container is necessary for
823// the quantum efficiency information.
824// If the bad pixel container is given all pixels which have the flag 'kUnsuitableRun' are ignored
825// in the calculation of the size average.
826//
827// Returns a TArrayF of dimension 2:
828// arr[0]: averaged conversion factors (default: -1.)
829// arr[1]: Error (rms) of averaged conversion factors (default: 0.)
830//
831TArrayF MCalibrationChargeCam::GetAveragedConvFADC2PhePerArea ( const MGeomCam &geom, const MCalibrationQECam &qecam,
832 const UInt_t ai, MBadPixelsCam *bad)
833{
834
835 const Int_t np = GetSize();
836
837 Double_t mean = 0.;
838 Double_t mean2 = 0.;
839 Int_t nr = 0;
840
841 MHCamera convcam(geom,"ConvFactors","Conversion Factors;Conv Factor [phe/FADC cnts];channels");
842
843 for (int i=0; i<np; i++)
844 {
845 if (bad && (*bad)[i].IsUnsuitable(MBadPixelsPix::kUnsuitableRun))
846 continue;
847
848 if (bad && (*bad)[i].IsUncalibrated(MBadPixelsPix::kDeviatingNumPhes))
849 continue;
850
851 const UInt_t aidx = geom[i].GetAidx();
852
853 if (ai != aidx)
854 continue;
855
856 const MCalibrationQEPix &qepix = (MCalibrationQEPix&)qecam[i];
857 if (!qepix.IsAverageQEFFactorAvailable())
858 continue;
859
860 const MCalibrationChargePix &pix = (MCalibrationChargePix&)(*this)[i];
861 const Float_t conv = pix.GetMeanConvFADC2Phe()/qepix.GetQECascadesFFactor()*MCalibrationQEPix::gkDefaultAverageQE;
862
863 mean += conv;
864 mean2 += conv*conv;
865 nr ++;
866
867 convcam.Fill(i,conv);
868 convcam.SetUsed(i);
869 }
870
871 Float_t mn = nr ? mean/nr : -1.;
872 Float_t sg = nr>1 ? TMath::Sqrt((mean2 - mean*mean/nr)/(nr-1)) : 0.;
873
874 const Int_t aidx = (Int_t)ai;
875
876 TH1D *h = convcam.ProjectionS(TArrayI(),TArrayI(1,&aidx),"_py",750);
877 h->SetDirectory(NULL);
878
879 TF1 *fit = NULL;
880
881 if (geom.InheritsFrom("MGeomCamMagic"))
882 {
883
884 fit = new TF1("fit","gaus",0.01,1.);
885
886 // Fix the ranges, as found by Nadia
887 if(aidx == 0)
888 {h->Fit("fit","REQ0", "",0.07,0.3);}
889 else
890 {h->Fit("fit","REQ0", "",0.15,1.0);}
891 }
892 else
893 {
894 h->Fit("gaus","Q0");
895 fit = h->GetFunction("gaus");
896 }
897
898 Float_t ci2 = fit->GetChisquare();
899 Float_t sigma = fit->GetParameter(2);
900
901 if (ci2 > 500. || sigma > sg)
902 {
903 if (geom.InheritsFrom("MGeomCamMagic"))
904 {
905 // Fix the ranges, as found by Nadia
906 if(aidx == 0)
907 {h->Fit("fit","REQ0", "",0.07,0.3);}
908 else
909 {h->Fit("fit","REQ0", "",0.15,1.0);}
910 }
911 else
912 {
913 h->Fit("gaus","MREQ0");
914 fit = h->GetFunction("gaus");
915 }
916
917 ci2 = fit->GetChisquare();
918 sigma = fit->GetParameter(2);
919 }
920
921 const Int_t ndf = fit->GetNDF();
922
923 if (ci2 < 500. && sigma < sg && ndf > 2)
924 {
925 mn = fit->GetParameter(1);
926 sg = sigma;
927 }
928
929 *fLog << inf << "Conversion Factors to equiv. photo-electrons area idx: " << aidx << ":" << endl;
930 *fLog << inf << "Mean: " << Form("%4.3f",mn)
931 << "+-" << Form("%4.3f",fit->GetParError(1))
932 << " Sigma: " << Form("%4.3f",sg) << "+-" << Form("%4.3f",fit->GetParError(2))
933 << " Chisquare: " << Form("%4.3f",fit->GetChisquare()) << " NDF : " << ndf << endl;
934
935 delete fit;
936 delete h;
937 gROOT->GetListOfFunctions()->Remove(fit);
938
939 TArrayF arr(2);
940 arr[0] = mn;
941 arr[1] = sg;
942
943 return arr;
944}
945
946// --------------------------------------------------------------------------
947//
948// Calculates the average conversion factor FADC counts to photons for camera sectors.
949// The geometry container is used to get the necessary
950// geometry information (area index). The MCalibrationQECam container is necessary for
951// the quantum efficiency information.
952// If the bad pixel container is given all pixels which have the flag 'kUnsuitableRun' are ignored
953// in the calculation of the size average.
954//
955// Returns a TArrayF of dimension 2:
956// arr[0]: averaged conversion factors (default: -1.)
957// arr[1]: Error (rms) of averaged conversion factors (default: 0.)
958//
959TArrayF MCalibrationChargeCam::GetAveragedConvFADC2PhotPerSector( const MGeomCam &geom, const MCalibrationQECam &qecam,
960 const UInt_t sec, MBadPixelsCam *bad)
961{
962 const Int_t np = GetSize();
963
964 Double_t mean = 0.;
965 Double_t mean2 = 0.;
966 Int_t nr = 0;
967
968 for (int i=0; i<np; i++)
969 {
970 if (bad && (*bad)[i].IsUnsuitable(MBadPixelsPix::kUnsuitableRun))
971 continue;
972
973 const UInt_t sector = geom[i].GetSector();
974
975 if (sector != sec)
976 continue;
977
978 const MCalibrationQEPix &qepix = (MCalibrationQEPix&)qecam[i];
979 if (!qepix.IsAverageQEFFactorAvailable())
980 continue;
981
982 const MCalibrationChargePix &pix = (MCalibrationChargePix&)(*this)[i];
983 const Float_t conv = pix.GetMeanConvFADC2Phe();
984 const Float_t qe = qepix.GetQECascadesFFactor();
985
986 mean += conv/qe;
987 mean2 += conv*conv/qe/qe;
988 nr ++;
989
990 }
991
992 TArrayF arr(2);
993 arr[0] = nr ? mean/nr : -1;
994 arr[1] = nr>1 ? TMath::Sqrt((mean2 - mean*mean/nr)/(nr-1)) : 0;
995 return arr;
996}
997
998// --------------------------------------------------------------------------
999//
1000// Calculates the average mean arrival times for pixel sizes.
1001// The geometry container is used to get the necessary
1002// geometry information (area index).
1003// If the bad pixel container is given all pixels which have the flag 'kUnsuitableRun' are ignored
1004// in the calculation of the size average.
1005//
1006// Returns a TArrayF of dimension 2:
1007// arr[0]: averaged mean arrival times (default: -1.)
1008// arr[1]: Error (rms) of averaged mean arrival times (default: 0.)
1009//
1010TArrayF MCalibrationChargeCam::GetAveragedArrivalTimeMeanPerArea(const MGeomCam &geom,
1011 const UInt_t ai, MBadPixelsCam *bad)
1012{
1013
1014 const Int_t np = GetSize();
1015
1016 Double_t mean = 0.;
1017 Double_t mean2 = 0.;
1018 Int_t nr = 0;
1019
1020 for (int i=0; i<np; i++)
1021 {
1022 if (bad && (*bad)[i].IsUnsuitable(MBadPixelsPix::kUnsuitableRun))
1023 continue;
1024
1025 const UInt_t aidx = geom[i].GetAidx();
1026
1027 if (ai != aidx)
1028 continue;
1029
1030 const MCalibrationChargePix &pix = (MCalibrationChargePix&)(*this)[i];
1031 const Float_t time = pix.GetAbsTimeMean();
1032
1033 mean += time ;
1034 mean2 += time*time;
1035 nr ++;
1036
1037 }
1038
1039 TArrayF arr(2);
1040 arr[0] = nr ? mean/nr : -1;
1041 arr[1] = nr>1 ? TMath::Sqrt((mean2 - mean*mean/nr)/(nr-1)) : 0;
1042 return arr;
1043}
1044
1045// --------------------------------------------------------------------------
1046//
1047// Calculates the average mean arrival times for camera sectors.
1048// The geometry container is used to get the necessary
1049// geometry information (area index).
1050// If the bad pixel container is given all pixels which have the flag 'kUnsuitableRun' are ignored
1051// in the calculation of the size average.
1052//
1053// Returns a TArrayF of dimension 2:
1054// arr[0]: averaged mean arrival times (default: -1.)
1055// arr[1]: Error (rms) of averaged mean arrival times (default: 0.)
1056//
1057TArrayF MCalibrationChargeCam::GetAveragedArrivalTimeMeanPerSector(const MGeomCam &geom,
1058 const UInt_t sec, MBadPixelsCam *bad)
1059{
1060 const Int_t np = GetSize();
1061
1062 Double_t mean = 0.;
1063 Double_t mean2 = 0.;
1064 Int_t nr = 0;
1065
1066 for (int i=0; i<np; i++)
1067 {
1068 if (bad && (*bad)[i].IsUnsuitable(MBadPixelsPix::kUnsuitableRun))
1069 continue;
1070
1071 const UInt_t sector = geom[i].GetSector();
1072
1073 if (sector != sec)
1074 continue;
1075
1076 const MCalibrationChargePix &pix = (MCalibrationChargePix&)(*this)[i];
1077 const Float_t time = pix.GetAbsTimeMean();
1078
1079 mean += time;
1080 mean2 += time*time;
1081 nr ++;
1082
1083 }
1084
1085 TArrayF arr(2);
1086 arr[0] = nr ? mean/nr : -1;
1087 arr[1] = nr>1 ? TMath::Sqrt((mean2 - mean*mean/nr)/(nr-1)) : 0;
1088 return arr;
1089}
1090
1091// --------------------------------------------------------------------------
1092//
1093// Calculates the average arrival time RMSs for pixel sizes.
1094// The geometry container is used to get the necessary
1095// geometry information (area index).
1096// If the bad pixel container is given all pixels which have the flag 'kUnsuitableRun' are ignored
1097// in the calculation of the size average.
1098//
1099// Returns a TArrayF of dimension 2:
1100// arr[0]: averaged arrival time RMSs (default: -1.)
1101// arr[1]: Error (rms) of averaged arrival time RMSs (default: 0.)
1102//
1103TArrayF MCalibrationChargeCam::GetAveragedArrivalTimeRmsPerArea ( const MGeomCam &geom,
1104 const UInt_t ai, MBadPixelsCam *bad)
1105{
1106
1107 const Int_t np = GetSize();
1108
1109 Double_t mean = 0.;
1110 Double_t mean2 = 0.;
1111 Int_t nr = 0;
1112
1113 for (int i=0; i<np; i++)
1114 {
1115 if (bad && (*bad)[i].IsUnsuitable(MBadPixelsPix::kUnsuitableRun))
1116 continue;
1117
1118 const UInt_t aidx = geom[i].GetAidx();
1119
1120 if (ai != aidx)
1121 continue;
1122
1123 const MCalibrationChargePix &pix = (MCalibrationChargePix&)(*this)[i];
1124 const Float_t rms = pix.GetAbsTimeRms();
1125
1126 mean += rms;
1127 mean2 += rms*rms;
1128 nr ++;
1129
1130 }
1131
1132 TArrayF arr(2);
1133 arr[0] = nr ? mean/nr : -1;
1134 arr[1] = nr>1 ? TMath::Sqrt((mean2 - mean*mean/nr)/(nr-1)) : 0;
1135 return arr;
1136}
1137
1138// --------------------------------------------------------------------------
1139//
1140// Calculates the average arrival time RMSs for camera sectors.
1141// The geometry container is used to get the necessary
1142// geometry information (area index).
1143// If the bad pixel container is given all pixels which have the flag 'kUnsuitableRun' are ignored
1144// in the calculation of the size average.
1145//
1146// Returns a TArrayF of dimension 2:
1147// arr[0]: averaged arrival time RMSs (default: -1.)
1148// arr[1]: Error (rms) of averaged arrival time RMSs (default: 0.)
1149//
1150TArrayF MCalibrationChargeCam::GetAveragedArrivalTimeRmsPerSector( const MGeomCam &geom, const UInt_t sec, MBadPixelsCam *bad)
1151{
1152 const Int_t np = GetSize();
1153
1154 Double_t mean = 0.;
1155 Double_t mean2 = 0.;
1156 Int_t nr = 0;
1157
1158 for (int i=0; i<np; i++)
1159 {
1160 if (bad && (*bad)[i].IsUnsuitable(MBadPixelsPix::kUnsuitableRun))
1161 continue;
1162
1163 const UInt_t sector = geom[i].GetSector();
1164
1165 if (sector != sec)
1166 continue;
1167
1168 const MCalibrationChargePix &pix = (MCalibrationChargePix&)(*this)[i];
1169 const Float_t rms = pix.GetAbsTimeRms();
1170
1171 mean += rms;
1172 mean2 += rms*rms;
1173 nr ++;
1174 }
1175
1176 TArrayF arr(2);
1177 arr[0] = nr ? mean/nr : -1;
1178 arr[1] = nr>1 ? TMath::Sqrt((mean2 - mean*mean/nr)/(nr-1)) : 0;
1179
1180 return arr;
1181}
Note: See TracBrowser for help on using the repository browser.