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

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