/* ======================================================================== *\
!
! *
! * This file is part of MARS, the MAGIC Analysis and Reconstruction
! * Software. It is distributed to you in the hope that it can be a useful
! * and timesaving tool in analysing Data of imaging Cerenkov telescopes.
! * It is distributed WITHOUT ANY WARRANTY.
! *
! * Permission to use, copy, modify and distribute this software and its
! * documentation for any purpose is hereby granted without fee,
! * provided that the above copyright notice appear in all copies and
! * that both that copyright notice and this permission notice appear
! * in supporting documentation. It is provided "as is" without express
! * or implied warranty.
! *
!
!
!   Author(s): Markus Gaug   11/2003 <mailto:markus@ifae.es>
!
!   Copyright: MAGIC Software Development, 2000-2004
!
!
\* ======================================================================== */
/////////////////////////////////////////////////////////////////////////////
//                                                               
// MCalibrationChargeCam                                               
//                                                               
// Storage container for charge calibration results from the signal distribution 
// fits (see MHCalibrationChargeCam and MHCalibrationChargePix), the calculation 
// of reduced sigmas and number of photo-electrons (this class) and conversion 
// factors sum FADC slices to photo-electrons (see MCalibrationChargeCalc)
//
// Individual pixels have to be cast when retrieved e.g.:
// MCalibrationChargePix &avpix = (MCalibrationChargePix&)(*fChargeCam)[i]
// 
// Averaged values over one whole area index (e.g. inner or outer pixels for 
// the MAGIC camera), can be retrieved via: 
// MCalibrationChargePix &avpix = (MCalibrationChargePix&)fChargeCam->GetAverageArea(i)
//
// Averaged values over one whole camera sector can be retrieved via: 
// MCalibrationChargePix &avpix = (MCalibrationChargePix&)fChargeCam->GetAverageSector(i)
//
// Note the averageing has been done on an event-by-event basis. Resulting 
// Sigma's of the Gauss fit have been multiplied with the square root of the number 
// of involved pixels in order to make a direct comparison possible with the mean of 
// sigmas. 
//
// Final numbers of uncalibrated or unreliable pixels can be retrieved via the commands:
// GetNumUncalibrated(aidx) and GetNumUnreliable(aidx) where aidx is the area index (0 for 
// inner and 1 for outer pixels in the MAGIC camera).
//
// The following "calibration" constants are used for the calibration of each pixel 
// (see MCalibrate):
//
// - MCalibrationQEPix::GetMeanConvFADC2Phe(): The mean conversion factor from the 
//   summed FADC slices to the number of photo-electrons (in first order independent 
//   on colour and intensity)
// - MCalibrationQEPix::GetMeanFFactorFADC2Phot(): The mean F-Factor of the total  
//   readout chain dividing the signal to noise of the incoming number of photons 
//   (= sqrt(number photons)) by the signal to noise of the outgoing summed FADC slices 
//   signal (= MCalibrationPix::GetMean() / MCalibrationChargePix::GetRSigma() )
//
// The following calibration constants can be retrieved directly from this class:
//
// - GetConversionFactorFFactor    ( Int_t idx, Float_t &mean, Float_t &err, Float_t &sigma );
//
// where: 
// - idx is the pixel software ID
// - "mean" is the mean conversion constant, to be multiplied with the retrieved signal 
//   in order to get a calibrated number of PHOTONS. 
// - "err" is the pure statistical uncertainty about the MEAN
// - "sigma", if mulitplied with the square root of signal, gives the approximate sigma of the 
//   retrieved mean number of incident Cherenkov photons.
//
// Note, Conversion is ALWAYS (included the F-Factor method!) from summed FADC slices to PHOTONS.
//
// See also: MCalibrationChargePix, MCalibrationChargeCalc, MCalibrationQECam
//           MHCalibrationChargePix, MHCalibrationChargeCam              
//           MCalibrationChargeBlindPix, MCalibrationChargePINDiode
//
/////////////////////////////////////////////////////////////////////////////
#include "MCalibrationChargeCam.h"

#include <TClonesArray.h>

