source: trunk/MagicSoft/Mars/mpedestal/MPedestalCam.cc@ 5332

Last change on this file since 5332 was 4794, checked in by gaug, 20 years ago
*** empty log message ***
File size: 16.8 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): Thomas Bretz 12/2000 <mailto:tbretz@uni-sw.gwdg.de>
19! Markus Gaug 02/2004 <mailto:markus@ifae.es>
20! Florian Goebel 06/2004 <mailto:fgoebel@mppmu.mpg.de>
21!
22! Copyright: MAGIC Software Development, 2000-2004
23!
24!
25\* ======================================================================== */
26
27/////////////////////////////////////////////////////////////////////////////
28// //
29// MPedestalCam //
30// //
31// Hold the Pedestal information for all pixels in the camera //
32// //
33/////////////////////////////////////////////////////////////////////////////
34#include "MPedestalCam.h"
35#include "MPedestalPix.h"
36
37#include <TArrayI.h>
38#include <TArrayF.h>
39#include <TArrayD.h>
40
41#include <TClonesArray.h>
42
43#include "MLog.h"
44#include "MLogManip.h"
45
46#include "MParList.h"
47
48#include "MGeomCam.h"
49#include "MGeomPix.h"
50
51#include "MBadPixelsCam.h"
52#include "MBadPixelsPix.h"
53
54ClassImp(MPedestalCam);
55
56using namespace std;
57
58// --------------------------------------------------------------------------
59//
60// Default constructor.
61//
62// Creates a TClonesArray of MPedestalPix containers, initialized to 1 entry, destinated
63// to hold one container per pixel. Later, a call to MPedestalCam::InitSize()
64// has to be performed (in MGeomApply).
65//
66// Creates a TClonesArray of MPedestalPix containers, initialized to 1 entry, destinated
67// to hold one container per pixel AREA. Later, a call to MPedestalCam::InitAreas()
68// has to be performed (in MGeomApply).
69//
70// Creates a TClonesArray of MPedestalPix containers, initialized to 1 entry, destinated
71// to hold one container per camera SECTOR. Later, a call to MPedestalCam::InitSectors()
72// has to be performed (in MGeomApply).
73//
74MPedestalCam::MPedestalCam(const char *name, const char *title)
75 : fTotalEntries(0)
76{
77 fName = name ? name : "MPedestalCam";
78 fTitle = title ? title : "Storage container for all Pedestal Information in the camera";
79
80 fArray = new TClonesArray("MPedestalPix", 1);
81 fAverageAreas = new TClonesArray("MPedestalPix", 1);
82 fAverageSectors = new TClonesArray("MPedestalPix", 1);
83}
84
85// --------------------------------------------------------------------------
86//
87// Deletes the following TClonesArray's of MPedestalPix containers (if exist):
88// - fArray
89// - fAverageAreas
90// - fAverageSectors
91//
92MPedestalCam::~MPedestalCam()
93{
94 delete fArray;
95 delete fAverageAreas;
96 delete fAverageSectors;
97}
98
99// --------------------------------------------------------------------------
100//
101// Set the size of the camera
102//
103void MPedestalCam::InitSize(const UInt_t i)
104{
105 fArray->ExpandCreate(i);
106}
107
108// -------------------------------------------------------------------
109//
110// Calls TClonesArray::ExpandCreate() for:
111// - fAverageAreas
112//
113void MPedestalCam::InitAverageAreas(const UInt_t i)
114{
115 fAverageAreas->ExpandCreate(i);
116}
117
118// -------------------------------------------------------------------
119//
120// Calls TClonesArray::ExpandCreate() for:
121// - fAverageSectors
122//
123void MPedestalCam::InitAverageSectors(const UInt_t i)
124{
125 fAverageSectors->ExpandCreate(i);
126}
127
128// -------------------------------------------------------------------
129//
130// Calls:
131// - InitSize()
132// - InitAverageAreas()
133// - InitAverageSectors()
134//
135void MPedestalCam::Init(const MGeomCam &geom)
136{
137 InitSize (geom.GetNumPixels() );
138 InitAverageAreas (geom.GetNumAreas() );
139 InitAverageSectors(geom.GetNumSectors());
140}
141
142// --------------------------------------------------------------------------
143//
144// This function returns the current size of the TClonesArray
145// independently if the MPedestalPix is filled with values or not.
146//
147// Get the size of the MPedestalCam
148//
149Int_t MPedestalCam::GetSize() const
150{
151 return fArray->GetEntriesFast();
152}
153
154// --------------------------------------------------------------------------
155//
156// Returns the current size of the TClonesArray fAverageAreas
157// independently if the MPedestalPix is filled with values or not.
158//
159const Int_t MPedestalCam::GetAverageAreas() const
160{
161 return fAverageAreas->GetEntriesFast();
162}
163
164// --------------------------------------------------------------------------
165//
166// Returns the current size of the TClonesArray fAverageSectors
167// independently if the MPedestalPix is filled with values or not.
168//
169const Int_t MPedestalCam::GetAverageSectors() const
170{
171 return fAverageSectors->GetEntriesFast();
172}
173
174// --------------------------------------------------------------------------
175//
176// Get i-th pixel (pixel number)
177//
178MPedestalPix &MPedestalCam::operator[](Int_t i)
179{
180 return *static_cast<MPedestalPix*>(fArray->UncheckedAt(i));
181}
182
183// --------------------------------------------------------------------------
184//
185// Get i-th pixel (pixel number)
186//
187const MPedestalPix &MPedestalCam::operator[](Int_t i) const
188{
189 return *static_cast<MPedestalPix*>(fArray->UncheckedAt(i));
190}
191
192// --------------------------------------------------------------------------
193//
194// Get i-th average pixel (area number)
195//
196MPedestalPix &MPedestalCam::GetAverageArea(UInt_t i)
197{
198 return *static_cast<MPedestalPix*>(fAverageAreas->UncheckedAt(i));
199}
200
201// --------------------------------------------------------------------------
202//
203// Get i-th average pixel (area number)
204//
205const MPedestalPix &MPedestalCam::GetAverageArea(UInt_t i) const
206{
207 return *static_cast<MPedestalPix*>(fAverageAreas->UncheckedAt(i));
208}
209
210// --------------------------------------------------------------------------
211//
212// Get i-th average pixel (sector number)
213//
214MPedestalPix &MPedestalCam::GetAverageSector(UInt_t i)
215{
216 return *static_cast<MPedestalPix*>(fAverageSectors->UncheckedAt(i));
217}
218
219// --------------------------------------------------------------------------
220//
221// Get i-th average pixel (sector number)
222//
223const MPedestalPix &MPedestalCam::GetAverageSector(UInt_t i) const
224{
225 return *static_cast<MPedestalPix*>(fAverageSectors->UncheckedAt(i));
226}
227
228// --------------------------------------
229//
230// Calls the ForEach macro for the TClonesArray fArray with the argument Clear()
231//
232// Loops over the fAverageAreas, calling the function Clear() for
233// every entry in fAverageAreas
234//
235// Loops over the fAverageSectors, calling the function Clear() for
236// every entry in fAverageSectors
237//
238void MPedestalCam::Clear(Option_t *o)
239{
240 fArray->ForEach(TObject, Clear)();
241
242 //
243 // another ForEach does not compile, thus have to do the loop ourselves:
244 //
245 for (Int_t i=0;i<GetAverageAreas();i++)
246 fAverageAreas[i].Clear();
247
248
249 //
250 // another ForEach does not compile, thus have to do the loop ourselves:
251 //
252 for (Int_t i=0;i<GetAverageSectors();i++)
253 fAverageSectors[i].Clear();
254
255 fTotalEntries = 0;
256}
257
258void MPedestalCam::Print(Option_t *o) const
259{
260 *fLog << all << GetDescriptor() << ":" << endl;
261 int id = 0;
262
263 TIter Next(fArray);
264 MPedestalPix *pix;
265 while ((pix=(MPedestalPix*)Next()))
266 {
267 id++;
268
269 if (!pix->IsValid())
270 continue;
271
272 *fLog << id-1 << ": ";
273 *fLog << pix->GetPedestal() << " " << pix->GetPedestalRms() << endl;
274 }
275}
276
277Float_t MPedestalCam::GetPedestalMin(const MGeomCam *geom) const
278{
279 if (fArray->GetEntries() <= 0)
280 return 50.;
281
282 Float_t minval = (*this)[0].GetPedestalRms();
283
284 for (Int_t i=1; i<fArray->GetEntries(); i++)
285 {
286 const MPedestalPix &pix = (*this)[i];
287
288 Float_t testval = pix.GetPedestalRms();
289
290 if (geom)
291 testval *= geom->GetPixRatio(i);
292
293 if (testval < minval)
294 minval = testval;
295 }
296 return minval;
297}
298
299Float_t MPedestalCam::GetPedestalMax(const MGeomCam *geom) const
300{
301 if (fArray->GetEntries() <= 0)
302 return 50.;
303
304 Float_t maxval = (*this)[0].GetPedestalRms();
305
306 for (Int_t i=1; i<fArray->GetEntries(); i++)
307 {
308 const MPedestalPix &pix = (*this)[i];
309
310 Float_t testval = pix.GetPedestalRms();
311
312 if (geom)
313 testval *= geom->GetPixRatio(i);
314
315 if (testval > maxval)
316 maxval = testval;
317 }
318 return maxval;
319}
320
321// --------------------------------------------------------------------------
322//
323// Calculates the average pedestal for pixel sizes.
324// The geometry container is used to get the necessary
325// geometry information (area index) If the bad pixel
326// container is given all pixels which have the flag 'kUnsuitableRun' are ignored
327// in the calculation of the size average.
328//
329// Returns a TArrayF of dimension 2:
330// arr[0]: averaged pedestal (default: -1.)
331// arr[1]: Error (rms) of averaged pedestal (default: 0.)
332//
333// ATTENTION: THE USER HAS TO DELETE THE RETURNED TARRAYF ACCORDINGLY
334//
335TArrayF *MPedestalCam::GetAveragedPedPerArea(const MGeomCam &geom, const UInt_t ai, MBadPixelsCam *bad)
336{
337
338 const Int_t np = GetSize();
339
340 Double_t mean = 0.;
341 Double_t mean2 = 0.;
342 Int_t nr = 0;
343
344 for (int i=0; i<np; i++)
345 {
346 if (bad && (*bad)[i].IsUnsuitable(MBadPixelsPix::kUnsuitableRun))
347 continue;
348
349 const UInt_t aidx = geom[i].GetAidx();
350
351 if (ai != aidx)
352 continue;
353
354 const MPedestalPix &pix = (*this)[i];
355
356 mean += pix.GetPedestal();
357 mean2 += pix.GetPedestal()*pix.GetPedestal();
358 nr ++;
359
360 }
361
362 TArrayF *arr = new TArrayF(2);
363 arr->AddAt(nr ? mean/nr : -1.,0);
364 arr->AddAt(nr>1 ? TMath::Sqrt((mean2 - mean*mean/nr)/(nr-1)) : 0. ,1);
365
366 return arr;
367}
368
369// --------------------------------------------------------------------------
370//
371// Calculates the average pedestal rms for pixel sizes.
372// The geometry container is used to get the necessary
373// geometry information (area index) If the bad pixel
374// container is given all pixels which have the flag 'kUnsuitableRun' are ignored
375// in the calculation of the size average.
376//
377// Returns a TArrayF of dimension 2:
378// arr[0]: averaged pedestal RMS (default: -1.)
379// arr[1]: Error (rms) of averaged pedestal RMS (default: 0.)
380//
381// ATTENTION: THE USER HAS TO DELETE THE RETURNED TARRAYF ACCORDINGLY
382//
383TArrayF *MPedestalCam::GetAveragedRmsPerArea(const MGeomCam &geom, const UInt_t ai, MBadPixelsCam *bad)
384{
385
386 const Int_t np = GetSize();
387
388 Double_t rms = 0.;
389 Double_t rms2 = 0.;
390 Int_t nr = 0;
391
392 for (int i=0; i<np; i++)
393 {
394 if (bad && (*bad)[i].IsUnsuitable(MBadPixelsPix::kUnsuitableRun))
395 continue;
396
397 const UInt_t aidx = geom[i].GetAidx();
398
399 if (ai != aidx)
400 continue;
401
402 const MPedestalPix &pix = (*this)[i];
403
404 rms += pix.GetPedestalRms();
405 rms2 += pix.GetPedestalRms()*pix.GetPedestalRms();
406 nr ++;
407 }
408
409 TArrayF *arr = new TArrayF(2);
410 arr->AddAt(nr ? rms/nr : -1.,0);
411 arr->AddAt(nr>1 ? TMath::Sqrt((rms2 - rms*rms/nr)/(nr-1)) : 0. ,1);
412
413 return arr;
414}
415
416// --------------------------------------------------------------------------
417//
418// Calculates the average pedestal for camera sectors.
419// The geometry container is used to get the necessary
420// geometry information (area index) If the bad pixel
421// container is given all pixels which have the flag 'kUnsuitableRun' are ignored
422// in the calculation of the size average.
423//
424// Returns a TArrayF of dimension 2:
425// arr[0]: averaged pedestal (default: -1.)
426// arr[1]: Error (rms) of averaged pedestal (default: 0.)
427//
428// ATTENTION: THE USER HAS TO DELETE THE RETURNED TARRAYF ACCORDINGLY
429//
430TArrayF *MPedestalCam::GetAveragedPedPerSector(const MGeomCam &geom, const UInt_t sec, MBadPixelsCam *bad)
431{
432
433 const Int_t np = GetSize();
434
435 Double_t mean = 0.;
436 Double_t mean2 = 0.;
437 Int_t nr = 0;
438
439 for (int i=0; i<np; i++)
440 {
441 if (bad && (*bad)[i].IsUnsuitable(MBadPixelsPix::kUnsuitableRun))
442 continue;
443
444 const UInt_t sector = geom[i].GetSector();
445
446 if (sec != sector)
447 continue;
448
449 const MPedestalPix &pix = (*this)[i];
450
451 mean += pix.GetPedestal();
452 mean2 += pix.GetPedestal()*pix.GetPedestal();
453 nr ++;
454
455 }
456
457 TArrayF *arr = new TArrayF(2);
458 arr->AddAt(nr ? mean/nr : -1.,0);
459 arr->AddAt(nr>1 ? TMath::Sqrt((mean2 - mean*mean/nr)/(nr-1)) : 0. ,1);
460
461 return arr;
462}
463
464// --------------------------------------------------------------------------
465//
466// Calculates the average pedestal rms for camera sectors.
467// The geometry container is used to get the necessary
468// geometry information (area index) If the bad pixel
469// container is given all pixels which have the flag 'kUnsuitableRun' are ignored
470// in the calculation of the size average.
471//
472// Returns a TArrayF of dimension 2:
473// arr[0]: averaged pedestal RMS (default: -1.)
474// arr[1]: Error (rms) of averaged pedestal RMS (default: 0.)
475//
476// ATTENTION: THE USER HAS TO DELETE THE RETURNED TARRAYF ACCORDINGLY
477//
478TArrayF *MPedestalCam::GetAveragedRmsPerSector(const MGeomCam &geom, const UInt_t sec, MBadPixelsCam *bad)
479{
480
481 const Int_t np = GetSize();
482
483 Double_t rms = 0.;
484 Double_t rms2 = 0.;
485 Int_t nr = 0;
486
487 for (int i=0; i<np; i++)
488 {
489 if (bad && (*bad)[i].IsUnsuitable(MBadPixelsPix::kUnsuitableRun))
490 continue;
491
492 const UInt_t sector = geom[i].GetSector();
493
494 if (sec != sector)
495 continue;
496
497 const MPedestalPix &pix = (*this)[i];
498
499 rms += pix.GetPedestalRms();
500 rms2 += pix.GetPedestalRms()*pix.GetPedestalRms();
501 nr ++;
502
503 }
504
505 TArrayF *arr = new TArrayF(2);
506 arr->AddAt(nr ? rms/nr : -1.,0);
507 arr->AddAt(nr>1 ? TMath::Sqrt((rms2 - rms*rms/nr)/(nr-1)) : 0. ,1);
508
509 return arr;
510
511}
512
513
514// --------------------------------------------------------------------------
515//
516// Calculates the avarage pedestal and pedestal rms for all sectors
517// and pixel sizes. The geometry container is used to get the necessary
518// geometry information (sector number, size index) If the bad pixel
519// container is given all pixels which have the flag 'kUnsuitableRun' are ignored
520// in the calculation of the sector and size average.
521//
522void MPedestalCam::ReCalc(const MGeomCam &geom, MBadPixelsCam *bad)
523{
524 const Int_t np = GetSize();
525 const Int_t ns = geom.GetNumSectors();
526 const Int_t na = geom.GetNumAreas();
527
528 TArrayI acnt(na);
529 TArrayI scnt(ns);
530 TArrayD asumx(na);
531 TArrayD ssumx(ns);
532 TArrayD asumr(na);
533 TArrayD ssumr(ns);
534
535 for (int i=0; i<np; i++)
536 {
537 if (bad && (*bad)[i].IsUnsuitable(MBadPixelsPix::kUnsuitableRun))
538 continue; //was: .IsBad()
539
540 // Create sums for areas and sectors
541 const UInt_t aidx = geom[i].GetAidx();
542 const UInt_t sect = geom[i].GetSector();
543
544 const MPedestalPix &pix = (*this)[i];
545
546 const UInt_t ne = pix.GetNumEvents();
547 const Float_t mean = pix.GetPedestal();
548 const Float_t rms = pix.GetPedestalRms();
549
550 asumx[aidx] += ne*mean;
551 asumr[aidx] += ne*rms;
552 acnt[aidx] += ne;
553
554 ssumx[sect] += ne*mean;
555 ssumr[sect] += ne*rms;
556 scnt[sect] += ne;
557 }
558
559 for (int i=0; i<ns; i++)
560 if (scnt[i]>0)
561 GetAverageSector(i).Set(ssumx[i]/scnt[i], ssumr[i]/scnt[i], 0, scnt[i]);
562 else
563 GetAverageSector(i).Clear();
564
565 for (int i=0; i<na; i++)
566 if (acnt[i]>0)
567 GetAverageArea(i).Set(asumx[i]/acnt[i], asumr[i]/acnt[i], 0, acnt[i]);
568 else
569 GetAverageArea(i).Clear();
570}
571
572Bool_t MPedestalCam::GetPixelContent(Double_t &val, Int_t idx, const MGeomCam &cam, Int_t type) const
573{
574 if (GetSize() <= idx)
575 return kFALSE;
576
577 if (!(*this)[idx].IsValid())
578 return kFALSE;
579
580 const Float_t ped = (*this)[idx].GetPedestal();
581 const Float_t rms = (*this)[idx].GetPedestalRms();
582
583 switch (type)
584 {
585 case 0:
586 val = ped;
587 break;
588 case 1:
589 val = fTotalEntries > 0 ?
590 rms/TMath::Sqrt((Float_t)fTotalEntries)
591 : (*this)[idx].GetPedestalError();
592 break;
593 case 2:
594 val = rms;
595 break;
596 case 3:
597 val = fTotalEntries > 0 ?
598 rms/TMath::Sqrt((Float_t)fTotalEntries*2.)
599 : (*this)[idx].GetPedestalRmsError();
600 break;
601 default:
602 return kFALSE;
603 }
604 return kTRUE;
605}
606
607void MPedestalCam::DrawPixelContent(Int_t idx) const
608{
609 *fLog << warn << "MPedestalCam::DrawPixelContent - not available." << endl;
610}
Note: See TracBrowser for help on using the repository browser.