Changeset 3636 for trunk/MagicSoft/Mars

Timestamp:

04/03/04 17:27:50 (21 years ago)

Author:

gaug

Message:

*** empty log message ***

Location:

trunk/MagicSoft/Mars

Files:

• Changelog (modified) (1 diff)
• mcalib/MHCalibrationCam.cc (modified) (6 diffs)
• mcalib/MHCalibrationCam.h (modified) (5 diffs)
• mcalib/MHCalibrationChargeCam.cc (modified) (22 diffs)
• mcalib/MHCalibrationChargeCam.h (modified) (4 diffs)
• mcalib/MHCalibrationChargeHiGainPix.cc (modified) (3 diffs)
• mcalib/MHCalibrationChargeHiGainPix.h (modified) (2 diffs)
• mcalib/MHCalibrationChargeLoGainPix.cc (modified) (3 diffs)
• mcalib/MHCalibrationChargeLoGainPix.h (modified) (2 diffs)
• mcalib/MHCalibrationChargePix.cc (modified) (10 diffs)
• mcalib/MHCalibrationChargePix.h (modified) (3 diffs)
• mcalib/MHCalibrationRelTimeCam.cc (modified) (7 diffs)
• mcalib/MHCalibrationRelTimeCam.h (modified) (1 diff)
• mcalib/MHCalibrationRelTimePix.cc (modified) (6 diffs)
• mcalib/MHCalibrationRelTimePix.h (modified) (2 diffs)
• mcalib/MHGausEvents.cc (modified) (22 diffs)
• mcalib/MHGausEvents.h (modified) (9 diffs)

trunk/MagicSoft/Mars/Changelog

@@ -18 +18 @@
 
                                                  -*-*- END OF LINE -*-*-
+
+ 2004/04/02: Markus Gaug
+   * mcalib/MHGausEvents.[h,cc]
+   * mcalib/MHCalibrationChargePix.[h,cc]
+   * mcalib/MHCalibrationChargeHiGainPix.[h,cc]
+   * mcalib/MHCalibrationChargeLoGainPix.[h,cc]
+   * mcalib/MHCalibrationRelTimePix.[h,cc]
+     - put fPixId, fPickup, fPickupLimt, CountPickup(), RepeatFit() and 
+       ChangeHistId() into MHGausEvents (before in the derived classes)
+     - put fChargeNbins, fChargeFirst, fChargeLast, 
+     - put fRelTimeNbins, fRelTimeFirst, fRelTimeLast together 
+       into MHGausEvents as fNbins, fFirst and fLast
+
+   * mcalib/MHCalibrationRelTimePix.[h,cc]
+     - remove Renormalization to time slices. Need to think about 
+       more direct way to implement
+
+   * mcalib/MHCalibrationCam.[h,cc]
+   * mcalib/MHCalibrationChargeCam.[h,cc]
+   * mcalib/MHCalibrationRelTimeCam.[h,cc]
+     - put most of the functionality into the base class MHCalibrationCam
+     - derived classes overload the functions SetupHists, ReInitHists, 
+       FillHists, FinalizeHists and FinalizeBadPixels. 
+     - functions FitHiGainArrays, FitLoGainArrays, FitHiGainHists, 
+       FitLoGainHists and InitHists can be used from base class. 
 
 

trunk/MagicSoft/Mars/mcalib/MHCalibrationCam.cc

@@ -53 +53 @@
 #include "MLogManip.h"
 
+#include "MCalibrationPix.h"
+#include "MCalibrationCam.h"
+
 #include "MHGausEvents.h"
 
+#include "MBadPixelsPix.h"
+#include "MBadPixelsCam.h"
+
 #include "MGeomCam.h"
+#include "MGeomPix.h"
+
+#include "MParList.h"
 
 ClassImp(MHCalibrationCam);
@@ -61 +70 @@
 using namespace std;
 
+const Int_t   MHCalibrationCam::fgAverageNbins    = 2000;
 const Int_t   MHCalibrationCam::fgPulserFrequency = 500;
 // --------------------------------------------------------------------------
 //
 // Default Constructor. 
+//
+// Sets:
+// - all pointers to NULL
 //
 // Initializes and sets owner of:
@@ -75 +88 @@
 //
 MHCalibrationCam::MHCalibrationCam(const char *name, const char *title)
+    :  fGeom(NULL), fBadPixels(NULL), fCam(NULL)
 {
     fName  = name  ? name  : "MHCalibrationCam";
@@ -97 +111 @@
     fAverageLoGainSectors->SetOwner();
 
+    SetAverageNbins();
     SetPulserFrequency();
 }
@@ -306 +321 @@
 // --------------------------------------------------------------------------
 //
+// Gets the pointers to:
+// - MGeomCam
+//
+// Calls SetupHists(const MParList *pList)
+//
+// Calls Delete-Function of:
+// - MHCalibrationCam::fHiGainArray, MHCalibrationCam::fLoGainArray
+// - MHCalibrationCam::fAverageHiGainAreas, MHCalibrationCam::fAverageLoGainAreas
+// - MHCalibrationCam::fAverageHiGainSectors, MHCalibrationCam::fAverageLoGainSectors
+//
+Bool_t MHCalibrationCam::SetupFill(const MParList *pList)
+{
+  
+  fGeom = (MGeomCam*)pList->FindObject("MGeomCam");
+  if (!fGeom)
+  {
+      *fLog << err << "MGeomCam not found... aborting." << endl;
+      return kFALSE;
+  }
+
+  fHiGainArray->Delete();
+  fLoGainArray->Delete();
+
+  fAverageHiGainAreas->Delete();
+  fAverageLoGainAreas->Delete();
+
+  fAverageHiGainSectors->Delete();
+  fAverageLoGainSectors->Delete();
+
+  return SetupHists(pList);
+}
+
+
+Bool_t MHCalibrationCam::SetupHists(const MParList *pList)
+{
+  return kTRUE;
+}
+
+// --------------------------------------------------------------------------
+//
+// Gets or creates the pointers to:
+// - MBadPixelsCam
+//
+// Searches pointer to:
+// - MArrivalTimeCam
+//
+// Initializes, if empty to MArrivalTimeCam::GetSize() for:
+// - MHCalibrationCam::fHiGainArray, MHCalibrationCam::fLoGainArray
+//
+// Initializes, if empty to MGeomCam::GetNumAreas() for:
+// - MHCalibrationCam::fAverageHiGainAreas, MHCalibrationCam::fAverageLoGainAreas
+//
+// Initializes, if empty to MGeomCam::GetNumSectors() for:
+// - MHCalibrationCam::fAverageHiGainSectors, MHCalibrationCam::fAverageLoGainSectors
+// 
+// Initializes TArray's to MGeomCam::GetNumAreas and MGeomCam::GetNumSectors, respectively
+// Fills with number of valid pixels (if !MBadPixelsPix::IsBad()):
+// - MHCalibrationCam::fAverageAreaNum[area index]
+// - MHCalibrationCam::fAverageSectorNum[area index]
+//
+// Calls InitializeHists() for every entry in:
+// - MHCalibrationCam::fHiGainArray
+// - MHCalibrationCam::fAverageHiGainAreas
+// - MHCalibrationCam::fAverageHiGainSectors
+//
+// Sets Titles and Names for the Histograms 
+// - MHCalibrationCam::fAverageHiGainAreas
+// - MHCalibrationCam::fAverageHiGainSectors
+// 
+Bool_t MHCalibrationCam::ReInit(MParList *pList)
+{
+
+  fBadPixels = (MBadPixelsCam*)pList->FindCreateObj("MBadPixelsCam");
+  if (!fBadPixels)
+      return kFALSE;
+
+
+  const Int_t npixels  = fGeom->GetNumPixels();
+  const Int_t nsectors = fGeom->GetNumSectors();
+  const Int_t nareas   = fGeom->GetNumAreas();
+
+  //
+  // The function TArrayF::Set() already sets all entries to 0.
+  //
+  fAverageAreaNum.        Set(nareas);
+  fAverageAreaSat.        Set(nareas);           
+  fAverageAreaSigma.      Set(nareas);      
+  fAverageAreaSigmaErr.   Set(nareas);   
+  fAverageAreaRelSigma.   Set(nareas);   
+  fAverageAreaRelSigmaErr.Set(nareas);
+  fAverageSectorNum.      Set(nsectors);
+
+  for (Int_t i=0; i<npixels; i++)
+    {
+
+      if ((*fBadPixels)[i].IsBad())
+        continue;
+
+      fAverageAreaNum  [(*fGeom)[i].GetAidx()  ]++;
+      fAverageSectorNum[(*fGeom)[i].GetSector()]++;
+    }
+
+  return ReInitHists(pList);
+}
+
+
+Bool_t MHCalibrationCam::ReInitHists(MParList *pList)
+{
+  return kTRUE;
+}
+
+
+
+//--------------------------------------------------------------------------------
+//
+// Retrieves from MGeomCam:
+// - number of pixels
+// - number of pixel areas
+// - number of sectors
+//
+// For all TObjArray's (including the averaged ones), the following steps are performed: 
+//
+// 1) Test size and return kFALSE if not matching
+// 2) 
+//
+Bool_t MHCalibrationCam::Fill(const MParContainer *par, const Stat_t w)
+{
+
+  const Int_t npixels  = fGeom->GetNumPixels();
+  const Int_t nareas   = fGeom->GetNumAreas();
+  const Int_t nsectors = fGeom->GetNumSectors();
+  
+  if (fHiGainArray->GetEntries() != npixels)
+    {
+      gLog << err << "ERROR - Size mismatch... abort." << endl;
+      return kFALSE;
+    }
+  
+  if (fLoGainArray->GetEntries() != npixels)
+    {
+      gLog << err << "ERROR - Size mismatch... abort." << endl;
+      return kFALSE;
+    }
+  
+  if (fAverageHiGainAreas->GetEntries() != nareas)
+    {
+      *fLog << err << "ERROR - Size mismatch in number of areas ... abort." << endl;
+      return kFALSE;
+    }
+
+  if (fAverageLoGainAreas->GetEntries() != nareas)
+    {
+      *fLog << err << "ERROR - Size mismatch in number of areas ... abort." << endl;
+      return kFALSE;
+    }
+
+  if (fAverageHiGainSectors->GetEntries() != nsectors)
+    {
+      *fLog << err << "ERROR - Size mismatch in number of sectors ... abort." << endl;
+      return kFALSE;
+    }
+
+  if (fAverageLoGainSectors->GetEntries() != nsectors)
+    {
+      *fLog << err << "ERROR - Size mismatch in number of sectors ... abort." << endl;
+      return kFALSE;
+    }
+
+  return FillHists(par, w);
+}
+
+Bool_t MHCalibrationCam::FillHists(const MParContainer *par, const Stat_t w)
+{
+  return kTRUE;
+}
+
+// --------------------------------------------------------------------------
+//
+// 1) FinalizeHists()
+// 2) FinalizeBadPixels()
+// 3) CalcAverageSigma()
+//
+Bool_t MHCalibrationCam::Finalize()
+{
+  if (!FinalizeHists())
+    return kFALSE;
+
+  FinalizeBadPixels();
+  CalcAverageSigma();
+
+  return kTRUE;
+}
+
+Bool_t MHCalibrationCam::FinalizeHists()
+{
+  return kTRUE;
+}
+
+void MHCalibrationCam::FinalizeBadPixels()
+{
+}
+
+
+// -------------------------------------------------------------
+//
+// If MBadPixelsPix::IsBad():
+// - calls MHGausEvents::SetExcluded()
+//
+// Calls:
+// - MHGausEvents::InitBins()
+// - MHGausEvents::ChangeHistId(i)
+// - MHGausEvents::SetEventFrequency(fPulserFrequency)
+// 
+void MHCalibrationCam::InitHists(MHGausEvents &hist, MBadPixelsPix &bad, const Int_t i)
+{
+
+  if (bad.IsBad())
+    hist.SetExcluded();
+
+  hist.InitBins();
+  hist.ChangeHistId(i);
+  hist.SetEventFrequency(fPulserFrequency);
+          
+}
+
+void MHCalibrationCam::FitHiGainArrays(MCalibrationCam &calcam, MBadPixelsCam &badcam,
+                                       MBadPixelsPix::UncalibratedType_t fittyp, 
+                                       MBadPixelsPix::UncalibratedType_t osctyp)
+{
+  
+  for (Int_t i=0; i<fHiGainArray->GetSize(); i++)
+    {
+      
+      MHGausEvents &hist = (*this)[i];
+      
+      if (hist.IsExcluded())
+        continue;
+      
+      MCalibrationPix &pix    = calcam[i];
+      MBadPixelsPix   &bad    = badcam[i];
+      
+      FitHiGainHists(hist,pix,bad,fittyp,osctyp);
+      
+    }
+
+  for (Int_t j=0; j<fAverageHiGainAreas->GetSize(); j++)
+    {
+      
+      MHGausEvents     &hist = GetAverageHiGainArea(j);      
+      MCalibrationPix  &pix  = calcam.GetAverageArea(j);
+      MBadPixelsPix    &bad  = calcam.GetAverageBadArea(j);        
+      
+      FitHiGainHists(hist,pix,bad,fittyp,osctyp);
+  }
+  
+
+  for (Int_t j=0; j<fAverageHiGainSectors->GetSize(); j++)
+    {
+      
+      MHGausEvents     &hist = GetAverageHiGainSector(j);      
+      MCalibrationPix  &pix  = calcam.GetAverageSector(j);
+      MBadPixelsPix    &bad  = calcam.GetAverageBadSector(j);        
+      
+      FitHiGainHists(hist,pix,bad,fittyp,osctyp);
+    }
+
+}
+
+void MHCalibrationCam::FitLoGainArrays(MCalibrationCam &calcam, MBadPixelsCam &badcam,
+                                            MBadPixelsPix::UncalibratedType_t fittyp, 
+                                            MBadPixelsPix::UncalibratedType_t osctyp)
+{
+  
+  for (Int_t i=0; i<fLoGainArray->GetSize(); i++)
+    {
+      
+      MHGausEvents &hist = (*this)(i);
+      
+      if (hist.IsExcluded())
+        continue;
+      
+      MCalibrationPix &pix    = calcam[i];
+      MBadPixelsPix   &bad    = badcam[i];
+      
+      FitLoGainHists(hist,pix,bad,fittyp,osctyp);
+      
+    }
+
+  for (Int_t j=0; j<fAverageLoGainAreas->GetSize(); j++)
+    {
+      
+      MHGausEvents     &hist = GetAverageLoGainArea(j);      
+      MCalibrationPix  &pix  = calcam.GetAverageArea(j);
+      MBadPixelsPix    &bad  = calcam.GetAverageBadArea(j);        
+      
+      FitLoGainHists(hist,pix,bad,fittyp,osctyp);
+  }
+  
+
+  for (Int_t j=0; j<fAverageLoGainSectors->GetSize(); j++)
+    {
+      
+      MHGausEvents     &hist = GetAverageLoGainSector(j);      
+      MCalibrationPix  &pix  = calcam.GetAverageSector(j);
+      MBadPixelsPix    &bad  = calcam.GetAverageBadSector(j);        
+      
+      FitLoGainHists(hist,pix,bad,fittyp,osctyp);
+    }
+}
+
+//------------------------------------------------------------
+//
+// For all averaged areas, the fitted sigma is multiplied with the square root of 
+// the number involved pixels
+//
+void MHCalibrationCam::CalcAverageSigma()
+{
+  
+  for (Int_t j=0; j<fGeom->GetNumAreas(); j++)
+    {
+  
+      MCalibrationPix &pix    = (*fCam).GetAverageArea(j);
+
+      fAverageAreaSigma[j]    = pix.GetSigma    () * TMath::Sqrt((Float_t)fAverageAreaNum[j]);
+      fAverageAreaSigmaErr[j] = pix.GetSigmaErr () * TMath::Sqrt((Float_t)fAverageAreaNum[j]);
+
+      pix.SetSigma   (fAverageAreaSigma[j]);
+      pix.SetSigmaErr(fAverageAreaSigmaErr[j]);
+
+      fAverageAreaRelSigma[j]   = fAverageAreaSigma[j] / pix.GetMean();
+      
+      Float_t relsigmaerr       =  fAverageAreaSigmaErr[j]*fAverageAreaSigmaErr[j] 
+                                / (fAverageAreaSigma[j]   *fAverageAreaSigma[j]   );
+      relsigmaerr               += pix.GetMeanErr()*pix.GetMeanErr() 
+                                / (pix.GetMean()   *pix.GetMean()   );
+      relsigmaerr               *= fAverageAreaRelSigma[j];
+      fAverageAreaRelSigmaErr[j] = TMath::Sqrt(relsigmaerr);
+    }
+}
+
+// ---------------------------------------------------------------------------
+//
+// Returns if the histogram is empty and sets the following flag:
+// - MBadPixelsPix::SetUnsuitable(MBadPixelsPix::kUnsuitableRun)
+//
+// Fits the histograms with a Gaussian, in case of failure 
+// calls MHGausEvents::RepeatFit(), in case of repeated failure 
+// calls MHGausEvents::BypassFit() and sets the following flags:
+// - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::UncalibratedType_t fittyp )
+// - MBadPixelsPix::SetUnsuitable(   MBadPixelsPix::kUnreliableRun   )
+// 
+// Creates the fourier spectrum and tests MHGausEvents::IsFourierSpectrumOK(). 
+// In case no, sets the following flags:
+// - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::UncalibratedType_t osctyp )
+// - MBadPixelsPix::SetUnsuitable(   MBadPixelsPix::kUnreliableRun     )
+//
+// Retrieves the results and store them in MCalibrationPix
+//
+void MHCalibrationCam::FitHiGainHists(MHGausEvents &hist, 
+                                      MCalibrationPix &pix, 
+                                      MBadPixelsPix &bad, 
+                                      MBadPixelsPix::UncalibratedType_t fittyp,
+                                      MBadPixelsPix::UncalibratedType_t osctyp)
+{
+
+  if (hist.IsEmpty())
+    return;
+  
+  //
+  // 2) Fit the Hi Gain histograms with a Gaussian
+  //
+  if (!hist.FitGaus())
+    //
+    // 3) In case of failure set the bit Fitted to false and take histogram means and RMS
+    //
+    if (!hist.RepeatFit())
+      {
+        hist.BypassFit();
+        bad.SetUncalibrated( fittyp );
+      }
+  
+  //
+  // 4) Check for oscillations
+  // 
+  hist.CreateFourierSpectrum();
+  
+  if (!hist.IsFourierSpectrumOK())
+    bad.SetUncalibrated( osctyp );
+  
+  //
+  // 5) Retrieve the results and store them in this class
+  //
+  pix.SetHiGainMean      ( hist.GetMean()      );
+  pix.SetHiGainMeanErr   ( hist.GetMeanErr()   );
+  pix.SetHiGainSigma     ( hist.GetSigma()     );
+  pix.SetHiGainSigmaErr  ( hist.GetSigmaErr()  );
+  pix.SetHiGainProb      ( hist.GetProb()      );
+  pix.SetHiGainNumPickup ( hist.GetPickup()    );
+  
+}
+
+
+// ---------------------------------------------------------------------------
+//
+// Returns if the histogram is empty and sets the following flag:
+// - MBadPixelsPix::SetUnsuitable(MBadPixelsPix::kUnsuitableRun)
+//
+// Fits the histograms with a Gaussian, in case of failure 
+// calls MHGausEvents::RepeatFit(), in case of repeated failure 
+// calls MHGausEvents::BypassFit() and sets the following flags:
+// - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::UncalibratedType_t fittyp )
+// - MBadPixelsPix::SetUnsuitable(   MBadPixelsPix::kUnreliableRun   )
+// 
+// Creates the fourier spectrum and tests MHGausEvents::IsFourierSpectrumOK(). 
+// In case no, sets the following flags:
+// - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::UncalibratedType_t osctyp )
+// - MBadPixelsPix::SetUnsuitable(   MBadPixelsPix::kUnreliableRun     )
+//
+// Retrieves the results and store them in MCalibrationPix
+//
+void MHCalibrationCam::FitLoGainHists(MHGausEvents &hist, 
+                                      MCalibrationPix &pix, 
+                                      MBadPixelsPix &bad, 
+                                      MBadPixelsPix::UncalibratedType_t fittyp,
+                                      MBadPixelsPix::UncalibratedType_t osctyp)
+{
+
+  if (hist.IsEmpty())
+      return;
+  
+
+  //
+  // 2) Fit the Hi Gain histograms with a Gaussian
+  //
+  if (!hist.FitGaus())
+    //
+    // 3) In case of failure set the bit Fitted to false and take histogram means and RMS
+    //
+    if (!hist.RepeatFit())
+      {
+        hist.BypassFit();
+        bad.SetUncalibrated( fittyp );
+      }
+  
+  //
+  // 4) Check for oscillations
+  // 
+  hist.CreateFourierSpectrum();
+  
+  if (!hist.IsFourierSpectrumOK())
+    bad.SetUncalibrated( osctyp );
+  
+  //
+  // 5) Retrieve the results and store them in this class
+  //
+  pix.SetLoGainMean      ( hist.GetMean()      );
+  pix.SetLoGainMeanErr   ( hist.GetMeanErr()   );
+  pix.SetLoGainSigma     ( hist.GetSigma()     );
+  pix.SetLoGainSigmaErr  ( hist.GetSigmaErr()  );
+  pix.SetLoGainProb      ( hist.GetProb()      );
+  pix.SetLoGainNumPickup ( hist.GetPickup()    );
+  
+}
+
+
+
+// --------------------------------------------------------------------------
+//
 // Dummy, needed by MCamEvent
 //
