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

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