Index: trunk/MagicSoft/Mars/manalysis/MHPedestalCam.cc
===================================================================
--- trunk/MagicSoft/Mars/manalysis/MHPedestalCam.cc	(revision 3642)
+++ trunk/MagicSoft/Mars/manalysis/MHPedestalCam.cc	(revision 3643)
@@ -28,6 +28,64 @@
 // MHPedestalCam                                                            //
 //                                                                         //
-// Hold the Pedestal information for all pixels in the camera              //
-//                                                                         //
+// Fills the extracted pedestals of MExtractedSignalCam into 
+// the MHGausEvents-classes MHPedestalPix for every:
+//
+// - Pixel, stored in the TObjArray's MHCalibrationCam::fHiGainArray  
+//   or MHCalibrationCam::fHiGainArray, respectively, depending if 
+//   MArrivalTimePix::IsLoGainUsed() is set.
+//
+// - Average pixel per AREA index (e.g. inner and outer for the MAGIC camera), 
+//   stored in the TObjArray's MHCalibrationCam::fAverageHiGainAreas and 
+//   MHCalibrationCam::fAverageHiGainAreas
+//
+// - Average pixel per camera SECTOR (e.g. sectors 1-6 for the MAGIC camera), 
+//   stored in the TObjArray's MHCalibrationCam::fAverageHiGainSectors 
+//   and MHCalibrationCam::fAverageHiGainSectors 
+//
+// Every pedestal is directly taken from MExtractedSignalCam, filled by the 
+// appropriate extractor. 
+//
+// The pedestals are 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 MHGausEvents::fProbLimit (default: 0.5%), 
+// the fit is declared valid.
+// Otherwise, the fit is repeated within ranges of the previous mean 
+// +- MHGausEvents::fPickupLimit (default: 5) sigma (see MHGausEvents::RepeatFit())
+// In case this does not make the fit valid, the histogram means and RMS's are 
+// taken directly (see MHGausEvents::BypassFit()).
+// 
+// Outliers of more than MHGausEvents::fPickupLimit (default: 5) sigmas 
+// from the mean are counted as Pickup events (stored in MHGausEvents::fPickup) 
+//
+// The number of summed FADC slices is taken to re-normalize 
+// the result to a pedestal per pixel with the formulas (see MHPedestalPix::Renorm()): 
+// - Mean Pedestal        / slice = Mean Pedestal        / Number slices
+// - Mean Pedestal Error  / slice = Mean Pedestal Error  / Number slices
+// - Sigma Pedestal       / slice = Sigma Pedestal       / Sqrt (Number slices)
+// - Sigma Pedestal Error / slice = Sigma Pedestal Error / Sqrt (Number slices)
+// 
+// The class also fills arrays with the signal vs. event number, creates a fourier 
+// spectrum (see MHGausEvents::CreateFourierSpectrum()) and investigates if the 
+// projected fourier components follow an exponential distribution. 
+// 
+// This same procedure is performed for the average pixels.
+//
+// The following results are written into MPedestalCam:
+//
+// - 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 "MHPedestalCam.h"
@@ -46,4 +104,7 @@
 #include "MPedestalPix.h"
 
+#include "MGeomCam.h"
+#include "MGeomPix.h"
+
 ClassImp(MHPedestalCam);
 
@@ -51,86 +112,62 @@
 // --------------------------------------------------------------------------
 //
-// Default constructor. Creates a MHPedestalPix object for each pixel
+// Default constructor. 
+//
+// Initializes:
+// - fExtractHiGainSlices to 1.
+// - fExtractLoGainSlices to 1.
 //
 MHPedestalCam::MHPedestalCam(const char *name, const char *title) 
-    : fExtractSlices(0.)
+    : fExtractHiGainSlices(1.), fExtractLoGainSlices(1.)
 {
     fName  = name  ? name  : "MHPedestalCam";
     fTitle = title ? title : "";
 
-    //
-    // 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();
-
-}
+}
+
+
 
 // --------------------------------------------------------------------------
 //
