source: trunk/MagicSoft/Mars/mcalib/MCalibrationBlindPix.cc@ 8462

Last change on this file since 8462 was 8106, checked in by tbretz, 19 years ago
*** empty log message ***
File size: 16.3 KB
Line 
1/* ======================================================================== *\
2! $Name: not supported by cvs2svn $:$Id: MCalibrationBlindPix.cc,v 1.15 2006-10-17 17:15:59 tbretz Exp $
3! --------------------------------------------------------------------------
4!
5! *
6! * This file is part of MARS, the MAGIC Analysis and Reconstruction
7! * Software. It is distributed to you in the hope that it can be a useful
8! * and timesaving tool in analysing Data of imaging Cerenkov telescopes.
9! * It is distributed WITHOUT ANY WARRANTY.
10! *
11! * Permission to use, copy, modify and distribute this software and its
12! * documentation for any purpose is hereby granted without fee,
13! * provided that the above copyright notice appear in all copies and
14! * that both that copyright notice and this permission notice appear
15! * in supporting documentation. It is provided "as is" without express
16! * or implied warranty.
17! *
18!
19!
20! Author(s): Markus Gaug 02/2004 <mailto:markus@ifae.es>
21!
22! Copyright: MAGIC Software Development, 2000-2004
23!
24!
25\* ======================================================================== */
26
27/////////////////////////////////////////////////////////////////////////////
28//
29// MCalibrationBlindPix
30//
31// Storage container of the fit results of the Blind Pixel signal
32// (from MHCalibrationChargeBlindPix).
33//
34// The Flux is calculated in photons per mm^2 in the camera plane.
35//
36// Currently, the following numbers are implemented:
37// - fArea: 100 mm^2
38// - Average QE of Blind Pixel:
39// fQEGreen: 0.154
40// fQEBlue : 0.226
41// fQEUV : 0.247
42// fQECT1 : 0.247
43// - Average QE Error of Blind Pixel:
44// fQEGreenErr: 0.015;
45// fQEBlueErr : 0.02;
46// fQEUVErr : 0.02;
47// fQECT1Err : 0.02;
48// - Attenuation factor Blind Pixel:
49// fAttGreen : 1.97;
50// fAttBlue : 1.96;
51// fAttUV : 1.95;
52// fAttCT1 : 1.95;
53//
54//
55/////////////////////////////////////////////////////////////////////////////
56#include "MCalibrationBlindPix.h"
57
58#include <TH1.h>
59
60#include "MLog.h"
61#include "MLogManip.h"
62
63ClassImp(MCalibrationBlindPix);
64
65using namespace std;
66
67const Float_t MCalibrationBlindPix::fgArea = 100;
68const Float_t MCalibrationBlindPix::fgAttGreen = 1.97;
69const Float_t MCalibrationBlindPix::fgAttBlue = 1.96;
70const Float_t MCalibrationBlindPix::fgAttUV = 1.95;
71const Float_t MCalibrationBlindPix::fgAttCT1 = 1.95;
72const Float_t MCalibrationBlindPix::fgAttErr = 0.01;
73const Float_t MCalibrationBlindPix::fgQEGreen = 0.154;
74const Float_t MCalibrationBlindPix::fgQEBlue = 0.226;
75const Float_t MCalibrationBlindPix::fgQEUV = 0.247;
76const Float_t MCalibrationBlindPix::fgQECT1 = 0.247;
77const Float_t MCalibrationBlindPix::fgQEErrGreen = 0.005;
78const Float_t MCalibrationBlindPix::fgQEErrBlue = 0.007;
79const Float_t MCalibrationBlindPix::fgQEErrUV = 0.01;
80const Float_t MCalibrationBlindPix::fgQEErrCT1 = 0.01;
81const Float_t MCalibrationBlindPix::fgCollEffGreen = 0.99;
82const Float_t MCalibrationBlindPix::fgCollEffBlue = 0.93;
83const Float_t MCalibrationBlindPix::fgCollEffUV = 0.90;
84const Float_t MCalibrationBlindPix::fgCollEffCT1 = 0.90;
85const Float_t MCalibrationBlindPix::fgCollEffErr = 0.05;
86
87// --------------------------------------------------------------------------
88//
89// Default Constructor.
90//
91// Calls:
92// - Clear()
93//
94// For backward-compatibility reasons, quantum eff., coll. eff. and att.
95// are intialized from the static members. This should, however, be
96// overwritten by a class deriving from MCalibrationBlindCam.
97//
98MCalibrationBlindPix::MCalibrationBlindPix(const char *name, const char *title)
99{
100
101 fName = name ? name : "MCalibrationBlindPix";
102 fTitle = title ? title : "Container of the fit results of the blind pixel";
103
104 Clear();
105
106 fArea = fgArea;
107 fAreaErr = 0.;
108
109 Float_t att[MCalibrationCam::gkNumPulserColors];
110
111 att [ MCalibrationCam::kGREEN ] = fgAttGreen;
112 att [ MCalibrationCam::kBLUE ] = fgAttBlue;
113 att [ MCalibrationCam::kUV ] = fgAttUV;
114 att [ MCalibrationCam::kCT1 ] = fgAttCT1;
115
116 SetAtt(MCalibrationCam::gkNumPulserColors,att);
117
118 Float_t atterr[MCalibrationCam::gkNumPulserColors];
119
120 atterr [ MCalibrationCam::kGREEN ] = fgAttErr;
121 atterr [ MCalibrationCam::kBLUE ] = fgAttErr;
122 atterr [ MCalibrationCam::kUV ] = fgAttErr;
123 atterr [ MCalibrationCam::kCT1 ] = fgAttErr;
124
125 SetAttErr(MCalibrationCam::gkNumPulserColors,atterr);
126
127 Float_t qe[MCalibrationCam::gkNumPulserColors];
128
129 qe [ MCalibrationCam::kGREEN ] = fgQEGreen;
130 qe [ MCalibrationCam::kBLUE ] = fgQEBlue;
131 qe [ MCalibrationCam::kUV ] = fgQEUV;
132 qe [ MCalibrationCam::kCT1 ] = fgQECT1;
133
134 SetQE(MCalibrationCam::gkNumPulserColors,qe);
135
136 Float_t qeerr[MCalibrationCam::gkNumPulserColors];
137
138 qeerr [ MCalibrationCam::kGREEN ] = fgQEErrGreen;
139 qeerr [ MCalibrationCam::kBLUE ] = fgQEErrBlue;
140 qeerr [ MCalibrationCam::kUV ] = fgQEErrUV;
141 qeerr [ MCalibrationCam::kCT1 ] = fgQEErrCT1;
142
143 SetQEErr(MCalibrationCam::gkNumPulserColors,qeerr);
144
145 Float_t colleff[MCalibrationCam::gkNumPulserColors];
146
147 colleff [ MCalibrationCam::kGREEN ] = fgCollEffGreen;
148 colleff [ MCalibrationCam::kBLUE ] = fgCollEffBlue;
149 colleff [ MCalibrationCam::kUV ] = fgCollEffUV;
150 colleff [ MCalibrationCam::kCT1 ] = fgCollEffCT1;
151
152 SetCollEff(MCalibrationCam::gkNumPulserColors,colleff);
153
154 Float_t collefferr[MCalibrationCam::gkNumPulserColors];
155
156 collefferr [ MCalibrationCam::kGREEN ] = fgCollEffErr;
157 collefferr [ MCalibrationCam::kBLUE ] = fgCollEffErr;
158 collefferr [ MCalibrationCam::kUV ] = fgCollEffErr;
159 collefferr [ MCalibrationCam::kCT1 ] = fgCollEffErr;
160
161 SetCollEffErr(MCalibrationCam::gkNumPulserColors,collefferr);
162}
163
164// ------------------------------------------------------------------------
165//
166// Sets:
167// - all flags to kFALSE
168// - all variables to -1.
169// - the fColor to MCalibrationCam::kNONE
170//
171// Calls:
172// - MCalibrationPix::Clear()
173//
174void MCalibrationBlindPix::Clear(Option_t *o)
175{
176
177 fFluxInsidePlexiglass = -1.;
178 fFluxInsidePlexiglassVar = -1.;
179 fLambda = -1.;
180 fLambdaCheck = -1.;
181 fLambdaVar = -1.;
182 fMu0 = -1.;
183 fMu0Err = -1.;
184 fMu1 = -1.;
185 fMu1Err = -1.;
186 fSigma0 = -1.;
187 fSigma0Err = -1.;
188 fSigma1 = -1.;
189 fSigma1Err = -1.;
190
191 SetOscillating ( kFALSE );
192 SetExcluded ( kFALSE );
193 SetChargeFitValid ( kFALSE );
194 SetPedestalFitOK ( kFALSE );
195 SetSinglePheFitOK ( kFALSE );
196 SetFluxInsidePlexiglassAvailable ( kFALSE );
197
198 SetColor(MCalibrationCam::kNONE);
199
200 MCalibrationPix::Clear();
201}
202
203void MCalibrationBlindPix::SetFluxInsidePlexiglassAvailable( const Bool_t b)
204{
205 b ? SETBIT(fFlags,kFluxInsidePlexiglassAvailable) : CLRBIT(fFlags,kFluxInsidePlexiglassAvailable);
206}
207
208
209// --------------------------------------------------------------------------
210//
211// Set the Oscillating Bit from outside
212//
213void MCalibrationBlindPix::SetOscillating( const Bool_t b)
214{
215 b ? SETBIT(fFlags,kOscillating) : CLRBIT(fFlags,kOscillating);
216}
217
218// -----------------------------------------------------
219//
220// copy 'constructor'
221//
222void MCalibrationBlindPix::Copy(TObject& object) const
223{
224
225 MCalibrationBlindPix &pix = (MCalibrationBlindPix&)object;
226
227 //
228 // Copy the data members
229 //
230 pix.fPixId = fPixId;
231 pix.fFlags = fFlags;
232 pix.fArea = fArea;
233 pix.fAreaErr = fAreaErr;
234 pix.fColor = fColor;
235
236 //
237 // Copy arrays
238 //
239 pix.fAtt = fAtt;
240 pix.fAttErr = fAttErr;
241 pix.fQE = fQE;
242 pix.fQEErr = fQEErr;
243 pix.fCollEff = fCollEff;
244 pix.fCollEffErr = fCollEffErr;
245
246 pix.fLambda = fLambda;
247 pix.fLambdaCheck = fLambdaCheck;
248 pix.fLambdaCheckErr = fLambdaCheckErr;
249 pix.fLambdaVar = fLambdaVar;
250 pix.fFluxInsidePlexiglass = fFluxInsidePlexiglass;
251 pix.fFluxInsidePlexiglassVar = fFluxInsidePlexiglassVar;
252 pix.fMu0 = fMu0;
253 pix.fMu0Err = fMu0Err;
254 pix.fMu1 = fMu1;
255 pix.fMu1Err = fMu1Err;
256 pix.fSigma0 = fSigma0;
257 pix.fSigma0Err = fSigma0Err;
258 pix.fSigma1 = fSigma1;
259 pix.fSigma1Err = fSigma1Err;
260
261}
262
263
264// --------------------------------------------------------------------------
265//
266// Set the ChargeFitValid Bit from outside
267//
268void MCalibrationBlindPix::SetChargeFitValid( const Bool_t b)
269{
270 b ? SETBIT(fFlags,kChargeFitValid) : CLRBIT(fFlags,kChargeFitValid);
271}
272
273// --------------------------------------------------------------------------
274//
275// Set the PedestalFitValid Bit from outside
276//
277void MCalibrationBlindPix::SetPedestalFitOK( const Bool_t b)
278{
279 b ? SETBIT(fFlags,kPedestalFitOK) : CLRBIT(fFlags,kPedestalFitOK);
280}
281
282// --------------------------------------------------------------------------
283//
284// Set the SinglePheFitValid Bit from outside
285//
286void MCalibrationBlindPix::SetSinglePheFitOK( const Bool_t b)
287{
288 b ? SETBIT(fFlags,kSinglePheFitOK) : CLRBIT(fFlags,kSinglePheFitOK);
289}
290
291// --------------------------------------------------------------------------
292//
293// Return -1 if fFluxInsidePlexiglassVar is smaller than 0.
294// Return square root of fFluxInsidePlexiglassVar
295//
296const Float_t MCalibrationBlindPix::GetFluxInsidePlexiglassErr() const
297{
298 if (fFluxInsidePlexiglassVar < 0.)
299 return -1.;
300
301 return TMath::Sqrt(fFluxInsidePlexiglassVar);
302}
303
304// --------------------------------------------------------------------------
305//
306// Return -1 if fFluxInsidePlexiglassVar is smaller than 0.
307// Return -1 if fFluxInsidePlexiglass is 0.
308// Return fFluxInsidePlexiglassVar / fFluxInsidePlexiglass^2
309//
310const Float_t MCalibrationBlindPix::GetFluxInsidePlexiglassRelVar() const
311{
312 if (fFluxInsidePlexiglassVar < 0.)
313 return -1.;
314
315 if (fFluxInsidePlexiglass == 0.)
316 return -1.;
317
318 return fFluxInsidePlexiglassVar / (fFluxInsidePlexiglass * fFluxInsidePlexiglass) ;
319}
320
321// --------------------------------------------------------------------------
322//
323// Return -1 if fLambdaVar is smaller than 0.
324// Return square root of fLambdaVar
325//
326const Float_t MCalibrationBlindPix::GetLambdaErr() const
327{
328 if (fLambdaVar < 0.)
329 return -1.;
330
331 return TMath::Sqrt(fLambdaVar);
332}
333
334// --------------------------------------------------------------------------
335//
336// Return -1 if fLambdaVar is smaller than 0.
337// Return -1 if fLambda is 0.
338// Return fLambdaVar / (fLambda * fLambda )
339//
340const Float_t MCalibrationBlindPix::GetLambdaRelVar() const
341{
342 if (fLambdaVar < 0.)
343 return -1.;
344
345 if (fLambda == 0.)
346 return -1.;
347
348 return fLambdaVar / fLambda / fLambda ;
349}
350
351// --------------------------------------------------------------------------
352//
353// Return fAreaErr^2 / (fArea^2 )
354//
355const Float_t MCalibrationBlindPix::GetAreaRelVar() const
356{
357 return fAreaErr*fAreaErr/fArea/fArea;
358}
359
360
361// --------------------------------------------------------------------------
362//
363// Return TMath::Power(10,fAtt[fColor])
364//
365const Float_t MCalibrationBlindPix::GetAtt() const
366{
367 return TMath::Power(10,fAtt[fColor]);
368}
369
370// --------------------------------------------------------------------------
371//
372// Return -1 if fAttErr[fColor] is smaller than 0.
373// Error of TMath::Power(10,fAtt[fColor]) = TMath::Power(10,fAtt[fColor])*ln(10.)*fAttErr[fColor]
374// Return fAttErr^2 / (fAtt^2 )
375//
376const Float_t MCalibrationBlindPix::GetAttRelVar() const
377{
378
379 const Float_t ferr = fAttErr[fColor];
380
381 if (ferr < 0.)
382 return -1.;
383
384 return ferr*ferr*2.3;
385}
386
387// --------------------------------------------------------------------------
388//
389// Return fQE[fColor]
390//
391const Float_t MCalibrationBlindPix::GetQE() const
392{
393 return fQE[fColor];
394}
395
396// --------------------------------------------------------------------------
397//
398// Return -1 if fQEErr[fColor] is smaller than 0.
399// Return fQEErr^2 / (fQE^2 )
400//
401const Float_t MCalibrationBlindPix::GetQERelVar() const
402{
403
404 if (fQEErr[fColor] < 0.)
405 return -1.;
406
407 return fQEErr[fColor]* fQEErr[fColor] / GetQE() / GetQE();
408}
409
410// --------------------------------------------------------------------------
411//
412// Return fCollEff[fColor]
413//
414const Float_t MCalibrationBlindPix::GetCollEff() const
415{
416 return fCollEff[fColor];
417}
418
419// --------------------------------------------------------------------------
420//
421// Return -1 if fCollEffErr[fColor] is smaller than 0.
422// Return fCollEffErr^2 / (fCollEff^2 )
423//
424const Float_t MCalibrationBlindPix::GetCollEffRelVar() const
425{
426
427 if (fCollEffErr[fColor] < 0.)
428 return -1.;
429
430 return fCollEffErr[fColor]* fCollEffErr[fColor] / GetCollEff() / GetCollEff();
431}
432
433// --------------------------------------------------------------------------
434//
435// Test bit kChargeFitValid
436//
437const Bool_t MCalibrationBlindPix::IsChargeFitValid() const
438{
439 return TESTBIT(fFlags,kChargeFitValid);
440}
441
442// --------------------------------------------------------------------------
443//
444// Test bit kOscillating
445//
446const Bool_t MCalibrationBlindPix::IsOscillating() const
447{
448 return TESTBIT(fFlags,kOscillating);
449}
450
451// --------------------------------------------------------------------------
452//
453// Test bit kPedestalFitValid
454//
455const Bool_t MCalibrationBlindPix::IsPedestalFitOK() const
456{
457 return TESTBIT(fFlags,kPedestalFitOK);
458}
459
460// --------------------------------------------------------------------------
461//
462// Test bit kSinglePheFitValid
463//
464const Bool_t MCalibrationBlindPix::IsSinglePheFitOK() const
465{
466 return TESTBIT(fFlags,kSinglePheFitOK);
467}
468
469// --------------------------------------------------------------------------
470//
471// Test bit kFluxInsidePlexiglassAvailable
472//
473const Bool_t MCalibrationBlindPix::IsFluxInsidePlexiglassAvailable() const
474{
475 return TESTBIT(fFlags,kFluxInsidePlexiglassAvailable);
476}
477
478
479// --------------------------------------------------------------------------
480//
481// Return kFALSE if IsChargeFitValid() is kFALSE
482//
483// Calculate fFluxInsidePlexiglass with the formula:
484// - fFluxInsidePlexiglass = fLambda
485// / GetCollEff()
486// / GetQE()
487// * GetAtt()
488// / fArea
489// - fFluxInsidePlexiglassVar = sqrt( fLambdaVar / ( fLambda * fLambda )
490// + GetQERelVar()
491// + GetCollEffRelVar()
492// + GetAttRelVar()
493// ) * fFluxInsidePlexiglass * * fFluxInsidePlexiglass
494//
495// If the fFluxInsidePlexiglass is smaller than 0., return kFALSE
496// If the Variance is smaller than 0., return kFALSE
497//
498// SetFluxInsidePlexiglassAvailable() and return kTRUE
499//
500Bool_t MCalibrationBlindPix::CalcFluxInsidePlexiglass()
501{
502
503 if (IsChargeFitValid())
504 return kFALSE;
505
506
507 //
508 // Start calculation of number of photons
509 // The blind pixel has exactly 100 mm^2 area (with negligible error),
510 //
511 fFluxInsidePlexiglass = fLambda / GetQE() * GetAtt() / GetCollEff() / fArea;
512
513 if (fFluxInsidePlexiglass < 0.)
514 return kFALSE;
515
516 fFluxInsidePlexiglassVar = GetLambdaRelVar() + GetQERelVar() + GetAttRelVar() + GetCollEffRelVar() + GetAreaRelVar();
517
518 //
519 // Finish calculation of errors -> convert from relative variance to absolute variance
520 //
521 fFluxInsidePlexiglassVar *= fFluxInsidePlexiglass * fFluxInsidePlexiglass;
522
523 if (fFluxInsidePlexiglassVar < 0.)
524 return kFALSE;
525
526 SetFluxInsidePlexiglassAvailable(kTRUE);
527
528 *fLog << inf << GetDescriptor()
529 << ": BlindPix #" << fPixId << ": Flux [ph/mm^2] inside Plexiglass: "
530 << Form("%5.3f +- %5.3f",fFluxInsidePlexiglass, GetFluxInsidePlexiglassErr()) << endl;
531
532 return kTRUE;
533}
534
535void MCalibrationBlindPix::Print(Option_t *opt) const
536{
537
538 *fLog << all << GetDescriptor()
539 << Form("BlindPixel: %3i" ,GetPixId())
540 << Form(" Lambda: %4.2f +- %4.2f",GetLambda(),GetLambdaErr())
541 << Form(" Mu0: %4.2f +- %4.2f" ,GetMu0(), GetMu0Err())
542 << Form(" Mu1: %4.2f +- %4.2f" ,GetMu1(), GetMu1Err())
543 << Form(" Sigma0: %4.2f +- %4.2f",GetSigma0(),GetSigma0Err())
544 << Form(" Sigma1: %4.2f +- %4.2f",GetSigma1(),GetSigma1Err())
545 << endl;
546 *fLog << all
547 << " Pedestal Fit OK? :" << IsPedestalFitOK()
548 << Form(" Lambda (Check): %4.2f +- %4.2f",GetLambdaCheck(),GetLambdaCheckErr())
549 << endl;
550 *fLog << all
551 << " Flux available? :" << IsFluxInsidePlexiglassAvailable()
552 << Form(" Flux: %4.2f +- %4.2f",GetFluxInsidePlexiglass(),GetFluxInsidePlexiglassErr())
553 << endl;
554}
Note: See TracBrowser for help on using the repository browser.