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

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