-// Delete the array conatining the pixel pedest information
-//
-MHPedestalCam::~MHPedestalCam()
-{
-  delete fArray;
-}
-
-// --------------------------------------------------------------------------
-//
-// Get i-th pixel (pixel number)
-//
-MHPedestalPix &MHPedestalCam::operator[](UInt_t i)
-{
-  return *(MHPedestalPix*)(fArray->At(i));
-}
-
-// --------------------------------------------------------------------------
-//
-// Get i-th pixel (pixel number)
-//
-const MHPedestalPix &MHPedestalCam::operator[](UInt_t i) const
-{
-  return *(MHPedestalPix*)(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 *MHPedestalCam::Clone(const char *) const
-{
-
-  const Int_t n = fArray->GetSize();
-  
-  //
-  // FIXME, this might be done faster and more elegant, by direct copy.
-  //
-  MHPedestalCam *cam = new MHPedestalCam;
-  
-  cam->fArray->Expand(n);
-  
-  for (int i=0; i<n; i++)
-    {
-      delete (*cam->fArray)[i];
-      (*cam->fArray)[i] = (*fArray)[i]->Clone();
-    }
-
-  return cam;
-}
-
-
-
-// --------------------------------------------------------------------------
-//
-// To setup the object we get the number of pixels from a MGeomCam object
-// in the Parameter list. 
-// MHPedestalPix sets its parameters to 0. (other than default which is -1.)
-//
-Bool_t MHPedestalCam::SetupFill(const MParList *pList)
-{
+// Searches pointer to:
+// - MPedestalCam
+// - MExtractedSignalCam
+//
+// Retrieves from MExtractedSignalCam:
+// - number of used High Gain FADC slices (MExtractedSignalCam::GetNumUsedHiGainFADCSlices())
+// - number of used Low  Gain FADC slices (MExtractedSignalCam::GetNumUsedLoGainFADCSlices())
+//
+// Initializes, if empty to MGeomCam::GetNumPixels():
+// - 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 MHCalibrationCam::InitPedHists() for every entry in:
+// - MHCalibrationCam::fHiGainArray, MHCalibrationCam::fLoGainArray
+// - MHCalibrationCam::fAverageHiGainAreas, MHCalibrationCam::fAverageLoGainAreas
+// - MHCalibrationCam::fAverageHiGainSectors, MHCalibrationCam::fAverageLoGainSectors
+//
+// 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 MHPedestalCam::ReInitHists(MParList *pList)
+{
+
+  MExtractedSignalCam *signal = (MExtractedSignalCam*)pList->FindObject("MExtractedSignalCam");
+  if (!signal)
+    {
+      gLog << err << "Cannot find MExtractedSignalCam... abort." << endl;
+      return kFALSE;
+    }
+  
 
   fPedestals = (MPedestalCam*)pList->FindObject("MPedestalCam");
@@ -142,12 +179,164 @@
     }
 
-  fArray->Delete();
+  Float_t sliceshi = signal->GetNumUsedHiGainFADCSlices();
+  Float_t sliceslo = signal->GetNumUsedLoGainFADCSlices();
+
+  if (sliceshi == 0.)
+    {
+      gLog << err << "Number of used signal slices in MExtractedSignalCam is zero  ... abort." 
+           << endl;
+      return kFALSE;
+    }
+
+  if (fExtractHiGainSlices != 0. && sliceshi != fExtractHiGainSlices )
+    {
+      gLog << err << "Number of used High Gain signal slices changed in MExtractedSignalCam  ... abort." 
+           << endl;
+      return kFALSE;
+    }
+
+  if (fExtractLoGainSlices != 0. && sliceslo != fExtractLoGainSlices )
+    {
+      gLog << err << "Number of used Low Gain signal slices changed in MExtractedSignalCam  ... abort." 
+           << endl;
+      return kFALSE;
+    }
+
+  fExtractHiGainSlices = sliceshi;
+  fExtractLoGainSlices = sliceslo;
+
+  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 MHPedestalPix;
+          InitPedHists((MHPedestalPix&)(*this)[i],i,fExtractHiGainSlices);
+      }
+  }
+
+  if (fLoGainArray->GetEntries()==0)
+  {
+      fLoGainArray->Expand(npixels);
+      for (Int_t i=0; i<npixels; i++)
+      {
+	  (*fLoGainArray)[i] = new MHPedestalPix;
+          InitPedHists((MHPedestalPix&)(*this)(i),i,fExtractLoGainSlices);
+      }
+  }
+
+  if (fAverageHiGainAreas->GetEntries()==0)
+  {
+    fAverageHiGainAreas->Expand(nareas);
+    
+    for (Int_t j=0; j<nareas; j++)
+      {
+        (*fAverageHiGainAreas)[j] = 
+          new MHPedestalPix("AverageHiGainArea",
+                                      "Average Pedestals area idx ");
+
+        InitPedHists((MHPedestalPix&)GetAverageHiGainArea(j),j,fExtractHiGainSlices);
+
+        GetAverageHiGainArea(j).GetHGausHist()->SetTitle("Pedestals 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 MHPedestalPix("AverageLoGainArea",
+                                      "Pedestals average Area idx ");
+
+        InitPedHists((MHPedestalPix&)GetAverageLoGainArea(j),j,fExtractLoGainSlices);
+
+        GetAverageLoGainArea(j).GetHGausHist()->SetTitle("Pedestals 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 MHPedestalPix("AverageHiGainSector",
+                                        "Pedestals average sector ");
+
+          InitPedHists((MHPedestalPix&)GetAverageHiGainSector(j),j,fExtractHiGainSlices);
+
+          GetAverageHiGainSector(j).GetHGausHist()->SetTitle("Pedestals average Sector ");
+          GetAverageHiGainSector(j).SetNbins(fAverageNbins);
+      }
+  }
+
+  if (fAverageLoGainSectors->GetEntries()==0)
+  {
+      fAverageLoGainSectors->Expand(nsectors);
+
+      for (Int_t j=0; j<nsectors; j++)
+      {
+	  (*fAverageLoGainSectors)[j] = 
+            new MHPedestalPix("AverageLoGainSector",
+                                        "Pedestals average sector ");
+
+          InitPedHists((MHPedestalPix&)GetAverageLoGainSector(j),j,fExtractLoGainSlices);
+          
+          GetAverageLoGainSector(j).GetHGausHist()->SetTitle("Pedestals average Sector ");
+          GetAverageLoGainSector(j).SetNbins(fAverageNbins);
+      }
+  }
+
   return kTRUE;
 
 }
 
