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

Last change on this file since 5029 was 5029, checked in by gaug, 20 years ago
*** empty log message ***
File size: 52.4 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 "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 <TPad.h>
144#include <TVirtualPad.h>
145#include <TCanvas.h>
146#include <TStyle.h>
147#include <TF1.h>
148#include <TLatex.h>
149#include <TLegend.h>
150#include <TGraph.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 or MCalibrationIntensityChargeCam
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// * Initialize an MHCalibrationChargePix
385// * Set Binning from fNbins, fFirst and fLast
386// * Set Binning of Abs Times histogram from fAbsNbins, fAbsFirst and fAbsLast
387// * Set Histgram names and titles from fHistName and fHistTitle
388// * Set Abs Times Histgram names and titles from fAbsHistName and fAbsHistTitle
389// * Set X-axis and Y-axis titles from fHistXTitle and fHistYTitle
390// * Set X-axis and Y-axis titles of Abs Times Histogram from fAbsHistXTitle and fAbsHistYTitle
391// * Call InitHists
392//
393//
394void MHCalibrationChargeCam::InitHiGainArrays(const Int_t npixels, const Int_t nareas, const Int_t nsectors)
395{
396
397 TH1F *h;
398
399 const Int_t higainsamples = fRunHeader->GetNumSamplesHiGain();
400
401 if (fHiGainArray->GetEntries()==0)
402 {
403 fHiGainArray->Expand(npixels);
404 for (Int_t i=0; i<npixels; i++)
405 {
406 (*fHiGainArray)[i] = new MHCalibrationChargePix(Form("%s%s",fHistName.Data(),"HiGainPix"),
407 Form("%s%s",fHistTitle.Data()," High Gain Pixel"));
408
409 MHCalibrationChargePix &pix = (MHCalibrationChargePix&)(*this)[i];
410
411 pix.SetNbins(fNbins);
412 pix.SetFirst(fFirst);
413 pix.SetLast (fLast);
414
415 pix.SetAbsTimeNbins(higainsamples);
416 pix.SetAbsTimeFirst(-0.5);
417 pix.SetAbsTimeLast(higainsamples-0.5);
418
419 h = pix.GetHGausHist();
420
421 h->SetName (Form("%s%s%s","H",fHistName.Data(),"HiGainPix"));
422 h->SetTitle(Form("%s%s",fHistTitle.Data()," High Gain Pixel "));
423 h->SetXTitle(fHistXTitle.Data());
424 h->SetYTitle(fHistYTitle.Data());
425
426 h = pix.GetHAbsTime();
427
428 h->SetName (Form("%s%s%s","H",fAbsHistName.Data(),"HiGainPix"));
429 h->SetTitle(Form("%s%s",fAbsHistTitle.Data()," High Gain Pixel "));
430 h->SetXTitle(fAbsHistXTitle.Data());
431 h->SetYTitle(fAbsHistYTitle.Data());
432
433 MBadPixelsPix &bad = fIntensBad ? (*fIntensBad)[i] : (*fBadPixels)[i];
434 InitHists((*this)[i],bad,i);
435 }
436 }
437
438
439 if (fAverageHiGainAreas->GetEntries()==0)
440 {
441 fAverageHiGainAreas->Expand(nareas);
442
443 for (Int_t j=0; j<nareas; j++)
444 {
445 (*fAverageHiGainAreas)[j] =
446 new MHCalibrationChargePix(Form("%s%s",fHistName.Data(),"HiGainArea"),
447 Form("%s%s",fHistTitle.Data()," High Gain Area Idx "));
448
449 MHCalibrationChargePix &pix = (MHCalibrationChargePix&)GetAverageHiGainArea(j);
450
451 pix.SetNbins(fNbins*(Int_t)TMath::Sqrt((Float_t)npixels/nareas));
452 pix.SetFirst(fFirst);
453 pix.SetLast (fLast);
454
455 pix.SetAbsTimeNbins(higainsamples);
456 pix.SetAbsTimeFirst(-0.5);
457 pix.SetAbsTimeLast(higainsamples-0.5);
458
459 h = pix.GetHGausHist();
460
461 h->SetName (Form("%s%s%s","H",fHistName.Data(),"HiGainArea"));
462 h->SetXTitle(fHistXTitle.Data());
463 h->SetYTitle(fHistYTitle.Data());
464
465 if (fGeom->InheritsFrom("MGeomCamMagic"))
466 {
467 h->SetTitle(Form("%s%s%s%s",fHistTitle.Data()," averaged on event-by-event basis ",
468 j==0 ? "Inner Pixels " : "Outer Pixels ","High Gain Runs: "));
469 pix.InitBins();
470 pix.SetEventFrequency(fPulserFrequency);
471 }
472 else
473 {
474 h->SetTitle(Form("%s%s",fHistTitle.Data(),
475 " averaged on event-by-event basis High Gain Area Idx "));
476 InitHists(pix,fIntensCam ? fIntensCam->GetAverageBadArea(j) : fCam->GetAverageBadArea(j),j);
477 }
478
479 h = pix.GetHAbsTime();
480
481 h->SetName (Form("%s%s%s","H",fAbsHistName.Data(),"HiGainArea"));
482 h->SetTitle(Form("%s%s",fAbsHistTitle.Data(),
483 " averaged on event-by-event basis High Gain Area Idx "));
484 h->SetXTitle(fAbsHistXTitle.Data());
485 h->SetYTitle(fAbsHistYTitle.Data());
486 }
487 }
488
489 if (fAverageHiGainSectors->GetEntries()==0)
490 {
491 fAverageHiGainSectors->Expand(nsectors);
492
493 for (Int_t j=0; j<nsectors; j++)
494 {
495 (*fAverageHiGainSectors)[j] =
496 new MHCalibrationChargePix(Form("%s%s",fHistName.Data(),"HiGainSector"),
497 Form("%s%s",fHistTitle.Data()," High Gain Sector "));
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 h = pix.GetHGausHist();
510
511 h->SetName (Form("%s%s%s","H",fHistName.Data(),"HiGainSector"));
512 h->SetTitle(Form("%s%s",fHistTitle.Data()," High Gain Sector "));
513 h->SetXTitle(fHistXTitle.Data());
514 h->SetYTitle(fHistYTitle.Data());
515
516 h = pix.GetHAbsTime();
517
518 h->SetName (Form("%s%s%s","H",fAbsHistName.Data(),"HiGainSector"));
519 h->SetTitle(Form("%s%s",fAbsHistTitle.Data(),
520 " averaged on event-by-event basis High Gain Area Sector "));
521 h->SetXTitle(fAbsHistXTitle.Data());
522 h->SetYTitle(fAbsHistYTitle.Data());
523
524 InitHists(pix,fIntensCam ? fIntensCam->GetAverageBadSector(j) : fCam->GetAverageBadSector(j),j);
525 }
526 }
527}
528
529//--------------------------------------------------------------------------------------
530//
531// Return, if IsLoGain() is kFALSE
532//
533// Retrieve:
534// - fRunHeader->GetNumSamplesHiGain();
535//
536// Initializes the Low Gain Arrays:
537//
538// - Expand fLoGainArrays to npixels
539// - Expand fAverageLoGainAreas to nareas
540// - Expand fAverageLoGainSectors to nsectors
541//
542// - For every entry in the expanded arrays:
543// * Initialize an MHCalibrationChargePix
544// * Set Binning from fNbins, fFirst and fLast
545// * Set Binning of Abs Times histogram from fAbsNbins, fAbsFirst and fAbsLast
546// * Set Histgram names and titles from fHistName and fHistTitle
547// * Set Abs Times Histgram names and titles from fAbsHistName and fAbsHistTitle
548// * Set X-axis and Y-axis titles from fHistXTitle and fHistYTitle
549// * Set X-axis and Y-axis titles of Abs Times Histogram from fAbsHistXTitle and fAbsHistYTitle
550// * Call InitHists
551//
552void MHCalibrationChargeCam::InitLoGainArrays(const Int_t npixels, const Int_t nareas, const Int_t nsectors)
553{
554
555 if (!IsLoGain())
556 return;
557
558 const Int_t logainsamples = fRunHeader->GetNumSamplesLoGain();
559
560 TH1F *h;
561
562 if (fLoGainArray->GetEntries()==0 )
563 {
564 fLoGainArray->Expand(npixels);
565
566 for (Int_t i=0; i<npixels; i++)
567 {
568 (*fLoGainArray)[i] =
569 new MHCalibrationChargePix(Form("%s%s",fHistName.Data(),"LoGainPix"),
570 Form("%s%s",fHistTitle.Data()," Low Gain Pixel"));
571
572 MHCalibrationChargePix &pix = (MHCalibrationChargePix&)(*this)(i);
573
574 pix.SetNbins(fLoGainNbins);
575 pix.SetFirst(fLoGainFirst);
576 pix.SetLast (fLoGainLast);
577
578 pix.SetAbsTimeNbins(logainsamples);
579 pix.SetAbsTimeFirst(-0.5);
580 pix.SetAbsTimeLast(logainsamples-0.5);
581
582 h = pix.GetHGausHist();
583
584 h->SetName (Form("%s%s%s","H",fHistName.Data(),"LoGainPix"));
585 h->SetTitle(Form("%s%s",fHistTitle.Data()," Low Gain Pixel "));
586 h->SetXTitle(fHistXTitle.Data());
587 h->SetYTitle(fHistYTitle.Data());
588
589 h = pix.GetHAbsTime();
590
591 h->SetName (Form("%s%s%s","H",fAbsHistName.Data(),"HiGainPix"));
592 h->SetTitle(Form("%s%s",fAbsHistTitle.Data()," High Gain Pixel "));
593 h->SetXTitle(fAbsHistXTitle.Data());
594 h->SetYTitle(fAbsHistYTitle.Data());
595
596 MBadPixelsPix &bad = fIntensBad ? (*fIntensBad)[i] : (*fBadPixels)[i];
597 InitHists(pix,bad,i);
598 }
599 }
600
601 if (fAverageLoGainAreas->GetEntries()==0)
602 {
603 fAverageLoGainAreas->Expand(nareas);
604
605 for (Int_t j=0; j<nareas; j++)
606 {
607 (*fAverageLoGainAreas)[j] =
608 new MHCalibrationChargePix(Form("%s%s",fHistName.Data(),"LoGainArea"),
609 Form("%s%s",fHistTitle.Data()," Low Gain Area Idx "));
610
611 MHCalibrationChargePix &pix = (MHCalibrationChargePix&)GetAverageLoGainArea(j);
612
613 pix.SetNbins(fLoGainNbins*(Int_t)TMath::Sqrt((Float_t)npixels/nareas));
614 pix.SetFirst(fLoGainFirst);
615 pix.SetLast (fLoGainLast);
616
617 pix.SetAbsTimeNbins(logainsamples);
618 pix.SetAbsTimeFirst(-0.5);
619 pix.SetAbsTimeLast(logainsamples-0.5);
620
621 h = pix.GetHGausHist();
622
623 h->SetName (Form("%s%s%s","H",fHistName.Data(),"LoGainArea"));
624 h->SetXTitle(fHistXTitle.Data());
625 h->SetYTitle(fHistYTitle.Data());
626
627
628 if (fGeom->InheritsFrom("MGeomCamMagic"))
629 {
630 h->SetTitle(Form("%s%s%s%s",fHistTitle.Data()," averaged on event-by-event basis ",
631 j==0 ? "Inner Pixels " : "Outer Pixels ","Low Gain Runs: "));
632 pix.InitBins();
633 pix.SetEventFrequency(fPulserFrequency);
634 }
635 else
636 {
637 h->SetTitle(Form("%s%s",fHistTitle.Data()," averaged on event-by-event basis Low Gain Area Idx "));
638 InitHists(pix,fIntensCam ? fIntensCam->GetAverageBadArea(j) : fCam->GetAverageBadArea(j),j);
639 }
640
641 h = pix.GetHAbsTime();
642
643 h->SetName (Form("%s%s%s","H",fAbsHistName.Data(),"LoGainArea"));
644 h->SetTitle(Form("%s%s",fAbsHistTitle.Data(),
645 " averaged on event-by-event basis Low Gain Area Idx "));
646 h->SetXTitle(fAbsHistXTitle.Data());
647 h->SetYTitle(fAbsHistYTitle.Data());
648
649 }
650 }
651
652
653 if (fAverageLoGainSectors->GetEntries()==0 && IsLoGain())
654 {
655 fAverageLoGainSectors->Expand(nsectors);
656
657 for (Int_t j=0; j<nsectors; j++)
658 {
659 (*fAverageLoGainSectors)[j] =
660 new MHCalibrationChargePix(Form("%s%s",fHistName.Data(),"LoGainSector"),
661 Form("%s%s",fHistTitle.Data()," Low Gain Sector "));
662
663 MHCalibrationChargePix &pix = (MHCalibrationChargePix&)GetAverageLoGainSector(j);
664
665 pix.SetNbins(fLoGainNbins*(Int_t)TMath::Sqrt((Float_t)npixels/nareas));
666 pix.SetFirst(fLoGainFirst);
667 pix.SetLast (fLoGainLast);
668
669 pix.SetAbsTimeNbins(logainsamples);
670 pix.SetAbsTimeFirst(-0.5);
671 pix.SetAbsTimeLast(logainsamples-0.5);
672
673 h = pix.GetHGausHist();
674
675 h->SetName (Form("%s%s%s","H",fHistName.Data(),"LoGainSector"));
676 h->SetTitle(Form("%s%s",fHistTitle.Data()," Low Gain Sector "));
677 h->SetXTitle(fHistXTitle.Data());
678 h->SetYTitle(fHistYTitle.Data());
679
680 h = pix.GetHAbsTime();
681
682 h->SetName (Form("%s%s%s","H",fAbsHistName.Data(),"LoGainSector"));
683 h->SetTitle(Form("%s%s",fAbsHistTitle.Data(),
684 " averaged on event-by-event basis Low Gain Area Sector "));
685 h->SetXTitle(fAbsHistXTitle.Data());
686 h->SetYTitle(fAbsHistYTitle.Data());
687
688 //
689 // Adapt the range for the case, the intense blue is used:
690 // FIXME: this is a nasty workaround, but for the moment necessary
691 // in order to avoid default memory space.
692 //
693 if (fGeom->InheritsFrom("MGeomCamMagic"))
694 {
695 if ( fColor == MCalibrationCam::kBLUE)
696 {
697 pix.SetFirst(-10.5);
698 pix.SetLast(999.5);
699 pix.SetNbins(3030);
700 }
701 }
702
703 InitHists(pix,fIntensCam ? fIntensCam->GetAverageBadSector(j) : fCam->GetAverageBadSector(j),j);
704 }
705 }
706}
707
708
709// --------------------------------------------------------------------------
710//
711// Retrieves from MExtractedSignalCam:
712// - first used LoGain FADC slice
713//
714// Retrieves from MGeomCam:
715// - number of pixels
716// - number of pixel areas
717// - number of sectors
718//
719// For all TObjArray's (including the averaged ones), the following steps are performed:
720//
721// 1) Fill Charges histograms (MHGausEvents::FillHistAndArray()) with:
722// - MExtractedSignalPix::GetExtractedSignalHiGain();
723// - MExtractedSignalPix::GetExtractedSignalLoGain();
724//
725// 2) Set number of saturated slices (MHCalibrationChargePix::AddSaturated()) with:
726// - MExtractedSignalPix::GetNumHiGainSaturated();
727// - MExtractedSignalPix::GetNumLoGainSaturated();
728//
729// 3) Fill AbsTime histograms (MHCalibrationChargePix::FillAbsTime()) with:
730// - MRawEvtPixelIter::GetIdxMaxHiGainSample();
731// - MRawEvtPixelIter::GetIdxMaxLoGainSample(first slice);
732//
733Bool_t MHCalibrationChargeCam::FillHists(const MParContainer *par, const Stat_t w)
734{
735
736 MExtractedSignalCam *signal = (MExtractedSignalCam*)par;
737 if (!signal)
738 {
739 *fLog << err << "No argument in MExtractedSignalCam::Fill... abort." << endl;
740 return kFALSE;
741 }
742
743 const UInt_t npixels = fGeom->GetNumPixels();
744 const UInt_t nareas = fGeom->GetNumAreas();
745 const UInt_t nsectors = fGeom->GetNumSectors();
746 const UInt_t lofirst = signal->GetFirstUsedSliceLoGain();
747
748 fSumhiarea .Reset();
749 fSumloarea .Reset();
750 fTimehiarea .Reset();
751 fTimeloarea .Reset();
752 fSumhisector.Reset();
753 fSumlosector.Reset();
754 fTimehisector.Reset();
755 fTimelosector.Reset();
756
757 fSathiarea .Reset();
758 fSatloarea .Reset();
759 fSathisector.Reset();
760 fSatlosector.Reset();
761
762 for (UInt_t i=0; i<npixels; i++)
763 {
764
765 MHCalibrationChargePix &histhi = (MHCalibrationChargePix&)(*this)[i];
766
767 if (histhi.IsExcluded())
768 continue;
769
770 const MExtractedSignalPix &pix = (*signal)[i];
771
772 const Float_t sumhi = pix.GetExtractedSignalHiGain();
773 const Int_t sathi = (Int_t)pix.GetNumHiGainSaturated();
774
775 if (IsOscillations())
776 histhi.FillHistAndArray(sumhi);
777 else
778 histhi.FillHist(sumhi);
779
780 histhi.AddSaturated(sathi);
781
782 const Int_t aidx = (*fGeom)[i].GetAidx();
783 const Int_t sector = (*fGeom)[i].GetSector();
784
785 fSumhiarea[aidx] += sumhi;
786 fSathiarea[aidx] += sathi;
787
788 fSumhisector[sector] += sumhi;
789 fSathisector[sector] += sathi;
790
791 if (IsLoGain())
792 {
793 MHCalibrationChargePix &histlo = (MHCalibrationChargePix&)(*this)(i);
794 const Float_t sumlo = pix.GetExtractedSignalLoGain();
795 const Int_t satlo = (Int_t)pix.GetNumLoGainSaturated();
796
797 if (IsOscillations())
798 histlo.FillHistAndArray(sumlo);
799 else
800 histlo.FillHist(sumlo);
801
802 histlo.AddSaturated(satlo);
803
804 fSumloarea[aidx] += sumlo;
805 fSatloarea[aidx] += satlo;
806 fSumlosector[sector] += sumlo;
807 fSatlosector[sector] += satlo;
808 }
809
810 }
811
812 MRawEvtPixelIter pixel(fRawEvt);
813 while (pixel.Next())
814 {
815
816 const UInt_t pixid = pixel.GetPixelId();
817
818 MHCalibrationChargePix &histhi = (MHCalibrationChargePix&)(*this)[pixid];
819
820 if (histhi.IsExcluded())
821 continue;
822
823 const Float_t timehi = (Float_t)pixel.GetIdxMaxHiGainSample();
824
825 histhi.FillAbsTime(timehi);
826
827 const Int_t aidx = (*fGeom)[pixid].GetAidx();
828 const Int_t sector = (*fGeom)[pixid].GetSector();
829
830 fTimehiarea [aidx] += timehi;
831 fTimehisector[sector] += timehi;
832
833 if (IsLoGain())
834 {
835 MHCalibrationChargePix &histlo = (MHCalibrationChargePix&)(*this)(pixid);
836
837 const Float_t timelo = (Float_t)pixel.GetIdxMaxLoGainSample(lofirst);
838 histlo.FillAbsTime(timelo);
839
840 fTimeloarea[aidx] += timelo;
841 fTimelosector[sector] += timelo;
842 }
843 }
844
845 for (UInt_t j=0; j<nareas; j++)
846 {
847
848 const Int_t npix = fAverageAreaNum[j];
849
850 if (npix == 0)
851 continue;
852
853 MHCalibrationChargePix &hipix = (MHCalibrationChargePix&)GetAverageHiGainArea(j);
854
855 if (IsOscillations())
856 hipix.FillHistAndArray(fSumhiarea [j]/npix);
857 else
858 hipix.FillHist(fSumhiarea[j]/npix);
859
860 hipix.AddSaturated ((Float_t)fSathiarea [j]/npix > 0.5 ? 1 : 0);
861 hipix.FillAbsTime (fTimehiarea[j]/npix);
862
863 if (IsLoGain())
864 {
865 MHCalibrationChargePix &lopix = (MHCalibrationChargePix&)GetAverageLoGainArea(j);
866 if (IsOscillations())
867 lopix.FillHistAndArray(fSumloarea [j]/npix);
868 else
869 lopix.FillHist(fSumloarea [j]/npix);
870 lopix.AddSaturated ((Float_t)fSatloarea [j]/npix > 0.5 ? 1 : 0);
871 lopix.FillAbsTime (fTimeloarea[j]/npix);
872 }
873 }
874
875 for (UInt_t j=0; j<nsectors; j++)
876 {
877
878 const Int_t npix = fAverageSectorNum[j];
879
880 if (npix == 0)
881 continue;
882
883 MHCalibrationChargePix &hipix = (MHCalibrationChargePix&)GetAverageHiGainSector(j);
884
885 if (IsOscillations())
886 hipix.FillHistAndArray(fSumhisector [j]/npix);
887 else
888 hipix.FillHist(fSumhisector [j]/npix);
889
890 hipix.AddSaturated ((Float_t)fSathisector[j]/npix > 0.5 ? 1 : 0);
891 hipix.FillAbsTime (fTimehisector[j]/npix);
892
893 if (IsLoGain())
894 {
895 MHCalibrationChargePix &lopix = (MHCalibrationChargePix&)GetAverageLoGainSector(j);
896
897 if (IsOscillations())
898 lopix.FillHistAndArray(fSumlosector [j]/npix);
899 else
900 lopix.FillHist(fSumlosector [j]/npix);
901
902 lopix.AddSaturated ((Float_t)fSatlosector[j]/npix > 0.5 ? 1 : 0);
903 lopix.FillAbsTime (fTimelosector[j]/npix);
904 }
905 }
906
907 return kTRUE;
908}
909
910// --------------------------------------------------------------------------
911//
912// For all TObjArray's (including the averaged ones), the following steps are performed:
913//
914// 1) Returns if the pixel is excluded.
915// 2) Tests saturation. In case yes, set the flag: MCalibrationPix::SetHiGainSaturation()
916// or the flag: MBadPixelsPix::SetUncalibrated( MBadPixelsPix::kLoGainSaturated )
917// 3) Store the absolute arrival times in the MCalibrationChargePix's. If flag
918// MCalibrationPix::IsHiGainSaturation() is set, the Low-Gain arrival times are stored,
919// otherwise the Hi-Gain ones.
920// 4) Calls to MHCalibrationCam::FitHiGainArrays() and MCalibrationCam::FitLoGainArrays()
921// with the flags:
922// - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::kHiGainNotFitted )
923// - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::kLoGainNotFitted )
924// - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::kHiGainOscillating )
925// - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::kLoGainOscillating )
926//
927Bool_t MHCalibrationChargeCam::FinalizeHists()
928{
929
930 *fLog << endl;
931
932 TH1F *h = NULL;
933
934 MCalibrationCam *chargecam = fIntensCam ? fIntensCam->GetCam() : fCam;
935 MBadPixelsCam *badcam = fIntensBad ? fIntensBad->GetCam() : fBadPixels;
936
937 for (Int_t i=0; i<fHiGainArray->GetSize(); i++)
938 {
939
940 MHCalibrationChargePix &histhi = (MHCalibrationChargePix&)(*this)[i];
941
942 if (histhi.IsExcluded())
943 continue;
944
945 MCalibrationChargePix &pix = (MCalibrationChargePix&)(*chargecam)[i] ;
946
947 if (histhi.GetSaturated() > fNumHiGainSaturationLimit*histhi.GetHGausHist()->GetEntries())
948 {
949 pix.SetHiGainSaturation();
950 if (IsOscillations())
951 histhi.CreateFourierSpectrum();
952 continue;
953 }
954
955 MBadPixelsPix &bad = (*badcam)[i];
956
957 h = histhi.GetHGausHist();
958
959 Stat_t overflow = h->GetBinContent(h->GetNbinsX()+1);
960 if (overflow > 0.1)
961 {
962 *fLog << warn << GetDescriptor()
963 << ": HiGain Histogram Overflow occurred " << overflow
964 << " times in pixel: " << i << " (without saturation!) " << endl;
965 bad.SetUncalibrated( MBadPixelsPix::kHiGainOverFlow );
966 }
967
968 overflow = h->GetBinContent(0);
969 if (overflow > 0.1)
970 {
971 *fLog << warn << GetDescriptor()
972 << ": HiGain Histogram Underflow occurred " << overflow
973 << " times in pixel: " << i << " (without saturation!) " << endl;
974 bad.SetUncalibrated( MBadPixelsPix::kHiGainOverFlow );
975 }
976
977 FinalizeAbsTimes(histhi, pix, bad, fFirstHiGain, fLastHiGain);
978 }
979
980 if (IsLoGain())
981 for (Int_t i=0; i<fLoGainArray->GetSize(); i++)
982 {
983
984 MHCalibrationChargePix &histlo = (MHCalibrationChargePix&)(*this)(i);
985 MBadPixelsPix &bad = (*badcam)[i];
986
987 if (histlo.IsExcluded())
988 continue;
989
990 if (histlo.GetSaturated() > fNumLoGainSaturationLimit*histlo.GetHGausHist()->GetEntries())
991 {
992 *fLog << warn << "Saturated Lo Gain histogram in pixel: " << i << endl;
993 bad.SetUncalibrated( MBadPixelsPix::kLoGainSaturation );
994 if (IsOscillations())
995 histlo.CreateFourierSpectrum();
996 continue;
997 }
998
999 h = histlo.GetHGausHist();
1000
1001 Stat_t overflow = h->GetBinContent(h->GetNbinsX()+1);
1002 if (overflow > 0.1)
1003 {
1004 *fLog << warn << GetDescriptor()
1005 << ": LoGain Histogram Overflow occurred " << overflow
1006 << " times in pixel: " << i << " (without saturation!) " << endl;
1007 bad.SetUncalibrated( MBadPixelsPix::kLoGainOverFlow );
1008 }
1009
1010 overflow = h->GetBinContent(0);
1011 if (overflow > 0.1)
1012 {
1013 *fLog << warn << GetDescriptor()
1014 << ": LoGain Histogram Underflow occurred " << overflow
1015 << " times in pixel: " << i << " (without saturation!) " << endl;
1016 bad.SetUncalibrated( MBadPixelsPix::kLoGainOverFlow );
1017 }
1018
1019 MCalibrationChargePix &pix = (MCalibrationChargePix&)(*chargecam)[i] ;
1020
1021 if (pix.IsHiGainSaturation())
1022 FinalizeAbsTimes(histlo, pix, bad, fFirstLoGain, fLastLoGain);
1023 }
1024
1025 for (Int_t j=0; j<fAverageHiGainAreas->GetSize(); j++)
1026 {
1027
1028 MHCalibrationChargePix &histhi = (MHCalibrationChargePix&)GetAverageHiGainArea(j);
1029 MCalibrationChargePix &pix = (MCalibrationChargePix&)chargecam->GetAverageArea(j);
1030
1031 if (histhi.GetSaturated() > fNumHiGainSaturationLimit*histhi.GetHGausHist()->GetEntries())
1032 {
1033 pix.SetHiGainSaturation();
1034 if (IsOscillations())
1035 histhi.CreateFourierSpectrum();
1036 continue;
1037 }
1038
1039 MBadPixelsPix &bad = chargecam->GetAverageBadArea(j);
1040 FinalizeAbsTimes(histhi, pix, bad, fFirstHiGain, fLastHiGain);
1041 }
1042
1043 if (IsLoGain())
1044 for (Int_t j=0; j<fAverageLoGainAreas->GetSize(); j++)
1045 {
1046
1047 MHCalibrationChargePix &histlo = (MHCalibrationChargePix&)GetAverageLoGainArea(j);
1048
1049 if (histlo.GetSaturated() > fNumLoGainSaturationLimit*histlo.GetHGausHist()->GetEntries())
1050 {
1051 *fLog << warn << "Saturated Lo Gain histogram in area idx: " << j << endl;
1052 histlo.CreateFourierSpectrum();
1053 continue;
1054 }
1055
1056 MCalibrationChargePix &pix = (MCalibrationChargePix&)chargecam->GetAverageArea(j);
1057
1058 if (pix.IsHiGainSaturation())
1059 {
1060 MBadPixelsPix &bad = chargecam->GetAverageBadArea(j);
1061 FinalizeAbsTimes(histlo, pix, bad, fFirstLoGain, fLastLoGain);
1062 }
1063
1064 }
1065
1066 for (Int_t j=0; j<fAverageHiGainSectors->GetSize(); j++)
1067 {
1068
1069 MHCalibrationChargePix &histhi = (MHCalibrationChargePix&)GetAverageHiGainSector(j);
1070 MCalibrationChargePix &pix = (MCalibrationChargePix&)chargecam->GetAverageSector(j);
1071
1072 if (histhi.GetSaturated() > fNumHiGainSaturationLimit*histhi.GetHGausHist()->GetEntries())
1073 {
1074 pix.SetHiGainSaturation();
1075 if (IsOscillations())
1076 histhi.CreateFourierSpectrum();
1077 continue;
1078 }
1079
1080 MBadPixelsPix &bad = chargecam->GetAverageBadSector(j);
1081 FinalizeAbsTimes(histhi, pix, bad, fFirstHiGain, fLastHiGain);
1082 }
1083
1084 if (IsLoGain())
1085 for (Int_t j=0; j<fAverageLoGainSectors->GetSize(); j++)
1086 {
1087
1088 MHCalibrationChargePix &histlo = (MHCalibrationChargePix&)GetAverageLoGainSector(j);
1089 MBadPixelsPix &bad = chargecam->GetAverageBadSector(j);
1090
1091 if (histlo.GetSaturated() > fNumLoGainSaturationLimit*histlo.GetHGausHist()->GetEntries())
1092 {
1093 *fLog << warn << "Saturated Lo Gain histogram in sector: " << j << endl;
1094 bad.SetUncalibrated( MBadPixelsPix::kLoGainSaturation );
1095 if (IsOscillations())
1096 histlo.CreateFourierSpectrum();
1097 continue;
1098 }
1099
1100 MCalibrationChargePix &pix = (MCalibrationChargePix&)chargecam->GetAverageSector(j);
1101
1102 if (pix.IsHiGainSaturation())
1103 FinalizeAbsTimes(histlo, pix, bad, fFirstLoGain, fLastLoGain);
1104 }
1105
1106 //
1107 // Perform the fitting for the High Gain (done in MHCalibrationCam)
1108 //
1109 FitHiGainArrays(*chargecam, *badcam,
1110 MBadPixelsPix::kHiGainNotFitted,
1111 MBadPixelsPix::kHiGainOscillating);
1112
1113 //
1114 // Perform the fitting for the Low Gain (done in MHCalibrationCam)
1115 //
1116 if (IsLoGain())
1117 FitLoGainArrays(*chargecam, *badcam,
1118 MBadPixelsPix::kLoGainNotFitted,
1119 MBadPixelsPix::kLoGainOscillating);
1120
1121 return kTRUE;
1122}
1123
1124// --------------------------------------------------------------------------------
1125//
1126// Fill the absolute time results into MCalibrationChargePix
1127//
1128// Check absolute time validity:
1129// - Mean arrival time is at least fTimeLowerLimit slices from the lower edge
1130// - Mean arrival time is at least fUpperLimit slices from the upper edge
1131//
1132void MHCalibrationChargeCam::FinalizeAbsTimes(MHCalibrationChargePix &hist, MCalibrationChargePix &pix, MBadPixelsPix &bad,
1133 Byte_t first, Byte_t last)
1134{
1135
1136 const Float_t mean = hist.GetAbsTimeMean();
1137 const Float_t rms = hist.GetAbsTimeRms();
1138
1139 pix.SetAbsTimeMean ( mean );
1140 pix.SetAbsTimeRms ( rms );
1141
1142 const Float_t lowerlimit = (Float_t)first + fTimeLowerLimit;
1143 const Float_t upperlimit = (Float_t)last + fTimeUpperLimit;
1144
1145 if ( mean < lowerlimit)
1146 {
1147 *fLog << warn << GetDescriptor()
1148 << Form("%s%3.1f%s%2.1f%s%4i",": Mean ArrivalTime: ",mean," smaller than ",fTimeLowerLimit,
1149 " FADC slices from lower edge in pixel ",hist.GetPixId()) << endl;
1150 bad.SetUncalibrated( MBadPixelsPix::kMeanTimeInFirstBin );
1151 }
1152
1153 if ( mean > upperlimit )
1154 {
1155 *fLog << warn << GetDescriptor()
1156 << Form("%s%3.1f%s%2.1f%s%4i",": Mean ArrivalTime: ",mean," greater than ",fTimeUpperLimit,
1157 " FADC slices from upper edge in pixel ",hist.GetPixId()) << endl;
1158 bad.SetUncalibrated( MBadPixelsPix::kMeanTimeInLast2Bins );
1159 }
1160}
1161
1162// --------------------------------------------------------------------------
1163//
1164// Sets all pixels to MBadPixelsPix::kUnsuitableRun, if following flags are set:
1165// - MBadPixelsPix::kLoGainSaturation
1166//
1167// Sets all pixels to MBadPixelsPix::kUnreliableRun, if following flags are set:
1168// - if MBadPixelsPix::kHiGainNotFitted and !MCalibrationPix::IsHiGainSaturation()
1169// - if MBadPixelsPix::kHiGainOscillating and !MCalibrationPix::IsHiGainSaturation()
1170// - if MBadPixelsPix::kLoGainNotFitted and MCalibrationPix::IsLoGainSaturation()
1171// - if MBadPixelsPix::kLoGainOscillating and MCalibrationPix::IsLoGainSaturation()
1172//
1173void MHCalibrationChargeCam::FinalizeBadPixels()
1174{
1175
1176 MBadPixelsCam *badcam = fIntensBad ? fIntensBad->GetCam() : fBadPixels;
1177 MCalibrationCam *chargecam = fIntensCam ? fIntensCam->GetCam() : fCam;
1178
1179 for (Int_t i=0; i<badcam->GetSize(); i++)
1180 {
1181
1182 MBadPixelsPix &bad = (*badcam)[i];
1183 MCalibrationPix &pix = (*chargecam)[i];
1184
1185 if (bad.IsUncalibrated( MBadPixelsPix::kHiGainNotFitted ))
1186 if (!pix.IsHiGainSaturation())
1187 bad.SetUnsuitable( MBadPixelsPix::kUnreliableRun );
1188
1189 if (bad.IsUncalibrated( MBadPixelsPix::kLoGainNotFitted ))
1190 if (pix.IsHiGainSaturation())
1191 bad.SetUnsuitable( MBadPixelsPix::kUnreliableRun );
1192
1193 if (bad.IsUncalibrated( MBadPixelsPix::kLoGainSaturation ))
1194 bad.SetUnsuitable( MBadPixelsPix::kUnsuitableRun );
1195
1196 if (IsOscillations())
1197 {
1198 if (bad.IsUncalibrated( MBadPixelsPix::kHiGainOscillating ))
1199 bad.SetUnsuitable( MBadPixelsPix::kUnreliableRun );
1200
1201 if (bad.IsUncalibrated( MBadPixelsPix::kLoGainOscillating ))
1202 if (pix.IsHiGainSaturation())
1203 bad.SetUnsuitable( MBadPixelsPix::kUnreliableRun );
1204 }
1205 }
1206
1207
1208
1209}
1210
1211// --------------------------------------------------------------------------
1212//
1213// Calls MHCalibrationPix::DrawClone() for pixel idx
1214//
1215void MHCalibrationChargeCam::DrawPixelContent(Int_t idx) const
1216{
1217 (*this)[idx].DrawClone();
1218}
1219
1220
1221// -----------------------------------------------------------------------------
1222//
1223// Default draw:
1224//
1225// Displays the averaged areas, both High Gain and Low Gain
1226//
1227// Calls the Draw of the fAverageHiGainAreas and fAverageLoGainAreas objects with options
1228//
1229void MHCalibrationChargeCam::Draw(const Option_t *opt)
1230{
1231
1232 const Int_t nareas = fAverageHiGainAreas->GetEntries();
1233 if (nareas == 0)
1234 return;
1235
1236 TString option(opt);
1237 option.ToLower();
1238
1239 if (!option.Contains("datacheck"))
1240 {
1241 MHCalibrationCam::Draw(opt);
1242 return;
1243 }
1244
1245 //
1246 // From here on , the datacheck - Draw
1247 //
1248 TVirtualPad *pad = gPad ? gPad : MH::MakeDefCanvas(this);
1249 pad->SetBorderMode(0);
1250 pad->Divide(1,nareas);
1251
1252 //
1253 // Loop over inner and outer pixels
1254 //
1255 for (Int_t i=0; i<nareas;i++)
1256 {
1257
1258 pad->cd(i+1);
1259
1260 MHCalibrationChargePix &hipix = (MHCalibrationChargePix&)GetAverageHiGainArea(i);
1261 //
1262 // Ask for Hi-Gain saturation
1263 //
1264 if (hipix.GetSaturated() > fNumHiGainSaturationLimit*hipix.GetHGausHist()->GetEntries() && IsLoGain())
1265 {
1266 MHCalibrationChargePix &lopix = (MHCalibrationChargePix&)GetAverageLoGainArea(i);
1267 DrawDataCheckPixel(lopix,i ? gkLoGainOuterRefLines : gkLoGainInnerRefLines);
1268 }
1269 else
1270 DrawDataCheckPixel(hipix,i ? gkHiGainOuterRefLines : gkHiGainInnerRefLines);
1271 }
1272}
1273
1274
1275// --------------------------------------------------------------------------
1276//
1277// Our own clone function is necessary since root 3.01/06 or Mars 0.4
1278// I don't know the reason.
1279//
1280// Creates new MHCalibrationChargeCam only for the Averaged Areas,
1281// the rest has to be retrieved directly, e.g. via:
1282// MHCalibrationChargeCam *cam = MParList::FindObject("MHCalibrationChargeCam");
1283// - cam->GetAverageSector(5).DrawClone();
1284// - (*cam)[100].DrawClone()
1285//
1286TObject *MHCalibrationChargeCam::Clone(const char *name) const
1287{
1288
1289 const Int_t navhi = fAverageHiGainAreas->GetEntries();
1290 // const Int_t nsehi = fAverageHiGainSectors->GetEntries();
1291 // const Int_t nselo = fAverageLoGainSectors->GetEntries();
1292
1293 //
1294 // FIXME, this might be done faster and more elegant, by direct copy.
1295 //
1296 MHCalibrationChargeCam *cam = new MHCalibrationChargeCam();
1297
1298 cam->fColor = fColor;
1299 cam->fPulserFrequency = fPulserFrequency;
1300 cam->fFlags = fFlags;
1301 cam->fRunNumbers = fRunNumbers;
1302 cam->fNbins = fNbins;
1303 cam->fFirst = fFirst;
1304 cam->fLast = fLast;
1305
1306 if (!IsAverageing())
1307 return cam;
1308
1309 cam->fAverageHiGainAreas->Expand(navhi);
1310
1311 for (int i=0; i<navhi; i++)
1312 (*cam->fAverageHiGainAreas) [i] = (*fAverageHiGainAreas) [i]->Clone();
1313 // cam->fAverageHiGainSectors->Expand(nsehi);
1314 // for (int i=0; i<nsehi; i++)
1315 // (*cam->fAverageHiGainSectors)[i] = (*fAverageHiGainSectors)[i]->Clone();
1316
1317 if (IsLoGain())
1318 {
1319 const Int_t navlo = fAverageLoGainAreas->GetEntries();
1320 cam->fAverageLoGainAreas->Expand(navlo);
1321
1322 for (int i=0; i<navlo; i++)
1323 (*cam->fAverageLoGainAreas) [i] = (*fAverageLoGainAreas) [i]->Clone();
1324 // cam->fAverageLoGainSectors->Expand(nselo);
1325 // for (int i=0; i<nselo; i++)
1326 // (*cam->fAverageLoGainSectors)[i] = (*fAverageLoGainSectors)[i]->Clone();
1327 }
1328
1329 cam->fAverageAreaNum = fAverageAreaNum;
1330 cam->fAverageAreaSat = fAverageAreaSat;
1331 cam->fAverageAreaSigma = fAverageAreaSigma;
1332 cam->fAverageAreaSigmaVar = fAverageAreaSigmaVar;
1333 cam->fAverageAreaRelSigma = fAverageAreaRelSigma;
1334 cam->fAverageAreaRelSigmaVar = fAverageAreaRelSigmaVar;
1335 cam->fAverageSectorNum = fAverageSectorNum;
1336
1337 return cam;
1338
1339}
1340
1341// -----------------------------------------------------------------------------
1342//
1343// Draw the average pixel for the datacheck:
1344//
1345// Displays the averaged areas, both High Gain and Low Gain
1346//
1347// Calls the Draw of the fAverageHiGainAreas and fAverageLoGainAreas objects with options
1348//
1349void MHCalibrationChargeCam::DrawDataCheckPixel(MHCalibrationChargePix &pix, const Float_t refline[])
1350{
1351
1352 TVirtualPad *newpad = gPad;
1353 newpad->Divide(1,2);
1354 newpad->cd(1);
1355
1356 gPad->SetTicks();
1357 if (!pix.IsEmpty() && !pix.IsOnlyOverflow() && !pix.IsOnlyUnderflow())
1358 gPad->SetLogy();
1359
1360 TH1F *hist = pix.GetHGausHist();
1361
1362 TH1F *null = new TH1F("Null",hist->GetTitle(),100,
1363 pix.GetFirst() > 0. ? pix.GetFirst() : 0.,
1364 pix.GetLast() > pix.GetFirst()
1365 ? ( pix.GetLast() > 450. ? 450. : pix.GetLast() )
1366 : pix.GetFirst()*2.);
1367
1368 null->SetMaximum(1.1*hist->GetMaximum());
1369 null->SetDirectory(NULL);
1370 null->SetBit(kCanDelete);
1371 null->SetStats(kFALSE);
1372 //
1373 // set the labels bigger
1374 //
1375 TAxis *xaxe = null->GetXaxis();
1376 TAxis *yaxe = null->GetYaxis();
1377 xaxe->CenterTitle();
1378 yaxe->CenterTitle();
1379 xaxe->SetTitleSize(0.07);
1380 yaxe->SetTitleSize(0.07);
1381 xaxe->SetTitleOffset(0.7);
1382 yaxe->SetTitleOffset(0.55);
1383 xaxe->SetLabelSize(0.06);
1384 yaxe->SetLabelSize(0.06);
1385 xaxe->SetTitle(hist->GetXaxis()->GetTitle());
1386 yaxe->SetTitle(hist->GetYaxis()->GetTitle());
1387
1388 null->Draw();
1389 hist->Draw("same");
1390
1391 gStyle->SetOptFit();
1392
1393 TF1 *fit = pix.GetFGausFit();
1394
1395 if (fit)
1396 {
1397 switch ( fColor )
1398 {
1399 case MCalibrationCam::kGREEN:
1400 fit->SetLineColor(kGreen);
1401 break;
1402 case MCalibrationCam::kBLUE:
1403 fit->SetLineColor(kBlue);
1404 break;
1405 case MCalibrationCam::kUV:
1406 fit->SetLineColor(106);
1407 break;
1408 case MCalibrationCam::kCT1:
1409 fit->SetLineColor(006);
1410 break;
1411 default:
1412 fit->SetLineColor(kRed);
1413 }
1414 fit->Draw("same");
1415 }
1416
1417 DisplayRefLines(null,refline);
1418
1419 newpad->cd(2);
1420 gPad->SetTicks();
1421
1422 TH1F *null2 = new TH1F("Null2",hist->GetTitle(),100,0.,pix.GetEvents()->GetSize()/pix.GetEventFrequency());
1423
1424 null2->SetMinimum(pix.GetMean()-10.*pix.GetSigma());
1425 null2->SetMaximum(pix.GetMean()+10.*pix.GetSigma());
1426 null2->SetDirectory(NULL);
1427 null2->SetBit(kCanDelete);
1428 null2->SetStats(kFALSE);
1429 //
1430 // set the labels bigger
1431 //
1432 TAxis *xaxe2 = null2->GetXaxis();
1433 TAxis *yaxe2 = null2->GetYaxis();
1434 xaxe2->CenterTitle();
1435 yaxe2->CenterTitle();
1436 xaxe2->SetTitleSize(0.07);
1437 yaxe2->SetTitleSize(0.07);
1438 xaxe2->SetTitleOffset(0.7);
1439 yaxe2->SetTitleOffset(0.55);
1440 xaxe2->SetLabelSize(0.06);
1441 yaxe2->SetLabelSize(0.06);
1442
1443 pix.CreateGraphEvents();
1444 TGraph *gr = pix.GetGraphEvents();
1445
1446 xaxe2->SetTitle(gr->GetXaxis()->GetTitle());
1447 yaxe2->SetTitle(gr->GetYaxis()->GetTitle());
1448
1449 null2->Draw();
1450
1451 pix.DrawEvents("same");
1452 return;
1453
1454}
1455
1456
1457void MHCalibrationChargeCam::DisplayRefLines(const TH1F *hist, const Float_t refline[]) const
1458{
1459
1460 TGraph *green1 = new TGraph(2);
1461 green1->SetPoint(0,refline[0],0.1);
1462 green1->SetPoint(1,refline[0],hist->GetMaximum());
1463 green1->SetBit(kCanDelete);
1464 green1->SetLineColor(kGreen);
1465 green1->SetLineStyle(2);
1466 green1->SetLineWidth(3);
1467 green1->Draw("L");
1468
1469 TGraph *green5 = new TGraph(2);
1470 green5->SetPoint(0,refline[6],0.1);
1471 green5->SetPoint(1,refline[6],hist->GetMaximum());
1472 green5->SetBit(kCanDelete);
1473 green5->SetLineColor(8);
1474 green5->SetLineStyle(2);
1475 green5->SetLineWidth(3);
1476 green5->Draw("L");
1477
1478 TGraph *blue1 = new TGraph(2);
1479 blue1->SetPoint(0,refline[1],0.1);
1480 blue1->SetPoint(1,refline[1],hist->GetMaximum());
1481 blue1->SetBit(kCanDelete);
1482 blue1->SetLineColor(227);
1483 blue1->SetLineStyle(2);
1484 blue1->SetLineWidth(3);
1485 blue1->Draw("L");
1486
1487 TGraph *blue5 = new TGraph(2);
1488 blue5->SetPoint(0,refline[2],0.1);
1489 blue5->SetPoint(1,refline[2],hist->GetMaximum());
1490 blue5->SetBit(kCanDelete);
1491 blue5->SetLineColor(68);
1492 blue5->SetLineStyle(2);
1493 blue5->SetLineWidth(3);
1494 blue5->Draw("L");
1495
1496 TGraph *blue10 = new TGraph(2);
1497 blue10->SetPoint(0,refline[3],0.1);
1498 blue10->SetPoint(1,refline[3],hist->GetMaximum());
1499 blue10->SetBit(kCanDelete);
1500 blue10->SetLineColor(4);
1501 blue10->SetLineStyle(2);
1502 blue10->SetLineWidth(3);
1503 blue10->Draw("L");
1504
1505 TGraph *uv10 = new TGraph(2);
1506 uv10->SetPoint(0,refline[4],0.1);
1507 uv10->SetPoint(1,refline[4],hist->GetMaximum());
1508 uv10->SetBit(kCanDelete);
1509 uv10->SetLineColor(106);
1510 uv10->SetLineStyle(2);
1511 uv10->SetLineWidth(3);
1512 uv10->Draw("L");
1513
1514 TGraph *ct1 = new TGraph(2);
1515 ct1->SetPoint(0,refline[5],0.1);
1516 ct1->SetPoint(1,refline[5],hist->GetMaximum());
1517 ct1->SetBit(kCanDelete);
1518 ct1->SetLineColor(6);
1519 ct1->SetLineStyle(2);
1520 ct1->SetLineWidth(3);
1521 ct1->Draw("L");
1522
1523 TLegend *leg = new TLegend(0.8,0.35,0.99,0.99);
1524 leg->SetBit(kCanDelete);
1525 leg->AddEntry(green1,"1 Led GREEN","l");
1526 leg->AddEntry(green5,"5 Leds GREEN","l");
1527 leg->AddEntry(blue1,"1 Led BLUE","l");
1528 leg->AddEntry(blue5,"5 Leds BLUE","l");
1529 leg->AddEntry(blue10,"10 Leds BLUE","l");
1530 leg->AddEntry(uv10,"10 Leds UV","l");
1531 leg->AddEntry(ct1,"CT1-Pulser","l");
1532
1533 leg->Draw();
1534}
1535
1536Int_t MHCalibrationChargeCam::ReadEnv(const TEnv &env, TString prefix, Bool_t print)
1537{
1538
1539 Bool_t rc = kFALSE;
1540
1541 if (MHCalibrationCam::ReadEnv(env,prefix,print))
1542 rc = kTRUE;
1543
1544 if (IsEnvDefined(env, prefix, "LoGainNbins", print))
1545 {
1546 SetLoGainNbins(GetEnvValue(env, prefix, "LoGainNbins", fLoGainNbins));
1547 rc = kTRUE;
1548 }
1549
1550 if (IsEnvDefined(env, prefix, "LoGainFirst", print))
1551 {
1552 SetLoGainFirst(GetEnvValue(env, prefix, "LoGainFirst", fLoGainFirst));
1553 rc = kTRUE;
1554 }
1555
1556 if (IsEnvDefined(env, prefix, "LoGainLast", print))
1557 {
1558 SetLoGainLast(GetEnvValue(env, prefix, "LoGainLast", fLoGainLast));
1559 rc = kTRUE;
1560 }
1561
1562 if (IsEnvDefined(env, prefix, "TimeLowerLimit", print))
1563 {
1564 SetTimeLowerLimit(GetEnvValue(env, prefix, "TimeLowerLimit", fTimeLowerLimit));
1565 rc = kTRUE;
1566 }
1567
1568 if (IsEnvDefined(env, prefix, "TimeUpperLimit", print))
1569 {
1570 SetTimeUpperLimit(GetEnvValue(env, prefix, "TimeUpperLimit", fTimeUpperLimit));
1571 rc = kTRUE;
1572 }
1573
1574
1575 return rc;
1576}
Note: See TracBrowser for help on using the repository browser.