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

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