#include "MLog.h"
#include "MLogManip.h"

#include "MGeomCam.h"
#include "MGeomPix.h"

#include "MBadPixelsCam.h"
#include "MBadPixelsPix.h"

#include "MCalibrationQECam.h"
#include "MCalibrationQEPix.h"

#include "MCalibrationChargePix.h"
#include "MCalibrationChargeBlindPix.h"
#include "MCalibrationChargePINDiode.h"

ClassImp(MCalibrationChargeCam);

using namespace std;
// --------------------------------------------------------------------------
//
// Default constructor. 
//
// Sets all pointers to 0
// 
// Creates a TClonesArray of MCalibrationChargePix containers, initialized to 1 entry, destinated 
// to hold one container per pixel. Later, a call to MCalibrationChargeCam::InitSize() 
// has to be performed (in MGeomApply). 
//
// Creates a TClonesArray of MCalibrationChargePix containers, initialized to 1 entry, destinated 
// to hold one container per pixel AREA. Later, a call to MCalibrationChargeCam::InitAreas() 
// has to be performed (in MGeomApply). 
//
// Creates a TClonesArray of MCalibrationChargePix containers, initialized to 1 entry, destinated
// to hold one container per camera SECTOR. Later, a call to MCalibrationChargeCam::InitSectors() 
// has to be performed (in MGeomApply). 
//
// Calls:
// - Clear()
//
MCalibrationChargeCam::MCalibrationChargeCam(const char *name, const char *title)
{
  fName  = name  ? name  : "MCalibrationChargeCam";
  fTitle = title ? title : "Storage container for the Calibration Information in the camera";
  
  fPixels           = new TClonesArray("MCalibrationChargePix",1);
  fAverageAreas     = new TClonesArray("MCalibrationChargePix",1);
  fAverageSectors   = new TClonesArray("MCalibrationChargePix",1);

  Clear();
}


// --------------------------------------
//
// Sets all variable to 0.
// Sets all flags to kFALSE
// Calls MCalibrationCam::Clear()
//
void MCalibrationChargeCam::Clear(Option_t *o)
{

  SetFFactorMethodValid    ( kFALSE );

  fNumPhotonsBlindPixelMethod    = 0.; 
  fNumPhotonsFFactorMethod       = 0.; 
  fNumPhotonsPINDiodeMethod      = 0.; 
  fNumPhotonsBlindPixelMethodErr = 0.; 
  fNumPhotonsFFactorMethodErr    = 0.;  
  fNumPhotonsPINDiodeMethodErr   = 0.;  

  MCalibrationCam::Clear();

  return;
}

// -----------------------------------------------
//
// Sets the kFFactorMethodValid bit from outside
//
void MCalibrationChargeCam::SetFFactorMethodValid(const Bool_t b)
{
    b ? SETBIT(fFlags, kFFactorMethodValid) : CLRBIT(fFlags, kFFactorMethodValid); 
}


// --------------------------------------------------------------------------
//
// Test bit kFFactorMethodValid
//
Bool_t  MCalibrationChargeCam::IsFFactorMethodValid()   const
{
  return TESTBIT(fFlags,kFFactorMethodValid);
}

