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

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