source: trunk/MagicSoft/Mars/mcalib/MCalibrationPINDiode.cc@ 3007

/* ======================================================================== *\
!
! *
! * 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-2001
!
!
\* ======================================================================== */

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

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

ClassImp(MCalibrationPINDiode);

using namespace std;
// --------------------------------------------------------------------------
//
// Default Constructor. 
//
MCalibrationPINDiode::MCalibrationPINDiode(const char *name, const char *title)
    : fChargeLimit(3.),
      fChargeErrLimit(0.),
      fChargeRelErrLimit(1.),
      fFlags(0)
{

  fName  = name  ? name  : "MCalibrationPINDiode";
  fTitle = title ? title : "Container of the MHCalibrationPINDiode and the fit results";

  fHist = new MHCalibrationPINDiode();

  if (!fHist)
    *fLog << warn << dbginf << " Could not create MHCalibrationPINDiode " << endl;
 
  Clear();

}

MCalibrationPINDiode::~MCalibrationPINDiode() 
{
  delete fHist;
}

// ------------------------------------------------------------------------
//
// Invalidate values
//
void MCalibrationPINDiode::Clear(Option_t *o)
{
  fHist->Reset();

  CLRBIT(fFlags, kExcluded);
  CLRBIT(fFlags, kChargeFitValid);
  CLRBIT(fFlags, kChargeFitValid);
  CLRBIT(fFlags, kFitted);

  fCharge                           =  -1.;
  fErrCharge                        =  -1.;
  fSigmaCharge                      =  -1.;
  fErrSigmaCharge                   =  -1.;
  fRSigmaSquare                     =  -1.;
  fChargeProb                       =  -1.;
  fPed                              =  -1.;
  fPedRms                           =  -1.;
  fTime                             =  -1.;
  fSigmaTime                        =  -1.;
  fTimeChiSquare                    =  -1.;

}


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

  fPed    = ped;    
  fPedRms = pedrms;
  
}

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


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

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

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

// --------------------------------------------------------------------------
//
// Set the Excluded Bit from outside 
//
void MCalibrationPINDiode::SetFitted(Bool_t b )
{ 
  b ?  SETBIT(fFlags, kFitted) : CLRBIT(fFlags, kFitted); 
}

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

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

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

Bool_t MCalibrationPINDiode::IsFitted() const 
{ 
  return TESTBIT(fFlags, kFitted);    
}

Bool_t MCalibrationPINDiode::FitCharge() 
{

  //
  // 1) Return if the charge distribution is already succesfully fitted  
  //    or if the histogram is empty
  //
  if (fHist->IsChargeFitOK() || fHist->IsEmpty())
    return kTRUE;

  //
  // 4) Fit the Lo Gain histograms with a Gaussian
  //
  if(fHist->FitCharge())
    {
      SETBIT(fFlags,kFitted);
    }
  else
    {
      *fLog << warn << "WARNING: Could not fit charges of PINDiode " << endl;
      //          
      // 5) In case of failure print out the fit results
      //
      //          fHist->PrintChargeFitResult();
      CLRBIT(fFlags,kFitted);
    }

  //
  // 6) Retrieve the results and store them in this class
  //
  fCharge         = fHist->GetChargeMean();
  fErrCharge      = fHist->GetChargeMeanErr(); 
  fSigmaCharge    = fHist->GetChargeSigma();
  fErrSigmaCharge = fHist->GetChargeSigmaErr();
  fChargeProb     = fHist->GetChargeProb();

  if (CheckChargeFitValidity())
    SETBIT(fFlags,kChargeFitValid);
  else
    {
      CLRBIT(fFlags,kChargeFitValid);
      return kFALSE;
    }

  return kTRUE;
  
}

//
// The check return kTRUE if:
//
// 1) PINDiode has a fitted charge greater than 5*PedRMS
// 2) PINDiode has a fit error greater than 0.
// 3) PINDiode has a fitted charge greater its charge error
// 4) PINDiode has a fit Probability greater than 0.0001 
// 5) PINDiode has a charge sigma bigger than its Pedestal RMS
// 
Bool_t MCalibrationPINDiode::CheckChargeFitValidity()
{

  if (TESTBIT(fFlags,kExcludeQualityCheck))
    return kTRUE;

  if (fCharge < fChargeLimit*GetPedRms())
    {
      *fLog << warn << "WARNING: Fitted Charge is smaller than "
            << fChargeLimit << " Pedestal RMS in PINDiode " << endl;
      return kFALSE;
    }
  
  if (fErrCharge < fChargeErrLimit) 
    {
      *fLog << warn << "WARNING: Error of Fitted Charge is smaller than "
            << fChargeErrLimit << " in PINDiode " << endl;
      return kFALSE;
    }
      
  if (fCharge < fChargeRelErrLimit*fErrCharge) 
    {
      *fLog << warn << "WARNING: Fitted Charge is smaller than "
            << fChargeRelErrLimit << "* its error in PINDiode " << endl;
      return kFALSE;
    }
      
  if (!fHist->IsChargeFitOK()) 
    {
      *fLog << warn << "WARNING: Probability of Fitted Charge too low 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;
}


// --------------------------------------------------------------------------
//
// 1) Fit the arrival times
// 2) Retrieve the results
// 3) Note that because of the low number of bins, the NDf is sometimes 0, so 
//    Root does not give a reasonable Probability, the Chisquare is more significant
//
Bool_t MCalibrationPINDiode::FitTime() 
{

  if(!fHist->FitTime())
    {
      *fLog << warn << "WARNING: Could not fit Hi Gain times of PIN Diode" << endl;
      fHist->PrintTimeFitResult();
      return kFALSE;
    }
  
  fTime          = fHist->GetRelTimeMean();
  fSigmaTime     = fHist->GetRelTimeSigma();
  fTimeProb      = fHist->GetRelTimeProb();
  
  if (CheckTimeFitValidity())
    SETBIT(fFlags,kTimeFitValid);
  else
    CLRBIT(fFlags,kTimeFitValid);

  return kTRUE;
}

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

  if (TESTBIT(fFlags,kExcludeQualityCheck))
    return kTRUE;

  return kTRUE;
}