// --------------------------------------------------------------------------
//
// Print first the well fitted pixels 
// and then the ones which are not FitValid
//
void MCalibrationChargeCam::Print(Option_t *o) const
{

  *fLog << all << GetDescriptor() << ":" << endl;
  int id = 0;
  
  *fLog << all << "Calibrated pixels:" << endl;
  *fLog << all << endl;

  TIter Next(fPixels);
  MCalibrationChargePix *pix;
  while ((pix=(MCalibrationChargePix*)Next()))
    {
      
      if (!pix->IsExcluded()) 
	{                            

          *fLog << all 
                << Form("%s%3i","Pixel: ",pix->GetPixId())
                << Form("%s%4.2f%s%4.2f","  Ped.Rms: ",pix->GetPedRms(),"+-",pix->GetPedRmsErr())
		<< Form("%s%4.2f%s%4.2f","  Charge: " ,pix->GetConvertedMean(),"+-",pix->GetConvertedSigma())
		<< Form("%s%4.2f%s%4.2f","  Red.Sigma: ",pix->GetConvertedRSigma(),"+-",pix->GetConvertedRSigmaErr())
		<< Form("%s%4.2f%s%4.2f","  Num.Phes: ",pix->GetPheFFactorMethod(),"+-",pix->GetPheFFactorMethodErr()) 
		<< Form("%s%4.2f%s%4.2f","  Conv.FADC2Phe: ",pix->GetMeanConvFADC2Phe(),"+-",pix->GetMeanConvFADC2PheErr())
                << " Saturated? :" << pix->IsHiGainSaturation() 
		<< endl;
          id++;
	}
    }
  
  *fLog << all << id << " pixels" << endl;
  id = 0;
  
   
  *fLog << all << endl;
  *fLog << all << "Excluded pixels:" << endl;
  *fLog << all << endl;
  
  id = 0;

  TIter Next4(fPixels);
  while ((pix=(MCalibrationChargePix*)Next4()))
  {
      if (pix->IsExcluded())
      {
        *fLog << all << pix->GetPixId() << " ";
        id++;

        if (!(id % 25))
          *fLog << endl;
      }
  }
  
  *fLog << endl;
  *fLog << all << id << " Excluded pixels " << endl;
  *fLog << endl;

  *fLog << all << endl;
  *fLog << all << "Averaged Areas:" << endl;
  *fLog << all << endl;

  TIter Next5(fAverageAreas);
  while ((pix=(MCalibrationChargePix*)Next5()))
  {
    *fLog << all << Form("%s%3i","Area Idx: ",pix->GetPixId())
          << "   Ped.  Rms: "            << pix->GetPedRms()        << " +- " << pix->GetPedRmsErr() 
          << "   Mean signal: "          << pix->GetMean()    << " +- " << pix->GetMeanErr() 
          << "   Sigma signal: "         << pix->GetSigma()    << " +- "<< pix->GetSigmaErr() 
          << "   Reduced Sigma: "        << pix->GetRSigma() 
          << "   Nr Phe's: "             << pix->GetPheFFactorMethod() 
          << endl;
  }

  *fLog << all << endl;
  *fLog << all << "Averaged Sectors:" << endl;
  *fLog << all << endl;

  TIter Next6(fAverageSectors);
  while ((pix=(MCalibrationChargePix*)Next6()))
  {
    *fLog << all << Form("%s%3i","Sector: ",pix->GetPixId())
          << "   Ped.  Rms: "            << pix->GetPedRms()        << " +- " << pix->GetPedRmsErr() 
          << "   Mean signal: "          << pix->GetMean()    << " +- " << pix->GetMeanErr() 
          << "   Sigma signal: "         << pix->GetSigma()    << " +- "<< pix->GetSigmaErr() 
          << "   Reduced Sigma: "        << pix->GetRSigma() 
          << "   Nr Phe's: "             << pix->GetPheFFactorMethod() 
          << endl;
  }
  *fLog << all << endl;
}


