source: trunk/MagicSoft/Mars/manalysis/MCalibrationCam.cc@ 2660

Last change on this file since 2660 was 2660, checked in by gaug, 22 years ago
*** empty log message ***
File size: 14.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!
19! Author(s): Markus Gaug 11/2003 <mailto:markus@ifae.es>
20!
21! Copyright: MAGIC Software Development, 2000-2001
22!
23!
24\* ======================================================================== */
25
26/////////////////////////////////////////////////////////////////////////////
27// //
28// MCalibrationCam //
29// //
30// Hold the Calibration information for all pixels in the camera //
31// //
32// The Classes MCalibrationPix's are stored in a TClonesArray //
33// The Class MCalibrationBlindPix and the MCalibrationPINDiode //
34// are stored in separate pointers //
35// //
36/////////////////////////////////////////////////////////////////////////////
37#include "MCalibrationCam.h"
38#include "MCalibrationPix.h"
39#include "MCalibrationBlindPix.h"
40#include "MCalibrationConfig.h"
41
42#include <TClonesArray.h>
43
44#include "MLog.h"
45#include "MLogManip.h"
46
47#include "MGeomCam.h"
48
49ClassImp(MCalibrationCam);
50
51using namespace std;
52// --------------------------------------------------------------------------
53//
54// Default constructor.
55// Creates a TClonesArray of MHCalibrationPix objects, initialized to 0 entries
56// Creates an MCalibrationBlindPix and an MCalibrationPINDiode
57//
58MCalibrationCam::MCalibrationCam(const char *name, const char *title)
59 : fMeanNrPhotAvailable(kFALSE),
60 fMeanNrPhotInnerPix(-1.),
61 fMeanNrPhotInnerPixErr(-1.)
62{
63 fName = name ? name : "MCalibrationCam";
64 fTitle = title ? title : "Storage container for the Calibration Information in the camera";
65
66 fPixels = new TClonesArray("MCalibrationPix",1);
67 fBlindPixel = new MCalibrationBlindPix();
68 fPINDiode = new MCalibrationPINDiode();
69}
70
71// --------------------------------------------------------------------------
72//
73// Delete the TClonesArray or MCalibrationPix's
74// Delete the MCalibrationPINDiode and the MCalibrationBlindPix
75//
76MCalibrationCam::~MCalibrationCam()
77{
78
79 //
80 // delete fPixels should delete all Objects stored inside
81 //
82 delete fPixels;
83 delete fBlindPixel;
84 delete fPINDiode;
85}
86
87// --------------------------------------------------------------------------
88//
89// This function return the size of the FILLED MCalibrationCam
90// It is NOT the size of the array fPixels !!!
91// Note that with every call to AddPixels, GetSize() might change
92//
93Int_t MCalibrationCam::GetSize() const
94{
95
96 //
97 // Here we get the number of "filled" fPixels!!
98 //
99 return fPixels->GetEntriesFast();
100}
101
102// --------------------------------------------------------------------------
103//
104// Get i-th pixel (pixel number)
105//
106MCalibrationPix &MCalibrationCam::operator[](Int_t i)
107{
108 return *static_cast<MCalibrationPix*>(fPixels->UncheckedAt(i));
109}
110
111// --------------------------------------------------------------------------
112//
113// Get i-th pixel (pixel number)
114//
115MCalibrationPix &MCalibrationCam::operator[](Int_t i) const
116{
117 return *static_cast<MCalibrationPix*>(fPixels->UncheckedAt(i));
118}
119
120// --------------------------------------------------------------------------
121//
122// Return a pointer to the pixel with the requested idx.
123// NULL if it doesn't exist.
124//
125MCalibrationPix *MCalibrationCam::GetCalibrationPix(Int_t idx) const
126{
127 if (idx<0)
128 return NULL;
129
130 if (!CheckBounds(idx))
131 return NULL;
132
133 return (MCalibrationPix*)fPixels->At(idx);
134}
135
136
137
138// --------------------------------------------------------------------------
139//
140// Check if position i is inside bounds
141//
142Bool_t MCalibrationCam::CheckBounds(Int_t i) const
143{
144 return i < fPixels->GetEntriesFast();
145}
146
147Bool_t MCalibrationCam::IsPixelUsed(Int_t idx) const
148{
149 if (!CheckBounds(idx))
150 return kFALSE;
151
152 return kTRUE;
153}
154
155Bool_t MCalibrationCam::IsPixelFitted(Int_t idx) const
156{
157 return ((*this)[idx].GetCharge() > 0. && (*this)[idx].GetErrCharge() > 0.);
158}
159
160
161
162void MCalibrationCam::Clear(Option_t *o)
163{
164 fPixels->ForEach(TObject, Clear)();
165}
166
167//
168// Perform the fits of the charges
169// If i=-1, then all pixels will be fitted
170// Otherwise only the one with index i
171//
172// The number of succesful fits is returned
173//
174UShort_t MCalibrationCam::FitCharge(Int_t i)
175{
176
177 UShort_t nsuccess = 0;
178
179 // invalidate i if it exceeds the number of entries in fPixels
180 if (i > fPixels->GetEntriesFast())
181 {
182 *fLog << warn << "Tried to fit pixel out of allowed range " << endl;
183 return 0;
184 }
185
186 if (i == -1) // loop over all events
187 {
188
189 TIter Next(fPixels);
190 MCalibrationPix *pix;
191 while ((pix=(MCalibrationPix*)Next()))
192 {
193 if (pix->FitCharge())
194 nsuccess++;
195 }
196 }
197 else // fit only the pixel with index i
198 {
199 if((*this)[i].FitCharge())
200 nsuccess++;
201 }
202
203 return nsuccess;
204
205}
206
207//
208// Perform the fits of the charges of all pixels
209// plus the blind pixel and the PIN Diode
210//
211// The number of succesful fits is returned
212//
213UShort_t MCalibrationCam::FitAllCharge()
214{
215
216 // FIXME: Once the blind pixel is fully working,
217 // there must be some penalty in case the
218 // fit does not succeed
219 //
220 UShort_t nsuccess = 0;
221
222 TIter Next(fPixels);
223 MCalibrationPix *pix;
224 while ((pix=(MCalibrationPix*)Next()))
225 {
226 if (pix->FitCharge())
227 nsuccess++;
228 }
229
230 if (fBlindPixel->FitCharge())
231 nsuccess++;
232
233 if (fPINDiode->FitCharge())
234 nsuccess++;
235
236 return nsuccess;
237
238}
239
240
241
242//
243// Perform the fits of the arrival times
244// If i=-1, then all pixels will be fitted
245// Otherwise only the one with index i
246//
247// The number of succesful fits is returned
248//
249UShort_t MCalibrationCam::FitTime(Int_t i)
250{
251
252 UShort_t nsuccess = 0;
253
254 // invalidate i if it exceeds the number of entries in fPixels
255 if (i > fPixels->GetEntriesFast())
256 {
257 *fLog << warn << "Tried to fit pixel out of allowed range " << endl;
258 return 0;
259 }
260
261 if (i == -1) // loop over all events
262 {
263
264 TIter Next(fPixels);
265 MCalibrationPix *pix;
266 while ((pix=(MCalibrationPix*)Next()))
267 {
268 if (pix->FitTime())
269 nsuccess++;
270 }
271 }
272 else // fit only the pixel with index i
273 {
274 if((*this)[i].FitTime())
275 nsuccess++;
276 }
277
278 return nsuccess;
279
280}
281
282
283//
284// Perform the fits of the times of all pixels
285// plus the blind pixel and the PIN Diode
286//
287// The number of succesful fits is returned
288//
289UShort_t MCalibrationCam::FitAllTime()
290{
291
292 // FIXME: Once the blind pixel is fully working,
293 // there must be some penalty in case the
294 // fit does not succeed
295 //
296 UShort_t nsuccess = 0;
297
298 TIter Next(fPixels);
299 MCalibrationPix *pix;
300 while ((pix=(MCalibrationPix*)Next()))
301 {
302 if (pix->FitTime())
303 nsuccess++;
304 }
305
306 if (fBlindPixel->FitTime())
307 nsuccess++;
308
309 if (fPINDiode->FitTime())
310 nsuccess++;
311
312 return nsuccess;
313
314}
315
316
317
318void MCalibrationCam::CutEdges()
319{
320
321 fBlindPixel->GetHist()->CutAllEdges();
322
323 TIter Next(fPixels);
324 MCalibrationPix *pix;
325 while ((pix=(MCalibrationPix*)Next()))
326 {
327 pix->GetHist()->CutAllEdges();
328 }
329
330 return;
331}
332
333// ---------------------------------------------------------------------
334//
335// This function simply allocates memory via the ROOT command:
336// (TObject**) TStorage::ReAlloc(fCont, newSize * sizeof(TObject*),
337// fSize * sizeof(TObject*));
338// newSize corresponds to size in our case
339// fSize is the old size (in most cases: 0)
340//
341void MCalibrationCam::InitSize(Int_t size)
342{
343
344 //
345 // check if we have already initialized to size
346 //
347 if (CheckBounds(size))
348 return;
349
350 //
351 // It is important to use Expand and NOT ExpandCreate, because
352 // we want to keep all pixel not used with a NULL pointer.
353 //
354 fPixels->ExpandCreate(size);
355
356}
357
358void MCalibrationCam::Print(Option_t *o) const
359{
360 *fLog << all << GetDescriptor() << ":" << endl;
361 int id = 0;
362
363 *fLog << "Succesfully calibrated pixels:" << endl;
364 *fLog << endl;
365
366 TIter Next(fPixels);
367 MCalibrationPix *pix;
368 while ((pix=(MCalibrationPix*)Next()))
369 {
370
371 if (pix->GetCharge() >= 0.)
372 {
373 *fLog << pix->GetPixId() << " Pedestals: " << pix->GetPed() << " +- " << pix->GetPedRms()
374 << " Reduced Charge: " << pix->GetCharge() << " +- "
375 << pix->GetSigmaCharge() << " Reduced Sigma: " << TMath::Sqrt(pix->GetRSigmaSquare()) << endl;
376 id++;
377 }
378 }
379
380 *fLog << id << " succesful pixels :-))" << endl;
381 id = 0;
382
383 *fLog << endl;
384 *fLog << "Pixels with errors:" << endl;
385 *fLog << endl;
386
387 TIter Next2(fPixels);
388 while ((pix=(MCalibrationPix*)Next2()))
389 {
390
391 if (pix->GetCharge() == -1.)
392 {
393 *fLog << pix->GetPixId() << " Pedestals: " << pix->GetPed() << " +- " << pix->GetPedRms()
394 << " Reduced Charge: " << pix->GetCharge() << " +- "
395 << pix->GetSigmaCharge() << " Reduced Sigma: " << TMath::Sqrt(pix->GetRSigmaSquare()) << endl;
396 id++;
397 }
398 }
399 *fLog << id << " pixels with errors :-((" << endl;
400
401}
402
403
404Bool_t MCalibrationCam::GetPixelContent(Double_t &val, Int_t idx, const MGeomCam &cam, Int_t type) const
405{
406
407 if (idx > GetSize())
408 return kFALSE;
409
410 switch (type)
411 {
412 case 0:
413 val = (*this)[idx].GetCharge();
414 break;
415 case 1:
416 val = (*this)[idx].GetErrCharge();
417 break;
418 case 2:
419 val = (*this)[idx].GetSigmaCharge();
420 break;
421 case 3:
422 val = (*this)[idx].GetErrSigmaCharge();
423 break;
424 case 4:
425 val = (*this)[idx].GetChargeProb();
426 break;
427 case 5:
428 val = (*this)[idx].GetTime();
429 break;
430 case 6:
431 val = (*this)[idx].GetSigmaTime();
432 break;
433 case 7:
434 val = (*this)[idx].GetTimeChiSquare();
435 break;
436 case 8:
437 val = (*this)[idx].GetPed();
438 break;
439 case 9:
440 val = (*this)[idx].GetPedRms();
441 break;
442 case 10:
443 if ((*this)[idx].GetRSigmaSquare() > 0.)
444 val = TMath::Sqrt((*this)[idx].GetRSigmaSquare());
445 else
446 val = -1.;
447 break;
448 case 15:
449 val = ((*this)[idx].GetSigmaCharge()/(*this)[idx].GetCharge())*
450 ((*this)[idx].GetSigmaCharge()/(*this)[idx].GetCharge());
451 break;
452 case 11:
453 val = (*this)[idx].GetPheFFactorMethod();
454 break;
455 case 12:
456 val = (*this)[idx].GetConversionFFactorMethod();
457 break;
458 case 13:
459 if (idx < 397)
460 val = (double)fMeanNrPhotInnerPix;
461 else
462 val = (double)fMeanNrPhotInnerPix*gkCalibrationOuterPixelArea;
463 break;
464 case 14:
465 if ((fMeanNrPhotInnerPix > 0. ) && ((*this)[idx].GetCharge() != -1.))
466 {
467 if (idx < 397)
468 val = fMeanNrPhotInnerPix / (*this)[idx].GetCharge();
469 else
470 val = fMeanNrPhotInnerPix*gkCalibrationOuterPixelArea / (*this)[idx].GetCharge();
471 }
472 else
473 {
474 val = -1.;
475 }
476 break;
477 default:
478 return kFALSE;
479 }
480 return val>=0;
481}
482
483void MCalibrationCam::DrawPixelContent(Int_t idx) const
484{
485
486 (*this)[idx].Draw();
487
488}
489
490Bool_t MCalibrationCam::CalcNrPhotInnerPixel()
491{
492 if (!fBlindPixel->IsFitOK())
493 return kFALSE;
494
495 const Float_t mean = fBlindPixel->GetLambda();
496 // const Float_t merr = fBlindPixel->GetErrLambda();
497
498 switch (fColor)
499 {
500 case kECGreen:
501 fMeanNrPhotInnerPix = (mean / gkCalibrationBlindPixelQEGreen) // real photons
502 *TMath::Power(10,gkCalibrationBlindPixelAttGreen) // correct for absorption
503 / gkCalibrationBlindPixelArea; // correct for area
504 break;
505 case kECBlue:
506 fMeanNrPhotInnerPix = (mean / gkCalibrationBlindPixelQEBlue )
507 *TMath::Power(10,gkCalibrationBlindPixelAttBlue)
508 / gkCalibrationBlindPixelArea;
509 break;
510 case kECUV:
511 fMeanNrPhotInnerPix = (mean / gkCalibrationBlindPixelQEUV )
512 *TMath::Power(10,gkCalibrationBlindPixelAttUV)
513 / gkCalibrationBlindPixelArea;
514 break;
515 case kECCT1:
516 default:
517 fMeanNrPhotInnerPix = (mean / gkCalibrationBlindPixelQECT1 )
518 *TMath::Power(10,gkCalibrationBlindPixelAttCT1)
519 / gkCalibrationBlindPixelArea;
520 break;
521 }
522
523 fMeanNrPhotAvailable = kTRUE;
524 return kTRUE;
525}
526
527
528Bool_t MCalibrationCam::GetConversionFactorBlindPixel(Int_t ipx, Float_t &mean, Float_t &err, Float_t &sigma)
529{
530
531 if (ipx < 0 || !IsPixelFitted(ipx))
532 return kFALSE;
533
534 if (!fMeanNrPhotAvailable)
535 if (!CalcNrPhotInnerPixel())
536 return kFALSE;
537
538 mean = fMeanNrPhotInnerPix / (*this)[ipx].GetCharge();
539
540 //
541 // Not yet ready , sorry
542 //
543 err = -1.;
544 sigma = -1.;
545
546 return kTRUE;
547}
548
549
550Bool_t MCalibrationCam::GetConversionFactorFFactor(Int_t ipx, Float_t &mean, Float_t &err, Float_t &sigma)
551{
552
553 if (ipx < 0 || !IsPixelFitted(ipx))
554 return kFALSE;
555
556 Float_t conv = (*this)[ipx].GetConversionFFactorMethod();
557
558 if (conv < 0.)
559 return kFALSE;
560
561 mean = conv;
562
563 //
564 // Not yet ready , sorry
565 //
566 err = -1.;
567 sigma = -1.;
568
569 return kTRUE;
570}
571
572
573Bool_t MCalibrationCam::GetConversionFactorPINDiode(Int_t ipx, Float_t &mean, Float_t &err, Float_t &sigma)
574{
575
576 return kTRUE;
577
578}
579
580Bool_t MCalibrationCam::GetConversionFactorCombined(Int_t ipx, Float_t &mean, Float_t &err, Float_t &sigma)
581{
582
583 return kTRUE;
584
585}
Note: See TracBrowser for help on using the repository browser.