@@ -343 +826 @@
   for (Int_t i=0; i<nareas;i++) 
     {
+
       pad->cd(2*(i+1)-1);
       GetAverageHiGainArea(i).Draw(opt);

trunk/MagicSoft/Mars/mcalib/MHCalibrationCam.h

@@ -21 +21 @@
 #endif
 
+#ifndef MARS_MBadPixelsPix
+#include "MBadPixelsPix.h"
+#endif
+
 class TText;
 class TArrayI;
 class TArrayF;
 class MHGausEvents;
+class MGeomCam;
+class MCalibrationCam;
+class MCalibrationPix;
+class MBadPixelsCam;
+class MBadPixelsPix;
 class MHCalibrationCam : public MH, public MCamEvent
 {
   
+private:
+
+  static const Int_t fgAverageNbins;     //! The default for fAverageNbins    (now set to: 2000)
+  static const Int_t fgPulserFrequency;  //! The default for fPulserFrequency (now set to: 500)
+
 protected:
 
-  static const Int_t fgPulserFrequency;  // The default for fPulserFrequency (now set to: 500)
+  Int_t   fAverageNbins;                // Number of bins for the average histograms
   Int_t   fPulserFrequency;             // Light pulser frequency
   
@@ -40 +54 @@
   TObjArray *fAverageLoGainSectors;     //-> Array of calibration pixels, one per camera sector
 
+  MGeomCam         *fGeom;              //!  Camera geometry
+  MBadPixelsCam    *fBadPixels;         //!  Bad Pixels storage container
+  MCalibrationCam  *fCam;               //!  Calibration Cam with the results
+
   TArrayI fAverageAreaNum;              // Number of pixels in average pixels per area
   TArrayI fAverageAreaSat;              // Number of saturated slices in average pixels per area
@@ -48 +66 @@
   TArrayI fAverageSectorNum;            // Number of pixels in average pixels per sector 
 
+  virtual Bool_t SetupHists(const MParList *pList);
+  virtual Bool_t ReInitHists(MParList *pList);  
+  virtual Bool_t FillHists(const MParContainer *par, const Stat_t w=1);
+  virtual Bool_t FinalizeHists();
+  virtual void   FinalizeBadPixels();
+  
+  void InitHists(MHGausEvents &hist, MBadPixelsPix &bad, const Int_t i);
+
+  void FitHiGainArrays(MCalibrationCam &calcam, MBadPixelsCam &badcam,
+                       MBadPixelsPix::UncalibratedType_t fittyp,
+                       MBadPixelsPix::UncalibratedType_t osctyp);
+  
+  void FitLoGainArrays(MCalibrationCam &calcam, MBadPixelsCam &badcam,
+                       MBadPixelsPix::UncalibratedType_t fittyp,
+                       MBadPixelsPix::UncalibratedType_t osctyp);
+  
+  void FitHiGainHists(MHGausEvents &hist, 
+                      MCalibrationPix &pix, 
+                      MBadPixelsPix &bad, 
+                      MBadPixelsPix::UncalibratedType_t fittyp,
+                      MBadPixelsPix::UncalibratedType_t osctyp);
+  
+  void FitLoGainHists(MHGausEvents &hist, 
+                      MCalibrationPix &pix, 
+                      MBadPixelsPix &bad, 
+                      MBadPixelsPix::UncalibratedType_t fittyp,
+                      MBadPixelsPix::UncalibratedType_t osctyp);
+
+  void CalcAverageSigma();
+  
   void DrawAverageSigma(Bool_t sat, Bool_t inner,
                         Float_t sigma, Float_t sigmaerr,
@@ -57 +105 @@
   ~MHCalibrationCam();
 
+  void SetAverageNbins(   const Int_t bins=fgAverageNbins ) { fAverageNbins = bins; }
   void SetPulserFrequency(const Int_t f=fgPulserFrequency) { fPulserFrequency = f; }
   
@@ -76 +125 @@
   MHGausEvents  &GetAverageLoGainSector(UInt_t i);
   const MHGausEvents   &GetAverageLoGainSector(UInt_t i)  const;
+
+  virtual Bool_t SetupFill(const MParList *pList);
+  virtual Bool_t ReInit   (      MParList *pList);
+  virtual Bool_t Fill     (const MParContainer *par, const Stat_t w=1);
+  virtual Bool_t Finalize ( );
 
   // Clone

trunk/MagicSoft/Mars/mcalib/MHCalibrationChargeCam.cc

@@ -127 +127 @@
 #include "MGeomPix.h"
 
+#include "MHGausEvents.h"
+
 #include "MBadPixelsCam.h"
 #include "MBadPixelsPix.h"
@@ -140 +142 @@
 using namespace std;
 
-const Int_t   MHCalibrationChargeCam::fgAverageNbins             = 4000;
 const Float_t MHCalibrationChargeCam::fgNumHiGainSaturationLimit = 0.01;
 const Float_t MHCalibrationChargeCam::fgNumLoGainSaturationLimit = 0.005;
@@ -150 +151 @@
 // - all pointers to NULL
 //
-// Initializes and sets owner of (done in MHCalibrationCam):
-// - fHiGainArray, fLoGainArray
-// - fAverageHiGainAreas, fAverageLoGainAreas
-// - fAverageHiGainSectors, fAverageLoGainSectors
-//
 // Initializes:
 // - fNumHiGainSaturationLimit to fgNumHiGainSaturationLimit
 // - fNumLoGainSaturationLimit to fgNumLoGainSaturationLimit
-// - fPulserFrequency to fgPulserFrequency (in MHCalibrationCam) 
 //
 MHCalibrationChargeCam::MHCalibrationChargeCam(const char *name, const char *title)
-    : fCam(NULL), fRawEvt(NULL), fGeom(NULL), fBadPixels(NULL)
+    : fRawEvt(NULL)
 {
     fName  = name  ? name  : "MHCalibrationChargeCam";
@@ -174 +169 @@
 // Gets the pointers to:
 // - MRawEvtData
-// - MGeomCam
-//
-// Calls Delete-Function of:
-// - MHCalibrationCam::fHiGainArray, MHCalibrationCam::fLoGainArray
-// - MHCalibrationCam::fAverageHiGainAreas, MHCalibrationCam::fAverageLoGainAreas
-// - MHCalibrationCam::fAverageHiGainSectors, MHCalibrationCam::fAverageLoGainSectors
-//
-Bool_t MHCalibrationChargeCam::SetupFill(const MParList *pList)
+//
+Bool_t MHCalibrationChargeCam::SetupHists(const MParList *pList)
 {
   
@@ -191 +180 @@
   }
 
-  fGeom = (MGeomCam*)pList->FindObject("MGeomCam");
-  if (!fGeom)
-  {
-      *fLog << err << "MGeomCam not found... aborting." << endl;
-      return kFALSE;
-  }
-
-  fHiGainArray->Delete();
-  fLoGainArray->Delete();
-
-  fAverageHiGainAreas->Delete();
-  fAverageLoGainAreas->Delete();
-
-  fAverageHiGainSectors->Delete();
-  fAverageLoGainSectors->Delete();
-
   return kTRUE;
 }
@@ -213 +186 @@
 //
 // Gets or creates the pointers to:
-// - MBadPixelsCam
 // - MCalibrationChargeCam
 //
@@ -219 +191 @@
 // - MExtractedSignalCam
 //
-// Initializes, if empty to MExtractedSignalCam::GetSize() for:
-// - MHCalibrationCam::fHiGainArray, MHCalibrationCam::fLoGainArray
-//
-// Initializes, if empty to MGeomCam::GetNumAreas() for:
-// - MHCalibrationCam::fAverageHiGainAreas, MHCalibrationCam::fAverageLoGainAreas
-//
-// Initializes, if empty to MGeomCam::GetNumSectors() for:
-// - MHCalibrationCam::fAverageHiGainSectors, MHCalibrationCam::fAverageLoGainSectors
-// 
-// Initializes TArrays to MGeomCam::GetNumAreas and MGeomCam::GetNumSectors, respectively
-// Fills with number of valid pixels (if !MBadPixelsPix::IsBad()):
-// - MHCalibrationCam::fAverageAreaNum[area index]
-// - MHCalibrationCam::fAverageSectorNum[area index]
-//
-// Calls InitializeHiGainHists() and InitializeLoGainHists() for every entry in:
+// Calls InitializeHists() for every entry in:
 // - MHCalibrationCam::fHiGainArray, MHCalibrationCam::fLoGainArray
 // - MHCalibrationCam::fAverageHiGainAreas, MHCalibrationCam::fAverageLoGainAreas
@@ -239 +197 @@
 //
 // Sets Titles and Names for the Charge Histograms and 
-// Sets number of bins to fAverageNbins for:
+// Sets number of bins to MHCalibrationCam::fAverageNbins for:
 // - MHCalibrationCam::fAverageHiGainAreas, MHCalibrationCam::fAverageLoGainAreas
 // - MHCalibrationCam::fAverageHiGainSectors, MHCalibrationCam::fAverageLoGainSectors
 // 
-Bool_t MHCalibrationChargeCam::ReInit(MParList *pList)
+Bool_t MHCalibrationChargeCam::ReInitHists(MParList *pList)
 {
 
-  fBadPixels = (MBadPixelsCam*)pList->FindCreateObj("MBadPixelsCam");
-  if (!fBadPixels)
-      return kFALSE;
-
-
-  fCam = (MCalibrationChargeCam*)pList->FindCreateObj("MCalibrationChargeCam");
+  fCam = (MCalibrationCam*)pList->FindCreateObj("MCalibrationChargeCam");
   if (!fCam)
       return kFALSE;
@@ -262 +215 @@
   }
 
-  const Int_t n        = signal->GetSize();
+  const Int_t npixels  = fGeom->GetNumPixels();
   const Int_t nsectors = fGeom->GetNumSectors();
   const Int_t nareas   = fGeom->GetNumAreas();
 
-  //
-  // The function TArrayF::Set() already sets all entries to 0.
-  //
-  fAverageAreaNum.        Set(nareas);
-  fAverageAreaSat.        Set(nareas);           
-  fAverageAreaSigma.      Set(nareas);      
-  fAverageAreaSigmaErr.   Set(nareas);   
-  fAverageAreaRelSigma.   Set(nareas);   
-  fAverageAreaRelSigmaErr.Set(nareas);
-  fAverageSectorNum.      Set(nsectors);
-
-  for (Int_t i=0; i<n; i++)
-    {
-      if ((*fBadPixels)[i].IsBad())
-        continue;
-      fAverageAreaNum  [(*fGeom)[i].GetAidx()  ]++;
-      fAverageSectorNum[(*fGeom)[i].GetSector()]++;
-    }
-
   if (fHiGainArray->GetEntries()==0)
   {
-      fHiGainArray->Expand(n);
-      for (Int_t i=0; i<n; i++)
+      fHiGainArray->Expand(npixels);
+      for (Int_t i=0; i<npixels; i++)
       {
           (*fHiGainArray)[i] = new MHCalibrationChargeHiGainPix;
-          InitializeHiGainHists((MHCalibrationChargePix&)(*this)[i],(*fBadPixels)[i],i);
+          InitHists((*this)[i],(*fBadPixels)[i],i);
       }
   }
 
- 
   if (fLoGainArray->GetEntries()==0)
   {
-      fLoGainArray->Expand(n);
-      
-      for (Int_t i=0; i<n; i++)
+      fLoGainArray->Expand(npixels);
+      
+      for (Int_t i=0; i<npixels; i++)
       {
           (*fLoGainArray)[i] = new MHCalibrationChargeLoGainPix;
-          InitializeLoGainHists((MHCalibrationChargePix&)(*this)(i),(*fBadPixels)[i],i);
+          InitHists((*this)(i),(*fBadPixels)[i],i);
       }
       
@@ -318 +251 @@
                                            "Average HiGain FADC sums area idx ");
 
+        InitHists(GetAverageHiGainArea(j),fCam->GetAverageBadArea(j),j);
+
         MHCalibrationChargePix &hist = (MHCalibrationChargePix&)GetAverageHiGainArea(j);
 
-        InitializeHiGainHists(hist,fCam->GetAverageBadArea(j),j);
-
         hist.GetHGausHist()->SetTitle("Summed FADC slices average HiGain Area Idx ");
-        hist.SetChargeNbins(fAverageNbins);
+        hist.SetNbins(fAverageNbins);
         hist.GetHAbsTime()->SetTitle("Absolute Arrival Time average HiGain Area Idx ");
       }
@@ -340 +273 @@
         MHCalibrationChargePix &hist = (MHCalibrationChargePix&)GetAverageLoGainArea(j);
 
-        InitializeLoGainHists(hist,fCam->GetAverageBadArea(j),j);
+        InitHists(hist,fCam->GetAverageBadArea(j),j);
 
         hist.GetHGausHist()->SetTitle("Summed FADC slices average LoGain Area Idx ");
-        hist.SetChargeNbins(fAverageNbins);
+        hist.SetNbins(fAverageNbins);
         hist.GetHAbsTime()->SetTitle("Absolute Arrival Time average LoGain Area Idx ");
 
@@ -361 +294 @@
           MHCalibrationChargePix &hist = (MHCalibrationChargePix&)GetAverageHiGainSector(j);
 
-          InitializeHiGainHists(hist,fCam->GetAverageBadSector(j),j);
+          InitHists(hist,fCam->GetAverageBadSector(j),j);
           
           hist.GetHGausHist()->SetTitle("Summed FADC slices average HiGain Sector ");
-          hist.SetChargeNbins(fAverageNbins);
+          hist.SetNbins(fAverageNbins);
           hist.GetHAbsTime()->SetTitle("Absolute Arrival Time average HiGain Sector ");
 
@@ -382 +315 @@
           MHCalibrationChargePix &hist = (MHCalibrationChargePix&)GetAverageLoGainSector(j);
 
-          InitializeLoGainHists(hist,fCam->GetAverageBadSector(j),j);
+          InitHists(hist,fCam->GetAverageBadSector(j),j);
           
           hist.GetHGausHist()->SetTitle("Summed FADC slices average LoGain Sector ");
-          hist.SetChargeNbins(fAverageNbins);
+          hist.SetNbins(fAverageNbins);
           hist.GetHAbsTime()->SetTitle("Absolute Arrival Time average LoGain Sector ");
 
@@ -394 +327 @@
 }
 
-// -------------------------------------------------------------
-//
-// If MBadPixelsPix::IsBad():
-// - calls MHCalibrationChargePix::SetExcluded()
-//
-// Calls:
-// - MHCalibrationChargePix::Init()
-// - MHCalibrationChargePix::ChangeHistId(i)
-// - MHGausEvents::SetEventFrequency(fPulserFrequency)
-// 
-void MHCalibrationChargeCam::InitializeHiGainHists(MHCalibrationChargePix &hist, MBadPixelsPix &bad, const Int_t i)
-{
-  
-  if (bad.IsBad())
-    {
-      *fLog << warn << "Excluded pixel: " << i << " from calibration " << endl;
-      hist.SetExcluded();
-    }
-  
-  hist.Init();
-  hist.ChangeHistId(i);
-  hist.SetEventFrequency(fPulserFrequency);
-          
-}
-
-// -------------------------------------------------------------
-//
-// If MBadPixelsPix::IsBad():
-// - calls MHCalibrationChargePix::SetExcluded()
-//
-// Calls:
-// - MHCalibrationChargePix::Init()
-// - MHCalibrationChargePix::ChangeHistId(i)
-// - MHGausEvents::SetEventFrequency(fPulserFrequency)
-// 
-void MHCalibrationChargeCam::InitializeLoGainHists(MHCalibrationChargePix &hist,MBadPixelsPix &bad, const Int_t i)
-{
-
-  if (bad.IsBad())
-    hist.SetExcluded();
-  
-  hist.Init();
-  hist.ChangeHistId(i);
-  hist.SetEventFrequency(fPulserFrequency);
-}
-
   
 // --------------------------------------------------------------------------
 //
 // Retrieves from MExtractedSignalCam:
+// - first used LoGain FADC slice
+//
+// Retrieves from MGeomCam:
 // - number of pixels
-// - first used LoGain FADC slice
-//
-// Retrieves from MGeomCam:
 // - number of pixel areas
 // - number of sectors
@@ -453 +340 @@
 // For all TObjArray's (including the averaged ones), the following steps are performed: 
 //
-// 1) Test size and return kFALSE if not matching
-//
-// 2) Fill Charges histograms (MHGausEvents::FillHistAndArray()) with:
+// 1) Fill Charges histograms (MHGausEvents::FillHistAndArray()) with:
 // - MExtractedSignalPix::GetExtractedSignalHiGain();
 // - MExtractedSignalPix::GetExtractedSignalLoGain();
 //
-// 3) Set number of saturated slices (MHCalibrationChargePix::SetSaturated()) with:
+// 2) Set number of saturated slices (MHCalibrationChargePix::SetSaturated()) with:
 // - MExtractedSignalPix::GetNumHiGainSaturated();
 // - MExtractedSignalPix::GetNumLoGainSaturated();
 //
-// 4) Fill AbsTime histograms (MHCalibrationChargePix::FillAbsTime()) with:
+// 3) Fill AbsTime histograms (MHCalibrationChargePix::FillAbsTime()) with:
 // - MRawEvtPixelIter::GetIdxMaxHiGainSample();       
 // - MRawEvtPixelIter::GetIdxMaxLoGainSample(first slice);
 //
-Bool_t MHCalibrationChargeCam::Fill(const MParContainer *par, const Stat_t w)
+Bool_t MHCalibrationChargeCam::FillHists(const MParContainer *par, const Stat_t w)
 {
 
@@ -477 +362 @@
     }
   
-  const Int_t npixels  = signal->GetSize();
-  const Int_t lofirst  = signal->GetFirstUsedSliceLoGain();
+  const Int_t npixels  = fGeom->GetNumPixels();
   const Int_t nareas   = fGeom->GetNumAreas();
   const Int_t nsectors = fGeom->GetNumSectors();
-
-  if (fHiGainArray->GetEntries() != npixels)
-    {
-      *fLog << err << "ERROR - Size mismatch in number of pixels ... abort." << endl;
-      return kFALSE;
-    }
-
-  if (fLoGainArray->GetEntries() != npixels)
-    {
-      *fLog << err << "ERROR - Size mismatch in number of pixels ... abort." << endl;
-      return kFALSE;
-    }
-
-  if (fAverageHiGainAreas->GetEntries() != nareas)
-    {
-      *fLog << err << "ERROR - Size mismatch in number of areas ... abort." << endl;
-      return kFALSE;
-    }
-
-  if (fAverageLoGainAreas->GetEntries() != nareas)
-    {
-      *fLog << err << "ERROR - Size mismatch in number of areas ... abort." << endl;
-      return kFALSE;
-    }
-
-  if (fAverageHiGainSectors->GetEntries() != nsectors)
-    {
-      *fLog << err << "ERROR - Size mismatch in number of sectors ... abort." << endl;
-      return kFALSE;
-    }
-
-  if (fAverageLoGainSectors->GetEntries() != nsectors)
-    {
-      *fLog << err << "ERROR - Size mismatch in number of sectors ... abort." << endl;
-      return kFALSE;
-    }
+  const Int_t lofirst  = signal->GetFirstUsedSliceLoGain();
 
   Float_t sumhiarea  [nareas],   sumloarea  [nareas],   timehiarea  [nareas],   timeloarea  [nareas];
@@ -526 +375 @@
     {
       sumhiarea  [j] = sumloarea  [j] = timehiarea  [j] =  timeloarea  [j] = 0.;
+      sathiarea  [j] = satloarea  [j] = 0;
+    }
+  
+  for (UInt_t j=0; j<nsectors; j++)
+    {
       sumhisector[j] = sumlosector[j] = timehisector[j] =  timelosector[j] = 0.;
-      sathiarea  [j] = satloarea  [j] = 0;
       sathisector[j] = satlosector[j] = 0;
     }
   
+
   for (Int_t i=0; i<npixels; i++)
     {
@@ -596 +450 @@
     }
   
-
   for (UInt_t j=0; j<nareas; j++)
     {
@@ -642 +495 @@
 // For all TObjArray's (including the averaged ones), the following steps are performed: 
 //
-// 1) Return if the pixel is excluded
-// 
-// 2) Call to FinalizeHiGainHists() and FinalizeLoGainHists()
-//
-// For all averaged areas, the fitted sigma is multiplied with the square root of 
-// the number involved pixels
-//
-Bool_t MHCalibrationChargeCam::Finalize()
+// 1) Returns if the pixel is excluded.
+// 2) Tests saturation. In case yes, set the flag: MCalibrationPix::SetHiGainSaturation()
+//    or the flag: MBadPixelsPix::SetUncalibrated( MBadPixelsPix::kLoGainSaturated )
+// 3) Store the absolute arrival times in the MCalibrationChargePix's. If flag 
+//    MCalibrationPix::IsHiGainSaturation() is set, the Low-Gain arrival times are stored, 
+//    otherwise the Hi-Gain ones.
+// 4) Calls to MHCalibrationCam::FitHiGainArrays() and MCalibrationCam::FitLoGainArrays() 
+//    and sets the flags (if  occurring):
+//    - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::kHiGainNotFitted )
+//    - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::kLoGainNotFitted )
+//    - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::kHiGainOscillating )
+//    - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::kLoGainOscillating )
+// 5) Retrieves the results and store them in MCalibrationChargePix
+//
+Bool_t MHCalibrationChargeCam::FinalizeHists()
 {
 
@@ -656 +516 @@
       
       MHCalibrationChargePix &histhi = (MHCalibrationChargePix&)(*this)[i];
-
+      MCalibrationChargePix  &pix    = (MCalibrationChargePix&)(*fCam)[i];
+      
       if (histhi.IsExcluded())
         continue;
       
-      MCalibrationChargePix  &pix    = (*fCam)[i];
+      if (histhi.GetSaturated() > fNumHiGainSaturationLimit*histhi.GetHGausHist()->GetEntries())
+        {
+          pix.SetHiGainSaturation();
+          histhi.CreateFourierSpectrum();
+          continue;
+        }
+
+      pix.SetAbsTimeMean ( histhi.GetAbsTimeMean());
+      pix.SetAbsTimeRms  ( histhi.GetAbsTimeRms() );
+    }
+
+  for (Int_t i=0; i<fLoGainArray->GetSize(); i++)
+    {
+      
+      MHCalibrationChargePix &histlo = (MHCalibrationChargePix&)(*this)(i);
       MBadPixelsPix          &bad    = (*fBadPixels)[i];
-      
-      FinalizeHiGainHists(histhi,pix,bad);
-      
-    }
-  
-  for (Int_t i=0; i<fLoGainArray->GetSize(); i++)
-    {
-      
-      MHCalibrationChargePix &histlo = (MHCalibrationChargePix&)(*this)(i);
 
       if (histlo.IsExcluded())
         continue;
       
-      MCalibrationChargePix        &pix    = (*fCam)[i];
-      MBadPixelsPix                &bad    = (*fBadPixels)[i];
-      
-      FinalizeLoGainHists(histlo,pix,bad);
-      
-    }
-  
+      if (histlo.GetSaturated() > fNumLoGainSaturationLimit*histlo.GetHGausHist()->GetEntries())
+        {
+          *fLog << warn << "Saturated Lo Gain histogram in pixel: " << i << endl;
+          bad.SetUncalibrated( MBadPixelsPix::kLoGainSaturation ); 
+          histlo.CreateFourierSpectrum();
+          continue;
+        }
+  
+      MCalibrationChargePix &pix    = (MCalibrationChargePix&)(*fCam)[i];
+      
+      if (pix.IsHiGainSaturation())
+        {
+          pix.SetAbsTimeMean     ( histlo.GetAbsTimeMean());
+          pix.SetAbsTimeRms      ( histlo.GetAbsTimeRms() );
+        }           
+    }
+
   for (Int_t j=0; j<fAverageHiGainAreas->GetSize(); j++)
     {
       
       MHCalibrationChargePix &histhi = (MHCalibrationChargePix&)GetAverageHiGainArea(j);      
-      MCalibrationChargePix  &pix    = fCam->GetAverageArea(j);
-      MBadPixelsPix          &bad    = fCam->GetAverageBadArea(j);        
-      
-      FinalizeHiGainHists(histhi,pix,bad);
+      MCalibrationChargePix  &pix    = (MCalibrationChargePix&)fCam->GetAverageArea(j);
+      
+      if (histhi.GetSaturated() > fNumHiGainSaturationLimit*histhi.GetHGausHist()->GetEntries())
+        {
+          pix.SetHiGainSaturation();
+          histhi.CreateFourierSpectrum();
+          continue;
+        }
+
+      pix.SetAbsTimeMean ( histhi.GetAbsTimeMean());
+      pix.SetAbsTimeRms  ( histhi.GetAbsTimeRms() );
     }
   
@@ -696 +579 @@
       
       MHCalibrationChargePix &histlo = (MHCalibrationChargePix&)GetAverageLoGainArea(j);      
-      MCalibrationChargePix  &pix    = fCam->GetAverageArea(j);
-      MBadPixelsPix          &bad    = fCam->GetAverageBadArea(j);        
-      
-      FinalizeLoGainHists(histlo,pix,bad);
+      MCalibrationChargePix  &pix    = (MCalibrationChargePix&)fCam->GetAverageArea(j);
+      
+      if (histlo.GetSaturated() > fNumLoGainSaturationLimit*histlo.GetHGausHist()->GetEntries())
+        {
+          *fLog << warn << "Saturated Lo Gain histogram in area idx: " << j << endl;
+          histlo.CreateFourierSpectrum();
+          continue;
+        }
+
+      if (pix.IsHiGainSaturation())
+        {
+          pix.SetAbsTimeMean  ( histlo.GetAbsTimeMean());
+          pix.SetAbsTimeRms   ( histlo.GetAbsTimeRms() );
+        }           
     }
 