// --------------------------------------------------------------------------
//
// The types are as follows:
// 
// Fitted values:
// ============== 
//
// 0: Fitted Charge                              (see MCalibrationPix::GetMean())
// 1: Error of fitted Charge                     (see MCalibrationPix::GetMeanErr())
// 2: Sigma of fitted Charge                     (see MCalibrationPix::GetSigma())
// 3: Error of Sigma of fitted Charge            (see MCalibrationPix::GetSigmaErr())
//
// Useful variables derived from the fit results:
// =============================================
//
// 4: Probability Gauss fit Charge distribution  (see MCalibrationPix::GetProb())
// 5: Reduced Sigma of fitted Charge             (see MCalibrationChargePix::GetRSigma())
// 6: Error Reduced Sigma of fitted Charge       (see MCalibrationChargePix::GetRSigmaErr())
// 7: Reduced Sigma per Charge                   (see MCalibrationChargePix::GetRSigmaPerCharge())
// 8: Error of Reduced Sigma per Charge          (see MCalibrationChargePix::GetRSigmaPerChargeErr())
//
// Results of the F-Factor calibration Method:
// ===========================================
//
//  9: Nr. Photo-electrons from F-Factor Method  (see MCalibrationChargePix::GetPheFFactorMethod())
// 10: Error Nr. Photo-el. from F-Factor Method  (see MCalibrationChargePix::GetPheFFactorMethodErr())
// 11: Conversion factor   from F-Factor Method  (see MCalibrationChargePix::GetMeanConvFADC2Phe()
// 12: Error conv. factor  from F-Factor Method  (see MCalibrationChargePix::GetMeanConvFADC2PheErr()
// 13: Overall F-Factor    from F-Factor Method  (see MCalibrationChargePix::GetMeanFFactorFADC2Phot()
// 14: Error F-Factor      from F-Factor Method  (see MCalibrationChargePix::GetMeanFFactorFADC2PhotErr()
// 15: Pixels valid calibration F-Factor-Method  (see MCalibrationChargePix::IsFFactorMethodValid())           
//
// Results of the Low-Gain vs. High-Gain Conversion:
// =================================================
//
// 16: Mean Signal Hi-Gain / Mean Signal Lo-Gain (see MCalibrationPix::GetHiLoMeansDivided())
// 17: Error Signal High-Gain / Signal Low Gain  (see MCalibrationPix::GetHiLoMeansDividedErr())
// 18: Sigma High-Gain / Sigma Low Gain          (see MCalibrationPix::GetHiLoSigmasDivided())
// 19: Error Sigma High-Gain / Sigma Low Gain    (see MCalibrationPix::GetHiLoSigmasDividedErr())
//                                                
// Localized defects:                             
// ==================
//
// 20: Excluded Pixels
// 21: Number of pickup   events in the Hi Gain  (see MCalibrationPix::GetHiGainNumPickup())
// 22: Number of pickup   events in the Lo Gain  (see MCalibrationPix::GetLoGainNumPickup())
// 23: Number of blackout events in the Hi Gain  (see MCalibrationPix::GetHiGainNumBlackout())
// 24: Number of blackout events in the Lo Gain  (see MCalibrationPix::GetLoGainNumBlackout())
//
// Other classifications of pixels:
// ================================
//
// 25: Pixels with saturated High-Gain           (see MCalibrationPix::IsHiGainSaturation())
//
// Calculated absolute arrival times (very low precision!):
// ========================================================
//
// 26: Absolute Arrival time of the signal       (see MCalibrationChargePix::GetAbsTimeMean())
// 27: RMS Ab.  Arrival time of the signal       (see MCalibrationChargePix::GetAbsTimeRms())
//
// Used Pedestals:
// ===============
//
// 28: Pedestal for entire signal extr. range    (see MCalibrationChargePix::Ped())
// 29: Error Pedestal entire signal extr. range  (see MCalibrationChargePix::PedErr())
// 30: Ped. RMS entire signal extraction range   (see MCalibrationChargePix::PedRms())
// 31: Error Ped. RMS entire signal extr. range  (see MCalibrationChargePix::PedRmsErr())
//
Bool_t MCalibrationChargeCam::GetPixelContent(Double_t &val, Int_t idx, const MGeomCam &cam, Int_t type) const
{

  if (idx > GetSize())
    return kFALSE;

  Float_t area = cam[idx].GetA();

 if (area == 0)
    return kFALSE;

 MCalibrationChargePix &pix = (MCalibrationChargePix&)(*this)[idx];

  switch (type)
    {
    case 0:
      if (pix.IsExcluded())
        return kFALSE;
      val = pix.GetMean();
      break;
    case 1:
      if (pix.IsExcluded())
        return kFALSE;
      val = pix.GetMeanErr();
      break;
    case 2:
      if (pix.IsExcluded())
        return kFALSE;
      val = pix.GetSigma();
      break;
    case 3:
      if (pix.IsExcluded())
        return kFALSE;
      val = pix.GetSigmaErr();
      break;
    case 4:
      if (pix.IsExcluded())
        return kFALSE;
      val = pix.GetProb();
      break;
    case 5:
      if (pix.IsExcluded())
        return kFALSE;
      if (pix.GetRSigma() == -1.)
	  return kFALSE;
      val = pix.GetRSigma();
      break;
    case 6:
      if (pix.IsExcluded())
        return kFALSE;    
      if (pix.GetRSigma() == -1.)
	  return kFALSE;
      val = pix.GetRSigmaErr();
      break;
    case 7:
      if (pix.IsExcluded())
        return kFALSE;
      val = pix.GetRSigmaPerCharge();
      break;
    case 8:
      if (pix.IsExcluded())
        return kFALSE;
      val =  pix.GetRSigmaPerChargeErr();
      break;
    case 9:
      if (pix.IsExcluded() || !pix.IsFFactorMethodValid())
        return kFALSE;
      val = pix.GetPheFFactorMethod();
      break;
    case 10:
      if (pix.IsExcluded() || !pix.IsFFactorMethodValid())
        return kFALSE;
      val = pix.GetPheFFactorMethodErr();
      break;
    case 11:
      if (pix.IsExcluded() || !pix.IsFFactorMethodValid())
        return kFALSE;
      val = pix.GetMeanConvFADC2Phe();
      break;
    case 12:
      if (pix.IsExcluded() || !pix.IsFFactorMethodValid())
        return kFALSE;
      val = pix.GetMeanConvFADC2PheErr();
      break;
    case 13:
      if (pix.IsExcluded() || !pix.IsFFactorMethodValid())
        return kFALSE;
      val = pix.GetMeanFFactorFADC2Phot();
      break;
    case 14:
      if (pix.IsExcluded() || !pix.IsFFactorMethodValid())
        return kFALSE;
      val = pix.GetMeanFFactorFADC2PhotErr();
      break;
    case 15:
      if (pix.IsExcluded())
        return kFALSE;
      if (pix.IsFFactorMethodValid())
        val = 1;
      else
        return kFALSE;
      break;
    case 16:
      if (pix.IsExcluded())
        return kFALSE;
      val = pix.GetHiLoMeansDivided();
      break;
    case 17:
      if (pix.IsExcluded())
        return kFALSE;
      val = pix.GetHiLoMeansDividedErr();
      break;
    case 18:
      if (pix.IsExcluded())
        return kFALSE;
      val = pix.GetHiLoSigmasDivided();
      break;
    case 19:
      if (pix.IsExcluded())
        return kFALSE;
      val = pix.GetHiLoSigmasDividedErr();
      break;
    case 20:
      if (pix.IsExcluded())
        val = 1.;
      else
        return kFALSE;
      break;
    case 21:
      if (pix.IsExcluded())
        return kFALSE;
      val = pix.GetHiGainNumPickup();
      break;
    case 22:
      if (pix.IsExcluded())
        return kFALSE;
      val = pix.GetLoGainNumPickup();
      break;
    case 23:
      if (pix.IsExcluded())
        return kFALSE;
      val = pix.GetHiGainNumBlackout();
      break;
    case 24:
      if (pix.IsExcluded())
        return kFALSE;
      val = pix.GetLoGainNumBlackout();
      break;
    case 25:
      if (pix.IsExcluded())
        return kFALSE;
      val = pix.IsHiGainSaturation();
      break;
    case 26:
      if (pix.IsExcluded())
        return kFALSE;
      val = pix.GetAbsTimeMean();
      break;
    case 27:
      if (pix.IsExcluded())
        return kFALSE;
      val = pix.GetAbsTimeRms();
      break;
    case 28:
      if (pix.IsExcluded())
        return kFALSE;
      val = pix.GetPed();
      break;
    case 29:
      if (pix.IsExcluded())
        return kFALSE;
      val = pix.GetPedErr();
      break;
    case 30:
      if (pix.IsExcluded())
        return kFALSE;
      val = pix.GetPedRms();
      break;
    case 31:
      if (pix.IsExcluded())
        return kFALSE;
      val = pix.GetPedErr()/2.;
      break;
    default:
      return kFALSE;
    }

  return val!=-1.;
}

