/* ======================================================================== *\
!
! *
! * 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   02/2004 <mailto:markus@ifae.es>
!
!   Copyright: MAGIC Software Development, 2000-2004
!
!
\* ======================================================================== */

/////////////////////////////////////////////////////////////////////////////
//                                                                         //
// MCalibrationChargePINDiode                                                            //
//                                                                         //
// This is the storage container to hold informations about the pedestal   //
// (offset) value of one Pixel (PMT).                                      //
//                                                                         //
/////////////////////////////////////////////////////////////////////////////
#include "MCalibrationChargePINDiode.h"

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

ClassImp(MCalibrationChargePINDiode);

using namespace std;

const Float_t MCalibrationChargePINDiode::fgChargeLimit        = 3.;
const Float_t MCalibrationChargePINDiode::fgChargeErrLimit     = 0.;    
const Float_t MCalibrationChargePINDiode::fgChargeRelErrLimit  = 1.; 

const Float_t MCalibrationChargePINDiode::fgConversionChargePhotons    = -1.; 
const Float_t MCalibrationChargePINDiode::fgConversionChargePhotonsErr = -1.; 
//
// Area of Inner Pixel w.r.t. PIN Diode (which is 1 cm2)
//
// Distance of PIN Diode to pulser D1:   1.5  +- 0.3 m
// Distance of Inner Pixel to pulser D2: 18.0 +- 0.5 m
//
//
//                 D1*D1
// conversion C = ------ = 0.0069
//                 D2*D2
//
// Delta C / C  = 2 * Sqrt( (Delta D1/D1)2 + (Delta D2/D2)2 )
// Delta C / C  = 0.4
// 
// C = 0.007 +- 0.003
//
const Float_t MCalibrationChargePINDiode::gkFluxCameravsPINDiode      = 0.007;
const Float_t MCalibrationChargePINDiode::gkFluxCameravsPINDiodeErr   = 0.003;
//
// Average QE of the PIN Diode
//
const Float_t MCalibrationChargePINDiode::gkPINDiodeQEGreen    = -1.0;
const Float_t MCalibrationChargePINDiode::gkPINDiodeQEBlue     = -1.0;
const Float_t MCalibrationChargePINDiode::gkPINDiodeQEUV       = -1.0;
const Float_t MCalibrationChargePINDiode::gkPINDiodeQECT1      = -1.0;
//
// Average QE of the PIN Diode
//
const Float_t MCalibrationChargePINDiode::gkPINDiodeQEGreenErr = -1.0;
const Float_t MCalibrationChargePINDiode::gkPINDiodeQEBlueErr  = -1.0;
const Float_t MCalibrationChargePINDiode::gkPINDiodeQEUVErr    = -1.0;
const Float_t MCalibrationChargePINDiode::gkPINDiodeQECT1Err   = -1.0;

const Float_t MCalibrationChargePINDiode::gkPINDiodeArea       = 100;
// --------------------------------------------------------------------------
//
// Default Constructor. 
//
MCalibrationChargePINDiode::MCalibrationChargePINDiode(const char *name, const char *title)
    : fFlags(0)
{

  fName  = name  ? name  : "MCalibrationChargePINDiode";
  fTitle = title ? title : "Container of the fit results of MHCalibrationChargePINDiode";

  Clear();

  SetChargeLimit();
  SetChargeErrLimit();
  SetChargeRelErrLimit();

  SetConversionChargePhotons();
  SetConversionChargePhotonsErr();

  SetExcluded(kFALSE);
  SetExcludeQualityCheck(kFALSE);
}

// ------------------------------------------------------------------------
//
// Invalidate values
//
void MCalibrationChargePINDiode::Clear(Option_t *o)
{

    SetChargeFitValid     ( kFALSE );
    SetTimeFitValid       ( kFALSE );
    SetMeanTimeInFirstBin ( kFALSE );
    SetMeanTimeInLastBin  ( kFALSE );
    
  fMeanCharge                       =  -1.;
  fMeanChargeErr                    =  -1.;
  fSigmaCharge                      =  -1.;
  fSigmaChargeErr                   =  -1.;
  fChargeProb                       =  -1.;
  fPed                              =  -1.;
  fPedRms                           =  -1.;

  fRmsChargeMean                    =  -1.;
  fRmsChargeMeanErr                 =  -1.;
  fRmsChargeSigma                   =  -1.;  
  fRmsChargeSigmaErr                =  -1.;

  fAbsTimeMean                      =  -1.;
  fAbsTimeRms                       =  -1.;

  fTimeLowerEdge                    =  -1.;
  fTimeUpperEdge                    =  -1.;

  fConvertedPhotons                 =  -1.;
  fConvertedPhotonsErr              =  -1.;

  fMeanFluxOutsidePlexiglass        =  -1.;  
  fMeanFluxErrOutsidePlexiglass     =  -1.;

}