@@ -706 +599 @@
       
       MHCalibrationChargePix &histhi = (MHCalibrationChargePix&)GetAverageHiGainSector(j);      
-      MCalibrationChargePix  &pix    = fCam->GetAverageSector(j);
+      MCalibrationChargePix  &pix    = (MCalibrationChargePix&)fCam->GetAverageSector(j);
+      
+      if (histhi.GetSaturated() > fNumHiGainSaturationLimit*histhi.GetHGausHist()->GetEntries())
+        {
+          pix.SetHiGainSaturation();
+          histhi.CreateFourierSpectrum();
+          continue;
+        }
+
+      pix.SetAbsTimeMean ( histhi.GetAbsTimeMean());
+      pix.SetAbsTimeRms  ( histhi.GetAbsTimeRms() );
+    }
+  
+  for (Int_t j=0; j<fAverageLoGainSectors->GetSize(); j++)
+    {
+      
+      MHCalibrationChargePix &histlo = (MHCalibrationChargePix&)GetAverageLoGainSector(j);      
+      MCalibrationChargePix  &pix    = (MCalibrationChargePix&)fCam->GetAverageSector(j);
       MBadPixelsPix          &bad    = fCam->GetAverageBadSector(j);        
       
-      FinalizeHiGainHists(histhi,pix,bad);
-    }
-  
-  for (Int_t j=0; j<fAverageLoGainSectors->GetSize(); j++)
-    {
-      
-      MHCalibrationChargePix &histlo = (MHCalibrationChargePix&)GetAverageLoGainSector(j);      
-      MCalibrationChargePix  &pix    = fCam->GetAverageSector(j);
-      MBadPixelsPix          &bad    = fCam->GetAverageBadSector(j);        
-      
-      FinalizeLoGainHists(histlo,pix,bad);
-    }
-  
-  for (Int_t j=0; j<fGeom->GetNumAreas();j++)
-    {
-      
-      MCalibrationChargePix &pix = fCam->GetAverageArea(j);
+      if (histlo.GetSaturated() > fNumLoGainSaturationLimit*histlo.GetHGausHist()->GetEntries())
+        {
+          *fLog << warn << "Saturated Lo Gain histogram in sector: " << j << endl;
+          bad.SetUncalibrated( MBadPixelsPix::kLoGainSaturation ); 
+          histlo.CreateFourierSpectrum();
+          continue;
+        }
 
       if (pix.IsHiGainSaturation())
-        fAverageAreaSat[j]++;
-
-      fAverageAreaSigma[j]    = pix.GetSigmaCharge    () * TMath::Sqrt((Float_t)fAverageAreaNum[j]);
-      fAverageAreaSigmaErr[j] = pix.GetSigmaChargeErr () * TMath::Sqrt((Float_t)fAverageAreaNum[j]);
-
-      pix.SetSigmaCharge   (fAverageAreaSigma[j]);
-      pix.SetSigmaChargeErr(fAverageAreaSigmaErr[j]);
-
-      fAverageAreaRelSigma[j]   = fAverageAreaSigma[j] / pix.GetMeanCharge();
-      
-      Float_t relsigmaerr       =  fAverageAreaSigmaErr[j]*fAverageAreaSigmaErr[j] 
-                                / (fAverageAreaSigma[j]   *fAverageAreaSigma[j]   );
-      relsigmaerr               += pix.GetMeanChargeErr()*pix.GetMeanChargeErr() 
-                                / (pix.GetMeanCharge()   *pix.GetMeanCharge()   );
-      relsigmaerr               *= fAverageAreaRelSigma[j];
-      fAverageAreaRelSigmaErr[j] = TMath::Sqrt(relsigmaerr);
-
-    }
-
+        {
+          pix.SetAbsTimeMean  ( histlo.GetAbsTimeMean());
+          pix.SetAbsTimeRms   ( histlo.GetAbsTimeRms() );
+        }           
+    }
+  
+  //
+  // Perform the fitting for the High Gain (done in MHCalibrationCam)
+  //
+  FitHiGainArrays((MCalibrationCam&)(*fCam),(*fBadPixels),
+                  MBadPixelsPix::kHiGainNotFitted,
+                  MBadPixelsPix::kHiGainOscillating);
+  //
+  // Perform the fitting for the Low Gain (done in MHCalibrationCam)
+  //
+  FitLoGainArrays((MCalibrationCam&)(*fCam),(*fBadPixels),
+                  MBadPixelsPix::kLoGainNotFitted,
+                  MBadPixelsPix::kLoGainOscillating);
+      
   return kTRUE;
 }
 
