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

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