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

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