-// ---------------------------------------------------------------------------
-//
-// Returns if the histogram is empty and sets the following flag:
-// - MBadPixelsPix::SetUnsuitable(MBadPixelsPix::kUnsuitableRun)
-//
-// Returns if the number of saturated slices (MHCalibrationChargePix::GetSaturated()) is 
-// higher than the allowed limit (fNumHiGainSaturationLimit*histogram-entries), creates 
-// the fourier spectrum and sets the following flags:
-// - MCalibrationChargePix::SetHiGainSaturation()
-//
-// Fits the histograms with a Gaussian, in case of failure 
-// calls MHCalibrationChargePix::RepeatFit(), in case of repeated failure 
-// calls MHGausEvents::BypassFit() and sets the following flags:
-// - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::kHiGainNotFitted )
-// - MBadPixelsPix::SetUnsuitable(   MBadPixelsPix::kUnreliableRun   )
-// 
-// Counts the number of pickup events
-//
-// Creates the fourier spectrum and tests MHGausEvents::IsFourierSpectrumOK(). 
-// In case no, sets the following flags:
-// - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::kHiGainOscillating )
-// - MBadPixelsPix::SetUnsuitable(   MBadPixelsPix::kUnreliableRun     )
-//
-// Retrieves the results and store them in MCalibrationChargePix
-//
-void MHCalibrationChargeCam::FinalizeHiGainHists(MHCalibrationChargePix &hist, 
-                                                 MCalibrationChargePix &pix, 
-                                                 MBadPixelsPix &bad)
+// --------------------------------------------------------------------------
+//
+// Takes the decisions under which conditions a pixel is declared:
+// MBadPixelsPix::kUnsuitableRun  or MBadPixelsPix::kUnreliableRun, namely:
+// * if MBadPixelsPix::kHiGainNotFitted   and !MCalibrationPix::IsHiGainSaturation()
+//   sets MBadPixelsPix::kUnreliableRun 
+// * if MBadPixelsPix::kHiGainOscillating and !MCalibrationPix::IsHiGainSaturation()
+//   sets MBadPixelsPix::kUnreliableRun 
+// * if MBadPixelsPix::kLoGainNotFitted   and  MCalibrationPix::IsLoGainSaturation()
+//   sets MBadPixelsPix::kUnreliableRun 
+// * if MBadPixelsPix::kLoGainOscillating and  MCalibrationPix::IsLoGainSaturation()
+//   sets MBadPixelsPix::kUnreliableRun 
+// * if MBadPixelsPix::kLoGainSaturation
+//   sets MBadPixelsPix::kUnsuitableRun 
+//
+void MHCalibrationChargeCam::FinalizeBadPixels()
 {
-
-    if (hist.IsEmpty())
-      {
-        *fLog << warn << "Empty Hi Gain histogram in pixel: " << pix.GetPixId() << endl;
-        bad.SetUnsuitable(MBadPixelsPix::kUnsuitableRun);        
-        return;
-      }
-    
-    if (hist.GetSaturated() > fNumHiGainSaturationLimit*hist.GetHGausHist()->GetEntries())
-    {
-      *fLog << warn << "Saturated Hi Gain histogram in pixel: " << pix.GetPixId() << endl;
-      pix.SetHiGainSaturation();
-      hist.CreateFourierSpectrum();
-      return;
-    }
-
-    //
-    // 2) Fit the Hi Gain histograms with a Gaussian
-    //
-    if (!hist.FitGaus())
-    //
-    // 3) In case of failure set the bit Fitted to false and take histogram means and RMS
-    //
-      if (!hist.RepeatFit())
-        {
-          hist.BypassFit();
-          *fLog << warn << "Hi Gain could not be fitted in pixel: " << pix.GetPixId() << endl;
-          bad.SetUncalibrated( MBadPixelsPix::kHiGainNotFitted );
-          bad.SetUnsuitable(   MBadPixelsPix::kUnreliableRun   );
-        }
-
-    //
-    // 4) Check for pickup
-    //
-    hist.CountPickup();
-    
-    //
-    // 5) Check for oscillations
-    // 
-    hist.CreateFourierSpectrum();
-    
-    if (!hist.IsFourierSpectrumOK())
-      {
-        *fLog << warn << "Oscillating Hi Gain in pixel: " << pix.GetPixId() << endl;
-        bad.SetUncalibrated( MBadPixelsPix::kHiGainOscillating );
-        bad.SetUnsuitable(   MBadPixelsPix::kUnreliableRun     );
-      }
-
-    //
-    // 6) Retrieve the results and store them in this class
-    //
-    pix.SetHiGainMeanCharge(     hist.GetMean()      );
-    pix.SetHiGainMeanChargeErr(  hist.GetMeanErr()   );
-    pix.SetHiGainSigmaCharge(    hist.GetSigma()     );
-    pix.SetHiGainSigmaChargeErr( hist.GetSigmaErr()  );
-    pix.SetHiGainChargeProb    ( hist.GetProb()      );
-    
-    pix.SetAbsTimeMean         ( hist.GetAbsTimeMean());
-    pix.SetAbsTimeRms          ( hist.GetAbsTimeRms() );
-    
-    pix.SetHiGainNumPickup     ( hist.GetPickup()     );
-    
+      
+  for (Int_t i=0; i<fBadPixels->GetSize(); i++)
+    {
+      
+      MBadPixelsPix    &bad    = (*fBadPixels)[i];
+      MCalibrationPix  &pix    = (*fCam)[i];
+
+      if (bad.IsUncalibrated( MBadPixelsPix::kHiGainNotFitted ))
+        if (!pix.IsHiGainSaturation())
+          bad.SetUnsuitable(   MBadPixelsPix::kUnreliableRun    );
+ 
+      if (bad.IsUncalibrated( MBadPixelsPix::kHiGainOscillating ))
+        bad.SetUnsuitable(   MBadPixelsPix::kUnreliableRun    );
+
+      if (bad.IsUncalibrated( MBadPixelsPix::kLoGainNotFitted ))
+        if (pix.IsHiGainSaturation())
+          bad.SetUnsuitable(   MBadPixelsPix::kUnreliableRun    );
+ 
+      if (bad.IsUncalibrated( MBadPixelsPix::kLoGainOscillating ))
+        if (pix.IsHiGainSaturation())
+          bad.SetUnsuitable(   MBadPixelsPix::kUnreliableRun    );
+
+      if (bad.IsUncalibrated( MBadPixelsPix::kLoGainSaturation ))
+          bad.SetUnsuitable(   MBadPixelsPix::kUnsuitableRun    );
+    }
 }
-
-// ---------------------------------------------------------------------------
-//
-// Returns if the histogram is empty. If yes and the flag MCalibrationChargePix::IsHiGainSaturation()
-// is set, sets the following flag:
-// - MBadPixelsPix::SetUnsuitable(MBadPixelsPix::kUnsuitableRun)
-//
-// Returns if the number of saturated slices (MHCalibrationChargePix::GetSaturated()) is 
-// higher than the allowed limit (fNumLoGainSaturationLimit*histogram-entries), creates 
-// the fourier spectrum and sets the following flags:
-// - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::kLoGainSaturation ); 
-// - MBadPixelsPix::SetUnsuitable(   MBadPixelsPix::kUnsuitableRun    );
-//
-// Fits the histograms with a Gaussian, in case of failure calls MHCalibrationChargePix::RepeatFit(), 
-// in case of repeated failure calls MHGausEvents::BypassFit() and sets the following flags:
-// - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::kLoGainNotFitted )
-//   If the flag MCalibrationChargePix::IsHiGainSaturation() is set, sets additionally:
-// - MBadPixelsPix::SetUnsuitable(   MBadPixelsPix::kUnreliableRun   )
-// 
-// Counts the number of pickup events
-//
-// Creates the fourier spectrum and tests MHGausEvents::IsFourierSpectrumOK(). 
-// In case no, sets the following flags:
-// - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::kLoGainOscillating )
-//   If the flag MCalibrationChargePix::IsHiGainSaturation() is set, sets additionally:
-// - MBadPixelsPix::SetUnsuitable(   MBadPixelsPix::kUnreliableRun     )
-//
-// Retrieves the results and store them in MCalibrationChargePix
-//
-void MHCalibrationChargeCam::FinalizeLoGainHists(MHCalibrationChargePix &hist, 
-                                                 MCalibrationChargePix &pix, 
-                                                 MBadPixelsPix &bad)
-{
-
-    if (hist.IsEmpty())
-      {
-        if (pix.IsHiGainSaturation())
-          bad.SetUnsuitable(MBadPixelsPix::kUnsuitableRun);
-        return;
-      }
-
-    if (hist.GetSaturated() > fNumLoGainSaturationLimit*hist.GetHGausHist()->GetEntries())
-    {
-      *fLog << warn << "Saturated Lo Gain histogram in pixel: " << pix.GetPixId() << endl;
-      bad.SetUncalibrated( MBadPixelsPix::kLoGainSaturation ); 
-      bad.SetUnsuitable(   MBadPixelsPix::kUnsuitableRun    );
-      hist.CreateFourierSpectrum();
-      return;
-    }
-
-    //
-    // 2) Fit the Lo Gain histograms with a Gaussian
-    //
-    if (!hist.FitGaus())
-    //
-    // 3) In case of failure set the bit kFitted to false and take histogram means and RMS
-    //
-      if (!hist.RepeatFit())
-        {
-          hist.BypassFit();
-          bad.SetUncalibrated( MBadPixelsPix::kLoGainNotFitted );
-          if (pix.IsHiGainSaturation())
-            bad.SetUnsuitable(MBadPixelsPix::kUnreliableRun);
-        }
-
-    //
-    // 4) Check for pickup
-    //
-    hist.CountPickup();
-    
-    //
-    // 5) Check for oscillations
-    // 
-    hist.CreateFourierSpectrum();
-    
-    if (!hist.IsFourierSpectrumOK())
-      {
-        bad.SetUncalibrated( MBadPixelsPix::kLoGainOscillating );
-        if (pix.IsHiGainSaturation())
-          bad.SetUnsuitable(MBadPixelsPix::kUnreliableRun);
-      }
-    //
-    // 6) Retrieve the results and store them in this class
-    //
-    pix.SetLoGainMeanCharge(     hist.GetMean()     );
-    pix.SetLoGainMeanChargeErr(  hist.GetMeanErr()  );
-    pix.SetLoGainSigmaCharge(    hist.GetSigma()    );
-    pix.SetLoGainSigmaChargeErr( hist.GetSigmaErr() );
-    pix.SetLoGainChargeProb    ( hist.GetProb()     );
-    
-    if (pix.IsHiGainSaturation())
-    {
-        pix.SetAbsTimeMean     ( hist.GetAbsTimeMean());
-        pix.SetAbsTimeRms      ( hist.GetAbsTimeRms() );
-    }       
-    
-    pix.SetLoGainNumPickup     (  hist.GetPickup()    );
-
-}    
 
 // --------------------------------------------------------------------------

trunk/MagicSoft/Mars/mcalib/MHCalibrationChargeCam.h

@@ -7 +7 @@
 
 class MRawEvtData;
-class MGeomCam;
-class MBadPixelsCam;
 class MBadPixelsPix;
 class MCalibrationChargeCam;
@@ -17 +15 @@
 private:
 
-  static const Int_t   fgAverageNbins;               // The default for fAverageNbins
   static const Float_t fgNumHiGainSaturationLimit;   // The default for fNumHiGainSaturationLimit
   static const Float_t fgNumLoGainSaturationLimit;   // The default for fNumLoGainSaturationLimit
   
-  Int_t   fAverageNbins;              // Number of bins for the average histograms
   Float_t fNumHiGainSaturationLimit;  // Rel. amount sat. higain FADC slices until pixel is called saturated
   Float_t fNumLoGainSaturationLimit;  // Rel. amount sat. logain FADC slices until pixel is called saturated
   
-  MCalibrationChargeCam  *fCam;       //!  Calibration Cam with the results
-  MRawEvtData            *fRawEvt;    //!  Raw event data 
-  MGeomCam               *fGeom;      //!  Camera geometry
-  MBadPixelsCam          *fBadPixels; //!  Bad Pixels storage container
+  MRawEvtData *fRawEvt;               //!  Raw event data 
 
-  void InitializeHiGainHists(MHCalibrationChargePix &hist, MBadPixelsPix &bad, const Int_t i);
-  void InitializeLoGainHists(MHCalibrationChargePix &hist, MBadPixelsPix &bad, const Int_t i);  
-  
+  Bool_t SetupHists(const MParList *pList);
+  Bool_t ReInitHists(MParList *pList);
+  Bool_t FillHists(const MParContainer *par, const Stat_t w=1);
+  Bool_t FinalizeHists();
+  void   FinalizeBadPixels();
+
   void FinalizeHiGainHists(MHCalibrationChargePix &hist, MCalibrationChargePix &pix, MBadPixelsPix &bad);
   void FinalizeLoGainHists(MHCalibrationChargePix &hist, MCalibrationChargePix &pix, MBadPixelsPix &bad);
@@ -41 +37 @@
   ~MHCalibrationChargeCam() {}
   
-
-  void SetAverageNbins(    const Int_t bins=fgAverageNbins  )          { fAverageNbins = bins; }
   void SetNumLoGainSaturationLimit( const Float_t lim=fgNumLoGainSaturationLimit) { fNumLoGainSaturationLimit = lim; }
   void SetNumHiGainSaturationLimit( const Float_t lim=fgNumHiGainSaturationLimit) { fNumHiGainSaturationLimit = lim; }
@@ -49 +43 @@
   Float_t GetNumLoGainSaturationLimit()      const  { return fNumLoGainSaturationLimit; }
 
-  Bool_t SetupFill(const MParList *pList);
-  Bool_t ReInit   (      MParList *pList);
-  Bool_t Fill     (const MParContainer *par, const Stat_t w=1);
-  Bool_t Finalize ( );
-
   Bool_t GetPixelContent ( Double_t &val, Int_t idx, const MGeomCam &cam, Int_t type=0) const;
   void   DrawPixelContent( Int_t num )  const;    
 
-  ClassDef(MHCalibrationChargeCam, 1)   // Histogram class for camera calibration 
+  ClassDef(MHCalibrationChargeCam, 1)   // Histogram class for camera charge calibration 
 };
 

trunk/MagicSoft/Mars/mcalib/MHCalibrationChargeHiGainPix.cc

@@ -60 +60 @@
 //
 // Sets: 