-// --------------------------------------------------------------------------
-//
-Bool_t MHPedestalCam::Fill(const MParContainer *par, const Stat_t w)
+
+// -------------------------------------------------------------
+//
+// If MBadPixelsPix::IsBad():
+// - calls MHGausEvents::SetExcluded()
+//
+// Calls:
+// - MHGausEvents::InitBins()
+// - MHGausEvents::ChangeHistId(i)
+// - MHGausEvents::SetEventFrequency(fPulserFrequency)
+// - MHPedestalPix::SetNSlices(nslices)
+//
+void MHPedestalCam::InitPedHists(MHPedestalPix &hist, const Int_t i, const Float_t nslices)
+{
+  
+  hist.InitBins();
+  hist.ChangeHistId(i);
+  hist.SetEventFrequency(fPulserFrequency);
+  hist.SetNSlices(nslices);
+}
+
+
+
+// -------------------------------------------------------------------------------
+//
+// Retrieves pointer to MExtractedSignalCam:
+//
+// Retrieves from MGeomCam:
+// - number of pixels
+// - number of pixel areas
+// - number of sectors
+//
+// Fills HiGain or LoGain histograms (MHGausEvents::FillHistAndArray()), respectively
+// with the signals MExtractedSignalCam::GetExtractedSignalHiGain() and 
+// MExtractedSignalCam::GetExtractedSignalLoGain(), respectively.
+//
+Bool_t MHPedestalCam::FillHists(const MParContainer *par, const Stat_t w)
 {
 
@@ -159,59 +348,73 @@
     }
   
-
-  Float_t slices = signal->GetNumUsedHiGainFADCSlices();
-
-  if (slices == 0.)
-    {
-      gLog << err << "Number of used signal slices in MExtractedSignalCam is zero  ... abort." 
-           << endl;
-      return kFALSE;
-    }
-
-  if (fExtractSlices != 0. && slices != fExtractSlices )
-    {
-      gLog << err << "Number of used signal slices changed in MExtractedSignalCam  ... abort." 
-           << endl;
-      return kFALSE;
-    }
-
-  fExtractSlices = slices;
-
-  const Int_t n = signal->GetSize();
-  
-  if (fArray->GetEntries()==0)
-    {
-      fArray->Expand(n);
-      
-      for (Int_t i=0; i<n; i++)
-        {
-          (*fArray)[i] = new MHPedestalPix;
-          MPedestalPix  &pix  = (*fPedestals)[i];
-          pix.InitUseHists();
-          MHPedestalPix &hist = (*this)[i];
-          hist.ChangeHistId(i);
-          hist.InitBins();
-        }
-    }
-  
-  if (fArray->GetEntries() != n)
-    {
-      gLog << err << "ERROR - Size mismatch... abort." << endl;
-      return kFALSE;
-    }
-  
-  for (Int_t i=0; i<n; i++)
-    {
+  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;
+    }
+  
+
+  for (Int_t i=0; i<npixels; i++)
+    {
+
+      MHGausEvents &histhi = (*this)[i];
+      MHGausEvents &histlo = (*this)(i);
+
+      if (histhi.IsExcluded())
+	continue;
 
       const MExtractedSignalPix &pix = (*signal)[i];
       
-      const Float_t signal = pix.GetExtractedSignalHiGain();
-
-      MHPedestalPix  &hist = (*this)[i];
-      //
-      // Don't fill signal per slice, we get completely screwed up 
-      // with the sigma. Better fill like it is and renorm later
-      //
-      hist.FillHistAndArray(signal);
+      const Float_t pedhi = pix.GetExtractedSignalHiGain();
+      const Float_t pedlo = pix.GetExtractedSignalLoGain();
+
+      const Int_t aidx   = (*fGeom)[i].GetAidx();
+      const Int_t sector = (*fGeom)[i].GetSector();
+
+      histhi.FillHistAndArray(pedhi) ;
+      sumareahi  [aidx]   += pedhi;
+      numareahi  [aidx]   ++;
+      sumsectorhi[sector] += pedhi;
+      numsectorhi[sector] ++;
+
+      histlo.FillHistAndArray(pedlo);
+      sumarealo  [aidx]   += pedlo;
+      numarealo  [aidx]   ++;
+      sumsectorlo[sector] += pedlo;
+      numsectorlo[sector] ++;
+
+    }
+  
+  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]);
     }
   
