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

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