-// - the default number for fChargeNbins  (fgChargeNbins)
-// - the default number for fChargeFirst  (fgChargeFirst)
-// - the default number for fChargeLast   (fgChargeLast)
+// - the default number for fNbins        (fgChargeNbins)
+// - the default number for fFirst        (fgChargeFirst)
+// - the default number for fLast         (fgChargeLast)
 // - the default number for fAbsTimeNbins (fgAbstTimeNbins)
 // - the default number for fAbsTimeFirst (fgAbsTimeFirst)
@@ -79 +79 @@
   fTitle = title ? title : "Fill the FADC sums of the HiGainPix events and perform the fits Pixel ";
   
-  SetChargeNbins();
-  SetChargeFirst();
-  SetChargeLast();
+  SetNbins ( fgChargeNbins );
+  SetFirst ( fgChargeFirst );
+  SetLast  ( fgChargeLast  );
   
   SetAbsTimeNbins();
@@ -87 +87 @@
   SetAbsTimeLast();
 
-  fHGausHist.SetName("HCalibrationChargeHiGainPix");
+  fHGausHist.SetName ("HCalibrationChargeHiGainPix");
   fHGausHist.SetTitle("Distribution of Summed Hi Gain FADC slices Pixel ");  
 
-  fHAbsTime.SetName("HAbsTimeHiGainPix");
+  fHAbsTime.SetName ("HAbsTimeHiGainPix");
   fHAbsTime.SetTitle("Distribution of Absolute Arrival Times Hi Gain Pixel ");  
 }

trunk/MagicSoft/Mars/mcalib/MHCalibrationChargeHiGainPix.h

@@ -12 +12 @@
 private:
 
-  static const Int_t   fgChargeNbins;     // Default for fChargeNBins  (now set to: 2000  )
-  static const Axis_t  fgChargeFirst;     // Default for fChargeFirst  (now set to: -0.5  )
-  static const Axis_t  fgChargeLast;      // Default for fChargeLast   (now set to: 1999.5)
-  static const Int_t   fgAbsTimeNbins;    // Default for fAbsTimeNbins (now set to: 20    )
-  static const Axis_t  fgAbsTimeFirst;    // Default for fAbsTimeFirst (now set to: -0.5  )
-  static const Axis_t  fgAbsTimeLast;     // Default for fAbsTimeLast  (now set to: 19.5  )
+  static const Int_t   fgChargeNbins;     // Default for MHGausEvents::fNbins (now set to: 2000  )
+  static const Axis_t  fgChargeFirst;     // Default for MHGausEvents::fFirst (now set to: -0.5  )
+  static const Axis_t  fgChargeLast;      // Default for MHGausEvents::fLast  (now set to: 1999.5)
+  static const Int_t   fgAbsTimeNbins;    // Default for fAbsTimeNbins        (now set to: 20    )
+  static const Axis_t  fgAbsTimeFirst;    // Default for fAbsTimeFirst        (now set to: -0.5  )
+  static const Axis_t  fgAbsTimeLast;     // Default for fAbsTimeLast         (now set to: 19.5  )
 
 public:
@@ -25 +25 @@
 
   // Setters
-  void SetChargeNbins(const Int_t  bins =fgChargeNbins)          { fChargeNbins = bins;     }
-  void SetChargeFirst(const Axis_t first=fgChargeFirst)          { fChargeFirst = first;    }
-  void SetChargeLast (const Axis_t last =fgChargeLast)           { fChargeLast  = last;     }
-  
   void SetAbsTimeNbins(const Int_t  bins =fgAbsTimeNbins)        { fAbsTimeNbins = bins;    }
   void SetAbsTimeFirst(const Axis_t first=fgAbsTimeFirst)        { fAbsTimeFirst = first;   }

trunk/MagicSoft/Mars/mcalib/MHCalibrationChargeLoGainPix.cc

@@ -60 +60 @@
 //
 // Sets: 
-// - the default number for fChargeNbins  (fgChargeNbins)
-// - the default number for fChargeFirst  (fgChargeFirst)
-// - the default number for fChargeLast   (fgChargeLast)
+// - the default number for fNbins        (fgChargeNbins)
+// - the default number for fFirst        (fgChargeFirst)
+// - the default number for fLast         (fgChargeLast)
 // - the default number for fAbsTimeNbins (fgAbstTimeNbins)
 // - the default number for fAbsTimeFirst (fgAbsTimeFirst)
@@ -79 +79 @@
   fTitle = title ? title : "Fill the FADC sums of the Low Gain events and perform the fits Pixel ";
   
-  SetChargeNbins();
-  SetChargeFirst();
-  SetChargeLast();
+  SetNbins ( fgChargeNbins );
+  SetFirst ( fgChargeFirst );
+  SetLast  ( fgChargeLast  );
   
   SetAbsTimeNbins();
@@ -87 +87 @@
   SetAbsTimeLast();
 
-  fHGausHist.SetName("HCalibrationChargeLoGainPix");
+  fHGausHist.SetName ("HCalibrationChargeLoGainPix");
   fHGausHist.SetTitle("Distribution of Summed Lo Gain FADC slices Pixel ");  
 
-  fHAbsTime.SetName("HAbsTimeLoGainPix");
+  fHAbsTime.SetName ("HAbsTimeLoGainPix");
   fHAbsTime.SetTitle("Distribution of Absolute Arrival Times Lo Gain Pixel ");  
 }

trunk/MagicSoft/Mars/mcalib/MHCalibrationChargeLoGainPix.h

@@ -11 +11 @@
 private:
 
-  static const Int_t   fgChargeNbins;      // Default for fChargeNBins  (now set to: 200   )
-  static const Axis_t  fgChargeFirst;      // Default for fChargeFirst  (now set to: -0.5  )
-  static const Axis_t  fgChargeLast;       // Default for fChargeLast   (now set to: 199.5 )
-  static const Int_t   fgAbsTimeNbins;     // Default for fAbsTimeNbins (now set to: 15    )
-  static const Axis_t  fgAbsTimeFirst;     // Default for fAbsTimeFirst (now set to: -0.5  )
-  static const Axis_t  fgAbsTimeLast;      // Default for fAbsTimeLast  (now set to: 14.5  )
+  static const Int_t   fgChargeNbins;      // Default for MHGausEvents::fNbins (now set to: 200   )
+  static const Axis_t  fgChargeFirst;      // Default for MHGausEvents::fFirst (now set to: -0.5  )
+  static const Axis_t  fgChargeLast;       // Default for MHGausEvents::fLast  (now set to: 199.5 )
+  static const Int_t   fgAbsTimeNbins;     // Default for fAbsTimeNbins        (now set to: 15    )
+  static const Axis_t  fgAbsTimeFirst;     // Default for fAbsTimeFirst        (now set to: -0.5  )
+  static const Axis_t  fgAbsTimeLast;      // Default for fAbsTimeLast         (now set to: 14.5  )
 
 public:
@@ -24 +24 @@
 
   // Setters
-  void SetChargeNbins(const Int_t  bins =fgChargeNbins)          { fChargeNbins = bins;     }
-  void SetChargeFirst(const Axis_t first=fgChargeFirst)          { fChargeFirst = first;    }
-  void SetChargeLast (const Axis_t last =fgChargeLast)           { fChargeLast  = last;     }
-  
   void SetAbsTimeNbins(const Int_t  bins =fgAbsTimeNbins)        { fAbsTimeNbins = bins;    }
   void SetAbsTimeFirst(const Axis_t first=fgAbsTimeFirst)        { fAbsTimeFirst = first;   }

trunk/MagicSoft/Mars/mcalib/MHCalibrationChargePix.cc

@@ -56 +56 @@
 const Axis_t  MHCalibrationChargePix::fgAbsTimeFirst    = -0.5;
 const Axis_t  MHCalibrationChargePix::fgAbsTimeLast     = 14.5;
-const Float_t MHCalibrationChargePix::fgPickupLimit     = 5.;
 // --------------------------------------------------------------------------
 //
@@ -62 +61 @@
 //
 // Sets: 
-// - the default number for fChargeNbins  (fgChargeNbins)
-// - the default number for fChargeFirst  (fgChargeFirst)
-// - the default number for fChargeLast   (fgChargeLast)
+// - the default number for fNbins        (fgChargeNbins)
+// - the default number for fFirst        (fgChargeFirst)
+// - the default number for fLast         (fgChargeLast)
 // - the default number for fAbsTimeNbins (fgAbstTimeNbins)
 // - the default number for fAbsTimeFirst (fgAbsTimeFirst)
 // - the default number for fAbsTimeLast  (fgAbsTimeLast)
-// - the default number for fPickupLimit  (fgPickupLimit)
 //
 // - the default name of the  fHGausHist ("HCalibrationCharge")
@@ -88 +86 @@
 //
 MHCalibrationChargePix::MHCalibrationChargePix(const char *name, const char *title)
-    : fHAbsTime(), fFlags(0)
+    : fHAbsTime()
 { 
   
@@ -94 +92 @@
   fTitle = title ? title : "Statistics of the FADC sums of calibration events";
 
-  SetChargeNbins();
-  SetChargeFirst();
-  SetChargeLast();
+  SetNbins ( fgChargeNbins );
+  SetFirst ( fgChargeFirst );
+  SetLast  ( fgChargeLast  );
 
   SetAbsTimeNbins();
   SetAbsTimeFirst();
   SetAbsTimeLast();
-
-  SetPickupLimit();
 
   fHGausHist.SetName("HCalibrationCharge");
@@ -124 +120 @@
 //
 // Sets:
-// - fHGausHist.SetBins(fChargeNbins,fChargeFirst,fChargeLast);
+// - fHGausHist.SetBins(fNbins,fFirst,fLast);
 // - fHAbsTime.SetBins(fAbsTimeNbins,fAbsTimeFirst,fAbsTimeLast);
 //
-void MHCalibrationChargePix::Init()
-{
-
-  fHGausHist.SetBins(fChargeNbins,fChargeFirst,fChargeLast);
+void MHCalibrationChargePix::InitBins()
+{
+  MHGausEvents::InitBins();
   fHAbsTime.SetBins(fAbsTimeNbins,fAbsTimeFirst,fAbsTimeLast);
 }
@@ -137 +132 @@
 //
 // Sets:
-// - all variables to 0., except fPixId to -1
-// - all flags to kFALSE
+// - all variables to 0.
 // 
 // Executes:
@@ -146 +140 @@
 {
 
-  fPixId     = -1;
   fSaturated = 0.;
-  fPickup    = 0.;
-
-  SetExcluded ( kFALSE );
-  
+
   MHGausEvents::Clear();
   return;
@@ -166 +156 @@
 // --------------------------------------------------------------------------
 //
-// Set fPixId to id
+// Calls MHGausEvents::ChangeHistId()
 //
 // Add id to names and titles of: 
-// - this
-// - fHGausHist
 // - fHAbsTime
 //
@@ -176 +164 @@
 {
 
-  fPixId = id;
-
-  fHGausHist.SetName(Form("%s%d", fHGausHist.GetName(), id));
-  fHGausHist.SetTitle(Form("%s%d", fHGausHist.GetTitle(), id));
-
-  fHAbsTime.SetName(Form("%s%d", fHAbsTime.GetName(), id));
+  MHGausEvents::ChangeHistId(id);
+
+  fHAbsTime.SetName (Form("%s%d", fHAbsTime.GetName(),  id));
   fHAbsTime.SetTitle(Form("%s%d", fHAbsTime.GetTitle(), id));
 
-  fName  = Form("%s%d", fName.Data(), id);
-  fTitle = Form("%s%d", fTitle.Data(), id);
-}
-
-// --------------------------------------------------------------------------
-//
-// Set Excluded bit from outside
-//
-void MHCalibrationChargePix::SetExcluded(const Bool_t b)
-{
-    b ? SETBIT(fFlags,kExcluded) : CLRBIT(fFlags,kExcluded);
-}
-
-const Bool_t MHCalibrationChargePix::IsExcluded()      const
-{
-  return TESTBIT(fFlags,kExcluded);
 }
 
@@ -317 +286 @@
 
 }
-
-// -----------------------------------------------------------------------------
-// 
-// Repeats the Gauss fit in a smaller range, defined by: 
-// 
-// min = GetMean() - fPickupLimit * GetSigma();
-// max = GetMean() + fPickupLimit * GetSigma();
-//
-Bool_t MHCalibrationChargePix::RepeatFit(const Option_t *option)
-{
-
-  //
-  // Get new fitting ranges
-  //
-  Axis_t rmin = GetMean() - fPickupLimit * GetSigma();
-  Axis_t rmax = GetMean() + fPickupLimit * GetSigma();
-
-  GetFGausFit()->SetRange(rmin,rmax);
-
-  GetHGausHist()->Fit(GetFGausFit(),option);
-
-  SetMean     ( GetFGausFit()->GetParameter(1) );
-  SetMeanErr  ( GetFGausFit()->GetParameter(2) );
-  SetSigma    ( GetFGausFit()->GetParError(1)  ); 
-  SetSigmaErr ( GetFGausFit()->GetParError(2)  ); 
-  SetProb     ( GetFGausFit()->GetProb()       );      
-
-  //
-  // The fit result is accepted under condition:
-  // 1) The results are not nan's
-  // 2) The NDF is not smaller than fNDFLimit (5)
-  // 3) The Probability is greater than fProbLimit (default 0.001 == 99.9%)
-  //
-  if (   TMath::IsNaN ( GetMean()     ) 
-      || TMath::IsNaN ( GetMeanErr()  )
-      || TMath::IsNaN ( GetProb()     )    
-      || TMath::IsNaN ( GetSigma()    )
-      || TMath::IsNaN ( GetSigmaErr() ) 
-      || GetFGausFit()->GetNDF() < fNDFLimit 
-      || GetProb() < fProbLimit )
-    return kFALSE;
-  
-  SetGausFitOK(kTRUE);
-  return kTRUE;
-
-}
-
-
-
-void MHCalibrationChargePix::CountPickup()
-{
-   fPickup  = (Int_t)GetHGausHist()->Integral(GetHGausHist()->GetXaxis()->FindBin(GetMean()+fPickupLimit*GetSigma()),
-                                               GetHGausHist()->GetXaxis()->GetLast(), 
-                                               "width");
-}
-
-
-
-
-
-

trunk/MagicSoft/Mars/mcalib/MHCalibrationChargePix.h

@@ -12 +12 @@
 private:
 
-  static const Int_t   fgChargeNbins;        // Default for fChargeNBins  (now set to: 2000  )
-  static const Axis_t  fgChargeFirst;        // Default for fChargeFirst  (now set to: -0.5  )
-  static const Axis_t  fgChargeLast;         // Default for fChargeLast   (now set to: 1999.5)
-  static const Int_t   fgAbsTimeNbins;       // Default for fAbsTimeNbins (now set to: 15    )
-  static const Axis_t  fgAbsTimeFirst;       // Default for fAbsTimeFirst (now set to: -0.5  )
-  static const Axis_t  fgAbsTimeLast;        // Default for fAbsTimeLast  (now set to: 14.5  )
-  static const Float_t fgPickupLimit;        // Default for fPickupLimit  (now set to:  5.   )
+  static const Int_t   fgChargeNbins;        // Default for MHGausEvents::fNBins (now set to: 2000  )
+  static const Axis_t  fgChargeFirst;        // Default for MHGausEvents::fFirst (now set to: -0.5  )
+  static const Axis_t  fgChargeLast;         // Default for MHGausEvents::fLast  (now set to: 1999.5)
+  static const Int_t   fgAbsTimeNbins;       // Default for fAbsTimeNbins        (now set to: 15    )
+  static const Axis_t  fgAbsTimeFirst;       // Default for fAbsTimeFirst        (now set to: -0.5  )
+  static const Axis_t  fgAbsTimeLast;        // Default for fAbsTimeLast         (now set to: 14.5  )
 
 protected:
 
-  Int_t    fPixId;         // The pixel ID
-
   TH1F     fHAbsTime;      // Histogram containing the absolute arrival times 
                           
-  Int_t    fChargeNbins;   // Number of  bins used for the fHGausHist
-  Axis_t   fChargeFirst;   // Lower bound bin used for the fHGausHist
-  Axis_t   fChargeLast;    // Upper bound bin used for the fHGausHist
-
   Int_t    fAbsTimeNbins;  // Number of  bins used for the fHAbsTime 
   Axis_t   fAbsTimeFirst;  // Lower bound bin used for the fHAbsTime
   Axis_t   fAbsTimeLast;   // Upper bound bin used for the fHAbsTime
 
-  Float_t  fPickupLimit;   // Upper number of sigmas from the fitted mean above which events are considered as pickup
-
   Float_t  fSaturated;     // Number of events classified as saturated
-  Float_t  fPickup;        // Number of events classified as pick-up
-
-  Byte_t   fFlags;         // Bit-field for the flags
-  enum     { kExcluded };  // Possible bits to be set
 
 public:
@@ -49 +36 @@
   virtual void Clear(Option_t *o="");
   virtual void Reset();  
-  virtual void Init();
-  virtual void ChangeHistId(Int_t i);
+  virtual void InitBins();
   
   // Setters 
-  virtual void SetChargeNbins(const Int_t  bins =fgChargeNbins)    { fChargeNbins = bins;   }
-  virtual void SetChargeFirst(const Axis_t first=fgChargeFirst)    { fChargeFirst = first;  }
-  virtual void SetChargeLast( const Axis_t last =fgChargeLast)     { fChargeLast  = last;   }
-  
   virtual void SetAbsTimeNbins(const Int_t  bins =fgAbsTimeNbins)  { fAbsTimeNbins = bins;  }
   virtual void SetAbsTimeFirst(const Axis_t first=fgAbsTimeFirst)  { fAbsTimeFirst = first; }
   virtual void SetAbsTimeLast( const Axis_t last =fgAbsTimeLast)   { fAbsTimeLast  = last;  }
 
-  virtual void SetPickupLimit( const Float_t  lim =fgPickupLimit)  { fPickupLimit  = lim;   }
-
   void SetSaturated          ( const Float_t f    )                { fSaturated += f;      }
-  void SetExcluded           ( const Bool_t b=kTRUE );
 
   // Getters
@@ -70 +49 @@
   const TH1F *GetHAbsTime()             const { return &fHAbsTime;  }
 
-  const Float_t  GetIntegral()          const;
-  
   const Float_t  GetAbsTimeMean(  )     const;
   const Float_t  GetAbsTimeRms()        const;
-
+  const Float_t  GetIntegral()          const;
   const Float_t  GetSaturated()         const { return fSaturated;   }
-  const Float_t  GetPickup()            const { return fPickup;      }
-
-  const Bool_t   IsExcluded()          const;
 
   // Fill histos
   Bool_t FillAbsTime(const Float_t t);
 
-  // Fits
-  Bool_t RepeatFit(const Option_t *option="RQ0");
-  
   // Draws
   virtual void Draw(Option_t *opt="");
 
   // Miscelleaneous
-  void CountPickup();
-
+  void ChangeHistId(Int_t id);
+  
   ClassDef(MHCalibrationChargePix, 1)     // Base Histogram class for a Charge Calibration Pixel
 };

trunk/MagicSoft/Mars/mcalib/MHCalibrationRelTimeCam.cc

@@ -22 +22 @@
 !
 \* ======================================================================== */
-
 /////////////////////////////////////////////////////////////////////////////
 //                                                                         //
-// MHCalibrationRelTimeCam                                                            //
+// MHCalibrationRelTimeCam                                                 //
 //                                                                         //
-// Hold the CalibrationRelTime information for all pixels in the camera              //
-//                                                                         //
+// Fills the extracted relative arrival times of MArrivalTimeCam into 
+// the MHGausEvents-classes MHCalibrationRelTimePix for every:
+//
+// - pixel, stored in the TObjArray's MHCalibrationCam::fHiGainArray 
+//
+// - average pixel per area index (e.g. inner and outer for the MAGIC camera), 
+//   stored in the TObjArray's MHCalibrationCam::fAverageHiGainAreas
+//
+// - average pixel per camera sector (e.g. sectors 1-6 for the MAGIC camera), 
+//   stored in the TObjArray's MHCalibrationCam::fAverageHiGainSectors 
+// 
+// Every time is filled into a histogram and an array, in order to perform 
+// a Fourier analysis (see MHGausEvents). The signals are moreover averaged on an 
+// event-by-event basis and written into the corresponding average pixels.
+//
+// The histograms are fitted to a Gaussian, mean and sigma with its errors 
+// and the fit probability are extracted. If none of these values are NaN's and 
+// if the probability is bigger than fProbLimit (default: 0.5%), the fit is valid.
+// Otherwise, the fit is repeated within ranges of the previous mean +- 5 sigma. 
+// In case this does not make the fit valid, the histogram means and RMS's are 
+// taken directly and the following flags are set:
+// - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::kRelTimeNotFitted ) and
+// - MBadPixelsPix::SetUnsuitable(   MBadPixelsPix::kUnreliableRun    ) 
+// 
+// Outliers of more than MHCalibrationChargePix::fPickupLimit (default: 5) sigmas 
+// from the mean are counted as PickUp events (stored in MHCalibrationRelTimePix::fPickup) 
+//
+// The class also fills arrays with the signal vs. event number, creates a fourier 
+// spectrum and investigates if the projected fourier components follow an exponential 
+// distribution. In case that the probability of the exponential fit is less than 
+// fProbLimit (default: 0.5%), the following flags are set:
+// - MBadPixelsPix::SetUncalibrated( MBadPixelsPix::kRelTimeOscillating ) and
+// - MBadPixelsPix::SetUnsuitable(   MBadPixelsPix::kUnreliableRun      )
+// 
+// This same procedure is performed for the average pixels.
+//
+// The following results are written into MCalibrationRelTimeCam:
+//
+// - MCalibrationPix::SetMean()
+// - MCalibrationPix::SetMeanErr()
+// - MCalibrationPix::SetSigma()
+// - MCalibrationPix::SetSigmaErr()
+// - MCalibrationPix::SetProb()
+// - MCalibrationPix::SetNumPickup()
+//
+// For all averaged areas, the fitted sigma is multiplied with the square root of 
+// the number involved pixels in order to be able to compare it to the average of 
+// sigmas in the camera.
+//                                                                         
 /////////////////////////////////////////////////////////////////////////////
 #include "MHCalibrationRelTimeCam.h"
+#include "MHCalibrationRelTimePix.h"
 
 #include "MLog.h"
@@ -37 +84 @@
 #include "MParList.h"
 
-#include "MHCalibrationRelTimePix.h"
+#include "MCalibrationRelTimeCam.h"
+#include "MCalibrationPix.h"
 
 #include "MArrivalTimeCam.h"
 #include "MArrivalTimePix.h"
 
+#include "MGeomCam.h"
+#include "MGeomPix.h"
+
+#include "MBadPixelsCam.h"
+#include "MBadPixelsPix.h"
+
 ClassImp(MHCalibrationRelTimeCam);
 
 using namespace std;
-const Float_t MHCalibrationRelTimeCam::fgTimeSliceWidth  = 3.3;
-const Int_t   MHCalibrationRelTimeCam::fgPulserFrequency = 500;
 // --------------------------------------------------------------------------
 //
-// Default constructor. Creates a MHCalibrationRelTimePix object for each pixel
+// Default Constructor. 
 //
 MHCalibrationRelTimeCam::MHCalibrationRelTimeCam(const char *name, const char *title) 
@@ -55 +107 @@
 
   fName  = name  ? name  : "MHCalibrationRelTimeCam";
-  fTitle = title ? title : "Histogram container for the relative time calibration of the camera";
-  
-  //
-  // loop over all Pixels and create two histograms
-  // one for the Low and one for the High gain
-  // connect all the histogram with the container fHist
-  //
-  fArray = new TObjArray;
-  fArray->SetOwner();
-  
-
-  SetTimeSliceWidth();
-  SetPulserFrequency();
+  fTitle = title ? title : "Histogram class for the relative time calibration of the camera";
+  
 }
 
 // --------------------------------------------------------------------------
 //
