source: trunk/MagicSoft/Mars/mcalib/MCalibrationChargeCalc.cc@ 4609

Last change on this file since 4609 was 4609, checked in by tbretz, 20 years ago
*** empty log message ***
File size: 64.0 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// MCalibrationChargeCalc
28//
29// Task to calculate the calibration conversion factors and quantum efficiencies
30// from the fit results to the summed FADC slice distributions delivered by
31// MCalibrationChargeCam, MCalibrationChargePix, MCalibrationChargeBlindPix and
32// MCalibrationChargePINDiode, calculated and filled by MHCalibrationChargeCam,
33// MHCalibrationChargePix, MHCalibrationChargeBlindPix and MHCalibrationChargePINDiode.
34//
35// PreProcess(): Initialize pointers to MCalibrationChargeCam, MCalibrationChargeBlindPix
36// MCalibrationChargePINDiode and MCalibrationQECam
37//
38// Initialize pulser light wavelength
39//
40// ReInit(): MCalibrationCam::InitSize(NumPixels) is called from MGeomApply (which allocates
41// memory in a TClonesArray of type MCalibrationChargePix)
42// Initializes pointer to MBadPixelsCam
43//
44// Process(): Nothing to be done, histograms getting filled by MHCalibrationChargeCam
45//
46// PostProcess(): - FinalizePedestals()
47// - FinalizeCharges()
48// - FinalizeFFactorMethod()
49// - FinalizeBadPixels()
50// - FinalizeBlindPixel()
51// - FinalizePINDiode()
52// - FinalizeFFactorQECam()
53// - FinalizeBlindPixelQECam()
54// - FinalizePINDiodeQECam()
55//
56// Input Containers:
57// MCalibrationChargeCam
58// MCalibrationChargeBlindPix
59// MCalibrationChargePINDiode
60// MCalibrationQECam
61// MPedestalCam
62// MBadPixelsCam
63// MGeomCam
64// MTime
65//
66// Output Containers:
67// MCalibrationChargeCam
68// MCalibrationChargeBlindPix
69// MCalibrationChargePINDiode
70// MCalibrationQECam
71// MBadPixelsCam
72//
73//
74// Preliminary description of the calibration in photons (email from 12/02/04)
75//
76// Why calibrating in photons:
77// ===========================
78//
79// At the Barcelona meeting in 2002, we decided to calibrate the camera in
80// photons. This for the following reasons:
81//
82// * The physical quantity arriving at the camera are photons. This is
83// the direct physical information from the air shower. The photons
84// have a flux and a spectrum.
85//
86// * The photon fluxes depend mostly on the shower energy (with
87// corrections deriving from the observation conditions), while the photon
88// spectra depend mostly on the observation conditions: zenith angle,
89// quality of the air, also the impact parameter of the shower.
90//
91// * The photomultiplier, in turn, has different response properties
92// (quantum efficiencies) for photons of different colour. (Moreover,
93// different pixels have slightly different quantum efficiencies).
94// The resulting number of photo-electrons is then amplified (linearly)
95// with respect to the photo-electron flux.
96//
97// * In the ideal case, one would like to disentagle the effects
98// of the observation conditions from the primary particle energy (which
99// one likes to measure). To do so, one needs:
100//
101// 1) A reliable calibration relating the FADC counts to the photo-electron
102// flux -> This is accomplished with the F-Factor method.
103//
104// 2) A reliable calibration of the wavelength-dependent quantum efficiency
105// -> This is accomplished with the combination of the three methods,
106// together with QE-measurements performed by David in order to do
107// the interpolation.
108//
109// 3) A reliable calibration of the observation conditions. This means:
110// - Tracing the atmospheric conditions -> LIDAR
111// - Tracing the observation zenith angle -> Drive System
112//
113// 4) Some knowlegde about the impact parameter:
114// - This is the only part which cannot be accomplished well with a
115// single telescope. We would thus need to convolute the spectrum
116// over the distribution of impact parameters.
117//
118//
119// How an ideal calibration would look like:
120// =========================================
121//
122// We know from the combined PIN-Diode and Blind-Pixel Method the response of
123// each pixel to well-measured light fluxes in three representative
124// wavelengths (green, blue, UV). We also know the response to these light
125// fluxes in photo-electrons. Thus, we can derive:
126//
127// - conversion factors to photo-electrons
128// - conversion factors to photons in three wavelengths.
129//
130// Together with David's measurements and some MC-simulation, we should be
131// able to derive tables for typical Cherenkov-photon spectra - convoluted
132// with the impact parameters and depending on the athmospheric conditions
133// and the zenith angle (the "outer parameters").
134//
135// From these tables we can create "calibration tables" containing some
136// effective quantum efficiency depending on these outer parameters and which
137// are different for each pixel.
138//
139// In an ideal MCalibrate, one would thus have to convert first the FADC
140// slices to Photo-electrons and then, depending on the outer parameters,
141// look up the effective quantum efficiency and get the mean number of
142// photons which is then used for the further analysis.
143//
144// How the (first) MAGIC calibration should look like:
145// ===================================================
146//
147// For the moment, we have only one reliable calibration method, although
148// with very large systematic errors. This is the F-Factor method. Knowing
149// that the light is uniform over the whole camera (which I would not at all
150// guarantee in the case of the CT1 pulser), one could in principle already
151// perform a relative calibration of the quantum efficiencies in the UV.
152// However, the spread in QE at UV is about 10-15% (according to the plot
153// that Abelardo sent around last time. The spread in photo-electrons is 15%
154// for the inner pixels, but much larger (40%) for the outer ones.
155//
156// I'm not sure if we can already say that we have measured the relative
157// difference in quantum efficiency for the inner pixels and produce a first
158// QE-table for each pixel. To so, I would rather check in other wavelengths
159// (which we can do in about one-two weeks when the optical transmission of
160// the calibration trigger is installed).
161//
162// Thus, for the moment being, I would join Thomas proposal to calibrate in
163// photo-electrons and apply one stupid average quantum efficiency for all
164// pixels. This keeping in mind that we will have much preciser information
165// in about one to two weeks.
166//
167//
168// What MCalibrate should calculate and what should be stored:
169// ===========================================================
170//
171// It is clear that in the end, MCerPhotEvt will store photons.
172// MCalibrationCam stores the conversionfactors to photo-electrons and also
173// some tables of how to apply the conversion to photons, given the outer
174// parameters. This is not yet implemented and not even discussed.
175//
176// To start, I would suggest that we define the "average quantum efficiency"
177// (maybe something like 25+-3%) and apply them equally to all
178// photo-electrons. Later, this average factor can be easily replaced by a
179// pixel-dependent factor and later by a (pixel-dependent) table.
180//
181//
182//
183//////////////////////////////////////////////////////////////////////////////
184#include "MCalibrationChargeCalc.h"
185
186#include <TSystem.h>
187#include <TH1.h>
188#include <TF1.h>
189
190#include "MLog.h"
191#include "MLogManip.h"
192
193#include "MParList.h"
194
195#include "MRawRunHeader.h"
196#include "MRawEvtPixelIter.h"
197
198#include "MGeomCam.h"
199#include "MGeomPix.h"
200#include "MHCamera.h"
201
202#include "MPedestalCam.h"
203#include "MPedestalPix.h"
204
205#include "MCalibrationChargeCam.h"
206#include "MCalibrationChargePix.h"
207#include "MCalibrationChargePINDiode.h"
208#include "MCalibrationChargeBlindPix.h"
209#include "MCalibrationChargeBlindCam.h"
210
211#include "MExtractedSignalCam.h"
212#include "MExtractedSignalPix.h"
213#include "MExtractedSignalBlindPixel.h"
214#include "MExtractedSignalPINDiode.h"
215
216#include "MBadPixelsCam.h"
217#include "MBadPixelsPix.h"
218
219#include "MCalibrationQECam.h"
220#include "MCalibrationQEPix.h"
221
222#include "MCalibrationCam.h"
223
224ClassImp(MCalibrationChargeCalc);
225
226using namespace std;
227
228const Float_t MCalibrationChargeCalc::fgChargeLimit = 2.5;
229const Float_t MCalibrationChargeCalc::fgChargeErrLimit = 0.;
230const Float_t MCalibrationChargeCalc::fgChargeRelErrLimit = 1.;
231const Float_t MCalibrationChargeCalc::fgLambdaErrLimit = 0.2;
232const Float_t MCalibrationChargeCalc::fgLambdaCheckLimit = 0.5;
233const Float_t MCalibrationChargeCalc::fgPheErrLimit = 4.5;
234const Float_t MCalibrationChargeCalc::fgFFactorErrLimit = 4.5;
235// --------------------------------------------------------------------------
236//
237// Default constructor.
238//
239// Sets the pointer to fQECam and fGeom to NULL
240//
241// Calls AddToBranchList for:
242// - MRawEvtData.fHiGainPixId
243// - MRawEvtData.fLoGainPixId
244// - MRawEvtData.fHiGainFadcSamples
245// - MRawEvtData.fLoGainFadcSamples
246//
247// Initializes:
248// - fChargeLimit to fgChargeLimit
249// - fChargeErrLimit to fgChargeErrLimit
250// - fChargeRelErrLimit to fgChargeRelErrLimit
251// - fFFactorErrLimit to fgFFactorErrLimit
252// - fLambdaCheckLimit to fgLambdaCheckLimit
253// - fLambdaErrLimit to fgLambdaErrLimit
254// - fPheErrLimit to fgPheErrLimit
255// - fPulserColor to MCalibrationCam::kCT1
256// - fOutputPath to "."
257// - fOutputFile to "ChargeCalibStat.txt"
258// - flag debug to kFALSE
259//
260// Calls:
261// - Clear()
262//
263MCalibrationChargeCalc::MCalibrationChargeCalc(const char *name, const char *title)
264 : fQECam(NULL), fGeom(NULL)
265{
266
267 fName = name ? name : "MCalibrationChargeCalc";
268 fTitle = title ? title : "Task to calculate the calibration constants and MCalibrationCam ";
269
270 AddToBranchList("MRawEvtData.fHiGainPixId");
271 AddToBranchList("MRawEvtData.fLoGainPixId");
272 AddToBranchList("MRawEvtData.fHiGainFadcSamples");
273 AddToBranchList("MRawEvtData.fLoGainFadcSamples");
274
275 SetChargeLimit ();
276 SetChargeErrLimit ();
277 SetChargeRelErrLimit ();
278 SetFFactorErrLimit ();
279 SetLambdaCheckLimit ();
280 SetLambdaErrLimit ();
281 SetPheErrLimit ();
282 SetOutputPath ();
283 SetOutputFile ();
284 SetDebug ( kFALSE );
285
286 Clear();
287
288}
289
290// --------------------------------------------------------------------------
291//
292// Sets:
293// - all variables to 0.,
294// - all flags to kFALSE
295// - all pointers to NULL
296// - the pulser colour to kNONE
297// - fBlindPixelFlags to 0
298// - fPINDiodeFlags to 0
299//
300void MCalibrationChargeCalc::Clear(const Option_t *o)
301{
302
303 fNumHiGainSamples = 0.;
304 fNumLoGainSamples = 0.;
305 fSqrtHiGainSamples = 0.;
306 fSqrtLoGainSamples = 0.;
307 fNumInnerFFactorMethodUsed = 0;
308
309 fBadPixels = NULL;
310 fCam = NULL;
311 fBlindPixel = NULL;
312 fBlindCam = NULL;
313 fPINDiode = NULL;
314 fPedestals = NULL;
315
316 SetPulserColor ( MCalibrationCam::kNONE );
317
318 fBlindPixelFlags.Set(0);
319 fPINDiodeFlags .Set(0);
320 fResultFlags .Set(0);
321}
322
323
324// -----------------------------------------------------------------------------------
325//
326// The following container are searched for and execution aborted if not in MParList:
327// - MPedestalCam
328//
329// The following containers are searched and created if they were not found:
330//
331// - MCalibrationQECam
332// - MBadPixelsCam
333//
334Int_t MCalibrationChargeCalc::PreProcess(MParList *pList)
335{
336
337 //
338 // Containers that have to be there.
339 //
340 fPedestals = (MPedestalCam*)pList->FindObject("MPedestalCam");
341 if (!fPedestals)
342 {
343 *fLog << err << "MPedestalCam not found... aborting" << endl;
344 return kFALSE;
345 }
346
347 //
348 // Containers that are created in case that they are not there.
349 //
350 fQECam = (MCalibrationQECam*)pList->FindCreateObj("MCalibrationQECam");
351 if (!fQECam)
352 {
353 *fLog << err << "Cannot find nor create MCalibrationQECam... aborting" << endl;
354 return kFALSE;
355 }
356
357 fBadPixels = (MBadPixelsCam*)pList->FindCreateObj("MBadPixelsCam");
358 if (!fBadPixels)
359 {
360 *fLog << err << "Could not find or create MBadPixelsCam ... aborting." << endl;
361 return kFALSE;
362 }
363
364
365 //
366 // Check the pulser colour --> FIXME: this solution is only valid until the arrival of the DM's
367 //
368 if (fPulserColor == MCalibrationCam::kNONE)
369 {
370 *fLog << endl;
371 *fLog << err << GetDescriptor()
372 << ": No Pulser colour has been chosen. Since the installation of the IFAE pulser box,"
373 << " you HAVE to provide the LEDs colour, otherwise there is no calibration. " << endl;
374 *fLog << "See e.g. the macro calibration.C " << endl;
375 return kFALSE;
376 }
377
378 return kTRUE;
379}
380
381
382// --------------------------------------------------------------------------
383//
384// Search for the following input containers and abort if not existing:
385// - MGeomCam
386// - MCalibrationChargeCam
387//
388// Search for the following input containers and give a warning if not existing:
389// - MCalibrationChargeBlindPix
390// - MCalibrationChargePINDiode
391//
392// It retrieves the following variables from MCalibrationChargeCam:
393//
394// - fNumHiGainSamples
395// - fNumLoGainSamples
396//
397// It defines the PixId of every pixel in:
398//
399// - MCalibrationChargeCam
400// - MCalibrationQECam
401//
402// It sets all pixels in excluded which have the flag fBadBixelsPix::IsBad() set in:
403//
404// - MCalibrationChargePix
405// - MCalibrationQEPix
406//
407// Sets the pulser colour and tests if it has not changed w.r.t. fPulserColor in:
408//
409// - MCalibrationChargeCam
410// - MCalibrationChargeBlindPix (if existing)
411// - MCalibrationChargePINDiode (if existing)
412//
413Bool_t MCalibrationChargeCalc::ReInit(MParList *pList )
414{
415
416 fGeom = (MGeomCam*)pList->FindObject("MGeomCam");
417 if (!fGeom)
418 {
419 *fLog << err << "No MGeomCam found... aborting." << endl;
420 return kFALSE;
421 }
422
423 fCam = (MCalibrationChargeCam*)pList->FindObject("MCalibrationChargeCam");
424 if (!fCam)
425 {
426 *fLog << err << "Cannot find MCalibrationChargeCam... aborting" << endl;
427 *fLog << err << "Maybe you forget to call an MFillH for the MHCalibrationChargeCam before..." << endl;
428 return kFALSE;
429 }
430
431 //
432 // Optional Containers
433 //
434 fBlindPixel = (MCalibrationChargeBlindPix*)pList->FindObject("MCalibrationChargeBlindPix");
435 if (!fBlindPixel)
436 {
437 fBlindCam = (MCalibrationChargeBlindCam*)pList->FindObject("MCalibrationChargeBlindCam");
438 if (!fBlindCam)
439 {
440 *fLog << endl;
441 *fLog << warn << GetDescriptor()
442 << ": MCalibrationChargeBlindPix nor MCalibrationChargeBlindCam "
443 << " found... no Blind Pixel method! " << endl;
444 }
445 }
446
447 fPINDiode = (MCalibrationChargePINDiode*)pList->FindObject("MCalibrationChargePINDiode");
448 if (!fPINDiode)
449 {
450 *fLog << endl;
451 *fLog << warn << GetDescriptor()
452 << ": MCalibrationChargePINDiode not found... no PIN Diode method! " << endl;
453 }
454
455
456 //
457 // Initialize the pulser colours
458 //
459 if (fCam->GetPulserColor() == MCalibrationCam::kNONE)
460 {
461 fCam->SetPulserColor( fPulserColor );
462
463 if (fBlindPixel)
464 fBlindPixel->SetColor( fPulserColor );
465
466 if (fBlindCam)
467 fBlindCam->SetColor( fPulserColor );
468
469 if (fPINDiode)
470 fPINDiode->SetColor( fPulserColor );
471 }
472
473 if (fPulserColor != fCam->GetPulserColor())
474 {
475 *fLog << err << GetDescriptor()
476 << ": Pulser colour has changed w.r.t. last file in MCalibrationChargeCam" << endl;
477 *fLog << err << "This feature is not yet implemented, sorry ... aborting " << endl;
478 return kFALSE;
479 }
480
481 if (fBlindPixel)
482 if (fPulserColor != fBlindPixel->GetColor())
483 {
484 *fLog << err << GetDescriptor()
485 << ": Pulser colour has changed w.r.t. last file in MCalibrationChargeBlindPix." << endl;
486 *fLog << err << "This feature is not yet implemented, sorry ... aborting " << endl;
487 return kFALSE;
488 }
489
490 if (fBlindCam)
491 if (fPulserColor != fBlindCam->GetColor())
492 {
493 *fLog << err << GetDescriptor()
494 << ": Pulser colour has changed w.r.t. last file in MCalibrationChargeBlindCam." << endl;
495 *fLog << err << "This feature is not yet implemented, sorry ... aborting " << endl;
496 return kFALSE;
497 }
498
499 if (fPINDiode)
500 if (fPulserColor != fPINDiode->GetColor())
501 {
502 *fLog << err << GetDescriptor()
503 << ": Pulser colour has changed w.r.t. last file in MCalibrationChargePINDiode." << endl;
504 *fLog << err << "This feature is not yet implemented, sorry ... aborting " << endl;
505 return kFALSE;
506 }
507
508
509 fNumHiGainSamples = fCam->GetNumHiGainFADCSlices();
510 fNumLoGainSamples = fCam->GetNumLoGainFADCSlices();
511
512 fSqrtHiGainSamples = TMath::Sqrt(fNumHiGainSamples);
513 fSqrtLoGainSamples = TMath::Sqrt(fNumLoGainSamples);
514
515 UInt_t npixels = fGeom->GetNumPixels();
516
517 for (UInt_t i=0; i<npixels; i++)
518 {
519
520 MCalibrationChargePix &pix = (MCalibrationChargePix&)(*fCam) [i];
521 MCalibrationQEPix &pqe = (MCalibrationQEPix&) (*fQECam)[i];
522 MBadPixelsPix &bad = (*fBadPixels)[i];
523
524 pix.SetPixId(i);
525 pqe.SetPixId(i);
526
527 if (bad.IsBad())
528 {
529 pix.SetExcluded();
530 pqe.SetExcluded();
531 continue;
532 }
533
534 if (IsDebug())
535 pix.SetDebug();
536 }
537
538 return kTRUE;
539}
540
541// ----------------------------------------------------------------------------------
542//
543// Nothing to be done in Process, but have a look at MHCalibrationChargeCam, instead
544//
545Int_t MCalibrationChargeCalc::Process()
546{
547 return kTRUE;
548}
549
550// -----------------------------------------------------------------------
551//
552// Return if number of executions is null.
553//
554// First loop over pixels, average areas and sectors, call:
555// - FinalizePedestals()
556// - FinalizeCharges()
557// for every entry. Count number of valid pixels in loop and return kFALSE
558// if there are none (the "Michele check").
559//
560// Call FinalizeBadPixels()
561//
562// Call FinalizeFFactorMethod() (second and third loop over pixels and areas)
563//
564// Call FinalizeBlindPixel()
565// Call FinalizePINDiode()
566//
567// Call FinalizeFFactorQECam() (fourth loop over pixels and areas)
568// Call FinalizeBlindPixelQECam() (fifth loop over pixels and areas)
569// Call FinalizePINDiodeQECam() (sixth loop over pixels and areas)
570//
571// Call FinalizeUnsuitablePixels()
572//
573// Call MParContainer::SetReadyToSave() for fCam, fQECam, fBadPixels and
574// fBlindPixel and fPINDiode if they exist
575//
576// Print out some statistics
577//
578Int_t MCalibrationChargeCalc::PostProcess()
579{
580
581 if (GetNumExecutions()==0)
582 return kFALSE;
583
584 *fLog << endl;
585
586 if (fPINDiode)
587 if (!fPINDiode->IsValid())
588 {
589 *fLog << warn << GetDescriptor()
590 << ": MCalibrationChargePINDiode is declared not valid... no PIN Diode method! " << endl;
591 fPINDiode = NULL;
592 }
593
594 if (fBlindPixel)
595 if (!fBlindPixel->IsValid())
596 {
597 *fLog << warn << GetDescriptor()
598 << ": MCalibrationChargeBlindPix is declared not valid... no Blind Pixel method! " << endl;
599 fBlindPixel = NULL;
600 }
601
602 if (fBlindCam)
603 if (!fBlindCam->IsValid())
604 {
605 *fLog << warn << GetDescriptor()
606 << ": MCalibrationChargeBlindCam is declared not valid... no Blind Pixel method! " << endl;
607 fBlindCam = NULL;
608 }
609
610 *fLog << endl;
611 //
612 // First loop over pixels, call FinalizePedestals and FinalizeCharges
613 //
614 Int_t nvalid = 0;
615
616 for (Int_t pixid=0; pixid<fPedestals->GetSize(); pixid++)
617 {
618
619 MCalibrationChargePix &pix = (MCalibrationChargePix&)(*fCam)[pixid];
620 //
621 // Check if the pixel has been excluded from the fits
622 //
623 if (pix.IsExcluded())
624 continue;
625
626 MPedestalPix &ped = (*fPedestals)[pixid];
627 MBadPixelsPix &bad = (*fBadPixels)[pixid];
628 const Int_t aidx = (*fGeom)[pixid].GetAidx();
629
630 FinalizePedestals(ped,pix,aidx);
631
632 if (FinalizeCharges(pix,bad,"pixel "))
633 nvalid++;
634 }
635
636 *fLog << endl;
637 //
638 // The Michele check ...
639 //
640 if (nvalid == 0)
641 {
642 *fLog << err << GetDescriptor() << ": All pixels have non-valid calibration. "
643 << "Did you forget to fill the histograms "
644 << "(filling MHCalibrationChargeCam from MExtractedSignalCam using MFillH) ? " << endl;
645 *fLog << err << GetDescriptor() << ": Or, maybe, you have used a pedestal run "
646 << "instead of a calibration run " << endl;
647 return kFALSE;
648 }
649
650 for (UInt_t aidx=0; aidx<fGeom->GetNumAreas(); aidx++)
651 {
652
653 const MPedestalPix &ped = fPedestals->GetAverageArea(aidx);
654 MCalibrationChargePix &pix = (MCalibrationChargePix&)fCam->GetAverageArea(aidx);
655
656 FinalizePedestals(ped,pix,aidx);
657 FinalizeCharges(pix, fCam->GetAverageBadArea(aidx),"area id");
658 }
659
660 *fLog << endl;
661
662 for (UInt_t sector=0; sector<fGeom->GetNumSectors(); sector++)
663 {
664
665 const MPedestalPix &ped = fPedestals->GetAverageSector(sector);
666
667 MCalibrationChargePix &pix = (MCalibrationChargePix&)fCam->GetAverageSector(sector);
668 FinalizePedestals(ped,pix, 0);
669 }
670
671 *fLog << endl;
672
673 //
674 // Finalize Bad Pixels
675 //
676 FinalizeBadPixels();
677
678 //
679 // Finalize F-Factor method
680 //
681 if (!FinalizeFFactorMethod())
682 {
683 *fLog << warn << "Could not calculate the photons flux from the F-Factor method " << endl;
684 fCam->SetFFactorMethodValid(kFALSE);
685 return kFALSE;
686 }
687 else
688 fCam->SetFFactorMethodValid(kTRUE);
689
690 *fLog << endl;
691 //
692 // Finalize Blind Pixel
693 //
694 if (FinalizeBlindPixel())
695 fQECam->SetBlindPixelMethodValid(kTRUE);
696 else
697 fQECam->SetBlindPixelMethodValid(kFALSE);
698
699 //
700 // Finalize PIN Diode
701 //
702 if (FinalizePINDiode())
703 fQECam->SetPINDiodeMethodValid(kTRUE);
704 else
705 fQECam->SetPINDiodeMethodValid(kFALSE);
706
707 //
708 // Finalize QE Cam
709 //
710 FinalizeFFactorQECam();
711 FinalizeBlindPixelQECam();
712 FinalizePINDiodeQECam();
713
714 //
715 // Re-direct the output to an ascii-file from now on:
716 //
717 MLog asciilog;
718 asciilog.SetOutputFile(GetOutputFile(),kTRUE);
719 SetLogStream(&asciilog);
720 //
721 // Finalize calibration statistics
722 //
723 FinalizeUnsuitablePixels();
724
725 fCam ->SetReadyToSave();
726 fQECam ->SetReadyToSave();
727 fBadPixels->SetReadyToSave();
728
729 if (fBlindPixel)
730 fBlindPixel->SetReadyToSave();
731 if (fBlindCam)
732 fBlindCam->SetReadyToSave();
733 if (fPINDiode)
734 fPINDiode->SetReadyToSave();
735
736 *fLog << inf << endl;
737 *fLog << GetDescriptor() << ": Fatal errors statistics:" << endl;
738
739 PrintUncalibrated(MBadPixelsPix::kChargeIsPedestal,
740 Form("%s%2.1f%s","Signal less than ",fChargeLimit," Pedestal RMS: "));
741 PrintUncalibrated(MBadPixelsPix::kChargeRelErrNotValid,
742 Form("%s%2.1f%s","Signal Error bigger than ",fChargeRelErrLimit," times Mean Signal: "));
743 PrintUncalibrated(MBadPixelsPix::kChargeSigmaNotValid,
744 "Signal Sigma smaller than Pedestal RMS: ");
745 PrintUncalibrated(MBadPixelsPix::kLoGainSaturation,
746 "Pixels with Low Gain Saturation: ");
747 PrintUncalibrated(MBadPixelsPix::kMeanTimeInFirstBin,
748 Form("%s%2.1f%s","Mean Abs. Arr. Time in First ",1.," Bin(s): "));
749 PrintUncalibrated(MBadPixelsPix::kMeanTimeInLast2Bins,
750 Form("%s%2.1f%s","Mean Abs. Arr. Time in Last ",2.," Bin(s): "));
751 PrintUncalibrated(MBadPixelsPix::kDeviatingNumPhes,
752 "Pixels with deviating number of phes: ");
753 PrintUncalibrated(MBadPixelsPix::kDeviatingFFactor,
754 "Pixels with deviating F-Factor: ");
755
756 *fLog << inf << endl;
757 *fLog << GetDescriptor() << ": Unreliable errors statistics:" << endl;
758
759 PrintUncalibrated(MBadPixelsPix::kHiGainOscillating,
760 "Pixels with changing Hi Gain signal over time: ");
761 PrintUncalibrated(MBadPixelsPix::kLoGainOscillating,
762 "Pixels with changing Lo Gain signal over time: ");
763 PrintUncalibrated(MBadPixelsPix::kHiGainNotFitted,
764 "Pixels with unsuccesful Gauss fit to the Hi Gain: ");
765 PrintUncalibrated(MBadPixelsPix::kLoGainNotFitted,
766 "Pixels with unsuccesful Gauss fit to the Lo Gain: ");
767
768 SetLogStream(&gLog);
769
770 return kTRUE;
771}
772
773// ----------------------------------------------------------------------------------
774//
775// Retrieves pedestal and pedestal RMS from MPedestalPix
776// Retrieves total entries from MPedestalCam
777// Sets pedestal*fNumHiGainSamples and pedestal*fNumLoGainSamples in MCalibrationChargePix
778// Sets pedRMS *fSqrtHiGainSamples and pedRMS *fSqrtLoGainSamples in MCalibrationChargePix
779//
780// If the flag MCalibrationPix::IsHiGainSaturation() is set, call also:
781// - MCalibrationChargePix::CalcLoGainPedestal()
782//
783void MCalibrationChargeCalc::FinalizePedestals(const MPedestalPix &ped, MCalibrationChargePix &cal, const Int_t aidx)
784{
785
786 //
787 // get the pedestals
788 //
789 const Float_t pedes = ped.GetPedestal();
790 const Float_t prms = ped.GetPedestalRms();
791 const Float_t num = TMath::Sqrt((Float_t)fPedestals->GetTotalEntries());
792
793
794 //
795 // set them in the calibration camera
796 //
797 if (cal.IsHiGainSaturation())
798 {
799 cal.SetPedestal(pedes* fNumLoGainSamples,
800 prms * fSqrtLoGainSamples,
801 prms * fNumLoGainSamples / num);
802 cal.CalcLoGainPedestal((Float_t)fNumLoGainSamples, aidx);
803 }
804 else
805 {
806 cal.SetPedestal(pedes* fNumHiGainSamples,
807 prms * fSqrtHiGainSamples,
808 prms * fNumHiGainSamples / num);
809 }
810
811}
812
813// ----------------------------------------------------------------------------------------------------
814//
815// Check fit results validity. Bad Pixels flags are set if:
816//
817// 1) Pixel has a mean smaller than fChargeLimit*PedRMS ( Flag: MBadPixelsPix::kChargeIsPedestal)
818// 2) Pixel has a mean error smaller than fChargeErrLimit ( Flag: MBadPixelsPix::kChargeErrNotValid)
819// 3) Pixel has mean smaller than fChargeRelVarLimit times its mean error
820// ( Flag: MBadPixelsPix::kChargeRelErrNotValid)
821// 4) Pixel has a sigma bigger than its Pedestal RMS ( Flag: MBadPixelsPix::kChargeSigmaNotValid )
822//
823// Further returns if flags: MBadPixelsPix::kUnsuitableRun is set
824//
825// Calls MCalibrationChargePix::CalcReducedSigma() and sets flag: MBadPixelsPix::kChargeIsPedestal
826// if not succesful.
827//
828// Calls MCalibrationChargePix::CalcFFactorMethod() and sets flag: MBadPixelsPix::kDeviatingNumPhes)
829// if not succesful.
830//
831Bool_t MCalibrationChargeCalc::FinalizeCharges(MCalibrationChargePix &cal, MBadPixelsPix &bad, const char* what)
832{
833
834 if (bad.IsUnsuitable(MBadPixelsPix::kUnsuitableRun))
835 return kFALSE;
836
837 if (cal.GetMean() < fChargeLimit*cal.GetPedRms())
838 {
839 *fLog << warn << GetDescriptor()
840 << Form(": Fitted Charge: %5.2f is smaller than %2.1f",cal.GetMean(),fChargeLimit)
841 << Form(" Pedestal RMS: %5.2f in %s%4i",cal.GetPedRms(),what,cal.GetPixId()) << endl;
842 bad.SetUncalibrated( MBadPixelsPix::kChargeIsPedestal);
843 }
844
845 if (cal.GetMean() < fChargeRelErrLimit*cal.GetMeanErr())
846 {
847 *fLog << warn << GetDescriptor()
848 << Form(": Fitted Charge: %4.2f is smaller than %2.1f",cal.GetMean(),fChargeRelErrLimit)
849 << Form(" times its error: %4.2f in %s%4i",cal.GetMeanErr(),what,cal.GetPixId()) << endl;
850 bad.SetUncalibrated( MBadPixelsPix::kChargeRelErrNotValid );
851 }
852
853 if (cal.GetSigma() < cal.GetPedRms())
854 {
855 *fLog << warn << GetDescriptor()
856 << Form(": Sigma of Fitted Charge: %6.2f is smaller than",cal.GetSigma())
857 << Form(" Ped. RMS: %5.2f in %s%4i",cal.GetPedRms(),what,cal.GetPixId()) << endl;
858 bad.SetUncalibrated( MBadPixelsPix::kChargeSigmaNotValid );
859 }
860
861 if (!cal.CalcReducedSigma())
862 {
863 *fLog << warn << GetDescriptor()
864 << Form(": Could not calculate the reduced sigma in %s: ",what)
865 << Form(" %4i",cal.GetPixId())
866 << endl;
867 bad.SetUncalibrated(MBadPixelsPix::kChargeIsPedestal);
868 return kFALSE;
869 }
870
871 if (!cal.CalcFFactorMethod())
872 {
873 *fLog << warn << GetDescriptor()
874 << Form(": Could not calculate the F-Factor in %s: ",what)
875 << Form(" %4i",cal.GetPixId())
876 << endl;
877 bad.SetUncalibrated(MBadPixelsPix::kDeviatingNumPhes);
878 return kFALSE;
879 }
880
881 return kTRUE;
882}
883
884
885
886// -----------------------------------------------------------------------------------
887//
888// Sets pixel to MBadPixelsPix::kUnsuitableRun, if one of the following flags is set:
889// - MBadPixelsPix::kChargeIsPedestal
890// - MBadPixelsPix::kChargeRelErrNotValid
891// - MBadPixelsPix::kChargeSigmaNotValid
892// - MBadPixelsPix::kMeanTimeInFirstBin
893// - MBadPixelsPix::kMeanTimeInLast2Bins
894// - MBadPixelsPix::kDeviatingNumPhes
895//
896// - Call MCalibrationPix::SetExcluded() for the bad pixels
897//
898void MCalibrationChargeCalc::FinalizeBadPixels()
899{
900
901 for (Int_t i=0; i<fBadPixels->GetSize(); i++)
902 {
903
904 MBadPixelsPix &bad = (*fBadPixels)[i];
905 MCalibrationPix &pix = (*fCam)[i];
906
907 if (bad.IsUncalibrated( MBadPixelsPix::kChargeIsPedestal))
908 bad.SetUnsuitable( MBadPixelsPix::kUnsuitableRun );
909
910 if (bad.IsUncalibrated( MBadPixelsPix::kChargeErrNotValid ))
911 bad.SetUnsuitable( MBadPixelsPix::kUnsuitableRun );
912
913 if (bad.IsUncalibrated( MBadPixelsPix::kChargeRelErrNotValid ))
914 bad.SetUnsuitable( MBadPixelsPix::kUnsuitableRun );
915
916 if (bad.IsUncalibrated( MBadPixelsPix::kChargeSigmaNotValid ))
917 bad.SetUnsuitable( MBadPixelsPix::kUnsuitableRun );
918
919 if (bad.IsUncalibrated( MBadPixelsPix::kMeanTimeInFirstBin ))
920 bad.SetUnsuitable( MBadPixelsPix::kUnsuitableRun );
921
922 if (bad.IsUncalibrated( MBadPixelsPix::kMeanTimeInLast2Bins ))
923 bad.SetUnsuitable( MBadPixelsPix::kUnsuitableRun );
924
925 if (bad.IsUncalibrated( MBadPixelsPix::kDeviatingNumPhes ))
926 bad.SetUnsuitable( MBadPixelsPix::kUnsuitableRun );
927
928 if (bad.IsUnsuitable( MBadPixelsPix::kUnsuitableRun ))
929 pix.SetExcluded();
930
931 }
932}
933
934// ------------------------------------------------------------------------
935//
936//
937// First loop: Calculate a mean and mean RMS of photo-electrons per area index
938// Include only pixels which are not MBadPixelsPix::kUnsuitableRun or
939// MBadPixelsPix::kUnreliableRun (see FinalizeBadPixels()) and set
940// MCalibrationChargePix::SetFFactorMethodValid(kFALSE) in that case.
941//
942// Second loop: Get weighted mean number of photo-electrons and its RMS including
943// only pixels with flag MCalibrationChargePix::IsFFactorMethodValid()
944// and further exclude those deviating by more than fPheErrLimit mean
945// sigmas from the mean (obtained in first loop). Set
946// MBadPixelsPix::kDeviatingNumPhes if excluded.
947//
948// Set weighted mean and variance of photo-electrons per area index in:
949// average area pixels of MCalibrationChargeCam (obtained from:
950// MCalibrationChargeCam::GetAverageArea() )
951//
952// Set weighted mean and variance of photo-electrons per sector in:
953// average sector pixels of MCalibrationChargeCam (obtained from:
954// MCalibrationChargeCam::GetAverageSector() )
955//
956Bool_t MCalibrationChargeCalc::FinalizeFFactorMethod()
957{
958
959 const UInt_t npixels = fGeom->GetNumPixels();
960 const UInt_t nareas = fGeom->GetNumAreas();
961 const UInt_t nsectors = fGeom->GetNumSectors();
962
963 Float_t lowlim [nareas];
964 Float_t upplim [nareas];
965 Double_t areavars [nareas];
966 Double_t areaweights[nareas], sectorweights [nsectors];
967 Double_t areaphes [nareas], sectorphes [nsectors];
968 Int_t numareavalid[nareas], numsectorvalid[nsectors];
969
970 memset(lowlim ,0, nareas * sizeof(Float_t));
971 memset(upplim ,0, nareas * sizeof(Float_t));
972 memset(areaphes ,0, nareas * sizeof(Double_t));
973 memset(areavars ,0, nareas * sizeof(Double_t));
974 memset(areaweights ,0, nareas * sizeof(Double_t));
975 memset(numareavalid ,0, nareas * sizeof(Int_t ));
976 memset(sectorweights ,0, nsectors * sizeof(Double_t));
977 memset(sectorphes ,0, nsectors * sizeof(Double_t));
978 memset(numsectorvalid,0, nsectors * sizeof(Int_t ));
979
980 //
981 // First loop: Get mean number of photo-electrons and the RMS
982 // The loop is only to recognize later pixels with very deviating numbers
983 //
984 MHCamera camphes(*fGeom,"Camphes","Phes in Camera");
985
986 for (UInt_t i=0; i<npixels; i++)
987 {
988
989 MCalibrationChargePix &pix = (MCalibrationChargePix&)(*fCam) [i];
990 MBadPixelsPix &bad = (*fBadPixels)[i];
991
992 if (!pix.IsFFactorMethodValid())
993 continue;
994
995 if (bad.IsUnsuitable(MBadPixelsPix::kUnsuitableRun))
996 {
997 pix.SetFFactorMethodValid(kFALSE);
998 continue;
999 }
1000
1001 // if (bad.IsUnsuitable(MBadPixelsPix::kUnreliableRun))
1002 // continue;
1003
1004 const Float_t nphe = pix.GetPheFFactorMethod();
1005 const Int_t aidx = (*fGeom)[i].GetAidx();
1006
1007 camphes.Fill(i,nphe);
1008 camphes.SetUsed(i);
1009
1010 areaphes [aidx] += nphe;
1011 areavars [aidx] += nphe*nphe;
1012 numareavalid[aidx] ++;
1013 }
1014
1015 for (UInt_t i=0; i<nareas; i++)
1016 {
1017 if (numareavalid[i] == 0)
1018 {
1019 *fLog << warn << GetDescriptor() << ": No pixels with valid number of photo-electrons found "
1020 << "in area index: " << i << endl;
1021 continue;
1022 }
1023
1024 if (numareavalid[i] == 1)
1025 areavars[i] = 0.;
1026 else if (numareavalid[i] == 0)
1027 {
1028 areaphes[i] = -1.;
1029 areaweights[i] = -1.;
1030 }
1031 else
1032 {
1033 areavars[i] = (areavars[i] - areaphes[i]*areaphes[i]/numareavalid[i]) / (numareavalid[i]-1);
1034 areaphes[i] = areaphes[i] / numareavalid[i];
1035 }
1036
1037 if (areavars[i] < 0.)
1038 {
1039 *fLog << warn << GetDescriptor() << ": No pixels with valid variance of photo-electrons found "
1040 << "in area index: " << i << endl;
1041 continue;
1042 }
1043
1044 lowlim [i] = areaphes[i] - fPheErrLimit*TMath::Sqrt(areavars[i]);
1045 upplim [i] = areaphes[i] + fPheErrLimit*TMath::Sqrt(areavars[i]);
1046
1047 TArrayI area(1);
1048 area[0] = i;
1049
1050 TH1D *hist = camphes.ProjectionS(TArrayI(),area,"_py",100);
1051 hist->Fit("gaus","Q");
1052 const Float_t mean = hist->GetFunction("gaus")->GetParameter(1);
1053 const Float_t sigma = hist->GetFunction("gaus")->GetParameter(2);
1054 const Int_t ndf = hist->GetFunction("gaus")->GetNDF();
1055
1056 if (IsDebug())
1057 camphes.DrawClone();
1058
1059 if (ndf < 2)
1060 {
1061 *fLog << warn << GetDescriptor() << ": Cannot use a Gauss fit to the number of photo-electrons "
1062 << "in the camera with area index: " << i << endl;
1063 *fLog << warn << GetDescriptor() << ": Number of dof.: " << ndf << " is smaller than 2 " << endl;
1064 *fLog << warn << GetDescriptor() << ": Will use the simple mean and rms " << endl;
1065 delete hist;
1066 continue;
1067 }
1068
1069 const Double_t prob = hist->GetFunction("gaus")->GetProb();
1070
1071 if (prob < 0.001)
1072 {
1073 *fLog << warn << GetDescriptor() << ": Cannot use a Gauss fit to the number of photo-electrons "
1074 << "in the camera with area index: " << i << endl;
1075 *fLog << warn << GetDescriptor() << ": Fit probability " << prob
1076 << " is smaller than 0.001 " << endl;
1077 *fLog << warn << GetDescriptor() << ": Will use the simple mean and rms " << endl;
1078 delete hist;
1079 continue;
1080 }
1081
1082 *fLog << inf << GetDescriptor() << ": Mean number of photo-electrons "
1083 << "with area idx " << i << ": "
1084 << Form("%7.2f+-%6.2f",mean,sigma) << endl;
1085
1086 lowlim [i] = mean - fPheErrLimit*sigma;
1087 upplim [i] = mean + fPheErrLimit*sigma;
1088
1089 delete hist;
1090 }
1091
1092 *fLog << endl;
1093
1094 memset(numareavalid,0,nareas*sizeof(Int_t));
1095 memset(areaphes ,0,nareas*sizeof(Double_t));
1096 memset(areavars ,0,nareas*sizeof(Double_t));
1097
1098 //
1099 // Second loop: Get mean number of photo-electrons and its RMS excluding
1100 // pixels deviating by more than fPheErrLimit sigma.
1101 // Set the conversion factor FADC counts to photo-electrons
1102 //
1103 for (UInt_t i=0; i<npixels; i++)
1104 {
1105
1106 MCalibrationChargePix &pix = (MCalibrationChargePix&)(*fCam)[i];
1107
1108 if (!pix.IsFFactorMethodValid())
1109 continue;
1110
1111 const Float_t nvar = pix.GetPheFFactorMethodVar();
1112
1113 if (nvar <= 0.)
1114 {
1115 pix.SetFFactorMethodValid(kFALSE);
1116 continue;
1117 }
1118
1119 MBadPixelsPix &bad = (*fBadPixels)[i];
1120
1121 const Int_t aidx = (*fGeom)[i].GetAidx();
1122 const Int_t sector = (*fGeom)[i].GetSector();
1123 const Float_t area = (*fGeom)[i].GetA();
1124 const Float_t nphe = pix.GetPheFFactorMethod();
1125
1126 if ( nphe < lowlim[aidx] || nphe > upplim[aidx] )
1127 {
1128 *fLog << warn << GetDescriptor() << ": Number of phes: "
1129 << Form("%7.2f out of %3.1f sigma limit: ",nphe,fPheErrLimit)
1130 << Form("[%7.2f,%7.2f] pixel%4i",lowlim[aidx],upplim[aidx],i) << endl;
1131 bad.SetUncalibrated( MBadPixelsPix::kDeviatingNumPhes );
1132 bad.SetUnsuitable ( MBadPixelsPix::kUnsuitableRun );
1133 pix.SetFFactorMethodValid(kFALSE);
1134 continue;
1135 }
1136
1137 areaweights [aidx] += nphe*nphe;
1138 areaphes [aidx] += nphe;
1139 numareavalid [aidx] ++;
1140
1141 if (aidx == 0)
1142 fNumInnerFFactorMethodUsed++;
1143
1144 sectorweights [sector] += nphe*nphe/area/area;
1145 sectorphes [sector] += nphe/area;
1146 numsectorvalid[sector] ++;
1147 }
1148
1149 *fLog << endl;
1150
1151 for (UInt_t aidx=0; aidx<nareas; aidx++)
1152 {
1153
1154 MCalibrationChargePix &apix = (MCalibrationChargePix&)fCam->GetAverageArea(aidx);
1155
1156 if (numareavalid[aidx] == 1)
1157 areaweights[aidx] = 0.;
1158 else if (numareavalid[aidx] == 0)
1159 {
1160 areaphes[aidx] = -1.;
1161 areaweights[aidx] = -1.;
1162 }
1163 else
1164 {
1165 areaweights[aidx] = (areaweights[aidx] - areaphes[aidx]*areaphes[aidx]/numareavalid[aidx])
1166 / (numareavalid[aidx]-1);
1167 areaphes[aidx] /= numareavalid[aidx];
1168 }
1169
1170 if (areaweights[aidx] < 0. || areaphes[aidx] <= 0.)
1171 {
1172 *fLog << warn << GetDescriptor()
1173 << ": Mean number phes from area index " << aidx << " could not be calculated: "
1174 << " Mean: " << areaphes[aidx]
1175 << " Variance: " << areaweights[aidx] << endl;
1176 apix.SetFFactorMethodValid(kFALSE);
1177 continue;
1178 }
1179
1180 *fLog << inf << "Average total number phes in area idx " << aidx << ": "
1181 << Form("%7.2f%s%6.2f",areaphes[aidx]," +- ",TMath::Sqrt(areaweights[aidx])) << endl;
1182
1183 apix.SetPheFFactorMethod ( areaphes[aidx] );
1184 apix.SetPheFFactorMethodVar( areaweights[aidx] / numareavalid[aidx] );
1185 apix.SetFFactorMethodValid ( kTRUE );
1186
1187 }
1188
1189 *fLog << endl;
1190
1191 for (UInt_t sector=0; sector<nsectors; sector++)
1192 {
1193
1194 if (numsectorvalid[sector] == 1)
1195 sectorweights[sector] = 0.;
1196 else if (numsectorvalid[sector] == 0)
1197 {
1198 sectorphes[sector] = -1.;
1199 sectorweights[sector] = -1.;
1200 }
1201 else
1202 {
1203 sectorweights[sector] = (sectorweights[sector]
1204 - sectorphes[sector]*sectorphes[sector]/numsectorvalid[sector]
1205 )
1206 / (numsectorvalid[sector]-1.);
1207 sectorphes[sector] /= numsectorvalid[sector];
1208 }
1209
1210 MCalibrationChargePix &spix = (MCalibrationChargePix&)fCam->GetAverageSector(sector);
1211
1212 if (sectorweights[sector] < 0. || sectorphes[sector] <= 0.)
1213 {
1214 *fLog << warn << GetDescriptor()
1215 <<": Mean number phes per area for sector " << sector << " could not be calculated: "
1216 << " Mean: " << sectorphes[sector]
1217 << " Variance: " << sectorweights[sector] << endl;
1218 spix.SetFFactorMethodValid(kFALSE);
1219 continue;
1220 }
1221
1222 *fLog << inf << "Average number phes per area in sector " << sector << ": "
1223 << Form("%5.2f+-%4.2f [phe/mm^2]",sectorphes[sector],TMath::Sqrt(sectorweights[sector]))
1224 << endl;
1225
1226 spix.SetPheFFactorMethod ( sectorphes[sector] );
1227 spix.SetPheFFactorMethodVar( sectorweights[sector] / numsectorvalid[sector]);
1228 spix.SetFFactorMethodValid ( kTRUE );
1229
1230 }
1231
1232 return kTRUE;
1233}
1234
1235
1236// ------------------------------------------------------------------------
1237//
1238// Returns kFALSE if pointer to MCalibrationChargeBlindPix is NULL
1239//
1240// The check returns kFALSE if:
1241//
1242// 1) fLambda and fLambdaCheck are separated relatively to each other by more than fLambdaCheckLimit
1243// 2) BlindPixel has an fLambdaErr greater than fLambdaErrLimit
1244//
1245// Calls:
1246// - MCalibrationChargeBlindPix::CalcFluxInsidePlexiglass()
1247//
1248Bool_t MCalibrationChargeCalc::FinalizeBlindPixel()
1249{
1250
1251 if (!fBlindPixel)
1252 return kFALSE;
1253
1254 const Float_t lambda = fBlindPixel->GetLambda();
1255 const Float_t lambdaerr = fBlindPixel->GetLambdaErr();
1256 const Float_t lambdacheck = fBlindPixel->GetLambdaCheck();
1257
1258 if (2.*(lambdacheck-lambda)/(lambdacheck+lambda) > fLambdaCheckLimit)
1259 {
1260 *fLog << warn << GetDescriptor()
1261 << Form("%s%4.2f%s%4.2f%s%4.2f%s",": Lambda: ",lambda," and Lambda-Check: ",
1262 lambdacheck," differ by more than ",fLambdaCheckLimit," in the Blind Pixel ")
1263 << endl;
1264 return kFALSE;
1265 }
1266
1267 if (lambdaerr > fLambdaErrLimit)
1268 {
1269 *fLog << warn << GetDescriptor()
1270 << Form("%s%4.2f%s%4.2f%s",": Error of Fitted Lambda: ",lambdaerr," is greater than ",
1271 fLambdaErrLimit," in Blind Pixel ") << endl;
1272 return kFALSE;
1273 }
1274
1275 if (!fBlindPixel->CalcFluxInsidePlexiglass())
1276 {
1277 *fLog << warn << "Could not calculate the flux of photons from the Blind Pixel, "
1278 << "will skip Blind Pixel Calibration " << endl;
1279 return kFALSE;
1280 }
1281
1282 return kTRUE;
1283}
1284
1285// ------------------------------------------------------------------------
1286//
1287// Returns kFALSE if pointer to MCalibrationChargePINDiode is NULL
1288//
1289// The check returns kFALSE if:
1290//
1291// 1) PINDiode has a fitted charge smaller than fChargeLimit*PedRMS
1292// 2) PINDiode has a fit error smaller than fChargeErrLimit
1293// 3) PINDiode has a fitted charge smaller its fChargeRelErrLimit times its charge error
1294// 4) PINDiode has a charge sigma smaller than its Pedestal RMS
1295//
1296// Calls:
1297// - MCalibrationChargePINDiode::CalcFluxOutsidePlexiglass()
1298//
1299Bool_t MCalibrationChargeCalc::FinalizePINDiode()
1300{
1301
1302 if (!fPINDiode)
1303 return kFALSE;
1304
1305 if (fPINDiode->GetMean() < fChargeLimit*fPINDiode->GetPedRms())
1306 {
1307 *fLog << warn << GetDescriptor() << ": Fitted Charge is smaller than "
1308 << fChargeLimit << " Pedestal RMS in PINDiode " << endl;
1309 return kFALSE;
1310 }
1311
1312 if (fPINDiode->GetMeanErr() < fChargeErrLimit)
1313 {
1314 *fLog << warn << GetDescriptor() << ": Error of Fitted Charge is smaller than "
1315 << fChargeErrLimit << " in PINDiode " << endl;
1316 return kFALSE;
1317 }
1318
1319 if (fPINDiode->GetMean() < fChargeRelErrLimit*fPINDiode->GetMeanErr())
1320 {
1321 *fLog << warn << GetDescriptor() << ": Fitted Charge is smaller than "
1322 << fChargeRelErrLimit << "* its error in PINDiode " << endl;
1323 return kFALSE;
1324 }
1325
1326 if (fPINDiode->GetSigma() < fPINDiode->GetPedRms())
1327 {
1328 *fLog << warn << GetDescriptor()
1329 << ": Sigma of Fitted Charge smaller than Pedestal RMS in PINDiode " << endl;
1330 return kFALSE;
1331 }
1332
1333
1334 if (!fPINDiode->CalcFluxOutsidePlexiglass())
1335 {
1336 *fLog << warn << "Could not calculate the flux of photons from the PIN Diode, "
1337 << "will skip PIN Diode Calibration " << endl;
1338 return kFALSE;
1339 }
1340
1341 return kTRUE;
1342}
1343
1344// ------------------------------------------------------------------------
1345//
1346// Calculate the average number of photons outside the plexiglass with the
1347// formula:
1348//
1349// av.Num.photons(area index) = av.Num.Phes(area index)
1350// / MCalibrationQEPix::GetDefaultQE(fPulserColor)
1351// / MCalibrationQEPix::GetPMTCollectionEff()
1352// / MCalibrationQEPix::GetLightGuidesEff(fPulserColor)
1353// / MCalibrationQECam::GetPlexiglassQE()
1354//
1355// Calculate the variance on the average number of photons assuming that the error on the
1356// Quantum efficiency is reduced by the number of used inner pixels, but the rest of the
1357// values keeps it ordinary error since it is systematic.
1358//
1359// Loop over pixels:
1360//
1361// - Continue, if not MCalibrationChargePix::IsFFactorMethodValid() and set:
1362// MCalibrationQEPix::SetFFactorMethodValid(kFALSE,fPulserColor)
1363//
1364// - Call MCalibrationChargePix::CalcMeanFFactor(av.Num.photons) and set:
1365// MCalibrationQEPix::SetFFactorMethodValid(kFALSE,fPulserColor) if not succesful
1366//
1367// - Calculate the quantum efficiency with the formula:
1368//
1369// QE = ( Num.Phes / av.Num.photons ) * MGeomCam::GetPixRatio()
1370//
1371// - Set QE in MCalibrationQEPix::SetQEFFactor ( QE, fPulserColor );
1372//
1373// - Set Variance of QE in MCalibrationQEPix::SetQEFFactorVar ( Variance, fPulserColor );
1374// - Set bit MCalibrationQEPix::SetFFactorMethodValid(kTRUE,fPulserColor)
1375//
1376// - Call MCalibrationQEPix::UpdateFFactorMethod()
1377//
1378void MCalibrationChargeCalc::FinalizeFFactorQECam()
1379{
1380
1381 if (fNumInnerFFactorMethodUsed < 2)
1382 {
1383 *fLog << warn << GetDescriptor()
1384 << ": Could not calculate F-Factor Method: Less than 2 inner pixels valid! " << endl;
1385 return;
1386 }
1387
1388 MCalibrationChargePix &avpix = (MCalibrationChargePix&)fCam->GetAverageArea(0);
1389 MCalibrationQEPix &qepix = (MCalibrationQEPix&) fQECam->GetAverageArea(0);
1390
1391 const Float_t avphotons = avpix.GetPheFFactorMethod()
1392 / qepix.GetDefaultQE(fPulserColor)
1393 / qepix.GetPMTCollectionEff()
1394 / qepix.GetLightGuidesEff(fPulserColor)
1395 / fQECam->GetPlexiglassQE();
1396
1397 const Float_t avphotrelvar = avpix.GetPheFFactorMethodRelVar()
1398 + qepix.GetDefaultQERelVar(fPulserColor) / fNumInnerFFactorMethodUsed
1399 + qepix.GetPMTCollectionEffRelVar()
1400 + qepix.GetLightGuidesEffRelVar(fPulserColor)
1401 + fQECam->GetPlexiglassQERelVar();
1402
1403 const UInt_t nareas = fGeom->GetNumAreas();
1404
1405 //
1406 // Set the results in the MCalibrationChargeCam
1407 //
1408 fCam->SetNumPhotonsFFactorMethod (avphotons);
1409 if (avphotrelvar > 0.)
1410 fCam->SetNumPhotonsFFactorMethodErr(TMath::Sqrt( avphotrelvar * avphotons * avphotons));
1411
1412 Float_t lowlim [nareas];
1413 Float_t upplim [nareas];
1414 Double_t avffactorphotons [nareas];
1415 Double_t avffactorphotvar [nareas];
1416 Int_t numffactor [nareas];
1417
1418 memset(lowlim ,0, nareas * sizeof(Float_t));
1419 memset(upplim ,0, nareas * sizeof(Float_t));
1420 memset(avffactorphotons,0, nareas * sizeof(Double_t));
1421 memset(avffactorphotvar,0, nareas * sizeof(Double_t));
1422 memset(numffactor ,0, nareas * sizeof(Int_t));
1423
1424 const UInt_t npixels = fGeom->GetNumPixels();
1425
1426 MHCamera camffactor(*fGeom,"Camffactor","F-Factor in Camera");
1427
1428 for (UInt_t i=0; i<npixels; i++)
1429 {
1430
1431 MCalibrationChargePix &pix = (MCalibrationChargePix&)(*fCam)[i];
1432 MCalibrationQEPix &qepix = (MCalibrationQEPix&) (*fQECam)[i];
1433 MBadPixelsPix &bad = (*fBadPixels)[i];
1434
1435 if (bad.IsUnsuitable(MBadPixelsPix::kUnsuitableRun))
1436 continue;
1437
1438 const Float_t photons = avphotons / fGeom->GetPixRatio(i);
1439 const Float_t qe = pix.GetPheFFactorMethod() / photons ;
1440
1441 if (!pix.CalcMeanFFactor( photons , avphotrelvar ))
1442 {
1443 (*fBadPixels)[i].SetUncalibrated( MBadPixelsPix::kDeviatingNumPhes );
1444 continue;
1445 }
1446
1447 const Float_t qerelvar = avphotrelvar + pix.GetPheFFactorMethodRelVar();
1448
1449 qepix.SetQEFFactor ( qe , fPulserColor );
1450 qepix.SetQEFFactorVar ( qerelvar*qe*qe, fPulserColor );
1451 qepix.SetFFactorMethodValid( kTRUE , fPulserColor );
1452
1453 if (!qepix.UpdateFFactorMethod())
1454 *fLog << warn << GetDescriptor()
1455 << ": Cannot update Quantum efficiencies with the F-Factor Method" << endl;
1456
1457 const Int_t aidx = (*fGeom)[i].GetAidx();
1458 const Float_t ffactor = pix.GetMeanFFactorFADC2Phot();
1459
1460 camffactor.Fill(i,ffactor);
1461 camffactor.SetUsed(i);
1462
1463 avffactorphotons[aidx] += ffactor;
1464 avffactorphotvar[aidx] += ffactor*ffactor;
1465 numffactor[aidx]++;
1466 }
1467
1468 for (UInt_t i=0; i<nareas; i++)
1469 {
1470
1471 if (numffactor[i] == 0)
1472 {
1473 *fLog << warn << GetDescriptor() << ": No pixels with valid total F-Factor found "
1474 << "in area index: " << i << endl;
1475 continue;
1476 }
1477
1478 avffactorphotvar[i] = (avffactorphotvar[i] - avffactorphotons[i]*avffactorphotons[i]/numffactor[i]) / (numffactor[i]-1.);
1479 avffactorphotons[i] = avffactorphotons[i] / numffactor[i];
1480
1481 if (avffactorphotvar[i] < 0.)
1482 {
1483 *fLog << warn << GetDescriptor() << ": No pixels with valid variance of total F-Factor found "
1484 << "in area index: " << i << endl;
1485 continue;
1486 }
1487
1488 lowlim [i] = 1.1; // Lowest known F-Factor of a PMT
1489 upplim [i] = avffactorphotons[i] + fFFactorErrLimit*TMath::Sqrt(avffactorphotvar[i]);
1490
1491 TArrayI area(1);
1492 area[0] = i;
1493
1494 TH1D *hist = camffactor.ProjectionS(TArrayI(),area,"_py",100);
1495 hist->Fit("gaus","Q");
1496 const Float_t mean = hist->GetFunction("gaus")->GetParameter(1);
1497 const Float_t sigma = hist->GetFunction("gaus")->GetParameter(2);
1498 const Int_t ndf = hist->GetFunction("gaus")->GetNDF();
1499
1500 if (IsDebug())
1501 camffactor.DrawClone();
1502
1503 if (ndf < 2)
1504 {
1505 *fLog << warn << GetDescriptor() << ": Cannot use a Gauss fit to the F-Factor "
1506 << "in the camera with area index: " << i << endl;
1507 *fLog << warn << GetDescriptor() << ": Number of dof.: " << ndf << " is smaller than 2 " << endl;
1508 *fLog << warn << GetDescriptor() << ": Will use the simple mean and rms " << endl;
1509 delete hist;
1510 continue;
1511 }
1512
1513 const Double_t prob = hist->GetFunction("gaus")->GetProb();
1514
1515 if (prob < 0.001)
1516 {
1517 *fLog << warn << GetDescriptor() << ": Cannot use a Gauss fit to the F-Factor "
1518 << "in the camera with area index: " << i << endl;
1519 *fLog << warn << GetDescriptor() << ": Fit probability " << prob
1520 << " is smaller than 0.001 " << endl;
1521 *fLog << warn << GetDescriptor() << ": Will use the simple mean and rms " << endl;
1522 delete hist;
1523 continue;
1524 }
1525
1526 *fLog << inf << endl;
1527 *fLog << inf << GetDescriptor() << ": Mean F-Factor "
1528 << "with area index #" << i << ": "
1529 << Form("%4.2f+-%4.2f",mean,sigma) << "fadc/ph" << endl;
1530
1531 lowlim [i] = 1.1;
1532 upplim [i] = mean + fFFactorErrLimit*sigma;
1533
1534 delete hist;
1535 }
1536
1537 *fLog << endl;
1538
1539 for (UInt_t i=0; i<npixels; i++)
1540 {
1541
1542 MCalibrationChargePix &pix = (MCalibrationChargePix&)(*fCam)[i];
1543 MBadPixelsPix &bad = (*fBadPixels)[i];
1544
1545 if (bad.IsUnsuitable(MBadPixelsPix::kUnsuitableRun))
1546 continue;
1547
1548 const Float_t ffactor = pix.GetMeanFFactorFADC2Phot();
1549 const Int_t aidx = (*fGeom)[i].GetAidx();
1550
1551 if ( ffactor < lowlim[aidx] || ffactor > upplim[aidx] )
1552 {
1553 *fLog << warn << GetDescriptor() << ": Overall F-Factor "
1554 << Form("%5.2f",ffactor) << " out of range ["
1555 << Form("%5.2f,%5.2f",lowlim[aidx],upplim[aidx]) << "] pixel " << i << endl;
1556
1557 bad.SetUncalibrated( MBadPixelsPix::kDeviatingFFactor );
1558 bad.SetUnsuitable ( MBadPixelsPix::kUnsuitableRun );
1559 }
1560 }
1561
1562 for (UInt_t i=0; i<npixels; i++)
1563 {
1564
1565 MCalibrationChargePix &pix = (MCalibrationChargePix&)(*fCam)[i];
1566 MCalibrationQEPix &qepix = (MCalibrationQEPix&) (*fQECam)[i];
1567 MBadPixelsPix &bad = (*fBadPixels)[i];
1568
1569 if (bad.IsUnsuitable(MBadPixelsPix::kUnsuitableRun))
1570 {
1571 qepix.SetFFactorMethodValid(kFALSE,fPulserColor);
1572 pix.SetFFactorMethodValid(kFALSE);
1573 pix.SetExcluded();
1574 continue;
1575 }
1576 }
1577}
1578
1579
1580// ------------------------------------------------------------------------
1581//
1582// Loop over pixels:
1583//
1584// - Continue, if not MCalibrationChargeBlindPix::IsFluxInsidePlexiglassAvailable() and set:
1585// MCalibrationQEPix::SetBlindPixelMethodValid(kFALSE,fPulserColor)
1586//
1587// - Calculate the quantum efficiency with the formula:
1588//
1589// QE = Num.Phes / MCalibrationChargeBlindPix::GetFluxInsidePlexiglass()
1590// / MGeomPix::GetA() * MCalibrationQECam::GetPlexiglassQE()
1591//
1592// - Set QE in MCalibrationQEPix::SetQEBlindPixel ( QE, fPulserColor );
1593// - Set Variance of QE in MCalibrationQEPix::SetQEBlindPixelVar ( Variance, fPulserColor );
1594// - Set bit MCalibrationQEPix::SetBlindPixelMethodValid(kTRUE,fPulserColor)
1595//
1596// - Call MCalibrationQEPix::UpdateBlindPixelMethod()
1597//
1598void MCalibrationChargeCalc::FinalizeBlindPixelQECam()
1599{
1600
1601 const UInt_t npixels = fGeom->GetNumPixels();
1602
1603 //
1604 // Set the results in the MCalibrationChargeCam
1605 //
1606 if (fBlindPixel)
1607 {
1608 if (fBlindPixel->IsFluxInsidePlexiglassAvailable())
1609 {
1610
1611 const Float_t photons = fBlindPixel->GetFluxInsidePlexiglass() * (*fGeom)[0].GetA()
1612 / fQECam->GetPlexiglassQE();
1613 fCam->SetNumPhotonsBlindPixelMethod(photons);
1614
1615 const Float_t photrelvar = fBlindPixel->GetFluxInsidePlexiglassRelVar()
1616 + fQECam->GetPlexiglassQERelVar();
1617 if (photrelvar > 0.)
1618 fCam->SetNumPhotonsBlindPixelMethodErr(TMath::Sqrt( photrelvar * photons * photons));
1619 }
1620 }
1621 //
1622 // With the knowledge of the overall photon flux, calculate the
1623 // quantum efficiencies after the Blind Pixel and PIN Diode method
1624 //
1625 for (UInt_t i=0; i<npixels; i++)
1626 {
1627
1628 MCalibrationQEPix &qepix = (MCalibrationQEPix&) (*fQECam)[i];
1629
1630 if (!fBlindPixel)
1631 {
1632 qepix.SetBlindPixelMethodValid(kFALSE, fPulserColor);
1633 continue;
1634 }
1635
1636 if (!fBlindPixel->IsFluxInsidePlexiglassAvailable())
1637 {
1638 qepix.SetBlindPixelMethodValid(kFALSE, fPulserColor);
1639 continue;
1640 }
1641
1642 MBadPixelsPix &bad = (*fBadPixels)[i];
1643 if (bad.IsUnsuitable (MBadPixelsPix::kUnsuitableRun))
1644 {
1645 qepix.SetBlindPixelMethodValid(kFALSE, fPulserColor);
1646 continue;
1647 }
1648
1649 MCalibrationChargePix &pix = (MCalibrationChargePix&)(*fCam)[i];
1650 MGeomPix &geo = (*fGeom)[i];
1651
1652 const Float_t qe = pix.GetPheFFactorMethod()
1653 / fBlindPixel->GetFluxInsidePlexiglass()
1654 / geo.GetA()
1655 * fQECam->GetPlexiglassQE();
1656
1657 const Float_t qerelvar = fBlindPixel->GetFluxInsidePlexiglassRelVar()
1658 + fQECam->GetPlexiglassQERelVar()
1659 + pix.GetPheFFactorMethodRelVar();
1660
1661 qepix.SetQEBlindPixel ( qe , fPulserColor );
1662 qepix.SetQEBlindPixelVar ( qerelvar*qe*qe, fPulserColor );
1663 qepix.SetBlindPixelMethodValid( kTRUE , fPulserColor );
1664
1665 if (!qepix.UpdateBlindPixelMethod())
1666 *fLog << warn << GetDescriptor()
1667 << ": Cannot update Quantum efficiencies with the Blind Pixel Method" << endl;
1668 }
1669}
1670
1671// ------------------------------------------------------------------------
1672//
1673// Loop over pixels:
1674//
1675// - Continue, if not MCalibrationChargePINDiode::IsFluxOutsidePlexiglassAvailable() and set:
1676// MCalibrationQEPix::SetPINDiodeMethodValid(kFALSE,fPulserColor)
1677//
1678// - Calculate the quantum efficiency with the formula:
1679//
1680// QE = Num.Phes / MCalibrationChargePINDiode::GetFluxOutsidePlexiglass() / MGeomPix::GetA()
1681//
1682// - Set QE in MCalibrationQEPix::SetQEPINDiode ( QE, fPulserColor );
1683// - Set Variance of QE in MCalibrationQEPix::SetQEPINDiodeVar ( Variance, fPulserColor );
1684// - Set bit MCalibrationQEPix::SetPINDiodeMethodValid(kTRUE,fPulserColor)
1685//
1686// - Call MCalibrationQEPix::UpdatePINDiodeMethod()
1687//
1688void MCalibrationChargeCalc::FinalizePINDiodeQECam()
1689{
1690
1691 const UInt_t npixels = fGeom->GetNumPixels();
1692
1693 //
1694 // With the knowledge of the overall photon flux, calculate the
1695 // quantum efficiencies after the PIN Diode method
1696 //
1697 for (UInt_t i=0; i<npixels; i++)
1698 {
1699
1700 MCalibrationQEPix &qepix = (MCalibrationQEPix&) (*fQECam)[i];
1701
1702 if (!fPINDiode)
1703 {
1704 qepix.SetPINDiodeMethodValid(kFALSE, fPulserColor);
1705 continue;
1706 }
1707
1708 if (!fPINDiode->IsFluxOutsidePlexiglassAvailable())
1709 {
1710 qepix.SetPINDiodeMethodValid(kFALSE, fPulserColor);
1711 continue;
1712 }
1713
1714 MBadPixelsPix &bad = (*fBadPixels)[i];
1715
1716 if (!bad.IsUnsuitable (MBadPixelsPix::kUnsuitableRun))
1717 {
1718 qepix.SetPINDiodeMethodValid(kFALSE, fPulserColor);
1719 continue;
1720 }
1721
1722 MCalibrationChargePix &pix = (MCalibrationChargePix&)(*fCam)[i];
1723 MGeomPix &geo = (*fGeom)[i];
1724
1725 const Float_t qe = pix.GetPheFFactorMethod()
1726 / fPINDiode->GetFluxOutsidePlexiglass()
1727 / geo.GetA();
1728
1729 const Float_t qerelvar = fPINDiode->GetFluxOutsidePlexiglassRelVar() + pix.GetPheFFactorMethodRelVar();
1730
1731 qepix.SetQEPINDiode ( qe , fPulserColor );
1732 qepix.SetQEPINDiodeVar ( qerelvar*qe*qe, fPulserColor );
1733 qepix.SetPINDiodeMethodValid( kTRUE , fPulserColor );
1734
1735 if (!qepix.UpdatePINDiodeMethod())
1736 *fLog << warn << GetDescriptor()
1737 << ": Cannot update Quantum efficiencies with the PIN Diode Method" << endl;
1738 }
1739}
1740
1741// -----------------------------------------------------------------------------------------------
1742//
1743// - Print out statistics about BadPixels of type UnsuitableType_t
1744// - store numbers of bad pixels of each type in fCam
1745//
1746void MCalibrationChargeCalc::FinalizeUnsuitablePixels()
1747{
1748
1749 *fLog << inf << endl;
1750 *fLog << GetDescriptor() << ": Charge Calibration status:" << endl;
1751 *fLog << dec << setfill(' ');
1752
1753 const Int_t nareas = fGeom->GetNumAreas();
1754
1755 Int_t counts[nareas];
1756 memset(counts,0,nareas*sizeof(Int_t));
1757
1758 for (Int_t i=0; i<fBadPixels->GetSize(); i++)
1759 {
1760 MBadPixelsPix &bad = (*fBadPixels)[i];
1761 if (!bad.IsBad())
1762 {
1763 const Int_t aidx = (*fGeom)[i].GetAidx();
1764 counts[aidx]++;
1765 }
1766 }
1767
1768 if (fGeom->InheritsFrom("MGeomCamMagic"))
1769 *fLog << " " << setw(7) << "Successfully calibrated Pixels: "
1770 << Form("%s%3i%s%3i","Inner: ",counts[0]," Outer: ",counts[1]) << endl;
1771
1772 memset(counts,0,nareas*sizeof(Int_t));
1773
1774 for (Int_t i=0; i<fBadPixels->GetSize(); i++)
1775 {
1776 MBadPixelsPix &bad = (*fBadPixels)[i];
1777 if (bad.IsUnsuitable(MBadPixelsPix::kUnsuitableRun))
1778 {
1779 const Int_t aidx = (*fGeom)[i].GetAidx();
1780 counts[aidx]++;
1781 }
1782 }
1783
1784 for (Int_t aidx=0; aidx<nareas; aidx++)
1785 fCam->SetNumUnsuitable(counts[aidx], aidx);
1786
1787 if (fGeom->InheritsFrom("MGeomCamMagic"))
1788 *fLog << " " << setw(7) << "Uncalibrated Pixels: "
1789 << Form("%s%3i%s%3i","Inner: ",counts[0]," Outer: ",counts[1]) << endl;
1790
1791 memset(counts,0,nareas*sizeof(Int_t));
1792
1793 for (Int_t i=0; i<fBadPixels->GetSize(); i++)
1794 {
1795 MBadPixelsPix &bad = (*fBadPixels)[i];
1796 if (bad.IsUnsuitable(MBadPixelsPix::kUnreliableRun))
1797 {
1798 const Int_t aidx = (*fGeom)[i].GetAidx();
1799 counts[aidx]++;
1800 }
1801 }
1802
1803 for (Int_t aidx=0; aidx<nareas; aidx++)
1804 fCam->SetNumUnreliable(counts[aidx], aidx);
1805
1806 *fLog << " " << setw(7) << "Unreliable Pixels: "
1807 << Form("%s%3i%s%3i","Inner: ",counts[0]," Outer: ",counts[1]) << endl;
1808
1809}
1810
1811// -----------------------------------------------------------------------------------------------
1812//
1813// Print out statistics about BadPixels of type UncalibratedType_t
1814//
1815void MCalibrationChargeCalc::PrintUncalibrated(MBadPixelsPix::UncalibratedType_t typ, const char *text) const
1816{
1817
1818 UInt_t countinner = 0;
1819 UInt_t countouter = 0;
1820
1821 for (Int_t i=0; i<fBadPixels->GetSize(); i++)
1822 {
1823 MBadPixelsPix &bad = (*fBadPixels)[i];
1824 if (bad.IsUncalibrated(typ))
1825 {
1826 if (fGeom->GetPixRatio(i) == 1.)
1827 countinner++;
1828 else
1829 countouter++;
1830 }
1831 }
1832
1833 *fLog << " " << setw(7) << text
1834 << Form("%s%3i%s%3i","Inner: ",countinner," Outer: ",countouter) << endl;
1835}
1836
1837// --------------------------------------------------------------------------
1838//
1839// Set the path for output file
1840//
1841void MCalibrationChargeCalc::SetOutputPath(TString path)
1842{
1843 fOutputPath = path;
1844 if (fOutputPath.EndsWith("/"))
1845 fOutputPath = fOutputPath(0, fOutputPath.Length()-1);
1846}
1847
1848void MCalibrationChargeCalc::SetOutputFile(TString file)
1849{
1850 fOutputFile = file;
1851}
1852
1853// --------------------------------------------------------------------------
1854//
1855// Get the output file
1856//
1857const char* MCalibrationChargeCalc::GetOutputFile()
1858{
1859 return Form("%s/%s", (const char*)fOutputPath, (const char*)fOutputFile);
1860}
1861
1862Int_t MCalibrationChargeCalc::ReadEnv(const TEnv &env, TString prefix, Bool_t print)
1863{
1864 Bool_t rc = kFALSE;
1865 if (IsEnvDefined(env, prefix, "ChargeLimit", print))
1866 {
1867 SetChargeLimit(GetEnvValue(env, prefix, "ChargeLimit", fChargeLimit));
1868 rc = kTRUE;
1869 }
1870
1871 if (IsEnvDefined(env, prefix, "ChargeErrLimit", print))
1872 {
1873 SetChargeErrLimit(GetEnvValue(env, prefix, "ChargeErrLimit", fChargeErrLimit));
1874 rc = kTRUE;
1875 }
1876
1877 if (IsEnvDefined(env, prefix, "ChargeRelErrLimit", print))
1878 {
1879 SetChargeRelErrLimit(GetEnvValue(env, prefix, "ChargeRelErrLimit", fChargeRelErrLimit));
1880 rc = kTRUE;
1881 }
1882 if (IsEnvDefined(env, prefix, "Debug", print))
1883 {
1884 SetDebug(GetEnvValue(env, prefix, "Debug", IsDebug()));
1885 rc = kTRUE;
1886 }
1887 if (IsEnvDefined(env, prefix, "FFactorErrLimit", print))
1888 {
1889 SetFFactorErrLimit(GetEnvValue(env, prefix, "FFactorErrLimit", fFFactorErrLimit));
1890 rc = kTRUE;
1891 }
1892 if (IsEnvDefined(env, prefix, "LambdaErrLimit", print))
1893 {
1894 SetLambdaErrLimit(GetEnvValue(env, prefix, "LambdaErrLimit", fLambdaErrLimit));
1895 rc = kTRUE;
1896 }
1897 if (IsEnvDefined(env, prefix, "LambdaCheckLimit", print))
1898 {
1899 SetLambdaCheckLimit(GetEnvValue(env, prefix, "LambdaCheckLimit", fLambdaCheckLimit));
1900 rc = kTRUE;
1901 }
1902 if (IsEnvDefined(env, prefix, "PheErrLimit", print))
1903 {
1904 SetPheErrLimit(GetEnvValue(env, prefix, "PheErrLimit", fPheErrLimit));
1905 rc = kTRUE;
1906 }
1907 // void SetPulserColor(const MCalibrationCam::PulserColor_t col) { fPulserColor = col; }
1908
1909 return rc;
1910}
Note: See TracBrowser for help on using the repository browser.