// --------------------------------------------------------------------------
//
// Calls MCalibrationChargePix::DrawClone()
//
void MCalibrationChargeCam::DrawPixelContent(Int_t idx) const
{
  MCalibrationChargePix &pix = (MCalibrationChargePix&)(*this)[idx];
  pix.DrawClone();
}



Bool_t MCalibrationChargeCam::GetConversionFactorFFactor(Int_t ipx, Float_t &mean, Float_t &err, Float_t &ffactor)
{

  MCalibrationChargePix &pix = (MCalibrationChargePix&)(*this)[ipx];

  Float_t conv = pix.GetMeanConvFADC2Phe();

  if (conv < 0.)
    return kFALSE;

  mean    = conv;
  err     = pix.GetMeanConvFADC2PheErr();
  ffactor = pix.GetMeanFFactorFADC2Phot();

  return kTRUE;
}


// --------------------------------------------------------------------------
//
// Calculates the average conversion factor FADC counts to photons for pixel sizes. 
// The geometry container is used to get the necessary
// geometry information (area index). The MCalibrationQECam container is necessary for 
// the quantum efficiency information. 
// If the bad pixel container is given all pixels which have the flag 'kUnsuitableRun' are ignored
// in the calculation of the size average.
//
// Returns a TArrayF of dimension 2: 
// arr[0]: averaged conversion factors (default: -1.)
// arr[1]: Error (rms) of averaged conversion factors (default: 0.)
//
// ATTENTION: THE USER HAS TO DELETE THE RETURNED TARRAYF ACCORDINGLY
//
TArrayF *MCalibrationChargeCam::GetAveragedConvFADC2PhotPerArea  ( const MGeomCam &geom, const MCalibrationQECam &qecam,
                                           const UInt_t ai,  MBadPixelsCam *bad)
{

  const Int_t np = GetSize();

  Double_t mean  = 0.;
  Double_t mean2 = 0.;
  Int_t    nr    = 0;

  for (int i=0; i<np; i++)
    {
      if (bad && (*bad)[i].IsUnsuitable(MBadPixelsPix::kUnsuitableRun))
        continue; 
      
      const UInt_t aidx = geom[i].GetAidx();
      
      if (ai != aidx)
        continue;

      const MCalibrationQEPix &qepix = (MCalibrationQEPix&)qecam[i];
      if (!qepix.IsAverageQECombinedAvailable())
        continue;
      
      const MCalibrationChargePix &pix = (MCalibrationChargePix&)(*this)[i];
      const Float_t conv = pix.GetMeanConvFADC2Phe();
      const Float_t qe   = qepix.GetQECascadesCombined();

      mean  += conv/qe;
      mean2 += conv*conv/qe/qe;
      nr    ++;
      
    }

  TArrayF *arr = new TArrayF(2);
  arr->AddAt(nr   ? mean/nr : -1.,0);
  arr->AddAt(nr>1 ? TMath::Sqrt((mean2 - mean*mean/nr)/(nr-1)) : 0. ,1);

  return arr;
}