-// Delete the array conatining the pixel pedest information
-//
-MHCalibrationRelTimeCam::~MHCalibrationRelTimeCam()
-{
-  delete fArray;
-}
-
-// --------------------------------------
-//
-void MHCalibrationRelTimeCam::Clear(Option_t *o)
-{
-
-  fArray->ForEach(TObject, Clear)();
-  return;
-}
-
-// --------------------------------------------------------------------------
-//
-// Get i-th pixel (pixel number)
-//
-MHCalibrationRelTimePix &MHCalibrationRelTimeCam::operator[](UInt_t i)
-{
-  return *(MHCalibrationRelTimePix*)(fArray->At(i));
-}
-
-// --------------------------------------------------------------------------
-//
-// Get i-th pixel (pixel number)
-//
-const MHCalibrationRelTimePix &MHCalibrationRelTimeCam::operator[](UInt_t i) const
-{
-  return *(MHCalibrationRelTimePix*)(fArray->At(i));
-}
-
-
-// --------------------------------------------------------------------------
-//
-// Our own clone function is necessary since root 3.01/06 or Mars 0.4
-// I don't know the reason
-//
-TObject *MHCalibrationRelTimeCam::Clone(const char *) const
-{
-
-  const Int_t n = fArray->GetSize();
-  
-  //
-  // FIXME, this might be done faster and more elegant, by direct copy.
-  //
-  MHCalibrationRelTimeCam *cam = new MHCalibrationRelTimeCam;
-  
-  cam->fArray->Expand(n);
-  
-  for (int i=0; i<n; i++)
-    {
-      delete (*cam->fArray)[i];
-      (*cam->fArray)[i] = (*fArray)[i]->Clone();
-    }
-
-  return cam;
-}
-
-
-
-// --------------------------------------------------------------------------
-//
-//
-Bool_t MHCalibrationRelTimeCam::SetupFill(const MParList *pList)
-{
+// Gets or creates the pointers to:
+// - MCalibrationRelTimeCam
+//
+// Searches pointer to:
+// - MArrivalTimeCam
+//
+// Initializes, if empty to MArrivalTimeCam::GetSize() for:
+// - MHCalibrationCam::fHiGainArray, MHCalibrationCam::fLoGainArray
+//
+// Initializes, if empty to MGeomCam::GetNumAreas() for:
+// - MHCalibrationCam::fAverageHiGainAreas, MHCalibrationCam::fAverageLoGainAreas
+//
+// Initializes, if empty to MGeomCam::GetNumSectors() for:
+// - MHCalibrationCam::fAverageHiGainSectors, MHCalibrationCam::fAverageLoGainSectors
+// 
+// Calls InitializeHists() for every entry in:
+// - MHCalibrationCam::fHiGainArray
+// - MHCalibrationCam::fAverageHiGainAreas
+// - MHCalibrationCam::fAverageHiGainSectors
+//
+// Sets Titles and Names for the Histograms 
+// - MHCalibrationCam::fAverageHiGainAreas
+// - MHCalibrationCam::fAverageHiGainSectors
+// 
+// Sets number of bins to MHCalibrationCam::fAverageNbins for:
+// - MHCalibrationCam::fAverageHiGainAreas, MHCalibrationCam::fAverageLoGainAreas
+// - MHCalibrationCam::fAverageHiGainSectors, MHCalibrationCam::fAverageLoGainSectors
+//
+Bool_t MHCalibrationRelTimeCam::ReInitHists(MParList *pList)
+{
+
+  fCam = (MCalibrationCam*)pList->FindCreateObj("MCalibrationRelTimeCam");
+  if (!fCam)
+      return kFALSE;
+
+  MArrivalTimeCam *signal = (MArrivalTimeCam*)pList->FindObject("MArrivalTimeCam");
+  if (!signal)
+  {
+      *fLog << err << "MArrivalTimeCam not found... abort." << endl;
+      return kFALSE;
+  }
+
+  const Int_t npixels  = fGeom->GetNumPixels();
+  const Int_t nsectors = fGeom->GetNumSectors();
+  const Int_t nareas   = fGeom->GetNumAreas();
+
+  if (fHiGainArray->GetEntries()==0)
+  {
+      fHiGainArray->Expand(npixels);
+      for (Int_t i=0; i<npixels; i++)
+      {
+          (*fHiGainArray)[i] = new MHCalibrationRelTimePix;
+          InitHists((*this)[i],(*fBadPixels)[i],i);
+      }
+  }
+
+  if (fLoGainArray->GetEntries()==0)
+  {
+      fLoGainArray->Expand(npixels);
+      for (Int_t i=0; i<npixels; i++)
+      {
+          (*fLoGainArray)[i] = new MHCalibrationRelTimePix;
+          InitHists((*this)(i),(*fBadPixels)[i],i);
+      }
+  }
+
+  if (fAverageHiGainAreas->GetEntries()==0)
+  {
+    fAverageHiGainAreas->Expand(nareas);
+    
+    for (Int_t j=0; j<nareas; j++)
+      {
+        (*fAverageHiGainAreas)[j] = 
+          new MHCalibrationRelTimePix("AverageHiGainArea",
+                                      "Average Relative Arrival Times area idx ");
+
+        InitHists(GetAverageHiGainArea(j),fCam->GetAverageBadArea(j),j);
+
+        GetAverageHiGainArea(j).GetHGausHist()->SetTitle("Relative Arrival Times average Area Idx ");
+        GetAverageHiGainArea(j).SetNbins(fAverageNbins);
+      }
+  }
+
+  if (fAverageLoGainAreas->GetEntries()==0)
+  {
+    fAverageLoGainAreas->Expand(nareas);
+    
+    for (Int_t j=0; j<nareas; j++)
+      {
+        (*fAverageLoGainAreas)[j] = 
+          new MHCalibrationRelTimePix("AverageLoGainArea",
+                                      "Relative Arrival Times average Area idx ");
+
+        InitHists(GetAverageLoGainArea(j),fCam->GetAverageBadArea(j),j);
+
+        GetAverageLoGainArea(j).GetHGausHist()->SetTitle("Relative Arrival Times average Area Idx ");
+        GetAverageLoGainArea(j).SetNbins(fAverageNbins);
+      }
+  }
+
+  if (fAverageHiGainSectors->GetEntries()==0)
+  {
+      fAverageHiGainSectors->Expand(nsectors);
+
+      for (Int_t j=0; j<nsectors; j++)
+      {
+          (*fAverageHiGainSectors)[j] = 
+            new MHCalibrationRelTimePix("AverageHiGainSector",
+                                        "Relative Arrival Times average sector ");
+
+          GetAverageHiGainSector(j).GetHGausHist()->SetTitle("Relative Arrival Times average Sector ");
+          GetAverageHiGainSector(j).SetNbins(fAverageNbins);
+      }
+  }
+
+  if (fAverageLoGainSectors->GetEntries()==0)
+  {
+      fAverageLoGainSectors->Expand(nsectors);
+
+      for (Int_t j=0; j<nsectors; j++)
+      {
+          (*fAverageLoGainSectors)[j] = 
+            new MHCalibrationRelTimePix("AverageLoGainSector",
+                                        "Relative Arrival Times average sector ");
+
+          InitHists(GetAverageLoGainSector(j),fCam->GetAverageBadSector(j),j);
+          
+          GetAverageLoGainSector(j).GetHGausHist()->SetTitle("Relative Arrival Times average Sector ");
+          GetAverageLoGainSector(j).SetNbins(fAverageNbins);
+      }
+  }
 
   return kTRUE;
-
-}
-
-// --------------------------------------------------------------------------
-//
-Bool_t MHCalibrationRelTimeCam::Fill(const MParContainer *par, const Stat_t w)
+}
+
+
+//--------------------------------------------------------------------------------
+//
+// Retrieves pointer to MArrivalTimeCam:
+//
+// Retrieves from MGeomCam:
+// - number of pixels
+// - number of pixel areas
+// - number of sectors
+//
+// Fill RelTimes histograms (MHGausEvents::FillHistAndArray()) with:
+// - MArrivalTimePix::GetArrivalTime(pixid) - MArrivalTimePix::GetArrivalTime(1);
+//   (i.e. the time difference between pixel i and pixel 1 (hardware number: 2)
+//
+Bool_t MHCalibrationRelTimeCam::FillHists(const MParContainer *par, const Stat_t w)
 {
 
@@ -157 +272 @@
     }
   
