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

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