// --------------------------------------------------------------------------
//
// Calculates the average conversion factor FADC counts to photons for camera sectors. 
// The geometry container is used to get the necessary
// geometry information (area index). The MCalibrationQECam container is necessary for 
// the quantum efficiency information. 
// If the bad pixel container is given all pixels which have the flag 'kUnsuitableRun' are ignored
// in the calculation of the size average.
//
// Returns a TArrayF of dimension 2: 
// arr[0]: averaged conversion factors (default: -1.)
// arr[1]: Error (rms) of averaged conversion factors (default: 0.)
//
// ATTENTION: THE USER HAS TO DELETE THE RETURNED TARRAYF ACCORDINGLY
//
TArrayF *MCalibrationChargeCam::GetAveragedConvFADC2PhotPerSector( const MGeomCam &geom, const MCalibrationQECam &qecam,
                                             const UInt_t sec, MBadPixelsCam *bad)
{
  const Int_t np = GetSize();

  Double_t mean  = 0.;
  Double_t mean2 = 0.;
  Int_t    nr    = 0;

  for (int i=0; i<np; i++)
    {
      if (bad && (*bad)[i].IsUnsuitable(MBadPixelsPix::kUnsuitableRun))
        continue; 
      
      const UInt_t sector = geom[i].GetSector();
      
      if (sector != sec)
        continue;

      const MCalibrationQEPix &qepix = (MCalibrationQEPix&)qecam[i];
      if (!qepix.IsAverageQECombinedAvailable())
        continue;
      
      const MCalibrationChargePix &pix = (MCalibrationChargePix&)(*this)[i];
      const Float_t conv = pix.GetMeanConvFADC2Phe();
      const Float_t qe = qepix.GetQECascadesCombined();

      mean  += conv/qe;
      mean2 += conv*conv/qe/qe;
      nr    ++;
      
    }

  TArrayF *arr = new TArrayF(2);
  arr->AddAt(nr   ? mean/nr : -1.,0);
  arr->AddAt(nr>1 ? TMath::Sqrt((mean2 - mean*mean/nr)/(nr-1)) : 0. ,1);

  return arr;
}  