-  const Int_t n = arrtime->GetSize();
-  
-  if (fArray->GetEntries()==0)
-    {
-      fArray->Expand(n);
+  const Int_t npixels  = fGeom->GetNumPixels();
+  const Int_t nareas   = fGeom->GetNumAreas();
+  const Int_t nsectors = fGeom->GetNumSectors();
+
+  Float_t sumareahi  [nareas],   sumarealo  [nareas];
+  Float_t sumsectorhi[nsectors], sumsectorlo[nsectors];
+  Int_t   numareahi  [nareas],   numarealo  [nareas];
+  Int_t   numsectorhi[nsectors], numsectorlo[nsectors];
+
+  for (UInt_t j=0; j<nareas; j++)
+    {
+      sumareahi[j] = sumarealo[j] = 0.;
+      numareahi[j] = numarealo[j] = 0;
+    }
+  
+  for (UInt_t j=0; j<nsectors; j++)
+    {
+      sumsectorhi[j] = sumsectorlo[j] = 0.;
+      numsectorhi[j] = numsectorlo[j] = 0;
+    }
+  
+  const MArrivalTimePix &refpix = (*arrtime)[1];
+  const Float_t reftime = refpix.IsLoGainUsed() 
+    ? refpix.GetArrivalTimeLoGain() : refpix.GetArrivalTimeHiGain();
+
+  for (Int_t i=0; i<npixels; i++)
+    {
+
+      MHGausEvents &histhi = (*this)[i];
+      MHGausEvents &histlo = (*this)(i);
+
+      if (histhi.IsExcluded())
+        continue;
+
+      const MArrivalTimePix &pix = (*arrtime)[i];
+      const Float_t time    = pix.IsLoGainUsed() 
+        ? pix.GetArrivalTimeLoGain() : pix.GetArrivalTimeHiGain();  
+      const Float_t reltime = time - reftime;
       
-      for (Int_t i=0; i<n; i++)
+      const Int_t aidx   = (*fGeom)[i].GetAidx();
+      const Int_t sector = (*fGeom)[i].GetSector();
+
+      if (pix.IsLoGainUsed())
+        { 
+          histlo.FillHistAndArray(reltime);
+          sumarealo  [aidx]   += reltime;
+          numarealo  [aidx]   ++;
+          sumsectorlo[sector] += reltime;
+          numsectorlo[sector] ++;
+        }
+      else
         {
-          (*fArray)[i] = new MHCalibrationRelTimePix;
-          MHCalibrationRelTimePix &hist = (*this)[i];
-          hist.ChangeHistId(i);
-          hist.Init();
-          hist.SetEventFrequency(fPulserFrequency);
+          histhi.FillHistAndArray(reltime) ;
+          sumareahi  [aidx]   += reltime;
+          numareahi  [aidx]   ++;
+          sumsectorhi[sector] += reltime;
+          numsectorhi[sector] ++;
         }
     }
   
-  if (fArray->GetEntries() != n)
-    {
-      gLog << err << "ERROR - Size mismatch... abort." << endl;
-      return kFALSE;
-    }
-  
-  const MArrivalTimePix &refpix = (*arrtime)[1];
-  const Float_t reftime = refpix.IsLoGainUsed() ? refpix.GetArrivalTimeLoGain() : refpix.GetArrivalTimeHiGain();
-
-  for (Int_t i=0; i<n; i++)
-    {
-        const MArrivalTimePix &pix = (*arrtime)[i];
-        const Float_t time = pix.IsLoGainUsed() ? pix.GetArrivalTimeLoGain() : pix.GetArrivalTimeHiGain();  
-
-        MHCalibrationRelTimePix  &hist = (*this)[i];
-        hist.FillHistAndArray(time - reftime);
-    }
-  
+  for (UInt_t j=0; j<nareas; j++)
+    {
+      MHGausEvents &histhi = GetAverageHiGainArea(j);
+      histhi.FillHistAndArray(numareahi[j] == 0 ? 0. : sumareahi[j]/numareahi[j]);
+
+      MHGausEvents &histlo = GetAverageLoGainArea(j);
+      histlo.FillHistAndArray(numarealo[j] == 0 ? 0. : sumarealo[j]/numarealo[j]);
+    }
+  
+  for (UInt_t j=0; j<nsectors; j++)
+    {
+      MHGausEvents &histhi = GetAverageHiGainSector(j);
+      histhi.FillHistAndArray(numsectorhi[j] == 0 ? 0. : sumsectorhi[j]/numsectorhi[j]);
+
+      MHGausEvents &histlo = GetAverageLoGainSector(j);
+      histlo.FillHistAndArray(numsectorlo[j] == 0 ? 0. : sumsectorlo[j]/numsectorlo[j]);
+    }
+
   return kTRUE;
 }
 
-Bool_t MHCalibrationRelTimeCam::Finalize()
-{
-
-  for (Int_t i=0; i<fArray->GetSize(); i++)
+// --------------------------------------------------------------------------
+//
+// Calls:
+// - MHCalibrationCam::FitHiGainArrays() with flags:
+//   MBadPixelsPix::kRelTimeNotFitted and MBadPixelsPix::kRelTimeOscillating
+// - MHCalibrationCam::FitLoGainArrays() with flags:
+//   MBadPixelsPix::kRelTimeNotFitted and MBadPixelsPix::kRelTimeOscillating
+// 
+Bool_t MHCalibrationRelTimeCam::FinalizeHists()
+{
+
+  FitHiGainArrays((MCalibrationCam&)(*fCam),*fBadPixels,
+                  MBadPixelsPix::kRelTimeNotFitted,
+                  MBadPixelsPix::kRelTimeOscillating);
+
+  FitLoGainArrays((MCalibrationCam&)(*fCam),*fBadPixels,
+                  MBadPixelsPix::kRelTimeNotFitted,
+                  MBadPixelsPix::kRelTimeOscillating);
+
+  return kTRUE;
+}
+
+// --------------------------------------------------------------------------
+//
+// Sets all pixels to MBadPixelsPix::kUnreliableRun, if following flags are set:
+// - MBadPixelsPix::kRelTimeNotFitted
+// - MBadPixelsPix::kRelTimeOscillating
+//
+void MHCalibrationRelTimeCam::FinalizeBadPixels()
+{
+
+  for (Int_t i=0; i<fBadPixels->GetSize(); i++)
     {
       
-      MHCalibrationRelTimePix &hist = (*this)[i];
+      MBadPixelsPix          &bad    = (*fBadPixels)[i];
+
+      if (bad.IsUncalibrated( MBadPixelsPix::kRelTimeNotFitted ))
+        bad.SetUnsuitable(   MBadPixelsPix::kUnreliableRun    );
+
+      if (bad.IsUncalibrated( MBadPixelsPix::kRelTimeOscillating))
+        bad.SetUnsuitable(   MBadPixelsPix::kUnreliableRun    );
       
-      //
-      // 1) Return if the charge distribution is already succesfully fitted
-      //    or if the histogram is empty
-      //
-      if (hist.IsGausFitOK() || hist.IsEmpty())
-        continue;
-      
-      //
-      // 2) Fit the Hi Gain histograms with a Gaussian
-      //
-      hist.FitGaus();
-      
-      //
-      // 3) If fit does not succeed , bypass the fit and take the histogram means and sigmas
-      //
-      if (!hist.IsGausFitOK())
-        hist.BypassFit();
-      
-      //
-      // 4) Create the fourier transform of the arrays 
-      //
-      hist.CreateFourierSpectrum();
-      
-      //
-      // 5) Renormalize to the real time in ns.
-      //
-      hist.Renorm(fTimeSliceWidth);
-      
-    }
-  return kTRUE;
+    }
 }
 
@@ -250 +414 @@
 // =============================================
 //
-// 4: Returned probability of Gauss fit to Charge distribution
+// 4: Returned probability of Gauss fit to RelTime distribution
 //
 // Localized defects:
