source: trunk/MagicSoft/Mars/mcalib/MCalibrationChargeCam.cc@ 3665

Last change on this file since 3665 was 3664, checked in by gaug, 22 years ago
*** empty log message ***
File size: 36.4 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 11/2003 <mailto:markus@ifae.es>
19!
20! Copyright: MAGIC Software Development, 2000-2004
21!
22!
23\* ======================================================================== */
24
25/////////////////////////////////////////////////////////////////////////////
26//
27// MCalibrationChargeCam
28//
29// Hold the whole Calibration results of the camera:
30//
31// 1) MCalibrationChargeCam initializes a TClonesArray whose elements are
32// pointers to MCalibrationChargePix Containers
33// 2) It initializes a pointer to an MCalibrationBlindPix container
34// 3) It initializes a pointer to an MCalibrationPINDiode container
35//
36//
37// The calculated values (types of GetPixelContent) are:
38//
39// --------------------------------------------------------------------------
40//
41// The types are as follows:
42//
43// Fitted values:
44// ==============
45//
46// 0: Fitted Charge
47// 1: Error of fitted Charge
48// 2: Sigma of fitted Charge
49// 3: Error of Sigma of fitted Charge
50//
51// Useful variables derived from the fit results:
52// =============================================
53//
54// 4: Returned probability of Gauss fit to Charge distribution
55// 5: Reduced Sigma of fitted Charge --> sqrt(sigma_Q^2 - PedRMS^2)
56// 6: Error Reduced Sigma of fitted Charge
57// 7: Reduced Sigma per Charge
58// 8: Error of Reduced Sigma per Charge
59//
60// Results of the different calibration methods:
61// =============================================
62//
63// 9: Number of Photo-electrons obtained with the F-Factor method
64// 10: Error on Number of Photo-electrons obtained with the F-Factor method
65// 11: Mean conversion factor obtained with the F-Factor method
66// 12: Error on the mean conversion factor obtained with the F-Factor method
67// 13: Overall F-Factor of the readout obtained with the F-Factor method
68// 14: Error on Overall F-Factor of the readout obtained with the F-Factor method
69// 15: Pixels with valid calibration by the F-Factor-Method
70// 16: Mean conversion factor obtained with the Blind Pixel method
71// 17: Error on the mean conversion factor obtained with the Blind Pixel method
72// 18: Overall F-Factor of the readout obtained with the Blind Pixel method
73// 19: Error on Overall F-Factor of the readout obtained with the Blind Pixel method
74// 20: Pixels with valid calibration by the Blind Pixel-Method
75// 21: Mean conversion factor obtained with the PIN Diode method
76// 22: Error on the mean conversion factor obtained with the PIN Diode method
77// 23: Overall F-Factor of the readout obtained with the PIN Diode method
78// 24: Error on Overall F-Factor of the readout obtained with the PIN Diode method
79// 25: Pixels with valid calibration by the PIN Diode-Method
80//
81// Localized defects:
82// ==================
83//
84// 26: Excluded Pixels
85// 27: Number of probable pickup events in the Hi Gain
86// 28: Number of probable pickup events in the Lo Gain
87//
88// Other classifications of pixels:
89// ================================
90//
91// 29: Pixels with saturated Hi-Gain
92//
93// Used Pedestals:
94// ===============
95//
96// 30: Mean Pedestal over the entire range of signal extraction
97// 31: Error on the Mean Pedestal over the entire range of signal extraction
98// 32: Pedestal RMS over the entire range of signal extraction
99// 33: Error on the Pedestal RMS over the entire range of signal extraction
100//
101// Calculated absolute arrival times (very low precision!):
102// ========================================================
103//
104// 34: Absolute Arrival time of the signal
105// 35: RMS of the Absolute Arrival time of the signal
106//
107/////////////////////////////////////////////////////////////////////////////
108#include "MCalibrationChargeCam.h"
109#include "MCalibrationCam.h"
110
111#include <TH2.h>
112#include <TCanvas.h>
113#include <TClonesArray.h>
114
115#include "MLog.h"
116#include "MLogManip.h"
117
118#include "MGeomCam.h"
119#include "MGeomPix.h"
120
121#include "MBadPixelsCam.h"
122#include "MBadPixelsPix.h"
123
124#include "MCalibrationChargePix.h"
125#include "MCalibrationChargeBlindPix.h"
126#include "MCalibrationChargePINDiode.h"
127
128ClassImp(MCalibrationChargeCam);
129
130using namespace std;
131
132const Float_t MCalibrationChargeCam::fgAverageQE = 0.18;
133const Float_t MCalibrationChargeCam::fgAverageQEErr = 0.02;
134const Float_t MCalibrationChargeCam::fgConvFFactorRelErrLimit = 0.35;
135const Float_t MCalibrationChargeCam::fgPheFFactorRelErrLimit = 5.;
136// --------------------------------------------------------------------------
137//
138// Default constructor.
139//
140// Creates a TClonesArray of MCalibrationPix containers, initialized to 1 entry
141// Later, a call to MCalibrationChargeCam::InitSize(Int_t size) has to be performed
142//
143// Creates an MCalibrationBlindPix container
144//
145MCalibrationChargeCam::MCalibrationChargeCam(const char *name, const char *title)
146 : fOffsets(NULL),
147 fSlopes(NULL),
148 fOffvsSlope(NULL)
149{
150 fName = name ? name : "MCalibrationChargeCam";
151 fTitle = title ? title : "Storage container for the Calibration Information in the camera";
152
153 fPixels = new TClonesArray("MCalibrationChargePix",1);
154 fAverageAreas = new TClonesArray("MCalibrationChargePix",1);
155 fAverageSectors = new TClonesArray("MCalibrationChargePix",1);
156
157 fFluxPhesInnerPixel = 0.;
158 fFluxPhesInnerPixelVar = 0.;
159 fFluxPhesOuterPixel = 0.;
160 fFluxPhesOuterPixelVar = 0.;
161
162 CLRBIT(fFlags,kBlindPixelMethodValid);
163 CLRBIT(fFlags,kFFactorMethodValid);
164 CLRBIT(fFlags,kPINDiodeMethodValid);
165
166 SetAverageQE();
167 SetConvFFactorRelErrLimit();
168 SetPheFFactorRelErrLimit();
169}
170
171// --------------------------------------------------------------------------
172//
173// Delete the TClonesArray of MCalibrationPix containers
174// Delete the MCalibrationPINDiode and the MCalibrationBlindPix
175//
176// Delete the histograms if they exist
177//
178MCalibrationChargeCam::~MCalibrationChargeCam()
179{
180
181
182 if (fOffsets)
183 delete fOffsets;
184 if (fSlopes)
185 delete fSlopes;
186 if (fOffvsSlope)
187 delete fOffvsSlope;
188
189}
190
191
192// --------------------------------------
193//
194void MCalibrationChargeCam::Clear(Option_t *o)
195{
196
197 fFluxPhesInnerPixel = 0.;
198 fFluxPhesInnerPixelVar = 0.;
199 fFluxPhesOuterPixel = 0.;
200 fFluxPhesOuterPixelVar = 0.;
201
202 CLRBIT(fFlags,kBlindPixelMethodValid);
203 CLRBIT(fFlags,kFFactorMethodValid);
204 CLRBIT(fFlags,kPINDiodeMethodValid);
205
206 MCalibrationCam::Clear();
207
208 return;
209}
210
211void MCalibrationChargeCam::SetFFactorMethodValid(const Bool_t b)
212{
213 b ? SETBIT(fFlags, kFFactorMethodValid) : CLRBIT(fFlags, kFFactorMethodValid);
214}
215
216void MCalibrationChargeCam::SetBlindPixelMethodValid(const Bool_t b)
217{
218 b ? SETBIT(fFlags, kBlindPixelMethodValid) : CLRBIT(fFlags, kBlindPixelMethodValid);
219}
220
221void MCalibrationChargeCam::SetPINDiodeMethodValid(const Bool_t b)
222{
223 b ? SETBIT(fFlags, kPINDiodeMethodValid) : CLRBIT(fFlags, kPINDiodeMethodValid);
224}
225
226Float_t MCalibrationChargeCam::GetFluxPhesInnerPixelErr() const
227{
228 if (fFluxPhesInnerPixelVar <= 0.)
229 return -1.;
230 return TMath::Sqrt(fFluxPhesInnerPixelVar);
231}
232
233Float_t MCalibrationChargeCam::GetFluxPhesOuterPixelErr() const
234{
235 if (fFluxPhesOuterPixelVar <= 0.)
236 return -1.;
237 return TMath::Sqrt(fFluxPhesOuterPixelVar);
238}
239
240Float_t MCalibrationChargeCam::GetFluxPhotonsInnerPixelErr() const
241{
242 if (fFluxPhotonsInnerPixelVar <= 0.)
243 return -1.;
244 return TMath::Sqrt(fFluxPhotonsInnerPixelVar);
245}
246
247Float_t MCalibrationChargeCam::GetFluxPhotonsOuterPixelErr() const
248{
249 if (fFluxPhotonsOuterPixelVar <= 0.)
250 return -1.;
251 return TMath::Sqrt(fFluxPhotonsOuterPixelVar);
252}
253
254Bool_t MCalibrationChargeCam::IsBlindPixelMethodValid() const
255{
256 return TESTBIT(fFlags,kBlindPixelMethodValid);
257}
258
259Bool_t MCalibrationChargeCam::IsPINDiodeMethodValid() const
260{
261 return TESTBIT(fFlags,kPINDiodeMethodValid);
262}
263
264
265// --------------------------------------------------------------------------
266//
267// Print first the well fitted pixels
268// and then the ones which are not FitValid
269//
270void MCalibrationChargeCam::Print(Option_t *o) const
271{
272
273 *fLog << all << GetDescriptor() << ":" << endl;
274 int id = 0;
275
276 *fLog << all << "Calibrated pixels:" << endl;
277 *fLog << all << endl;
278
279 TIter Next(fPixels);
280 MCalibrationChargePix *pix;
281 while ((pix=(MCalibrationChargePix*)Next()))
282 {
283
284 if (!pix->IsExcluded())
285 {
286
287 *fLog << all << "Pix " << pix->GetPixId()
288 << ": Ped. Rms: " << pix->GetPedRms() << " +- " << pix->GetPedRmsErr()
289 << " Mean signal: " << pix->GetMean() << " +- " << pix->GetSigma()
290 << " Reduced Sigma: " << pix->GetRSigma()
291 << " Nr Phe's: " << pix->GetPheFFactorMethod()
292 << " Saturated? :" << pix->IsHiGainSaturation()
293 << endl;
294 id++;
295 }
296 }
297
298 *fLog << all << id << " pixels" << endl;
299 id = 0;
300
301
302 *fLog << all << endl;
303 *fLog << all << "Excluded pixels:" << endl;
304 *fLog << all << endl;
305
306 id = 0;
307
308 TIter Next4(fPixels);
309 while ((pix=(MCalibrationChargePix*)Next4()))
310 {
311 if (pix->IsExcluded())
312 {
313 *fLog << all << pix->GetPixId() << endl;
314 id++;
315 }
316 }
317 *fLog << all << id << " Excluded pixels " << endl;
318 *fLog << endl;
319
320 TIter Next5(fAverageAreas);
321 while ((pix=(MCalibrationChargePix*)Next5()))
322 {
323 *fLog << all << "Average Area:"
324 << " Ped. Rms: " << pix->GetPedRms() << " +- " << pix->GetPedRmsErr()
325 << " Mean signal: " << pix->GetMean() << " +- " << pix->GetMeanErr()
326 << " Sigma signal: " << pix->GetSigma() << " +- "<< pix->GetSigmaErr()
327 << " Reduced Sigma: " << pix->GetRSigma()
328 << " Nr Phe's: " << pix->GetPheFFactorMethod()
329 << endl;
330 }
331
332 TIter Next6(fAverageSectors);
333 while ((pix=(MCalibrationChargePix*)Next5()))
334 {
335 *fLog << all << "Average Sector:"
336 << " Ped. Rms: " << pix->GetPedRms() << " +- " << pix->GetPedRmsErr()
337 << " Mean signal: " << pix->GetMean() << " +- " << pix->GetMeanErr()
338 << " Sigma signal: " << pix->GetSigma() << " +- "<< pix->GetSigmaErr()
339 << " Reduced Sigma: " << pix->GetRSigma()
340 << " Nr Phe's: " << pix->GetPheFFactorMethod()
341 << endl;
342 }
343
344}
345
346
347// --------------------------------------------------------------------------
348//
349// The types are as follows:
350//
351// Fitted values:
352// ==============
353//
354// 0: Fitted Charge
355// 1: Error of fitted Charge
356// 2: Sigma of fitted Charge
357// 3: Error of Sigma of fitted Charge
358//
359// Useful variables derived from the fit results:
360// =============================================
361//
362// 4: Returned probability of Gauss fit to Charge distribution
363// 5: Reduced Sigma of fitted Charge --> sqrt(sigma_Q^2 - PedRMS^2)
364// 6: Error Reduced Sigma of fitted Charge
365// 7: Reduced Sigma per Charge
366// 8: Error of Reduced Sigma per Charge
367//
368// Useful variables derived from the fit results:
369// =============================================
370//
371// 4: Returned probability of Gauss fit to Charge distribution
372// 5: Reduced Sigma of fitted Charge --> sqrt(sigma_Q^2 - PedRMS^2)
373// 6: Error Reduced Sigma of fitted Charge
374// 7: Reduced Sigma per Charge
375// 8: Error of Reduced Sigma per Charge
376//
377// Results of the different calibration methods:
378// =============================================
379//
380// 9: Number of Photo-electrons obtained with the F-Factor method
381// 10: Error on Number of Photo-electrons obtained with the F-Factor method
382// 11: Mean conversion factor obtained with the F-Factor method
383// 12: Error on the mean conversion factor obtained with the F-Factor method
384// 13: Overall F-Factor of the readout obtained with the F-Factor method
385// 14: Error on Overall F-Factor of the readout obtained with the F-Factor method
386// 15: Pixels with valid calibration by the F-Factor-Method
387// 16: Mean conversion factor obtained with the Blind Pixel method
388// 17: Error on the mean conversion factor obtained with the Blind Pixel method
389// 18: Overall F-Factor of the readout obtained with the Blind Pixel method
390// 19: Error on Overall F-Factor of the readout obtained with the Blind Pixel method
391// 20: Pixels with valid calibration by the Blind Pixel-Method
392// 21: Mean conversion factor obtained with the PIN Diode method
393// 22: Error on the mean conversion factor obtained with the PIN Diode method
394// 23: Overall F-Factor of the readout obtained with the PIN Diode method
395// 24: Error on Overall F-Factor of the readout obtained with the PIN Diode method
396// 25: Pixels with valid calibration by the PIN Diode-Method
397//
398// Localized defects:
399// ==================
400//
401// 26: Excluded Pixels
402// 27: Number of probable pickup events in the Hi Gain
403// 28: Number of probable pickup events in the Lo Gain
404//
405// Other classifications of pixels:
406// ================================
407//
408// 29: Pixels with saturated Hi-Gain
409//
410// Used Pedestals:
411// ===============
412//
413// 30: Mean Pedestal over the entire range of signal extraction
414// 31: Error on the Mean Pedestal over the entire range of signal extraction
415// 32: Pedestal RMS over the entire range of signal extraction
416// 33: Error on the Pedestal RMS over the entire range of signal extraction
417//
418// Calculated absolute arrival times (very low precision!):
419// ========================================================
420//
421// 34: Absolute Arrival time of the signal
422// 35: RMS of the Absolute Arrival time of the signal
423//
424Bool_t MCalibrationChargeCam::GetPixelContent(Double_t &val, Int_t idx, const MGeomCam &cam, Int_t type) const
425{
426
427 if (idx > GetSize())
428 return kFALSE;
429
430 Float_t area = cam[idx].GetA();
431
432 if (area == 0)
433 return kFALSE;
434
435 MCalibrationChargePix &pix = (MCalibrationChargePix&)(*this)[idx];
436
437 switch (type)
438 {
439 case 0:
440 if (pix.IsExcluded())
441 return kFALSE;
442 val = pix.GetMean();
443 break;
444 case 1:
445 if (pix.IsExcluded())
446 return kFALSE;
447 val = pix.GetMeanErr();
448 break;
449 case 2:
450 if (pix.IsExcluded())
451 return kFALSE;
452 val = pix.GetSigma();
453 break;
454 case 3:
455 if (pix.IsExcluded())
456 return kFALSE;
457 val = pix.GetSigmaErr();
458 break;
459 case 4:
460 if (pix.IsExcluded())
461 return kFALSE;
462 val = pix.GetProb();
463 break;
464 case 5:
465 if (pix.IsExcluded())
466 return kFALSE;
467 if (pix.GetRSigma() == -1.)
468 return kFALSE;
469 val = pix.GetRSigma();
470 break;
471 case 6:
472 if (pix.IsExcluded())
473 return kFALSE;
474 if (pix.GetRSigma() == -1.)
475 return kFALSE;
476 val = pix.GetRSigmaErr();
477 break;
478 case 7:
479 if (pix.IsExcluded())
480 return kFALSE;
481 if (pix.GetRSigma() == -1.)
482 return kFALSE;
483 if (pix.GetMean() == 0.)
484 return kFALSE;
485 val = pix.GetRSigma() / pix.GetMean();
486 break;
487 case 8:
488 if (pix.IsExcluded())
489 return kFALSE;
490 if (pix.GetRSigma() <= 0.)
491 return kFALSE;
492 if (pix.GetMean() <= 0.)
493 return kFALSE;
494 if (pix.GetRSigmaErr() <= 0.)
495 return kFALSE;
496 if (pix.GetMeanErr() <= 0.)
497 return kFALSE;
498 // relative error Rsigma square
499 val = pix.GetRSigmaErr()* pix.GetRSigmaErr()
500 / (pix.GetRSigma() * pix.GetRSigma() );
501 // relative error square
502 val += pix.GetMeanErr() * pix.GetMeanErr()
503 / (pix.GetMean() * pix.GetMean() );
504 // calculate relative error out of squares
505 val = TMath::Sqrt(val) ;
506 // multiply with value to get absolute error
507 val *= pix.GetRSigma() / pix.GetMean();
508 break;
509 case 9:
510 if (pix.IsExcluded() || !pix.IsFFactorMethodValid())
511 return kFALSE;
512 val = pix.GetPheFFactorMethod();
513 break;
514 case 10:
515 if (pix.IsExcluded() || !pix.IsFFactorMethodValid())
516 return kFALSE;
517 val = pix.GetPheFFactorMethodErr();
518 break;
519 case 11:
520 if (pix.IsExcluded() || !pix.IsFFactorMethodValid())
521 return kFALSE;
522 val = pix.GetMeanConversionFFactorMethod();
523 break;
524 case 12:
525 if (pix.IsExcluded() || !pix.IsFFactorMethodValid())
526 return kFALSE;
527 val = pix.GetConversionFFactorMethodErr();
528 break;
529 case 13:
530 if (pix.IsExcluded() || !pix.IsFFactorMethodValid())
531 return kFALSE;
532 val = pix.GetTotalFFactorFFactorMethod();
533 break;
534 case 14:
535 if (pix.IsExcluded() || !pix.IsFFactorMethodValid())
536 return kFALSE;
537 val = pix.GetTotalFFactorFFactorMethodErr();
538 break;
539 case 15:
540 if (pix.IsExcluded())
541 return kFALSE;
542 if (pix.IsFFactorMethodValid())
543 val = 1;
544 else
545 return kFALSE;
546 break;
547 case 16:
548 if (pix.IsExcluded() || !pix.IsBlindPixelMethodValid())
549 return kFALSE;
550 val = pix.GetMeanConversionBlindPixelMethod();
551 break;
552 case 17:
553 if (pix.IsExcluded() || !pix.IsBlindPixelMethodValid())
554 return kFALSE;
555 val = pix.GetConversionBlindPixelMethodErr();
556 break;
557 case 18:
558 if (pix.IsExcluded() || !pix.IsBlindPixelMethodValid())
559 return kFALSE;
560 val = pix.GetTotalFFactorBlindPixelMethod();
561 break;
562 case 19:
563 if (pix.IsExcluded() || !pix.IsBlindPixelMethodValid())
564 return kFALSE;
565 val = pix.GetTotalFFactorBlindPixelMethodErr();
566 break;
567 case 20:
568 if (pix.IsExcluded())
569 return kFALSE;
570 if (pix.IsBlindPixelMethodValid())
571 val = 1;
572 else
573 return kFALSE;
574 break;
575 case 21:
576 if (pix.IsExcluded() || !pix.IsPINDiodeMethodValid())
577 return kFALSE;
578 val = pix.GetMeanConversionPINDiodeMethod();
579 break;
580 case 22:
581 if (pix.IsExcluded() || !pix.IsPINDiodeMethodValid())
582 return kFALSE;
583 val = pix.GetConversionPINDiodeMethodErr();
584 break;
585 case 23:
586 if (pix.IsExcluded() || !pix.IsPINDiodeMethodValid())
587 return kFALSE;
588 val = pix.GetTotalFFactorPINDiodeMethod();
589 break;
590 case 24:
591 if (pix.IsExcluded() || !pix.IsPINDiodeMethodValid())
592 return kFALSE;
593 val = pix.GetTotalFFactorPINDiodeMethodErr();
594 break;
595 case 25:
596 if (pix.IsExcluded())
597 return kFALSE;
598 if (pix.IsPINDiodeMethodValid())
599 val = 1;
600 else
601 return kFALSE;
602 break;
603 case 26:
604 if (pix.IsExcluded())
605 val = 1.;
606 else
607 return kFALSE;
608 break;
609 case 27:
610 if (pix.IsExcluded())
611 return kFALSE;
612 val = pix.GetHiGainNumPickup();
613 break;
614 case 28:
615 if (pix.IsExcluded())
616 return kFALSE;
617 val = pix.GetLoGainNumPickup();
618 break;
619 case 29:
620 if (pix.IsExcluded())
621 return kFALSE;
622 val = pix.IsHiGainSaturation();
623 break;
624 case 30:
625 if (pix.IsExcluded())
626 return kFALSE;
627 val = pix.GetPed();
628 break;
629 case 31:
630 if (pix.IsExcluded())
631 return kFALSE;
632 val = pix.GetPedErr();
633 break;
634 case 32:
635 if (pix.IsExcluded())
636 return kFALSE;
637 val = pix.GetPedRms();
638 break;
639 case 33:
640 if (pix.IsExcluded())
641 return kFALSE;
642 val = pix.GetPedErr()/2.;
643 break;
644 case 34:
645 if (pix.IsExcluded())
646 return kFALSE;
647 val = pix.GetAbsTimeMean();
648 break;
649 case 35:
650 if (pix.IsExcluded())
651 return kFALSE;
652 val = pix.GetAbsTimeRms();
653 break;
654 default:
655 return kFALSE;
656 }
657
658 return val!=-1.;
659
660}
661
662// --------------------------------------------------------------------------
663//
664// What MHCamera needs in order to draw an individual pixel in the camera
665//
666void MCalibrationChargeCam::DrawPixelContent(Int_t idx) const
667{
668 MCalibrationChargePix &pix = (MCalibrationChargePix&)(*this)[idx];
669 pix.DrawClone();
670}
671
672//
673// Calculate the weighted mean of the phe's of all inner and outer pixels, respectively.
674// Bad pixels are excluded from the calculation. Two loops are performed to exclude pixels
675// which are fPheFFactorRelLimit sigmas from the mean.
676//
677Bool_t MCalibrationChargeCam::CalcFluxPhotonsFFactorMethod(const MGeomCam &geom, MBadPixelsCam &bad)
678{
679
680 const Float_t avQERelVar = fAverageQEVar / (fAverageQE * fAverageQE );
681
682 Float_t sumweightsinner = 0.;
683 Float_t sumphesinner = 0.;
684 Float_t sumweightsouter = 0.;
685 Float_t sumphesouter = 0.;
686 Int_t validinner = 0;
687 Int_t validouter = 0;
688
689 TIter Next(fPixels);
690 MCalibrationChargePix *pix;
691 while ((pix=(MCalibrationChargePix*)Next()))
692 {
693
694 if (!pix->IsFFactorMethodValid())
695 continue;
696
697 const Int_t idx = pix->GetPixId();
698
699 if(!bad[idx].IsCalibrationResultOK())
700 continue;
701
702 const Float_t nphe = pix->GetPheFFactorMethod();
703 const Float_t npheerr = pix->GetPheFFactorMethodErr();
704 const Float_t ratio = geom.GetPixRatio(idx);
705
706 if (npheerr > 0.)
707 {
708 const Float_t weight = nphe*nphe/npheerr/npheerr;
709 //
710 // first the inner pixels:
711 //
712 if (ratio == 1.)
713 {
714 sumweightsinner += weight;
715 sumphesinner += weight*nphe;
716 validinner++;
717 }
718 else
719 {
720 //
721 // now the outers
722 //
723 sumweightsouter += weight;
724 sumphesouter += weight*nphe;
725 validouter++;
726 }
727 } /* if npheerr != 0 */
728 } /* while ((pix=(MCalibrationChargePix*)Next())) */
729
730 if (sumweightsinner <= 0. || sumphesinner <= 0.)
731 {
732 *fLog << warn << " Mean number of phe's from inner pixels cannot be calculated: "
733 << " Sum of weights: " << sumweightsinner
734 << " Sum of weighted phes: " << sumphesinner << endl;
735 return kFALSE;
736 }
737 else
738 {
739 fFluxPhesInnerPixel = sumphesinner/sumweightsinner;
740 fFluxPhesInnerPixelVar = (1./sumweightsinner)*fFluxPhesInnerPixel*fFluxPhesInnerPixel;
741 }
742
743 if (sumweightsouter <= 0. || sumphesouter <= 0.)
744 {
745 *fLog << warn << " Mean number of phe's from outer pixels cannot be calculated: "
746 << " Sum of weights or sum of weighted phes is 0. " << endl;
747 }
748 else
749 {
750 fFluxPhesOuterPixel = sumphesouter/sumweightsouter;
751 fFluxPhesOuterPixelVar = (1./sumweightsouter)*fFluxPhesOuterPixel*fFluxPhesOuterPixel;
752 }
753
754 Float_t meanFluxPhotonsRelVar = fFluxPhesInnerPixelVar
755 / (fFluxPhesInnerPixel * fFluxPhesInnerPixel);
756
757 fFluxPhotonsInnerPixel = fFluxPhesInnerPixel/fAverageQE;
758 fFluxPhotonsInnerPixelVar = (meanFluxPhotonsRelVar + avQERelVar)
759 * fFluxPhotonsInnerPixel * fFluxPhotonsInnerPixel;
760
761 fFluxPhotonsOuterPixel = 4. *fFluxPhotonsInnerPixel;
762 fFluxPhotonsOuterPixelVar = 16.*fFluxPhotonsInnerPixelVar;
763
764 *fLog << inf << " Mean number of photo-electrons from inner pixels (F-Factor Method): "
765 << fFluxPhesInnerPixel << " +- " << GetFluxPhesInnerPixelErr() << endl;
766
767 *fLog << inf << " Mean number of photons from inner pixels (F-Factor Method): "
768 << fFluxPhotonsInnerPixel << " +- " << GetFluxPhotonsInnerPixelErr() << endl;
769
770 //
771 // Here starts the second loop discarting pixels out of the range:
772 //
773 const Float_t innervar = (Float_t)validinner*fPheFFactorRelVarLimit*fFluxPhesInnerPixelVar;
774 const Float_t outervar = (Float_t)validouter*fPheFFactorRelVarLimit*fFluxPhesOuterPixelVar;
775
776 Float_t innererr;
777 Float_t outererr;
778
779 if (innervar > 0.)
780 innererr = TMath::Sqrt(innervar);
781 else
782 {
783 *fLog << warn << GetDescriptor() << ": Variance of mean number of photo-electrons for inner pixels "
784 << " is smaller than 0. " << endl;
785 SetFFactorMethodValid(kFALSE);
786 return kFALSE;
787 }
788
789 if (outervar > 0.)
790 outererr = TMath::Sqrt(outervar);
791 else
792 {
793 *fLog << warn << GetDescriptor() << ": Variance of mean number of photo-electrons for outer pixels "
794 << " is smaller than 0. " << endl;
795 SetFFactorMethodValid(kFALSE);
796 return kFALSE;
797 }
798
799 const Float_t lowerpheinnerlimit = fFluxPhesInnerPixel - innererr;
800 const Float_t upperpheinnerlimit = fFluxPhesInnerPixel + innererr;
801
802 const Float_t lowerpheouterlimit = fFluxPhesOuterPixel - outererr;
803 const Float_t upperpheouterlimit = fFluxPhesOuterPixel + outererr;
804
805 sumweightsinner = 0.;
806 sumphesinner = 0.;
807 sumweightsouter = 0.;
808 sumphesouter = 0.;
809
810 TIter Next2(fPixels);
811 MCalibrationChargePix *pix2;
812 while ((pix2=(MCalibrationChargePix*)Next2()))
813 {
814
815 if (!pix2->IsFFactorMethodValid())
816 continue;
817
818 const Int_t idx = pix2->GetPixId();
819
820 if(!bad[idx].IsCalibrationResultOK())
821 continue;
822
823 const Float_t nphe = pix2->GetPheFFactorMethod();
824 const Float_t npheerr = pix2->GetPheFFactorMethodErr();
825 const Float_t ratio = geom.GetPixRatio(idx);
826
827 if (npheerr > 0.)
828 {
829 const Float_t weight = nphe*nphe/npheerr/npheerr;
830 //
831 // first the inner pixels:
832 //
833 if (ratio == 1.)
834 {
835
836 if (nphe < lowerpheinnerlimit || nphe > upperpheinnerlimit)
837 {
838 bad[idx].SetUncalibrated(MBadPixelsPix::kDeviatingNumPhes);
839 bad[idx].SetUnsuitable( MBadPixelsPix::kUnreliableRun);
840 continue;
841 }
842
843 sumweightsinner += weight;
844 sumphesinner += weight*nphe;
845 }
846 else
847 {
848 //
849 // now the outers
850 //
851 if (nphe < lowerpheouterlimit || nphe > upperpheouterlimit)
852 {
853 bad[idx].SetUncalibrated(MBadPixelsPix::kDeviatingNumPhes);
854 bad[idx].SetUnsuitable( MBadPixelsPix::kUnreliableRun);
855 continue;
856 }
857
858 sumweightsouter += weight;
859 sumphesouter += weight*nphe;
860 }
861 } /* if npheerr != 0 */
862 } /* while ((pix2=(MCalibrationChargePix*)Next2())) */
863
864 if (sumweightsinner <= 0. || sumphesinner <= 0.)
865 {
866 *fLog << warn << " Mean number of phe's from inner pixels cannot be calculated: "
867 << " Sum of weights: " << sumweightsinner
868 << " Sum of weighted phes: " << sumphesinner << endl;
869 return kFALSE;
870 }
871 else
872 {
873 fFluxPhesInnerPixel = sumphesinner/sumweightsinner;
874 fFluxPhesInnerPixelVar = (1./sumweightsinner)*fFluxPhesInnerPixel*fFluxPhesInnerPixel;
875
876 }
877
878 if (sumweightsouter <= 0. || sumphesouter <= 0.)
879 {
880 *fLog << warn << " Mean number of phe's from outer pixels cannot be calculated: "
881 << " Sum of weights or sum of weighted phes is 0. " << endl;
882 }
883 else
884 {
885 fFluxPhesOuterPixel = sumphesouter/sumweightsouter;
886 fFluxPhesOuterPixelVar = (1./sumweightsouter)*fFluxPhesOuterPixel*fFluxPhesOuterPixel;
887 }
888
889 meanFluxPhotonsRelVar = fFluxPhesInnerPixelVar
890 / (fFluxPhesInnerPixel * fFluxPhesInnerPixel);
891
892 fFluxPhotonsInnerPixel = fFluxPhesInnerPixel/fAverageQE;
893 fFluxPhotonsInnerPixelVar = (meanFluxPhotonsRelVar + avQERelVar)
894 * fFluxPhotonsInnerPixel * fFluxPhotonsInnerPixel;
895
896 fFluxPhotonsOuterPixel = 4. *fFluxPhotonsInnerPixel;
897 fFluxPhotonsOuterPixelVar = 16.*fFluxPhotonsInnerPixelVar;
898
899 *fLog << inf << " Mean number of photo-electrons from inner pixels (F-Factor Method): "
900 << fFluxPhesInnerPixel << " +- " << GetFluxPhesInnerPixelErr() << endl;
901
902 *fLog << inf << " Mean number of photons from inner pixels (F-Factor Method): "
903 << fFluxPhotonsInnerPixel << " +- " << GetFluxPhotonsInnerPixelErr() << endl;
904
905 return kTRUE;
906}
907
908void MCalibrationChargeCam::ApplyFFactorCalibration(const MGeomCam &geom, const MBadPixelsCam &bad)
909{
910
911 const Float_t meanphotRelVar = fFluxPhotonsInnerPixelVar
912 /( fFluxPhotonsInnerPixel * fFluxPhotonsInnerPixel );
913
914 TIter Next(fPixels);
915 MCalibrationChargePix *pix;
916 while ((pix=(MCalibrationChargePix*)Next()))
917 {
918
919 const Int_t idx = pix->GetPixId();
920
921 if (!(bad[idx].IsCalibrationResultOK()))
922 {
923 pix->SetFFactorMethodValid(kFALSE);
924 continue;
925 }
926
927 if (!(pix->IsFFactorMethodValid()))
928 continue;
929
930 const Float_t ratio = geom.GetPixRatio(idx);
931 //
932 // Calculate the conversion factor between PHOTONS and FADC counts
933 //
934 // Nphot = Nphe / avQE
935 // conv = Nphot / FADC counts
936 //
937 Float_t conv;
938
939 if (ratio == 1.)
940 conv = fFluxPhotonsInnerPixel / pix->GetMean();
941 else
942 conv = fFluxPhotonsOuterPixel / pix->GetMean();
943
944 if (conv <= 0.)
945 {
946 pix->SetFFactorMethodValid(kFALSE);
947 continue;
948 }
949
950 const Float_t chargeRelVar = pix->GetMeanErr() * pix->GetMeanErr()
951 / ( pix->GetMean() * pix->GetMean());
952 const Float_t rsigmaRelVar = pix->GetRSigmaErr() * pix->GetRSigmaErr()
953 / (pix->GetRSigma() * pix->GetRSigma()) ;
954
955 const Float_t convrelvar = meanphotRelVar + chargeRelVar;
956
957 if (convrelvar > fConvFFactorRelVarLimit)
958 {
959 *fLog << warn << GetDescriptor() << ": Conversion Factor F-Factor Method Rel. Variance: "
960 << convrelvar << " above limit of: " << fConvFFactorRelVarLimit
961 << " in pixel: " << idx << endl;
962 pix->SetFFactorMethodValid(kFALSE);
963 continue;
964 }
965
966 //
967 // Calculate the Total F-Factor of the camera (in photons)
968 //
969 const Float_t totalFFactor = (pix->GetRSigma()/pix->GetMean())
970 *TMath::Sqrt(fFluxPhotonsInnerPixel);
971
972 //
973 // Calculate the error of the Total F-Factor of the camera ( in photons )
974 //
975 const Float_t totalFFactorVar = rsigmaRelVar + chargeRelVar + meanphotRelVar;
976
977 if (convrelvar > 0. && conv > 0.)
978 pix->SetConversionFFactorMethod(conv,
979 TMath::Sqrt(convrelvar)*conv,
980 totalFFactor*TMath::Sqrt(conv));
981
982 pix->SetTotalFFactorFFactorMethod( totalFFactor );
983
984 if (totalFFactorVar > 0.)
985 pix->SetTotalFFactorFFactorMethodErr(TMath::Sqrt(totalFFactorVar));
986
987 pix->SetFFactorMethodValid();
988 }
989}
990
991
992
993void MCalibrationChargeCam::ApplyBlindPixelCalibration(const MGeomCam &geom, const MBadPixelsCam &bad, const MCalibrationChargeBlindPix &blindpix)
994{
995
996 Float_t flux = blindpix.GetFluxInsidePlexiglass();
997 Float_t fluxerr = blindpix.GetFluxInsidePlexiglassErr();
998
999 TIter Next(fPixels);
1000 MCalibrationChargePix *pix;
1001 while ((pix=(MCalibrationChargePix*)Next()))
1002 {
1003
1004 const Int_t idx = pix->GetPixId();
1005
1006 if (!(bad[idx].IsCalibrationResultOK()))
1007 {
1008 pix->SetBlindPixelMethodValid(kFALSE);
1009 continue;
1010 }
1011
1012 const Float_t charge = pix->GetMean();
1013 const Float_t area = geom[idx].GetA();
1014 const Float_t chargeerr = pix->GetMeanErr();
1015
1016 const Float_t nphot = flux * area;
1017 const Float_t nphoterr = fluxerr * area;
1018 const Float_t conversion = nphot/charge;
1019 Float_t conversionerr;
1020
1021 conversionerr = nphoterr/charge
1022 * nphoterr/charge ;
1023 conversionerr += chargeerr/charge
1024 * chargeerr/charge
1025 * conversion*conversion;
1026 conversionerr = TMath::Sqrt(conversionerr);
1027
1028 const Float_t conversionsigma = 0.;
1029
1030 pix->SetConversionBlindPixelMethod(conversion, conversionerr, conversionsigma);
1031
1032 if (conversionerr/conversion < 0.1)
1033 pix->SetBlindPixelMethodValid();
1034 }
1035}
1036
1037void MCalibrationChargeCam::ApplyPINDiodeCalibration(const MGeomCam &geom, const MBadPixelsCam &bad, const MCalibrationChargePINDiode &pindiode)
1038{
1039
1040 Float_t flux = pindiode.GetFluxOutsidePlexiglass();
1041 Float_t fluxerr = pindiode.GetFluxOutsidePlexiglassErr();
1042
1043 TIter Next(fPixels);
1044 MCalibrationChargePix *pix;
1045 while ((pix=(MCalibrationChargePix*)Next()))
1046 {
1047
1048 const Int_t idx = pix->GetPixId();
1049
1050 if (!(bad[idx].IsCalibrationResultOK()))
1051 {
1052 pix->SetPINDiodeMethodValid(kFALSE);
1053 continue;
1054 }
1055
1056 const Float_t charge = pix->GetMean();
1057 const Float_t area = geom[idx].GetA();
1058 const Float_t chargeerr = pix->GetMeanErr();
1059
1060 const Float_t nphot = flux * area;
1061 const Float_t nphoterr = fluxerr * area;
1062 const Float_t conversion = nphot/charge;
1063
1064 Float_t conversionerr;
1065
1066 conversionerr = nphoterr/charge
1067 * nphoterr/charge ;
1068 conversionerr += chargeerr/charge
1069 * chargeerr/charge
1070 * conversion*conversion;
1071 if (conversionerr > 0.)
1072 conversionerr = TMath::Sqrt(conversionerr);
1073
1074 const Float_t conversionsigma = 0.;
1075
1076 pix->SetConversionPINDiodeMethod(conversion, conversionerr, conversionsigma);
1077
1078 if (conversionerr/conversion < 0.1)
1079 pix->SetPINDiodeMethodValid();
1080
1081 }
1082}
1083
1084
1085Bool_t MCalibrationChargeCam::GetConversionFactorBlindPixel(Int_t ipx, Float_t &mean, Float_t &err, Float_t &sigma)
1086{
1087
1088 MCalibrationChargePix &pix = (MCalibrationChargePix&)(*this)[ipx];
1089
1090 mean = pix.GetMeanConversionBlindPixelMethod();
1091 err = pix.GetConversionBlindPixelMethodErr();
1092 sigma = pix.GetSigmaConversionBlindPixelMethod();
1093
1094 return kTRUE;
1095}
1096
1097
1098Bool_t MCalibrationChargeCam::GetConversionFactorFFactor(Int_t ipx, Float_t &mean, Float_t &err, Float_t &sigma)
1099{
1100
1101 MCalibrationChargePix &pix = (MCalibrationChargePix&)(*this)[ipx];
1102
1103 Float_t conv = pix.GetMeanConversionFFactorMethod();
1104
1105 if (conv < 0.)
1106 return kFALSE;
1107
1108 mean = conv;
1109 err = pix.GetConversionFFactorMethodErr();
1110 sigma = pix.GetSigmaConversionFFactorMethod();
1111
1112 return kTRUE;
1113}
1114
1115
1116//-----------------------------------------------------------------------------------
1117//
1118// Calculates the conversion factor between the integral of FADCs slices
1119// (as defined in the signal extractor MExtractSignal.cc)
1120// and the number of photons reaching the plexiglass for one Inner Pixel
1121//
1122Bool_t MCalibrationChargeCam::GetConversionFactorPINDiode(Int_t ipx, Float_t &mean, Float_t &err, Float_t &sigma)
1123{
1124
1125 MCalibrationChargePix &pix = (MCalibrationChargePix&)(*this)[ipx];
1126
1127 mean = pix.GetMeanConversionPINDiodeMethod();
1128 err = pix.GetConversionPINDiodeMethodErr();
1129 sigma = pix.GetSigmaConversionPINDiodeMethod();
1130
1131 return kFALSE;
1132
1133}
1134
1135//-----------------------------------------------------------------------------------
1136//
1137// Calculates the best combination of the three used methods possible
1138// between the integral of FADCs slices
1139// (as defined in the signal extractor MExtractSignal.cc)
1140// and the number of photons reaching one Inner Pixel.
1141// The procedure is not yet defined.
1142//
1143Bool_t MCalibrationChargeCam::GetConversionFactorCombined(Int_t ipx, Float_t &mean, Float_t &err, Float_t &sigma)
1144{
1145 return kFALSE;
1146
1147}
1148
1149/*
1150void MCalibrationChargeCam::DrawHiLoFits()
1151{
1152
1153 if (!fOffsets)
1154 fOffsets = new TH1D("pp","Offsets of the HiGain LoGain Fit",100,-600.,400.);
1155 if (!fSlopes)
1156 fSlopes = new TH1D("mm","Slopes of the HiGain LoGain Fit",100,-2.,2.);
1157 if (!fOffvsSlope)
1158 fOffvsSlope = new TH2D("aa","Slopes vs Offsets of the HiGain LoGain Fit",100,-600.,400.,100,-2.,2.);
1159
1160 TIter Next(fPixels);
1161 MCalibrationPix *pix;
1162 MHCalibrationPixel *hist;
1163 while ((pix=(MCalibrationPix*)Next()))
1164 {
1165 hist = pix->GetHist();
1166 hist->FitHiGainvsLoGain();
1167 fOffsets->Fill(hist->GetOffset(),1.);
1168 fSlopes->Fill(hist->GetSlope(),1.);
1169 fOffvsSlope->Fill(hist->GetOffset(),hist->GetSlope(),1.);
1170 }
1171
1172 TCanvas *c1 = new TCanvas();
1173
1174 c1->Divide(1,3);
1175 c1->cd(1);
1176 fOffsets->Draw();
1177 gPad->Modified();
1178 gPad->Update();
1179
1180 c1->cd(2);
1181 fSlopes->Draw();
1182 gPad->Modified();
1183 gPad->Update();
1184
1185 c1->cd(3);
1186 fOffvsSlope->Draw("col1");
1187 gPad->Modified();
1188 gPad->Update();
1189}
1190
1191*/
Note: See TracBrowser for help on using the repository browser.