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

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