@@ -261 +425 @@
 {
 
-  if (fArray->GetSize() <= idx)
+  if (fHiGainArray->GetSize() <= idx)
     return kFALSE;
 
+  MHCalibrationRelTimePix &hist = (MHCalibrationRelTimePix&)(*this)[idx];      
+
   switch (type)
     {
     case 0:
-      val = (*this)[idx].GetMean();
+      val = hist.GetMean();
       break;
     case 1:
-      val = (*this)[idx].GetMeanErr();
+      val = hist.GetMeanErr();
       break;
     case 2:
-      val = (*this)[idx].GetSigma();
+      val = hist.GetSigma();
       break;
     case 3:
-      val = (*this)[idx].GetSigmaErr();
+      val = hist.GetSigmaErr();
       break;
     case 4:
-      val = (*this)[idx].GetProb();
+      val = hist.GetProb();
       break;
     case 5:
-      if (!(*this)[idx].IsGausFitOK())
+      if (!hist.IsGausFitOK())
         val = 1.;
       break;
     case 6:
-      if (!(*this)[idx].IsFourierSpectrumOK())
+      if (!hist.IsFourierSpectrumOK())
         val = 1.;
       break;
@@ -297 +463 @@
 void MHCalibrationRelTimeCam::DrawPixelContent(Int_t idx) const
 {
-  (*this)[idx].DrawClone();
-}
+  const MHCalibrationRelTimePix &hist = (MHCalibrationRelTimePix&)(*this)[idx];      
+  hist.DrawClone();
+}

trunk/MagicSoft/Mars/mcalib/MHCalibrationRelTimeCam.h

@@ -2 +2 @@
 #define MARS_MHCalibrationRelTimeCam
 
-#ifndef ROOT_TObjArray
-#include <TObjArray.h>
+#ifndef MARS_MHCalibrationCam
+#include "MHCalibrationCam.h"
 #endif
 
-#ifndef MARS_MH
-#include "MH.h"
-#endif
-#ifndef MARS_MCamEvent
-#include "MCamEvent.h"
-#endif
-
-class MHCalibrationRelTimePix;
-class MHCalibrationRelTimeCam : public MH, public MCamEvent
+class MGeomCam;
+class MHCalibrationRelTimeCam : public MHCalibrationCam
 {
 
 private:
 
-  static const Float_t fgTimeSliceWidth;       // Default for fTimeSliceWidth
-  static const Int_t   fgPulserFrequency;      // Default for fPulserFrequency
-
-  Float_t fTimeSliceWidth;                    // FADC slice time width
-  Int_t   fPulserFrequency;                   // The pulser frequency
+  Bool_t ReInitHists(MParList *pList);
+  Bool_t FillHists(const MParContainer *par, const Stat_t w=1);
+  Bool_t FinalizeHists();
+  void   FinalizeBadPixels();
   
-  TObjArray  *fArray;       //-> List of MHCalibrationRelTimePix's
-
 public:
 
   MHCalibrationRelTimeCam(const char *name=NULL, const char *title=NULL);
-  ~MHCalibrationRelTimeCam();
+  ~MHCalibrationRelTimeCam() {}
 
-  void Clear(Option_t *o="");
-  
-  MHCalibrationRelTimePix &operator[](UInt_t i);
-  const MHCalibrationRelTimePix &operator[](UInt_t i) const;
-  
-  Bool_t SetupFill(const MParList *pList);
-  Bool_t Fill(const MParContainer *par, const Stat_t w=1);
-  Bool_t Finalize();
-
-  // Setters
-  void SetTimeSliceWidth( const Float_t width=fgTimeSliceWidth)  {  fTimeSliceWidth = width; }
-  void SetPulserFrequency( const Int_t  f=fgPulserFrequency)     { fPulserFrequency = f;     }
-  
-  TObject *Clone(const char *) const;
-  
   Bool_t GetPixelContent(Double_t &val, Int_t idx, const MGeomCam &cam, Int_t type=0) const;
   void DrawPixelContent(Int_t idx) const;

trunk/MagicSoft/Mars/mcalib/MHCalibrationRelTimePix.cc

@@ -40 +40 @@
 using namespace std;
 //
-const Int_t   MHCalibrationRelTimePix::fgChargeNbins    = 900;
-const Axis_t  MHCalibrationRelTimePix::fgChargeFirst    = -13.;
-const Axis_t  MHCalibrationRelTimePix::fgChargeLast     =  13.;
+const Int_t   MHCalibrationRelTimePix::fgRelTimeNbins    = 900;
+const Axis_t  MHCalibrationRelTimePix::fgRelTimeFirst    = -13.;
+const Axis_t  MHCalibrationRelTimePix::fgRelTimeLast     =  13.;
 // --------------------------------------------------------------------------
 //
@@ -48 +48 @@
 //
 // Sets: 
-// - the default number for fChargeNbins  (fgChargeNbins)
-// - the default number for fChargeFirst  (fgChargeFirst)
-// - the default number for fChargeLast   (fgChargeLast)
+// - the default number for fNbins  (fgRelTimeNbins)
+// - the default number for fFirst  (fgRelTimeFirst)
+// - the default number for fLast   (fgRelTimeLast)
 //
 // - the default name of the  fHGausHist ("HCalibrationRelTime")
@@ -64 +64 @@
 //
 MHCalibrationRelTimePix::MHCalibrationRelTimePix(const char *name, const char *title) 
-    : fPixId(-1)
 { 
 
@@ -70 +69 @@
   fTitle = title ? title : "Histogrammed Calibration Relative Arrival Time events";
 
-  SetChargeNbins();
-  SetChargeFirst();
-  SetChargeLast();
+  SetNbins ( fgRelTimeNbins );
+  SetFirst ( fgRelTimeFirst );
+  SetLast  ( fgRelTimeLast  );
 
   // Create a large number of bins, later we will rebin
@@ -96 +95 @@
 
   fPixId = -1;
+
   MHGausEvents::Clear();
   return;
@@ -109 +109 @@
 }
 
-// --------------------------------------------------------------------------
-//
-// Sets:
-// - fHGausHist.SetBins(fChargeNbins,fChargeFirst,fChargeLast);
-//
-void MHCalibrationRelTimePix::Init()
-{
 
-  fHGausHist.SetBins(fChargeNbins,fChargeFirst,fChargeLast);
-
-}
-
-// --------------------------------------------------------------------------
-//
-// - Set fPixId to id
-//
-// Add id to names and titles of:
-// - fHGausHist
-//
-void MHCalibrationRelTimePix::ChangeHistId(Int_t id)
-{
-
-  fPixId = id;
-
-  fHGausHist.SetName(  Form("%s%d", fHGausHist.GetName(),  id));
-  fHGausHist.SetTitle( Form("%s%d", fHGausHist.GetTitle(), id));
-
-}
-
-
-// ----------------------------------------------------------------------
-//
-// Renorm the results from FADC slices to times. 
-// The parameters slicewidth is the inverse of the FADC frequency and has the unit ns.
-//
-void MHCalibrationRelTimePix::Renorm(const Float_t slicewidth)
-{
-
-  SetMean(     GetMean()    * slicewidth  );
-  SetMeanErr(  GetMeanErr() * slicewidth  );
-  SetSigma(    GetSigma()   * slicewidth  );
-  SetSigmaErr( GetSigmaErr()* slicewidth  );
-  
-}
-

trunk/MagicSoft/Mars/mcalib/MHCalibrationRelTimePix.h

@@ -11 +11 @@
 private:
 
-  static const Int_t   fgChargeNbins;      // Default for fChargeNBins  (now set to: 900   )
-  static const Axis_t  fgChargeFirst;      // Default for fChargeFirst  (now set to: -13.5 )
-  static const Axis_t  fgChargeLast;       // Default for fChargeLast   (now set to:  13.5 )
-
-  Int_t   fChargeNbins;                   // Number of  bins used for the fHGausHist
-  Axis_t  fChargeFirst;                   // Lower bound bin used for the fHGausHist
-  Axis_t  fChargeLast;                    // Upper bound bin used for the fHGausHist
-
-  Int_t fPixId;                           // The pixel ID
+  static const Int_t   fgRelTimeNbins;      //! Default for MHGausEvents::fNBins  (now set to: 900   )
+  static const Axis_t  fgRelTimeFirst;      //! Default for MHGausEvents::fFirst  (now set to: -13.5 )
+  static const Axis_t  fgRelTimeLast;       //! Default for MHGausEvents::fLast   (now set to:  13.5 )
 
 public:
@@ -26 +20 @@
   ~MHCalibrationRelTimePix() {}
   
-
   void Clear(Option_t *o="");
   void Reset();
-  void Init();
-  
-  // Setters
-  void SetChargeNbins(const Int_t  bins =fgChargeNbins)    { fChargeNbins = bins; }
-  void SetChargeFirst(const Axis_t first=fgChargeFirst)    { fChargeFirst = first; }
-  void SetChargeLast( const Axis_t last =fgChargeLast)     { fChargeLast  = last; }
-
-  // Others
-  void ChangeHistId(Int_t i);
-  void Renorm(const Float_t slicewidth);
   
   ClassDef(MHCalibrationRelTimePix, 1)     // Histogram class for Relative Time Calibration pixel

trunk/MagicSoft/Mars/mcalib/MHGausEvents.cc

@@ -40 +40 @@
 //     It is created with a default name and title and resides in directory NULL.
 //     - Any class deriving from MHGausEvents needs to apply a binning to fHGausHist
-//       (e.g. with the function TH1F::SetBins(..) )
+//       (e.g. by setting the variables fNbins, fFirst, fLast and calling the function 
+//       InitBins() or by directly calling fHGausHist.SetBins(..) )
 //     - The histogram is filled with the functions FillHist() or FillHistAndArray(). 
 //     - The histogram can be fitted with the function FitGaus(). This involves stripping 
@@ -48 +49 @@
 //       The NDF is compared to fNDFLimit and a check is made whether results are NaNs. 
 //       Accordingly, a flag IsGausFitOK() is set.
+//     - One can repeat the fit within a given amount of sigmas from the previous mean 
+//       (specified by the variable fPickupLimit) with the function RepeatFit()
+//     - One can completely skip the fitting to set mean, sigma and its errors directly 
+//       from the histograms with the function BypassFit()
 // 
 //  2) The TArrayF fEvents:
@@ -96 +101 @@
 const Float_t  MHGausEvents::fgProbLimit            = 0.001;
 const Int_t    MHGausEvents::fgNDFLimit             = 2;
+const Float_t  MHGausEvents::fgPickupLimit          = 5.;
 const Int_t    MHGausEvents::fgPowerProbabilityBins = 20;
 const Int_t    MHGausEvents::fgBinsAfterStripping   = 40;
@@ -105 +111 @@
 // - the default Probability Limit for fProbLimit  (fgProbLimit)
 // - the default NDF Limit for fNDFLimit  (fgNDFLimit)
+// - the default number for fPickupLimit  (fgPickupLimit)
 // - the default number of bins after stripping for fBinsAfterStipping (fgBinsAfterStipping)
 // - the default name of the fHGausHist   ("HGausHist")
 // - the default title of the fHGausHist  ("Histogram of Events with Gaussian Distribution")
 // - the default directory of the fHGausHist (NULL)
+// - the default for fNbins (100)
+// - the default for fFirst (0.)
+// - the default for fLast  (100.)
 //
 // Initializes:
@@ -121 +131 @@
       fPowerSpectrum(NULL),
       fGraphEvents(NULL), fGraphPowerSpectrum(NULL),
+      fNbins(100), fFirst(0.), fLast(100.), fPixId(-1),
       fHGausHist(), fEvents(0),
       fFGausFit(NULL), fFExpFit(NULL)
@@ -133 +144 @@
   SetProbLimit();
   SetNDFLimit();
+  SetPickupLimit();
   SetBinsAfterStripping();
 
@@ -183 +195 @@
 // --------------------------------------------------------------------------
 //
+// Default Clear(), can be overloaded.
+//
 // Sets:
 // - all other pointers to NULL
-// - all variables to 0., except fEventFrequency
+// - all variables to 0., except fPixId to -1 and keep fEventFrequency
 // - all flags to kFALSE
 // 
@@ -197 +211 @@
   SetExpFitOK         ( kFALSE );
   SetFourierSpectrumOK( kFALSE );
+  SetExcluded         ( kFALSE );
 
   fMean              = 0.;
@@ -205 +220 @@
   
   fCurrentSize       = 0;
+  fPixId             = -1;
 
   if (fHPowerProbability)
@@ -247 +263 @@
 
 // --------------------------------------------------------------------------
+//
+// Default Reset(), can be overloaded.
 //
 // Executes:
@@ -264 +282 @@
 // --------------------------------------------------------------------------
 //
+// Default InitBins, can be overloaded.
+//
+// Executes:
+// - fHGausHist.SetBins(fNbins,fFirst,fLast)
+//
+void MHGausEvents::InitBins()
+{
+  fHGausHist.SetBins(fNbins,fFirst,fLast);
+}
+
+// --------------------------------------------------------------------------
+//
 // Executes:
 // - FillArray()
@@ -303 +333 @@
 }
 
+// --------------------------------------------------------------------------
+//
+// Set Excluded bit from outside
+//
+void MHGausEvents::SetExcluded(const Bool_t b)
+{
+    b ? SETBIT(fFlags,kExcluded) : CLRBIT(fFlags,kExcluded);
+}
+
 
 const Double_t MHGausEvents::GetChiSquare()  const 
@@ -366 +405 @@
 {
   return TESTBIT(fFlags,kExpFitOK);
+}
+
+const Bool_t MHGausEvents::IsExcluded() const
+{
+  return TESTBIT(fFlags,kExcluded);
 }
 
@@ -409 +453 @@
   if (fFExpFit)
     return;
+
+  if (fEvents.GetSize() < 8)
+    {
+      *fLog << warn << "Cannot create Fourier spectrum in pixel: " << fPixId 
+            << ". Number of events smaller than 8 " << endl;
+      return;
+    }
+  
 
   //
@@ -424 +476 @@
   // This cuts only the non-used zero's, but MFFT will later cut the rest
   MArray::StripZeros(fEvents);
+
+  if (fEvents.GetSize() < 8)
+    {
+      *fLog << warn << "Cannot create Fourier spectrum. " << endl;
+      *fLog << warn << "Number of events (after stripping of zeros) is smaller than 8 " 
+            << "in pixel: " << fPixId << endl;
+      return;
+    }
 
   MFFT fourier;
@@ -506 +566 @@
   if (!fFGausFit) 
     {
-    *fLog << warn << dbginf << "WARNING: Could not create fit function for Gauss fit" << endl;
+    *fLog << warn << dbginf << "WARNING: Could not create fit function for Gauss fit " 
+          << "in pixel: " << fPixId << endl;
     return kFALSE;
     }
@@ -528 +589 @@
   // The fit result is accepted under condition:
   // 1) The results are not nan's
-  // 2) The NDF is not smaller than fNDFLimit (5)
-  // 3) The Probability is greater than fProbLimit (default 0.001 == 99.9%)
+  // 2) The NDF is not smaller than fNDFLimit (default: fgNDFLimit)
+  // 3) The Probability is greater than fProbLimit (default: fgProbLimit)
   //
   if (   TMath::IsNaN(fMean) 
@@ -545 +606 @@
 
 // -----------------------------------------------------------------------------
+//
+// If flag IsGausFitOK() is set (histogram already successfully fitted), 
+// returns kTRUE
+// 
+// If both fMean and fSigma are still zero, call FitGaus() 
+// 
+// Repeats the Gauss fit in a smaller range, defined by: 
+// 
+// min = GetMean() - fPickupLimit * GetSigma();
+// max = GetMean() + fPickupLimit * GetSigma();
+//
+// The fit results are retrieved and stored in class-own variables.  
+//
+// A flag IsGausFitOK() is set according to whether the fit probability 
+// is smaller or bigger than fProbLimit, whether the NDF is bigger than 
+// fNDFLimit and whether results are NaNs.
+//
+Bool_t MHGausEvents::RepeatFit(const Option_t *option)
+{
+
+  if (IsGausFitOK())
+    return kTRUE;
+
+  if ((fMean == 0.) && (fSigma == 0.))
+    return FitGaus();
+
+  //
+  // Get new fitting ranges
+  //
+  Axis_t rmin = fMean - fPickupLimit * fSigma;
+  Axis_t rmax = fMean + fPickupLimit * fSigma;
+
+  Axis_t hmin = fHGausHist.GetBinCenter(fHGausHist.GetXaxis()->GetFirst());
+  Axis_t hmax = fHGausHist.GetBinCenter(fHGausHist.GetXaxis()->GetLast()) ;
+
+  fFGausFit->SetRange(hmin < rmin ? rmin : hmin , hmax > rmax ? rmax : hmax);
+
+  fHGausHist.Fit(fFGausFit,option);
+
+  fMean     = fFGausFit->GetParameter(1);
+  fMeanErr  = fFGausFit->GetParameter(2);
+  fSigma    = fFGausFit->GetParError(1) ; 
+  fSigmaErr = fFGausFit->GetParError(2) ; 
+  fProb     = fFGausFit->GetProb()      ;      
+
+  //
+  // The fit result is accepted under condition:
+  // 1) The results are not nan's
+  // 2) The NDF is not smaller than fNDFLimit (default: fgNDFLimit)
+  // 3) The Probability is greater than fProbLimit (default: fgProbLimit)
+  //
+  if (   TMath::IsNaN ( fMean     ) 
+      || TMath::IsNaN ( fMeanErr  )
+      || TMath::IsNaN ( fProb     )    
+      || TMath::IsNaN ( fSigma    )
+      || TMath::IsNaN ( fSigmaErr ) 
+      || fFGausFit->GetNDF() < fNDFLimit 
+      || fProb < fProbLimit )
+    return kFALSE;
+  
+  SetGausFitOK(kTRUE);
+  return kTRUE;
+
+}
+
+// -----------------------------------------------------------------------------
 // 
 // Bypasses the Gauss fit by taking mean and RMS from the histogram
@@ -559 +686 @@
   if (entries <= 0.)
     {
-      *fLog << warn << GetDescriptor() << ": Cannot bypass fit. Number of entries smaller or equal 0" << endl;
+      *fLog << warn << GetDescriptor() 
+            << ": Cannot bypass fit. Number of entries smaller or equal 0 in pixel: " << fPixId << endl;
       return;
     }
   
-  SetMean     ( fHGausHist.GetMean() );
-  SetMeanErr  ( fHGausHist.GetRMS() / TMath::Sqrt(entries) );
-  SetSigma    ( fHGausHist.GetRMS() );
-  SetSigmaErr ( fHGausHist.GetRMS() / TMath::Sqrt(entries) / 2. );
+  fMean     = fHGausHist.GetMean();
+  fMeanErr  = fHGausHist.GetRMS() / TMath::Sqrt(entries);
+  fSigma    = fHGausHist.GetRMS() ;
+  fSigmaErr = fHGausHist.GetRMS() / TMath::Sqrt(entries) / 2.;
 }
 
@@ -577 +705 @@
   
   *fLog << all                                                        << endl;
-  *fLog << all << "Results of the Gauss Fit: "                        << endl;
+  *fLog << all << "Results of the Gauss Fit in pixel: " << fPixId     << endl;
   *fLog << all << "Mean: "             << GetMean()                   << endl;
   *fLog << all << "Sigma: "            << GetSigma()                  << endl;
@@ -585 +713 @@
   *fLog << all                                                        << endl;
   
+}
+
+// --------------------------------------------------------------------------
+//
+// - Set fPixId to id
+//
+// Add id to names and titles of:
+// - fHGausHist
+//
+void MHGausEvents::ChangeHistId(Int_t id)
+{
+
+  fPixId = id;
+
+  fHGausHist.SetName(  Form("%s%d", fHGausHist.GetName(),  id));
+  fHGausHist.SetTitle( Form("%s%d", fHGausHist.GetTitle(), id));
+
+  fName  = Form("%s%d", fName.Data(),  id);
+  fTitle = Form("%s%d", fTitle.Data(), id);
+
 }
 
@@ -766 +914 @@
 }
 
-
+const Double_t MHGausEvents::GetPickup() const 
+{
+  
+  if ((fMean == 0.) && (fSigma == 0.))
+    return -1.;
+
+  const Int_t first = fHGausHist.GetXaxis()->FindBin(fMean+fPickupLimit*fSigma);
+  const Int_t last  = fHGausHist.GetXaxis()->GetLast();
+
+  if (first >= last)
+    return 0.;
+  
+  return fHGausHist.Integral(first, last, "width");
+}
+
+

trunk/MagicSoft/Mars/mcalib/MHGausEvents.h

@@ -20 +20 @@
 private:
 
-  const static Float_t  fgProbLimit;            // Default for fProbLimit (now 0.001)
-  const static Int_t    fgNDFLimit;             // Default for fNDFLimit  (now 2)
-  const static Int_t    fgPowerProbabilityBins; // Default for fPowerProbabilityBins (now 20)
-  const static Int_t    fgBinsAfterStripping;   // Default for fBinsAfterStripping (now 40)
-
+  const static Float_t  fgProbLimit;            //! Default for fProbLimit (now set to: 0.001)
+  const static Int_t    fgNDFLimit;             //! Default for fNDFLimit  (now set to: 2)
+  const static Float_t  fgPickupLimit;          //! Default for fPickupLimit (now set to: 5. )
+  const static Int_t    fgPowerProbabilityBins; //! Default for fPowerProbabilityBins (now set to: 20)
+  const static Int_t    fgBinsAfterStripping;   //! Default for fBinsAfterStripping (now set to: 40)
+  
   Int_t    fPowerProbabilityBins;      // Bins for the projected power spectrum
   Int_t    fBinsAfterStripping;        // Bins for the Gauss Histogram after stripping off the zeros at both ends
@@ -41 +42 @@
   Double_t fProb;                      //  Probability of the Gauss fit (derived from Chi-Square and NDF
 
-  enum { kGausFitOK, kExpFitOK, kFourierSpectrumOK }; // Bits to hold information about fit results 
+  enum { kGausFitOK, kExpFitOK, kFourierSpectrumOK,
+         kExcluded };                  //  Bits to hold information about fit results 
   
-  Byte_t fFlags;                       //   Byte to hold the bits fit result bits
+  Byte_t fFlags;                       //  Byte to hold the bits fit result bits
   
-  Int_t fCurrentSize;                  //   Current size of the array fEvents
-  
+  Int_t fCurrentSize;                  //  Current size of the array fEvents
+
 protected:
 
+  Int_t   fNbins;                      // Number histogram bins for fHGausHist (used by Init())
+  Axis_t  fFirst;                      // Lower histogram edge  for fHGausHist (used by Init()) 
+  Axis_t  fLast;                       // Upper histogram edge  for fHGausHist (used by Init()) 
+  Int_t   fPixId;                      //  Pixel ID 
+  
   TH1F    fHGausHist;                  // Histogram which should hold the Gaussian distribution
   TArrayF fEvents;                     // Array which holds the entries of GausHist
@@ -57 +64 @@
   Float_t  fProbLimit;                 // Probability limit for judgement if fit is OK 
   Int_t    fNDFLimit;                  // NDF limit for judgement if fit is OK
+  Float_t  fPickupLimit;               // Upper number of sigmas from the mean until events are considered as pickup
 
   Float_t *CreateEventXaxis(Int_t n);  //   Create an x-axis for the Event TGraphs
@@ -75 +83 @@
   virtual void Clear(Option_t *o="");
   virtual void Reset();  
-
+  virtual void InitBins();
+  
   // Setters
   void  SetEventFrequency(const Float_t f)   { fEventFrequency = f; }
@@ -84 +93 @@
   void  SetSigmaErr( const Double_t d )   { fSigmaErr = d;   }
   void  SetProb    ( const Double_t d )   { fProb     = d;   }
+  void  SetPixId    ( const Int_t i    )   { fPixId    = i;   }
+
+  void  SetNbins    ( const Int_t i    )   { fNbins    = i;   }  
+  void  SetFirst   ( const Double_t d )   { fFirst    = d;   }
+  void  SetLast    ( const Double_t d )   { fLast     = d;   }
   
+  void  SetNDFLimit(  const Int_t   lim=fgNDFLimit  ) {  fNDFLimit = lim;  }  
+  void  SetPickupLimit( const Float_t  lim =fgPickupLimit)  { fPickupLimit  = lim;   }
   void  SetProbLimit( const Float_t lim=fgProbLimit ) {  fProbLimit = lim; }
-  void  SetNDFLimit(  const Int_t   lim=fgNDFLimit  ) {  fNDFLimit = lim;  }  
 
   // Setters ONLY for MC:
-  void  SetGausFitOK(         const Bool_t b );
-  void  SetExpFitOK(          const Bool_t b );
-  void  SetFourierSpectrumOK( const Bool_t b );
+  void  SetGausFitOK(         const Bool_t b=kTRUE );
+  void  SetExpFitOK(          const Bool_t b=kTRUE );
+  void  SetFourierSpectrumOK( const Bool_t b=kTRUE );
+  void  SetExcluded         ( const Bool_t b=kTRUE );  
 
   // Getters
@@ -101 +117 @@
   const Double_t GetProb()          const { return fProb;      }
   const Int_t    GetNdf()           const;
+  const Double_t GetPickup()        const;
+  const Int_t    GetPixId()         const { return fPixId;     }
 
   const Double_t GetSlope()        const;
@@ -136 +154 @@
   const Bool_t IsEmpty()              const;
   const Bool_t IsFourierSpectrumOK()  const;
+  const Bool_t IsExcluded         ()  const;
 
   // Fill
@@ -145 +164 @@
   Bool_t FitGaus(Option_t *option="RQ0",
                  const Double_t xmin=0., const Double_t xmax=0.); // Fit the histogram HGausHist with a Gaussian
-  void BypassFit();   // Take mean and RMS from the histogram
-  
+  Bool_t RepeatFit(const Option_t *option="RQ0");      // Repeat fit within limits defined by fPickupLimit
+  void BypassFit();                                    // Take mean and RMS from the histogram
   
   // Draws
@@ -155 +174 @@
   
   // Miscelleaneous
+  virtual void ChangeHistId(Int_t id);   // Changes names and titles of the histogram
   void CreateFourierSpectrum();             // Create the fourier spectrum out of fEvents
   void CreateGraphEvents();                 // Create the TGraph fGraphEvents of fEvents