source: trunk/MagicSoft/Mars/mcalib/MCalibrationChargePix.cc@ 3394

Last change on this file since 3394 was 3356, checked in by gaug, 21 years ago
*** empty log message ***
File size: 27.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// //
27// MCalibrationChargePix //
28// //
29// Storage container to hold informations about the calibration values //
30// values of one Pixel (PMT). //
31// //
32// The following values are initialized to meaningful values:
33//
34// - The Electronic Rms to 1.5 per FADC slice
35// - The uncertainty about the Electronic RMS to 0.3 per slice
36// - The F-Factor is assumed to have been measured in Munich to 1.13 - 1.17.
37// with the Munich definition of the F-Factor, thus:
38// F = Sigma(Out)/Mean(Out) * Mean(In)/Sigma(In)
39// Mean F-Factor = 1.15
40// Error F-Factor = 0.02
41//
42// - Average QE: (email David Paneque, 14.2.04):
43//
44// The conversion factor that comes purely from QE folded to a Cherenkov
45// spectrum has to be multiplied by:
46// * Plexiglass window -->> 0.96 X 0.96
47// * PMT photoelectron collection efficiency -->> 0.9
48// * Light guides efficiency -->> 0.94
49//
50// Concerning the light guides effiency estimation... Daniel Ferenc
51// is preparing some work (simulations) to estimate it. Yet so far, he has
52// been busy with other stuff, and this work is still UNfinished.
53//
54// The estimation I did comes from:
55// 1) Reflectivity of light guide walls is 85 % (aluminum)
56// 2) At ZERO degree light incidence, 37% of the light hits such walls
57// (0.15X37%= 5.6% of light lost)
58// 3) When increasing the light incidence angle, more and more light hits
59// the walls.
60//
61// However, the loses due to larger amount of photons hitting the walls is more
62// or less counteracted by the fact that more and more photon trajectories cross
63// the PMT photocathode twice, increasing the effective sensitivity of the PMT.
64//
65// Jurgen Gebauer did some quick measurements about this issue. I attach a
66// plot. You can see that the angular dependence is (more or less) in agreement
67// with a CosTheta function (below 20-25 degrees),
68// which is the variation of teh entrance window cross section. So, in
69// first approximation, no loses when increasing light incidence angle;
70// and therefore, the factor 0.94.
71//
72// So, summarizing... I would propose the following conversion factors
73// (while working with CT1 cal box) in order to get the final number of photons
74// from the detected measured size in ADC counts.
75//
76// Nph = ADC * FmethodConversionFactor / ConvPhe-PhFactor
77//
78// FmethodConversionFactor ; measured for individual pmts
79//
80// ConvPhe-PhFactor = 0.98 * 0.23 * 0.90 * 0.94 * 0.96 * 0.96 = 0.18
81//
82// I would not apply any smearing of this factor (which we have in nature),
83// since we might be applying it to PMTs in the totally wrong direction.
84//
85//
86/////////////////////////////////////////////////////////////////////////////
87#include "MCalibrationChargePix.h"
88
89#include "MLog.h"
90#include "MLogManip.h"
91
92ClassImp(MCalibrationChargePix);
93
94using namespace std;
95
96const Float_t MCalibrationChargePix::gkElectronicPedRms = 1.5;
97const Float_t MCalibrationChargePix::gkElectronicPedRmsErr = 0.3;
98const Float_t MCalibrationChargePix::gkFFactor = 1.15;
99const Float_t MCalibrationChargePix::gkFFactorErr = 0.02;
100
101const Float_t MCalibrationChargePix::gkAverageQE = 0.18;
102const Float_t MCalibrationChargePix::gkAverageQEErr = 0.02;
103const Float_t MCalibrationChargePix::gkConversionHiLo = 10.;
104const Float_t MCalibrationChargePix::gkConversionHiLoErr = 2.5;
105
106const Float_t MCalibrationChargePix::fgChargeLimit = 3.;
107const Float_t MCalibrationChargePix::fgChargeErrLimit = 0.;
108const Float_t MCalibrationChargePix::fgChargeRelErrLimit = 1.;
109const Float_t MCalibrationChargePix::fgTimeLimit = 1.5;
110const Float_t MCalibrationChargePix::fgTimeErrLimit = 3.;
111const Float_t MCalibrationChargePix::fgConvFFactorRelErrLimit = 0.25;
112// --------------------------------------------------------------------------
113//
114// Default Constructor:
115//
116MCalibrationChargePix::MCalibrationChargePix(const char *name, const char *title)
117 : fPixId(-1),
118 fFlags(0)
119{
120
121 fName = name ? name : "MCalibrationChargePix";
122 fTitle = title ? title : "Container of the fit results of MHCalibrationChargePixs ";
123
124 Clear();
125
126 //
127 // At the moment, we don't have a database, yet,
128 // so we get it from the configuration file
129 //
130 SetConversionHiLo();
131 SetConversionHiLoErr();
132
133 SetAverageQE();
134 SetChargeLimit();
135 SetChargeErrLimit();
136
137 SetChargeRelErrLimit();
138 SetTimeLimit();
139 SetTimeErrLimit();
140 SetConvFFactorRelErrLimit();
141}
142
143// ------------------------------------------------------------------------
144//
145// Invalidate values
146//
147void MCalibrationChargePix::Clear(Option_t *o)
148{
149
150 SetHiGainSaturation ( kFALSE );
151 SetLoGainSaturation ( kFALSE );
152 SetHiGainFitted ( kFALSE );
153 SetLoGainFitted ( kFALSE );
154 SetHiGainOscillating ( kFALSE );
155 SetLoGainOscillating ( kFALSE );
156 SetExcluded ( kFALSE );
157 SetChargeValid ( kFALSE );
158 SetBlindPixelMethodValid ( kFALSE );
159 SetFFactorMethodValid ( kFALSE );
160 SetPINDiodeMethodValid ( kFALSE );
161 SetCombinedMethodValid ( kFALSE );
162
163 SetMeanTimeInFirstBin ( kFALSE );
164 SetMeanTimeInLastBin ( kFALSE );
165
166 fHiGainMeanCharge = -1.;
167 fHiGainMeanChargeErr = -1.;
168 fHiGainSigmaCharge = -1.;
169 fHiGainSigmaChargeErr = -1.;
170 fHiGainChargeProb = -1.;
171
172 fLoGainMeanCharge = -1.;
173 fLoGainMeanChargeErr = -1.;
174 fLoGainSigmaCharge = -1.;
175 fLoGainSigmaChargeErr = -1.;
176 fLoGainChargeProb = -1.;
177
178 fRSigmaCharge = -1.;
179 fRSigmaChargeErr = -1.;
180
181 fHiGainNumPickup = -1;
182 fLoGainNumPickup = -1;
183
184 fNumLoGainSamples = -1.;
185
186 fPed = -1.;
187 fPedRms = -1.;
188 fPedErr = -1.;
189
190 fLoGainPedRms = -1.;
191 fLoGainPedRmsErr = -1.;
192
193 fTimeFirstHiGain = 0 ;
194 fTimeLastHiGain = 0 ;
195 fTimeFirstLoGain = 0 ;
196 fTimeLastLoGain = 0 ;
197
198 fAbsTimeMean = -1.;
199 fAbsTimeRms = -1.;
200
201 fPheFFactorMethod = -1.;
202 fPheFFactorMethodErr = -1.;
203
204 fMeanConversionFFactorMethod = -1.;
205 fMeanConversionBlindPixelMethod = -1.;
206 fMeanConversionPINDiodeMethod = -1.;
207 fMeanConversionCombinedMethod = -1.;
208
209 fConversionFFactorMethodErr = -1.;
210 fConversionBlindPixelMethodErr = -1.;
211 fConversionPINDiodeMethodErr = -1.;
212 fConversionCombinedMethodErr = -1.;
213
214 fSigmaConversionFFactorMethod = -1.;
215 fSigmaConversionBlindPixelMethod = -1.;
216 fSigmaConversionPINDiodeMethod = -1.;
217 fSigmaConversionCombinedMethod = -1.;
218
219 fTotalFFactorFFactorMethod = -1.;
220 fTotalFFactorBlindPixelMethod = -1.;
221 fTotalFFactorPINDiodeMethod = -1.;
222 fTotalFFactorCombinedMethod = -1.;
223
224}
225
226
227void MCalibrationChargePix::DefinePixId(Int_t i)
228{
229 fPixId = i;
230}
231
232
233// --------------------------------------------------------------------------
234//
235// Set the pedestals from outside
236//
237void MCalibrationChargePix::SetPedestal(const Float_t ped, const Float_t pedrms, const Float_t pederr)
238{
239
240 fPed = ped;
241 fPedRms = pedrms;
242 fPedErr = pederr;
243}
244
245
246void MCalibrationChargePix::SetMeanCharge( const Float_t f )
247{
248 if (IsHiGainSaturation())
249 fLoGainMeanCharge = f;
250 else
251 fHiGainMeanCharge = f;
252}
253
254void MCalibrationChargePix::SetMeanChargeErr( const Float_t f )
255{
256 if (IsHiGainSaturation())
257 fLoGainMeanChargeErr = f;
258 else
259 fHiGainMeanChargeErr = f;
260
261}
262
263void MCalibrationChargePix::SetSigmaCharge( const Float_t f )
264{
265 if (IsHiGainSaturation())
266 fLoGainSigmaCharge = f;
267 else
268 fHiGainSigmaCharge = f;
269}
270
271
272void MCalibrationChargePix::SetSigmaChargeErr( const Float_t f )
273{
274 if (IsHiGainSaturation())
275 fLoGainSigmaChargeErr = f;
276 else
277 fHiGainSigmaChargeErr = f;
278
279}
280
281
282
283// --------------------------------------------------------------------------
284//
285// Set the conversion factors from outside (only for MC)
286//
287void MCalibrationChargePix::SetConversionFFactorMethod(Float_t c, Float_t err, Float_t sig)
288{
289 fMeanConversionFFactorMethod = c;
290 fConversionFFactorMethodErr = err;
291 fSigmaConversionFFactorMethod = sig;
292}
293
294// --------------------------------------------------------------------------
295//
296// Set the conversion factors from outside (only for MC)
297//
298void MCalibrationChargePix::SetConversionCombinedMethod(Float_t c, Float_t err, Float_t sig)
299{
300 fMeanConversionCombinedMethod = c;
301 fConversionCombinedMethodErr = err;
302 fSigmaConversionCombinedMethod = sig;
303}
304
305
306// --------------------------------------------------------------------------
307//
308// Set the conversion factors from outside (only for MC)
309//
310void MCalibrationChargePix::SetConversionBlindPixelMethod(Float_t c, Float_t err, Float_t sig)
311{
312 fMeanConversionBlindPixelMethod = c;
313 fConversionBlindPixelMethodErr = err;
314 fSigmaConversionBlindPixelMethod = sig;
315}
316
317// --------------------------------------------------------------------------
318//
319// Set the conversion factors from outside (only for MC)
320//
321void MCalibrationChargePix::SetConversionPINDiodeMethod(Float_t c, Float_t err, Float_t sig)
322{
323 fMeanConversionPINDiodeMethod = c ;
324 fConversionPINDiodeMethodErr = err;
325 fSigmaConversionPINDiodeMethod = sig;
326}
327
328// --------------------------------------------------------------------------
329//
330// Set the Hi Gain Saturation Bit from outside
331//
332void MCalibrationChargePix::SetHiGainSaturation(Bool_t b)
333{
334 b ? SETBIT(fFlags, kHiGainSaturation) : CLRBIT(fFlags, kHiGainSaturation);
335}
336
337// --------------------------------------------------------------------------
338//
339// Set the Lo Gain Saturation Bit from outside
340//
341void MCalibrationChargePix::SetLoGainSaturation(Bool_t b)
342{
343 b ? SETBIT(fFlags, kLoGainSaturation) : CLRBIT(fFlags, kLoGainSaturation);
344}
345
346// --------------------------------------------------------------------------
347//
348// Set the Excluded Bit from outside
349//
350void MCalibrationChargePix::SetExcluded(Bool_t b )
351{
352 b ? SETBIT(fFlags, kExcluded) : CLRBIT(fFlags, kExcluded);
353}
354
355// --------------------------------------------------------------------------
356//
357// Set the Excluded Bit from outside
358//
359void MCalibrationChargePix::SetChargeValid(Bool_t b )
360{
361 b ? SETBIT(fFlags, kChargeValid) : CLRBIT(fFlags, kChargeValid);
362}
363
364// --------------------------------------------------------------------------
365//
366// Set the Fitted Bit from outside
367//
368void MCalibrationChargePix::SetHiGainFitted(Bool_t b )
369{
370 b ? SETBIT(fFlags, kHiGainFitted) : CLRBIT(fFlags, kHiGainFitted);
371}
372
373// --------------------------------------------------------------------------
374//
375// Set the Fitted Bit from outside
376//
377void MCalibrationChargePix::SetLoGainFitted(const Bool_t b )
378{
379 b ? SETBIT(fFlags, kLoGainFitted) : CLRBIT(fFlags, kLoGainFitted);
380}
381
382// --------------------------------------------------------------------------
383//
384// Set the Ocillating Bit from outside
385//
386void MCalibrationChargePix::SetHiGainOscillating(const Bool_t b )
387{
388 b ? SETBIT(fFlags, kHiGainOscillating) : CLRBIT(fFlags, kHiGainOscillating);
389}
390
391// --------------------------------------------------------------------------
392//
393// Set the Ocillating Bit from outside
394//
395void MCalibrationChargePix::SetLoGainOscillating(const Bool_t b )
396{
397 b ? SETBIT(fFlags, kLoGainOscillating) : CLRBIT(fFlags, kLoGainOscillating);
398}
399
400// --------------------------------------------------------------------------
401//
402// Set the Excluded Bit from outside
403//
404void MCalibrationChargePix::SetBlindPixelMethodValid(const Bool_t b )
405{
406 b ? SETBIT(fFlags, kBlindPixelMethodValid) : CLRBIT(fFlags, kBlindPixelMethodValid);
407}
408
409// --------------------------------------------------------------------------
410//
411// Set the Excluded Bit from outside
412//
413void MCalibrationChargePix::SetFFactorMethodValid(const Bool_t b )
414{
415 b ? SETBIT(fFlags, kFFactorMethodValid) : CLRBIT(fFlags, kFFactorMethodValid);
416}
417
418// --------------------------------------------------------------------------
419//
420// Set the Excluded Bit from outside
421//
422void MCalibrationChargePix::SetPINDiodeMethodValid(const Bool_t b )
423{
424 b ? SETBIT(fFlags, kPINDiodeMethodValid) : CLRBIT(fFlags, kPINDiodeMethodValid);
425}
426
427// --------------------------------------------------------------------------
428//
429// Set the Excluded Bit from outside
430//
431void MCalibrationChargePix::SetCombinedMethodValid(const Bool_t b )
432{
433 b ? SETBIT(fFlags, kCombinedMethodValid) : CLRBIT(fFlags, kCombinedMethodValid);
434}
435
436void MCalibrationChargePix::SetAbsTimeBordersHiGain(const Byte_t f, const Byte_t l)
437{
438
439 fTimeFirstHiGain = f;
440 fTimeLastHiGain = l;
441
442}
443
444void MCalibrationChargePix::SetAbsTimeBordersLoGain(const Byte_t f, const Byte_t l)
445{
446
447 fTimeFirstLoGain = f;
448 fTimeLastLoGain = l;
449
450}
451
452void MCalibrationChargePix::SetMeanTimeInFirstBin(const Bool_t b)
453{
454 b ? SETBIT(fFlags,kMeanTimeInFirstBin) : CLRBIT(fFlags,kMeanTimeInFirstBin);
455}
456
457void MCalibrationChargePix::SetMeanTimeInLastBin(const Bool_t b)
458{
459 b ? SETBIT(fFlags,kMeanTimeInLastBin) : CLRBIT(fFlags,kMeanTimeInLastBin);
460}
461
462Float_t MCalibrationChargePix::GetMeanCharge() const
463{
464 return IsHiGainSaturation() ? fLoGainMeanCharge : fHiGainMeanCharge ;
465}
466
467Float_t MCalibrationChargePix::GetMeanChargeErr() const
468{
469 return IsHiGainSaturation() ? fLoGainMeanChargeErr : fHiGainMeanChargeErr ;
470}
471
472Float_t MCalibrationChargePix::GetChargeProb() const
473{
474 return IsHiGainSaturation() ? fLoGainChargeProb : fHiGainChargeProb ;
475}
476
477Float_t MCalibrationChargePix::GetSigmaCharge() const
478{
479 return IsHiGainSaturation() ? fLoGainSigmaCharge : fHiGainSigmaCharge ;
480}
481
482Float_t MCalibrationChargePix::GetSigmaChargeErr() const
483{
484 return IsHiGainSaturation() ? fLoGainSigmaChargeErr : fHiGainSigmaChargeErr ;
485}
486
487Bool_t MCalibrationChargePix::IsFitted() const
488{
489 return IsHiGainSaturation() ? IsLoGainFitted() : IsHiGainFitted();
490}
491
492Bool_t MCalibrationChargePix::IsOscillating() const
493{
494 return IsHiGainSaturation() ? IsLoGainOscillating() : IsHiGainOscillating();
495}
496
497Bool_t MCalibrationChargePix::IsMeanTimeInFirstBin() const
498{
499 return TESTBIT(fFlags,kMeanTimeInFirstBin);
500}
501
502Bool_t MCalibrationChargePix::IsMeanTimeInLastBin() const
503{
504 return TESTBIT(fFlags,kMeanTimeInLastBin);
505}
506
507
508
509Bool_t MCalibrationChargePix::IsExcluded() const
510{
511 return TESTBIT(fFlags,kExcluded);
512}
513
514Bool_t MCalibrationChargePix::IsHiGainSaturation() const
515{
516 return TESTBIT(fFlags,kHiGainSaturation);
517}
518
519Bool_t MCalibrationChargePix::IsLoGainSaturation() const
520{
521 return TESTBIT(fFlags,kLoGainSaturation);
522}
523
524Bool_t MCalibrationChargePix::IsChargeValid() const
525{
526 return TESTBIT(fFlags, kChargeValid);
527}
528
529Bool_t MCalibrationChargePix::IsHiGainFitted() const
530{
531 return TESTBIT(fFlags, kHiGainFitted);
532}
533
534Bool_t MCalibrationChargePix::IsLoGainFitted() const
535{
536 return TESTBIT(fFlags, kLoGainFitted);
537}
538
539Bool_t MCalibrationChargePix::IsHiGainOscillating() const
540{
541 return TESTBIT(fFlags, kHiGainOscillating);
542}
543
544Bool_t MCalibrationChargePix::IsLoGainOscillating() const
545{
546 return TESTBIT(fFlags, kLoGainOscillating);
547}
548
549Bool_t MCalibrationChargePix::IsBlindPixelMethodValid() const
550{
551 return TESTBIT(fFlags, kBlindPixelMethodValid);
552}
553
554Bool_t MCalibrationChargePix::IsFFactorMethodValid() const
555{
556 return TESTBIT(fFlags, kFFactorMethodValid);
557}
558
559Bool_t MCalibrationChargePix::IsPINDiodeMethodValid() const
560{
561 return TESTBIT(fFlags, kPINDiodeMethodValid);
562}
563
564Bool_t MCalibrationChargePix::IsCombinedMethodValid() const
565{
566 return TESTBIT(fFlags, kCombinedMethodValid);
567}
568
569
570//
571// The check return kTRUE if:
572//
573// 1) Pixel has a fitted charge greater than fChargeLimit*PedRMS
574// 2) Pixel has a fit error greater than fChargeErrLimit
575// 3) Pixel has a fitted charge greater its fChargeRelErrLimit times its charge error
576// 4) Pixel has a charge sigma bigger than its Pedestal RMS
577//
578Bool_t MCalibrationChargePix::CheckChargeValidity()
579{
580
581 if (!IsFitted())
582 return kFALSE;
583
584 if (GetMeanCharge() < fChargeLimit*GetPedRms())
585 {
586 *fLog << warn << "WARNING: Fitted Charge is smaller than "
587 << fChargeLimit << " Pedestal RMS in Pixel " << fPixId << endl;
588 return kFALSE;
589 }
590
591 if (GetMeanChargeErr() < fChargeErrLimit)
592 {
593 *fLog << warn << "WARNING: Error of Fitted Charge is smaller than "
594 << fChargeErrLimit << " in Pixel " << fPixId << endl;
595 return kFALSE;
596 }
597
598 if (GetMeanCharge() < fChargeRelErrLimit*GetMeanChargeErr())
599 {
600 *fLog << warn << "WARNING: Fitted Charge is smaller than "
601 << fChargeRelErrLimit << "* its error in Pixel " << fPixId << endl;
602 return kFALSE;
603 }
604
605 if (GetSigmaCharge() < GetPedRms())
606 {
607 *fLog << warn << "WARNING: Sigma of Fitted Charge smaller than Pedestal RMS in Pixel " << fPixId << endl;
608 return kFALSE;
609 }
610
611 SetChargeValid();
612 return kTRUE;
613}
614
615//
616// The check returns kTRUE if:
617//
618// The mean arrival time is at least 1.0 slices from the used edge slices
619//
620Bool_t MCalibrationChargePix::CheckTimeValidity()
621{
622
623 const Float_t loweredge = IsHiGainSaturation() ? fTimeFirstHiGain : fTimeFirstLoGain;
624 const Float_t upperedge = IsHiGainSaturation() ? fTimeLastHiGain : fTimeLastLoGain;
625
626 if ( fAbsTimeMean < loweredge+1.)
627 {
628 *fLog << warn << "WARNING: Mean ArrivalTime in first extraction bin of the Pixel " << fPixId << endl;
629 SetMeanTimeInFirstBin();
630 return kFALSE;
631 }
632
633 if ( fAbsTimeMean > upperedge-1.)
634 {
635 *fLog << warn << "WARNING: Mean ArrivalTime in last extraction bin of the Pixel " << fPixId << endl;
636 SetMeanTimeInLastBin();
637 return kFALSE;
638 }
639
640 return kTRUE;
641}
642
643void MCalibrationChargePix::CalcLoGainPed()
644{
645
646 Float_t pedRmsSquare = fPedRms * fPedRms;
647 Float_t pedRmsSquareErrSquare = fPedErr * fPedErr * pedRmsSquare; // fPedRmsErr = fPedErr/2.
648 Float_t pedRmsSquareErr;
649
650 //
651 // We do not know the Lo Gain Pedestal RMS, so we have to retrieve it
652 // from the HI GAIN (all calculation per slice up to now):
653 //
654 // We extract the pure NSB contribution:
655 //
656 const Float_t elecRmsSquare = fElectronicPedRms * fElectronicPedRms;
657 const Float_t elecRmsSquareErrSquare = 4.*fElectronicPedRmsErr * fElectronicPedRmsErr * elecRmsSquare;
658
659 Float_t nsbSquare = pedRmsSquare - elecRmsSquare;
660 Float_t nsbSquareRelErrSquare = (pedRmsSquareErrSquare + elecRmsSquareErrSquare)
661 / (nsbSquare * nsbSquare) ;
662
663 if (nsbSquare < 0.)
664 nsbSquare = 0.;
665
666 //
667 // Now, we divide the NSB by the conversion factor and
668 // add it quadratically to the electronic noise
669 //
670 const Float_t conversionSquare = fConversionHiLo * fConversionHiLo;
671 const Float_t conversionSquareRelErrSquare = 4.*fConversionHiLoErr * fConversionHiLoErr / conversionSquare;
672
673 const Float_t convertedNsbSquare = nsbSquare / conversionSquare;
674 const Float_t convertedNsbSquareErrSquare = (nsbSquareRelErrSquare + conversionSquareRelErrSquare)
675 * convertedNsbSquare * convertedNsbSquare;
676
677 pedRmsSquare = convertedNsbSquare + elecRmsSquare;
678 pedRmsSquareErr = TMath::Sqrt( convertedNsbSquareErrSquare + elecRmsSquareErrSquare );
679
680 //
681 // Correct also for the conversion to Hi-Gain:
682 //
683 fLoGainPedRms = TMath::Sqrt(pedRmsSquare) * fConversionHiLo;
684 fLoGainPedRmsErr = 0.5 * pedRmsSquareErr / fLoGainPedRms * fConversionHiLo;
685}
686
687Bool_t MCalibrationChargePix::CalcReducedSigma()
688{
689
690 const Float_t sigmaSquare = GetSigmaCharge() * GetSigmaCharge();
691 const Float_t sigmaSquareErrSquare = 4.*GetSigmaChargeErr()* GetSigmaChargeErr() * sigmaSquare;
692
693 Float_t pedRmsSquare;
694 Float_t pedRmsSquareErrSquare;
695
696 if (IsHiGainSaturation())
697 {
698 CalcLoGainPed();
699
700 pedRmsSquare = fLoGainPedRms * fLoGainPedRms;
701 pedRmsSquareErrSquare = 4.* fLoGainPedRmsErr * fLoGainPedRmsErr * pedRmsSquare;
702 } /* if (HiGainSaturation) */
703 else
704 {
705
706 pedRmsSquare = fPedRms * fPedRms;
707 pedRmsSquareErrSquare = fPedErr * fPedErr * pedRmsSquare; // fPedRmsErr = fPedErr/2.
708 }
709 //
710 // Calculate the reduced sigmas
711 //
712 const Float_t rsigmachargesquare = sigmaSquare - pedRmsSquare;
713
714 if (rsigmachargesquare <= 0.)
715 {
716 *fLog << warn
717 << "WARNING: Cannot calculate the reduced sigma: smaller than 0 in pixel "
718 << fPixId << endl;
719 return kFALSE;
720 }
721
722 fRSigmaCharge = TMath::Sqrt(rsigmachargesquare);
723 fRSigmaChargeErr = TMath::Sqrt(sigmaSquareErrSquare + pedRmsSquareErrSquare) / 2. / fRSigmaCharge;
724
725 return kTRUE;
726}
727
728//
729// Calculate the number of photo-electrons after the F-Factor method
730// Calculate the errors of the F-Factor method
731//
732Bool_t MCalibrationChargePix::CalcFFactorMethod()
733{
734
735 if (!IsChargeValid())
736 return kFALSE;
737
738 if (fRSigmaCharge < 0.)
739 return kFALSE;
740
741 //
742 // Square all variables in order to avoid applications of square root
743 //
744 // First the relative error squares
745 //
746 const Float_t chargeSquare = GetMeanCharge() * GetMeanCharge();
747 const Float_t chargeSquareRelErrSquare = 4.* GetMeanChargeErr() * GetMeanChargeErr() / chargeSquare;
748
749 const Float_t chargeRelErrSquare = GetMeanChargeErr() * GetMeanChargeErr()
750 / (GetMeanCharge() * GetMeanCharge());
751
752 const Float_t ffactorsquare = gkFFactor * gkFFactor;
753 const Float_t ffactorsquareRelErrSquare = 4.*gkFFactorErr * gkFFactorErr / ffactorsquare;
754
755 const Float_t avQERelErrSquare = fAverageQEErr * fAverageQEErr / fAverageQE / fAverageQE;
756
757 const Float_t avQEFFactor = TMath::Sqrt( ( 1. - fAverageQE ) / fAverageQE );
758 const Float_t avQEFFactorErr = 1./ ( 2. * avQEFFactor ) * fAverageQEErr
759 / ( fAverageQE * fAverageQE );
760
761 const Float_t avQEFFactorRelErrSquare = avQEFFactorErr * avQEFFactorErr
762 / ( avQEFFactor * avQEFFactor) ;
763
764 const Float_t rsigmaSquare = fRSigmaCharge * fRSigmaCharge;
765 const Float_t rsigmaSquareRelErrSquare = 4.* fRSigmaChargeErr * fRSigmaChargeErr / rsigmaSquare;
766
767
768 //
769 // Calculate the number of phe's from the F-Factor method
770 // (independent on Hi Gain or Lo Gain)
771 //
772 fPheFFactorMethod = ffactorsquare * chargeSquare / rsigmaSquare;
773
774 if (fPheFFactorMethod < 0.)
775 return kFALSE;
776
777 //
778 // Calculate the Error of Nphe
779 //
780 const Float_t pheFFactorRelErrSquare = ffactorsquareRelErrSquare
781 + chargeSquareRelErrSquare
782 + rsigmaSquareRelErrSquare;
783 fPheFFactorMethodErr = TMath::Sqrt(pheFFactorRelErrSquare) * fPheFFactorMethod;
784
785 //
786 // Calculate the conversion factor between PHOTONS and FADC counts
787 //
788 // Nphot = Nphe / avQE
789 // conv = Nphot / FADC counts
790 //
791 fMeanConversionFFactorMethod = fPheFFactorMethod / fAverageQE / GetMeanCharge();
792
793
794 //
795 // Calculate the error of the mean conversion factor between PHOTONS and FADC counts
796 //
797 const Float_t convRelErrSquare = ( pheFFactorRelErrSquare + chargeRelErrSquare + avQERelErrSquare);
798
799 if (convRelErrSquare < 0.)
800 return kFALSE;
801
802
803 const Float_t convRelErr = TMath::Sqrt(convRelErrSquare);
804 fConversionFFactorMethodErr = convRelErr * fMeanConversionFFactorMethod;
805
806 if (convRelErr < fConvFFactorRelErrLimit)
807 SetFFactorMethodValid();
808
809 //
810 // Calculate the Total F-Factor of the camera (in photons)
811 //
812 fTotalFFactorFFactorMethod = (fRSigmaCharge/GetMeanCharge())*TMath::Sqrt(fPheFFactorMethod);
813 fTotalFFactorFFactorMethod *= avQEFFactor;
814
815 //
816 // Calculate the error of the Total F-Factor of the camera ( in photons )
817 //
818 const Float_t rsigmaChargeRelErrSquare = fRSigmaChargeErr * fRSigmaChargeErr
819 / (fRSigmaCharge * fRSigmaCharge) ;
820
821 fTotalFFactorErrFFactorMethod = TMath::Sqrt( rsigmaChargeRelErrSquare
822 + chargeRelErrSquare
823 + pheFFactorRelErrSquare
824 + avQEFFactorRelErrSquare );
825 fTotalFFactorErrFFactorMethod *= fTotalFFactorFFactorMethod;
826
827 //
828 // Calculate the sigma of the conversion from FADC counts to photons
829 //
830 fSigmaConversionFFactorMethod = GetTotalFFactorFFactorMethod()*TMath::Sqrt(fMeanConversionFFactorMethod);
831
832 return kTRUE;
833}
834
835void MCalibrationChargePix::ApplyLoGainConversion()
836{
837
838 const Float_t chargeRelErrSquare = fLoGainMeanChargeErr * fLoGainMeanChargeErr
839 /( fLoGainMeanCharge * fLoGainMeanCharge );
840 const Float_t sigmaRelErrSquare = fLoGainSigmaChargeErr * fLoGainSigmaChargeErr
841 /( fLoGainSigmaCharge * fLoGainSigmaCharge );
842 const Float_t conversionRelErrSquare = fConversionHiLoErr * fConversionHiLoErr
843 /( fConversionHiLo * fConversionHiLo );
844
845 fLoGainMeanCharge *= fConversionHiLo;
846 fLoGainMeanChargeErr = TMath::Sqrt(chargeRelErrSquare + conversionRelErrSquare) * fLoGainMeanCharge;
847
848 fLoGainSigmaCharge *= fConversionHiLo;
849 fLoGainSigmaChargeErr = TMath::Sqrt(sigmaRelErrSquare + conversionRelErrSquare) * fLoGainSigmaCharge;
850
851 fElectronicPedRms = gkElectronicPedRms * TMath::Sqrt(fNumLoGainSamples);
852 fElectronicPedRmsErr = gkElectronicPedRmsErr * TMath::Sqrt(fNumLoGainSamples);
853
854}
855
Note: See TracBrowser for help on using the repository browser.