source: trunk/MagicSoft/Mars/mhcalib/MHCalibrationChargeCam.cc@ 5729

Last change on this file since 5729 was 5687, checked in by gaug, 20 years ago
*** empty log message ***
File size: 48.6 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 MHCalibrationPix-classes
29// MHCalibrationChargeHiGainPix and MHCalibrationChargeLoGainPix for every:
30//
31// - Pixel, stored in the TOrdCollection'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 TOrdCollection'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 TOrdCollection'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// +- MHCalibrationPix::fPickupLimit (default: 5) sigma (see MHCalibrationPix::RepeatFit())
57// In case this does not make the fit valid, the histogram means and RMS's are
58// taken directly (see MHCalibrationPix::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 MHCalibrationPix::fPickupLimit (default: 5) sigmas
64// from the mean are counted as Pickup events (stored in MHCalibrationPix::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 "MHCalibrationChargePix.h"
120#include "MHCalibrationPix.h"
121
122#include "MCalibrationIntensityCam.h"
123#include "MCalibrationChargeCam.h"
124#include "MCalibrationChargePix.h"
125
126#include "MGeomCam.h"
127#include "MGeomPix.h"
128
129#include "MBadPixelsIntensityCam.h"
130#include "MBadPixelsCam.h"
131#include "MBadPixelsPix.h"
132
133#include "MRawEvtData.h"
134#include "MRawRunHeader.h"
135#include "MRawEvtPixelIter.h"
136
137#include "MExtractedSignalCam.h"
138#include "MExtractedSignalPix.h"
139
140#include "MArrayI.h"
141#include "MArrayD.h"
142
143#include <TOrdCollection.h>
144#include <TPad.h>
145#include <TVirtualPad.h>
146#include <TCanvas.h>
147#include <TStyle.h>
148#include <TF1.h>
149#include <TLatex.h>
150#include <TLegend.h>
151#include <TGraph.h>
152
153ClassImp(MHCalibrationChargeCam);
154
155using namespace std;
156
157const Int_t MHCalibrationChargeCam::fgChargeHiGainNbins = 500;
158const Axis_t MHCalibrationChargeCam::fgChargeHiGainFirst = -100.5;
159const Axis_t MHCalibrationChargeCam::fgChargeHiGainLast = 1899.5;
160const Int_t MHCalibrationChargeCam::fgChargeLoGainNbins = 620;
161const Axis_t MHCalibrationChargeCam::fgChargeLoGainFirst = -350.5;
162const Axis_t MHCalibrationChargeCam::fgChargeLoGainLast = 2049.5;
163const Float_t MHCalibrationChargeCam::fgProbLimit = 0.00001;
164const TString MHCalibrationChargeCam::gsHistName = "Charge";
165const TString MHCalibrationChargeCam::gsHistTitle = "Signals";
166const TString MHCalibrationChargeCam::gsHistXTitle = "Signal [FADC counts]";
167const TString MHCalibrationChargeCam::gsHistYTitle = "Nr. events";
168const TString MHCalibrationChargeCam::gsAbsHistName = "AbsTime";
169const TString MHCalibrationChargeCam::gsAbsHistTitle = "Abs. Arr. Times";
170const TString MHCalibrationChargeCam::gsAbsHistXTitle = "Time [FADC slices]";
171const TString MHCalibrationChargeCam::gsAbsHistYTitle = "Nr. events";
172const Float_t MHCalibrationChargeCam::fgNumHiGainSaturationLimit = 0.01;
173const Float_t MHCalibrationChargeCam::fgNumLoGainSaturationLimit = 0.005;
174const Float_t MHCalibrationChargeCam::fgTimeLowerLimit = 1.;
175const Float_t MHCalibrationChargeCam::fgTimeUpperLimit = 2.;
176// 1Led Green, 1 LED blue, 5 LEDs blue, 10 LEDs blue, 10 LEDs UV, CT1, 5Leds Green
177const Float_t MHCalibrationChargeCam::gkHiGainInnerRefLines[7] = { 245., 323. , 1065., 1467., 180., 211. , 533.5};
178const Float_t MHCalibrationChargeCam::gkHiGainOuterRefLines[7] = { 217., 307.5, 932. , 1405., 167., 183.5, 405.5};
179const Float_t MHCalibrationChargeCam::gkLoGainInnerRefLines[7] = { 20.8, 28.0 , 121. , 200.2, 16.5, 13.5 , 41.7 };
180const Float_t MHCalibrationChargeCam::gkLoGainOuterRefLines[7] = { 18.9, 26.0 , 108.3, 198. , 14.0, 11. , 42. };
181// --------------------------------------------------------------------------
182//
183// Default Constructor.
184//
185// Sets:
186// - all pointers to NULL
187//
188// Initializes:
189// - fNumHiGainSaturationLimit to fgNumHiGainSaturationLimit
190// - fNumLoGainSaturationLimit to fgNumLoGainSaturationLimit
191// - fTimeLowerLimit to fgTimeLowerLimit
192// - fTimeUpperLimit to fgTimeUpperLimit
193//
194// - fNbins to fgChargeHiGainNbins
195// - fFirst to fgChargeHiGainFirst
196// - fLast to fgChargeHiGainLast
197//
198// - fLoGainNbins to fgChargeLoGainNbins
199// - fLoGainFirst to fgChargeLoGainFirst
200// - fLoGainLast to fgChargeLoGainLast
201//
202// - fHistName to gsHistName
203// - fHistTitle to gsHistTitle
204// - fHistXTitle to gsHistXTitle
205// - fHistYTitle to gsHistYTitle
206//
207// - fAbsHistName to gsAbsHistName
208// - fAbsHistTitle to gsAbsHistTitle
209// - fAbsHistXTitle to gsAbsHistXTitle
210// - fAbsHistYTitle to gsAbsHistYTitle
211//
212MHCalibrationChargeCam::MHCalibrationChargeCam(const char *name, const char *title)
213 : fRawEvt(NULL)
214{
215
216 fName = name ? name : "MHCalibrationChargeCam";
217 fTitle = title ? title : "Class to fill the calibration histograms ";
218
219 SetNumHiGainSaturationLimit(fgNumHiGainSaturationLimit);
220 SetNumLoGainSaturationLimit(fgNumLoGainSaturationLimit);
221 SetTimeLowerLimit();
222 SetTimeUpperLimit();
223
224 SetNbins(fgChargeHiGainNbins);
225 SetFirst(fgChargeHiGainFirst);
226 SetLast (fgChargeHiGainLast );
227
228 SetProbLimit(fgProbLimit);
229
230 SetLoGainNbins(fgChargeLoGainNbins);
231 SetLoGainFirst(fgChargeLoGainFirst);
232 SetLoGainLast (fgChargeLoGainLast );
233
234 SetHistName (gsHistName .Data());
235 SetHistTitle (gsHistTitle .Data());
236 SetHistXTitle(gsHistXTitle.Data());
237 SetHistYTitle(gsHistYTitle.Data());
238
239 SetAbsHistName (gsAbsHistName .Data());
240 SetAbsHistTitle (gsAbsHistTitle .Data());
241 SetAbsHistXTitle(gsAbsHistXTitle.Data());
242 SetAbsHistYTitle(gsAbsHistYTitle.Data());
243}
244
245// --------------------------------------------------------------------------
246//
247// Creates new MHCalibrationChargeCam only with the averaged areas:
248// the rest has to be retrieved directly, e.g. via:
249// MHCalibrationChargeCam *cam = MParList::FindObject("MHCalibrationChargeCam");
250// - cam->GetAverageSector(5).DrawClone();
251// - (*cam)[100].DrawClone()
252//
253TObject *MHCalibrationChargeCam::Clone(const char *) const
254{
255
256 MHCalibrationChargeCam *cam = new MHCalibrationChargeCam();
257
258 //
259 // Copy the data members
260 //
261 cam->fColor = fColor;
262 cam->fRunNumbers = fRunNumbers;
263 cam->fPulserFrequency = fPulserFrequency;
264 cam->fFlags = fFlags;
265 cam->fNbins = fNbins;
266 cam->fFirst = fFirst;
267 cam->fLast = fLast;
268 cam->fLoGainNbins = fLoGainNbins;
269 cam->fLoGainFirst = fLoGainFirst;
270 cam->fLoGainLast = fLoGainLast;
271
272 //
273 // Copy the MArrays
274 //
275 cam->fAverageAreaRelSigma = fAverageAreaRelSigma;
276 cam->fAverageAreaRelSigmaVar = fAverageAreaRelSigmaVar;
277 cam->fAverageAreaSat = fAverageAreaSat;
278 cam->fAverageAreaSigma = fAverageAreaSigma;
279 cam->fAverageAreaSigmaVar = fAverageAreaSigmaVar;
280 cam->fAverageAreaNum = fAverageAreaNum;
281 cam->fAverageSectorNum = fAverageSectorNum;
282
283 if (!IsAverageing())
284 return cam;
285
286 const Int_t navhi = fAverageHiGainAreas->GetSize();
287
288 for (int i=0; i<navhi; i++)
289 cam->fAverageHiGainAreas->AddAt(GetAverageHiGainArea(i).Clone(),i);
290
291 if (IsLoGain())
292 {
293
294 const Int_t navlo = fAverageLoGainAreas->GetSize();
295 for (int i=0; i<navlo; i++)
296 cam->fAverageLoGainAreas->AddAt(GetAverageLoGainArea(i).Clone(),i);
297
298 }
299
300 return cam;
301}
302
303// --------------------------------------------------------------------------
304//
305// Gets the pointers to:
306// - MRawEvtData
307//
308Bool_t MHCalibrationChargeCam::SetupHists(const MParList *pList)
309{
310
311 fRawEvt = (MRawEvtData*)pList->FindObject("MRawEvtData");
312 if (!fRawEvt)
313 {
314 *fLog << err << dbginf << "MRawEvtData not found... aborting." << endl;
315 return kFALSE;
316 }
317
318 return kTRUE;
319}
320
321// --------------------------------------------------------------------------
322//
323// Gets or creates the pointers to:
324// - MExtractedSignalCam
325// - MCalibrationChargeCam or MCalibrationIntensityChargeCam
326// - MBadPixelsCam
327//
328// Initializes the number of used FADC slices from MExtractedSignalCam
329// into MCalibrationChargeCam and test for changes in that variable
330//
331// Calls:
332// - InitHiGainArrays()
333// - InitLoGainArrays()
334//
335// Sets:
336// - fSumhiarea to nareas
337// - fSumloarea to nareas
338// - fTimehiarea to nareas
339// - fTimeloarea to nareas
340// - fSumhisector to nsectors
341// - fSumlosector to nsectors
342// - fTimehisector to nsectors
343// - fTimelosector to nsectors
344// - fSathiarea to nareas
345// - fSatloarea to nareas
346// - fSathisector to nsectors
347// - fSatlosector to nsectors
348//
349Bool_t MHCalibrationChargeCam::ReInitHists(MParList *pList)
350{
351
352 MExtractedSignalCam *signal =
353 (MExtractedSignalCam*)pList->FindObject(AddSerialNumber("MExtractedSignalCam"));
354 if (!signal)
355 {
356 *fLog << err << "MExtractedSignalCam not found... abort." << endl;
357 return kFALSE;
358 }
359
360 fIntensCam = (MCalibrationIntensityCam*)pList->FindObject(AddSerialNumber("MCalibrationIntensityChargeCam"));
361 if (fIntensCam)
362 *fLog << inf << "Found MCalibrationIntensityChargeCam ... " << endl;
363 else
364 {
365 fCam = (MCalibrationCam*)pList->FindObject(AddSerialNumber("MCalibrationChargeCam"));
366 if (!fCam)
367 {
368 fCam = (MCalibrationCam*)pList->FindCreateObj(AddSerialNumber("MCalibrationChargeCam"));
369 if (!fCam)
370 {
371 *fLog << err << "Cannot find nor create MCalibrationChargeCam ... abort." << endl;
372 return kFALSE;
373 }
374 fCam->Init(*fGeom);
375 }
376 }
377
378 fFirstHiGain = signal->GetFirstUsedSliceHiGain();
379 fLastHiGain = signal->GetLastUsedSliceHiGain();
380 fFirstLoGain = signal->GetFirstUsedSliceLoGain();
381 fLastLoGain = signal->GetLastUsedSliceLoGain();
382
383 /*
384 const Float_t numhigain = signal->GetNumUsedHiGainFADCSlices();
385 const Float_t numlogain = signal->GetNumUsedLoGainFADCSlices();
386
387 if (fCam)
388 {
389 if (fCam->GetNumHiGainFADCSlices() == 0.)
390 fCam->SetNumHiGainFADCSlices ( numhigain );
391 else if (fCam->GetNumHiGainFADCSlices() != numhigain)
392 {
393 *fLog << err << GetDescriptor()
394 << ": Number of High Gain FADC extraction slices has changed, abort..." << endl;
395 return kFALSE;
396 }
397
398 if (fCam->GetNumLoGainFADCSlices() == 0.)
399 fCam->SetNumLoGainFADCSlices ( numlogain );
400 else if (fCam->GetNumLoGainFADCSlices() != numlogain)
401 {
402 *fLog << err << GetDescriptor()
403 << ": Number of Low Gain FADC extraction slices has changes, abort..." << endl;
404 return kFALSE;
405 }
406 }
407 */
408
409 const Int_t npixels = fGeom->GetNumPixels();
410 const Int_t nsectors = fGeom->GetNumSectors();
411 const Int_t nareas = fGeom->GetNumAreas();
412
413 //
414 // In case of the intense blue, double the range
415 //
416 if (fGeom->InheritsFrom("MGeomCamMagic"))
417 if ( fColor == MCalibrationCam::kBLUE)
418 SetLoGainLast(2.*fLoGainLast - fLoGainFirst);
419
420 InitHiGainArrays(npixels,nareas,nsectors);
421 InitLoGainArrays(npixels,nareas,nsectors);
422
423 fSumhiarea .Set(nareas);
424 fSumloarea .Set(nareas);
425 fTimehiarea .Set(nareas);
426 fTimeloarea .Set(nareas);
427 fSumhisector.Set(nsectors);
428 fSumlosector.Set(nsectors);
429 fTimehisector.Set(nsectors);
430 fTimelosector.Set(nsectors);
431
432 fSathiarea .Set(nareas);
433 fSatloarea .Set(nareas);
434 fSathisector.Set(nsectors);
435 fSatlosector.Set(nsectors);
436
437 return kTRUE;
438}
439
440// --------------------------------------------------------------------------
441//
442// Retrieve:
443// - fRunHeader->GetNumSamplesHiGain();
444//
445// Initializes the High Gain Arrays:
446//
447// - For every entry in the expanded arrays:
448// * Initialize an MHCalibrationChargePix
449// * Set Binning from fNbins, fFirst and fLast
450// * Set Binning of Abs Times histogram from fAbsNbins, fAbsFirst and fAbsLast
451// * Set Histgram names and titles from fHistName and fHistTitle
452// * Set Abs Times Histgram names and titles from fAbsHistName and fAbsHistTitle
453// * Set X-axis and Y-axis titles from fHistXTitle and fHistYTitle
454// * Set X-axis and Y-axis titles of Abs Times Histogram from fAbsHistXTitle and fAbsHistYTitle
455// * Call InitHists
456//
457//
458void MHCalibrationChargeCam::InitHiGainArrays(const Int_t npixels, const Int_t nareas, const Int_t nsectors)
459{
460
461 TH1F *h;
462
463 const Int_t higainsamples = fRunHeader->GetNumSamplesHiGain();
464
465 if (fHiGainArray->GetSize()==0)
466 {
467 for (Int_t i=0; i<npixels; i++)
468 {
469 fHiGainArray->AddAt(new MHCalibrationChargePix(Form("%sHiGainPix%04d",fHistName.Data(),i),
470 Form("%s High Gain Pixel%04d",fHistTitle.Data(),i)),i);
471
472 MHCalibrationChargePix &pix = (MHCalibrationChargePix&)(*this)[i];
473
474 pix.SetNbins(fNbins);
475 pix.SetFirst(fFirst);
476 pix.SetLast (fLast);
477
478 pix.SetProbLimit(fProbLimit);
479
480 pix.SetAbsTimeNbins(higainsamples);
481 pix.SetAbsTimeFirst(-0.5);
482 pix.SetAbsTimeLast(higainsamples-0.5);
483
484 InitHists(pix,fIntensBad ? (*fIntensBad)[i] : (*fBadPixels)[i],i);
485
486 h = pix.GetHAbsTime();
487
488 h->SetName (Form("H%sHiGainPix%04d",fAbsHistName.Data(),i));
489 h->SetTitle(Form("%s High Gain Pixel %04d",fAbsHistTitle.Data(),i));
490 h->SetXTitle(fAbsHistXTitle.Data());
491 h->SetYTitle(fAbsHistYTitle.Data());
492 }
493 }
494
495
496 if (fAverageHiGainAreas->GetSize()==0)
497 {
498 for (Int_t j=0; j<nareas; j++)
499 {
500 fAverageHiGainAreas->AddAt(new MHCalibrationChargePix(Form("%sHiGainArea%d",fHistName.Data(),j),
501 Form("%s High Gain Area Idx %d",fHistTitle.Data(),j)),j);
502
503 MHCalibrationChargePix &pix = (MHCalibrationChargePix&)GetAverageHiGainArea(j);
504
505 pix.SetNbins(fNbins*(Int_t)TMath::Sqrt((Float_t)npixels/nareas));
506 pix.SetFirst(fFirst);
507 pix.SetLast (fLast);
508
509 pix.SetAbsTimeNbins(higainsamples);
510 pix.SetAbsTimeFirst(-0.5);
511 pix.SetAbsTimeLast(higainsamples-0.5);
512
513 InitHists(pix,fIntensCam ? fIntensCam->GetAverageBadArea(j) : fCam->GetAverageBadArea(j),j);
514
515 h = pix.GetHAbsTime();
516
517 h->SetName (Form("H%sHiGainArea%d",fAbsHistName.Data(),j));
518 h->SetTitle(Form("%s%s%d",fAbsHistTitle.Data(),
519 " averaged on event-by-event basis High Gain Area Idx ",j));
520 h->SetXTitle(fAbsHistXTitle.Data());
521 h->SetYTitle(fAbsHistYTitle.Data());
522 }
523 }
524
525 if (fAverageHiGainSectors->GetSize()==0)
526 {
527 for (Int_t j=0; j<nsectors; j++)
528 {
529 fAverageHiGainSectors->AddAt(new MHCalibrationChargePix(Form("%sHiGainSector%02d",fHistName.Data(),j),
530 Form("%s High Gain Sector %02d",fHistTitle.Data(),j)),j);
531
532 MHCalibrationChargePix &pix = (MHCalibrationChargePix&)GetAverageHiGainSector(j);
533
534 pix.SetNbins(fNbins*(Int_t)TMath::Sqrt((Float_t)npixels/nareas));
535 pix.SetFirst(fFirst);
536 pix.SetLast (fLast);
537
538 pix.SetAbsTimeNbins(higainsamples);
539 pix.SetAbsTimeFirst(-0.5);
540 pix.SetAbsTimeLast(higainsamples-0.5);
541
542 InitHists(pix,fIntensCam ? fIntensCam->GetAverageBadSector(j) : fCam->GetAverageBadSector(j),j);
543
544 h = pix.GetHAbsTime();
545
546 h->SetName (Form("H%sHiGainSector%02d",fAbsHistName.Data(),j));
547 h->SetTitle(Form("%s%s%02d",fAbsHistTitle.Data(),
548 " averaged on event-by-event basis High Gain Area Sector ",j));
549 h->SetXTitle(fAbsHistXTitle.Data());
550 h->SetYTitle(fAbsHistYTitle.Data());
551 }
552 }
553}
554
555//--------------------------------------------------------------------------------------
556//
557// Return, if IsLoGain() is kFALSE
558//
559// Retrieve:
560// - fRunHeader->GetNumSamplesHiGain();
561//
562// Initializes the Low Gain Arrays:
563//
564// - For every entry in the expanded arrays:
565// * Initialize an MHCalibrationChargePix
566// * Set Binning from fNbins, fFirst and fLast
567// * Set Binning of Abs Times histogram from fAbsNbins, fAbsFirst and fAbsLast
568// * Set Histgram names and titles from fHistName and fHistTitle
569// * Set Abs Times Histgram names and titles from fAbsHistName and fAbsHistTitle
570// * Set X-axis and Y-axis titles from fHistXTitle and fHistYTitle
571// * Set X-axis and Y-axis titles of Abs Times Histogram from fAbsHistXTitle and fAbsHistYTitle
572// * Call InitHists
573//
574void MHCalibrationChargeCam::InitLoGainArrays(const Int_t npixels, const Int_t nareas, const Int_t nsectors)
575{
576
577 if (!IsLoGain())
578 return;
579
580 const Int_t logainsamples = fRunHeader->GetNumSamplesLoGain();
581
582 TH1F *h;
583
584 if (fLoGainArray->GetSize()==0 )
585 {
586 for (Int_t i=0; i<npixels; i++)
587 {
588 fLoGainArray->AddAt(new MHCalibrationChargePix(Form("%sLoGainPix%04d",fHistName.Data(),i),
589 Form("%s Low Gain Pixel %04d",fHistTitle.Data(),i)),i);
590
591 MHCalibrationChargePix &pix = (MHCalibrationChargePix&)(*this)(i);
592
593 pix.SetNbins(fLoGainNbins);
594 pix.SetFirst(fLoGainFirst);
595 pix.SetLast (fLoGainLast);
596
597 pix.SetProbLimit(fProbLimit);
598
599 pix.SetAbsTimeNbins(logainsamples);
600 pix.SetAbsTimeFirst(-0.5);
601 pix.SetAbsTimeLast(logainsamples-0.5);
602
603 InitHists(pix,fIntensBad ? (*fIntensBad)[i] : (*fBadPixels)[i],i);
604
605 h = pix.GetHAbsTime();
606
607 h->SetName (Form("H%sLoGainPix%04d",fAbsHistName.Data(),i));
608 h->SetTitle(Form("%s Low Gain Pixel %04d",fAbsHistTitle.Data(),i));
609 h->SetXTitle(fAbsHistXTitle.Data());
610 h->SetYTitle(fAbsHistYTitle.Data());
611 }
612 }
613
614 if (fAverageLoGainAreas->GetSize()==0)
615 {
616 for (Int_t j=0; j<nareas; j++)
617 {
618 fAverageLoGainAreas->AddAt(new MHCalibrationChargePix(Form("%sLoGainArea%d",fHistName.Data(),j),
619 Form("%s Low Gain Area Idx %d",fHistTitle.Data(),j)),j);
620
621 MHCalibrationChargePix &pix = (MHCalibrationChargePix&)GetAverageLoGainArea(j);
622
623 pix.SetNbins(fLoGainNbins*(Int_t)TMath::Sqrt((Float_t)npixels/nareas));
624 pix.SetFirst(fLoGainFirst);
625 pix.SetLast (fLoGainLast);
626
627 pix.SetAbsTimeNbins(logainsamples);
628 pix.SetAbsTimeFirst(-0.5);
629 pix.SetAbsTimeLast(logainsamples-0.5);
630
631 InitHists(pix,fIntensCam ? fIntensCam->GetAverageBadArea(j) : fCam->GetAverageBadArea(j),j);
632
633 h = pix.GetHAbsTime();
634
635 h->SetName (Form("H%sLoGainArea%02d",fAbsHistName.Data(),j));
636 h->SetTitle(Form("%s%s%02d",fAbsHistTitle.Data(),
637 " averaged on event-by-event basis Low Gain Area Idx ",j));
638 h->SetXTitle(fAbsHistXTitle.Data());
639 h->SetYTitle(fAbsHistYTitle.Data());
640 }
641 }
642
643
644 if (fAverageLoGainSectors->GetSize()==0 && IsLoGain())
645 {
646 for (Int_t j=0; j<nsectors; j++)
647 {
648 fAverageLoGainSectors->AddAt(new MHCalibrationChargePix(Form("%sLoGainSector%02d",fHistName.Data(),j),
649 Form("%s Low Gain Sector %02d",fHistTitle.Data(),j)),j);
650
651 MHCalibrationChargePix &pix = (MHCalibrationChargePix&)GetAverageLoGainSector(j);
652
653 pix.SetNbins(fLoGainNbins*(Int_t)TMath::Sqrt((Float_t)npixels/nareas));
654 pix.SetFirst(fLoGainFirst);
655 pix.SetLast (fLoGainLast);
656
657 pix.SetAbsTimeNbins(logainsamples);
658 pix.SetAbsTimeFirst(-0.5);
659 pix.SetAbsTimeLast(logainsamples-0.5);
660
661 InitHists(pix,fIntensCam ? fIntensCam->GetAverageBadSector(j) : fCam->GetAverageBadSector(j),j);
662
663 h = pix.GetHAbsTime();
664
665 h->SetName (Form("H%sLoGainSector%02d",fAbsHistName.Data(),j));
666 h->SetTitle(Form("%s%s%02d",fAbsHistTitle.Data(),
667 " averaged on event-by-event basis Low Gain Area Sector ",j));
668 h->SetXTitle(fAbsHistXTitle.Data());
669 h->SetYTitle(fAbsHistYTitle.Data());
670 }
671 }
672}
673
674
675// --------------------------------------------------------------------------
676//
677// Retrieves from MExtractedSignalCam:
678// - first used LoGain FADC slice
679//
680// Retrieves from MGeomCam:
681// - number of pixels
682// - number of pixel areas
683// - number of sectors
684//
685// For all TOrdCollection's (including the averaged ones), the following steps are performed:
686//
687// 1) Fill Charges histograms (MHGausEvents::FillHistAndArray()) with:
688// - MExtractedSignalPix::GetExtractedSignalHiGain();
689// - MExtractedSignalPix::GetExtractedSignalLoGain();
690//
691// 2) Set number of saturated slices (MHCalibrationChargePix::AddSaturated()) with:
692// - MExtractedSignalPix::GetNumHiGainSaturated();
693// - MExtractedSignalPix::GetNumLoGainSaturated();
694//
695// 3) Fill AbsTime histograms (MHCalibrationChargePix::FillAbsTime()) with:
696// - MRawEvtPixelIter::GetIdxMaxHiGainSample();
697// - MRawEvtPixelIter::GetIdxMaxLoGainSample(first slice);
698//
699Bool_t MHCalibrationChargeCam::FillHists(const MParContainer *par, const Stat_t w)
700{
701
702 MExtractedSignalCam *signal = (MExtractedSignalCam*)par;
703 if (!signal)
704 {
705 *fLog << err << "No argument in MExtractedSignalCam::Fill... abort." << endl;
706 return kFALSE;
707 }
708
709 const UInt_t npixels = fGeom->GetNumPixels();
710 const UInt_t nareas = fGeom->GetNumAreas();
711 const UInt_t nsectors = fGeom->GetNumSectors();
712 const UInt_t lofirst = signal->GetFirstUsedSliceLoGain();
713
714 fSumhiarea .Reset();
715 fSumloarea .Reset();
716 fTimehiarea .Reset();
717 fTimeloarea .Reset();
718 fSumhisector.Reset();
719 fSumlosector.Reset();
720 fTimehisector.Reset();
721 fTimelosector.Reset();
722
723 fSathiarea .Reset();
724 fSatloarea .Reset();
725 fSathisector.Reset();
726 fSatlosector.Reset();
727
728 for (UInt_t i=0; i<npixels; i++)
729 {
730
731 MHCalibrationChargePix &histhi = (MHCalibrationChargePix&)(*this)[i];
732
733 if (histhi.IsExcluded())
734 continue;
735
736 const MExtractedSignalPix &pix = (*signal)[i];
737
738 const Float_t sumhi = pix.GetExtractedSignalHiGain();
739 const Int_t sathi = (Int_t)pix.GetNumHiGainSaturated();
740
741 if (IsOscillations())
742 histhi.FillHistAndArray(sumhi);
743 else
744 histhi.FillHist(sumhi);
745
746 histhi.AddSaturated(sathi);
747
748 const Int_t aidx = (*fGeom)[i].GetAidx();
749 const Int_t sector = (*fGeom)[i].GetSector();
750
751 fSumhiarea[aidx] += sumhi;
752 fSathiarea[aidx] += sathi;
753
754 fSumhisector[sector] += sumhi;
755 fSathisector[sector] += sathi;
756
757 if (IsLoGain())
758 {
759 MHCalibrationChargePix &histlo = (MHCalibrationChargePix&)(*this)(i);
760 const Float_t sumlo = pix.GetExtractedSignalLoGain();
761 const Int_t satlo = (Int_t)pix.GetNumLoGainSaturated();
762
763 if (IsOscillations())
764 histlo.FillHistAndArray(sumlo);
765 else
766 histlo.FillHist(sumlo);
767
768 histlo.AddSaturated(satlo);
769
770 fSumloarea[aidx] += sumlo;
771 fSatloarea[aidx] += satlo;
772 fSumlosector[sector] += sumlo;
773 fSatlosector[sector] += satlo;
774 }
775
776 }
777
778 MRawEvtPixelIter pixel(fRawEvt);
779 while (pixel.Next())
780 {
781
782 const UInt_t pixid = pixel.GetPixelId();
783
784 MHCalibrationChargePix &histhi = (MHCalibrationChargePix&)(*this)[pixid];
785
786 if (histhi.IsExcluded())
787 continue;
788
789 const Float_t timehi = (Float_t)pixel.GetIdxMaxHiGainSample();
790
791 histhi.FillAbsTime(timehi);
792
793 const Int_t aidx = (*fGeom)[pixid].GetAidx();
794 const Int_t sector = (*fGeom)[pixid].GetSector();
795
796 fTimehiarea [aidx] += timehi;
797 fTimehisector[sector] += timehi;
798
799 if (IsLoGain())
800 {
801 MHCalibrationChargePix &histlo = (MHCalibrationChargePix&)(*this)(pixid);
802
803 const Float_t timelo = (Float_t)pixel.GetIdxMaxLoGainSample(lofirst);
804 histlo.FillAbsTime(timelo);
805
806 fTimeloarea[aidx] += timelo;
807 fTimelosector[sector] += timelo;
808 }
809 }
810
811 for (UInt_t j=0; j<nareas; j++)
812 {
813
814 const Int_t npix = fAverageAreaNum[j];
815
816 if (npix == 0)
817 continue;
818
819 MHCalibrationChargePix &hipix = (MHCalibrationChargePix&)GetAverageHiGainArea(j);
820
821
822 if (IsOscillations())
823 hipix.FillHistAndArray(fSumhiarea [j]/npix);
824 else
825 hipix.FillHist(fSumhiarea[j]/npix);
826
827 hipix.AddSaturated ((Float_t)fSathiarea [j]/npix > 0.5 ? 1 : 0);
828 hipix.FillAbsTime (fTimehiarea[j]/npix);
829
830 if (IsLoGain())
831 {
832 MHCalibrationChargePix &lopix = (MHCalibrationChargePix&)GetAverageLoGainArea(j);
833 if (IsOscillations())
834 lopix.FillHistAndArray(fSumloarea [j]/npix);
835 else
836 lopix.FillHist(fSumloarea [j]/npix);
837 lopix.AddSaturated ((Float_t)fSatloarea [j]/npix > 0.5 ? 1 : 0);
838 lopix.FillAbsTime (fTimeloarea[j]/npix);
839 }
840 }
841
842 for (UInt_t j=0; j<nsectors; j++)
843 {
844
845 const Int_t npix = fAverageSectorNum[j];
846
847 if (npix == 0)
848 continue;
849
850 MHCalibrationChargePix &hipix = (MHCalibrationChargePix&)GetAverageHiGainSector(j);
851
852 if (IsOscillations())
853 hipix.FillHistAndArray(fSumhisector [j]/npix);
854 else
855 hipix.FillHist(fSumhisector [j]/npix);
856
857 hipix.AddSaturated ((Float_t)fSathisector[j]/npix > 0.5 ? 1 : 0);
858 hipix.FillAbsTime (fTimehisector[j]/npix);
859
860 if (IsLoGain())
861 {
862 MHCalibrationChargePix &lopix = (MHCalibrationChargePix&)GetAverageLoGainSector(j);
863
864 if (IsOscillations())
865 lopix.FillHistAndArray(fSumlosector [j]/npix);
866 else
867 lopix.FillHist(fSumlosector [j]/npix);
868
869 lopix.AddSaturated ((Float_t)fSatlosector[j]/npix > 0.5 ? 1 : 0);
870 lopix.FillAbsTime (fTimelosector[j]/npix);
871 }
872 }
873
874 return kTRUE;
875}
876
877// --------------------------------------------------------------------------
878//
879// For all TOrdCollection's (including the averaged ones), the following steps are performed:
880//
881// 1) Returns if the pixel is excluded.
882// 2) Tests saturation. In case yes, set the flag: MCalibrationPix::SetHiGainSaturation()
883// or the flag: MBadPixelsPix::SetUncalibrated( MBadPixelsPix::kLoGainSaturated )
884// 3) Store the absolute arrival times in the MCalibrationChargePix's. If flag
885// MCalibrationPix::IsHiGainSaturation() is set, the Low-Gain arrival times are stored,
886// otherwise the Hi-Gain ones.
887// 4) Calls to MHCalibrationCam::FitHiGainArrays() and MCalibrationCam::FitLoGainArrays()
888// with the flags:
889// - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::kHiGainNotFitted )
890// - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::kLoGainNotFitted )
891// - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::kHiGainOscillating )
892// - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::kLoGainOscillating )
893//
894Bool_t MHCalibrationChargeCam::FinalizeHists()
895{
896
897 *fLog << endl;
898
899 TH1F *h = NULL;
900
901 MCalibrationCam *chargecam = fIntensCam ? fIntensCam->GetCam() : fCam;
902 MBadPixelsCam *badcam = fIntensBad ? fIntensBad->GetCam() : fBadPixels;
903
904 for (Int_t i=0; i<fHiGainArray->GetSize(); i++)
905 {
906
907 MHCalibrationChargePix &histhi = (MHCalibrationChargePix&)(*this)[i];
908
909 if (histhi.IsExcluded())
910 continue;
911
912 MCalibrationChargePix &pix = (MCalibrationChargePix&)(*chargecam)[i] ;
913
914 if (histhi.GetSaturated() > fNumHiGainSaturationLimit*histhi.GetHGausHist()->GetEntries())
915 {
916 pix.SetHiGainSaturation();
917 if (IsOscillations())
918 histhi.CreateFourierSpectrum();
919 continue;
920 }
921
922 MBadPixelsPix &bad = (*badcam)[i];
923
924 h = histhi.GetHGausHist();
925
926 Stat_t overflow = h->GetBinContent(h->GetNbinsX()+1);
927 if (overflow > 0.1)
928 {
929 *fLog << warn
930 << "HiGain Hist-overflow " << overflow
931 << " times in pix: " << i << " (w/o saturation!) " << endl;
932 bad.SetUncalibrated( MBadPixelsPix::kHiGainOverFlow );
933 }
934
935 overflow = h->GetBinContent(0);
936 if (overflow > 0.1)
937 {
938 *fLog << warn
939 << "HiGain Hist-underflow " << overflow
940 << " times in pix: " << i << " (w/o saturation!) " << endl;
941 bad.SetUncalibrated( MBadPixelsPix::kHiGainOverFlow );
942 }
943
944 FinalizeAbsTimes(histhi, pix, bad, fFirstHiGain, fLastHiGain);
945 }
946
947 if (IsLoGain())
948 for (Int_t i=0; i<fLoGainArray->GetSize(); i++)
949 {
950
951 MHCalibrationChargePix &histlo = (MHCalibrationChargePix&)(*this)(i);
952 MBadPixelsPix &bad = (*badcam)[i];
953
954 if (histlo.IsExcluded())
955 continue;
956
957 if (histlo.GetSaturated() > fNumLoGainSaturationLimit*histlo.GetHGausHist()->GetEntries())
958 {
959 *fLog << warn << "Saturated Lo Gain histogram in pixel: " << i << endl;
960 bad.SetUncalibrated( MBadPixelsPix::kLoGainSaturation );
961 if (IsOscillations())
962 histlo.CreateFourierSpectrum();
963 continue;
964 }
965
966 h = histlo.GetHGausHist();
967
968 Stat_t overflow = h->GetBinContent(h->GetNbinsX()+1);
969 if (overflow > 0.1)
970 {
971 *fLog << warn
972 << "LoGain Hist-overflow " << overflow
973 << " times in pix: " << i << " (w/o saturation!) " << endl;
974 bad.SetUncalibrated( MBadPixelsPix::kLoGainOverFlow );
975 }
976
977 overflow = h->GetBinContent(0);
978 if (overflow > 0.1)
979 {
980 *fLog << warn
981 << "LoGain Hist-underflow " << overflow
982 << " times in pix: " << i << " (w/o saturation!) " << endl;
983 bad.SetUncalibrated( MBadPixelsPix::kLoGainOverFlow );
984 }
985
986 MCalibrationChargePix &pix = (MCalibrationChargePix&)(*chargecam)[i] ;
987
988 if (pix.IsHiGainSaturation())
989 FinalizeAbsTimes(histlo, pix, bad, fFirstLoGain, fLastLoGain);
990 }
991
992 for (Int_t j=0; j<fAverageHiGainAreas->GetSize(); j++)
993 {
994
995 MHCalibrationChargePix &histhi = (MHCalibrationChargePix&)GetAverageHiGainArea(j);
996 MCalibrationChargePix &pix = (MCalibrationChargePix&)chargecam->GetAverageArea(j);
997
998 if (histhi.GetSaturated() > fNumHiGainSaturationLimit*histhi.GetHGausHist()->GetEntries())
999 {
1000 pix.SetHiGainSaturation();
1001 if (IsOscillations())
1002 histhi.CreateFourierSpectrum();
1003 continue;
1004 }
1005
1006 MBadPixelsPix &bad = chargecam->GetAverageBadArea(j);
1007 FinalizeAbsTimes(histhi, pix, bad, fFirstHiGain, fLastHiGain);
1008 }
1009
1010 if (IsLoGain())
1011 for (Int_t j=0; j<fAverageLoGainAreas->GetSize(); j++)
1012 {
1013
1014 MHCalibrationChargePix &histlo = (MHCalibrationChargePix&)GetAverageLoGainArea(j);
1015
1016 if (histlo.GetSaturated() > fNumLoGainSaturationLimit*histlo.GetHGausHist()->GetEntries())
1017 {
1018 *fLog << warn << "Saturated Lo Gain histogram in area idx: " << j << endl;
1019 histlo.CreateFourierSpectrum();
1020 continue;
1021 }
1022
1023 MCalibrationChargePix &pix = (MCalibrationChargePix&)chargecam->GetAverageArea(j);
1024
1025 if (pix.IsHiGainSaturation())
1026 {
1027 MBadPixelsPix &bad = chargecam->GetAverageBadArea(j);
1028 FinalizeAbsTimes(histlo, pix, bad, fFirstLoGain, fLastLoGain);
1029 }
1030
1031 }
1032
1033 for (Int_t j=0; j<fAverageHiGainSectors->GetSize(); j++)
1034 {
1035
1036 MHCalibrationChargePix &histhi = (MHCalibrationChargePix&)GetAverageHiGainSector(j);
1037 MCalibrationChargePix &pix = (MCalibrationChargePix&)chargecam->GetAverageSector(j);
1038
1039 if (histhi.GetSaturated() > fNumHiGainSaturationLimit*histhi.GetHGausHist()->GetEntries())
1040 {
1041 pix.SetHiGainSaturation();
1042 if (IsOscillations())
1043 histhi.CreateFourierSpectrum();
1044 continue;
1045 }
1046
1047 MBadPixelsPix &bad = chargecam->GetAverageBadSector(j);
1048 FinalizeAbsTimes(histhi, pix, bad, fFirstHiGain, fLastHiGain);
1049 }
1050
1051 if (IsLoGain())
1052 for (Int_t j=0; j<fAverageLoGainSectors->GetSize(); j++)
1053 {
1054
1055 MHCalibrationChargePix &histlo = (MHCalibrationChargePix&)GetAverageLoGainSector(j);
1056 MBadPixelsPix &bad = chargecam->GetAverageBadSector(j);
1057
1058 if (histlo.GetSaturated() > fNumLoGainSaturationLimit*histlo.GetHGausHist()->GetEntries())
1059 {
1060 *fLog << warn << "Saturated Lo Gain histogram in sector: " << j << endl;
1061 bad.SetUncalibrated( MBadPixelsPix::kLoGainSaturation );
1062 if (IsOscillations())
1063 histlo.CreateFourierSpectrum();
1064 continue;
1065 }
1066
1067 MCalibrationChargePix &pix = (MCalibrationChargePix&)chargecam->GetAverageSector(j);
1068
1069 if (pix.IsHiGainSaturation())
1070 FinalizeAbsTimes(histlo, pix, bad, fFirstLoGain, fLastLoGain);
1071 }
1072
1073 //
1074 // Perform the fitting for the High Gain (done in MHCalibrationCam)
1075 //
1076 FitHiGainArrays(*chargecam, *badcam,
1077 MBadPixelsPix::kHiGainNotFitted,
1078 MBadPixelsPix::kHiGainOscillating);
1079
1080 //
1081 // Perform the fitting for the Low Gain (done in MHCalibrationCam)
1082 //
1083 if (IsLoGain())
1084 FitLoGainArrays(*chargecam, *badcam,
1085 MBadPixelsPix::kLoGainNotFitted,
1086 MBadPixelsPix::kLoGainOscillating);
1087
1088
1089 return kTRUE;
1090}
1091
1092// --------------------------------------------------------------------------------
1093//
1094// Fill the absolute time results into MCalibrationChargePix
1095//
1096// Check absolute time validity:
1097// - Mean arrival time is at least fTimeLowerLimit slices from the lower edge
1098// - Mean arrival time is at least fUpperLimit slices from the upper edge
1099//
1100void MHCalibrationChargeCam::FinalizeAbsTimes(MHCalibrationChargePix &hist, MCalibrationChargePix &pix, MBadPixelsPix &bad,
1101 Byte_t first, Byte_t last)
1102{
1103
1104 const Float_t mean = hist.GetAbsTimeMean();
1105 const Float_t rms = hist.GetAbsTimeRms();
1106
1107 pix.SetAbsTimeMean ( mean );
1108 pix.SetAbsTimeRms ( rms );
1109
1110 const Float_t lowerlimit = (Float_t)first + fTimeLowerLimit;
1111 const Float_t upperlimit = (Float_t)last - fTimeUpperLimit;
1112
1113 if ( mean < lowerlimit)
1114 {
1115 *fLog << warn
1116 << Form("Mean ArrivalTime: %3.1f < %2.1f slices from lower edge: %2i in pixel %s",
1117 mean,fTimeLowerLimit,(Int_t)first,hist.GetName()) << endl;
1118 bad.SetUncalibrated( MBadPixelsPix::kMeanTimeInFirstBin );
1119 }
1120
1121 if ( mean > upperlimit )
1122 {
1123 *fLog << warn
1124 << Form("Mean ArrivalTime: %3.1f > %2.1f slices from upper edge: %2i in pixel %s",
1125 mean,fTimeUpperLimit,(Int_t)last,hist.GetName()) << endl;
1126 bad.SetUncalibrated( MBadPixelsPix::kMeanTimeInLast2Bins );
1127 }
1128}
1129
1130// --------------------------------------------------------------------------
1131//
1132// Sets all pixels to MBadPixelsPix::kUnsuitableRun, if following flags are set:
1133// - MBadPixelsPix::kLoGainSaturation
1134//
1135// Sets all pixels to MBadPixelsPix::kUnreliableRun, if following flags are set:
1136// - if MBadPixelsPix::kHiGainNotFitted and !MCalibrationPix::IsHiGainSaturation()
1137// - if MBadPixelsPix::kHiGainOscillating and !MCalibrationPix::IsHiGainSaturation()
1138// - if MBadPixelsPix::kLoGainNotFitted and MCalibrationPix::IsLoGainSaturation()
1139// - if MBadPixelsPix::kLoGainOscillating and MCalibrationPix::IsLoGainSaturation()
1140//
1141void MHCalibrationChargeCam::FinalizeBadPixels()
1142{
1143
1144 MBadPixelsCam *badcam = fIntensBad ? fIntensBad->GetCam() : fBadPixels;
1145 MCalibrationCam *chargecam = fIntensCam ? fIntensCam->GetCam() : fCam;
1146
1147 for (Int_t i=0; i<badcam->GetSize(); i++)
1148 {
1149
1150 MBadPixelsPix &bad = (*badcam)[i];
1151 MCalibrationPix &pix = (*chargecam)[i];
1152
1153 if (bad.IsUncalibrated( MBadPixelsPix::kHiGainNotFitted ))
1154 if (!pix.IsHiGainSaturation())
1155 bad.SetUnsuitable( MBadPixelsPix::kUnreliableRun );
1156
1157 if (bad.IsUncalibrated( MBadPixelsPix::kLoGainNotFitted ))
1158 if (pix.IsHiGainSaturation())
1159 bad.SetUnsuitable( MBadPixelsPix::kUnreliableRun );
1160
1161 if (bad.IsUncalibrated( MBadPixelsPix::kLoGainSaturation ))
1162 bad.SetUnsuitable( MBadPixelsPix::kUnsuitableRun );
1163
1164 if (IsOscillations())
1165 {
1166 if (bad.IsUncalibrated( MBadPixelsPix::kHiGainOscillating ))
1167 bad.SetUnsuitable( MBadPixelsPix::kUnreliableRun );
1168
1169 if (bad.IsUncalibrated( MBadPixelsPix::kLoGainOscillating ))
1170 if (pix.IsHiGainSaturation())
1171 bad.SetUnsuitable( MBadPixelsPix::kUnreliableRun );
1172 }
1173 }
1174
1175
1176
1177}
1178
1179// --------------------------------------------------------------------------
1180//
1181// Calls MHCalibrationPix::DrawClone() for pixel idx
1182//
1183void MHCalibrationChargeCam::DrawPixelContent(Int_t idx) const
1184{
1185 (*this)[idx].DrawClone();
1186}
1187
1188
1189// -----------------------------------------------------------------------------
1190//
1191// Default draw:
1192//
1193// Displays the averaged areas, both High Gain and Low Gain
1194//
1195// Calls the Draw of the fAverageHiGainAreas and fAverageLoGainAreas objects with options
1196//
1197void MHCalibrationChargeCam::Draw(const Option_t *opt)
1198{
1199
1200 const Int_t nareas = fAverageHiGainAreas->GetSize();
1201 if (nareas == 0)
1202 return;
1203
1204 TString option(opt);
1205 option.ToLower();
1206
1207 if (!option.Contains("datacheck"))
1208 {
1209 MHCalibrationCam::Draw(opt);
1210 return;
1211 }
1212
1213 //
1214 // From here on , the datacheck - Draw
1215 //
1216 TVirtualPad *pad = gPad ? gPad : MH::MakeDefCanvas(this);
1217 pad->SetBorderMode(0);
1218 pad->Divide(1,nareas);
1219
1220 //
1221 // Loop over inner and outer pixels
1222 //
1223 for (Int_t i=0; i<nareas;i++)
1224 {
1225
1226 pad->cd(i+1);
1227
1228 MHCalibrationChargePix &hipix = (MHCalibrationChargePix&)GetAverageHiGainArea(i);
1229 //
1230 // Ask for Hi-Gain saturation
1231 //
1232 if (hipix.GetSaturated() > fNumHiGainSaturationLimit*hipix.GetHGausHist()->GetEntries() && IsLoGain())
1233 {
1234 MHCalibrationChargePix &lopix = (MHCalibrationChargePix&)GetAverageLoGainArea(i);
1235 DrawDataCheckPixel(lopix,i ? gkLoGainOuterRefLines : gkLoGainInnerRefLines);
1236 }
1237 else
1238 DrawDataCheckPixel(hipix,i ? gkHiGainOuterRefLines : gkHiGainInnerRefLines);
1239
1240 }
1241}
1242
1243// -----------------------------------------------------------------------------
1244//
1245// Draw the average pixel for the datacheck:
1246//
1247// Displays the averaged areas, both High Gain and Low Gain
1248//
1249// Calls the Draw of the fAverageHiGainAreas and fAverageLoGainAreas objects with options
1250//
1251void MHCalibrationChargeCam::DrawDataCheckPixel(MHCalibrationChargePix &pix, const Float_t refline[])
1252{
1253
1254 TVirtualPad *newpad = gPad;
1255 newpad->Divide(1,2);
1256 newpad->cd(1);
1257
1258 gPad->SetTicks();
1259 if (!pix.IsEmpty() && !pix.IsOnlyOverflow() && !pix.IsOnlyUnderflow())
1260 gPad->SetLogy();
1261
1262 TH1F *hist = pix.GetHGausHist();
1263
1264 TH1F *null = new TH1F("Null",hist->GetTitle(),100,
1265 pix.GetFirst() > 0. ? pix.GetFirst() : 0.,
1266 pix.GetLast() > pix.GetFirst()
1267 ? ( pix.GetLast() > 450.
1268 ? ( fColor == MCalibrationCam::kBLUE ? 800. : 450. )
1269 : pix.GetLast() )
1270 : pix.GetFirst()*2.);
1271
1272 null->SetMaximum(1.1*hist->GetMaximum());
1273 null->SetDirectory(NULL);
1274 null->SetBit(kCanDelete);
1275 null->SetStats(kFALSE);
1276 //
1277 // set the labels bigger
1278 //
1279 TAxis *xaxe = null->GetXaxis();
1280 TAxis *yaxe = null->GetYaxis();
1281 xaxe->CenterTitle();
1282 yaxe->CenterTitle();
1283 xaxe->SetTitleSize(0.07);
1284 yaxe->SetTitleSize(0.07);
1285 xaxe->SetTitleOffset(0.7);
1286 yaxe->SetTitleOffset(0.55);
1287 xaxe->SetLabelSize(0.06);
1288 yaxe->SetLabelSize(0.06);
1289 xaxe->SetTitle(hist->GetXaxis()->GetTitle());
1290 yaxe->SetTitle(hist->GetYaxis()->GetTitle());
1291
1292 null->Draw();
1293 hist->Draw("same");
1294
1295 gStyle->SetOptFit();
1296
1297 TF1 *fit = pix.GetFGausFit();
1298
1299 if (fit)
1300 {
1301 switch ( fColor )
1302 {
1303 case MCalibrationCam::kGREEN:
1304 fit->SetLineColor(kGreen);
1305 break;
1306 case MCalibrationCam::kBLUE:
1307 fit->SetLineColor(kBlue);
1308 break;
1309 case MCalibrationCam::kUV:
1310 fit->SetLineColor(106);
1311 break;
1312 case MCalibrationCam::kCT1:
1313 fit->SetLineColor(006);
1314 break;
1315 default:
1316 fit->SetLineColor(kRed);
1317 }
1318 fit->Draw("same");
1319 }
1320
1321 DisplayRefLines(null,refline);
1322
1323 newpad->cd(2);
1324 gPad->SetTicks();
1325
1326 TH1F *null2 = new TH1F("Null2",hist->GetTitle(),100,0.,pix.GetEvents()->GetSize()/pix.GetEventFrequency());
1327
1328 null2->SetMinimum(pix.GetMean()-10.*pix.GetSigma());
1329 null2->SetMaximum(pix.GetMean()+10.*pix.GetSigma());
1330 null2->SetDirectory(NULL);
1331 null2->SetBit(kCanDelete);
1332 null2->SetStats(kFALSE);
1333 //
1334 // set the labels bigger
1335 //
1336 TAxis *xaxe2 = null2->GetXaxis();
1337 TAxis *yaxe2 = null2->GetYaxis();
1338 xaxe2->CenterTitle();
1339 yaxe2->CenterTitle();
1340 xaxe2->SetTitleSize(0.07);
1341 yaxe2->SetTitleSize(0.07);
1342 xaxe2->SetTitleOffset(0.7);
1343 yaxe2->SetTitleOffset(0.55);
1344 xaxe2->SetLabelSize(0.06);
1345 yaxe2->SetLabelSize(0.06);
1346
1347 pix.CreateGraphEvents();
1348 TGraph *gr = pix.GetGraphEvents();
1349
1350 xaxe2->SetTitle(gr->GetXaxis()->GetTitle());
1351 yaxe2->SetTitle(gr->GetYaxis()->GetTitle());
1352
1353 null2->Draw();
1354
1355 pix.DrawEvents("same");
1356
1357 // newpad->cd(3);
1358 // pix.DrawPowerSpectrum(*newpad,4);
1359
1360 return;
1361
1362}
1363
1364
1365void MHCalibrationChargeCam::DisplayRefLines(const TH1F *hist, const Float_t refline[]) const
1366{
1367
1368 TGraph *green1 = new TGraph(2);
1369 green1->SetPoint(0,refline[0],0.1);
1370 green1->SetPoint(1,refline[0],hist->GetMaximum());
1371 green1->SetBit(kCanDelete);
1372 green1->SetLineColor(kGreen);
1373 green1->SetLineStyle(2);
1374 green1->SetLineWidth(3);
1375 green1->Draw("L");
1376
1377 TGraph *green5 = new TGraph(2);
1378 green5->SetPoint(0,refline[6],0.1);
1379 green5->SetPoint(1,refline[6],hist->GetMaximum());
1380 green5->SetBit(kCanDelete);
1381 green5->SetLineColor(8);
1382 green5->SetLineStyle(2);
1383 green5->SetLineWidth(3);
1384 green5->Draw("L");
1385
1386 TGraph *blue1 = new TGraph(2);
1387 blue1->SetPoint(0,refline[1],0.1);
1388 blue1->SetPoint(1,refline[1],hist->GetMaximum());
1389 blue1->SetBit(kCanDelete);
1390 blue1->SetLineColor(227);
1391 blue1->SetLineStyle(2);
1392 blue1->SetLineWidth(3);
1393 blue1->Draw("L");
1394
1395 TGraph *blue5 = new TGraph(2);
1396 blue5->SetPoint(0,refline[2],0.1);
1397 blue5->SetPoint(1,refline[2],hist->GetMaximum());
1398 blue5->SetBit(kCanDelete);
1399 blue5->SetLineColor(68);
1400 blue5->SetLineStyle(2);
1401 blue5->SetLineWidth(3);
1402 blue5->Draw("L");
1403
1404 TGraph *blue10 = new TGraph(2);
1405 blue10->SetPoint(0,refline[3],0.1);
1406 blue10->SetPoint(1,refline[3],hist->GetMaximum());
1407 blue10->SetBit(kCanDelete);
1408 blue10->SetLineColor(4);
1409 blue10->SetLineStyle(2);
1410 blue10->SetLineWidth(3);
1411 blue10->Draw("L");
1412
1413 TGraph *uv10 = new TGraph(2);
1414 uv10->SetPoint(0,refline[4],0.1);
1415 uv10->SetPoint(1,refline[4],hist->GetMaximum());
1416 uv10->SetBit(kCanDelete);
1417 uv10->SetLineColor(106);
1418 uv10->SetLineStyle(2);
1419 uv10->SetLineWidth(3);
1420 uv10->Draw("L");
1421
1422 TGraph *ct1 = new TGraph(2);
1423 ct1->SetPoint(0,refline[5],0.1);
1424 ct1->SetPoint(1,refline[5],hist->GetMaximum());
1425 ct1->SetBit(kCanDelete);
1426 ct1->SetLineColor(6);
1427 ct1->SetLineStyle(2);
1428 ct1->SetLineWidth(3);
1429 ct1->Draw("L");
1430
1431 TLegend *leg = new TLegend(0.8,0.35,0.99,0.99);
1432 leg->SetBit(kCanDelete);
1433 leg->AddEntry(green1,"1 Led GREEN","l");
1434 leg->AddEntry(green5,"5 Leds GREEN","l");
1435 leg->AddEntry(blue1,"1 Led BLUE","l");
1436 leg->AddEntry(blue5,"5 Leds BLUE","l");
1437 leg->AddEntry(blue10,"10 Leds BLUE","l");
1438 leg->AddEntry(uv10,"10 Leds UV","l");
1439 leg->AddEntry(ct1,"CT1-Pulser","l");
1440
1441 leg->Draw();
1442}
1443
1444Int_t MHCalibrationChargeCam::ReadEnv(const TEnv &env, TString prefix, Bool_t print)
1445{
1446
1447 Bool_t rc = kFALSE;
1448
1449 if (MHCalibrationCam::ReadEnv(env,prefix,print))
1450 rc = kTRUE;
1451
1452 if (IsEnvDefined(env, prefix, "LoGainNbins", print))
1453 {
1454 SetLoGainNbins(GetEnvValue(env, prefix, "LoGainNbins", fLoGainNbins));
1455 rc = kTRUE;
1456 }
1457
1458 if (IsEnvDefined(env, prefix, "LoGainFirst", print))
1459 {
1460 SetLoGainFirst(GetEnvValue(env, prefix, "LoGainFirst", fLoGainFirst));
1461 rc = kTRUE;
1462 }
1463
1464 if (IsEnvDefined(env, prefix, "LoGainLast", print))
1465 {
1466 SetLoGainLast(GetEnvValue(env, prefix, "LoGainLast", fLoGainLast));
1467 rc = kTRUE;
1468 }
1469
1470 if (IsEnvDefined(env, prefix, "TimeLowerLimit", print))
1471 {
1472 SetTimeLowerLimit(GetEnvValue(env, prefix, "TimeLowerLimit", fTimeLowerLimit));
1473 rc = kTRUE;
1474 }
1475
1476 if (IsEnvDefined(env, prefix, "TimeUpperLimit", print))
1477 {
1478 SetTimeUpperLimit(GetEnvValue(env, prefix, "TimeUpperLimit", fTimeUpperLimit));
1479 rc = kTRUE;
1480 }
1481
1482
1483 return rc;
1484}
Note: See TracBrowser for help on using the repository browser.