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

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