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

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