source: trunk/MagicSoft/Mars/mbadpixels/MBadPixelsCam.cc@ 5563

Last change on this file since 5563 was 5517, checked in by gaug, 21 years ago
*** empty log message ***
File size: 19.4 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): Thomas Bretz 1/2004 <mailto:tbretz@astro.uni-wuerzburg.de>
19! Markus Gaug 3/2004 <mailto:markus@ifae.es>
20!
21! Copyright: MAGIC Software Development, 2000-2004
22!
23!
24\* ======================================================================== */
25
26/////////////////////////////////////////////////////////////////////////////
27//
28// MBadPixelsCam
29//
30//
31// Storage container to store bad pixel of the camera...
32//
33/////////////////////////////////////////////////////////////////////////////
34#include "MBadPixelsCam.h"
35
36#include <iostream>
37
38#include <TClonesArray.h>
39
40#include "MLog.h"
41#include "MLogManip.h"
42
43#include "MBadPixelsPix.h"
44
45#include "MGeomCam.h"
46#include "MGeomPix.h"
47
48ClassImp(MBadPixelsCam);
49
50using namespace std;
51
52// --------------------------------------------------------------------------
53//
54// Default constructor.
55//
56MBadPixelsCam::MBadPixelsCam(const char *name, const char *title)
57{
58 fName = name ? name : "MBadPixelsCam";
59 fTitle = title ? title : "Storage container to store bad pixel information";
60
61 fArray = new TClonesArray("MBadPixelsPix", 1);
62}
63
64MBadPixelsCam::MBadPixelsCam(const MBadPixelsCam &cam)
65{
66 fName = "MBadPixelsCam";
67 fTitle = "Storage container to store bad pixel information";
68
69 fArray = new TClonesArray("MBadPixelsPix", 1);
70 cam.Copy(*this);
71}
72
73// --------------------------------------------------------------------------
74//
75// Delete the array conatining the bad pixel information
76//
77MBadPixelsCam::~MBadPixelsCam()
78{
79 delete fArray;
80}
81
82// --------------------------------------------------------------------------
83//
84// Set the size of the camera
85//
86void MBadPixelsCam::InitSize(const UInt_t i)
87{
88 fArray->ExpandCreate(i);
89}
90
91// --------------------------------------------------------------------------
92//
93// Get the size of the MBadPixelsCam
94//
95Int_t MBadPixelsCam::GetSize() const
96{
97 return fArray->GetEntriesFast();
98}
99
100// --------------------------------------------------------------------------
101//
102// Copy 'constructor'
103//
104void MBadPixelsCam::Copy(TObject &object) const
105{
106 MBadPixelsCam &cam = (MBadPixelsCam&)object;
107
108 const Int_t n = GetSize();
109
110 if (n==0)
111 return;
112
113 cam.InitSize(n);
114 for (int i=0; i<n; i++)
115 (*this)[i].Copy(cam[i]);
116}
117
118// --------------------------------------------------------------------------
119//
120// Get i-th pixel (pixel number)
121//
122MBadPixelsPix &MBadPixelsCam::operator[](Int_t i)
123{
124 return *static_cast<MBadPixelsPix*>(fArray->UncheckedAt(i));
125}
126
127// --------------------------------------------------------------------------
128//
129// Get i-th pixel (pixel number)
130//
131const MBadPixelsPix &MBadPixelsCam::operator[](Int_t i) const
132{
133 return *static_cast<MBadPixelsPix*>(fArray->UncheckedAt(i));
134}
135
136// --------------------------------------------------------------------------
137//
138// Merge two MBadPixelsCam together, see also MBadPixelsPix::Merge
139//
140void MBadPixelsCam::Merge(const MBadPixelsCam &cam)
141{
142 const Int_t n = cam.GetSize();
143 if (n==0)
144 {
145 *fLog << warn << "MBadPixelsCam::Merge: Container empty." << endl;
146 return;
147 }
148
149 if (GetSize()==0)
150 InitSize(n);
151
152 if (n!=GetSize())
153 {
154 *fLog << warn << "MBadPixelsCam::Merge: Size mismatch... ignored." << endl;
155 return;
156 }
157
158 for (int i=0; i<n; i++)
159 (*this)[i].Merge(cam[i]);
160}
161
162// --------------------------------------------------------------------------
163//
164// Clear the contents of all bad pixels (set=0 means Ok)
165//
166void MBadPixelsCam::Clear(Option_t *o)
167{
168 fArray->ForEach(TObject, Clear)();
169}
170
171// --------------------------------------------------------------------------
172//
173// Reset event depending bits
174//
175void MBadPixelsCam::Reset()
176{
177 fArray->ForEach(MParContainer, Reset)();
178}
179
180// --------------------------------------------------------------------------
181//
182// Calculate the number of pixels with the given type-flags.
183//
184// The second argument aidx is the area index (see MGeomCam, MGeomPix)
185// The default (or any value less than 0) means: all
186//
187// Returns -1 if the geometry doesn't match.
188//
189Short_t MBadPixelsCam::GetNumUnsuitable(MBadPixelsPix::UnsuitableType_t type, const MGeomCam *geom, Int_t aidx) const
190{
191 const UInt_t n = GetSize();
192
193 if (aidx>=0 && geom->GetNumPixels()!=n)
194 {
195 *fLog << err << GetDescriptor() << "ERROR - Geometry (" << geom->ClassName() << ") size mismatch!" << endl;
196 return -1;
197 }
198
199 Short_t rc = 0;
200 for (UInt_t i=0; i<n; i++)
201 {
202 if (aidx>=0 && (*geom)[i].GetAidx()!=aidx)
203 continue;
204
205 if ((*this)[i].IsUnsuitable(type))
206 rc++;
207 }
208 return rc;
209}
210
211// --------------------------------------------------------------------------
212//
213// Counts the number of neighbors matching NOT UnsuitableType type
214//
215Short_t MBadPixelsCam::GetNumSuitableNeighbors(MBadPixelsPix::UnsuitableType_t type, const MGeomPix &pix) const
216{
217 const Int_t n2 = pix.GetNumNeighbors();
218
219 Int_t cnt=0;
220 for (int j=0; j<n2; j++)
221 {
222 const Int_t id2 = pix.GetNeighbor(j);
223 if (!(*this)[id2].IsUnsuitable(type))
224 cnt++;
225 }
226
227 return cnt;
228}
229
230// --------------------------------------------------------------------------
231//
232// Calculate the number of pixels which are - under no circumstances -
233// interpolatable, called isolated. This means that a pixel (its own status
234// doesn't matter) has less than two reliable neighbor pixels.
235//
236// The second argument aidx is the area index (see MGeomCam, MGeomPix)
237// The default (or any value less than 0) means: all
238//
239// Returns -1 if the geometry doesn't match.
240//
241Short_t MBadPixelsCam::GetNumIsolated(MBadPixelsPix::UnsuitableType_t type, const MGeomCam &geom, Int_t aidx) const
242{
243 const Int_t n = geom.GetNumPixels();
244
245 if (n!=GetSize())
246 {
247 *fLog << err << GetDescriptor() << "ERROR - Geometry (" << geom.ClassName() << ") size mismatch!" << endl;
248 return -1;
249 }
250
251 Short_t rc = 0;
252 for (int i=0; i<n; i++)
253 {
254 const MGeomPix &pix = geom[i];
255 if (aidx>=0 && pix.GetAidx()!=aidx)
256 continue;
257
258 if (GetNumSuitableNeighbors(type, pix)<2)
259 rc++;
260 }
261 return rc;
262}
263
264// --------------------------------------------------------------------------
265//
266// This is a helper function which calculates the size of a single cluster
267// by iterative calling.
268//
269// If a pixel matches the criterias the counter is increased by 1 and
270// the function is called for all its neighbors. If
271//
272// The second argument aidx is the area index (see MGeomCam, MGeomPix)
273// The default (or any value less than 0) means: all
274//
275// Returns -1 if the geometry doesn't match.
276//
277Short_t MBadPixelsCam::GetNumMaxCluster(MBadPixelsPix::UnsuitableType_t type, TObjArray &list, Int_t idx, Int_t aidx) const
278{
279 const MGeomPix *pix = (MGeomPix*)list[idx];
280 if (!pix)
281 return 0;
282
283 if (GetNumSuitableNeighbors(type, *pix)>1)
284 return 0;
285
286 if (aidx>=0 && pix->GetAidx()!=aidx)
287 return 1;
288
289 list.RemoveAt(idx);
290
291 Short_t cnt = 1;
292 const Int_t n = pix->GetNumNeighbors();
293 for (int i=0; i<n; i++)
294 cnt += GetNumMaxCluster(type, list, pix->GetNeighbor(i), aidx);
295
296 return cnt;
297}
298
299// --------------------------------------------------------------------------
300//
301// Returns the size of the biggest cluster with the given USuitableType
302// type and the given area index.
303//
304// The second argument aidx is the area index (see MGeomCam, MGeomPix)
305// The default (or any value less than 0) means: all
306//
307// Returns -1 if the geometry doesn't match.
308//
309Short_t MBadPixelsCam::GetNumMaxCluster(MBadPixelsPix::UnsuitableType_t type, const MGeomCam &geom, Int_t aidx) const
310{
311 const Int_t n = geom.GetNumPixels();
312
313 if (n!=GetSize())
314 {
315 *fLog << err << GetDescriptor() << "ERROR - Geometry (" << geom.ClassName() << ") size mismatch!" << endl;
316 return -1;
317 }
318
319 TObjArray list(n);
320 for (int i=0; i<n; i++)
321 list.AddAt(&geom[i], i);
322
323 Short_t max = 0;
324 for (int i=0; i<n; i++)
325 max = TMath::Max(GetNumMaxCluster(type, list, i, aidx), max);
326
327 return max;
328}
329
330// --------------------------------------------------------------------------
331//
332// Print the contents of all bad pixels
333//
334void MBadPixelsCam::Print(Option_t *o) const
335{
336
337 *fLog << all << GetDescriptor() << ":" << endl;
338 *fLog << "Pixels without problems:" << endl;
339
340 Int_t count = 0;
341
342 for (Int_t i=0; i<GetSize(); i++)
343 {
344 if (!(*this)[i].IsUnsuitable(MBadPixelsPix::kUnsuitableRun))
345 {
346 *fLog << i << " ";
347 count ++;
348 }
349
350 if (count == 0)
351 continue;
352
353 if (!(count % 25))
354 *fLog << endl;
355 }
356 *fLog << count << " normal pixels" << endl;
357 *fLog << endl;
358
359 *fLog << "Pixels unsuited for the whole run:" << endl;
360
361 for (Int_t i=0; i<GetSize(); i++)
362 {
363 if ((*this)[i].IsUnsuitable(MBadPixelsPix::kUnsuitableRun))
364 {
365 *fLog << i << " ";
366 count ++;
367 }
368
369 if (count == 0)
370 continue;
371
372 if (!(count % 25))
373 *fLog << endl;
374 }
375
376 *fLog << count << " unsuited pixels :-(" << endl;
377 *fLog << endl;
378
379 count = 0;
380
381 *fLog << all << "Pixels unreliable for the whole run:" << endl;
382
383 for (Int_t i=0; i<GetSize(); i++)
384 {
385 if ((*this)[i].IsUnsuitable(MBadPixelsPix::kUnreliableRun))
386 {
387 *fLog << i << " ";
388 count ++;
389 }
390
391 if (count == 0)
392 continue;
393
394 if (!(count % 25))
395 *fLog << endl;
396 }
397
398 *fLog << count << " unreliable pixels :-(" << endl;
399 *fLog << endl;
400 *fLog << endl;
401 *fLog << all << "Unsuited pixels statistics:" << endl;
402 *fLog << endl;
403
404 PrintBadPixels(MBadPixelsPix::kChargeIsPedestal,"Signal smaller 3 Pedestal RMS");
405 PrintBadPixels(MBadPixelsPix::kChargeSigmaNotValid,"Signal smaller 3 Pedestal RMS");
406 PrintBadPixels(MBadPixelsPix::kChargeRelErrNotValid,"Signal Rel. error too large");
407 PrintBadPixels(MBadPixelsPix::kLoGainSaturation,"Low Gain Saturation");
408 PrintBadPixels(MBadPixelsPix::kMeanTimeInFirstBin,"Mean Arr. Time In First Extraction Bin");
409 PrintBadPixels(MBadPixelsPix::kMeanTimeInLast2Bins,"Mean Arr. Time In Last 2 Extraction Bins");
410 PrintBadPixels(MBadPixelsPix::kDeviatingNumPhes,"Deviating Number of Photo-electrons");
411 PrintBadPixels(MBadPixelsPix::kDeviatingNumPhots,"Deviating Number of Photons");
412 PrintBadPixels(MBadPixelsPix::kHiGainOverFlow,"High-Gain Histogram Overflow");
413 PrintBadPixels(MBadPixelsPix::kLoGainOverFlow,"Low-Gain Histogram Overflow");
414
415 *fLog << endl;
416 *fLog << all << "Unreliable pixels statistics:" << endl;
417 *fLog << endl;
418
419 PrintBadPixels(MBadPixelsPix::kChargeSigmaNotValid,"Signal Sigma smaller Pedestal RMS");
420 PrintBadPixels(MBadPixelsPix::kHiGainNotFitted ,"High Gain Signals could not be fitted");
421 PrintBadPixels(MBadPixelsPix::kLoGainNotFitted ,"Low Gain Signals could not be fitted");
422 PrintBadPixels(MBadPixelsPix::kRelTimeNotFitted ,"Relative Arr. Times could not be fitted");
423 PrintBadPixels(MBadPixelsPix::kHiGainOscillating ,"High Gain Signals Oscillation");
424 PrintBadPixels(MBadPixelsPix::kLoGainOscillating ,"Low Gain Signals Oscillation");
425 PrintBadPixels(MBadPixelsPix::kRelTimeOscillating ,"Relative Arr. Times Oscillation");
426 PrintBadPixels(MBadPixelsPix::kDeviatingFFactor ,"Deviating global F-Factor");
427}
428
429void MBadPixelsCam::PrintBadPixels( MBadPixelsPix::UncalibratedType_t typ, const char *text) const
430{
431
432
433 *fLog << "Pixels with " << text << ": " << endl;
434 UInt_t count = 0;
435
436 for (Int_t i=0; i<GetSize(); i++)
437 {
438 if ((*this)[i].IsUncalibrated(typ))
439 {
440 *fLog << i << " ";
441 count++;
442 }
443
444 if (count == 0)
445 continue;
446
447 if (!(count % 25))
448 *fLog << endl;
449 }
450
451 *fLog << Form("%3i",count) << " pixels in total " << endl;
452 *fLog << endl;
453}
454
455// --------------------------------------------------------------------------
456//
457// Read from an ascii file of the format:
458// pixel1 pixel2 pixel3 pixel4
459// while pixel1,2,3,4 are the pixel indices used in the software.
460//
461// To read the pixels corresponding to a given run you can give run!=0
462// and a file of the format:
463// 1234: 17 193 292
464// 1235: 17 193 292 293
465//
466void MBadPixelsCam::AsciiRead(istream &fin, UInt_t run=0)
467{
468 Int_t len;
469 TString str;
470
471 while (1)
472 {
473 str.ReadLine(fin);
474 if (!fin)
475 return;
476
477 Int_t r;
478
479 const Int_t n = sscanf(str.Data(), " %d : %n", &r, &len);
480 if (n!=1)
481 return;
482
483 if (run==0 || run && (UInt_t)r==run)
484 break;
485 }
486
487 str.Remove(0, len);
488
489 while (1)
490 {
491 Int_t idx;
492 const Int_t n = sscanf(str.Data(), " %d %n", &idx, &len);
493
494 if (n!=1)
495 break;
496
497 str.Remove(0, len);
498
499 if (idx>=GetSize())
500 InitSize(idx+1);
501
502 (*this)[idx].SetUnsuitable(MBadPixelsPix::kUnsuitableRun);
503 }
504}
505
506// --------------------------------------------------------------------------
507//
508// Write the information into an ascii file. If a run-number is given the
509// run-number will lead the line.
510//
511Bool_t MBadPixelsCam::AsciiWrite(ostream &fout, UInt_t run=0) const
512{
513 if (run)
514 fout << run << ":";
515
516 for (int i=0; i<GetSize(); i++)
517 if ((*this)[i].IsUnsuitable(MBadPixelsPix::kUnsuitableRun))
518 fout << " " << i;
519
520 if (run && GetSize())
521 fout << endl;
522
523 return kTRUE;
524}
525
526// --------------------------------------------------------------------------
527//
528// The types are the following:
529// 0: MBadPixelsPix::GetInfo()[0]
530// 1: MBadPixelsPix::IsUnsuitable(MBadPixelsPix::kUnsuitableRun)
531// 2: MBadPixelsPix::IsUnsuitable(MBadPixelsPix::kUnsuitableEvt)
532// 3: MBadPixelsPix::IsUnsuitable(MBadPixelsPix::kUnreliableRun)
533// 4: MBadPixelsPix::IsHiGainBad()
534// 5: MBadPixelsPix::IsLoGainBad()
535// 6: MBadPixelsPix::GetUnsuitableCalibration()
536// 7: MBadPixelsPix::GetUnreliableCalibration()
537// 8: MBadPixelsPix::IsUncalibrated(MBadPixelsPix::kHiGainNotFitted)
538// 9: MBadPixelsPix::IsUncalibrated(MBadPixelsPix::kLoGainNotFitted)
539// 10: MBadPixelsPix::IsUncalibrated(MBadPixelsPix::kHiGainOscillating)
540// 11: MBadPixelsPix::IsUncalibrated(MBadPixelsPix::kLoGainOscillating)
541// 12: MBadPixelsPix::IsUncalibrated(MBadPixelsPix::kLoGainSaturation )
542// 13: MBadPixelsPix::IsUncalibrated(MBadPixelsPix::kChargeIsPedestal )
543// 14: MBadPixelsPix::IsUncalibrated(MBadPixelsPix::kChargeErrNotValid)
544// 15: MBadPixelsPix::IsUncalibrated(MBadPixelsPix::kChargeRelErrNotValid)
545// 16: MBadPixelsPix::IsUncalibrated(MBadPixelsPix::kChargeSigmaNotValid )
546// 17: MBadPixelsPix::IsUncalibrated(MBadPixelsPix::kMeanTimeInFirstBin )
547// 18: MBadPixelsPix::IsUncalibrated(MBadPixelsPix::kMeanTimeInLast2Bins )
548// 19: MBadPixelsPix::IsUncalibrated(MBadPixelsPix::kDeviatingNumPhes )
549// 20: MBadPixelsPix::IsUncalibrated(MBadPixelsPix::kRelTimeNotFitted )
550// 21: MBadPixelsPix::IsUncalibrated(MBadPixelsPix::kRelTimeOscillating )
551// 22: MBadPixelsPix::IsUncalibrated(MBadPixelsPix::kDeviatingNumPhots )
552// 23: MBadPixelsPix::IsUncalibrated(MBadPixelsPix::kHiGainOverFlow )
553// 24: MBadPixelsPix::IsUncalibrated(MBadPixelsPix::kLoGainOverFlow )
554//
555Bool_t MBadPixelsCam::GetPixelContent(Double_t &val, Int_t idx, const MGeomCam &cam, Int_t type) const
556{
557
558 if (idx >= GetSize())
559 return kFALSE;
560
561 switch (type)
562 {
563 case 0:
564 return (*this)[idx].GetInfo()[0];
565 case 1:
566 if (!(*this)[idx].IsUnsuitable(MBadPixelsPix::kUnsuitableRun))
567 return kFALSE;
568 val = 1;
569 break;
570 case 2:
571 if (!(*this)[idx].IsUnsuitable(MBadPixelsPix::kUnsuitableEvt))
572 return kFALSE;
573 val = 1;
574 break;
575 case 3:
576 if (!(*this)[idx].IsUnsuitable(MBadPixelsPix::kUnreliableRun))
577 return kFALSE;
578 val = 1;
579 break;
580 case 4:
581 if (!(*this)[idx].IsHiGainBad())
582 return kFALSE;
583 val = 1;
584 break;
585 case 5:
586 if (!(*this)[idx].IsLoGainBad())
587 return kFALSE;
588 val = 1;
589 break;
590 case 6:
591 if ((*this)[idx].GetUnsuitableCalibration())
592 val = (*this)[idx].GetUnsuitableCalibration();
593 else
594 return kFALSE;
595 break;
596 case 7:
597 if ((*this)[idx].GetUnreliableCalibration())
598 val = (*this)[idx].GetUnreliableCalibration();
599 else
600 return kFALSE;
601 break;
602 case 8:
603 if (!(*this)[idx].IsUncalibrated(MBadPixelsPix::kHiGainNotFitted))
604 return kFALSE;
605 val = 1;
606 break;
607 case 9:
608 if (!(*this)[idx].IsUncalibrated(MBadPixelsPix::kLoGainNotFitted))
609 return kFALSE;
610 val = 1;
611 break;
612 case 10:
613 if (!(*this)[idx].IsUncalibrated(MBadPixelsPix::kHiGainOscillating))
614 return kFALSE;
615 val = 1;
616 break;
617 case 11:
618 if (!(*this)[idx].IsUncalibrated(MBadPixelsPix::kLoGainOscillating))
619 return kFALSE;
620 val = 1;
621 break;
622 case 12:
623 if (!(*this)[idx].IsUncalibrated(MBadPixelsPix::kLoGainSaturation ))
624 return kFALSE;
625 val = 1;
626 break;
627 case 13:
628 if (!(*this)[idx].IsUncalibrated(MBadPixelsPix::kChargeIsPedestal ))
629 return kFALSE;
630 val = 1;
631 break;
632 case 14:
633 if (!(*this)[idx].IsUncalibrated(MBadPixelsPix::kChargeErrNotValid))
634 return kFALSE;
635 val = 1;
636 break;
637 case 15:
638 if (!(*this)[idx].IsUncalibrated(MBadPixelsPix::kChargeRelErrNotValid))
639 return kFALSE;
640 val = 1;
641 break;
642 case 16:
643 if (!(*this)[idx].IsUncalibrated(MBadPixelsPix::kChargeSigmaNotValid ))
644 return kFALSE;
645 val = 1;
646 break;
647 case 17:
648 if (!(*this)[idx].IsUncalibrated(MBadPixelsPix::kMeanTimeInFirstBin ))
649 return kFALSE;
650 val = 1;
651 break;
652 case 18:
653 if (!(*this)[idx].IsUncalibrated(MBadPixelsPix::kMeanTimeInLast2Bins ))
654 return kFALSE;
655 val = 1;
656 break;
657 case 19:
658 if (!(*this)[idx].IsUncalibrated(MBadPixelsPix::kDeviatingNumPhes ))
659 return kFALSE;
660 val = 1;
661 break;
662 case 20:
663 if (!(*this)[idx].IsUncalibrated(MBadPixelsPix::kRelTimeNotFitted ))
664 return kFALSE;
665 val = 1;
666 break;
667 case 21:
668 if (!(*this)[idx].IsUncalibrated(MBadPixelsPix::kRelTimeOscillating))
669 return kFALSE;
670 val = 1;
671 break;
672 case 22:
673 if (!(*this)[idx].IsUncalibrated(MBadPixelsPix::kDeviatingNumPhots))
674 return kFALSE;
675 val = 1;
676 break;
677 case 23:
678 if (!(*this)[idx].IsUncalibrated(MBadPixelsPix::kHiGainOverFlow ))
679 return kFALSE;
680 val = 1;
681 break;
682 case 24:
683 if (!(*this)[idx].IsUncalibrated(MBadPixelsPix::kLoGainOverFlow ))
684 return kFALSE;
685 val = 1;
686 break;
687 default:
688 return kFALSE;
689 }
690
691 return kTRUE;
692}
693
694void MBadPixelsCam::DrawPixelContent(Int_t num) const
695{
696 *fLog << warn << "MBadPixelsCam::DrawPixelContent - not available." << endl;
697}
Note: See TracBrowser for help on using the repository browser.