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

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