source: trunk/MagicSoft/Mars/mcalib/MCalibrationRelTimeCalc.cc@ 8040

Last change on this file since 8040 was 8015, checked in by tbretz, 18 years ago
*** empty log message ***
File size: 17.3 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): Markus Gaug 04/2004 <mailto:markus@ifae.es>
19!
20! Copyright: MAGIC Software Development, 2000-2004
21!
22!
23\* ======================================================================== */
24
25//////////////////////////////////////////////////////////////////////////////
26//
27// MCalibrationRelTimeCalc
28//
29// Task to finalize the relative time calibration obtained
30// from the fit results to the summed FADC slice distributions delivered by
31// MCalibrationRelTimeCam, calculated and filled by MHCalibrationRelTimeCam,
32//
33// PreProcess(): Initialize pointers to MCalibrationRelTimeCam,
34//
35// ReInit(): Initializes pointer to MBadPixelsCam
36//
37// Process(): Nothing to be done, histograms getting filled by MHCalibrationChargeCam
38//
39// PostProcess(): - FinalizeRelTimes()
40// - FinalizeBadPixels()
41//
42// Class Version 2:
43// + Byte_t fCheckFlags; // Bit-field to hold the possible check flags
44//
45//
46// Input Containers:
47// MCalibrationRelTimeCam
48// MBadPixelsCam
49// MGeomCam
50//
51// Output Containers:
52// MCalibrationRelTimeCam
53// MBadPixelsCam
54//
55//
56//////////////////////////////////////////////////////////////////////////////
57#include "MCalibrationRelTimeCalc.h"
58
59#include "MLog.h"
60#include "MLogManip.h"
61
62#include "MParList.h"
63
64#include "MStatusDisplay.h"
65
66#include "MGeomCam.h"
67#include "MGeomPix.h"
68
69#include "MCalibrationIntensityRelTimeCam.h"
70#include "MCalibrationRelTimeCam.h"
71#include "MCalibrationRelTimePix.h"
72
73#include "MBadPixelsIntensityCam.h"
74#include "MBadPixelsCam.h"
75#include "MBadPixelsPix.h"
76
77ClassImp(MCalibrationRelTimeCalc);
78
79using namespace std;
80
81const Float_t MCalibrationRelTimeCalc::fgRelTimeResolutionLimit = 1.0;
82
83// --------------------------------------------------------------------------
84//
85// Default constructor.
86//
87// Sets all pointers to NULL
88//
89// Initializes:
90// - fRelTimeResolutionLimit to fgRelTimeResolutionimit
91// - fOutputPath to "."
92// - fOutputFile to "TimeCalibStat.txt"
93//
94// Calls:
95// - Clear()
96//
97MCalibrationRelTimeCalc::MCalibrationRelTimeCalc(const char *name, const char *title)
98 : fGeom(NULL), fFlags(0)
99{
100
101 fName = name ? name : "MCalibrationRelTimeCalc";
102 fTitle = title ? title : "Task to finalize the relative time calibration";
103
104 SetCheckFitResults ( kFALSE );
105 SetCheckDeviatingBehavior( kFALSE );
106 SetCheckHistOverflow ( kFALSE );
107 SetCheckOscillations ( kFALSE );
108
109 SetRelTimeResolutionLimit();
110 SetOutputPath();
111 SetOutputFile("");
112
113 Clear();
114
115}
116
117// --------------------------------------------------------------------------
118//
119// Sets:
120// - all flags to kFALSE
121// - all pointers to NULL
122//
123void MCalibrationRelTimeCalc::Clear(const Option_t *o)
124{
125
126 fIntensBad = NULL;
127 fBadPixels = NULL;
128 fCam = NULL;
129 fIntensCam = NULL;
130
131}
132
133// --------------------------------------------------------------------------
134//
135// Search for the following input containers and abort if not existing:
136// - MGeomCam
137// - MCalibrationIntensityRelTimeCam or MCalibrationRelTimeCam
138// - MBadPixelsIntensityCam or MBadPixelsCam
139//
140// It defines the PixId of every pixel in:
141//
142// - MCalibrationRelTimeCam
143//
144// It sets all pixels in excluded which have the flag fBadBixelsPix::IsBad() set in:
145//
146// - MCalibrationRelTimePix
147//
148Bool_t MCalibrationRelTimeCalc::ReInit(MParList *pList )
149{
150
151 fGeom = (MGeomCam*)pList->FindObject("MGeomCam");
152 if (!fGeom)
153 {
154 *fLog << err << "No MGeomCam found... aborting." << endl;
155 return kFALSE;
156 }
157
158 fIntensBad = (MBadPixelsIntensityCam*)pList->FindObject(AddSerialNumber("MBadPixelsIntensityCam"));
159 if (fIntensBad)
160 *fLog << inf << "Found MBadPixelsIntensityCam ... " << endl;
161 else
162 {
163 fBadPixels = (MBadPixelsCam*)pList->FindObject(AddSerialNumber("MBadPixelsCam"));
164 if (!fBadPixels)
165 {
166 *fLog << err << "Cannot find MBadPixelsCam ... abort." << endl;
167 return kFALSE;
168 }
169 }
170
171 fIntensCam = (MCalibrationIntensityRelTimeCam*)pList->FindObject(AddSerialNumber("MCalibrationIntensityRelTimeCam"));
172 if (fIntensCam)
173 *fLog << inf << "Found MCalibrationIntensityRelTimeCam ... " << endl;
174 else
175 {
176 fCam = (MCalibrationRelTimeCam*)pList->FindObject(AddSerialNumber("MCalibrationRelTimeCam"));
177 if (!fCam)
178 {
179 *fLog << err << "Cannot find MCalibrationRelTimeCam ... abort." << endl;
180 *fLog << err << "Maybe you forget to call an MFillH for the MHCalibrationRelTimeCam before..." << endl;
181 return kFALSE;
182 }
183 }
184
185 UInt_t npixels = fGeom->GetNumPixels();
186
187 MCalibrationRelTimeCam *relcam = fIntensCam
188 ? (MCalibrationRelTimeCam*)fIntensCam->GetCam() : fCam;
189 MBadPixelsCam *badcam = fIntensBad
190 ? (MBadPixelsCam*) fIntensBad->GetCam() : fBadPixels;
191
192 for (UInt_t i=0; i<npixels; i++)
193 {
194
195 MCalibrationRelTimePix &pix = (MCalibrationRelTimePix&)(*relcam)[i];
196 MBadPixelsPix &bad = (*badcam)[i];
197
198 if (bad.IsBad())
199 {
200 pix.SetExcluded();
201 continue;
202 }
203
204 if (IsDebug())
205 pix.SetDebug();
206 }
207
208 return kTRUE;
209}
210
211// -----------------------------------------------------------------------
212//
213// Return if number of executions is null.
214//
215Int_t MCalibrationRelTimeCalc::PostProcess()
216{
217
218 if (GetNumExecutions()==0)
219 return kTRUE;
220
221 return Finalize();
222}
223
224// -----------------------------------------------------------------------
225//
226// First loop over pixels, average areas and sectors, call:
227// - FinalizeRelTimes()
228// for every entry. Count number of valid pixels in loop and return kFALSE
229// if there are none (the "Michele check").
230//
231// Call FinalizeBadPixels()
232//
233// Call MParContainer::SetReadyToSave() for fCam
234//
235// Print out some statistics
236//
237Int_t MCalibrationRelTimeCalc::Finalize()
238{
239
240 //
241 // First loop over pixels, call FinalizePedestals and FinalizeRelTimes
242 //
243 FinalizeRelTimes();
244
245 //
246 // Finalize Bad Pixels
247 //
248 FinalizeBadPixels();
249
250 //
251 // Re-direct the output to an ascii-file from now on:
252 //
253 MLog *asciilog = fOutputFile.IsNull() ? 0 : new MLog;
254 if (asciilog)
255 {
256 asciilog->SetOutputFile(GetOutputFile(),kTRUE);
257 SetLogStream(asciilog);
258 }
259
260 //
261 // Finalize calibration statistics
262 //
263 FinalizeUnsuitablePixels();
264
265 if (fIntensCam)
266 fIntensCam->SetReadyToSave();
267 else
268 fCam ->SetReadyToSave();
269
270 if (fIntensBad)
271 fIntensBad->SetReadyToSave();
272 else
273 fBadPixels->SetReadyToSave();
274
275 *fLog << inf << endl;
276 *fLog << GetDescriptor() << ": Errors statistics:" << endl;
277
278 PrintUncalibrated(MBadPixelsPix::kDeviatingTimeResolution,
279 Form("%s%2.1f%s","Time resol. less than ", fRelTimeResolutionLimit, " FADC sl. from Mean:"));
280 PrintUncalibrated(MBadPixelsPix::kRelTimeOscillating,
281 "Pixels with changing Rel. Times over time:");
282 PrintUncalibrated(MBadPixelsPix::kRelTimeNotFitted,
283 "Pixels with unsuccesful Gauss fit to the times:");
284
285 if (asciilog)
286 {
287 SetLogStream(&gLog);
288 delete asciilog;
289 }
290
291 return kTRUE;
292}
293
294
295// ----------------------------------------------------------------------------------------------------
296//
297//
298// First loop: Calculate a mean and mean RMS of time resolution per area index
299// Include only pixels which are not MBadPixelsPix::kUnsuitableRun or
300// MBadPixelsPix::kUnreliableRun (see FinalizeBadPixels())
301//
302// Second loop: Exclude those deviating by more than fRelTimeResolutionLimit FADC slices
303// from the mean (obtained in first loop). Set
304// MBadPixelsPix::kDeviatingTimeResolution if excluded.
305//
306void MCalibrationRelTimeCalc::FinalizeRelTimes()
307{
308
309 MCalibrationRelTimeCam *relcam = fIntensCam
310 ? (MCalibrationRelTimeCam*)fIntensCam->GetCam() : fCam;
311 MBadPixelsCam *badcam = fIntensBad
312 ? (MBadPixelsCam*) fIntensBad->GetCam() : fBadPixels;
313
314 const UInt_t npixels = fGeom->GetNumPixels();
315 const UInt_t nareas = fGeom->GetNumAreas();
316
317 TArrayF lowlim (nareas);
318 TArrayF upplim (nareas);
319 TArrayF areasum (nareas);
320 // Float_t areasum2 [nareas];
321 TArrayI numareavalid (nareas);
322 TArrayI useunreliable(nareas);
323
324 //
325 // Apero loop: Count number of unreliable pixels:
326 //
327 for (UInt_t i=0; i<npixels; i++)
328 {
329 MBadPixelsPix &bad = (*badcam)[i];
330 const Int_t aidx = (*fGeom)[i].GetAidx();
331
332 if (bad.IsUnsuitable(MBadPixelsPix::kUnsuitableRun))
333 continue;
334
335 if (bad.IsUnsuitable(MBadPixelsPix::kUnreliableRun))
336 continue;
337
338 numareavalid[aidx] ++;
339 }
340
341 for (UInt_t aidx=0; aidx<nareas; aidx++)
342 if (numareavalid[aidx] < 100)
343 useunreliable[aidx] = 1;
344
345 numareavalid.Reset();
346 //
347 // First loop: Get mean time resolution the RMS
348 // The loop is only to recognize later pixels with very deviating numbers
349 //
350 for (UInt_t i=0; i<npixels; i++)
351 {
352
353 MCalibrationRelTimePix &pix = (MCalibrationRelTimePix&)(*relcam)[i];
354 MBadPixelsPix &bad = (*badcam)[i];
355
356 if (pix.IsExcluded())
357 continue;
358
359 if (bad.IsUnsuitable(MBadPixelsPix::kUnsuitableRun))
360 continue;
361
362 const Int_t aidx = (*fGeom)[i].GetAidx();
363
364 if (!useunreliable[aidx])
365 if (bad.IsUnsuitable(MBadPixelsPix::kUnreliableRun))
366 continue;
367
368 const Float_t res = pix.GetTimePrecision();
369
370 areasum [aidx] += res;
371 // areasum2 [aidx] += res*res;
372 numareavalid[aidx] ++;
373 }
374
375
376 for (UInt_t aidx=0; aidx<nareas; aidx++)
377 {
378 if (numareavalid[aidx] < 20)
379 {
380 *fLog << warn << GetDescriptor() << ": Less than 20 pixels with valid time resolution found "
381 << "in area index: " << aidx << endl;
382 continue;
383 }
384
385 // Calculate the rms out of sum2:
386 /*
387 areasum2[aidx] = (areasum2[aidx] - areasum[aidx]*areasum[aidx]/numareavalid[aidx]);
388 areasum2[aidx] /= (numareavalid[aidx]-1.);
389 */
390 areasum [aidx] /= numareavalid[aidx];
391 lowlim [aidx] = 0.;
392 upplim [aidx] = areasum [aidx] + fRelTimeResolutionLimit;
393
394 }
395 *fLog << endl;
396
397
398 for (UInt_t i=0; i<npixels; i++)
399 {
400
401 MCalibrationRelTimePix &pix = (MCalibrationRelTimePix&)(*relcam)[i];
402 MBadPixelsPix &bad = (*badcam)[i];
403
404 if (pix.IsExcluded())
405 continue;
406
407 if (bad.IsUnsuitable(MBadPixelsPix::kUnsuitableRun))
408 continue;
409
410 const Float_t res = pix.GetTimePrecision();
411 const Int_t aidx = (*fGeom)[i].GetAidx();
412
413 if ( res < lowlim[aidx] || res > upplim[aidx] )
414 {
415 *fLog << warn << "Deviating time resolution: "
416 << Form("%4.2f",res) << " out of range ["
417 << Form("%4.2f,%4.2f",lowlim[aidx],upplim[aidx]) << "] in pixel " << i << endl;
418 bad.SetUncalibrated( MBadPixelsPix::kDeviatingTimeResolution);
419 pix.SetExcluded();
420 }
421 }
422}
423
424
425// -----------------------------------------------------------------------------------
426//
427// Sets pixel to MBadPixelsPix::kUnsuitableRun, if one of the following flags is set:
428// - MBadPixelsPix::kRelTimeIsPedestal
429// - MBadPixelsPix::kRelTimeErrNotValid
430// - MBadPixelsPix::kRelTimeRelErrNotValid
431//
432// - Call MCalibrationPix::SetExcluded() for the bad pixels
433//
434void MCalibrationRelTimeCalc::FinalizeBadPixels()
435{
436
437 MCalibrationRelTimeCam *relcam = fIntensCam
438 ? (MCalibrationRelTimeCam*)fIntensCam->GetCam() : fCam;
439 MBadPixelsCam *badcam = fIntensBad
440 ? (MBadPixelsCam*) fIntensBad->GetCam() : fBadPixels;
441
442 for (Int_t i=0; i<badcam->GetSize(); i++)
443 {
444
445 MBadPixelsPix &bad = (*badcam)[i];
446 MCalibrationRelTimePix &pix = (MCalibrationRelTimePix&)(*relcam)[i];
447
448 if (IsCheckDeviatingBehavior())
449 if (bad.IsUncalibrated(MBadPixelsPix::kDeviatingTimeResolution))
450 bad.SetUnsuitable(MBadPixelsPix::kUnreliableRun);
451
452 if (IsCheckFitResults())
453 if (bad.IsUncalibrated(MBadPixelsPix::kRelTimeNotFitted))
454 bad.SetUnsuitable(MBadPixelsPix::kUnreliableRun);
455
456 if (IsCheckOscillations())
457 if (bad.IsUncalibrated(MBadPixelsPix::kRelTimeOscillating))
458 bad.SetUnsuitable(MBadPixelsPix::kUnreliableRun);
459
460 if (bad.IsUnsuitable(MBadPixelsPix::kUnsuitableRun))
461 pix.SetExcluded();
462 }
463
464}
465
466
467// -----------------------------------------------------------------------------------------------
468//
469// - Print out statistics about BadPixels of type UnsuitableType_t
470// - store numbers of bad pixels of each type in fIntensCam or fCam
471//
472void MCalibrationRelTimeCalc::FinalizeUnsuitablePixels()
473{
474 *fLog << inf << endl;
475 *fLog << GetDescriptor() << ": Rel. Times Calibration status:" << endl;
476 *fLog << dec;
477
478 MCalibrationRelTimeCam *relcam = fIntensCam ? (MCalibrationRelTimeCam*)fIntensCam->GetCam() : fCam;
479 const MBadPixelsCam *badcam = fIntensBad ? (MBadPixelsCam*)fIntensBad->GetCam() : fBadPixels;
480
481 const Int_t nareas = fGeom->GetNumAreas();
482
483 TArrayI unsuit(nareas);
484 TArrayI unrel(nareas);
485
486 for (int aidx=0; aidx<nareas; aidx++)
487 {
488 unsuit[aidx] += badcam->GetNumUnsuitable(MBadPixelsPix::kUnsuitableRun, fGeom, aidx);
489 unrel[aidx] += badcam->GetNumUnsuitable(MBadPixelsPix::kUnreliableRun, fGeom, aidx);
490 relcam->SetNumUnsuitable(unsuit[aidx], aidx);
491 relcam->SetNumUnreliable(unrel[aidx], aidx);
492 }
493
494 if (fGeom->InheritsFrom("MGeomCamMagic"))
495 {
496 PrintUncalibrated("Uncalibrated Pixels:", unsuit[0], unsuit[1]);
497 PrintUncalibrated("Unreliable Pixels:", unrel[0], unrel[1]);
498 }
499}
500
501// -----------------------------------------------------------------------------------------------
502//
503// Print out statistics about BadPixels of type UncalibratedType_t
504//
505void MCalibrationRelTimeCalc::PrintUncalibrated(MBadPixelsPix::UncalibratedType_t typ, const char *text) const
506{
507 const MBadPixelsCam *badcam = fIntensBad ? (MBadPixelsCam*)fIntensBad->GetCam() : fBadPixels;
508
509 UInt_t countinner = 0;
510 UInt_t countouter = 0;
511 for (Int_t i=0; i<badcam->GetSize(); i++)
512 {
513 if ((*badcam)[i].IsUncalibrated(typ))
514 {
515 if (fGeom->GetPixRatio(i) == 1.)
516 countinner++;
517 else
518 countouter++;
519 }
520 }
521
522 PrintUncalibrated(text, countinner, countouter);
523}
524
525void MCalibrationRelTimeCalc::PrintUncalibrated(const char *text, Int_t in, Int_t out) const
526{
527 *fLog << " " << setfill(' ') << setw(56) << setiosflags(ios::left) << text;
528 *fLog << " Inner: " << Form("%3i", in);
529 *fLog << " Outer: " << Form("%3i", out);
530 *fLog << resetiosflags(ios::left) << endl;
531
532}
533
534// --------------------------------------------------------------------------
535//
536// Set the path for output file
537//
538void MCalibrationRelTimeCalc::SetOutputPath(TString path)
539{
540 fOutputPath = path;
541 if (fOutputPath.EndsWith("/"))
542 fOutputPath = fOutputPath(0, fOutputPath.Length()-1);
543}
544
545// --------------------------------------------------------------------------
546//
547// Get the output file
548//
549const char* MCalibrationRelTimeCalc::GetOutputFile()
550{
551 return Form("%s/%s", (const char*)fOutputPath, (const char*)fOutputFile);
552}
553
554
555// --------------------------------------------------------------------------
556//
557// MCalibrationRelTimeCam.CheckFitResults: Yes
558// MCalibrationRelTimeCam.CheckDeviatingBehavior: Yes
559// MCalibrationRelTimeCam.CheckHistOverflow: Yes
560// MCalibrationRelTimeCam.CheckOscillations: Yes
561//
562Int_t MCalibrationRelTimeCalc::ReadEnv(const TEnv &env, TString prefix, Bool_t print)
563{
564 Bool_t rc = kFALSE;
565
566 if (IsEnvDefined(env, prefix, "CheckFitResults", print))
567 {
568 SetCheckFitResults(GetEnvValue(env, prefix, "CheckFitResults", IsCheckFitResults()));
569 rc = kTRUE;
570 }
571 if (IsEnvDefined(env, prefix, "CheckDeviatingBehavior", print))
572 {
573 SetCheckDeviatingBehavior(GetEnvValue(env, prefix, "CheckDeviatingBehavior", IsCheckDeviatingBehavior()));
574 rc = kTRUE;
575 }
576 if (IsEnvDefined(env, prefix, "CheckHistOverflow", print))
577 {
578 SetCheckHistOverflow(GetEnvValue(env, prefix, "CheckHistOverflow", IsCheckHistOverflow()));
579 rc = kTRUE;
580 }
581
582 if (IsEnvDefined(env, prefix, "CheckOscillations", print))
583 {
584 SetCheckOscillations(GetEnvValue(env, prefix, "CheckOscillations", IsCheckOscillations()));
585 rc = kTRUE;
586 }
587 if (IsEnvDefined(env, prefix, "RelTimeResolutionLimit", print))
588 {
589 SetRelTimeResolutionLimit(GetEnvValue(env, prefix, "RelTimeResolutionLimit", fRelTimeResolutionLimit));
590 rc = kTRUE;
591 }
592
593 return rc;
594}
Note: See TracBrowser for help on using the repository browser.