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

Last change on this file since 3034 was 3014, checked in by gaug, 21 years ago
*** empty log message ***
File size: 11.0 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//
145const MHPedestalPixel &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 // fHArray->ForEach(TObject, Clear)();
158
159 fTotalEntries = 0;
160 fExtractSlices = 0.;
161}
162
163// --------------------------------------------------------------------------
164//
165// To setup the object we get the number of pixels from a MGeomCam object
166// in the Parameter list.
167// MPedestalPix sets its parameters to 0. (other than default which is -1.)
168//
169Bool_t MPedestalCam::SetupFill(const MParList *pList)
170{
171
172 fArray->ForEach(MPedestalPix, InitUseHists)();
173 fHArray->Delete();
174
175 return kTRUE;
176
177}
178
179// --------------------------------------------------------------------------
180//
181Bool_t MPedestalCam::Fill(const MParContainer *par, const Stat_t w)
182{
183
184 MExtractedSignalCam *signal = (MExtractedSignalCam*)par;
185 if (!signal)
186 {
187 gLog << err << "No argument in MExtractedSignalCam::Fill... abort." << endl;
188 return kFALSE;
189 }
190
191 Float_t slices = (Float_t)signal->GetNumUsedFADCSlices();
192
193 if (slices == 0.)
194 {
195 gLog << err << "Number of used signal slices in MExtractedSignalCam is zero ... abort."
196 << endl;
197 return kFALSE;
198 }
199
200 if (fExtractSlices != 0. && slices != fExtractSlices )
201 {
202 gLog << err << "Number of used signal slices changed in MExtractedSignalCam ... abort."
203 << endl;
204 return kFALSE;
205 }
206
207 fExtractSlices = slices;
208
209 const Int_t n = signal->GetSize();
210
211 if (fHArray->GetEntries()==0)
212 {
213 fHArray->Expand(n);
214
215 for (Int_t i=0; i<n; i++)
216 {
217 (*fHArray)[i] = new MHPedestalPixel;
218 MHPedestalPixel &hist = (*this)(i);
219 hist.ChangeHistId(i);
220 }
221 }
222
223 if (fHArray->GetEntries() != n)
224 {
225 gLog << err << "ERROR - Size mismatch... abort." << endl;
226 return kFALSE;
227 }
228
229 for (Int_t i=0; i<n; i++)
230 {
231 const MExtractedSignalPix &pix = (*signal)[i];
232
233 const Float_t sig = pix.GetExtractedSignalHiGain();
234 const Float_t sigPerSlice = sig/fExtractSlices;
235
236 MHPedestalPixel &hist = (*this)(i);
237
238 hist.FillCharge(sigPerSlice);
239 hist.FillChargevsN(sigPerSlice);
240 }
241
242 return kTRUE;
243}
244
245Bool_t MPedestalCam::Finalize()
246{
247 for (Int_t i=0; i<fHArray->GetSize(); i++)
248 {
249 MHPedestalPixel &hist = (*this)(i);
250
251 //
252 // 1) Return if the charge distribution is already succesfully fitted
253 // or if the histogram is empty
254 //
255 if (hist.IsFitOK() || hist.IsEmpty())
256 continue;
257
258 hist.CutAllEdges();
259
260 //
261 // 2) Fit the Hi Gain histograms with a Gaussian
262 //
263 hist.FitCharge();
264 }
265 return kTRUE;
266}
267
268
269void MPedestalCam::Print(Option_t *o) const
270{
271 *fLog << all << GetDescriptor() << ":" << endl;
272 int id = 0;
273
274 TIter Next(fArray);
275 MPedestalPix *pix;
276 while ((pix=(MPedestalPix*)Next()))
277 {
278 id++;
279
280 if (!pix->IsValid())
281 continue;
282
283 *fLog << id-1 << ": ";
284 *fLog << pix->GetPedestal() << " " << pix->GetPedestalRms() << endl;
285 }
286}
287
288Float_t MPedestalCam::GetPedestalMin(const MGeomCam *geom) const
289{
290 if (fArray->GetEntries() <= 0)
291 return 50.;
292
293 Float_t minval = (*this)[0].GetPedestalRms();
294
295 for (Int_t i=1; i<fArray->GetEntries(); i++)
296 {
297 const MPedestalPix &pix = (*this)[i];
298
299 Float_t testval = pix.GetPedestalRms();
300
301 if (geom)
302 testval *= geom->GetPixRatio(i);
303
304 if (testval < minval)
305 minval = testval;
306 }
307 return minval;
308}
309
310Float_t MPedestalCam::GetPedestalMax(const MGeomCam *geom) const
311{
312 if (fArray->GetEntries() <= 0)
313 return 50.;
314
315 Float_t maxval = (*this)[0].GetPedestalRms();
316
317 for (Int_t i=1; i<fArray->GetEntries(); i++)
318 {
319 const MPedestalPix &pix = (*this)[i];
320
321 Float_t testval = pix.GetPedestalRms();
322
323 if (geom)
324 testval *= geom->GetPixRatio(i);
325
326 if (testval > maxval)
327 maxval = testval;
328 }
329 return maxval;
330}
331
332Bool_t MPedestalCam::GetPixelContent(Double_t &val, Int_t idx, const MGeomCam &cam, Int_t type) const
333{
334
335 if (!(*this)[idx].IsValid())
336 return kFALSE;
337
338 if ((UInt_t)idx == gkBlindPixelId)
339 return kFALSE;
340
341 const Float_t ped = (*this)[idx].GetPedestal();
342 const Float_t rms = (*this)[idx].GetPedestalRms();
343
344 const Float_t pederr = rms/TMath::Sqrt((Float_t)fTotalEntries);
345 const Float_t rmserr = rms/TMath::Sqrt((Float_t)fTotalEntries/2.);
346
347 const MHPedestalPixel &hist = (*this)(idx);
348
349 const Float_t mean = hist.GetChargeMean();
350 const Float_t meanerr = hist.GetChargeMeanErr() * TMath::Sqrt((Float_t)fExtractSlices);
351 const Float_t sigma = hist.GetChargeSigma() * TMath::Sqrt((Float_t)fExtractSlices);
352 const Float_t sigmaerr = hist.GetChargeSigmaErr()* TMath::Sqrt((Float_t)fExtractSlices);
353 const Float_t prob = hist.GetChargeProb();
354
355 switch (type)
356 {
357 case 0:
358 val = ped;
359 break;
360 case 1:
361 val = pederr;
362 break;
363 case 2:
364 val = rms;
365 break;
366 case 3:
367 val = rmserr;
368 break;
369 case 4:
370 // if ((*this)[idx].IsFitValid())
371 val = mean;
372 // else
373 // return kFALSE;
374 break;
375 case 5:
376 // if ((*this)[idx].IsFitValid())
377 val = meanerr;
378 // else
379 // return kFALSE;
380 break;
381 case 6:
382 // if ((*this)[idx].IsFitValid())
383 val = sigma;
384 // else
385 // return kFALSE;
386 break;
387 case 7:
388 // if ((*this)[idx].IsFitValid())
389 val = sigmaerr;
390 // else
391 // return kFALSE;
392 break;
393 case 8:
394 // if ((*this)[idx].IsFitValid())
395 val = prob;
396 // else
397 // return kFALSE;
398 break;
399 case 9:
400 // if ((*this)[idx].IsFitValid())
401 val = 2.*(ped-mean)/(ped+mean);
402 // else
403 // return kFALSE;
404 break;
405 case 10:
406 val = TMath::Sqrt((pederr*pederr + meanerr*meanerr) * (ped*ped + mean*mean))
407 *2./(ped+mean)/(ped+mean);
408 break;
409 case 11:
410 val = 2.*(pederr - meanerr)/(pederr + meanerr);
411 break;
412 case 12:
413 val = 2.*(rms-sigma)/(rms+sigma);
414 break;
415 case 13:
416 val = TMath::Sqrt((rmserr*rmserr + sigmaerr*sigmaerr) * (rms*rms + sigma*sigma))
417 *2./(rms+sigma)/(rms+sigma);
418 break;
419 case 14:
420 // if ((*this)[idx].IsFitValid())
421 val = 2.*(rmserr - sigmaerr)/(rmserr + sigmaerr);
422 // else
423 // return kFALSE;
424 break;
425 default:
426 return kFALSE;
427 }
428 return kTRUE;
429}
430
431void MPedestalCam::DrawPixelContent(Int_t idx) const
432{
433 (*this)[idx].Draw();
434}
Note: See TracBrowser for help on using the repository browser.