// --------------------------------------------------------------------------
//
// Calculates the average mean arrival times for pixel sizes. 
// The geometry container is used to get the necessary
// geometry information (area index). 
// If the bad pixel container is given all pixels which have the flag 'kUnsuitableRun' are ignored
// in the calculation of the size average.
//
// Returns a TArrayF of dimension 2: 
// arr[0]: averaged mean arrival times (default: -1.)
// arr[1]: Error (rms) of averaged mean arrival times (default: 0.)
//
// ATTENTION: THE USER HAS TO DELETE THE RETURNED TARRAYF ACCORDINGLY
//
TArrayF *MCalibrationChargeCam::GetAveragedArrivalTimeMeanPerArea  ( const MGeomCam &geom, 
                                           const UInt_t ai,  MBadPixelsCam *bad)
{

  const Int_t np = GetSize();

  Double_t mean  = 0.;
  Double_t mean2 = 0.;
  Int_t    nr    = 0;

  for (int i=0; i<np; i++)
    {
      if (bad && (*bad)[i].IsUnsuitable(MBadPixelsPix::kUnsuitableRun))
        continue; 
      
      const UInt_t aidx = geom[i].GetAidx();
      
      if (ai != aidx)
        continue;

      const MCalibrationChargePix &pix = (MCalibrationChargePix&)(*this)[i];
      const Float_t time = pix.GetAbsTimeMean();

      mean  += time   ;
      mean2 += time*time;
      nr    ++;
      
    }

  TArrayF *arr = new TArrayF(2);
  arr->AddAt(nr   ? mean/nr : -1.,0);
  arr->AddAt(nr>1 ? TMath::Sqrt((mean2 - mean*mean/nr)/(nr-1)) : 0. ,1);

  return arr;
}

// --------------------------------------------------------------------------
//
// Calculates the average mean arrival times for camera sectors. 
// The geometry container is used to get the necessary
// geometry information (area index).
// If the bad pixel container is given all pixels which have the flag 'kUnsuitableRun' are ignored
// in the calculation of the size average.
//
// Returns a TArrayF of dimension 2: 
// arr[0]: averaged mean arrival times (default: -1.)
// arr[1]: Error (rms) of averaged mean arrival times (default: 0.)
//
// ATTENTION: THE USER HAS TO DELETE THE RETURNED TARRAYF ACCORDINGLY
//
TArrayF *MCalibrationChargeCam::GetAveragedArrivalTimeMeanPerSector( const MGeomCam &geom, 
                                             const UInt_t sec, MBadPixelsCam *bad)
{
  const Int_t np = GetSize();

  Double_t mean  = 0.;
  Double_t mean2 = 0.;
  Int_t    nr    = 0;

  for (int i=0; i<np; i++)
    {
      if (bad && (*bad)[i].IsUnsuitable(MBadPixelsPix::kUnsuitableRun))
        continue; 
      
      const UInt_t sector = geom[i].GetSector();
      
      if (sector != sec)
        continue;

      const MCalibrationChargePix &pix = (MCalibrationChargePix&)(*this)[i];
      const Float_t time = pix.GetAbsTimeMean();

      mean  += time;
      mean2 += time*time;
      nr    ++;
      
    }

  TArrayF *arr = new TArrayF(2);
  arr->AddAt(nr   ? mean/nr : -1.,0);
  arr->AddAt(nr>1 ? TMath::Sqrt((mean2 - mean*mean/nr)/(nr-1)) : 0. ,1);

  return arr;
}  

