source: trunk/MagicSoft/Mars/mcalib/MHCalibrationChargeCam.cc@ 4323

Last change on this file since 4323 was 4315, checked in by gaug, 21 years ago
*** empty log message ***
File size: 31.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// MHCalibrationChargeCam
27//
28// Fills the extracted signals of MExtractedSignalCam into the MHGausEvents-classes
29// MHCalibrationChargeHiGainPix and MHCalibrationChargeLoGainPix for every:
30//
31// - Pixel, stored in the TObjArray's MHCalibrationCam::fHiGainArray and
32// MHCalibrationCam::fLoGainArray
33//
34// - Average pixel per AREA index (e.g. inner and outer for the MAGIC camera),
35// stored in the TObjArray's MHCalibrationCam::fAverageHiGainAreas and
36// MHCalibrationCam::fAverageLoGainAreas
37//
38// - Average pixel per camera SECTOR (e.g. sectors 1-6 for the MAGIC camera),
39// stored in the TObjArray's MHCalibrationCam::fAverageHiGainSectors and
40// MHCalibrationCam::fAverageLoGainSectors
41//
42// Every signal is taken from MExtractedSignalCam and filled into a histogram and
43// an array, in order to perform a Fourier analysis (see MHGausEvents).
44// The signals are moreover averaged on an event-by-event basis and written into
45// the corresponding average pixels.
46//
47// Additionally, the (FADC slice) position of the maximum is stored in an Absolute
48// Arrival Time histogram. This histogram serves for a rough cross-check if the
49// signal does not lie at or outside the edges of the extraction window.
50//
51// The Charge histograms are fitted to a Gaussian, mean and sigma with its errors
52// and the fit probability are extracted. If none of these values are NaN's and
53// if the probability is bigger than MHGausEvents::fProbLimit (default: 0.5%),
54// the fit is declared valid.
55// Otherwise, the fit is repeated within ranges of the previous mean
56// +- MHGausEvents::fPickupLimit (default: 5) sigma (see MHGausEvents::RepeatFit())
57// In case this does not make the fit valid, the histogram means and RMS's are
58// taken directly (see MHGausEvents::BypassFit()) and the following flags are set:
59// - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::kHiGainNotFitted ) or
60// - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::kLoGainNotFitted ) and
61// - MBadPixelsPix::SetUnsuitable( MBadPixelsPix::kUnreliableRun )
62//
63// Outliers of more than MHGausEvents::fPickupLimit (default: 5) sigmas
64// from the mean are counted as Pickup events (stored in MHGausEvents::fPickup)
65//
66// Unless more than fNumHiGainSaturationLimit (default: 1%) of the overall FADC
67// slices show saturation, the following flag is set:
68// - MCalibrationChargePix::SetHiGainSaturation();
69// In that case, the calibration constants are derived from the low-gain results.
70//
71// If more than fNumLoGainSaturationLimit (default: 1%) of the overall
72// low-gain FADC slices saturate, the following flags are set:
73// - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::kLoGainSaturation ) and
74// - MBadPixelsPix::SetUnsuitable( MBadPixelsPix::kUnsuitableRun )
75//
76// The class also fills arrays with the signal vs. event number, creates a fourier
77// spectrum and investigates if the projected fourier components follow an exponential
78// distribution. In case that the probability of the exponential fit is less than
79// MHGausEvents::fProbLimit (default: 0.5%), the following flags are set:
80// - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::kHiGainOscillating ) or
81// - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::kLoGainOscillating ) and
82// - MBadPixelsPix::SetUnsuitable( MBadPixelsPix::kUnreliableRun )
83//
84// This same procedure is performed for the average pixels.
85//
86// The following results are written into MCalibrationChargeCam:
87//
88// - MCalibrationPix::SetHiGainSaturation()
89// - MCalibrationPix::SetHiGainMean()
90// - MCalibrationPix::SetHiGainMeanErr()
91// - MCalibrationPix::SetHiGainSigma()
92// - MCalibrationPix::SetHiGainSigmaErr()
93// - MCalibrationPix::SetHiGainProb()
94// - MCalibrationPix::SetHiGainNumPickup()
95//
96// - MCalibrationPix::SetLoGainMean()
97// - MCalibrationPix::SetLoGainMeanErr()
98// - MCalibrationPix::SetLoGainSigma()
99// - MCalibrationPix::SetLoGainSigmaErr()
100// - MCalibrationPix::SetLoGainProb()
101// - MCalibrationPix::SetLoGainNumPickup()
102//
103// - MCalibrationChargePix::SetAbsTimeMean()
104// - MCalibrationChargePix::SetAbsTimeRms()
105//
106// For all averaged areas, the fitted sigma is multiplied with the square root of
107// the number involved pixels in order to be able to compare it to the average of
108// sigmas in the camera.
109//
110/////////////////////////////////////////////////////////////////////////////
111#include "MHCalibrationChargeCam.h"
112#include "MHCalibrationCam.h"
113
114#include "MLog.h"
115#include "MLogManip.h"
116
117#include "MParList.h"
118
119#include "MHCalibrationChargeHiGainPix.h"
120#include "MHCalibrationChargeLoGainPix.h"
121#include "MHCalibrationChargePix.h"
122
123#include "MCalibrationChargeCam.h"
124#include "MCalibrationChargePix.h"
125
126#include "MGeomCam.h"
127#include "MGeomPix.h"
128
129#include "MHGausEvents.h"
130
131#include "MBadPixelsCam.h"
132#include "MBadPixelsPix.h"
133
134#include "MRawEvtData.h"
135#include "MRawEvtPixelIter.h"
136
137#include "MExtractedSignalCam.h"
138#include "MExtractedSignalPix.h"
139
140#include <TPad.h>
141#include <TVirtualPad.h>
142#include <TCanvas.h>
143
144ClassImp(MHCalibrationChargeCam);
145
146using namespace std;
147
148const Float_t MHCalibrationChargeCam::fgNumHiGainSaturationLimit = 0.01;
149const Float_t MHCalibrationChargeCam::fgNumLoGainSaturationLimit = 0.005;
150const Float_t MHCalibrationChargeCam::fgTimeLowerLimit = 1.;
151const Float_t MHCalibrationChargeCam::fgTimeUpperLimit = 2.;
152// --------------------------------------------------------------------------
153//
154// Default Constructor.
155//
156// Sets:
157// - all pointers to NULL
158//
159// Initializes:
160// - fNumHiGainSaturationLimit to fgNumHiGainSaturationLimit
161// - fNumLoGainSaturationLimit to fgNumLoGainSaturationLimit
162// - fTimeLowerLimit to fgTimeLowerLimit
163// - fTimeUpperLimit to fgTimeUpperLimit
164//
165MHCalibrationChargeCam::MHCalibrationChargeCam(const char *name, const char *title)
166 : fRawEvt(NULL)
167{
168 fName = name ? name : "MHCalibrationChargeCam";
169 fTitle = title ? title : "Class to fill the calibration histograms ";
170
171 SetNumHiGainSaturationLimit(fgNumHiGainSaturationLimit);
172 SetNumLoGainSaturationLimit(fgNumLoGainSaturationLimit);
173 SetTimeLowerLimit();
174 SetTimeUpperLimit();
175}
176
177// --------------------------------------------------------------------------
178//
179// Gets the pointers to:
180// - MRawEvtData
181//
182Bool_t MHCalibrationChargeCam::SetupHists(const MParList *pList)
183{
184
185 fRawEvt = (MRawEvtData*)pList->FindObject("MRawEvtData");
186 if (!fRawEvt)
187 {
188 *fLog << err << dbginf << "MRawEvtData not found... aborting." << endl;
189 return kFALSE;
190 }
191
192 return kTRUE;
193}
194
195// --------------------------------------------------------------------------
196//
197// Gets or creates the pointers to:
198// - MExtractedSignalCam
199// - MCalibrationChargeCam
200// - MBadPixelsCam
201//
202// Initializes the number of used FADC slices from MExtractedSignalCam
203// into MCalibrationChargeCam and test for changes in that variable
204//
205// Initializes, if empty to MGeomCam::GetNumPixels():
206// - MHCalibrationCam::fHiGainArray, MHCalibrationCam::fLoGainArray
207//
208// Initializes, if empty to MGeomCam::GetNumAreas() for:
209// - MHCalibrationCam::fAverageHiGainAreas, MHCalibrationCam::fAverageLoGainAreas
210//
211// Initializes, if empty to MGeomCam::GetNumSectors() for:
212// - MHCalibrationCam::fAverageHiGainSectors, MHCalibrationCam::fAverageLoGainSectors
213//
214// Calls MHCalibrationCam::InitHists() for every entry in:
215// - MHCalibrationCam::fHiGainArray, MHCalibrationCam::fLoGainArray
216// - MHCalibrationCam::fAverageHiGainAreas, MHCalibrationCam::fAverageLoGainAreas
217// - MHCalibrationCam::fAverageHiGainSectors, MHCalibrationCam::fAverageLoGainSectors
218//
219// Sets Titles and Names for the Charge Histograms:
220// - MHCalibrationCam::fAverageHiGainAreas
221// - MHCalibrationCam::fAverageHiGainSectors
222//
223// Sets number of bins to MHCalibrationCam::fAverageNbins for:
224// - MHCalibrationCam::fAverageHiGainAreas, MHCalibrationCam::fAverageLoGainAreas
225// - MHCalibrationCam::fAverageHiGainSectors, MHCalibrationCam::fAverageLoGainSectors
226//
227Bool_t MHCalibrationChargeCam::ReInitHists(MParList *pList)
228{
229
230 MExtractedSignalCam *signal = (MExtractedSignalCam*)pList->FindObject("MExtractedSignalCam");
231 if (!signal)
232 {
233 *fLog << err << "MExtractedSignalCam not found... abort." << endl;
234 return kFALSE;
235 }
236
237 fCam = (MCalibrationCam*)pList->FindObject("MCalibrationChargeCam");
238 if (!fCam)
239 {
240 fCam = (MCalibrationCam*)pList->FindCreateObj(AddSerialNumber("MCalibrationChargeCam"));
241 if (!fCam)
242 {
243 gLog << err << "Cannot find nor create MCalibrationChargeCam ... abort." << endl;
244 return kFALSE;
245 }
246 else
247 fCam->Init(*fGeom);
248 }
249
250 fFirstHiGain = signal->GetFirstUsedSliceHiGain();
251 fLastHiGain = signal->GetLastUsedSliceHiGain();
252 fFirstLoGain = signal->GetFirstUsedSliceLoGain();
253 fLastLoGain = signal->GetLastUsedSliceLoGain();
254
255 const Float_t numhigain = signal->GetNumUsedHiGainFADCSlices();
256 const Float_t numlogain = signal->GetNumUsedLoGainFADCSlices();
257
258 if (fCam->GetNumHiGainFADCSlices() == 0.)
259 fCam->SetNumHiGainFADCSlices ( numhigain );
260 else if (fCam->GetNumHiGainFADCSlices() != numhigain)
261 {
262 *fLog << err << GetDescriptor()
263 << ": Number of High Gain FADC extraction slices has changed, abort..." << endl;
264 return kFALSE;
265 }
266
267 if (fCam->GetNumLoGainFADCSlices() == 0.)
268 fCam->SetNumLoGainFADCSlices ( numlogain );
269 else if (fCam->GetNumLoGainFADCSlices() != numlogain)
270 {
271 *fLog << err << GetDescriptor()
272 << ": Number of Low Gain FADC extraction slices has changes, abort..." << endl;
273 return kFALSE;
274 }
275
276 const Int_t npixels = fGeom->GetNumPixels();
277 const Int_t nsectors = fGeom->GetNumSectors();
278 const Int_t nareas = fGeom->GetNumAreas();
279
280 if (fHiGainArray->GetEntries()==0)
281 {
282 fHiGainArray->Expand(npixels);
283 for (Int_t i=0; i<npixels; i++)
284 {
285 (*fHiGainArray)[i] = new MHCalibrationChargeHiGainPix;
286 InitHists((*this)[i],(*fBadPixels)[i],i);
287 }
288 }
289
290 if (fLoGainArray->GetEntries()==0)
291 {
292 fLoGainArray->Expand(npixels);
293
294 for (Int_t i=0; i<npixels; i++)
295 {
296 (*fLoGainArray)[i] = new MHCalibrationChargeLoGainPix;
297 InitHists((*this)(i),(*fBadPixels)[i],i);
298 }
299
300 }
301
302 if (fAverageHiGainAreas->GetEntries()==0)
303 {
304 fAverageHiGainAreas->Expand(nareas);
305
306 for (Int_t j=0; j<nareas; j++)
307 {
308 (*fAverageHiGainAreas)[j] =
309 new MHCalibrationChargeHiGainPix("AverageHiGainArea",
310 "Average HiGain FADC sums area idx ");
311
312 MHCalibrationChargePix &hist = (MHCalibrationChargePix&)GetAverageHiGainArea(j);
313
314 hist.GetHGausHist()->SetTitle("Summed FADC slices average HiGain Area Idx ");
315 hist.SetNbins(fAverageNbins);
316 hist.SetLast(2.*hist.GetLast());
317 hist.GetHAbsTime()->SetTitle("Absolute Arrival Time average HiGain Area Idx ");
318
319 InitHists(hist,fCam->GetAverageBadArea(j),j);
320
321 }
322 }
323
324
325 if (fAverageLoGainAreas->GetEntries()==0)
326 {
327 fAverageLoGainAreas->Expand(nareas);
328
329 for (Int_t j=0; j<nareas; j++)
330 {
331 (*fAverageLoGainAreas)[j] =
332 new MHCalibrationChargeLoGainPix("AverageLoGainArea",
333 "Average LoGain FADC sums of pixel area idx ");
334
335 MHCalibrationChargePix &hist = (MHCalibrationChargePix&)GetAverageLoGainArea(j);
336
337 hist.GetHGausHist()->SetTitle("Summed FADC slices average LoGain Area Idx ");
338 hist.SetNbins(fAverageNbins);
339 hist.GetHAbsTime()->SetTitle("Absolute Arrival Time average LoGain Area Idx ");
340
341 InitHists(hist,fCam->GetAverageBadArea(j),j);
342
343 }
344 }
345
346 if (fAverageHiGainSectors->GetEntries()==0)
347 {
348 fAverageHiGainSectors->Expand(nsectors);
349
350 for (Int_t j=0; j<nsectors; j++)
351 {
352 (*fAverageHiGainSectors)[j] =
353 new MHCalibrationChargeHiGainPix("AverageHiGainSector",
354 "Average HiGain FADC sums of pixel sector ");
355
356 MHCalibrationChargePix &hist = (MHCalibrationChargePix&)GetAverageHiGainSector(j);
357
358 hist.GetHGausHist()->SetTitle("Summed FADC slices average HiGain Sector ");
359 hist.SetNbins(fAverageNbins);
360 hist.SetLast (2.*hist.GetLast());
361 hist.GetHAbsTime()->SetTitle("Absolute Arrival Time average HiGain Sector ");
362
363 InitHists(hist,fCam->GetAverageBadSector(j),j);
364
365 }
366 }
367
368 if (fAverageLoGainSectors->GetEntries()==0)
369 {
370 fAverageLoGainSectors->Expand(nsectors);
371
372 for (Int_t j=0; j<nsectors; j++)
373 {
374 (*fAverageLoGainSectors)[j] =
375 new MHCalibrationChargeLoGainPix("AverageLoGainSector",
376 "Average LoGain FADC sums of pixel sector ");
377
378 MHCalibrationChargePix &hist = (MHCalibrationChargePix&)GetAverageLoGainSector(j);
379
380 hist.GetHGausHist()->SetTitle("Summed FADC slices average LoGain Sector ");
381 hist.SetNbins(fAverageNbins);
382 hist.GetHAbsTime()->SetTitle("Absolute Arrival Time average LoGain Sector ");
383
384 InitHists(hist,fCam->GetAverageBadSector(j),j);
385
386 }
387 }
388
389 return kTRUE;
390}
391
392
393// --------------------------------------------------------------------------
394//
395// Retrieves from MExtractedSignalCam:
396// - first used LoGain FADC slice
397//
398// Retrieves from MGeomCam:
399// - number of pixels
400// - number of pixel areas
401// - number of sectors
402//
403// For all TObjArray's (including the averaged ones), the following steps are performed:
404//
405// 1) Fill Charges histograms (MHGausEvents::FillHistAndArray()) with:
406// - MExtractedSignalPix::GetExtractedSignalHiGain();
407// - MExtractedSignalPix::GetExtractedSignalLoGain();
408//
409// 2) Set number of saturated slices (MHCalibrationChargePix::SetSaturated()) with:
410// - MExtractedSignalPix::GetNumHiGainSaturated();
411// - MExtractedSignalPix::GetNumLoGainSaturated();
412//
413// 3) Fill AbsTime histograms (MHCalibrationChargePix::FillAbsTime()) with:
414// - MRawEvtPixelIter::GetIdxMaxHiGainSample();
415// - MRawEvtPixelIter::GetIdxMaxLoGainSample(first slice);
416//
417Bool_t MHCalibrationChargeCam::FillHists(const MParContainer *par, const Stat_t w)
418{
419
420 MExtractedSignalCam *signal = (MExtractedSignalCam*)par;
421 if (!signal)
422 {
423 *fLog << err << "No argument in MExtractedSignalCam::Fill... abort." << endl;
424 return kFALSE;
425 }
426
427 const UInt_t npixels = fGeom->GetNumPixels();
428 const UInt_t nareas = fGeom->GetNumAreas();
429 const UInt_t nsectors = fGeom->GetNumSectors();
430 const UInt_t lofirst = signal->GetFirstUsedSliceLoGain();
431
432 Float_t sumhiarea [nareas], sumloarea [nareas], timehiarea [nareas], timeloarea [nareas];
433 Float_t sumhisector[nsectors], sumlosector[nsectors], timehisector[nsectors], timelosector[nsectors];
434 Int_t sathiarea [nareas], satloarea [nareas];
435 Int_t sathisector[nsectors], satlosector[nsectors];
436
437 memset(sumhiarea, 0, nareas * sizeof(Float_t));
438 memset(sumloarea, 0, nareas * sizeof(Float_t));
439 memset(timehiarea, 0, nareas * sizeof(Float_t));
440 memset(timeloarea, 0, nareas * sizeof(Float_t));
441 memset(sathiarea, 0, nareas * sizeof(Int_t ));
442 memset(satloarea, 0, nareas * sizeof(Int_t ));
443 memset(sumhisector, 0, nsectors*sizeof(Float_t));
444 memset(sumlosector, 0, nsectors*sizeof(Float_t));
445 memset(timehisector,0, nsectors*sizeof(Float_t));
446 memset(timelosector,0, nsectors*sizeof(Float_t));
447 memset(sathisector, 0, nsectors*sizeof(Int_t ));
448 memset(satlosector, 0, nsectors*sizeof(Int_t ));
449
450 for (UInt_t i=0; i<npixels; i++)
451 {
452
453 MHCalibrationChargePix &histhi = (MHCalibrationChargePix&)(*this)[i];
454 MHCalibrationChargePix &histlo = (MHCalibrationChargePix&)(*this)(i);
455
456 if (histhi.IsExcluded())
457 continue;
458
459 const MExtractedSignalPix &pix = (*signal)[i];
460
461 const Float_t sumhi = pix.GetExtractedSignalHiGain();
462 const Float_t sumlo = pix.GetExtractedSignalLoGain();
463
464 histhi.FillHistAndArray(sumhi);
465 histlo.FillHistAndArray(sumlo);
466
467 const Int_t sathi = (Int_t)pix.GetNumHiGainSaturated();
468 const Int_t satlo = (Int_t)pix.GetNumLoGainSaturated();
469
470 histhi.SetSaturated(sathi);
471 histlo.SetSaturated(satlo);
472
473 const Int_t aidx = (*fGeom)[i].GetAidx();
474 const Int_t sector = (*fGeom)[i].GetSector();
475
476 sumhiarea[aidx] += sumhi;
477 sumloarea[aidx] += sumlo;
478 sathiarea[aidx] += sathi;
479 satloarea[aidx] += satlo;
480
481 sumhisector[sector] += sumhi;
482 sumlosector[sector] += sumlo;
483 sathisector[sector] += sathi;
484 satlosector[sector] += satlo;
485 }
486
487 MRawEvtPixelIter pixel(fRawEvt);
488 while (pixel.Next())
489 {
490
491 const UInt_t pixid = pixel.GetPixelId();
492
493 MHCalibrationChargePix &histhi = (MHCalibrationChargePix&)(*this)[pixid];
494 MHCalibrationChargePix &histlo = (MHCalibrationChargePix&)(*this)(pixid);
495
496 if (histhi.IsExcluded())
497 continue;
498
499 const Float_t timehi = (Float_t)pixel.GetIdxMaxHiGainSample();
500 const Float_t timelo = (Float_t)pixel.GetIdxMaxLoGainSample(lofirst);
501
502 histhi.FillAbsTime(timehi);
503 histlo.FillAbsTime(timelo);
504
505 const Int_t aidx = (*fGeom)[pixid].GetAidx();
506 const Int_t sector = (*fGeom)[pixid].GetSector();
507
508 timehiarea[aidx] += timehi;
509 timeloarea[aidx] += timelo;
510
511 timehisector[sector] += timehi;
512 timelosector[sector] += timelo;
513 }
514
515 for (UInt_t j=0; j<nareas; j++)
516 {
517
518 const Int_t npix = fAverageAreaNum[j];
519
520 MHCalibrationChargePix &hipix = (MHCalibrationChargePix&)GetAverageHiGainArea(j);
521 MHCalibrationChargePix &lopix = (MHCalibrationChargePix&)GetAverageLoGainArea(j);
522
523 hipix.FillHistAndArray(sumhiarea[j]/npix);
524 lopix.FillHistAndArray(sumloarea[j]/npix);
525
526 hipix.SetSaturated((Float_t)sathiarea[j]/npix);
527 lopix.SetSaturated((Float_t)satloarea[j]/npix);
528
529 hipix.FillAbsTime(timehiarea[j]/npix);
530 lopix.FillAbsTime(timeloarea[j]/npix);
531
532 }
533
534 for (UInt_t j=0; j<nsectors; j++)
535 {
536
537 const Int_t npix = fAverageSectorNum[j];
538
539 MHCalibrationChargePix &hipix = (MHCalibrationChargePix&)GetAverageHiGainSector(j);
540 MHCalibrationChargePix &lopix = (MHCalibrationChargePix&)GetAverageLoGainSector(j);
541
542 hipix.FillHistAndArray(sumhisector[j]/npix);
543 lopix.FillHistAndArray(sumlosector[j]/npix);
544
545 hipix.SetSaturated((Float_t)sathisector[j]/npix);
546 lopix.SetSaturated((Float_t)satlosector[j]/npix);
547
548 hipix.FillAbsTime(timehisector[j]/npix);
549 lopix.FillAbsTime(timelosector[j]/npix);
550
551 }
552
553 return kTRUE;
554}
555
556// --------------------------------------------------------------------------
557//
558// For all TObjArray's (including the averaged ones), the following steps are performed:
559//
560// 1) Returns if the pixel is excluded.
561// 2) Tests saturation. In case yes, set the flag: MCalibrationPix::SetHiGainSaturation()
562// or the flag: MBadPixelsPix::SetUncalibrated( MBadPixelsPix::kLoGainSaturated )
563// 3) Store the absolute arrival times in the MCalibrationChargePix's. If flag
564// MCalibrationPix::IsHiGainSaturation() is set, the Low-Gain arrival times are stored,
565// otherwise the Hi-Gain ones.
566// 4) Calls to MHCalibrationCam::FitHiGainArrays() and MCalibrationCam::FitLoGainArrays()
567// with the flags:
568// - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::kHiGainNotFitted )
569// - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::kLoGainNotFitted )
570// - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::kHiGainOscillating )
571// - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::kLoGainOscillating )
572//
573Bool_t MHCalibrationChargeCam::FinalizeHists()
574{
575
576 for (Int_t i=0; i<fHiGainArray->GetSize(); i++)
577 {
578
579 MHCalibrationChargePix &histhi = (MHCalibrationChargePix&)(*this)[i];
580 MCalibrationChargePix &pix = (MCalibrationChargePix&)(*fCam)[i];
581 MBadPixelsPix &bad = (*fBadPixels)[i];
582
583 if (histhi.IsExcluded())
584 continue;
585
586 if (histhi.GetSaturated() > fNumHiGainSaturationLimit*histhi.GetHGausHist()->GetEntries())
587 {
588 pix.SetHiGainSaturation();
589 histhi.CreateFourierSpectrum();
590 continue;
591 }
592
593 FinalizeAbsTimes(histhi, pix, bad, fFirstHiGain, fLastHiGain);
594 }
595
596 for (Int_t i=0; i<fLoGainArray->GetSize(); i++)
597 {
598
599 MHCalibrationChargePix &histlo = (MHCalibrationChargePix&)(*this)(i);
600 MBadPixelsPix &bad = (*fBadPixels)[i];
601
602 if (histlo.IsExcluded())
603 continue;
604
605 if (histlo.GetSaturated() > fNumLoGainSaturationLimit*histlo.GetHGausHist()->GetEntries())
606 {
607 *fLog << warn << "Saturated Lo Gain histogram in pixel: " << i << endl;
608 bad.SetUncalibrated( MBadPixelsPix::kLoGainSaturation );
609 histlo.CreateFourierSpectrum();
610 continue;
611 }
612
613 MCalibrationChargePix &pix = (MCalibrationChargePix&)(*fCam)[i];
614
615 if (pix.IsHiGainSaturation())
616 FinalizeAbsTimes(histlo, pix, bad, fFirstLoGain, fLastLoGain);
617 }
618
619 for (Int_t j=0; j<fAverageHiGainAreas->GetSize(); j++)
620 {
621
622 MHCalibrationChargePix &histhi = (MHCalibrationChargePix&)GetAverageHiGainArea(j);
623 MCalibrationChargePix &pix = (MCalibrationChargePix&)fCam->GetAverageArea(j);
624 MBadPixelsPix &bad = fCam->GetAverageBadArea(j);
625
626 if (histhi.GetSaturated() > fNumHiGainSaturationLimit*histhi.GetHGausHist()->GetEntries())
627 {
628 pix.SetHiGainSaturation();
629 histhi.CreateFourierSpectrum();
630 continue;
631 }
632
633 FinalizeAbsTimes(histhi, pix, bad, fFirstHiGain, fLastHiGain);
634 }
635
636 for (Int_t j=0; j<fAverageLoGainAreas->GetSize(); j++)
637 {
638
639 MHCalibrationChargePix &histlo = (MHCalibrationChargePix&)GetAverageLoGainArea(j);
640 MCalibrationChargePix &pix = (MCalibrationChargePix&)fCam->GetAverageArea(j);
641 MBadPixelsPix &bad = fCam->GetAverageBadArea(j);
642
643 if (histlo.GetSaturated() > fNumLoGainSaturationLimit*histlo.GetHGausHist()->GetEntries())
644 {
645 *fLog << warn << "Saturated Lo Gain histogram in area idx: " << j << endl;
646 histlo.CreateFourierSpectrum();
647 continue;
648 }
649
650 if (pix.IsHiGainSaturation())
651 FinalizeAbsTimes(histlo, pix, bad, fFirstLoGain, fLastLoGain);
652 }
653
654 for (Int_t j=0; j<fAverageHiGainSectors->GetSize(); j++)
655 {
656
657 MHCalibrationChargePix &histhi = (MHCalibrationChargePix&)GetAverageHiGainSector(j);
658 MCalibrationChargePix &pix = (MCalibrationChargePix&)fCam->GetAverageSector(j);
659 MBadPixelsPix &bad = fCam->GetAverageBadSector(j);
660
661 if (histhi.GetSaturated() > fNumHiGainSaturationLimit*histhi.GetHGausHist()->GetEntries())
662 {
663 pix.SetHiGainSaturation();
664 histhi.CreateFourierSpectrum();
665 continue;
666 }
667
668 FinalizeAbsTimes(histhi, pix, bad, fFirstHiGain, fLastHiGain);
669 }
670
671 for (Int_t j=0; j<fAverageLoGainSectors->GetSize(); j++)
672 {
673
674 MHCalibrationChargePix &histlo = (MHCalibrationChargePix&)GetAverageLoGainSector(j);
675 MCalibrationChargePix &pix = (MCalibrationChargePix&)fCam->GetAverageSector(j);
676 MBadPixelsPix &bad = fCam->GetAverageBadSector(j);
677
678 if (histlo.GetSaturated() > fNumLoGainSaturationLimit*histlo.GetHGausHist()->GetEntries())
679 {
680 *fLog << warn << "Saturated Lo Gain histogram in sector: " << j << endl;
681 bad.SetUncalibrated( MBadPixelsPix::kLoGainSaturation );
682 histlo.CreateFourierSpectrum();
683 continue;
684 }
685
686 if (pix.IsHiGainSaturation())
687 FinalizeAbsTimes(histlo, pix, bad, fFirstLoGain, fLastLoGain);
688 }
689
690 //
691 // Perform the fitting for the High Gain (done in MHCalibrationCam)
692 //
693 FitHiGainArrays((MCalibrationCam&)(*fCam),(*fBadPixels),
694 MBadPixelsPix::kHiGainNotFitted,
695 MBadPixelsPix::kHiGainOscillating);
696 //
697 // Perform the fitting for the Low Gain (done in MHCalibrationCam)
698 //
699 FitLoGainArrays((MCalibrationCam&)(*fCam),(*fBadPixels),
700 MBadPixelsPix::kLoGainNotFitted,
701 MBadPixelsPix::kLoGainOscillating);
702
703 return kTRUE;
704}
705
706// --------------------------------------------------------------------------------
707//
708// Fill the absolute time results into MCalibrationChargePix
709//
710// Check absolute time validity:
711// - Mean arrival time is at least fTimeLowerLimit slices from the lower edge
712// - Mean arrival time is at least fUpperLimit slices from the upper edge
713//
714void MHCalibrationChargeCam::FinalizeAbsTimes(MHCalibrationChargePix &hist, MCalibrationChargePix &pix, MBadPixelsPix &bad,
715 Byte_t first, Byte_t last)
716{
717
718 const Float_t mean = hist.GetAbsTimeMean();
719 const Float_t rms = hist.GetAbsTimeRms();
720
721 pix.SetAbsTimeMean ( mean );
722 pix.SetAbsTimeRms ( rms );
723
724 const Float_t lowerlimit = (Float_t)first + fTimeLowerLimit;
725 const Float_t upperlimit = (Float_t)last + fTimeUpperLimit;
726
727 if ( mean < lowerlimit)
728 {
729 *fLog << warn << GetDescriptor()
730 << Form("%s%3.1f%s%2.1f%s%4i",": Mean ArrivalTime: ",mean," smaller than ",fTimeLowerLimit,
731 " FADC slices from lower edge in pixel ",hist.GetPixId()) << endl;
732 bad.SetUncalibrated( MBadPixelsPix::kMeanTimeInFirstBin );
733 }
734
735 if ( mean > upperlimit )
736 {
737 *fLog << warn << GetDescriptor()
738 << Form("%s%3.1f%s%2.1f%s%4i",": Mean ArrivalTime: ",mean," greater than ",fTimeUpperLimit,
739 " FADC slices from upper edge in pixel ",hist.GetPixId()) << endl;
740 bad.SetUncalibrated( MBadPixelsPix::kMeanTimeInLast2Bins );
741 }
742}
743
744// --------------------------------------------------------------------------
745//
746// Sets all pixels to MBadPixelsPix::kUnsuitableRun, if following flags are set:
747// - MBadPixelsPix::kLoGainSaturation
748//
749// Sets all pixels to MBadPixelsPix::kUnreliableRun, if following flags are set:
750// - if MBadPixelsPix::kHiGainNotFitted and !MCalibrationPix::IsHiGainSaturation()
751// - if MBadPixelsPix::kHiGainOscillating and !MCalibrationPix::IsHiGainSaturation()
752// - if MBadPixelsPix::kLoGainNotFitted and MCalibrationPix::IsLoGainSaturation()
753// - if MBadPixelsPix::kLoGainOscillating and MCalibrationPix::IsLoGainSaturation()
754//
755void MHCalibrationChargeCam::FinalizeBadPixels()
756{
757
758 for (Int_t i=0; i<fBadPixels->GetSize(); i++)
759 {
760
761 MBadPixelsPix &bad = (*fBadPixels)[i];
762 MCalibrationPix &pix = (*fCam)[i];
763
764 if (bad.IsUncalibrated( MBadPixelsPix::kHiGainNotFitted ))
765 if (!pix.IsHiGainSaturation())
766 bad.SetUnsuitable( MBadPixelsPix::kUnreliableRun );
767
768 if (bad.IsUncalibrated( MBadPixelsPix::kHiGainOscillating ))
769 bad.SetUnsuitable( MBadPixelsPix::kUnreliableRun );
770
771 if (bad.IsUncalibrated( MBadPixelsPix::kLoGainNotFitted ))
772 if (pix.IsHiGainSaturation())
773 bad.SetUnsuitable( MBadPixelsPix::kUnreliableRun );
774
775 if (bad.IsUncalibrated( MBadPixelsPix::kLoGainOscillating ))
776 if (pix.IsHiGainSaturation())
777 bad.SetUnsuitable( MBadPixelsPix::kUnreliableRun );
778
779 if (bad.IsUncalibrated( MBadPixelsPix::kLoGainSaturation ))
780 bad.SetUnsuitable( MBadPixelsPix::kUnsuitableRun );
781 }
782}
783
784// --------------------------------------------------------------------------
785//
786// Dummy, needed by MCamEvent
787//
788Bool_t MHCalibrationChargeCam::GetPixelContent(Double_t &val, Int_t idx, const MGeomCam &cam, Int_t type) const
789{
790 return kTRUE;
791}
792
793// --------------------------------------------------------------------------
794//
795// Calls MHGausEvents::DrawClone() for pixel idx
796//
797void MHCalibrationChargeCam::DrawPixelContent(Int_t idx) const
798{
799 (*this)[idx].DrawClone();
800}
801
802
803// -----------------------------------------------------------------------------
804//
805// Default draw:
806//
807// Displays the averaged areas, both High Gain and Low Gain
808//
809// Calls the Draw of the fAverageHiGainAreas and fAverageLoGainAreas objects with options
810//
811void MHCalibrationChargeCam::Draw(const Option_t *opt)
812{
813
814 const Int_t nareas = fAverageHiGainAreas->GetEntries();
815 if (nareas == 0)
816 return;
817
818 TString option(opt);
819 option.ToLower();
820
821 if (!option.Contains("datacheck"))
822 {
823 MHCalibrationCam::Draw(opt);
824 return;
825 }
826
827 TVirtualPad *pad = gPad ? gPad : MH::MakeDefCanvas(this);
828 pad->SetBorderMode(0);
829
830 pad->Divide(2,nareas);
831
832 for (Int_t i=0; i<nareas;i++)
833 {
834 pad->cd(2*(i+1)-1);
835 GetAverageHiGainArea(i).Draw("datacheck");
836 pad->cd(2*(i+1));
837 GetAverageLoGainArea(i).Draw("datacheck");
838 }
839}
840
841
842// --------------------------------------------------------------------------
843//
844// Our own clone function is necessary since root 3.01/06 or Mars 0.4
845// I don't know the reason.
846//
847// Creates new MHCalibrationCam
848//
849TObject *MHCalibrationChargeCam::Clone(const char *) const
850{
851
852 const Int_t navhi = fAverageHiGainAreas->GetEntries();
853 const Int_t navlo = fAverageLoGainAreas->GetEntries();
854 const Int_t nsehi = fAverageHiGainSectors->GetEntries();
855 const Int_t nselo = fAverageLoGainSectors->GetEntries();
856
857 //
858 // FIXME, this might be done faster and more elegant, by direct copy.
859 //
860 MHCalibrationChargeCam *cam = new MHCalibrationChargeCam();
861
862 cam->fAverageHiGainAreas->Expand(navhi);
863 cam->fAverageLoGainAreas->Expand(navlo);
864 cam->fAverageHiGainSectors->Expand(nsehi);
865 cam->fAverageLoGainSectors->Expand(nselo);
866
867 cam->fAverageHiGainAreas->Expand(navhi);
868 cam->fAverageLoGainAreas->Expand(navlo);
869 cam->fAverageHiGainSectors->Expand(nsehi);
870 cam->fAverageLoGainSectors->Expand(nselo);
871
872 for (int i=0; i<navhi; i++)
873 {
874 // delete (*cam->fAverageHiGainAreas)[i];
875 (*cam->fAverageHiGainAreas)[i] = (*fAverageHiGainAreas)[i]->Clone();
876 }
877 for (int i=0; i<navlo; i++)
878 {
879 // delete (*cam->fAverageLoGainAreas)[i];
880 (*cam->fAverageLoGainAreas)[i] = (*fAverageLoGainAreas)[i]->Clone();
881 }
882 for (int i=0; i<nsehi; i++)
883 {
884 // delete (*cam->fAverageHiGainSectors)[i];
885 (*cam->fAverageHiGainSectors)[i] = (*fAverageHiGainSectors)[i]->Clone();
886 }
887 for (int i=0; i<nselo; i++)
888 {
889 // delete (*cam->fAverageLoGainSectors)[i];
890 (*cam->fAverageLoGainSectors)[i] = (*fAverageLoGainSectors)[i]->Clone();
891 }
892
893 cam->fAverageAreaNum = fAverageAreaNum;
894 cam->fAverageAreaSat = fAverageAreaSat;
895 cam->fAverageAreaSigma = fAverageAreaSigma;
896 cam->fAverageAreaSigmaVar = fAverageAreaSigmaVar;
897 cam->fAverageAreaRelSigma = fAverageAreaRelSigma;
898 cam->fAverageAreaRelSigmaVar = fAverageAreaRelSigmaVar;
899 cam->fAverageSectorNum = fAverageSectorNum;
900 cam->fRunNumbers = fRunNumbers;
901
902 cam->fPulserFrequency = fPulserFrequency;
903 cam->fAverageNbins = fAverageNbins;
904
905 return cam;
906
907}
908
Note: See TracBrowser for help on using the repository browser.