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

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