// --------------------------------------------------------------------------
//
// Calculates the average arrival time RMSs for pixel sizes. 
// The geometry container is used to get the necessary
// geometry information (area index). 
// If the bad pixel container is given all pixels which have the flag 'kUnsuitableRun' are ignored
// in the calculation of the size average.
//
// Returns a TArrayF of dimension 2: 
// arr[0]: averaged arrival time RMSs (default: -1.)
// arr[1]: Error (rms) of averaged arrival time RMSs (default: 0.)
//
// ATTENTION: THE USER HAS TO DELETE THE RETURNED TARRAYF ACCORDINGLY
//
TArrayF *MCalibrationChargeCam::GetAveragedArrivalTimeRmsPerArea  ( const MGeomCam &geom, 
                                           const UInt_t ai,  MBadPixelsCam *bad)
{

  const Int_t np = GetSize();

  Double_t mean  = 0.;
  Double_t mean2 = 0.;
  Int_t    nr    = 0;

  for (int i=0; i<np; i++)
    {
      if (bad && (*bad)[i].IsUnsuitable(MBadPixelsPix::kUnsuitableRun))
        continue; 
      
      const UInt_t aidx = geom[i].GetAidx();
      
      if (ai != aidx)
        continue;

      const MCalibrationChargePix &pix = (MCalibrationChargePix&)(*this)[i];
      const Float_t rms = pix.GetAbsTimeRms();

      mean  += rms;
      mean2 += rms*rms;
      nr    ++;
      
    }

  TArrayF *arr = new TArrayF(2);
  arr->AddAt(nr   ? mean/nr : -1.,0);
  arr->AddAt(nr>1 ? TMath::Sqrt((mean2 - mean*mean/nr)/(nr-1)) : 0. ,1);

  return arr;
}

// --------------------------------------------------------------------------
//
// Calculates the average arrival time RMSs for camera sectors. 
// The geometry container is used to get the necessary
// geometry information (area index).
// If the bad pixel container is given all pixels which have the flag 'kUnsuitableRun' are ignored
// in the calculation of the size average.
//
// Returns a TArrayF of dimension 2: 
// arr[0]: averaged arrival time RMSs (default: -1.)
// arr[1]: Error (rms) of averaged arrival time RMSs (default: 0.)
//
// ATTENTION: THE USER HAS TO DELETE THE RETURNED TARRAYF ACCORDINGLY
//
TArrayF *MCalibrationChargeCam::GetAveragedArrivalTimeRmsPerSector( const MGeomCam &geom, const UInt_t sec, MBadPixelsCam *bad)
{
  const Int_t np = GetSize();

  Double_t mean  = 0.;
  Double_t mean2 = 0.;
  Int_t    nr    = 0;

  for (int i=0; i<np; i++)
    {
      if (bad && (*bad)[i].IsUnsuitable(MBadPixelsPix::kUnsuitableRun))
        continue; 
      
      const UInt_t sector = geom[i].GetSector();
      
      if (sector != sec)
        continue;

      const MCalibrationChargePix &pix = (MCalibrationChargePix&)(*this)[i];
      const Float_t rms = pix.GetAbsTimeRms();

      mean  += rms;
      mean2 += rms*rms;
      nr    ++;
    }

  TArrayF *arr = new TArrayF(2);
  arr->AddAt(nr   ? mean/nr : -1.,0);
  arr->AddAt(nr>1 ? TMath::Sqrt((mean2 - mean*mean/nr)/(nr-1)) : 0. ,1);

  return arr;
}  
