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

Last change on this file since 3699 was 3684, checked in by gaug, 21 years ago
*** empty log message ***
File size: 26.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// MCalibrationChargePix
26//
27// Storage container of the calibrated Charge of one pixel.
28//
29// The following values are initialized to meaningful values:
30//
31// - The Electronic Rms to 1.5 per FADC slice
32// - The uncertainty about the Electronic RMS to 0.3 per slice
33// - The F-Factor is assumed to have been measured in Munich to 1.13 - 1.17.
34// with the Munich definition of the F-Factor, thus:
35// F = Sigma(Out)/Mean(Out) * Mean(In)/Sigma(In)
36// Mean F-Factor (gkFFactor) = 1.15
37// Error F-Factor (gkFFactorErr) = 0.02
38//
39// The following variables are calculated inside this class:
40// - fLoGainPedRmsSquare and fLoGainPedRmsSquareVar (see CalcLoGainPedestal())
41// - fRSigmaSquare and fRSigmaSquareVar (see CalcReducedSigma() )
42// - fPheFFactorMethod and fPheFFactorMethodVar (see CalcFFactorMethod() )
43//
44// The following variables are set by MHCalibrationChargeCam:
45// - fAbsTimeMean and fAbsTimeRms
46// - all variables in MCalibrationPix
47//
48// The following variables are set by MCalibrationChargeCalc:
49// - fPed, fPedVar and fPedRms
50// - fMeanConvFADC2Phe
51// - fConvFADC2PheVar
52// - fSigmaConvFADC2Phe
53// - fTotalFFactorFFactorMethod
54// - fTotalFFactorFFactorMethodVar
55//
56// The following variables are not yet implemented:
57// - fConversionHiLo and fConversionHiLoVar (now set fixed to 10. +- 2.5)
58//
59// Error of all variables are calculated by error-propagation. Note that internally,
60// all error variables contain Variances in order to save the CPU-intensive square rooting
61//
62// Low-Gain variables are stored internally unconverted, i.e. directly from the summed
63// FADC slices extraction results, but can be retrieved converted to High-Gain amplifications
64// by calls to: GetConvertedLoGainMean() or GetConvertedLoGainSigma()
65//
66// See also: MCalibrationChargeCam, MCalibrationChargeCalc,
67// MHCalibrationChargeCam, MHCalibrationChargePix
68//
69/////////////////////////////////////////////////////////////////////////////
70#include "MCalibrationChargePix.h"
71
72#include "MLog.h"
73#include "MLogManip.h"
74
75#include "MBadPixelsPix.h"
76
77ClassImp(MCalibrationChargePix);
78
79using namespace std;
80
81const Float_t MCalibrationChargePix::gkElectronicPedRms = 1.5;
82const Float_t MCalibrationChargePix::gkElectronicPedRmsErr = 0.3;
83const Float_t MCalibrationChargePix::gkFFactor = 1.15;
84const Float_t MCalibrationChargePix::gkFFactorErr = 0.02;
85
86const Float_t MCalibrationChargePix::fgConversionHiLo = 10.;
87const Float_t MCalibrationChargePix::fgConversionHiLoErr = 2.5;
88const Float_t MCalibrationChargePix::fgPheFFactorMethodLimit = 5.;
89const Float_t MCalibrationChargePix::fgConvFFactorRelErrLimit = 0.35;
90// --------------------------------------------------------------------------
91//
92// Default Constructor:
93//
94// Sets:
95// - fCalibFlags to 0
96// - fConversionHiLo to fgConversionHiLo
97// - fConversionHiLoVar to square of fgConversionHiLoErr
98// - fConvFFactorRelErrLimit to fgConvFFactorRelErrLimit*fgConvFFactorRelErrLimit
99// - fPheFFactorLimit to fgPheFFactorLimit
100//
101// Calls:
102// - Clear()
103//
104MCalibrationChargePix::MCalibrationChargePix(const char *name, const char *title)
105 : fCalibFlags(0)
106{
107
108 fName = name ? name : "MCalibrationChargePix";
109 fTitle = title ? title : "Container of the fit results of MHCalibrationChargePixs ";
110
111 //
112 // At the moment, we don't have a database, yet,
113 // so we get it from the configuration file
114 //
115 SetConversionHiLo();
116 SetConversionHiLoErr();
117
118 SetPheFFactorMethodLimit();
119 SetConvFFactorRelErrLimit();
120
121 Clear();
122}
123
124// ------------------------------------------------------------------------
125//
126// Sets:
127// - all flags to kFALSE
128// - all variables to -1.
129//
130// Calls:
131// - MCalibrationPix::Clear()
132//
133void MCalibrationChargePix::Clear(Option_t *o)
134{
135
136 SetFFactorMethodValid ( kFALSE );
137
138 fRSigmaSquare = -1.;
139 fRSigmaSquareVar = -1.;
140
141 fPed = -1.;
142 fPedRms = -1.;
143 fPedVar = -1.;
144
145 fLoGainPedRmsSquare = -1.;
146 fLoGainPedRmsSquareVar = -1.;
147
148 fAbsTimeMean = -1.;
149 fAbsTimeRms = -1.;
150
151 fPheFFactorMethod = -1.;
152 fPheFFactorMethodVar = -1.;
153
154 fMeanConvFADC2Phe = -1.;
155 fMeanConvFADC2PheVar = -1.;
156 fMeanFFactorFADC2Phot = -1.;
157 fMeanFFactorFADC2PhotVar = -1.;
158
159 MCalibrationPix::Clear();
160}
161
162
163// --------------------------------------------------------------------------
164//
165// Set F-Factor Method Validity Bit from outside
166//
167void MCalibrationChargePix::SetFFactorMethodValid(const Bool_t b )
168{
169 b ? SETBIT(fCalibFlags, kFFactorMethodValid) : CLRBIT(fCalibFlags, kFFactorMethodValid);
170}
171
172// --------------------------------------------------------------------------
173//
174// Set pedestals from outside (done by MCalibrationChargeCalc)
175//
176void MCalibrationChargePix::SetPedestal(const Float_t ped, const Float_t pedrms, const Float_t pederr)
177{
178
179 fPed = ped;
180 fPedRms = pedrms;
181 fPedVar = pederr*pederr;
182}
183
184// -------------------------------------------------------------------------------
185//
186// Get the conversion Error Hi-Gain to Low-Gain:
187// - If fConversionHiLoVar is smaller than 0 (i.e. has not yet been set), return -1.
188//
189Float_t MCalibrationChargePix::GetConversionHiLoErr() const
190{
191 if (fConversionHiLoVar < 0.)
192 return -1.;
193
194 return TMath::Sqrt(fConversionHiLoVar);
195}
196
197// --------------------------------------------------------------------------
198//
199// Get the relative variance of the conversion factor between higain and logain:
200// - If fConversionHiLo is 0, return -1.
201// - If fConversionHiLoVar is smaller than 0, return -1.
202// - Else returns: fConversionHiLoVar / fConversionHiLo^2
203//
204const Float_t MCalibrationChargePix::GetConversionHiLoRelVar() const
205{
206
207 if (fConversionHiLoVar < 0.)
208 return -1.;
209
210 if (fConversionHiLo == 0.)
211 return -1.;
212
213 return fConversionHiLoVar / (fConversionHiLo * fConversionHiLo);
214}
215
216
217// --------------------------------------------------------------------------
218//
219// Get the relative variance of the conversion factor between higain and logain:
220// - If gkFFactor is 0, return -1.
221// - If gkFFactorErr is smaller than 0, return -1.
222// - Else returns: gkFFactorErr^2 / gkFFactor*^2
223//
224const Float_t MCalibrationChargePix::GetFFactorRelVar() const
225{
226
227 if (gkFFactorErr < 0.)
228 return -1.;
229
230 if (gkFFactor == 0.)
231 return -1.;
232
233 return gkFFactorErr * gkFFactorErr / (gkFFactor * gkFFactor);
234}
235
236
237//
238// Get the Error of the Mean pedestals:
239// Returns square root of fPedVar
240//
241Float_t MCalibrationChargePix::GetPedErr() const
242{
243 return TMath::Sqrt(fPedVar);
244}
245
246// --------------------------------------------------------------------------
247//
248// Get the pedestals RMS:
249// - Test bit kHiGainSaturation:
250// If yes, return square root of fLoGainPedRmsSquare (if greater than 0, otherwise -1.),
251// If no, return fPedRms
252//
253Float_t MCalibrationChargePix::GetPedRms() const
254{
255
256 if (IsHiGainSaturation())
257 if (fLoGainPedRmsSquare < 0.)
258 return -1.;
259 else
260 return TMath::Sqrt(fLoGainPedRmsSquare);
261
262 return fPedRms;
263}
264
265// --------------------------------------------------------------------------
266//
267// Get the Error of the pedestals RMS:
268// - Test bit kHiGainSaturation:
269// If yes, return square root of (0.25*fLoGainPedRmsSquareVar/ fLoGainPedRmsSquare) (if greater than 0, otherwise -1.)
270// If no , return square root of (fPedVar) (if greater than 0, otherwise -1.), divided by 2.
271//
272Float_t MCalibrationChargePix::GetPedRmsErr() const
273{
274 if (IsHiGainSaturation())
275 if (fLoGainPedRmsSquareVar < 0.)
276 return -1.;
277 else
278 return TMath::Sqrt(0.25*fLoGainPedRmsSquareVar/fLoGainPedRmsSquare);
279 else
280 if (fPedVar < 0.)
281 return -1.;
282 else
283 return TMath::Sqrt(fPedVar)/2.;
284}
285
286
287// --------------------------------------------------------------------------
288//
289// Get the Low Gain Mean converted to High Gain amplification:
290// Returns fLoGainMean multiplied with fConversionHiLo
291//
292Float_t MCalibrationChargePix::GetConvertedLoGainMean() const
293{
294 return fLoGainMean * fConversionHiLo;
295}
296
297// --------------------------------------------------------------------------
298//
299// Get the Error of the converted Low Gain Mean:
300//
301// Returns -1 if the variable fLoGainMean or fLoGainMeanVar are smaller than 0.
302//
303// Returns the square root of the quadratic sum of the relative variances of
304// the fLoGainMean and fConversionHiLo, mulitplied with GetConvertedLoGainMean()
305//
306Float_t MCalibrationChargePix::GetConvertedLoGainMeanErr() const
307{
308
309 const Float_t logainrelvar = GetLoGainMeanRelVar();
310
311 if (logainrelvar < 0.)
312 return -1.;
313
314 return TMath::Sqrt(logainrelvar + GetConversionHiLoRelVar()) * GetConvertedLoGainMean();
315}
316
317// --------------------------------------------------------------------------
318//
319// Get the Low Gain Sigma converted to High Gain amplification:
320// Returns fLoGainSigma multiplied with fConversionHiLo
321//
322Float_t MCalibrationChargePix::GetConvertedLoGainSigma() const
323{
324 return fLoGainSigma * fConversionHiLo;
325}
326
327// --------------------------------------------------------------------------
328//
329// Get the Error of the converted Low Gain Sigma:
330//
331// Returns -1 if the variable fLoGainSigma or fLoGainSigmaVar are smaller than 0.
332//
333// Returns the square root of the quadratic sum of the relative variances of
334// the fLoGainSigma and fConversionHiLo, mulitplied with GetConvertedLoGainSigma()
335//
336Float_t MCalibrationChargePix::GetConvertedLoGainSigmaErr() const
337{
338
339 if (fLoGainSigmaVar < 0.)
340 return -1.;
341
342 if (fLoGainSigma < 0.)
343 return -1.;
344
345 const Float_t sigmaRelVar = fLoGainSigmaVar
346 /( fLoGainSigma * fLoGainSigma );
347
348 return TMath::Sqrt(sigmaRelVar+GetConversionHiLoRelVar()) * GetConvertedLoGainSigma();
349}
350
351// --------------------------------------------------------------------------
352//
353// Get the reduced Sigma:
354// - If fRSigmaSquare is smaller than 0 (i.e. has not yet been set), return -1.
355// - Test bit kHiGainSaturation:
356// If yes, return square root of fRSigmaSquare, multiplied with fConversionHiLo,
357// If no , return square root of fRSigmaSquare
358//
359Float_t MCalibrationChargePix::GetRSigma() const
360{
361 if (fRSigmaSquare < 0)
362 return -1;
363
364 const Float_t rsigma = TMath::Sqrt(fRSigmaSquare);
365
366 return IsHiGainSaturation() ? rsigma*fConversionHiLo : rsigma ;
367}
368
369// --------------------------------------------------------------------------
370//
371// Get the error of the reduced Sigma:
372// - If fRSigmaSquareVar is smaller than 0 (i.e. has not yet been set), return -1.
373// - Calculate the absolute variance of the reduced sigma with the formula:
374// reduced sigma variance = 0.25 * fRSigmaSquareVar / fRSigmaSquare
375// - Test bit kHiGainSaturation:
376// If yes, returns the square root of the quadratic sum of the relative variances of the
377// reduced sigma and fConversionHiLo, mulitplied with GetRSigma()
378// Else returns the square root of rel. (0.25*fRSigmaSquareVar / fRSigmaSquare)
379//
380Float_t MCalibrationChargePix::GetRSigmaErr() const
381{
382
383 if (fRSigmaSquareVar < 0)
384 return -1;
385
386 //
387 // SigmaSquareVar = 4. * Sigma * Sigma * Var(sigma)
388 // ==> Var(sigma) = 0.25 * SigmaSquareVar / (Sigma * Sigma)
389 //
390 const Float_t rsigmaVar = 0.25 * fRSigmaSquareVar / fRSigmaSquare;
391
392 if (IsHiGainSaturation())
393 return TMath::Sqrt(rsigmaVar/fRSigmaSquare + GetConversionHiLoRelVar()) * GetRSigma();
394 else
395 return TMath::Sqrt(rsigmaVar);
396
397}
398
399// --------------------------------------------------------------------------
400//
401// Get the reduced Sigma Square:
402// - If fRSigmaSquare is smaller than 0 (i.e. has not yet been set), return -1.
403// - Test bit kHiGainSaturation:
404// If yes, return fRSigmaSquare, multiplied with fConversionHiLo^2,
405// If no , return fRSigmaSquare
406//
407Float_t MCalibrationChargePix::GetRSigmaSquare() const
408{
409 if (fRSigmaSquare < 0)
410 return -1;
411
412 return IsHiGainSaturation() ? fRSigmaSquare*fConversionHiLo*fConversionHiLo : fRSigmaSquare ;
413}
414
415// --------------------------------------------------------------------------
416//
417// Get the relative variance of the reduced Sigma:
418// - If fRSigmaSquareVar is smaller than 0 (i.e. has not yet been set), return -1.
419// - Calculate the relative variance of the reduced sigma squares with the formula:
420// reduced sigma rel. variance = 0.25 * fRSigmaSquareVar / fRSigmaSquare / fRSigmaSquare
421// - Test bit kHiGainSaturation:
422// If yes, returns the sum of the relative variances of the reduced sigma and fConversionHiLo
423// Else returns the relative variance of the reduced sigma
424//
425Float_t MCalibrationChargePix::GetRSigmaRelVar() const
426{
427
428 if (fRSigmaSquareVar < 0)
429 return -1;
430
431 //
432 // SigmaSquareVar = 4. * Sigma * Sigma * Var(sigma)
433 // ==> Var(sigma) = 0.25 * SigmaSquareVar / (Sigma * Sigma)
434 //
435 const Float_t rsigmaRelVar = 0.25 * fRSigmaSquareVar / ( fRSigmaSquare * fRSigmaSquare );
436
437 if (IsHiGainSaturation())
438 return rsigmaRelVar + GetConversionHiLoRelVar();
439 else
440 return rsigmaRelVar;
441}
442
443// --------------------------------------------------------------------------
444//
445// Get the error on the number of photo-electrons (F-Factor Method):
446// - If fPheFFactorMethodVar is smaller than 0 (i.e. has not yet been set), return -1.
447// - Else returns the square root of fPheFFactorMethodVar
448//
449Float_t MCalibrationChargePix::GetPheFFactorMethodErr() const
450{
451 if (fPheFFactorMethodVar < 0.)
452 return -1.;
453 return TMath::Sqrt(fPheFFactorMethodVar);
454}
455
456// --------------------------------------------------------------------------
457//
458// Get the error on the mean total F-Factor of the signal readout (F-Factor Method):
459// - If fMeanFFactorFADC2PhotVar is smaller than 0 (i.e. has not yet been set), return -1.
460// - Else returns the square root of fMeanFFactorFADC2PhotVar
461//
462Float_t MCalibrationChargePix::GetMeanFFactorFADC2PhotErr() const
463{
464 if (fMeanFFactorFADC2PhotVar < 0.)
465 return -1.;
466 return TMath::Sqrt(fMeanFFactorFADC2PhotVar);
467}
468
469// --------------------------------------------------------------------------
470//
471// Get the relative variance on the number of photo-electrons (F-Factor Method):
472// - If fPheFFactorMethodVar is smaller than 0 (i.e. has not yet been set), return -1.
473// - If fPheFFactorMethod is 0, return -1.
474// - Else returns fPheFFactorMethodVar / fPheFFactorMethod^2
475//
476Float_t MCalibrationChargePix::GetPheFFactorMethodRelVar() const
477{
478 if (fPheFFactorMethodVar < 0.)
479 return -1.;
480 if (fPheFFactorMethod == 0.)
481 return -1.;
482
483 return fPheFFactorMethodVar / (fPheFFactorMethod * fPheFFactorMethod);
484}
485
486
487// --------------------------------------------------------------------------
488//
489// Get the error on the mean conversion factor (FFactor  Method):
490// - If fMeanConvFADC2PheVar is smaller than 0 (i.e. has not yet been set), return -1.
491// - Else returns the square root of fMeanConvFADC2PheVar
492//
493Float_t MCalibrationChargePix::GetMeanConvFADC2PheErr() const
494{
495 if (fMeanConvFADC2PheVar < 0.)
496 return -1.;
497 return TMath::Sqrt(fMeanConvFADC2PheVar);
498}
499
500// --------------------------------------------------------------------------
501//
502// Test bit kFFactorMethodValid
503//
504Bool_t MCalibrationChargePix::IsFFactorMethodValid() const
505{
506 return TESTBIT(fCalibFlags, kFFactorMethodValid);
507}
508
509
510// ----------------------------------------------------------------------------
511//
512// - If fSigma is smaller than 0 (i.e. has not yet been set), return kFALSE
513// - If fPedRms is smaller than 0 (i.e. has not yet been set), return kFALSE
514//
515// Calculate the reduced sigma of the low-Gain FADC slices:
516// - Test bit IsHiGainSaturation() for the Sigma:
517// If yes, take fLoGainSigma and fLoGainSigmaVar
518// If no , take fHiGainSigma and fHiGainSigmaVar
519//
520// - Test bit IsHiGainSaturation() for the pedRMS:
521// If yes, take fLoGainPedRmsSquare and fLoGainPedRmsSquareVar
522// If no , take fPedRms and fPedVar
523//
524// - Calculate the reduced sigma with the formula:
525// fRSigmaSquare = Sigma*Sigma - pedRMS*pedRMS
526//
527// - If fRSigmaSquare is smaller than 0, give a warning and return kFALSE
528//
529// - Calculate the variance of the reduced sigma with the formula:
530// fRSigmaSquareVar = 4.* (sigmaVar*Sigma*Sigma + pedRmsVar*pedRMS*pedRMS)
531//
532// A back-transformation to the corr. amplification factor of the High-Gain is done
533// in GetRSigma() and GetRSigmaErr()
534//
535Bool_t MCalibrationChargePix::CalcReducedSigma()
536{
537
538 if (GetSigma() < 0.)
539 return kFALSE;
540
541 if (GetPedRms() < 0.)
542 return kFALSE;
543
544 const Float_t sigma = IsHiGainSaturation() ? fLoGainSigma : fHiGainSigma ;
545 const Float_t sigmavar = IsHiGainSaturation() ? fLoGainSigmaVar : fHiGainSigmaVar;
546 const Float_t pedRmsSquare = IsHiGainSaturation() ? fLoGainPedRmsSquare : fPedRms*fPedRms;
547 const Float_t pedRmsSquareVar = IsHiGainSaturation() ? fLoGainPedRmsSquareVar : 0.25*fPedVar*pedRmsSquare;
548
549 const Float_t sigmaSquare = sigma * sigma;
550 const Float_t sigmaSquareVar = 4. * sigmavar * sigmaSquare;
551
552 //
553 // Calculate the reduced sigmas
554 //
555 fRSigmaSquare = sigmaSquare - pedRmsSquare;
556 if (fRSigmaSquare <= 0.)
557 {
558 *fLog << warn
559 << "WARNING: Cannot calculate the reduced sigma: smaller than 0 in pixel "
560 << fPixId << endl;
561 return kFALSE;
562 }
563
564 fRSigmaSquareVar = 4. * (sigmaSquareVar + pedRmsSquareVar);
565
566 return kTRUE;
567}
568
569// ------------------------------------------------------------------
570//
571// If fRSigmaSquare is smaller than 0 (i.e. has not yet been set),
572// set kFFactorMethodValid to kFALSE and return kFALSE
573//
574// Calculate the number of photo-electrons with the F-Factor method:
575// - Test bit IsHiGainSaturation() for the Mean Sum of FADC slices:
576// If yes, take fLoGainMean and fLoGainMeanVar
577// If no , take fHiGainMean and fHiGainMeanVar
578//
579// - Test bit IsHiGainSaturation() for the pedRMS:
580// If yes, take fLoGainPedRmsSquare and fLoGainPedRmsSquareVar
581// If no , take fPedRms and fPedVar
582//
583// - Calculate the number of photo-electrons with the formula:
584// fPheFFactorMethod = gkFFactor*gkFFactor * Mean * Mean / fRSigmaSquare
585//
586// - Calculate the Variance on the photo-electrons with the formula:
587// fPheFFactorMethodVar = ( 4. * gkFFactorErr * gkFFactorErr / ( gkFFactor * gkFFactor )
588// + 4. * Mean Var. / ( Mean * Mean )
589// + fRSigmaSquareVar / fRSigmaSquare
590// ) * fPheFFactorMethod * fPheFFactorMethod
591//
592// - If fPheFFactorMethod is less than fPheFFactorMethodLimit,
593// set kFFactorMethodValid to kFALSE and return kFALSE
594// else: Set kFFactorMethodValid to kTRUE and return kTRUE
595//
596Bool_t MCalibrationChargePix::CalcFFactorMethod()
597{
598
599 if (fRSigmaSquare < 0.)
600 return kFALSE;
601
602 //
603 // Square all variables in order to avoid applications of square root
604 //
605 const Float_t meanSquare = GetMean() * GetMean();
606 const Float_t meanSquareRelVar = 4.* GetMeanRelVar();
607
608 const Float_t ffactorsquare = gkFFactor * gkFFactor;
609 const Float_t ffactorsquareRelVar = 4.* GetFFactorRelVar();
610
611 const Float_t rsigmaSquareRelVar = fRSigmaSquareVar / fRSigmaSquare / fRSigmaSquare;
612 //
613 // Calculate the number of phe's from the F-Factor method
614 // (independent on Hi Gain or Lo Gain)
615 //
616 fPheFFactorMethod = ffactorsquare * meanSquare / fRSigmaSquare;
617
618 if (fPheFFactorMethod < fPheFFactorMethodLimit)
619 return kFALSE;
620
621 //
622 // Calculate the Error of Nphe
623 //
624 const Float_t pheRelVar = ffactorsquareRelVar + meanSquareRelVar + rsigmaSquareRelVar;
625 fPheFFactorMethodVar = pheRelVar * fPheFFactorMethod * fPheFFactorMethod;
626
627 if (fPheFFactorMethodVar < 0. )
628 return kFALSE;
629
630 fMeanConvFADC2Phe = fPheFFactorMethod / GetMean();
631
632 if (fMeanConvFADC2Phe < 0. )
633 return kFALSE;
634
635 //
636 // In the calculation of the number of phe's one mean square has already been used.
637 // Now, we divide by another mean, so one mean calcels out, we cannot directly propagate
638 // the errors, but have to take account of this cancellation:
639 //
640 const Float_t convrelvar = ffactorsquareRelVar + GetMeanRelVar() + rsigmaSquareRelVar;
641
642 if (convrelvar > fConvFFactorRelVarLimit || convrelvar < 0.)
643 {
644 *fLog << warn << GetDescriptor() << ": Conversion F-Factor Method Rel. Variance: "
645 << convrelvar << " above limits of: [0," << Form("%3.2f",fConvFFactorRelVarLimit)
646 << "] in pixel: " << fPixId << endl;
647 return kFALSE;
648 }
649
650 fMeanConvFADC2PheVar = convrelvar * fMeanConvFADC2Phe * fMeanConvFADC2Phe;
651
652 SetFFactorMethodValid(kTRUE);
653 return kTRUE;
654}
655
656// ----------------------------------------------------------------------------------
657//
658// If photflux is smaller or equal 0, return kFALSE
659//
660// Calculate the total F-Factor with the formula:
661// fMeanFFactorFADC2Phot = Sqrt ( fRSigmaSquare ) / GetMean() * sqrt( nphotons )
662//
663// Calculate the error of the total F-Factor
664//
665Bool_t MCalibrationChargePix::CalcMeanFFactor( const Float_t nphotons, const Float_t nphotonsrelvar )
666{
667
668 if (nphotons <= 0.)
669 {
670 *fLog << warn << GetDescriptor() << ": Assumed photon flux is smaller or equal 0." << endl;
671 return kFALSE;
672 }
673
674 if (nphotonsrelvar < 0.)
675 {
676 *fLog << warn << GetDescriptor() << ": Assumed photon flux variance is smaller than 0." << endl;
677 return kFALSE;
678 }
679
680 fMeanFFactorFADC2Phot = TMath::Sqrt(fRSigmaSquare * nphotons) / GetMean() ;
681
682 if (fMeanFFactorFADC2Phot < 0.)
683 {
684 *fLog << warn << GetDescriptor() << ": F-Factor photons to FADC counts smaller than 0." << endl;
685 return kFALSE;
686 }
687
688 const Float_t ffactorrelvar = 0.25 * fRSigmaSquareVar / ( fRSigmaSquare * fRSigmaSquare)
689 + GetMeanRelVar()
690 + 0.25 * nphotonsrelvar;
691
692 fMeanFFactorFADC2PhotVar = ffactorrelvar * fMeanFFactorFADC2Phot * fMeanFFactorFADC2Phot;
693
694 return kTRUE;
695}
696
697
698// ----------------------------------------------------------------------------
699//
700// - If fPed is smaller than 0 (i.e. has not yet been set), return.
701// - If fPedVar is smaller than 0 (i.e. has not yet been set), return.
702//
703// Calculate the electronic pedestal RMS with the formula:
704// - elec. pedestal = gkElectronicPedRms * sqrt(logainsamples)
705//
706// Calculate the night sky background ped. RMS contribution ("NSB") in the high-gain
707// from the high gain Pedestal RMS with the formula:
708// - HiGain NSB square = fPedRms * fPedRms - elec.ped.* elec.ped.
709// - Var(HiGain NSB square) = fPedVar * fPedRms * fPedRms + 4.*elecPedRmsVar * elec.ped.* elec.ped.
710//
711// If HiGain NSB square is smaller than 0., set it to zero. (but not the error!)
712//
713// Convert the NSB ped. RMS contribution to the low-gain with the formula:
714// - LoGain NSB square = - HiGain NSB square / (fConversionHiLo*fConversionHiLo)
715// - Var(LoGain NSB square) = ( Var(HiGain NSB square) / (HiGain NSB square * HiGain NSB square)
716// + GetConversionHiLoRelVar()
717// ) * LoGain NSB square * LoGain NSB square
718//
719// - Low Gain Ped RMS Square = LoGain NSB square + elec.ped. square
720// Var (Low Gain Ped RMS Square) = Var(LoGain NSB square) + Var(elec.ped. square)
721//
722void MCalibrationChargePix::CalcLoGainPedestal(Float_t logainsamples)
723{
724
725 if (fPedRms < 0.)
726 return;
727
728 if (fPedVar < 0.)
729 return;
730
731 const Float_t elecPedRms = gkElectronicPedRms * TMath::Sqrt(logainsamples);
732 const Float_t elecPedRmsVar = gkElectronicPedRmsErr * gkElectronicPedRmsErr * logainsamples;
733
734 Float_t pedRmsSquare = fPedRms * fPedRms;
735 Float_t pedRmsSquareVar = fPedVar * pedRmsSquare; // fPedRmsErr = fPedErr/2.
736
737 //
738 // We do not know the Lo Gain Pedestal RMS, so we have to retrieve it
739 // from the HI GAIN (all calculation per slice up to now):
740 //
741 // We extract the pure NSB contribution:
742 //
743 const Float_t elecRmsSquare = elecPedRms * elecPedRms;
744 const Float_t elecRmsSquareVar = 4.*elecPedRmsVar * elecRmsSquare;
745
746 Float_t higainNsbSquare = pedRmsSquare - elecRmsSquare;
747 Float_t higainNsbSquareRelVar = (pedRmsSquareVar + elecRmsSquareVar)
748 / (higainNsbSquare * higainNsbSquare) ;
749
750 if (higainNsbSquare < 0.)
751 higainNsbSquare = 0.;
752
753 //
754 // Now, we divide the NSB by the conversion factor and
755 // add it quadratically to the electronic noise
756 //
757 const Float_t conversionSquare = fConversionHiLo * fConversionHiLo;
758 const Float_t conversionSquareRelVar = 4.* GetConversionHiLoRelVar();
759
760 const Float_t logainNsbSquare = higainNsbSquare / conversionSquare;
761 const Float_t logainNsbSquareVar = ( higainNsbSquareRelVar + conversionSquareRelVar )
762 * logainNsbSquare * logainNsbSquare;
763
764 fLoGainPedRmsSquare = logainNsbSquare + elecRmsSquare;
765 fLoGainPedRmsSquareVar = logainNsbSquareVar + elecRmsSquareVar;
766}
767
Note: See TracBrowser for help on using the repository browser.