source: trunk/MagicSoft/Mars/mcalib/MCalibrationPix.cc@ 3149

Last change on this file since 3149 was 3147, checked in by gaug, 22 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 11/2003 <mailto:markus@ifae.es>
19!
20! Copyright: MAGIC Software Development, 2000-2001
21!
22!
23\* ======================================================================== */
24
25/////////////////////////////////////////////////////////////////////////////
26// //
27// MCalibrationPix //
28// //
29// This is the storage container to hold informations about the pedestal //
30// (offset) value 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/////////////////////////////////////////////////////////////////////////////
43#include "MCalibrationPix.h"
44#include "MCalibrationConfig.h"
45
46#include "MLog.h"
47#include "MLogManip.h"
48
49ClassImp(MCalibrationPix);
50
51using namespace std;
52
53const Float_t MCalibrationPix::gkElectronicPedRms = 1.5;
54const Float_t MCalibrationPix::gkErrElectronicPedRms = 0.3;
55const Float_t MCalibrationPix::gkFFactor = 1.15;
56const Float_t MCalibrationPix::gkFFactorError = 0.02;
57const Float_t MCalibrationPix::gkChargeLimit = 3.;
58const Float_t MCalibrationPix::gkChargeErrLimit = 0.;
59const Float_t MCalibrationPix::gkChargeRelErrLimit = 1.;
60const Float_t MCalibrationPix::gkTimeLimit = 1.5;
61const Float_t MCalibrationPix::gkTimeErrLimit = 3.;
62const Float_t MCalibrationPix::gkConvFFactorRelErrorLimit = 0.1;
63
64const Float_t MCalibrationPix::gkAverageQE = 0.20;
65const Float_t MCalibrationPix::gkAverageQEErr = 0.02;
66const Float_t MCalibrationPix::gkConversionHiLo = 10.;
67const Float_t MCalibrationPix::gkConversionHiLoError = 2.5;
68// --------------------------------------------------------------------------
69//
70// Default Constructor:
71//
72MCalibrationPix::MCalibrationPix(const char *name, const char *title)
73 : fPixId(-1),
74 fFlags(0)
75{
76
77 fName = name ? name : "MCalibrationPixel";
78 fTitle = title ? title : "Container of the MHCalibrationPixels and the fit results";
79
80 //
81 // At the moment, we don't have a database, yet,
82 // so we get it from the configuration file
83 //
84 fConversionHiLo = gkConversionHiLo;
85 fConversionHiLoError = gkConversionHiLoError;
86
87 fHist = new MHCalibrationPixel("MHCalibrationPixel","Calibration Histograms Pixel ");
88
89 if (!fHist)
90 *fLog << warn << dbginf << " Could not create MHCalibrationPixel " << endl;
91
92 Clear();
93
94 SetAverageQE();
95
96}
97
98MCalibrationPix::~MCalibrationPix()
99{
100 delete fHist;
101}
102
103
104// ------------------------------------------------------------------------
105//
106// Invalidate values
107//
108void MCalibrationPix::Clear(Option_t *o)
109{
110
111 fHist->Reset();
112
113 CLRBIT(fFlags, kHiGainSaturation);
114 CLRBIT(fFlags, kExcluded);
115 CLRBIT(fFlags, kExcludeQualityCheck);
116 CLRBIT(fFlags, kChargeValid);
117 CLRBIT(fFlags, kFitted);
118 CLRBIT(fFlags, kOscillating);
119 CLRBIT(fFlags, kBlindPixelMethodValid);
120 CLRBIT(fFlags, kFFactorMethodValid);
121 CLRBIT(fFlags, kPINDiodeMethodValid);
122 CLRBIT(fFlags, kCombinedMethodValid);
123
124 fCharge = -1.;
125 fErrCharge = -1.;
126 fSigmaCharge = -1.;
127 fErrSigmaCharge = -1.;
128 fRSigmaCharge = -1.;
129 fErrRSigmaCharge = -1.;
130
131 fChargeProb = -1.;
132 fPed = -1.;
133 fPedRms = -1.;
134 fErrPedRms = 0.;
135
136 fNumHiGainSamples = -1.;
137 fNumLoGainSamples = -1.;
138
139 fTimeFirstHiGain = 0 ;
140 fTimeLastHiGain = 0 ;
141 fTimeFirstLoGain = 0 ;
142 fTimeLastLoGain = 0 ;
143
144 fAbsTimeMean = -1.;
145 fAbsTimeRms = -1.;
146
147 fPheFFactorMethod = -1.;
148 fPheFFactorMethodError = -1.;
149
150 fMeanConversionFFactorMethod = -1.;
151 fMeanConversionBlindPixelMethod = -1.;
152 fMeanConversionPINDiodeMethod = -1.;
153 fMeanConversionCombinedMethod = -1.;
154
155 fErrorConversionFFactorMethod = -1.;
156 fErrorConversionBlindPixelMethod = -1.;
157 fErrorConversionPINDiodeMethod = -1.;
158 fErrorConversionCombinedMethod = -1.;
159
160 fSigmaConversionFFactorMethod = -1.;
161 fSigmaConversionBlindPixelMethod = -1.;
162 fSigmaConversionPINDiodeMethod = -1.;
163 fSigmaConversionCombinedMethod = -1.;
164
165 fFactorCalculated = kFALSE;
166
167}
168
169
170void MCalibrationPix::DefinePixId(Int_t i)
171{
172
173 fPixId = i;
174 fHist->ChangeHistId(i);
175
176}
177
178
179// --------------------------------------------------------------------------
180//
181// Set the pedestals from outside
182//
183void MCalibrationPix::SetPedestal(const Float_t ped, const Float_t pedrms,
184 const Float_t higainsamp, const Float_t logainsamp )
185{
186
187 fPed = ped;
188 fPedRms = pedrms;
189
190 fNumHiGainSamples = higainsamp;
191 fNumLoGainSamples = logainsamp;
192
193}
194
195// --------------------------------------------------------------------------
196//
197// Set the conversion factors from outside (only for MC)
198//
199void MCalibrationPix::SetConversionFFactorMethod(Float_t c, Float_t err, Float_t sig)
200{
201 fMeanConversionFFactorMethod = c;
202 fErrorConversionFFactorMethod = err;
203 fSigmaConversionFFactorMethod = sig;
204}
205
206// --------------------------------------------------------------------------
207//
208// Set the conversion factors from outside (only for MC)
209//
210void MCalibrationPix::SetConversionCombinedMethod(Float_t c, Float_t err, Float_t sig)
211{
212 fMeanConversionCombinedMethod = c;
213 fErrorConversionCombinedMethod = err;
214 fSigmaConversionCombinedMethod = sig;
215}
216
217
218// --------------------------------------------------------------------------
219//
220// Set the conversion factors from outside (only for MC)
221//
222void MCalibrationPix::SetConversionBlindPixelMethod(Float_t c, Float_t err, Float_t sig)
223{
224 fMeanConversionBlindPixelMethod = c;
225 fErrorConversionBlindPixelMethod = err;
226 fSigmaConversionBlindPixelMethod = sig;
227}
228
229// --------------------------------------------------------------------------
230//
231// Set the conversion factors from outside (only for MC)
232//
233void MCalibrationPix::SetConversionPINDiodeMethod(Float_t c, Float_t err, Float_t sig)
234{
235 fMeanConversionPINDiodeMethod = c ;
236 fErrorConversionPINDiodeMethod = err;
237 fSigmaConversionPINDiodeMethod = sig;
238}
239
240// --------------------------------------------------------------------------
241//
242// Set the Hi Gain Saturation Bit from outside (only for MC)
243//
244void MCalibrationPix::SetHiGainSaturation(Bool_t b)
245{
246
247 if (b)
248 {
249 SETBIT(fFlags, kHiGainSaturation);
250 fHist->SetUseLoGain(1);
251 }
252 else
253 {
254 CLRBIT(fFlags, kHiGainSaturation);
255 fHist->SetUseLoGain(0);
256 }
257}
258
259// --------------------------------------------------------------------------
260//
261// Set the Excluded Bit from outside
262//
263void MCalibrationPix::SetExcluded(Bool_t b )
264{
265 b ? SETBIT(fFlags, kExcluded) : CLRBIT(fFlags, kExcluded);
266}
267
268
269// --------------------------------------------------------------------------
270//
271// Set the Excluded Bit from outside
272//
273void MCalibrationPix::SetExcludeQualityCheck(Bool_t b )
274{
275 b ? SETBIT(fFlags, kExcludeQualityCheck) : CLRBIT(fFlags, kExcludeQualityCheck);
276}
277
278// --------------------------------------------------------------------------
279//
280// Set the Excluded Bit from outside
281//
282void MCalibrationPix::SetChargeValid(Bool_t b )
283{
284 b ? SETBIT(fFlags, kChargeValid) : CLRBIT(fFlags, kChargeValid);
285}
286
287// --------------------------------------------------------------------------
288//
289// Set the Excluded Bit from outside
290//
291void MCalibrationPix::SetFitted(Bool_t b )
292{
293 b ? SETBIT(fFlags, kFitted) : CLRBIT(fFlags, kFitted);
294}
295
296// --------------------------------------------------------------------------
297//
298// Set the Excluded Bit from outside
299//
300void MCalibrationPix::SetOscillating(Bool_t b )
301{
302 b ? SETBIT(fFlags, kOscillating) : CLRBIT(fFlags, kOscillating);
303}
304
305// --------------------------------------------------------------------------
306//
307// Set the Excluded Bit from outside
308//
309void MCalibrationPix::SetBlindPixelMethodValid(Bool_t b )
310{
311 b ? SETBIT(fFlags, kBlindPixelMethodValid) : CLRBIT(fFlags, kBlindPixelMethodValid);
312}
313
314// --------------------------------------------------------------------------
315//
316// Set the Excluded Bit from outside
317//
318void MCalibrationPix::SetFFactorMethodValid(Bool_t b )
319{
320 b ? SETBIT(fFlags, kFFactorMethodValid) : CLRBIT(fFlags, kFFactorMethodValid);
321}
322
323// --------------------------------------------------------------------------
324//
325// Set the Excluded Bit from outside
326//
327void MCalibrationPix::SetPINDiodeMethodValid(Bool_t b )
328{
329 b ? SETBIT(fFlags, kPINDiodeMethodValid) : CLRBIT(fFlags, kPINDiodeMethodValid);
330}
331
332void MCalibrationPix::SetAbsTimeBordersHiGain(Byte_t f, Byte_t l)
333{
334
335 fTimeFirstHiGain = f;
336 fTimeLastHiGain = l;
337
338}
339
340void MCalibrationPix::SetAbsTimeBordersLoGain(Byte_t f, Byte_t l)
341{
342
343 fTimeFirstLoGain = f;
344 fTimeLastLoGain = l;
345
346}
347
348Float_t MCalibrationPix::GetPheFFactorMethod()
349{
350
351 if (!fFactorCalculated)
352 CalcFFactorMethod();
353
354 return fPheFFactorMethod;
355
356}
357
358Float_t MCalibrationPix::GetPheFFactorMethodError()
359{
360
361 if (!fFactorCalculated)
362 CalcFFactorMethod();
363
364 return fPheFFactorMethodError;
365
366}
367
368
369Float_t MCalibrationPix::GetTotalFFactorFFactorMethod()
370{
371 if (!fFactorCalculated)
372 CalcFFactorMethod();
373
374 if (fPheFFactorMethod > 0)
375 return (fRSigmaCharge/fCharge)*TMath::Sqrt(fPheFFactorMethod);
376 else
377 return -1.;
378}
379
380Float_t MCalibrationPix::GetTotalFFactorErrorFFactorMethod()
381{
382
383 if (!fFactorCalculated)
384 CalcFFactorMethod();
385
386 const Float_t rsigmaChargeRelErrSquare = fErrRSigmaCharge * fErrRSigmaCharge
387 / (fRSigmaCharge * fRSigmaCharge) ;
388 const Float_t rChargeRelErrSquare = fErrCharge * fErrCharge
389 / (fCharge * fCharge) ;
390 const Float_t rPheRelErrSquare = fPheFFactorMethodError * fPheFFactorMethodError
391 / (fPheFFactorMethod * fPheFFactorMethod) ;
392
393 return TMath::Sqrt(rsigmaChargeRelErrSquare+rChargeRelErrSquare+rPheRelErrSquare);
394}
395
396
397Float_t MCalibrationPix::GetTotalFFactorBlindPixelMethod()
398{
399 return 1.;
400}
401
402Float_t MCalibrationPix::GetTotalFFactorErrorBlindPixelMethod()
403{
404
405 return 1.;
406}
407
408Float_t MCalibrationPix::GetTotalFFactorPINDiodeMethod()
409{
410 return 1.;
411}
412
413Float_t MCalibrationPix::GetTotalFFactorErrorPINDiodeMethod()
414{
415
416 return 1.;
417}
418
419Float_t MCalibrationPix::GetTotalFFactorCombinedMethod()
420{
421 return 1.;
422}
423
424Float_t MCalibrationPix::GetTotalFFactorErrorCombinedMethod()
425{
426
427 return 1.;
428}
429
430//
431// FIXME: This is a preliminary solution, the qe shall be
432// calibrated itself!
433//
434Float_t MCalibrationPix::GetMeanConversionFFactorMethod()
435{
436
437 if (!fFactorCalculated)
438 CalcFFactorMethod();
439
440 return fMeanConversionFFactorMethod/fAverageQE;
441
442}
443
444Float_t MCalibrationPix::GetErrorConversionFFactorMethod()
445{
446
447 if (!fFactorCalculated)
448 CalcFFactorMethod();
449
450 Float_t var = fErrorConversionFFactorMethod*fErrorConversionFFactorMethod
451 / (fMeanConversionFFactorMethod * fMeanConversionFFactorMethod);
452
453 var += fAverageQEErr * fAverageQEErr
454 / (fAverageQE * fAverageQE);
455
456 if (var > 0)
457 return -1.;
458
459 var = TMath::Sqrt(var);
460
461 return var*GetMeanConversionFFactorMethod();
462
463}
464
465Float_t MCalibrationPix::GetSigmaConversionFFactorMethod()
466{
467
468 if (!fFactorCalculated)
469 CalcFFactorMethod();
470
471 return fSigmaConversionFFactorMethod;
472
473}
474
475
476Bool_t MCalibrationPix::IsExcluded() const
477 {
478 return TESTBIT(fFlags,kExcluded);
479 }
480
481Bool_t MCalibrationPix::IsExcludeQualityCheck() const
482 {
483 return TESTBIT(fFlags,kExcludeQualityCheck);
484 }
485
486Bool_t MCalibrationPix::IsHiGainSaturation() const
487 {
488 return TESTBIT(fFlags,kHiGainSaturation);
489 }
490
491Bool_t MCalibrationPix::IsChargeValid() const
492{
493 return TESTBIT(fFlags, kChargeValid);
494}
495
496Bool_t MCalibrationPix::IsFitted() const
497{
498 return TESTBIT(fFlags, kFitted);
499}
500
501Bool_t MCalibrationPix::IsOscillating()
502{
503
504 if (TESTBIT(fFlags, kOscillating))
505 return kTRUE;
506
507 if (fHist->CheckOscillations())
508 {
509 SETBIT(fFlags,kOscillating);
510 return kTRUE;
511 }
512
513 return kFALSE;
514}
515
516Bool_t MCalibrationPix::IsBlindPixelMethodValid() const
517{
518 return TESTBIT(fFlags, kBlindPixelMethodValid);
519}
520
521Bool_t MCalibrationPix::IsFFactorMethodValid()
522{
523
524 if (!fFactorCalculated)
525 CalcFFactorMethod();
526
527 return TESTBIT(fFlags, kFFactorMethodValid);
528}
529
530Bool_t MCalibrationPix::IsPINDiodeMethodValid() const
531{
532 return TESTBIT(fFlags, kPINDiodeMethodValid);
533}
534
535Bool_t MCalibrationPix::IsCombinedMethodValid()
536{
537 return TESTBIT(fFlags, kCombinedMethodValid);
538}
539
540
541// --------------------------------------------------------------------------
542//
543// 1) Return if the charge distribution is already succesfully fitted
544// or if the histogram is empty
545// 2) Set a lower Fit range according to 1.5 Pedestal RMS in order to avoid
546// possible remaining cosmics to spoil the fit.
547// 3) Decide if the LoGain Histogram is fitted or the HiGain Histogram
548// 4) Fit the histograms with a Gaussian
549// 5) In case of failure set the bit kFitted to false
550// 6) Retrieve the results and store them in this class
551// 7) Calculate the number of photo-electrons after the F-Factor method
552// 8) Calculate the errors of the F-Factor method
553//
554// The fits are declared valid (fFitValid = kTRUE), if:
555//
556// 1) Pixel has a fitted charge greater than 3*PedRMS
557// 2) Pixel has a fit error greater than 0.
558// 3) Pixel has a fit Probability greater than 0.0001
559// 4) Pixel has a charge sigma bigger than its Pedestal RMS
560// 5) If FitTimes is used,
561// the mean arrival time is at least 1.0 slices from the used edge slices
562// (this stage is only performed in the times fit)
563//
564// If the histogram is empty, all values are set to -1.
565//
566// The conversion factor after the F-Factor method is declared valid, if:
567//
568// 1) fFitValid is kTRUE
569// 2) Conversion Factor is bigger than 0.
570// 3) The error of the conversion factor is smaller than 10%
571//
572Bool_t MCalibrationPix::FitCharge()
573{
574
575 //
576 // 1) Return if the charge distribution is already succesfully fitted
577 // or if the histogram is empty
578 //
579 if (fHist->IsChargeFitOK() || fHist->IsEmpty())
580 return kTRUE;
581
582 //
583 // 2) Set a lower Fit range according to 1.5 Pedestal RMS in order to avoid
584 // possible remaining cosmics to spoil the fit.
585 //
586 // if (fPed && fPedRms)
587 // fHist->SetLowerFitRange(1.5*fPedRms);
588 // else
589 // *fLog << warn << "WARNING: Cannot set lower fit range: Pedestals not available" << endl;
590
591 //
592 // 3) Decide if the LoGain Histogram is fitted or the HiGain Histogram
593 //
594 if (fHist->UseLoGain())
595 SetHiGainSaturation();
596
597 //
598 // 4) Fit the Lo Gain histograms with a Gaussian
599 //
600 if (fHist->FitCharge())
601 SETBIT(fFlags,kFitted);
602 else
603 {
604 *fLog << warn << "WARNING: Could not fit charges of pixel " << fPixId << endl;
605 //
606 // 5) In case of failure set the bit kFitted to false
607 //
608 CLRBIT(fFlags,kFitted);
609 }
610
611 //
612 // 6) Retrieve the results and store them in this class
613 // If fFitted is false, we get the means and RMS of the histogram!!
614 //
615 fCharge = fHist->GetChargeMean();
616 fErrCharge = fHist->GetChargeMeanErr();
617 fSigmaCharge = fHist->GetChargeSigma();
618 fErrSigmaCharge = fHist->GetChargeSigmaErr();
619 fChargeProb = fHist->GetChargeProb();
620
621
622 fAbsTimeMean = fHist->GetAbsTimeMean();
623 fAbsTimeMeanErr = fHist->GetAbsTimeMeanErr();
624 fAbsTimeRms = fHist->GetAbsTimeRms();
625
626 if (CheckTimeFitValidity())
627 SETBIT(fFlags,kTimeFitValid);
628 else
629 CLRBIT(fFlags,kTimeFitValid);
630
631 //
632 // Calculate the conversion factors
633 //
634 if (IsHiGainSaturation())
635 ApplyLoGainConversion();
636
637 if (CheckChargeValidity())
638 SETBIT(fFlags,kChargeValid);
639 else
640 {
641 CLRBIT(fFlags,kChargeValid);
642 return kFALSE;
643 }
644
645 return kTRUE;
646
647}
648
649//
650// Calculate the number of photo-electrons after the F-Factor method
651// Calculate the errors of the F-Factor method
652//
653Bool_t MCalibrationPix::CalcFFactorMethod()
654{
655
656 if ( (fCharge == -1.)
657 || (fErrCharge < 0.)
658 || (fSigmaCharge < 0.)
659 || (fPedRms < 0.) )
660 {
661 *fLog << warn << GetDescriptor() << "Cannot calculate the FFactor Method! "
662 << "Some of the needed parameters are not available ";
663 CLRBIT(fFlags,kFFactorMethodValid);
664 return kFALSE;
665 }
666
667 //
668 // Square all variables in order to avoid applications of square root
669 //
670 // First the relative error squares
671 //
672 const Float_t chargeSquare = fCharge* fCharge;
673 const Float_t chargeSquareRelErrSquare = 4.*fErrCharge*fErrCharge / chargeSquare;
674
675 const Float_t ffactorsquare = gkFFactor * gkFFactor;
676 const Float_t ffactorsquareRelErrSquare = 4.*gkFFactorError * gkFFactorError / ffactorsquare;
677 //
678 // Now the absolute error squares
679 //
680 const Float_t sigmaSquare = fSigmaCharge*fSigmaCharge;
681 const Float_t sigmaSquareErrSquare = 4.*fErrSigmaCharge*fErrSigmaCharge * sigmaSquare;
682
683 Float_t pedRmsSquare = fPedRms* fPedRms;
684 Float_t pedRmsSquareErrSquare = 4.*fErrPedRms*fErrPedRms * pedRmsSquare;
685
686 if (!IsHiGainSaturation())
687 { /* HiGain */
688
689 pedRmsSquare *= fNumHiGainSamples;
690 pedRmsSquareErrSquare *= fNumHiGainSamples*fNumHiGainSamples;
691 }
692 else
693 { /* LoGain */
694
695 //
696 // We do not know the Lo Gain Pedestal RMS, so we have to retrieve it
697 // from the HI GAIN (all calculation per slice up to now):
698 //
699 // We extract the pure NSB contribution:
700 //
701 const Float_t elecRmsSquare = gkElectronicPedRms *gkElectronicPedRms;
702 const Float_t elecRmsSquareErrSquare = 4.*gkErrElectronicPedRms*gkErrElectronicPedRms * elecRmsSquare;
703
704 Float_t nsbSquare = pedRmsSquare - elecRmsSquare;
705 Float_t nsbSquareRelErrSquare = (pedRmsSquareErrSquare + elecRmsSquareErrSquare)
706 / (nsbSquare * nsbSquare) ;
707
708 if (nsbSquare < 0.)
709 nsbSquare = 0.;
710
711 //
712 // Now, we divide the NSB by the conversion factor and
713 // add it quadratically to the electronic noise
714 //
715 const Float_t conversionSquare = fConversionHiLo *fConversionHiLo;
716 const Float_t conversionSquareRelErrSquare = 4.*fConversionHiLoError*fConversionHiLoError/conversionSquare;
717
718 const Float_t convertedNsbSquare = nsbSquare / conversionSquare;
719 const Float_t convertedNsbSquareErrSquare = (nsbSquareRelErrSquare + conversionSquareRelErrSquare)
720 * convertedNsbSquare * convertedNsbSquare;
721
722 pedRmsSquare = convertedNsbSquare + elecRmsSquare;
723 pedRmsSquareErrSquare = convertedNsbSquareErrSquare + elecRmsSquareErrSquare;
724
725 //
726 // Now, correct for the number of used FADC slices in the LoGain:
727 //
728 pedRmsSquare *= fNumLoGainSamples;
729 pedRmsSquareErrSquare *= fNumLoGainSamples*fNumLoGainSamples;
730 //
731 // Correct also for the conversion to Hi-Gain:
732 //
733 pedRmsSquare *= fConversionHiLo*fConversionHiLo;
734 pedRmsSquareErrSquare *= fConversionHiLo*fConversionHiLo*fConversionHiLo*fConversionHiLo;
735
736 } /* if (HiGainSaturation) */
737
738 //
739 // Calculate the reduced sigmas
740 //
741 const Float_t rsigmachargesquare = sigmaSquare - pedRmsSquare;
742 if (rsigmachargesquare <= 0.)
743 {
744 *fLog << warn
745 << "WARNING: Cannot apply F-Factor calibration: Reduced Sigma smaller than 0 in pixel "
746 << fPixId << endl;
747 CLRBIT(fFlags,kFFactorMethodValid);
748 fFactorCalculated = kTRUE;
749 return kFALSE;
750 }
751
752 const Float_t rSigmaSquareRelErrSquare = (sigmaSquareErrSquare + pedRmsSquareErrSquare)
753 / (rsigmachargesquare * rsigmachargesquare) ;
754
755 fRSigmaCharge = TMath::Sqrt(rsigmachargesquare);
756 fErrRSigmaCharge = TMath::Sqrt(sigmaSquareErrSquare + pedRmsSquareErrSquare);
757
758
759 //
760 // Calculate the number of phe's from the F-Factor method
761 // (independent on Hi Gain or Lo Gain)
762 //
763 fPheFFactorMethod = ffactorsquare * chargeSquare / rsigmachargesquare;
764
765 const Float_t pheFFactorRelErrSquare = ffactorsquareRelErrSquare
766 + chargeSquareRelErrSquare
767 + rSigmaSquareRelErrSquare ;
768
769 fPheFFactorMethodError = TMath::Sqrt(pheFFactorRelErrSquare) * fPheFFactorMethod;
770
771 const Float_t chargeRelErrSquare = fErrCharge*fErrCharge / (fCharge * fCharge);
772
773 fMeanConversionFFactorMethod = fPheFFactorMethod / fCharge ;
774 fErrorConversionFFactorMethod = ( pheFFactorRelErrSquare + chargeRelErrSquare )
775 * fMeanConversionFFactorMethod * fMeanConversionFFactorMethod;
776
777 const Float_t convrelerror = fErrorConversionFFactorMethod/fMeanConversionFFactorMethod;
778
779 if ( (fMeanConversionFFactorMethod > 0.) && (convrelerror < gkConvFFactorRelErrorLimit))
780 SETBIT(fFlags,kFFactorMethodValid);
781
782 fFactorCalculated = kTRUE;
783
784 fSigmaConversionFFactorMethod = GetTotalFFactorFFactorMethod()*TMath::Sqrt(fMeanConversionFFactorMethod);
785
786 return kTRUE;
787}
788
789
790//
791// The check returns kTRUE if:
792//
793// 0) Pixel has BIT fitted set:
794// This means:
795// a) No result is a nan
796// b) The NDF is not smaller than fNDFLimit (5)
797// c) The Probability is greater than gkProbLimit (default 0.001 == 99.9%)
798// 1) Pixel has a fitted charge greater than 3*PedRMS
799// 2) Pixel has a fit error greater than 0.
800// 3) Pixel has a fitted charge greater its charge error
801// 4) Pixel has a fit Probability greater than 0.0001
802// 5) Pixel has a charge sigma bigger than its Pedestal RMS
803//
804Bool_t MCalibrationPix::CheckChargeValidity()
805{
806
807 if (!IsFitted())
808 return kFALSE;
809
810 if (IsExcludeQualityCheck())
811 return kTRUE;
812
813 Float_t pedestal;
814
815 if (!IsHiGainSaturation()) /* higain */
816 pedestal = GetPedRms()*TMath::Sqrt(fNumHiGainSamples);
817 else /* logain */
818 pedestal = GetPedRms()*TMath::Sqrt(fNumLoGainSamples);
819
820
821 if (fCharge < gkChargeLimit*pedestal)
822 {
823 *fLog << warn << "WARNING: Fitted Charge is smaller than "
824 << gkChargeLimit << " Pedestal RMS in Pixel " << fPixId << endl;
825 return kFALSE;
826 }
827
828 if (fErrCharge < gkChargeErrLimit)
829 {
830 *fLog << warn << "WARNING: Error of Fitted Charge is smaller than "
831 << gkChargeErrLimit << " in Pixel " << fPixId << endl;
832 return kFALSE;
833 }
834
835 if (fCharge < gkChargeRelErrLimit*fErrCharge)
836 {
837 *fLog << warn << "WARNING: Fitted Charge is smaller than "
838 << gkChargeRelErrLimit << "* its error in Pixel " << fPixId << endl;
839 return kFALSE;
840 }
841
842 if (!fHist->IsChargeFitOK())
843 {
844 *fLog << warn << "WARNING: Probability of Fitted Charge too low in Pixel "
845 << fPixId << endl;
846 return kFALSE;
847 }
848
849 if (fSigmaCharge < pedestal)
850 {
851 *fLog << warn << "WARNING: Sigma of Fitted Charge smaller than Pedestal RMS in Pixel "
852 << fPixId << endl;
853 return kFALSE;
854 }
855 return kTRUE;
856}
857
858//
859// The check return kTRUE if:
860//
861// 0) No value is nan
862// 1) Pixel has a fitted rel. time smaller than 3*FADC slices
863// 2) Pixel has a fit error greater than 0.
864// 4) Pixel has a fit Probability greater than 0.001
865// 5) The absolute arrival time is at least 1.0 slices from the used edge slices
866//
867Bool_t MCalibrationPix::CheckTimeFitValidity()
868{
869
870
871 if (IsExcludeQualityCheck())
872 return kTRUE;
873
874 if (IsHiGainSaturation())
875 {
876
877 if (fAbsTimeMean < (Float_t)fTimeFirstLoGain+1)
878 {
879 *fLog << warn
880 << "WARNING: Some absolute times smaller than limit in Pixel "
881 << fPixId << " time: " << fAbsTimeMean
882 << " Limit: " << (Float_t)fTimeFirstLoGain+1. << endl;
883 return kFALSE;
884 }
885
886 if (fAbsTimeMean > (Float_t)fTimeLastLoGain-1)
887 {
888 *fLog << warn
889 << "WARNING: Some absolute times bigger than limit in Pixel "
890 << fPixId << " time: " << fAbsTimeMean
891 << " Limit: " << (Float_t)fTimeLastLoGain-1. << endl;
892 return kFALSE;
893 }
894
895 }
896 else
897 {
898
899 if (fAbsTimeMean < (Float_t)fTimeFirstHiGain+1.)
900 {
901 *fLog << warn
902 << "WARNING: Some absolute times smaller than limit in Pixel "
903 << fPixId << " time: " << fAbsTimeMean
904 << " Limit: " << (Float_t)fTimeFirstHiGain+1. << endl;
905 // return kFALSE;
906 }
907
908 if (fAbsTimeMean > (Float_t)fTimeLastHiGain-1.)
909 {
910 *fLog << warn
911 << "WARNING: Some absolute times bigger than limit in Pixel "
912 << fPixId << " time: " << fAbsTimeMean
913 << " Limit: " << (Float_t)fTimeLastHiGain-1. << endl;
914 // return kFALSE;
915 }
916
917 }
918
919
920
921 return kTRUE;
922}
923
924
925void MCalibrationPix::CheckOscillations()
926{
927 fHist->CheckOscillations();
928}
929
930void MCalibrationPix::ApplyLoGainConversion()
931{
932
933 const Float_t chargeRelErrSquare = fErrCharge*fErrCharge
934 /( fCharge * fCharge);
935 const Float_t sigmaRelErrSquare = fErrSigmaCharge*fErrSigmaCharge
936 /( fSigmaCharge * fSigmaCharge);
937 const Float_t conversionRelErrSquare = fConversionHiLoError*fConversionHiLoError
938 /(fConversionHiLo * fConversionHiLo);
939
940 fCharge *= fConversionHiLo;
941 fErrCharge = TMath::Sqrt(chargeRelErrSquare + conversionRelErrSquare) * fCharge;
942
943 fSigmaCharge *= fConversionHiLo;
944 fErrSigmaCharge = TMath::Sqrt(sigmaRelErrSquare + conversionRelErrSquare) * fSigmaCharge;
945
946}
Note: See TracBrowser for help on using the repository browser.