source: trunk/MagicSoft/Mars/mcalib/MHCalibrationChargeCam.cc@ 4900

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