// --------------------------------------------------------------------------
//
// Set the pedestals from outside
//
void MCalibrationChargePINDiode::SetPedestal(Float_t ped, Float_t pedrms)
{

  fPed    = ped;    
  fPedRms = pedrms;
  
}

// --------------------------------------------------------------------------
//
// Set the Excluded Bit from outside 
//
void MCalibrationChargePINDiode::SetExcluded(Bool_t b )
{ 
  b ?  SETBIT(fFlags, kExcluded) : CLRBIT(fFlags, kExcluded); 
}


// --------------------------------------------------------------------------
//
// Set the Excluded Bit from outside 
//
void MCalibrationChargePINDiode::SetExcludeQualityCheck(Bool_t b )
{ 
  b ?  SETBIT(fFlags, kExcludeQualityCheck) : CLRBIT(fFlags, kExcludeQualityCheck); 
}

// --------------------------------------------------------------------------
//
// Set the Excluded Bit from outside 
//
void MCalibrationChargePINDiode::SetChargeFitValid(Bool_t b )    
{ 
  b ?  SETBIT(fFlags, kChargeFitValid) : CLRBIT(fFlags, kChargeFitValid); 
}

// --------------------------------------------------------------------------
//
// Set the Excluded Bit from outside 
//
void MCalibrationChargePINDiode::SetTimeFitValid(Bool_t b )    
{ 
  b ?  SETBIT(fFlags, kTimeFitValid) : CLRBIT(fFlags, kTimeFitValid); 
}

void MCalibrationChargePINDiode::SetFluxOutsidePlexiglassAvailable (const Bool_t b)
{
  b ?  SETBIT(fFlags, kFluxOutsidePlexiglassAvailable) : CLRBIT(fFlags, kFluxOutsidePlexiglassAvailable); 
}

void MCalibrationChargePINDiode::SetMeanTimeInFirstBin(const Bool_t b)
{
  b ? SETBIT(fFlags,kMeanTimeInFirstBin) : CLRBIT(fFlags,kMeanTimeInFirstBin);
}

void MCalibrationChargePINDiode::SetMeanTimeInLastBin(const Bool_t b)
{
  b ? SETBIT(fFlags,kMeanTimeInLastBin) : CLRBIT(fFlags,kMeanTimeInLastBin);
}

Bool_t MCalibrationChargePINDiode::IsMeanTimeInFirstBin() const
{
  return TESTBIT(fFlags,kMeanTimeInFirstBin);
}

Bool_t MCalibrationChargePINDiode::IsMeanTimeInLastBin() const
{
  return TESTBIT(fFlags,kMeanTimeInLastBin);
}

Bool_t MCalibrationChargePINDiode::IsExcluded()       const
{ 
   return TESTBIT(fFlags,kExcluded);  
}

Bool_t MCalibrationChargePINDiode::IsChargeFitValid() const 
{
  return TESTBIT(fFlags, kChargeFitValid);  
}

Bool_t MCalibrationChargePINDiode::IsTimeFitValid()   const 
{
  return TESTBIT(fFlags, kTimeFitValid);  
}

//
// The check return kTRUE if:
//
// 1) PINDiode has a fitted charge greater than fChargeLimit*PedRMS
// 2) PINDiode has a fit error greater than fChargeErrLimit
// 3) PINDiode has a fitted charge greater its fChargeRelErrLimit times its charge error
// 4) PINDiode has a charge sigma bigger than its Pedestal RMS
// 
Bool_t MCalibrationChargePINDiode::CheckChargeFitValidity()
{


  if (fMeanCharge < fChargeLimit*GetPedRms())
    {
      *fLog << warn << "WARNING: Fitted Charge is smaller than "
            << fChargeLimit << " Pedestal RMS in PINDiode " << endl;
      return kFALSE;
    }
  
  if (fMeanChargeErr < fChargeErrLimit) 
    {
      *fLog << warn << "WARNING: Error of Fitted Charge is smaller than "
            << fChargeErrLimit << " in PINDiode " << endl;
      return kFALSE;
    }
      
  if (fMeanCharge < fChargeRelErrLimit*fMeanChargeErr) 
    {
      *fLog << warn << "WARNING: Fitted Charge is smaller than "
            << fChargeRelErrLimit << "* its error in PINDiode " << endl;
      return kFALSE;
    }
      
  if (fSigmaCharge < GetPedRms())
    {
      *fLog << warn << "WARNING: Sigma of Fitted Charge smaller than Pedestal RMS in PINDiode " << endl;
      return kFALSE;
    }
  return kTRUE;
}