@@ -219,30 +422,25 @@
 }
 
-Bool_t MHPedestalCam::Finalize()
-{
-    for (Int_t i=0; i<fArray->GetSize(); i++)
-    {
-
-        MHPedestalPix &hist = (*this)[i];
-
-        //
-        // 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();
-        hist.Renorm(fExtractSlices);
-        hist.CreateFourierSpectrum();
-
-    }
-    return kTRUE;
-}
-
 // --------------------------------------------------------------------------
+//
+// Calls:
+// - MHCalibrationCam::FitHiGainArrays() with Bad Pixels flags 0
+// - MHCalibrationCam::FitLoGainArrays() with Bad Pixels flags 0
+// 
+Bool_t MHPedestalCam::FinalizeHists()
+{
+
+  FitHiGainArrays((MCalibrationCam&)(*fCam),*fBadPixels, 
+                  MBadPixelsPix::kHiGainNotFitted,
+                  MBadPixelsPix::kHiGainOscillating);
+  FitLoGainArrays((MCalibrationCam&)(*fCam),*fBadPixels,
+                  MBadPixelsPix::kLoGainNotFitted,
+                  MBadPixelsPix::kLoGainOscillating);
+
+  return kTRUE;
+  
+}
+
+// ------------------------------------------------------------------
 //
 // The types are as follows:
@@ -277,5 +475,5 @@
 {
 
-  if (fArray->GetSize() <= idx)
+  if (fHiGainArray->GetSize() <= idx)
     return kFALSE;
 
@@ -344,4 +542,8 @@
 }
 
+// --------------------------------------------------------------------------
+//
+// Calls MHGausEvents::DrawClone() for pixel idx
+//
 void MHPedestalCam::DrawPixelContent(Int_t idx) const
 {
Index: trunk/MagicSoft/Mars/manalysis/MHPedestalCam.h
===================================================================
--- trunk/MagicSoft/Mars/manalysis/MHPedestalCam.h	(revision 3642)
+++ trunk/MagicSoft/Mars/manalysis/MHPedestalCam.h	(revision 3643)
@@ -32,5 +32,5 @@
   void DrawPixelContent(Int_t idx) const;
   
-  ClassDef(MHPedestalCam, 1)	// Storage Container for all pedestal information of the camera
+  ClassDef(MHPedestalCam, 1)	// Histogram class for Charge Camera Pedestals 
 };
 
Index: trunk/MagicSoft/Mars/manalysis/MHPedestalPix.h
===================================================================
--- trunk/MagicSoft/Mars/manalysis/MHPedestalPix.h	(revision 3642)
+++ trunk/MagicSoft/Mars/manalysis/MHPedestalPix.h	(revision 3643)
@@ -11,5 +11,5 @@
 private:
 
-  static const Int_t   fgChargeNbins;        // Default for fNBins          (now set to: 450  )
+  static const Int_t   fgChargeNbins;        // Default for fNbins          (now set to: 450  )
   static const Axis_t  fgChargeFirst;        // Default for fFirst          (now set to: -0.5  )
   static const Axis_t  fgChargeLast;         // Default for fLast           (now set to: 449.5)
@@ -31,5 +31,5 @@
   void Renorm();  
 
-  ClassDef(MHPedestalPix, 1)     // Histograms for each calibrated pixel
+  ClassDef(MHPedestalPix, 1)     // Histogram class for Charge Pedestal Pixel
 };