source: trunk/MagicSoft/Mars/manalysis/MPedestalCam.cc@ 3051

Last change on this file since 3051 was 3048, checked in by gaug, 22 years ago
*** empty log message ***
File size: 12.1 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!
21! Copyright: MAGIC Software Development, 2000-2004
22!
23!
24\* ======================================================================== */
25
26/////////////////////////////////////////////////////////////////////////////
27// //
28// MPedestalCam //
29// //
30// Hold the Pedestal information for all pixels in the camera //
31// //
32/////////////////////////////////////////////////////////////////////////////
33#include "MPedestalCam.h"
34#include "MPedestalPix.h"
35
36#include <TClonesArray.h>
37
38#include "MLog.h"
39#include "MLogManip.h"
40
41#include "MParList.h"
42
43#include "MHExtractedSignalPix.h"
44
45#include "MExtractedSignalCam.h"
46#include "MExtractedSignalPix.h"
47
48#include "MGeomCam.h"
49
50ClassImp(MPedestalCam);
51
52using namespace std;
53
54const UInt_t MPedestalCam::gkBlindPixelId = 559;
55// --------------------------------------------------------------------------
56//
57// Default constructor. Creates a MPedestalPix object for each pixel
58//
59MPedestalCam::MPedestalCam(const char *name, const char *title)
60 : fTotalEntries(0), fExtractSlices(0.)
61{
62 fName = name ? name : "MPedestalCam";
63 fTitle = title ? title : "Storage container for all Pedestal Information in the camera";
64
65 fArray = new TClonesArray("MPedestalPix", 1);
66
67 //
68 // loop over all Pixels and create two histograms
69 // one for the Low and one for the High gain
70 // connect all the histogram with the container fHist
71 //
72 fHArray = new TObjArray;
73 fHArray->SetOwner();
74
75}
76
77// --------------------------------------------------------------------------
78//
79// Delete the array conatining the pixel pedest information
80//
81MPedestalCam::~MPedestalCam()
82{
83 delete fArray;
84 delete fHArray;
85}
86
87// --------------------------------------------------------------------------
88//
89// Set the size of the camera
90//
91void MPedestalCam::InitSize(const UInt_t i)
92{
93 fArray->ExpandCreate(i);
94}
95
96// --------------------------------------------------------------------------
97//
98// This function returns the current size of the TClonesArray
99// independently if the MPedestalPix is filled with values or not.
100//
101// Get the size of the MPedestalCam
102//
103Int_t MPedestalCam::GetSize() const
104{
105 return fArray->GetEntriesFast();
106}
107
108Int_t MPedestalCam::GetHistSize() const
109{
110 return fHArray->GetEntriesFast();
111}
112
113
114// --------------------------------------------------------------------------
115//
116// Get i-th pixel (pixel number)
117//
118MPedestalPix &MPedestalCam::operator[](Int_t i)
119{
120 return *static_cast<MPedestalPix*>(fArray->UncheckedAt(i));
121}
122
123// --------------------------------------------------------------------------
124//
125// Get i-th pixel (pixel number)
126//
127MPedestalPix &MPedestalCam::operator[](Int_t i) const
128{
129 return *static_cast<MPedestalPix*>(fArray->UncheckedAt(i));
130}
131
132// --------------------------------------------------------------------------
133//
134// Get i-th pixel (pixel number)
135//
136MHPedestalPixel &MPedestalCam::operator()(UInt_t i)
137{
138 return *static_cast<MHPedestalPixel*>(fHArray->UncheckedAt(i));
139}
140
141// --------------------------------------------------------------------------
142//
143// Get i-th pixel (pixel number)
144//
145MHPedestalPixel &MPedestalCam::operator()(UInt_t i) const
146{
147 return *static_cast<MHPedestalPixel*>(fHArray->UncheckedAt(i));
148}
149
150
151// -------------------------------------------------------------------------
152//
153void MPedestalCam::Clear(Option_t *o)
154{
155
156 fArray->ForEach(TObject, Clear)();
157
158 fTotalEntries = 0;
159 fExtractSlices = 0.;
160}
161
162
163// --------------------------------------------------------------------------
164//
165// Our own clone function is necessary since root 3.01/06 or Mars 0.4
166// I don't know the reason
167//
168TObject *MPedestalCam::Clone(const char *) const
169{
170
171 const Int_t n1 = fArray->GetSize();
172 const Int_t n2 = fHArray->GetSize();
173
174 //
175 // FIXME, this might be done faster and more elegant, by direct copy.
176 //
177 MPedestalCam *cam = new MPedestalCam;
178
179 cam->fArray->Expand(n1);
180 cam->fHArray->Expand(n2);
181
182 for (int i=0; i<n1; i++)
183 {
184 delete (*cam->fArray)[i];
185 (*cam->fArray)[i] = (*fArray)[i]->Clone();
186 }
187
188 for (int i=0; i<n2; i++)
189 {
190 delete (*cam->fHArray)[i];
191 (*cam->fHArray)[i] = (*fHArray)[i]->Clone();
192 }
193 return cam;
194}
195
196
197
198// --------------------------------------------------------------------------
199//
200// To setup the object we get the number of pixels from a MGeomCam object
201// in the Parameter list.
202// MPedestalPix sets its parameters to 0. (other than default which is -1.)
203//
204Bool_t MPedestalCam::SetupFill(const MParList *pList)
205{
206
207 fHArray->Delete();
208
209 return kTRUE;
210
211}
212
213// --------------------------------------------------------------------------
214//
215Bool_t MPedestalCam::Fill(const MParContainer *par, const Stat_t w)
216{
217
218 MExtractedSignalCam *signal = (MExtractedSignalCam*)par;
219 if (!signal)
220 {
221 gLog << err << "No argument in MExtractedSignalCam::Fill... abort." << endl;
222 return kFALSE;
223 }
224
225
226 Float_t slices = (Float_t)signal->GetNumUsedFADCSlices();
227
228 if (slices == 0.)
229 {
230 gLog << err << "Number of used signal slices in MExtractedSignalCam is zero ... abort."
231 << endl;
232 return kFALSE;
233 }
234
235 if (fExtractSlices != 0. && slices != fExtractSlices )
236 {
237 gLog << err << "Number of used signal slices changed in MExtractedSignalCam ... abort."
238 << endl;
239 return kFALSE;
240 }
241
242 fExtractSlices = slices;
243
244 const Int_t n = signal->GetSize();
245
246 if (fHArray->GetEntries()==0)
247 {
248 fHArray->Expand(n);
249
250 for (Int_t i=0; i<n; i++)
251 {
252 (*fHArray)[i] = new MHPedestalPixel;
253 MHPedestalPixel &hist = (*this)(i);
254 hist.ChangeHistId(i);
255 MPedestalPix &pix = (*this)[i];
256 pix.InitUseHists();
257 }
258 }
259
260 if (fHArray->GetEntries() != n)
261 {
262 gLog << err << "ERROR - Size mismatch... abort." << endl;
263 return kFALSE;
264 }
265
266 for (Int_t i=0; i<n; i++)
267 {
268
269 const MExtractedSignalPix &pix = (*signal)[i];
270
271 const Float_t sig = pix.GetExtractedSignalHiGain();
272
273 MHPedestalPixel &hist = (*this)(i);
274 //
275 // Don't fill signal per slice, we get completely screwed up
276 // with the sigma. Better fill like it is and renorm later
277 //
278 // const Float_t sigPerSlice = sig/fExtractSlices;
279 // hist.FillCharge(sigPerSlice);
280 // hist.FillChargevsN(sigPerSlice);
281 const Float_t signal = sig;
282 hist.FillCharge(signal);
283 hist.FillChargevsN(signal);
284 }
285
286 return kTRUE;
287}
288
289Bool_t MPedestalCam::Finalize()
290{
291 for (Int_t i=0; i<fHArray->GetSize(); i++)
292 {
293
294 MHPedestalPixel &hist = (*this)(i);
295
296 //
297 // 1) Return if the charge distribution is already succesfully fitted
298 // or if the histogram is empty
299 //
300 if (hist.IsFitOK() || hist.IsEmpty())
301 continue;
302
303 hist.CutAllEdges();
304
305 //
306 // 2) Fit the Hi Gain histograms with a Gaussian
307 //
308 hist.FitCharge();
309 hist.Renorm(fExtractSlices);
310
311 }
312 return kTRUE;
313}
314
315
316void MPedestalCam::Print(Option_t *o) const
317{
318 *fLog << all << GetDescriptor() << ":" << endl;
319 int id = 0;
320
321 TIter Next(fArray);
322 MPedestalPix *pix;
323 while ((pix=(MPedestalPix*)Next()))
324 {
325 id++;
326
327 if (!pix->IsValid())
328 continue;
329
330 *fLog << id-1 << ": ";
331 *fLog << pix->GetPedestal() << " " << pix->GetPedestalRms() << endl;
332 }
333}
334
335Float_t MPedestalCam::GetPedestalMin(const MGeomCam *geom) const
336{
337 if (fArray->GetEntries() <= 0)
338 return 50.;
339
340 Float_t minval = (*this)[0].GetPedestalRms();
341
342 for (Int_t i=1; i<fArray->GetEntries(); i++)
343 {
344 const MPedestalPix &pix = (*this)[i];
345
346 Float_t testval = pix.GetPedestalRms();
347
348 if (geom)
349 testval *= geom->GetPixRatio(i);
350
351 if (testval < minval)
352 minval = testval;
353 }
354 return minval;
355}
356
357Float_t MPedestalCam::GetPedestalMax(const MGeomCam *geom) const
358{
359 if (fArray->GetEntries() <= 0)
360 return 50.;
361
362 Float_t maxval = (*this)[0].GetPedestalRms();
363
364 for (Int_t i=1; i<fArray->GetEntries(); i++)
365 {
366 const MPedestalPix &pix = (*this)[i];
367
368 Float_t testval = pix.GetPedestalRms();
369
370 if (geom)
371 testval *= geom->GetPixRatio(i);
372
373 if (testval > maxval)
374 maxval = testval;
375 }
376 return maxval;
377}
378
379Bool_t MPedestalCam::GetPixelContent(Double_t &val, Int_t idx, const MGeomCam &cam, Int_t type) const
380{
381
382 if (GetSize() <= idx)
383 return kFALSE;
384
385 if (!(*this)[idx].IsValid())
386 return kFALSE;
387
388 if ((UInt_t)idx == gkBlindPixelId)
389 return kFALSE;
390
391 const Float_t ped = (*this)[idx].GetPedestal();
392 const Float_t rms = (*this)[idx].GetPedestalRms();
393
394 const Float_t pederr = rms/TMath::Sqrt((Float_t)fTotalEntries);
395 const Float_t rmserr = rms/TMath::Sqrt((Float_t)fTotalEntries*2.);
396
397 Float_t mean = 0.;
398 Float_t meanerr = 0.;
399 Float_t sigma = 0.;
400 Float_t sigmaerr = 0.;
401 Float_t prob = 0.;
402
403 if (type > 3)
404 if (GetHistSize() != 0)
405 {
406 if (!(*this)(idx).IsFitOK())
407 return kFALSE;
408
409 const MHPedestalPixel &hist = (*this)(idx);
410 mean = hist.GetChargeMean();
411 meanerr = hist.GetChargeMeanErr();
412 sigma = hist.GetChargeSigma() ;
413 sigmaerr = hist.GetChargeSigmaErr();
414 prob = hist.GetChargeProb();
415 }
416 else
417 return kFALSE;
418
419 switch (type)
420 {
421 case 0:
422 val = ped;
423 break;
424 case 1:
425 val = pederr;
426 break;
427 case 2:
428 val = rms;
429 break;
430 case 3:
431 val = rmserr;
432 break;
433 case 4:
434 if (!(*this)(idx).IsFitOK())
435 return kFALSE;
436 val = mean;
437 break;
438 case 5:
439 if (!(*this)(idx).IsFitOK())
440 return kFALSE;
441 val = meanerr;
442 break;
443 case 6:
444 if (!(*this)(idx).IsFitOK())
445 return kFALSE;
446 val = sigma;
447 break;
448 case 7:
449 if (!(*this)(idx).IsFitOK())
450 return kFALSE;
451 val = sigmaerr;
452 break;
453 case 8:
454 if (!(*this)(idx).IsFitOK())
455 return kFALSE;
456 val = prob;
457 break;
458 case 9:
459 if (!(*this)(idx).IsFitOK())
460 return kFALSE;
461 val = 2.*(ped-mean)/(ped+mean);
462 break;
463 case 10:
464 if (!(*this)(idx).IsFitOK())
465 return kFALSE;
466 val = TMath::Sqrt((pederr*pederr + meanerr*meanerr) * (ped*ped + mean*mean))
467 *2./(ped+mean)/(ped+mean);
468 break;
469 case 11:
470 if (!(*this)(idx).IsFitOK())
471 return kFALSE;
472 val = 2.*(pederr - meanerr)/(pederr + meanerr);
473 break;
474 case 12:
475 if (!(*this)(idx).IsFitOK())
476 return kFALSE;
477 val = 2.*(sigma-rms)/(sigma+rms);
478 break;
479 case 13:
480 if (!(*this)(idx).IsFitOK())
481 return kFALSE;
482 val = TMath::Sqrt((rmserr*rmserr + sigmaerr*sigmaerr) * (rms*rms + sigma*sigma))
483 *2./(rms+sigma)/(rms+sigma);
484 break;
485 case 14:
486 if (!(*this)(idx).IsFitOK())
487 return kFALSE;
488 val = 2.*(sigmaerr - rmserr)/(sigmaerr + rmserr);
489 break;
490 default:
491 return kFALSE;
492 }
493 return kTRUE;
494}
495
496void MPedestalCam::DrawPixelContent(Int_t idx) const
497{
498 (*this)(idx).Draw();
499}
Note: See TracBrowser for help on using the repository browser.