Context Navigation

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

Visit:

Last change on this file since 3936 was 3936, checked in by gaug, 21 years ago
* empty log message *
File size: 52.5 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
189	#include "MLog.h"
190	#include "MLogManip.h"
191
192	#include "MParList.h"
193
194	#include "MRawRunHeader.h"
195	#include "MRawEvtPixelIter.h"
196
197	#include "MGeomCam.h"
198	#include "MGeomPix.h"
199
200	#include "MPedestalCam.h"
201	#include "MPedestalPix.h"
202
203	#include "MCalibrationChargeCam.h"
204	#include "MCalibrationChargePix.h"
205	#include "MCalibrationChargePINDiode.h"
206	#include "MCalibrationChargeBlindPix.h"
207
208	#include "MExtractedSignalCam.h"
209	#include "MExtractedSignalPix.h"
210	#include "MExtractedSignalBlindPixel.h"
211	#include "MExtractedSignalPINDiode.h"
212
213	#include "MBadPixelsCam.h"
214	#include "MBadPixelsPix.h"
215
216	#include "MCalibrationQECam.h"
217	#include "MCalibrationQEPix.h"
218
219	#include "MCalibrationCam.h"
220
221	ClassImp(MCalibrationChargeCalc);
222
223	using namespace std;
224
225	const Float_t MCalibrationChargeCalc::fgChargeLimit = 3.;
226	const Float_t MCalibrationChargeCalc::fgChargeErrLimit = 0.;
227	const Float_t MCalibrationChargeCalc::fgChargeRelErrLimit = 1.;
228	const Float_t MCalibrationChargeCalc::fgLambdaErrLimit = 0.2;
229	const Float_t MCalibrationChargeCalc::fgLambdaCheckLimit = 0.2;
230	const Float_t MCalibrationChargeCalc::fgPheErrLimit = 4.;
231	const Float_t MCalibrationChargeCalc::fgFFactorErrLimit = 3.;
232	// --------------------------------------------------------------------------
233	//
234	// Default constructor.
235	//
236	// Sets all pointers to NULL
237	//
238	// Calls AddToBranchList for:
239	// - MRawEvtData.fHiGainPixId
240	// - MRawEvtData.fLoGainPixId
241	// - MRawEvtData.fHiGainFadcSamples
242	// - MRawEvtData.fLoGainFadcSamples
243	//
244	// Initializes:
245	// - fChargeLimit to fgChargeLimit
246	// - fChargeErrLimit to fgChargeErrLimit
247	// - fChargeRelErrLimit to fgChargeRelErrLimit
248	// - fFFactorErrLimit to fgFFactorErrLimit
249	// - fLambdaCheckLimit to fgLambdaCheckLimit
250	// - fLambdaErrLimit to fgLambdaErrLimit
251	// - fPheErrLimit to fgPheErrLimit
252	// - fPulserColor to MCalibrationCam::kCT1
253	//
254	// Calls:
255	// - Clear()
256	//
257	MCalibrationChargeCalc::MCalibrationChargeCalc(const char name, const char title)
258	: fBadPixels(NULL), fCam(NULL), fBlindPixel(NULL), fPINDiode(NULL),
259	fQECam(NULL), fGeom(NULL), fPedestals(NULL), fEvtTime(NULL)
260	{
261
262	fName = name ? name : "MCalibrationChargeCalc";
263	fTitle = title ? title : "Task to calculate the calibration constants and MCalibrationCam ";
264
265	AddToBranchList("MRawEvtData.fHiGainPixId");
266	AddToBranchList("MRawEvtData.fLoGainPixId");
267	AddToBranchList("MRawEvtData.fHiGainFadcSamples");
268	AddToBranchList("MRawEvtData.fLoGainFadcSamples");
269
270	SetChargeLimit();
271	SetChargeErrLimit();
272	SetChargeRelErrLimit();
273	SetFFactorErrLimit();
274	SetLambdaCheckLimit();
275	SetLambdaErrLimit();
276	SetPheErrLimit();
277	SetPulserColor(MCalibrationCam::kNONE);
278
279	Clear();
280
281	}
282
283	// --------------------------------------------------------------------------
284	//
285	// Sets:
286	// - all variables to 0.,
287	// - all flags to kFALSE
288	//
289	void MCalibrationChargeCalc::Clear(const Option_t *o)
290	{
291
292	fNumHiGainSamples = 0.;
293	fNumLoGainSamples = 0.;
294	fSqrtHiGainSamples = 0.;
295	fSqrtLoGainSamples = 0.;
296	SkipHiLoGainCalibration( kFALSE );
297	}
298
299
300	// -----------------------------------------------------------------------------------
301	//
302	// The following container are searched for and execution aborted if not in MParList:
303	// - MPedestalCam
304	//
305	// The following containers are searched and created if they were not found:
306	//
307	// - MCalibrationQECam
308	// - MBadPixelsCam
309	//
310	// The following output containers are only searched, but not created. If they
311	// cannot be found, the corresponding calibration part is only skipped.
312	//
313	// - MTime
314	//
315	Int_t MCalibrationChargeCalc::PreProcess(MParList *pList)
316	{
317
318	//
319	// Containers that have to be there.
320	//
321	fPedestals = (MPedestalCam*)pList->FindObject("MPedestalCam");
322	if (!fPedestals)
323	{
324	*fLog << err << "MPedestalCam not found... aborting" << endl;
325	return kFALSE;
326	}
327
328	//
329	// Containers that are created in case that they are not there.
330	//
331	fQECam = (MCalibrationQECam*)pList->FindCreateObj("MCalibrationQECam");
332	if (!fQECam)
333	{
334	*fLog << err << "Cannot find nor create MCalibrationQECam... aborting" << endl;
335	return kFALSE;
336	}
337
338	fBadPixels = (MBadPixelsCam*)pList->FindCreateObj("MBadPixelsCam");
339	if (!fBadPixels)
340	{
341	*fLog << err << "Could not find or create MBadPixelsCam ... aborting." << endl;
342	return kFALSE;
343	}
344
345
346	fEvtTime = (MTime*)pList->FindObject("MTime");
347
348	//
349	// Check the pulser colour --> FIXME: this solution is only valid until the arrival of the DM's
350	//
351	if (fPulserColor == MCalibrationCam::kNONE)
352	{
353	*fLog << endl;
354	*fLog << err << GetDescriptor()
355	<< ": No Pulser colour has been chosen. Since the installation of the IFAE pulser box,"
356	<< " you HAVE to provide the LEDs colour, otherwise there is no calibration. " << endl;
357	*fLog << "See e.g. the macro calibration.C " << endl;
358	return kFALSE;
359	}
360
361	return kTRUE;
362	}
363
364
365	// --------------------------------------------------------------------------
366	//
367	// Search for the following input containers and abort if not existing:
368	// - MGeomCam
369	// - MCalibrationChargeCam
370	//
371	// Search for the following input containers and give a warning if not existing:
372	// - MCalibrationChargeBlindPix
373	// - MCalibrationChargePINDiode
374	//
375	// It retrieves the following variables from MCalibrationChargeCam:
376	//
377	// - fNumHiGainSamples
378	// - fNumLoGainSamples
379	//
380	// It defines the PixId of every pixel in:
381	//
382	// - MCalibrationChargeCam
383	// - MCalibrationQECam
384	//
385	// It sets all pixels in excluded which have the flag fBadBixelsPix::IsBad() set in:
386	//
387	// - MCalibrationChargePix
388	// - MCalibrationQEPix
389	//
390	// Sets the pulser colour and tests if it has not changed w.r.t. fPulserColor in:
391	//
392	// - MCalibrationChargeCam
393	// - MCalibrationChargeBlindPix (if existing)
394	// - MCalibrationChargePINDiode (if existing)
395	//
396	Bool_t MCalibrationChargeCalc::ReInit(MParList *pList )
397	{
398
399	fGeom = (MGeomCam*)pList->FindObject("MGeomCam");
400	if (!fGeom)
401	{
402	*fLog << err << "No MGeomCam found... aborting." << endl;
403	return kFALSE;
404	}
405
406	fCam = (MCalibrationChargeCam*)pList->FindObject("MCalibrationChargeCam");
407	if (!fCam)
408	{
409	*fLog << err << "Cannot find MCalibrationChargeCam... aborting" << endl;
410	*fLog << err << "Maybe you forget to call an MFillH for the MHCalibrationChargeCam before..." << endl;
411	return kFALSE;
412	}
413
414	//
415	// Optional Containers
416	//
417	fBlindPixel = (MCalibrationChargeBlindPix*)pList->FindObject("MCalibrationChargeBlindPix");
418	if (!fBlindPixel)
419	*fLog << warn << GetDescriptor()
420	<< ": MCalibrationChargeBlindPix not found... no blind pixel method! " << endl;
421
422	fPINDiode = (MCalibrationChargePINDiode*)pList->FindObject("MCalibrationChargePINDiode");
423	if (!fPINDiode)
424	*fLog << warn << GetDescriptor()
425	<< "MCalibrationChargePINDiode not found... no PIN Diode method! " << endl;
426
427
428	//
429	// Initialize the pulser colours
430	//
431	if (fCam->GetPulserColor() == MCalibrationCam::kNONE)
432	{
433	fCam->SetPulserColor( fPulserColor );
434
435	if (fBlindPixel)
436	fBlindPixel->SetColor( fPulserColor );
437
438	if (fPINDiode)
439	fPINDiode->SetColor( fPulserColor );
440	}
441
442	if (fPulserColor != fCam->GetPulserColor())
443	{
444	*fLog << err << GetDescriptor()
445	<< ": Pulser colour has changed w.r.t. last file in MCalibrationChargeCam" << endl;
446	*fLog << err << "This feature is not yet implemented, sorry ... aborting " << endl;
447	return kFALSE;
448	}
449
450	if (fBlindPixel)
451	if (fPulserColor != fBlindPixel->GetColor())
452	{
453	*fLog << err << GetDescriptor()
454	<< ": Pulser colour has changed w.r.t. last file in MCalibrationChargeBlindPix." << endl;
455	*fLog << err << "This feature is not yet implemented, sorry ... aborting " << endl;
456	return kFALSE;
457	}
458
459	if (fPINDiode)
460	if (fPulserColor != fPINDiode->GetColor())
461	{
462	*fLog << err << GetDescriptor()
463	<< ": Pulser colour has changed w.r.t. last file in MCalibrationChargePINDiode." << endl;
464	*fLog << err << "This feature is not yet implemented, sorry ... aborting " << endl;
465	return kFALSE;
466	}
467
468
469	fNumHiGainSamples = fCam->GetNumHiGainFADCSlices();
470	fNumLoGainSamples = fCam->GetNumLoGainFADCSlices();
471	fSqrtHiGainSamples = TMath::Sqrt(fNumHiGainSamples);
472	fSqrtLoGainSamples = TMath::Sqrt(fNumLoGainSamples);
473
474	UInt_t npixels = fGeom->GetNumPixels();
475
476	for (UInt_t i=0; i<npixels; i++)
477	{
478
479	MCalibrationChargePix &pix = (MCalibrationChargePix&)(*fCam) [i];
480	MCalibrationQEPix &pqe = (MCalibrationQEPix&) (*fQECam)[i];
481	MBadPixelsPix &bad = (*fBadPixels)[i];
482
483	pix.SetPixId(i);
484	pqe.SetPixId(i);
485
486	if (bad.IsBad())
487	{
488	pix.SetExcluded();
489	pqe.SetExcluded();
490	continue;
491	}
492
493	}
494
495	return kTRUE;
496	}
497
498	// ----------------------------------------------------------------------------------
499	//
500	// Nothing to be done in Process, but have a look at MHCalibrationChargeCam, instead
501	//
502	Int_t MCalibrationChargeCalc::Process()
503	{
504	return kTRUE;
505	}
506
507	// -----------------------------------------------------------------------
508	//
509	// Return if number of executions is null.
510	//
511	// First loop over pixels, average areas and sectors, call:
512	// - FinalizePedestals()
513	// - FinalizeCharges()
514	// for every entry. Count number of valid pixels in loop and return kFALSE
515	// if there are none (the "Michele check").
516	//
517	// Call FinalizeBadPixels()
518	//
519	// Call FinalizeFFactorMethod() (second and third loop over pixels and areas)
520	//
521	// Call FinalizeBlindPixel()
522	// Call FinalizePINDiode()
523	//
524	// Call FinalizeFFactorQECam() (fourth loop over pixels and areas)
525	// Call FinalizeBlindPixelQECam() (fifth loop over pixels and areas)
526	// Call FinalizePINDiodeQECam() (sixth loop over pixels and areas)
527	//
528	// Call MParContainer::SetReadyToSave() for fCam, fQECam, fBadPixels and
529	// fBlindPixel and fPINDiode if they exist
530	//
531	// Print out some statistics
532	//
533	Int_t MCalibrationChargeCalc::PostProcess()
534	{
535
536	if (GetNumExecutions()==0)
537	return kFALSE;
538
539	if (fPINDiode)
540	if (!fPINDiode->IsValid())
541	{
542	*fLog << warn << GetDescriptor()
543	<< ": MCalibrationChargePINDiode is declared not valid... no PIN Diode method! " << endl;
544	fPINDiode = NULL;
545	}
546
547	if (fBlindPixel)
548	if (!fBlindPixel->IsValid())
549	{
550	*fLog << warn << GetDescriptor()
551	<< ": MCalibrationChargeBlindPix is declared not valid... no Blind Pixel method! " << endl;
552	fBlindPixel = NULL;
553	}
554
555	//
556	// First loop over pixels, call FinalizePedestals and FinalizeCharges
557	//
558	Int_t nvalid = 0;
559
560	for (Int_t pixid=0; pixid<fPedestals->GetSize(); pixid++)
561	{
562
563	MCalibrationChargePix &pix = (MCalibrationChargePix&)(*fCam)[pixid];
564	//
565	// Check if the pixel has been excluded from the fits
566	//
567	if (pix.IsExcluded())
568	continue;
569
570	MPedestalPix &ped = (*fPedestals)[pixid];
571	MBadPixelsPix &bad = (*fBadPixels)[pixid];
572
573	FinalizePedestals(ped,pix);
574
575	if (FinalizeCharges(pix,bad))
576	nvalid++;
577	}
578
579	//
580	// The Michele check ...
581	//
582	if (nvalid == 0)
583	{
584	*fLog << err << GetDescriptor() << ": All pixels have non-valid calibration. "
585	<< "Did you forget to fill the histograms "
586	<< "(filling MHCalibrationChargeCam from MExtractedSignalCam using MFillH) ? " << endl;
587	*fLog << err << GetDescriptor() << ": Or, maybe, you have used a pedestal run "
588	<< "instead of a calibration run " << endl;
589	return kFALSE;
590	}
591
592	for (UInt_t aidx=0; aidx<fGeom->GetNumAreas(); aidx++)
593	{
594
595	const MPedestalPix &ped = fPedestals->GetAverageArea(aidx);
596	MCalibrationChargePix &pix = (MCalibrationChargePix&)fCam->GetAverageArea(aidx);
597
598	FinalizePedestals(ped,pix);
599	FinalizeCharges(pix, fCam->GetAverageBadArea(aidx));
600	}
601
602	for (UInt_t sector=0; sector<fGeom->GetNumSectors(); sector++)
603	{
604
605	const MPedestalPix &ped = fPedestals->GetAverageSector(sector);
606	MCalibrationChargePix &pix = (MCalibrationChargePix&)fCam->GetAverageSector(sector);
607
608	FinalizePedestals(ped,pix);
609	FinalizeCharges(pix, fCam->GetAverageBadSector(sector));
610	}
611
612	//
613	// Finalize Bad Pixels
614	//
615	FinalizeBadPixels();
616
617	//
618	// Finalize F-Factor method
619	//
620	if (!FinalizeFFactorMethod())
621	{
622	*fLog << warn << "Could not calculate the photons flux from the F-Factor method " << endl;
623	fCam->SetFFactorMethodValid(kFALSE);
624	return kFALSE;
625	}
626	else
627	fCam->SetFFactorMethodValid(kTRUE);
628
629	//
630	// Finalize Blind Pixel
631	//
632	if (FinalizeBlindPixel())
633	fQECam->SetBlindPixelMethodValid(kTRUE);
634	else
635	fQECam->SetBlindPixelMethodValid(kFALSE);
636
637	//
638	// Finalize PIN Diode
639	//
640	if (FinalizePINDiode())
641	fQECam->SetPINDiodeMethodValid(kTRUE);
642	else
643	fQECam->SetPINDiodeMethodValid(kFALSE);
644
645	//
646	// Finalize QE Cam
647	//
648	FinalizeFFactorQECam();
649	FinalizeBlindPixelQECam();
650	FinalizePINDiodeQECam();
651
652	//
653	// Finalize calibration statistics
654	//
655	FinalizeUnsuitablePixels();
656
657	fCam ->SetReadyToSave();
658	fQECam ->SetReadyToSave();
659	fBadPixels->SetReadyToSave();
660
661	if (fBlindPixel)
662	fBlindPixel->SetReadyToSave();
663	if (fPINDiode)
664	fPINDiode->SetReadyToSave();
665
666	*fLog << inf << endl;
667	*fLog << GetDescriptor() << ": Errors statistics:" << endl;
668
669	PrintUncalibrated(MBadPixelsPix::kChargeIsPedestal,
670	Form("%s%2.1f%s","Signal less than ",fChargeLimit," Pedestal RMS: "));
671	PrintUncalibrated(MBadPixelsPix::kChargeErrNotValid,
672	Form("%s%2.1f%s","Signal Error smaller than ",fChargeErrLimit,": "));
673	PrintUncalibrated(MBadPixelsPix::kChargeRelErrNotValid,
674	Form("%s%2.1f%s","Signal Error bigger than ",fChargeRelErrLimit," times Mean Signal: "));
675	PrintUncalibrated(MBadPixelsPix::kChargeSigmaNotValid,
676	"Signal Sigma smaller than Pedestal RMS: ");
677	PrintUncalibrated(MBadPixelsPix::kLoGainSaturation,
678	"Pixels with Low Gain Saturation: ");
679	PrintUncalibrated(MBadPixelsPix::kMeanTimeInFirstBin,
680	Form("%s%2.1f%s","Mean Abs. Arr. Time in First ",1.," Bin(s): "));
681	PrintUncalibrated(MBadPixelsPix::kMeanTimeInLast2Bins,
682	Form("%s%2.1f%s","Mean Abs. Arr. Time in Last ",2.," Bin(s): "));
683	PrintUncalibrated(MBadPixelsPix::kHiGainOscillating,
684	"Pixels with changing Hi Gain signal over time: ");
685	PrintUncalibrated(MBadPixelsPix::kLoGainOscillating,
686	"Pixels with changing Lo Gain signal over time: ");
687	PrintUncalibrated(MBadPixelsPix::kDeviatingNumPhes,
688	"Pixels with deviating number of phes: ");
689	PrintUncalibrated(MBadPixelsPix::kHiGainNotFitted,
690	"Pixels with unsuccesful Gauss fit to the Hi Gain: ");
691	PrintUncalibrated(MBadPixelsPix::kLoGainNotFitted,
692	"Pixels with unsuccesful Gauss fit to the Lo Gain: ");
693
694	return kTRUE;
695	}
696
697	// ----------------------------------------------------------------------------------
698	//
699	// Retrieves pedestal and pedestal RMS from MPedestalPix
700	// Retrieves total entries from MPedestalCam
701	// Sets pedestalfNumHiGainSamples and pedestalfNumLoGainSamples in MCalibrationChargePix
702	// Sets pedRMS fSqrtHiGainSamples and pedRMS fSqrtLoGainSamples in MCalibrationChargePix
703	//
704	// If the flag MCalibrationPix::IsHiGainSaturation() is set, call also:
705	// - MCalibrationChargePix::CalcLoGainPedestal()
706	//
707	void MCalibrationChargeCalc::FinalizePedestals(const MPedestalPix &ped, MCalibrationChargePix &cal)
708	{
709
710	//
711	// get the pedestals
712	//
713	const Float_t pedes = ped.GetPedestal();
714	const Float_t prms = ped.GetPedestalRms();
715	const Float_t num = TMath::Sqrt((Float_t)fPedestals->GetTotalEntries());
716
717	//
718	// set them in the calibration camera
719	//
720	if (cal.IsHiGainSaturation())
721	{
722	cal.SetPedestal(pedes* fNumLoGainSamples,
723	prms * fSqrtLoGainSamples,
724	prms * fNumLoGainSamples / num);
725	cal.CalcLoGainPedestal((Float_t)fNumLoGainSamples);
726	}
727	else
728	{
729	cal.SetPedestal(pedes* fNumHiGainSamples,
730	prms * fSqrtHiGainSamples,
731	prms * fNumHiGainSamples / num);
732	}
733
734	}
735
736	// ----------------------------------------------------------------------------------------------------
737	//
738	// Check fit results validity. Bad Pixels flags are set if:
739	//
740	// 1) Pixel has a mean smaller than fChargeLimit*PedRMS ( Flag: MBadPixelsPix::kChargeIsPedestal)
741	// 2) Pixel has a mean error smaller than fChargeErrLimit ( Flag: MBadPixelsPix::kChargeErrNotValid)
742	// 3) Pixel has mean smaller than fChargeRelVarLimit times its mean error
743	// ( Flag: MBadPixelsPix::kChargeRelErrNotValid)
744	// 4) Pixel has a sigma bigger than its Pedestal RMS ( Flag: MBadPixelsPix::kChargeSigmaNotValid )
745	//
746	// Further returns if flags: MBadPixelsPix::kUnsuitableRun is set
747	//
748	// Calls MCalibrationChargePix::CalcReducedSigma() and sets flag: MBadPixelsPix::kChargeIsPedestal
749	// if not succesful.
750	//
751	// Calls MCalibrationChargePix::CalcFFactorMethod() and sets flag: MBadPixelsPix::kDeviatingNumPhes)
752	// if not succesful.
753	//
754	Bool_t MCalibrationChargeCalc::FinalizeCharges(MCalibrationChargePix &cal, MBadPixelsPix &bad)
755	{
756
757	if (cal.GetMean() < fChargeLimit*cal.GetPedRms())
758	{
759	*fLog << warn << GetDescriptor() << ": Fitted Charge: " << cal.GetMean() << " is smaller than "
760	<< fChargeLimit << " Pedestal RMS: " << cal.GetPedRms() << " in Pixel " << cal.GetPixId() << endl;
761	bad.SetUncalibrated( MBadPixelsPix::kChargeIsPedestal);
762	}
763
764	if (cal.GetMeanErr() < fChargeErrLimit)
765	{
766	*fLog << warn << GetDescriptor() << ": Error of Fitted Charge: " << cal.GetMeanErr()
767	<< " is smaller than " << fChargeErrLimit << " in Pixel " << cal.GetPixId() << endl;
768	bad.SetUncalibrated( MBadPixelsPix::kChargeErrNotValid );
769	}
770
771	if (cal.GetMean() < fChargeRelErrLimit*cal.GetMeanErr())
772	{
773	*fLog << warn << GetDescriptor() << ": Fitted Charge: " << cal.GetMean() << " is smaller than "
774	<< fChargeRelErrLimit << "* its error: " << cal.GetMeanErr()
775	<< " in Pixel " << cal.GetPixId() << endl;
776	bad.SetUncalibrated( MBadPixelsPix::kChargeRelErrNotValid );
777	}
778
779	if (cal.GetSigma() < cal.GetPedRms())
780	{
781	*fLog << warn << GetDescriptor() << ": Sigma of Fitted Charge: " << cal.GetSigma()
782	<< " smaller than Pedestal RMS: " << cal.GetPedRms() << " in Pixel " << cal.GetPixId() << endl;
783	bad.SetUncalibrated( MBadPixelsPix::kChargeSigmaNotValid );
784	}
785
786	if (bad.IsUnsuitable(MBadPixelsPix::kUnsuitableRun))
787	return kFALSE;
788
789	if (!cal.CalcReducedSigma())
790	{
791	*fLog << warn << GetDescriptor()
792	<< ": Could not calculate reduced sigmas of pixel: " << cal.GetPixId() << endl;
793	bad.SetUncalibrated(MBadPixelsPix::kChargeIsPedestal);
794	return kFALSE;
795	}
796
797	if (!cal.CalcFFactorMethod())
798	{
799	*fLog << warn << GetDescriptor()
800	<< ": Could not calculate F-Factor of pixel: " << cal.GetPixId() << endl;
801	bad.SetUncalibrated(MBadPixelsPix::kDeviatingNumPhes);
802	return kFALSE;
803	}
804
805	return kTRUE;
806	}
807
808
809
810	// -----------------------------------------------------------------------------------
811	//
812	// Sets pixel to MBadPixelsPix::kUnsuitableRun, if one of the following flags is set:
813	// - MBadPixelsPix::kChargeIsPedestal
814	// - MBadPixelsPix::kChargeErrNotValid
815	// - MBadPixelsPix::kChargeRelErrNotValid
816	// - MBadPixelsPix::kChargeSigmaNotValid
817	// - MBadPixelsPix::kMeanTimeInFirstBin
818	// - MBadPixelsPix::kMeanTimeInLast2Bins
819	// - MBadPixelsPix::kDeviatingNumPhes
820	//
821	// - Call MCalibrationPix::SetExcluded() for the bad pixels
822	//
823	void MCalibrationChargeCalc::FinalizeBadPixels()
824	{
825
826	for (Int_t i=0; i<fBadPixels->GetSize(); i++)
827	{
828
829	MBadPixelsPix &bad = (*fBadPixels)[i];
830	MCalibrationPix &pix = (*fCam)[i];
831
832	if (bad.IsUncalibrated( MBadPixelsPix::kChargeIsPedestal))
833	bad.SetUnsuitable( MBadPixelsPix::kUnsuitableRun );
834
835	if (bad.IsUncalibrated( MBadPixelsPix::kChargeErrNotValid ))
836	bad.SetUnsuitable( MBadPixelsPix::kUnsuitableRun );
837
838	if (bad.IsUncalibrated( MBadPixelsPix::kChargeRelErrNotValid ))
839	bad.SetUnsuitable( MBadPixelsPix::kUnsuitableRun );
840
841	if (bad.IsUncalibrated( MBadPixelsPix::kChargeSigmaNotValid ))
842	bad.SetUnsuitable( MBadPixelsPix::kUnsuitableRun );
843
844	if (bad.IsUncalibrated( MBadPixelsPix::kMeanTimeInFirstBin ))
845	bad.SetUnsuitable( MBadPixelsPix::kUnsuitableRun );
846
847	if (bad.IsUncalibrated( MBadPixelsPix::kMeanTimeInLast2Bins ))
848	bad.SetUnsuitable( MBadPixelsPix::kUnsuitableRun );
849
850	if (bad.IsUncalibrated( MBadPixelsPix::kDeviatingNumPhes ))
851	bad.SetUnsuitable( MBadPixelsPix::kUnsuitableRun );
852	// bad.SetUnsuitable( MBadPixelsPix::kUnreliableRun );
853
854	if (bad.IsUnsuitable( MBadPixelsPix::kUnsuitableRun ))
855	pix.SetExcluded();
856
857	}
858	}
859
860	// ------------------------------------------------------------------------
861	//
862	//
863	// First loop: Calculate a mean and mean RMS of photo-electrons per area index
864	// Include only pixels which are not MBadPixelsPix::kUnsuitableRun or
865	// MBadPixelsPix::kUnreliableRun (see FinalizeBadPixels()) and set
866	// MCalibrationChargePix::SetFFactorMethodValid(kFALSE) in that case.
867	//
868	// Second loop: Get weighted mean number of photo-electrons and its RMS including
869	// only pixels with flag MCalibrationChargePix::IsFFactorMethodValid()
870	// and further exclude those deviating by more than fPheErrLimit mean
871	// sigmas from the mean (obtained in first loop). Set
872	// MBadPixelsPix::kDeviatingNumPhes if excluded.
873	//
874	// Set weighted mean and variance of photo-electrons per area index in:
875	// average area pixels of MCalibrationChargeCam (obtained from:
876	// MCalibrationChargeCam::GetAverageArea() )
877	//
878	// Set weighted mean and variance of photo-electrons per sector in:
879	// average sector pixels of MCalibrationChargeCam (obtained from:
880	// MCalibrationChargeCam::GetAverageSector() )
881	//
882	Bool_t MCalibrationChargeCalc::FinalizeFFactorMethod()
883	{
884
885	const UInt_t npixels = fGeom->GetNumPixels();
886	const UInt_t nareas = fGeom->GetNumAreas();
887	const UInt_t nsectors = fGeom->GetNumSectors();
888
889	Float_t lowlim [nareas];
890	Float_t upplim [nareas];
891	Float_t areavars [nareas];
892	Float_t areaweights [nareas], sectorweights [nsectors];
893	Float_t areaphes [nareas], sectorphes [nsectors];
894	Int_t numareavalid[nareas], numsectorvalid[nsectors];
895
896	memset(lowlim ,0, nareas * sizeof(Float_t));
897	memset(upplim ,0, nareas * sizeof(Float_t));
898	memset(areaphes ,0, nareas * sizeof(Float_t));
899	memset(areavars ,0, nareas * sizeof(Float_t));
900	memset(areaweights ,0, nareas * sizeof(Float_t));
901	memset(numareavalid ,0, nareas * sizeof(Int_t ));
902	memset(sectorweights ,0, nsectors * sizeof(Float_t));
903	memset(sectorphes ,0, nsectors * sizeof(Float_t));
904	memset(numsectorvalid,0, nsectors * sizeof(Int_t ));
905
906	//
907	// First loop: Get mean number of photo-electrons and the RMS
908	// The loop is only to recognize later pixels with very deviating numbers
909	//
910	for (UInt_t i=0; i<npixels; i++)
911	{
912
913	MCalibrationChargePix &pix = (MCalibrationChargePix&)(*fCam) [i];
914	MBadPixelsPix &bad = (*fBadPixels)[i];
915
916	if (!pix.IsFFactorMethodValid())
917	continue;
918
919	if (bad.IsUnsuitable(MBadPixelsPix::kUnsuitableRun))
920	{
921	pix.SetFFactorMethodValid(kFALSE);
922	continue;
923	}
924
925	if (bad.IsUnsuitable(MBadPixelsPix::kUnreliableRun))
926	continue;
927
928	const Float_t nphe = pix.GetPheFFactorMethod();
929	const Float_t nvar = pix.GetPheFFactorMethodVar();
930	const Int_t aidx = (*fGeom)[i].GetAidx();
931
932	if (nvar > 0.)
933	{
934	areaphes [aidx] += nphe;
935	areavars [aidx] += nvar;
936	numareavalid[aidx] ++;
937	}
938	}
939
940	for (UInt_t i=0; i<nareas; i++)
941	{
942	if (numareavalid[i] == 0)
943	{
944	*fLog << warn << GetDescriptor() << ": No pixels with valid number of photo-electrons found "
945	<< "in area index: " << i << endl;
946	continue;
947	}
948
949	areaphes[i] = areaphes[i] / numareavalid[i];
950	areavars[i] = areavars[i] / numareavalid[i];
951	lowlim [i] = areaphes[i] - fPheErrLimit*TMath::Sqrt(areavars[i]);
952	upplim [i] = areaphes[i] + fPheErrLimit*TMath::Sqrt(areavars[i]);
953	}
954
955	memset(numareavalid,0,nareas*sizeof(Int_t));
956	memset(areaphes ,0,nareas*sizeof(Int_t));
957	memset(areavars ,0,nareas*sizeof(Int_t));
958
959	//
960	// Second loop: Get weighted mean number of photo-electrons and its RMS excluding
961	// pixels deviating by more than fPheErrLimit sigma.
962	// Set the conversion factor FADC counts to photo-electrons
963	//
964	for (UInt_t i=0; i<npixels; i++)
965	{
966
967	MCalibrationChargePix &pix = (MCalibrationChargePix&)(*fCam)[i];
968
969	if (!pix.IsFFactorMethodValid())
970	continue;
971
972	const Float_t nvar = pix.GetPheFFactorMethodVar();
973
974	if (nvar <= 0.)
975	{
976	pix.SetFFactorMethodValid(kFALSE);
977	continue;
978	}
979
980	MBadPixelsPix &bad = (*fBadPixels)[i];
981
982	const Int_t aidx = (*fGeom)[i].GetAidx();
983	const Int_t sector = (*fGeom)[i].GetSector();
984	const Float_t nphe = pix.GetPheFFactorMethod();
985
986	if ( nphe < lowlim[aidx] \|\| nphe > upplim[aidx] )
987	{
988	*fLog << warn << GetDescriptor() << ": Deviating number of photo-electrons: "
989	<< Form("%4.2f",nphe) << " out of accepted limits: ["
990	<< Form("%4.2f%s%4.2f",lowlim[aidx],",",upplim[aidx]) << "] in pixel " << i << endl;
991	bad.SetUncalibrated( MBadPixelsPix::kDeviatingNumPhes );
992	bad.SetUnsuitable ( MBadPixelsPix::kUnsuitableRun );
993	pix.SetExcluded();
994	// bad.SetUnsuitable ( MBadPixelsPix::kUnreliableRun );
995	continue;
996	}
997
998	const Float_t weight = 1./nvar;
999
1000	areaweights [aidx] += weight;
1001	areaphes [aidx] += weight*nphe;
1002	numareavalid [aidx] ++;
1003	sectorweights [sector] += weight;
1004	sectorphes [sector] += weight*nphe;
1005	numsectorvalid[sector] ++;
1006	}
1007
1008	for (UInt_t aidx=0; aidx<nareas; aidx++)
1009	{
1010
1011	MCalibrationChargePix &apix = (MCalibrationChargePix&)fCam->GetAverageArea(aidx);
1012
1013	if (areaweights[aidx] <= 0. \|\| areaphes[aidx] <= 0.)
1014	{
1015	*fLog << warn << " Mean number of phe's from area index " << aidx << " cannot be calculated: "
1016	<< " Sum of weights: " << areaweights[aidx]
1017	<< " Sum of weighted phes: " << areaphes[aidx] << endl;
1018	apix.SetFFactorMethodValid(kFALSE);
1019	continue;
1020	}
1021
1022	*fLog << inf << "Replacing number photo-electrons of average area idx " << aidx << ": "
1023	<< Form("%5.3f%s%5.3f",apix.GetPheFFactorMethod()," +- ",apix.GetPheFFactorMethodErr()) << endl;
1024	*fLog << inf << " by average number of photo-electrons from area idx " << aidx << ": "
1025	<< Form("%5.3f%s%5.3f",areaphes[aidx] / areaweights[aidx]," +- ",
1026	TMath::Sqrt(1./areaweights[aidx])) << endl;
1027
1028	apix.SetPheFFactorMethod ( areaphes[aidx]/ areaweights[aidx] );
1029	apix.SetPheFFactorMethodVar( 1. / areaweights[aidx] );
1030	apix.SetFFactorMethodValid ( kTRUE );
1031
1032	}
1033
1034	for (UInt_t sector=0; sector<nsectors; sector++)
1035	{
1036
1037	MCalibrationChargePix &spix = (MCalibrationChargePix&)fCam->GetAverageSector(sector);
1038
1039	if (sectorweights[sector] <= 0. \|\| sectorphes[sector] <= 0.)
1040	{
1041	*fLog << warn << " Mean number of phe's from sector " << sector << " cannot be calculated: "
1042	<< " Sum of weights: " << sectorweights[sector]
1043	<< " Sum of weighted phes: " << sectorphes[sector] << endl;
1044	spix.SetFFactorMethodValid(kFALSE);
1045	continue;
1046	}
1047
1048	*fLog << inf << "Replacing number photo-electrons of average sector " << sector << ": "
1049	<< Form("%5.3f%s%5.3f",spix.GetPheFFactorMethod()," +- ",spix.GetPheFFactorMethodErr()) << endl;
1050	*fLog << inf << " by average number photo-electrons from sector " << sector << ": "
1051	<< Form("%5.3f%s%5.3f",sectorphes[sector]/ sectorweights[sector]," +- ",
1052	TMath::Sqrt(1./sectorweights[sector])) << endl;
1053
1054	spix.SetPheFFactorMethod ( sectorphes[sector]/ sectorweights[sector] );
1055	spix.SetPheFFactorMethodVar( 1. / sectorweights[sector] );
1056	spix.SetFFactorMethodValid ( kTRUE );
1057
1058	}
1059
1060	return kTRUE;
1061	}
1062
1063
1064	// ------------------------------------------------------------------------
1065	//
1066	// Returns kFALSE if pointer to MCalibrationChargeBlindPix is NULL
1067	//
1068	// The check returns kFALSE if:
1069	//
1070	// 1) fLambda and fLambdaCheck are separated relatively to each other by more than fLambdaCheckLimit
1071	// 2) BlindPixel has an fLambdaErr greater than fLambdaErrLimit
1072	//
1073	// Calls:
1074	// - MCalibrationChargeBlindPix::CalcFluxInsidePlexiglass()
1075	//
1076	Bool_t MCalibrationChargeCalc::FinalizeBlindPixel()
1077	{
1078
1079	if (!fBlindPixel)
1080	return kFALSE;
1081
1082	const Float_t lambda = fBlindPixel->GetLambda();
1083	const Float_t lambdaerr = fBlindPixel->GetLambdaErr();
1084	const Float_t lambdacheck = fBlindPixel->GetLambdaCheck();
1085
1086	if (2.*(lambdacheck-lambda)/(lambdacheck+lambda) < fLambdaCheckLimit)
1087	{
1088	*fLog << warn << GetDescriptor() << ": Lambda and Lambda-Check differ by more than "
1089	<< fLambdaCheckLimit << " in the Blind Pixel " << endl;
1090	return kFALSE;
1091	}
1092
1093	if (lambdaerr < fLambdaErrLimit)
1094	{
1095	*fLog << warn << GetDescriptor() << ": Error of Fitted Lambda is greater than "
1096	<< fLambdaErrLimit << " in Blind Pixel " << endl;
1097	return kFALSE;
1098	}
1099
1100	if (!fBlindPixel->CalcFluxInsidePlexiglass())
1101	{
1102	*fLog << warn << "Could not calculate the flux of photons from the Blind Pixel, "
1103	<< "will skip Blind Pixel Calibration " << endl;
1104	return kFALSE;
1105	}
1106
1107	return kTRUE;
1108	}
1109
1110	// ------------------------------------------------------------------------
1111	//
1112	// Returns kFALSE if pointer to MCalibrationChargePINDiode is NULL
1113	//
1114	// The check returns kFALSE if:
1115	//
1116	// 1) PINDiode has a fitted charge smaller than fChargeLimit*PedRMS
1117	// 2) PINDiode has a fit error smaller than fChargeErrLimit
1118	// 3) PINDiode has a fitted charge smaller its fChargeRelErrLimit times its charge error
1119	// 4) PINDiode has a charge sigma smaller than its Pedestal RMS
1120	//
1121	// Calls:
1122	// - MCalibrationChargePINDiode::CalcFluxOutsidePlexiglass()
1123	//
1124	Bool_t MCalibrationChargeCalc::FinalizePINDiode()
1125	{
1126
1127	if (!fPINDiode)
1128	return kFALSE;
1129
1130	if (fPINDiode->GetMean() < fChargeLimit*fPINDiode->GetPedRms())
1131	{
1132	*fLog << warn << GetDescriptor() << ": Fitted Charge is smaller than "
1133	<< fChargeLimit << " Pedestal RMS in PINDiode " << endl;
1134	return kFALSE;
1135	}
1136
1137	if (fPINDiode->GetMeanErr() < fChargeErrLimit)
1138	{
1139	*fLog << warn << GetDescriptor() << ": Error of Fitted Charge is smaller than "
1140	<< fChargeErrLimit << " in PINDiode " << endl;
1141	return kFALSE;
1142	}
1143
1144	if (fPINDiode->GetMean() < fChargeRelErrLimit*fPINDiode->GetMeanErr())
1145	{
1146	*fLog << warn << GetDescriptor() << ": Fitted Charge is smaller than "
1147	<< fChargeRelErrLimit << "* its error in PINDiode " << endl;
1148	return kFALSE;
1149	}
1150
1151	if (fPINDiode->GetSigma() < fPINDiode->GetPedRms())
1152	{
1153	*fLog << warn << GetDescriptor()
1154	<< ": Sigma of Fitted Charge smaller than Pedestal RMS in PINDiode " << endl;
1155	return kFALSE;
1156	}
1157
1158
1159	if (!fPINDiode->CalcFluxOutsidePlexiglass())
1160	{
1161	*fLog << warn << "Could not calculate the flux of photons from the PIN Diode, "
1162	<< "will skip PIN Diode Calibration " << endl;
1163	return kFALSE;
1164	}
1165
1166	return kTRUE;
1167	}
1168
1169	// ------------------------------------------------------------------------
1170	//
1171	// Calculate the average number of photons outside the plexiglass with the
1172	// formula:
1173	//
1174	// av.Num.photons(area index) = av.Num.Phes(area index)
1175	// / MCalibrationQEPix::GetDefaultQE(fPulserColor)
1176	// / MCalibrationQECam::GetPlexiglassQE()
1177	//
1178	// Calculate the variance on the average number of photons.
1179	//
1180	// Loop over pixels:
1181	//
1182	// - Continue, if not MCalibrationChargePix::IsFFactorMethodValid() and set:
1183	// MCalibrationQEPix::SetFFactorMethodValid(kFALSE,fPulserColor)
1184	//
1185	// - Call MCalibrationChargePix::CalcMeanFFactor(av.Num.photons) and set:
1186	// MCalibrationQEPix::SetFFactorMethodValid(kFALSE,fPulserColor) if not succesful
1187	//
1188	// - Calculate the quantum efficiency with the formula:
1189	//
1190	// QE = ( Num.Phes / av.Num.photons ) * MGeomCam::GetPixRatio()
1191	//
1192	// - Set QE in MCalibrationQEPix::SetQEFFactor ( QE, fPulserColor );
1193	// - Set Variance of QE in MCalibrationQEPix::SetQEFFactorVar ( Variance, fPulserColor );
1194	// - Set bit MCalibrationQEPix::SetFFactorMethodValid(kTRUE,fPulserColor)
1195	//
1196	// - Call MCalibrationQEPix::UpdateFFactorMethod()
1197	//
1198	void MCalibrationChargeCalc::FinalizeFFactorQECam()
1199	{
1200
1201	MCalibrationChargePix &avpix = (MCalibrationChargePix&)fCam->GetAverageArea(0);
1202	MCalibrationQEPix &qepix = (MCalibrationQEPix&) fQECam->GetAverageArea(0);
1203
1204	const Float_t avphotons = avpix.GetPheFFactorMethod()
1205	/ qepix.GetQEFFactor(fPulserColor)
1206	/ fQECam->GetPlexiglassQE();
1207
1208	const Float_t avphotrelvar = avpix.GetPheFFactorMethodRelVar()
1209	+ qepix.GetQEFFactorRelVar(fPulserColor)
1210	+ fQECam->GetPlexiglassQERelVar();
1211
1212	const UInt_t nareas = fGeom->GetNumAreas();
1213
1214	Float_t lowlim [nareas];
1215	Float_t upplim [nareas];
1216	Float_t avffactorphotons [nareas];
1217	Float_t avffactorphotvar [nareas];
1218	Int_t numffactor [nareas];
1219
1220	memset(lowlim ,0, nareas * sizeof(Float_t));
1221	memset(upplim ,0, nareas * sizeof(Float_t));
1222	memset(avffactorphotons,0, nareas * sizeof(Float_t));
1223	memset(avffactorphotvar,0, nareas * sizeof(Float_t));
1224	memset(numffactor ,0, nareas * sizeof(Int_t));
1225
1226	const UInt_t npixels = fGeom->GetNumPixels();
1227
1228	for (UInt_t i=0; i<npixels; i++)
1229	{
1230
1231	MCalibrationChargePix &pix = (MCalibrationChargePix&)(*fCam)[i];
1232	MCalibrationQEPix &qepix = (MCalibrationQEPix&) (*fQECam)[i];
1233
1234	if (!pix.IsFFactorMethodValid())
1235	{
1236	qepix.SetFFactorMethodValid(kFALSE,fPulserColor);
1237	continue;
1238	}
1239
1240	const Float_t photons = avphotons / fGeom->GetPixRatio(i);
1241	const Float_t qe = pix.GetPheFFactorMethod() / photons ;
1242
1243	if (!pix.CalcMeanFFactor( photons , avphotrelvar ))
1244	{
1245	pix.SetFFactorMethodValid(kFALSE);
1246	qepix.SetFFactorMethodValid(kFALSE, fPulserColor);
1247	(*fBadPixels)[i].SetUncalibrated( MBadPixelsPix::kDeviatingNumPhes );
1248	}
1249
1250	const Float_t qerelvar = avphotrelvar + pix.GetPheFFactorMethodRelVar();
1251
1252	qepix.SetQEFFactor ( qe , fPulserColor );
1253	qepix.SetQEFFactorVar ( qerelvarqeqe, fPulserColor );
1254	qepix.SetFFactorMethodValid( kTRUE , fPulserColor );
1255
1256	if (!qepix.UpdateFFactorMethod())
1257	*fLog << warn << GetDescriptor()
1258	<< ": Cannot update Quantum efficiencies with the F-Factor Method" << endl;
1259
1260	const Int_t aidx = (*fGeom)[i].GetAidx();
1261
1262	avffactorphotons[aidx] += pix.GetMeanFFactorFADC2Phot();
1263	avffactorphotvar[aidx] += pix.GetMeanFFactorFADC2PhotVar();
1264	numffactor[aidx]++;
1265	}
1266
1267	for (UInt_t i=0; i<nareas; i++)
1268	{
1269	avffactorphotons[i] /= numffactor[i];
1270	avffactorphotvar[i] /= numffactor[i];
1271	lowlim [i] = avffactorphotons[i] - fFFactorErrLimit*TMath::Sqrt(avffactorphotvar[i]);
1272	upplim [i] = avffactorphotons[i] + fFFactorErrLimit*TMath::Sqrt(avffactorphotvar[i]);
1273	}
1274
1275	for (UInt_t i=0; i<npixels; i++)
1276	{
1277
1278	MCalibrationChargePix &pix = (MCalibrationChargePix&)(*fCam)[i];
1279
1280	if (!pix.IsFFactorMethodValid())
1281	continue;
1282
1283	const Float_t ffactor = pix.GetMeanFFactorFADC2Phot();
1284	MBadPixelsPix &bad = (*fBadPixels)[i];
1285
1286	const Int_t aidx = (*fGeom)[i].GetAidx();
1287
1288	if ( ffactor < lowlim[aidx] \|\| ffactor > upplim[aidx] )
1289	{
1290	*fLog << warn << GetDescriptor() << ": Deviating F-Factor: "
1291	<< Form("%4.2f",ffactor) << " out of accepted limits: ["
1292	<< Form("%4.2f%s%4.2f",lowlim[aidx],",",upplim[aidx]) << "] in pixel " << i << endl;
1293	bad.SetUncalibrated( MBadPixelsPix::kDeviatingFFactor );
1294	bad.SetUnsuitable ( MBadPixelsPix::kUnsuitableRun );
1295	pix.SetExcluded();
1296	}
1297	}
1298	}
1299
1300
1301	// ------------------------------------------------------------------------
1302	//
1303	// Loop over pixels:
1304	//
1305	// - Continue, if not MCalibrationChargeBlindPix::IsFluxInsidePlexiglassAvailable() and set:
1306	// MCalibrationQEPix::SetBlindPixelMethodValid(kFALSE,fPulserColor)
1307	//
1308	// - Calculate the quantum efficiency with the formula:
1309	//
1310	// QE = Num.Phes / MCalibrationChargeBlindPix::GetFluxInsidePlexiglass()
1311	// / MGeomPix::GetA() * MCalibrationQECam::GetPlexiglassQE()
1312	//
1313	// - Set QE in MCalibrationQEPix::SetQEBlindPixel ( QE, fPulserColor );
1314	// - Set Variance of QE in MCalibrationQEPix::SetQEBlindPixelVar ( Variance, fPulserColor );
1315	// - Set bit MCalibrationQEPix::SetBlindPixelMethodValid(kTRUE,fPulserColor)
1316	//
1317	// - Call MCalibrationQEPix::UpdateBlindPixelMethod()
1318	//
1319	void MCalibrationChargeCalc::FinalizeBlindPixelQECam()
1320	{
1321
1322	const UInt_t npixels = fGeom->GetNumPixels();
1323
1324	//
1325	// With the knowledge of the overall photon flux, calculate the
1326	// quantum efficiencies after the Blind Pixel and PIN Diode method
1327	//
1328	for (UInt_t i=0; i<npixels; i++)
1329	{
1330
1331	MCalibrationQEPix &qepix = (MCalibrationQEPix&) (*fQECam)[i];
1332
1333	if (!fBlindPixel)
1334	{
1335	qepix.SetBlindPixelMethodValid(kFALSE, fPulserColor);
1336	continue;
1337	}
1338
1339	if (!fBlindPixel->IsFluxInsidePlexiglassAvailable())
1340	{
1341	qepix.SetBlindPixelMethodValid(kFALSE, fPulserColor);
1342	continue;
1343	}
1344
1345	MBadPixelsPix &bad = (*fBadPixels)[i];
1346
1347	if (!bad.IsUnsuitable (MBadPixelsPix::kUnsuitableRun))
1348	{
1349	qepix.SetBlindPixelMethodValid(kFALSE, fPulserColor);
1350	continue;
1351	}
1352
1353	MCalibrationChargePix &pix = (MCalibrationChargePix&)(*fCam)[i];
1354	MGeomPix &geo = (*fGeom)[i];
1355
1356	const Float_t qe = pix.GetPheFFactorMethod()
1357	/ fBlindPixel->GetFluxInsidePlexiglass()
1358	/ geo.GetA()
1359	* fQECam->GetPlexiglassQE();
1360
1361	const Float_t qerelvar = fBlindPixel->GetFluxInsidePlexiglassRelVar()
1362	+ fQECam->GetPlexiglassQERelVar()
1363	+ pix.GetPheFFactorMethodRelVar();
1364
1365	qepix.SetQEBlindPixel ( qe , fPulserColor );
1366	qepix.SetQEBlindPixelVar ( qerelvarqeqe, fPulserColor );
1367	qepix.UpdateBlindPixelMethod();
1368	}
1369	}
1370
1371	// ------------------------------------------------------------------------
1372	//
1373	// Loop over pixels:
1374	//
1375	// - Continue, if not MCalibrationChargePINDiode::IsFluxOutsidePlexiglassAvailable() and set:
1376	// MCalibrationQEPix::SetPINDiodeMethodValid(kFALSE,fPulserColor)
1377	//
1378	// - Calculate the quantum efficiency with the formula:
1379	//
1380	// QE = Num.Phes / MCalibrationChargePINDiode::GetFluxOutsidePlexiglass() / MGeomPix::GetA()
1381	//
1382	// - Set QE in MCalibrationQEPix::SetQEPINDiode ( QE, fPulserColor );
1383	// - Set Variance of QE in MCalibrationQEPix::SetQEPINDiodeVar ( Variance, fPulserColor );
1384	// - Set bit MCalibrationQEPix::SetPINDiodeMethodValid(kTRUE,fPulserColor)
1385	//
1386	// - Call MCalibrationQEPix::UpdatePINDiodeMethod()
1387	//
1388	void MCalibrationChargeCalc::FinalizePINDiodeQECam()
1389	{
1390
1391	const UInt_t npixels = fGeom->GetNumPixels();
1392
1393	//
1394	// With the knowledge of the overall photon flux, calculate the
1395	// quantum efficiencies after the PIN Diode method
1396	//
1397	for (UInt_t i=0; i<npixels; i++)
1398	{
1399
1400	MCalibrationQEPix &qepix = (MCalibrationQEPix&) (*fQECam)[i];
1401
1402	if (!fPINDiode)
1403	{
1404	qepix.SetPINDiodeMethodValid(kFALSE, fPulserColor);
1405	continue;
1406	}
1407
1408	if (!fPINDiode->IsFluxOutsidePlexiglassAvailable())
1409	{
1410	qepix.SetPINDiodeMethodValid(kFALSE, fPulserColor);
1411	continue;
1412	}
1413
1414	MBadPixelsPix &bad = (*fBadPixels)[i];
1415
1416	if (!bad.IsUnsuitable (MBadPixelsPix::kUnsuitableRun))
1417	{
1418	qepix.SetPINDiodeMethodValid(kFALSE, fPulserColor);
1419	continue;
1420	}
1421
1422	MCalibrationChargePix &pix = (MCalibrationChargePix&)(*fCam)[i];
1423	MGeomPix &geo = (*fGeom)[i];
1424
1425	const Float_t qe = pix.GetPheFFactorMethod()
1426	/ fPINDiode->GetFluxOutsidePlexiglass()
1427	/ geo.GetA();
1428
1429	const Float_t qerelvar = fPINDiode->GetFluxOutsidePlexiglassRelVar() + pix.GetPheFFactorMethodRelVar();
1430
1431	qepix.SetQEPINDiode ( qe , fPulserColor );
1432	qepix.SetQEPINDiodeVar ( qerelvarqeqe, fPulserColor );
1433	qepix.UpdateBlindPixelMethod();
1434	}
1435	}
1436
1437	// -----------------------------------------------------------------------------------------------
1438	//
1439	// - Print out statistics about BadPixels of type UnsuitableType_t
1440	// - store numbers of bad pixels of each type in fCam
1441	//
1442	void MCalibrationChargeCalc::FinalizeUnsuitablePixels()
1443	{
1444
1445	*fLog << inf << endl;
1446	*fLog << GetDescriptor() << ": Charge Calibration status:" << endl;
1447	*fLog << dec << setfill(' ');
1448
1449	const Int_t nareas = fGeom->GetNumAreas();
1450
1451	Int_t counts[nareas];
1452	memset(counts,0,nareas*sizeof(Int_t));
1453
1454	for (Int_t i=0; i<fBadPixels->GetSize(); i++)
1455	{
1456	MBadPixelsPix &bad = (*fBadPixels)[i];
1457	if (!bad.IsBad())
1458	{
1459	const Int_t aidx = (*fGeom)[i].GetAidx();
1460	counts[aidx]++;
1461	}
1462	}
1463
1464	if (fGeom->InheritsFrom("MGeomCamMagic"))
1465	*fLog << " " << setw(7) << "Successfully calibrated Pixels: "
1466	<< Form("%s%3i%s%3i","Inner: ",counts[0]," Outer: ",counts[1]) << endl;
1467
1468	memset(counts,0,nareas*sizeof(Int_t));
1469
1470	for (Int_t i=0; i<fBadPixels->GetSize(); i++)
1471	{
1472	MBadPixelsPix &bad = (*fBadPixels)[i];
1473	if (bad.IsUnsuitable(MBadPixelsPix::kUnsuitableRun))
1474	{
1475	const Int_t aidx = (*fGeom)[i].GetAidx();
1476	counts[aidx]++;
1477	}
1478	}
1479
1480	for (Int_t aidx=0; aidx<nareas; aidx++)
1481	fCam->SetNumUnsuitable(counts[aidx], aidx);
1482
1483	if (fGeom->InheritsFrom("MGeomCamMagic"))
1484	*fLog << " " << setw(7) << "Uncalibrated Pixels: "
1485	<< Form("%s%3i%s%3i","Inner: ",counts[0]," Outer: ",counts[1]) << endl;
1486
1487	memset(counts,0,nareas*sizeof(Int_t));
1488
1489	for (Int_t i=0; i<fBadPixels->GetSize(); i++)
1490	{
1491	MBadPixelsPix &bad = (*fBadPixels)[i];
1492	if (bad.IsUnsuitable(MBadPixelsPix::kUnreliableRun))
1493	{
1494	const Int_t aidx = (*fGeom)[i].GetAidx();
1495	counts[aidx]++;
1496	}
1497	}
1498
1499	for (Int_t aidx=0; aidx<nareas; aidx++)
1500	fCam->SetNumUnreliable(counts[aidx], aidx);
1501
1502	*fLog << " " << setw(7) << "Unreliable Pixels: "
1503	<< Form("%s%3i%s%3i","Inner: ",counts[0]," Outer: ",counts[1]) << endl;
1504
1505	}
1506
1507	// -----------------------------------------------------------------------------------------------
1508	//
1509	// Print out statistics about BadPixels of type UncalibratedType_t
1510	//
1511	void MCalibrationChargeCalc::PrintUncalibrated(MBadPixelsPix::UncalibratedType_t typ, const char *text) const
1512	{
1513
1514	UInt_t countinner = 0;
1515	UInt_t countouter = 0;
1516	for (Int_t i=0; i<fBadPixels->GetSize(); i++)
1517	{
1518	MBadPixelsPix &bad = (*fBadPixels)[i];
1519	if (bad.IsUncalibrated(typ))
1520	{
1521	if (fGeom->GetPixRatio(i) == 1.)
1522	countinner++;
1523	else
1524	countouter++;
1525	}
1526	}
1527
1528	*fLog << " " << setw(7) << text
1529	<< Form("%s%3i%s%3i","Inner: ",countinner," Outer: ",countouter) << endl;
1530	}
1531

Note: See TracBrowser for help on using the repository browser.

Download in other formats: