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

Last change on this file since 8439 was 8417, checked in by tbretz, 18 years ago
*** empty log message ***
File size: 48.4 KB
Line 
1/* ======================================================================== *\
2! $Name: not supported by cvs2svn $:$Id: MHCalibrationChargeCam.cc,v 1.56 2007-04-18 15:33:56 tbretz Exp $
3! --------------------------------------------------------------------------
4!
5! *
6! * This file is part of MARS, the MAGIC Analysis and Reconstruction
7! * Software. It is distributed to you in the hope that it can be a useful
8! * and timesaving tool in analysing Data of imaging Cerenkov telescopes.
9! * It is distributed WITHOUT ANY WARRANTY.
10! *
11! * Permission to use, copy, modify and distribute this software and its
12! * documentation for any purpose is hereby granted without fee,
13! * provided that the above copyright notice appear in all copies and
14! * that both that copyright notice and this permission notice appear
15! * in supporting documentation. It is provided "as is" without express
16! * or implied warranty.
17! *
18!
19!
20! Author(s): Markus Gaug 02/2004 <mailto:markus@ifae.es>
21! Author(s): Thomas Bretz <mailto:tbretz@astro.uni-wuerzburg.de>
22!
23! Copyright: MAGIC Software Development, 2000-2007
24!
25!
26\* ======================================================================== */
27/////////////////////////////////////////////////////////////////////////////
28//
29// MHCalibrationChargeCam
30//
31// Fills the extracted signals of MExtractedSignalCam into the MHCalibrationPix-classes
32// MHCalibrationChargeHiGainPix and MHCalibrationChargeLoGainPix for every:
33//
34// - Pixel, stored in the TOrdCollection's MHCalibrationCam::fHiGainArray and
35// MHCalibrationCam::fLoGainArray
36//
37// - Average pixel per AREA index (e.g. inner and outer for the MAGIC camera),
38// stored in the TOrdCollection's MHCalibrationCam::fAverageHiGainAreas and
39// MHCalibrationCam::fAverageLoGainAreas
40//
41// - Average pixel per camera SECTOR (e.g. sectors 1-6 for the MAGIC camera),
42// stored in the TOrdCollection's MHCalibrationCam::fAverageHiGainSectors and
43// MHCalibrationCam::fAverageLoGainSectors
44//
45// Every signal is taken from MExtractedSignalCam and filled into a histogram and
46// an array, in order to perform a Fourier analysis (see MHGausEvents).
47// The signals are moreover averaged on an event-by-event basis and written into
48// the corresponding average pixels.
49//
50// Additionally, the (FADC slice) position of the maximum is stored in an Absolute
51// Arrival Time histogram. This histogram serves for a rough cross-check if the
52// signal does not lie at or outside the edges of the extraction window.
53//
54// The Charge histograms are fitted to a Gaussian, mean and sigma with its errors
55// and the fit probability are extracted. If none of these values are NaN's and
56// if the probability is bigger than MHGausEvents::fProbLimit (default: 0.5%),
57// the fit is declared valid.
58// Otherwise, the fit is repeated within ranges of the previous mean
59// +- MHCalibrationPix::fPickupLimit (default: 5) sigma (see MHCalibrationPix::RepeatFit())
60// In case this does not make the fit valid, the histogram means and RMS's are
61// taken directly (see MHCalibrationPix::BypassFit()) and the following flags are set:
62// - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::kHiGainNotFitted ) or
63// - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::kLoGainNotFitted ) and
64// - MBadPixelsPix::SetUnsuitable( MBadPixelsPix::kUnreliableRun )
65//
66// Outliers of more than MHCalibrationPix::fPickupLimit (default: 5) sigmas
67// from the mean are counted as Pickup events (stored in MHCalibrationPix::fPickup)
68//
69// If more than fNumHiGainSaturationLimit (default: 15%) of the overall FADC
70// slices show saturation, the following flag is set:
71// - MCalibrationChargePix::SetHiGainSaturation();
72// In that case, the calibration constants are derived from the low-gain results.
73//
74// If more than fNumLoGainSaturationLimit (default: 1%) of the overall
75// low-gain FADC slices saturate, the following flags are set:
76// - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::kLoGainSaturation ) and
77// - MBadPixelsPix::SetUnsuitable( MBadPixelsPix::kUnsuitableRun )
78//
79// The class also fills arrays with the signal vs. event number, creates a fourier
80// spectrum and investigates if the projected fourier components follow an exponential
81// distribution. In case that the probability of the exponential fit is less than
82// MHGausEvents::fProbLimit (default: 0.5%), the following flags are set:
83// - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::kHiGainOscillating ) or
84// - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::kLoGainOscillating ) and
85// - MBadPixelsPix::SetUnsuitable( MBadPixelsPix::kUnreliableRun )
86//
87// This same procedure is performed for the average pixels.
88//
89// The following results are written into MCalibrationChargeCam:
90//
91// - MCalibrationPix::SetHiGainSaturation()
92// - MCalibrationPix::SetHiGainMean()
93// - MCalibrationPix::SetHiGainMeanErr()
94// - MCalibrationPix::SetHiGainSigma()
95// - MCalibrationPix::SetHiGainSigmaErr()
96// - MCalibrationPix::SetHiGainProb()
97// - MCalibrationPix::SetHiGainNumPickup()
98//
99// - MCalibrationPix::SetLoGainMean()
100// - MCalibrationPix::SetLoGainMeanErr()
101// - MCalibrationPix::SetLoGainSigma()
102// - MCalibrationPix::SetLoGainSigmaErr()
103// - MCalibrationPix::SetLoGainProb()
104// - MCalibrationPix::SetLoGainNumPickup()
105//
106// - MCalibrationChargePix::SetAbsTimeMean()
107// - MCalibrationChargePix::SetAbsTimeRms()
108//
109// For all averaged areas, the fitted sigma is multiplied with the square root of
110// the number involved pixels in order to be able to compare it to the average of
111// sigmas in the camera.
112//
113// Class Version 2:
114// + Float_t fNumLoGainBlackoutLimit; // Rel. amount blackout logain events until pixel is declared unsuitable
115//
116// Class Version 3:
117// - Byte_t fFirstHiGain; // First used slice High Gain
118// - Byte_t fLastHiGain; // Last used slice High Gain
119// - Byte_t fFirstLoGain; // First used slice Low Gain
120// - Byte_t fLastLoGain; // Last used slice Low Gain
121//
122/////////////////////////////////////////////////////////////////////////////
123#include "MHCalibrationChargeCam.h"
124
125#include <TOrdCollection.h>
126#include <TPad.h>
127#include <TVirtualPad.h>
128#include <TCanvas.h>
129#include <TStyle.h>
130#include <TF1.h>
131#include <TLine.h>
132#include <TLatex.h>
133#include <TLegend.h>
134#include <TGraph.h>
135#include <TEnv.h>
136
137#include "MLog.h"
138#include "MLogManip.h"
139
140#include "MParList.h"
141
142#include "MHCalibrationChargePix.h"
143#include "MHCalibrationPix.h"
144
145#include "MCalibrationIntensityCam.h"
146#include "MCalibrationChargeCam.h"
147#include "MCalibrationChargePix.h"
148
149#include "MGeomCam.h"
150#include "MGeomPix.h"
151
152#include "MBadPixelsCam.h"
153#include "MBadPixelsPix.h"
154
155#include "MRawRunHeader.h"
156#include "MRawEvtPixelIter.h"
157
158#include "MPedestalSubtractedEvt.h"
159
160#include "MExtractedSignalCam.h"
161#include "MExtractedSignalPix.h"
162
163ClassImp(MHCalibrationChargeCam);
164
165using namespace std;
166
167const Int_t MHCalibrationChargeCam::fgChargeHiGainNbins = 500;
168const Axis_t MHCalibrationChargeCam::fgChargeHiGainFirst = -98.;
169const Axis_t MHCalibrationChargeCam::fgChargeHiGainLast = 1902.;
170const Int_t MHCalibrationChargeCam::fgChargeLoGainNbins = 500;
171const Axis_t MHCalibrationChargeCam::fgChargeLoGainFirst = -99.;
172const Axis_t MHCalibrationChargeCam::fgChargeLoGainLast = 901.;
173const Float_t MHCalibrationChargeCam::fgProbLimit = 0.00000001;
174const TString MHCalibrationChargeCam::gsHistName = "Charge";
175const TString MHCalibrationChargeCam::gsHistTitle = "Signals";
176const TString MHCalibrationChargeCam::gsHistXTitle = "Signal [FADC counts]";
177const TString MHCalibrationChargeCam::gsHistYTitle = "Nr. events";
178const TString MHCalibrationChargeCam::gsAbsHistName = "AbsTime";
179const TString MHCalibrationChargeCam::gsAbsHistTitle = "Abs. Arr. Times";
180const TString MHCalibrationChargeCam::gsAbsHistXTitle = "Time [FADC slices]";
181const TString MHCalibrationChargeCam::gsAbsHistYTitle = "Nr. events";
182const Float_t MHCalibrationChargeCam::fgNumHiGainSaturationLimit = 0.15;
183const Float_t MHCalibrationChargeCam::fgNumLoGainSaturationLimit = 0.005;
184const Float_t MHCalibrationChargeCam::fgNumLoGainBlackoutLimit = 0.05;
185const Float_t MHCalibrationChargeCam::fgLoGainBlackoutLimit = 3.5;
186const Float_t MHCalibrationChargeCam::fgLoGainPickupLimit = 3.5;
187const Float_t MHCalibrationChargeCam::fgTimeLowerLimit = 1.;
188const Float_t MHCalibrationChargeCam::fgTimeUpperLimit = 3.;
189const TString MHCalibrationChargeCam::fgReferenceFile = "mjobs/calibrationref.rc";
190
191// --------------------------------------------------------------------------
192//
193// Default Constructor.
194//
195// Sets:
196// - all pointers to NULL
197//
198// Initializes:
199// - fNumHiGainSaturationLimit to fgNumHiGainSaturationLimit
200// - fNumLoGainSaturationLimit to fgNumLoGainSaturationLimit
201// - fTimeLowerLimit to fgTimeLowerLimit
202// - fTimeUpperLimit to fgTimeUpperLimit
203// - fNumLoGainBlackoutLimit to fgNumLoGainBlackoutLimit
204//
205// - fNbins to fgChargeHiGainNbins
206// - fFirst to fgChargeHiGainFirst
207// - fLast to fgChargeHiGainLast
208//
209// - fLoGainNbins to fgChargeLoGainNbins
210// - fLoGainFirst to fgChargeLoGainFirst
211// - fLoGainLast to fgChargeLoGainLast
212//
213// - fHistName to gsHistName
214// - fHistTitle to gsHistTitle
215// - fHistXTitle to gsHistXTitle
216// - fHistYTitle to gsHistYTitle
217//
218// - fAbsHistName to gsAbsHistName
219// - fAbsHistTitle to gsAbsHistTitle
220// - fAbsHistXTitle to gsAbsHistXTitle
221// - fAbsHistYTitle to gsAbsHistYTitle
222//
223MHCalibrationChargeCam::MHCalibrationChargeCam(const char *name, const char *title)
224 : fRawEvt(NULL), fSignal(NULL)
225{
226
227 fName = name ? name : "MHCalibrationChargeCam";
228 fTitle = title ? title : "Class to fill the calibration histograms ";
229
230 SetNumHiGainSaturationLimit(fgNumHiGainSaturationLimit);
231 SetNumLoGainSaturationLimit(fgNumLoGainSaturationLimit);
232
233 SetNumLoGainBlackoutLimit (fgNumLoGainBlackoutLimit);
234
235 SetTimeLowerLimit();
236 SetTimeUpperLimit();
237
238 SetBinning(fgChargeHiGainNbins, fgChargeHiGainFirst, fgChargeHiGainLast );
239 SetBinningLoGain(fgChargeLoGainNbins, fgChargeLoGainFirst, fgChargeLoGainLast );
240
241 SetProbLimit(fgProbLimit);
242
243 SetHistName (gsHistName .Data());
244 SetHistTitle (gsHistTitle .Data());
245 SetHistXTitle(gsHistXTitle.Data());
246 SetHistYTitle(gsHistYTitle.Data());
247
248 SetAbsHistName (gsAbsHistName .Data());
249 SetAbsHistTitle (gsAbsHistTitle .Data());
250 SetAbsHistXTitle(gsAbsHistXTitle.Data());
251 SetAbsHistYTitle(gsAbsHistYTitle.Data());
252
253 SetReferenceFile();
254
255 fInnerRefCharge = 278.;
256 fOuterRefCharge = 282.;
257}
258
259// --------------------------------------------------------------------------
260//
261// Creates new MHCalibrationChargeCam only with the averaged areas:
262// the rest has to be retrieved directly, e.g. via:
263// MHCalibrationChargeCam *cam = MParList::FindObject("MHCalibrationChargeCam");
264// - cam->GetAverageSector(5).DrawClone();
265// - (*cam)[100].DrawClone()
266//
267TObject *MHCalibrationChargeCam::Clone(const char *) const
268{
269
270 MHCalibrationChargeCam *cam = new MHCalibrationChargeCam();
271
272 //
273 // Copy the data members
274 //
275 cam->fColor = fColor;
276 cam->fRunNumbers = fRunNumbers;
277 cam->fPulserFrequency = fPulserFrequency;
278 cam->fFlags = fFlags;
279 cam->fNbins = fNbins;
280 cam->fFirst = fFirst;
281 cam->fLast = fLast;
282 cam->fLoGainNbins = fLoGainNbins;
283 cam->fLoGainFirst = fLoGainFirst;
284 cam->fLoGainLast = fLoGainLast;
285 cam->fReferenceFile = fReferenceFile;
286 cam->fInnerRefCharge = fInnerRefCharge;
287 cam->fOuterRefCharge = fOuterRefCharge;
288
289 //
290 // Copy the MArrays
291 //
292 cam->fAverageAreaRelSigma = fAverageAreaRelSigma;
293 cam->fAverageAreaRelSigmaVar = fAverageAreaRelSigmaVar;
294 cam->fAverageAreaSat = fAverageAreaSat;
295 cam->fAverageAreaSigma = fAverageAreaSigma;
296 cam->fAverageAreaSigmaVar = fAverageAreaSigmaVar;
297 cam->fAverageAreaNum = fAverageAreaNum;
298 cam->fAverageSectorNum = fAverageSectorNum;
299
300 if (!IsAverageing())
301 return cam;
302
303 const Int_t navhi = fAverageHiGainAreas->GetSize();
304
305 for (int i=0; i<navhi; i++)
306 cam->fAverageHiGainAreas->AddAt(GetAverageHiGainArea(i).Clone(),i);
307
308 if (IsLoGain())
309 {
310
311 const Int_t navlo = fAverageLoGainAreas->GetSize();
312 for (int i=0; i<navlo; i++)
313 cam->fAverageLoGainAreas->AddAt(GetAverageLoGainArea(i).Clone(),i);
314
315 }
316
317 return cam;
318}
319
320// --------------------------------------------------------------------------
321//
322// Gets the pointers to:
323// - MRawEvtData
324//
325Bool_t MHCalibrationChargeCam::SetupHists(const MParList *pList)
326{
327 fRawEvt = (MPedestalSubtractedEvt*)pList->FindObject("MPedestalSubtractedEvt");
328 if (!fRawEvt)
329 {
330 *fLog << err << dbginf << "MPedestalSubtractedEvt not found... aborting." << endl;
331 return kFALSE;
332 }
333
334 return kTRUE;
335}
336
337// --------------------------------------------------------------------------
338//
339// Gets or creates the pointers to:
340// - MExtractedSignalCam
341// - MCalibrationChargeCam or MCalibrationIntensityChargeCam
342// - MBadPixelsCam
343//
344// Initializes the number of used FADC slices from MExtractedSignalCam
345// into MCalibrationChargeCam and test for changes in that variable
346//
347// Calls:
348// - InitHiGainArrays()
349// - InitLoGainArrays()
350//
351// Sets:
352// - fSumhiarea to nareas
353// - fSumloarea to nareas
354// - fTimehiarea to nareas
355// - fTimeloarea to nareas
356// - fSumhisector to nsectors
357// - fSumlosector to nsectors
358// - fTimehisector to nsectors
359// - fTimelosector to nsectors
360// - fSathiarea to nareas
361// - fSatloarea to nareas
362// - fSathisector to nsectors
363// - fSatlosector to nsectors
364//
365Bool_t MHCalibrationChargeCam::ReInitHists(MParList *pList)
366{
367
368 fSignal = (MExtractedSignalCam*)pList->FindObject(AddSerialNumber("MExtractedSignalCam"));
369 if (!fSignal)
370 {
371 *fLog << err << "MExtractedSignalCam not found... abort." << endl;
372 return kFALSE;
373 }
374
375 if (!InitCams(pList,"Charge"))
376 return kFALSE;
377
378 const Int_t npixels = fGeom->GetNumPixels();
379 const Int_t nsectors = fGeom->GetNumSectors();
380 const Int_t nareas = fGeom->GetNumAreas();
381
382 //
383 // In case of the intense blue, double the range
384 //
385 if (fGeom->InheritsFrom("MGeomCamMagic"))
386 if ( fColor == MCalibrationCam::kBLUE)
387 SetBinning(fLoGainNbins, fLoGainFirst, 2.*fLoGainLast - fLoGainFirst);
388
389 InitHiGainArrays(npixels,nareas,nsectors);
390 InitLoGainArrays(npixels,nareas,nsectors);
391
392 fSumhiarea .Set(nareas);
393 fSumloarea .Set(nareas);
394 fTimehiarea .Set(nareas);
395 fTimeloarea .Set(nareas);
396 fSumhisector.Set(nsectors);
397 fSumlosector.Set(nsectors);
398 fTimehisector.Set(nsectors);
399 fTimelosector.Set(nsectors);
400
401 fSathiarea .Set(nareas);
402 fSatloarea .Set(nareas);
403 fSathisector.Set(nsectors);
404 fSatlosector.Set(nsectors);
405
406 return kTRUE;
407}
408
409// --------------------------------------------------------------------------
410//
411// Retrieve:
412// - fRunHeader->GetNumSamplesHiGain();
413//
414// Initializes the High Gain Arrays:
415//
416// - For every entry in the expanded arrays:
417// * Initialize an MHCalibrationChargePix
418// * Set Binning from fNbins, fFirst and fLast
419// * Set Binning of Abs Times histogram from fAbsNbins, fAbsFirst and fAbsLast
420// * Set Histgram names and titles from fHistName and fHistTitle
421// * Set Abs Times Histgram names and titles from fAbsHistName and fAbsHistTitle
422// * Set X-axis and Y-axis titles from fHistXTitle and fHistYTitle
423// * Set X-axis and Y-axis titles of Abs Times Histogram from fAbsHistXTitle and fAbsHistYTitle
424// * Call InitHists
425//
426//
427void MHCalibrationChargeCam::InitHiGainArrays(const Int_t npixels, const Int_t nareas, const Int_t nsectors)
428{
429
430 TH1F *h;
431
432 const Int_t higainsamples = fRunHeader->GetNumSamplesHiGain();
433
434 if (fHiGainArray->GetSize()==0)
435 {
436 for (Int_t i=0; i<npixels; i++)
437 {
438 fHiGainArray->AddAt(new MHCalibrationChargePix(Form("%sHiGainPix%04d",fHistName.Data(),i),
439 Form("%s High Gain Pixel%04d",fHistTitle.Data(),i)),i);
440
441 MHCalibrationChargePix &pix = (MHCalibrationChargePix&)(*this)[i];
442
443 pix.SetBinning(fNbins, fFirst, fLast);
444 pix.SetProbLimit(fProbLimit);
445
446 pix.SetBinningAbsTime(higainsamples, -0.5, higainsamples-0.5);
447
448 InitHists(pix,(*fBadPixels)[i], i);
449
450 h = pix.GetHAbsTime();
451
452 h->SetName (Form("H%sHiGainPix%04d",fAbsHistName.Data(),i));
453 h->SetTitle(Form("%s High Gain Pixel %04d",fAbsHistTitle.Data(),i));
454 h->SetXTitle(fAbsHistXTitle.Data());
455 h->SetYTitle(fAbsHistYTitle.Data());
456 }
457 }
458
459
460 if (fAverageHiGainAreas->GetSize()==0)
461 {
462 for (Int_t j=0; j<nareas; j++)
463 {
464 fAverageHiGainAreas->AddAt(new MHCalibrationChargePix(Form("%sHiGainArea%d",fHistName.Data(),j),
465 Form("%s High Gain Area Idx %d",fHistTitle.Data(),j)),j);
466
467 MHCalibrationChargePix &pix = (MHCalibrationChargePix&)GetAverageHiGainArea(j);
468
469 pix.SetBinning(fNbins*(Int_t)TMath::Sqrt((Float_t)npixels/nareas),
470 fFirst, fLast);
471
472 pix.SetBinningAbsTime(higainsamples, -0.5, higainsamples-0.5);
473
474 InitHists(pix,fIntensCam ? fIntensCam->GetAverageBadArea(j) : fCam->GetAverageBadArea(j),j);
475
476 h = pix.GetHAbsTime();
477
478 h->SetName (Form("H%sHiGainArea%d",fAbsHistName.Data(),j));
479 h->SetTitle(Form("%s%s%d",fAbsHistTitle.Data(),
480 " averaged on event-by-event basis High Gain Area Idx ",j));
481 h->SetXTitle(fAbsHistXTitle.Data());
482 h->SetYTitle(fAbsHistYTitle.Data());
483 }
484 }
485
486 if (fAverageHiGainSectors->GetSize()==0)
487 {
488 for (Int_t j=0; j<nsectors; j++)
489 {
490 fAverageHiGainSectors->AddAt(new MHCalibrationChargePix(Form("%sHiGainSector%02d",fHistName.Data(),j),
491 Form("%s High Gain Sector %02d",fHistTitle.Data(),j)),j);
492
493 MHCalibrationChargePix &pix = (MHCalibrationChargePix&)GetAverageHiGainSector(j);
494
495 pix.SetBinning(fNbins*(Int_t)TMath::Sqrt((Float_t)npixels/nareas),
496 fFirst, fLast);
497
498 pix.SetBinningAbsTime(higainsamples, -0.5, higainsamples-0.5);
499
500 InitHists(pix,fIntensCam ? fIntensCam->GetAverageBadSector(j) : fCam->GetAverageBadSector(j),j);
501
502 h = pix.GetHAbsTime();
503
504 h->SetName (Form("H%sHiGainSector%02d",fAbsHistName.Data(),j));
505 h->SetTitle(Form("%s%s%02d",fAbsHistTitle.Data(),
506 " averaged on event-by-event basis High Gain Area Sector ",j));
507 h->SetXTitle(fAbsHistXTitle.Data());
508 h->SetYTitle(fAbsHistYTitle.Data());
509 }
510 }
511}
512
513//--------------------------------------------------------------------------------------
514//
515// Return, if IsLoGain() is kFALSE
516//
517// Retrieve:
518// - fRunHeader->GetNumSamplesHiGain();
519//
520// Initializes the Low Gain Arrays:
521//
522// - For every entry in the expanded arrays:
523// * Initialize an MHCalibrationChargePix
524// * Set Binning from fNbins, fFirst and fLast
525// * Set Binning of Abs Times histogram from fAbsNbins, fAbsFirst and fAbsLast
526// * Set Histgram names and titles from fHistName and fHistTitle
527// * Set Abs Times Histgram names and titles from fAbsHistName and fAbsHistTitle
528// * Set X-axis and Y-axis titles from fHistXTitle and fHistYTitle
529// * Set X-axis and Y-axis titles of Abs Times Histogram from fAbsHistXTitle and fAbsHistYTitle
530// * Call InitHists
531//
532void MHCalibrationChargeCam::InitLoGainArrays(const Int_t npixels, const Int_t nareas, const Int_t nsectors)
533{
534
535 if (!IsLoGain())
536 return;
537
538
539 const Int_t logainsamples = fRunHeader->GetNumSamplesLoGain();
540
541 TH1F *h;
542
543 if (fLoGainArray->GetSize()==0 )
544 {
545 for (Int_t i=0; i<npixels; i++)
546 {
547 fLoGainArray->AddAt(new MHCalibrationChargePix(Form("%sLoGainPix%04d",fHistName.Data(),i),
548 Form("%s Low Gain Pixel %04d",fHistTitle.Data(),i)),i);
549
550 MHCalibrationChargePix &pix = (MHCalibrationChargePix&)(*this)(i);
551
552 pix.SetBinning(fLoGainNbins, fLoGainFirst, fLoGainLast);
553 pix.SetProbLimit(fProbLimit);
554
555 pix.SetBinningAbsTime(logainsamples, -0.5, logainsamples-0.5);
556 pix.SetPickupLimit(fgLoGainPickupLimit);
557 pix.SetBlackoutLimit(fgLoGainBlackoutLimit);
558
559 InitHists(pix,(*fBadPixels)[i], i);
560
561 h = pix.GetHAbsTime();
562
563 h->SetName (Form("H%sLoGainPix%04d",fAbsHistName.Data(),i));
564 h->SetTitle(Form("%s Low Gain Pixel %04d",fAbsHistTitle.Data(),i));
565 h->SetXTitle(fAbsHistXTitle.Data());
566 h->SetYTitle(fAbsHistYTitle.Data());
567 }
568 }
569
570 if (fAverageLoGainAreas->GetSize()==0)
571 {
572 for (Int_t j=0; j<nareas; j++)
573 {
574 fAverageLoGainAreas->AddAt(new MHCalibrationChargePix(Form("%sLoGainArea%d",fHistName.Data(),j),
575 Form("%s Low Gain Area Idx %d",fHistTitle.Data(),j)),j);
576
577 MHCalibrationChargePix &pix = (MHCalibrationChargePix&)GetAverageLoGainArea(j);
578
579 pix.SetBinning(fLoGainNbins*(Int_t)TMath::Sqrt((Float_t)npixels/nareas),
580 fLoGainFirst, fLoGainLast);
581
582 pix.SetBinningAbsTime(logainsamples, -0.5, logainsamples-0.5);
583
584 InitHists(pix,fIntensCam ? fIntensCam->GetAverageBadArea(j) : fCam->GetAverageBadArea(j),j);
585
586 h = pix.GetHAbsTime();
587
588 h->SetName (Form("H%sLoGainArea%02d",fAbsHistName.Data(),j));
589 h->SetTitle(Form("%s%s%02d",fAbsHistTitle.Data(),
590 " averaged on event-by-event basis Low Gain Area Idx ",j));
591 h->SetXTitle(fAbsHistXTitle.Data());
592 h->SetYTitle(fAbsHistYTitle.Data());
593 }
594 }
595
596
597 if (fAverageLoGainSectors->GetSize()==0 && IsLoGain())
598 {
599 for (Int_t j=0; j<nsectors; j++)
600 {
601 fAverageLoGainSectors->AddAt(new MHCalibrationChargePix(Form("%sLoGainSector%02d",fHistName.Data(),j),
602 Form("%s Low Gain Sector %02d",fHistTitle.Data(),j)),j);
603
604 MHCalibrationChargePix &pix = (MHCalibrationChargePix&)GetAverageLoGainSector(j);
605
606 pix.SetBinning(fLoGainNbins*(Int_t)TMath::Sqrt((Float_t)npixels/nareas),
607 fLoGainFirst, fLoGainLast);
608
609 pix.SetBinningAbsTime(logainsamples, -0.5, logainsamples-0.5);
610
611 InitHists(pix,fIntensCam ? fIntensCam->GetAverageBadSector(j) : fCam->GetAverageBadSector(j),j);
612
613 h = pix.GetHAbsTime();
614
615 h->SetName (Form("H%sLoGainSector%02d",fAbsHistName.Data(),j));
616 h->SetTitle(Form("%s%s%02d",fAbsHistTitle.Data(),
617 " averaged on event-by-event basis Low Gain Area Sector ",j));
618 h->SetXTitle(fAbsHistXTitle.Data());
619 h->SetYTitle(fAbsHistYTitle.Data());
620 }
621 }
622}
623
624
625// --------------------------------------------------------------------------
626//
627// Retrieves from MExtractedSignalCam:
628// - first used LoGain FADC slice
629//
630// Retrieves from MGeomCam:
631// - number of pixels
632// - number of pixel areas
633// - number of sectors
634//
635// For all TOrdCollection's (including the averaged ones), the following steps are performed:
636//
637// 1) Fill Charges histograms (MHGausEvents::FillHistAndArray()) with:
638// - MExtractedSignalPix::GetExtractedSignalHiGain();
639// - MExtractedSignalPix::GetExtractedSignalLoGain();
640//
641// 2) Set number of saturated slices (MHCalibrationChargePix::AddSaturated()) with:
642// - MExtractedSignalPix::IsHiGainSaturated();
643// - MExtractedSignalPix::GetNumLoGainSaturated();
644//
645// 3) Fill AbsTime histograms (MHCalibrationChargePix::FillAbsTime()) with:
646// - MRawEvtPixelIter::GetIdxMaxHiGainSample();
647// - MRawEvtPixelIter::GetIdxMaxLoGainSample(first slice);
648//
649Bool_t MHCalibrationChargeCam::FillHists(const MParContainer *par, const Stat_t w)
650{
651
652 MExtractedSignalCam *signal = (MExtractedSignalCam*)par;
653 if (!signal)
654 {
655 *fLog << err << "No argument in MExtractedSignalCam::Fill... abort." << endl;
656 return kFALSE;
657 }
658
659 const UInt_t npixels = fGeom->GetNumPixels();
660 const UInt_t nareas = fGeom->GetNumAreas();
661 const UInt_t nsectors = fGeom->GetNumSectors();
662
663 const Int_t hifirst = fSignal->GetFirstUsedSliceHiGain();
664 const Int_t hilast = fSignal->GetLastUsedSliceHiGain();
665 const Int_t lofirst = fSignal->GetFirstUsedSliceLoGain();
666 const Int_t lolast = fSignal->GetLastUsedSliceLoGain();
667
668 fSumhiarea .Reset();
669 fSumloarea .Reset();
670 fTimehiarea .Reset();
671 fTimeloarea .Reset();
672 fSumhisector.Reset();
673 fSumlosector.Reset();
674 fTimehisector.Reset();
675 fTimelosector.Reset();
676
677 fSathiarea .Reset();
678 fSatloarea .Reset();
679 fSathisector.Reset();
680 fSatlosector.Reset();
681
682 for (UInt_t i=0; i<npixels; i++)
683 {
684 MHCalibrationChargePix &histhi = (MHCalibrationChargePix&)(*this)[i];
685 if (histhi.IsExcluded())
686 continue;
687
688 const Int_t aidx = (*fGeom)[i].GetAidx();
689 const Int_t sector = (*fGeom)[i].GetSector();
690
691 const MExtractedSignalPix &pix = (*signal)[i];
692
693 const Int_t sathi = pix.IsHiGainSaturated() ? 1 : 0;
694 histhi.AddSaturated(sathi);
695 if (sathi)
696 {
697 fSathiarea[aidx]++;
698 fSathisector[sector]++;
699 }
700
701 if (pix.IsHiGainValid())
702 {
703 const Float_t sumhi = pix.GetExtractedSignalHiGain();
704 if (IsOscillations())
705 histhi.FillHistAndArray(sumhi);
706 else
707 histhi.FillHist(sumhi);
708
709 fSumhiarea[aidx] += sumhi;
710 fSumhisector[sector] += sumhi;
711 }
712
713 if (!IsLoGain())
714 continue;
715
716 MHCalibrationChargePix &histlo = (MHCalibrationChargePix&)(*this)(i);
717
718 const Int_t satlo = pix.IsLoGainSaturated() ? 1 : 0;
719 histlo.AddSaturated(satlo);
720 if (satlo)
721 {
722 fSatloarea[aidx]++;
723 fSatlosector[sector]++;
724 }
725
726 /*
727 // Previously:
728 fSatloarea[aidx] += pix.GetNumHiGainSaturated();
729 fSatlosector[sector] += pix.GetNumHiGainSaturated();
730 */
731
732 if (pix.IsLoGainValid())
733 {
734 const Float_t sumlo = pix.GetExtractedSignalLoGain();
735
736 if (IsOscillations())
737 histlo.FillHistAndArray(sumlo);
738 else
739 histlo.FillHist(sumlo);
740
741 fSumloarea[aidx] += sumlo;
742 fSumlosector[sector] += sumlo;
743 }
744 }
745
746 const Int_t npix = fRawEvt->GetNumPixels();
747 const Int_t nhi = fRunHeader->GetNumSamplesHiGain();
748 for (int pixid=0; pixid<npix; pixid++)
749 {
750 MHCalibrationChargePix &histhi = (MHCalibrationChargePix&)(*this)[pixid];
751 if (histhi.IsExcluded())
752 continue;
753
754 const Float_t timehi = fRawEvt->GetMax(pixid, hifirst, hilast)+hifirst;
755 histhi.FillAbsTime(timehi);
756
757 const Int_t aidx = (*fGeom)[pixid].GetAidx();
758 const Int_t sector = (*fGeom)[pixid].GetSector();
759
760 fTimehiarea [aidx] += timehi;
761 fTimehisector[sector] += timehi;
762
763 if (IsLoGain())
764 {
765 MHCalibrationChargePix &histlo = (MHCalibrationChargePix&)(*this)(pixid);
766
767 const Float_t timelo = fRawEvt->GetMax(pixid, nhi+lofirst, nhi+lolast)+lofirst;
768 histlo.FillAbsTime(timelo);
769
770 fTimeloarea[aidx] += timelo;
771 fTimelosector[sector] += timelo;
772 }
773 }
774
775 for (UInt_t j=0; j<nareas; j++)
776 {
777
778 const Int_t npix = fAverageAreaNum[j];
779
780 if (npix == 0)
781 continue;
782
783 MHCalibrationChargePix &hipix = (MHCalibrationChargePix&)GetAverageHiGainArea(j);
784
785 if (IsOscillations())
786 hipix.FillHistAndArray(fSumhiarea [j]/npix);
787 else
788 hipix.FillHist(fSumhiarea[j]/npix);
789
790 hipix.AddSaturated ((Float_t)fSathiarea [j]/npix > 0.5 ? 1 : 0);
791 hipix.FillAbsTime (fTimehiarea[j]/npix);
792
793 if (IsLoGain())
794 {
795 MHCalibrationChargePix &lopix = (MHCalibrationChargePix&)GetAverageLoGainArea(j);
796 if (IsOscillations())
797 lopix.FillHistAndArray(fSumloarea [j]/npix);
798 else
799 lopix.FillHist(fSumloarea [j]/npix);
800 lopix.AddSaturated ((Float_t)fSatloarea [j]/npix > 0.5 ? 1 : 0);
801 lopix.FillAbsTime (fTimeloarea[j]/npix);
802 }
803 }
804
805 for (UInt_t j=0; j<nsectors; j++)
806 {
807
808 const Int_t npix = fAverageSectorNum[j];
809
810 if (npix == 0)
811 continue;
812
813 MHCalibrationChargePix &hipix = (MHCalibrationChargePix&)GetAverageHiGainSector(j);
814
815 if (IsOscillations())
816 hipix.FillHistAndArray(fSumhisector [j]/npix);
817 else
818 hipix.FillHist(fSumhisector [j]/npix);
819
820 hipix.AddSaturated ((Float_t)fSathisector[j]/npix > 0.5 ? 1 : 0);
821 hipix.FillAbsTime (fTimehisector[j]/npix);
822
823 if (IsLoGain())
824 {
825 MHCalibrationChargePix &lopix = (MHCalibrationChargePix&)GetAverageLoGainSector(j);
826
827 if (IsOscillations())
828 lopix.FillHistAndArray(fSumlosector [j]/npix);
829 else
830 lopix.FillHist(fSumlosector [j]/npix);
831
832 lopix.AddSaturated ((Float_t)fSatlosector[j]/npix > 0.5 ? 1 : 0);
833 lopix.FillAbsTime (fTimelosector[j]/npix);
834 }
835 }
836
837 return kTRUE;
838}
839
840// --------------------------------------------------------------------------
841//
842// For all TOrdCollection's (including the averaged ones), the following steps are performed:
843//
844// 1) Returns if the pixel is excluded.
845// 2) Tests saturation. In case yes, set the flag: MCalibrationPix::SetHiGainSaturation()
846// or the flag: MBadPixelsPix::SetUncalibrated( MBadPixelsPix::kLoGainSaturated )
847// 3) Store the absolute arrival times in the MCalibrationChargePix's. If flag
848// MCalibrationPix::IsHiGainSaturation() is set, the Low-Gain arrival times are stored,
849// otherwise the Hi-Gain ones.
850// 4) Calls to MHCalibrationCam::FitHiGainArrays() and MCalibrationCam::FitLoGainArrays()
851// with the flags:
852// - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::kHiGainNotFitted )
853// - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::kLoGainNotFitted )
854// - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::kHiGainOscillating )
855// - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::kLoGainOscillating )
856//
857Bool_t MHCalibrationChargeCam::FinalizeHists()
858{
859
860 *fLog << endl;
861
862 TH1F *h = NULL;
863
864 MCalibrationCam *chargecam = fIntensCam ? fIntensCam->GetCam() : fCam;
865
866 const Int_t hifirst = fSignal->GetFirstUsedSliceHiGain();
867 const Int_t hilast = fSignal->GetLastUsedSliceHiGain();
868 const Int_t lofirst = fSignal->GetFirstUsedSliceLoGain();
869 const Int_t lolast = fSignal->GetLastUsedSliceLoGain();
870
871 for (Int_t i=0; i<fHiGainArray->GetSize(); i++)
872 {
873
874 MHCalibrationChargePix &histhi = (MHCalibrationChargePix&)(*this)[i];
875
876 if (histhi.IsExcluded())
877 continue;
878
879 MCalibrationChargePix &pix = (MCalibrationChargePix&)(*chargecam)[i] ;
880
881 const Int_t numsat = histhi.GetSaturated();
882
883 pix.SetNumSaturated(numsat);
884
885 if (numsat > fNumHiGainSaturationLimit*histhi.GetHGausHist()->GetEntries())
886 {
887 pix.SetHiGainSaturation();
888 if (IsOscillations())
889 histhi.CreateFourierSpectrum();
890 continue;
891 }
892
893 MBadPixelsPix &bad = (*fBadPixels)[i];
894
895 h = histhi.GetHGausHist();
896
897 Stat_t overflow = h->GetBinContent(h->GetNbinsX()+1);
898 if (overflow > fOverflowLimit*histhi.GetHGausHist()->GetEntries())
899 {
900 *fLog << warn
901 << "HiGain Hist-overflow " << overflow
902 << " times in: " << histhi.GetName() << " (w/o saturation!) " << endl;
903 bad.SetUncalibrated( MBadPixelsPix::kHiGainOverFlow );
904 }
905
906 overflow = h->GetBinContent(0);
907 if (overflow > fOverflowLimit*histhi.GetHGausHist()->GetEntries())
908 {
909 *fLog << warn
910 << "HiGain Hist-underflow " << overflow
911 << " times in pix: " << histhi.GetName() << " (w/o saturation!) " << endl;
912 bad.SetUncalibrated( MBadPixelsPix::kHiGainOverFlow );
913 }
914
915 FinalizeAbsTimes(histhi, pix, bad, hifirst, hilast);
916 }
917
918 if (IsLoGain())
919 for (Int_t i=0; i<fLoGainArray->GetSize(); i++)
920 {
921
922 MHCalibrationChargePix &histlo = (MHCalibrationChargePix&)(*this)(i);
923
924 if (histlo.IsExcluded())
925 continue;
926
927 MBadPixelsPix &bad = (*fBadPixels)[i];
928 if (histlo.GetSaturated() > fNumLoGainSaturationLimit*histlo.GetHGausHist()->GetEntries())
929 {
930 *fLog << warn << "Pixel " << setw(4) << i << ": More than " << Form("%.1f%%", fNumLoGainSaturationLimit*100) << " lo-gains saturated." << endl;
931 bad.SetUncalibrated( MBadPixelsPix::kLoGainSaturation );
932 if (IsOscillations())
933 histlo.CreateFourierSpectrum();
934 continue;
935 }
936
937 MCalibrationChargePix &pix = (MCalibrationChargePix&)(*chargecam)[i] ;
938 if (!pix.IsHiGainSaturation())
939 continue;
940
941 h = histlo.GetHGausHist();
942
943 Stat_t overflow = h->GetBinContent(h->GetNbinsX()+1);
944 if (overflow > fOverflowLimit*histlo.GetHGausHist()->GetEntries())
945 {
946 *fLog << warn
947 << "LoGain Hist-overflow " << overflow
948 << " times in: " << histlo.GetName() << " (w/o saturation!) " << endl;
949 bad.SetUncalibrated( MBadPixelsPix::kLoGainOverFlow );
950 }
951
952 overflow = h->GetBinContent(0);
953 if (overflow > fOverflowLimit*histlo.GetHGausHist()->GetEntries())
954 {
955 *fLog << warn
956 << "LoGain Hist-underflow " << overflow
957 << " times in: " << histlo.GetName() << " (w/o saturation!) " << endl;
958 bad.SetUncalibrated( MBadPixelsPix::kLoGainOverFlow );
959 }
960
961 FinalizeAbsTimes(histlo, pix, bad, lofirst, lolast);
962 }
963
964 for (Int_t j=0; j<fAverageHiGainAreas->GetSize(); j++)
965 {
966
967 MHCalibrationChargePix &histhi = (MHCalibrationChargePix&)GetAverageHiGainArea(j);
968 MCalibrationChargePix &pix = (MCalibrationChargePix&)chargecam->GetAverageArea(j);
969
970 if (histhi.GetSaturated() > fNumHiGainSaturationLimit*histhi.GetHGausHist()->GetEntries())
971 {
972 pix.SetHiGainSaturation();
973 if (IsOscillations())
974 histhi.CreateFourierSpectrum();
975 continue;
976 }
977
978 MBadPixelsPix &bad = chargecam->GetAverageBadArea(j);
979 FinalizeAbsTimes(histhi, pix, bad, hifirst, hilast);
980 }
981
982 if (IsLoGain())
983 for (Int_t j=0; j<fAverageLoGainAreas->GetSize(); j++)
984 {
985
986 MHCalibrationChargePix &histlo = (MHCalibrationChargePix&)GetAverageLoGainArea(j);
987
988 if (histlo.GetSaturated() > fNumLoGainSaturationLimit*histlo.GetHGausHist()->GetEntries())
989 {
990 *fLog << warn << "Area " << setw(4) << j << ": More than " << Form("%.1f%%", fNumLoGainSaturationLimit*100) << " lo-gains saturated." << endl;
991 if (IsOscillations())
992 histlo.CreateFourierSpectrum();
993 continue;
994 }
995
996 MCalibrationChargePix &pix = (MCalibrationChargePix&)chargecam->GetAverageArea(j);
997
998 if (pix.IsHiGainSaturation())
999 {
1000 MBadPixelsPix &bad = chargecam->GetAverageBadArea(j);
1001 FinalizeAbsTimes(histlo, pix, bad, lofirst, lolast);
1002 }
1003
1004 }
1005
1006 for (Int_t j=0; j<fAverageHiGainSectors->GetSize(); j++)
1007 {
1008
1009 MHCalibrationChargePix &histhi = (MHCalibrationChargePix&)GetAverageHiGainSector(j);
1010 MCalibrationChargePix &pix = (MCalibrationChargePix&)chargecam->GetAverageSector(j);
1011
1012 if (histhi.GetSaturated() > fNumHiGainSaturationLimit*histhi.GetHGausHist()->GetEntries())
1013 {
1014 pix.SetHiGainSaturation();
1015 if (IsOscillations())
1016 histhi.CreateFourierSpectrum();
1017 continue;
1018 }
1019
1020 MBadPixelsPix &bad = chargecam->GetAverageBadSector(j);
1021 FinalizeAbsTimes(histhi, pix, bad, hifirst, hilast);
1022 }
1023
1024 if (IsLoGain())
1025 for (Int_t j=0; j<fAverageLoGainSectors->GetSize(); j++)
1026 {
1027
1028 MHCalibrationChargePix &histlo = (MHCalibrationChargePix&)GetAverageLoGainSector(j);
1029 MBadPixelsPix &bad = chargecam->GetAverageBadSector(j);
1030
1031 if (histlo.GetSaturated() > fNumLoGainSaturationLimit*histlo.GetHGausHist()->GetEntries())
1032 {
1033 *fLog << warn << "Sector " << setw(4) << j << ": More than " << Form("%.1f%%", fNumLoGainSaturationLimit*100) << " lo-gains saturated." << endl;
1034 bad.SetUncalibrated( MBadPixelsPix::kLoGainSaturation );
1035 if (IsOscillations())
1036 histlo.CreateFourierSpectrum();
1037 continue;
1038 }
1039
1040 MCalibrationChargePix &pix = (MCalibrationChargePix&)chargecam->GetAverageSector(j);
1041
1042 if (pix.IsHiGainSaturation())
1043 FinalizeAbsTimes(histlo, pix, bad, lofirst, lolast);
1044 }
1045
1046 //
1047 // Perform the fitting for the High Gain (done in MHCalibrationCam)
1048 //
1049 FitHiGainArrays(*chargecam, *fBadPixels,
1050 MBadPixelsPix::kHiGainNotFitted,
1051 MBadPixelsPix::kHiGainOscillating);
1052
1053 //
1054 // Perform the fitting for the Low Gain (done in MHCalibrationCam)
1055 //
1056 if (IsLoGain())
1057 FitLoGainArrays(*chargecam, *fBadPixels,
1058 MBadPixelsPix::kLoGainNotFitted,
1059 MBadPixelsPix::kLoGainOscillating);
1060
1061 //
1062 // Check for pixels which have the high-gain saturated, but the low-gain
1063 // switch not applied in sufficient cases. Have to exclude these pixels,
1064 // although they not abnormal. They simply cannot be calibrated...
1065 //
1066 for (Int_t i=0; i<fLoGainArray->GetSize(); i++)
1067 {
1068
1069 MHCalibrationChargePix &histlo = (MHCalibrationChargePix&)(*this)(i);
1070 if (histlo.IsExcluded())
1071 continue;
1072
1073 MCalibrationChargePix &pix = (MCalibrationChargePix&)(*chargecam)[i] ;
1074 if (!pix.IsHiGainSaturation())
1075 continue;
1076
1077 //
1078 // Now,treat only low-gain calibrated pixels:
1079 //
1080 const Double_t lim = fNumLoGainBlackoutLimit*histlo.GetHGausHist()->GetEntries();
1081 if (histlo.GetBlackout() <= lim)
1082 continue;
1083
1084 MBadPixelsPix &bad = (*fBadPixels)[i];
1085 bad.SetUnsuitable(MBadPixelsPix::kUnsuitableRun);
1086 bad.SetUncalibrated(MBadPixelsPix::kLoGainBlackout);
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 Int_t first, Int_t last)
1102{
1103 const Float_t mean = hist.GetAbsTimeMean();
1104 const Float_t rms = hist.GetAbsTimeRms();
1105
1106 pix.SetAbsTimeMean(mean);
1107 pix.SetAbsTimeRms(rms);
1108
1109 const Float_t lowerlimit = (Float_t)first+1;// + fTimeLowerLimit;
1110 const Float_t upperlimit = (Float_t)last -1;// - fTimeUpperLimit;
1111
1112 // FIXME: instead of checking whether the maximum is in the first or
1113 // last extracted slice we should check whether the extractor
1114 // was able to properly extract the signal!!!
1115
1116 if (mean<lowerlimit)
1117 {
1118 *fLog << warn << hist.GetName() << ": Mean Arr.Time: "
1119 << Form("%4.1f < %4.1f (first used fadc sl+1)", mean, lowerlimit) << endl;
1120 bad.SetUncalibrated(MBadPixelsPix::kMeanTimeInFirstBin);
1121 }
1122
1123 if (mean>upperlimit)
1124 {
1125 *fLog << warn << hist.GetName() << ": Mean Arr.Time: "
1126 << Form("%4.1f > %4.1f (last used fadc sl-1)", mean, upperlimit) << endl;
1127 bad.SetUncalibrated(MBadPixelsPix::kMeanTimeInLast2Bins);
1128 }
1129}
1130
1131// --------------------------------------------------------------------------
1132//
1133// Sets all pixels to MBadPixelsPix::kUnsuitableRun, if following flags are set:
1134// - MBadPixelsPix::kLoGainSaturation
1135//
1136// Sets all pixels to MBadPixelsPix::kUnreliableRun, if following flags are set:
1137// - if MBadPixelsPix::kHiGainNotFitted and !MCalibrationPix::IsHiGainSaturation()
1138// - if MBadPixelsPix::kHiGainOscillating and !MCalibrationPix::IsHiGainSaturation()
1139// - if MBadPixelsPix::kLoGainNotFitted and MCalibrationPix::IsLoGainSaturation()
1140// - if MBadPixelsPix::kLoGainOscillating and MCalibrationPix::IsLoGainSaturation()
1141//
1142void MHCalibrationChargeCam::FinalizeBadPixels()
1143{
1144
1145 MCalibrationCam *chargecam = fIntensCam ? fIntensCam->GetCam() : fCam;
1146
1147 for (Int_t i=0; i<fBadPixels->GetSize(); i++)
1148 {
1149
1150 MBadPixelsPix &bad = (*fBadPixels)[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 if (!pix.IsHiGainSaturation())
1168 bad.SetUnsuitable(MBadPixelsPix::kUnreliableRun);
1169
1170 if (bad.IsUncalibrated(MBadPixelsPix::kLoGainOscillating))
1171 if (pix.IsHiGainSaturation())
1172 bad.SetUnsuitable(MBadPixelsPix::kUnreliableRun);
1173 }
1174 }
1175}
1176
1177// --------------------------------------------------------------------------
1178//
1179// Calls MHCalibrationPix::DrawClone() for pixel idx
1180//
1181void MHCalibrationChargeCam::DrawPixelContent(Int_t idx) const
1182{
1183 (*this)[idx].DrawClone();
1184}
1185
1186
1187// -----------------------------------------------------------------------------
1188//
1189// Default draw:
1190//
1191// Displays the averaged areas, both High Gain and Low Gain
1192//
1193// Calls the Draw of the fAverageHiGainAreas and fAverageLoGainAreas objects with options
1194//
1195void MHCalibrationChargeCam::Draw(const Option_t *opt)
1196{
1197
1198 const Int_t nareas = fAverageHiGainAreas->GetSize();
1199 if (nareas == 0)
1200 return;
1201
1202 TString option(opt);
1203 option.ToLower();
1204
1205 if (!option.Contains("datacheck"))
1206 {
1207 MHCalibrationCam::Draw(opt);
1208 return;
1209 }
1210
1211 //
1212 // From here on , the datacheck - Draw
1213 //
1214 TVirtualPad *pad = gPad ? gPad : MH::MakeDefCanvas(this);
1215 pad->SetBorderMode(0);
1216 pad->Divide(1,nareas);
1217
1218 //
1219 // Loop over inner and outer pixels
1220 //
1221 for (Int_t i=0; i<nareas;i++)
1222 {
1223
1224 pad->cd(i+1);
1225
1226 MHCalibrationChargePix &hipix = (MHCalibrationChargePix&)GetAverageHiGainArea(i);
1227 //
1228 // Ask for Hi-Gain saturation
1229 //
1230 if (hipix.GetSaturated() > fNumHiGainSaturationLimit*hipix.GetHGausHist()->GetEntries() && IsLoGain())
1231 {
1232 MHCalibrationChargePix &lopix = (MHCalibrationChargePix&)GetAverageLoGainArea(i);
1233 DrawDataCheckPixel(lopix,i ? fOuterRefCharge : fInnerRefCharge);
1234 }
1235 else
1236 DrawDataCheckPixel(hipix,i ? fOuterRefCharge : fInnerRefCharge);
1237
1238 }
1239}
1240
1241// -----------------------------------------------------------------------------
1242//
1243// Draw the average pixel for the datacheck:
1244//
1245// Displays the averaged areas, both High Gain and Low Gain
1246//
1247// Calls the Draw of the fAverageHiGainAreas and fAverageLoGainAreas objects with options
1248//
1249void MHCalibrationChargeCam::DrawDataCheckPixel(MHCalibrationChargePix &pix, const Float_t refline)
1250{
1251 if (pix.IsEmpty())
1252 return;
1253
1254 TVirtualPad *pad = gPad;
1255 pad->Divide(1,2, 1e-10, 1e-10);
1256 pad->cd(1);
1257
1258 gPad->SetBorderMode(0);
1259 gPad->SetTicks();
1260 if (!pix.IsEmpty() && !pix.IsOnlyOverflow() && !pix.IsOnlyUnderflow())
1261 gPad->SetLogy();
1262
1263 TH1F *hist = pix.GetHGausHist();
1264 TAxis *xaxe = hist->GetXaxis();
1265 TAxis *yaxe = hist->GetYaxis();
1266 xaxe->CenterTitle();
1267 yaxe->CenterTitle();
1268 xaxe->SetTitleSize(0.07);
1269 yaxe->SetTitleSize(0.07);
1270 xaxe->SetTitleOffset(0.7);
1271 yaxe->SetTitleOffset(0.55);
1272 xaxe->SetLabelSize(0.06);
1273 yaxe->SetLabelSize(0.06);
1274 xaxe->SetTitle(hist->GetXaxis()->GetTitle());
1275 yaxe->SetTitle(hist->GetYaxis()->GetTitle());
1276 xaxe->SetRange(hist->GetMaximumBin()-30, hist->GetMaximumBin()+30);
1277
1278 gStyle->SetOptFit();
1279
1280 hist->Draw();
1281
1282 TF1 *fit = pix.GetFGausFit();
1283 if (fit)
1284 {
1285 switch (fColor)
1286 {
1287 case MCalibrationCam::kGREEN:
1288 fit->SetLineColor(kGreen);
1289 break;
1290 case MCalibrationCam::kBLUE:
1291 fit->SetLineColor(kBlue);
1292 break;
1293 case MCalibrationCam::kUV:
1294 fit->SetLineColor(106);
1295 break;
1296 case MCalibrationCam::kCT1:
1297 fit->SetLineColor(006);
1298 break;
1299 default:
1300 fit->SetLineColor(kRed);
1301 }
1302 fit->Draw("same");
1303 }
1304
1305 DisplayRefLines(hist,refline);
1306
1307 pad->cd(2);
1308 gPad->SetBorderMode(0);
1309 gPad->SetTicks();
1310
1311 pix.CreateGraphEvents();
1312 TGraph *gr = pix.GetGraphEvents();
1313 if (gr)
1314 {
1315 TH1F *null2 = gr->GetHistogram();
1316
1317 null2->SetMinimum(pix.GetMean()-10.*pix.GetSigma());
1318 null2->SetMaximum(pix.GetMean()+10.*pix.GetSigma());
1319 null2->SetStats(kFALSE);
1320
1321 //
1322 // set the labels bigger
1323 //
1324 TAxis *xaxe2 = null2->GetXaxis();
1325 TAxis *yaxe2 = null2->GetYaxis();
1326 xaxe2->CenterTitle();
1327 yaxe2->CenterTitle();
1328 xaxe2->SetTitleSize(0.07);
1329 yaxe2->SetTitleSize(0.07);
1330 xaxe2->SetTitleOffset(0.7);
1331 yaxe2->SetTitleOffset(0.55);
1332 xaxe2->SetLabelSize(0.06);
1333 yaxe2->SetLabelSize(0.06);
1334 xaxe2->SetRangeUser(0, pix.GetEvents()->GetSize()/pix.GetEventFrequency());
1335 }
1336
1337 pix.DrawEvents();
1338
1339 // newpad->cd(3);
1340 // pix.DrawPowerSpectrum(*newpad,4);
1341}
1342
1343void MHCalibrationChargeCam::DisplayRefLines(const TH1F *hist, const Float_t refline) const
1344{
1345 TLine *line = new TLine(refline, 0, refline, hist->GetMaximum());
1346 line->SetLineColor(106);
1347 line->SetLineStyle(2);
1348 line->SetLineWidth(3);
1349 line->SetBit(kCanDelete);
1350 line->Draw();
1351
1352 TLegend *leg = new TLegend(0.8,0.01,0.99,0.45);
1353 leg->AddEntry(line, "10 Leds UV", "l");
1354 leg->SetBit(kCanDelete);
1355 leg->Draw();
1356}
1357
1358Int_t MHCalibrationChargeCam::ReadEnv(const TEnv &env, TString prefix, Bool_t print)
1359{
1360
1361 Int_t rc = MHCalibrationCam::ReadEnv(env,prefix,print);
1362 if (rc==kERROR)
1363 return kERROR;
1364
1365 if (IsEnvDefined(env, prefix, "HiGainNbins", print))
1366 {
1367 fNbins = GetEnvValue(env, prefix, "HiGainNbins", fNbins);
1368 rc = kTRUE;
1369 }
1370
1371 if (IsEnvDefined(env, prefix, "HiGainFirst", print))
1372 {
1373 fFirst = GetEnvValue(env, prefix, "HiGainFirst", fFirst);
1374 rc = kTRUE;
1375 }
1376
1377 if (IsEnvDefined(env, prefix, "HiGainLast", print))
1378 {
1379 fLast = GetEnvValue(env, prefix, "HiGainLast", fLast);
1380 rc = kTRUE;
1381 }
1382
1383 if (IsEnvDefined(env, prefix, "LoGainNbins", print))
1384 {
1385 fLoGainNbins = GetEnvValue(env, prefix, "LoGainNbins", fLoGainNbins);
1386 rc = kTRUE;
1387 }
1388
1389 if (IsEnvDefined(env, prefix, "LoGainFirst", print))
1390 {
1391 fLoGainFirst = GetEnvValue(env, prefix, "LoGainFirst", fLoGainFirst);
1392 rc = kTRUE;
1393 }
1394
1395 if (IsEnvDefined(env, prefix, "LoGainLast", print))
1396 {
1397 fLoGainLast = GetEnvValue(env, prefix, "LoGainLast", fLoGainLast);
1398 rc = kTRUE;
1399 }
1400
1401 if (IsEnvDefined(env, prefix, "TimeLowerLimit", print))
1402 {
1403 SetTimeLowerLimit(GetEnvValue(env, prefix, "TimeLowerLimit", fTimeLowerLimit));
1404 rc = kTRUE;
1405 }
1406
1407 if (IsEnvDefined(env, prefix, "TimeUpperLimit", print))
1408 {
1409 SetTimeUpperLimit(GetEnvValue(env, prefix, "TimeUpperLimit", fTimeUpperLimit));
1410 rc = kTRUE;
1411 }
1412
1413 if (IsEnvDefined(env, prefix, "ReferenceFile", print))
1414 {
1415 SetReferenceFile(GetEnvValue(env,prefix,"ReferenceFile",fReferenceFile.Data()));
1416 rc = kTRUE;
1417 }
1418
1419 if (IsEnvDefined(env, prefix, "NumHiGainSaturationLimit", print))
1420 {
1421 SetNumHiGainSaturationLimit(GetEnvValue(env, prefix, "NumHiGainSaturationLimit", fNumHiGainSaturationLimit));
1422 rc = kTRUE;
1423 }
1424
1425 if (IsEnvDefined(env, prefix, "NumLoGainSaturationLimit", print))
1426 {
1427 SetNumLoGainSaturationLimit(GetEnvValue(env, prefix, "NumLoGainSaturationLimit", fNumLoGainSaturationLimit));
1428 rc = kTRUE;
1429 }
1430
1431 if (IsEnvDefined(env, prefix, "NumLoGainBlackoutLimit", print))
1432 {
1433 SetNumLoGainBlackoutLimit(GetEnvValue(env, prefix, "NumLoGainBlackoutLimit", fNumLoGainBlackoutLimit));
1434 rc = kTRUE;
1435 }
1436
1437 TEnv refenv(fReferenceFile);
1438
1439 fInnerRefCharge = refenv.GetValue("InnerRefCharge",fInnerRefCharge);
1440 fOuterRefCharge = refenv.GetValue("OuterRefCharge",fOuterRefCharge);
1441
1442 return rc;
1443}
Note: See TracBrowser for help on using the repository browser.