//
// The check returns kTRUE if:
//
// The mean arrival time is at least 1.0 slices from the used edge slices 
//
Bool_t MCalibrationChargePINDiode::CheckTimeFitValidity()
{

  if ( fAbsTimeMean < fTimeLowerEdge+1.)
    {
      *fLog << warn << "WARNING: Mean ArrivalTime in first extraction bin of the PINDiode " << endl;
      SetMeanTimeInFirstBin();
      return kFALSE;
    }

  if ( fAbsTimeMean > fTimeUpperEdge-1.)
    {
      *fLog << warn << "WARNING: Mean ArrivalTime in last extraction bin of the PINDiode " << endl;
      SetMeanTimeInLastBin();
      return kFALSE;
    }

  return kTRUE;
}

Bool_t MCalibrationChargePINDiode::CalcFluxOutsidePlexiglass()
{

  if (IsChargeFitValid())
    return kFALSE;
  
  // Start calculation of number of photons per mm^2 on the camera
  fMeanFluxOutsidePlexiglass  = fConvertedPhotons * gkPINDiodeArea;
  // Correct for the distance between camera and PIN Diode and for different areas.
  fMeanFluxOutsidePlexiglass *= gkFluxCameravsPINDiode;

  // Start calculation of number of photons relative Variance (!!)
  fMeanFluxErrOutsidePlexiglass  = fConvertedPhotonsErr * fConvertedPhotonsErr 
                                 / fConvertedPhotons    / fConvertedPhotons    ;
  fMeanFluxErrOutsidePlexiglass += gkFluxCameravsPINDiodeErr*gkFluxCameravsPINDiodeErr
                                 / gkFluxCameravsPINDiode/gkFluxCameravsPINDiode;
  
  switch (fColor)
    {
    case kGREEN:
      fMeanFluxOutsidePlexiglass    /= gkPINDiodeQEGreen;
      fMeanFluxErrOutsidePlexiglass += gkPINDiodeQEGreenErr*gkPINDiodeQEGreenErr
                                     / gkPINDiodeQEGreen/gkPINDiodeQEGreen;
      break;
    case kBLUE:
      fMeanFluxOutsidePlexiglass    /= gkPINDiodeQEBlue;
      fMeanFluxErrOutsidePlexiglass += gkPINDiodeQEBlueErr*gkPINDiodeQEBlueErr
                                     / gkPINDiodeQEBlue/gkPINDiodeQEBlue;
      break; 
    case kUV:
      fMeanFluxOutsidePlexiglass    /= gkPINDiodeQEUV;
      fMeanFluxErrOutsidePlexiglass += gkPINDiodeQEUVErr*gkPINDiodeQEUVErr
                                     / gkPINDiodeQEUV/gkPINDiodeQEUV;
      break;
    case kCT1:
    default:
      fMeanFluxOutsidePlexiglass    /= gkPINDiodeQECT1;
      fMeanFluxErrOutsidePlexiglass += gkPINDiodeQECT1Err*gkPINDiodeQECT1Err
                                     / gkPINDiodeQECT1/gkPINDiodeQECT1;
      break;
    }


  *fLog << inf << endl;
  *fLog << inf << " Mean Photon flux [ph/mm^2] outside Plexiglass: " 
        << fMeanFluxOutsidePlexiglass << endl;

  if (fMeanFluxOutsidePlexiglass < 0.)
      return kFALSE;

  if (fMeanFluxErrOutsidePlexiglass < 0.)
      return kFALSE;

  SetFluxOutsidePlexiglassAvailable();  

  // Finish calculation of errors -> convert from relative variance to absolute error
  fMeanFluxErrOutsidePlexiglass = TMath::Sqrt(fMeanFluxErrOutsidePlexiglass);
  fMeanFluxErrOutsidePlexiglass *= fMeanFluxOutsidePlexiglass;

  *fLog << inf << " Error on Photon flux [ph/mm^2] outside Plexiglass: " 
        << fMeanFluxErrOutsidePlexiglass << endl;
  *fLog << inf << endl;

  return kTRUE;
}



