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

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