Context Navigation

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

Visit:

Last change on this file since 3481 was 3479, checked in by gaug, 21 years ago
* empty log message *
File size: 39.3 KB

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

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

Download in other formats: