source: trunk/MagicSoft/Mars/manalysis/MCalibrationCam.cc@ 2541

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

/////////////////////////////////////////////////////////////////////////////
//                                                                         //
// MCalibrationCam                                                         //
//                                                                         //
// Hold the Calibration information for all pixels in the camera           //
//                                                                         //
// The Classes MCalibrationPix's are stored in a TClonesArray              //
// The Class MCalibrationBlindPix and the MCalibrationPINDiode             // 
//           are stored in separate pointers                               //
//                                                                         //
/////////////////////////////////////////////////////////////////////////////
#include "MCalibrationCam.h"
#include "MCalibrationPix.h"
#include "MCalibrationBlindPix.h"
#include "MCalibrationConfig.h"

#include <TClonesArray.h>

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

#include "MGeomCam.h"

ClassImp(MCalibrationCam);

using namespace std;
// --------------------------------------------------------------------------
//
// Default constructor. 
// Creates a TClonesArray of MHCalibrationPix objects, initialized to 0 entries
// Creates an MCalibrationBlindPix and an MCalibrationPINDiode
//
MCalibrationCam::MCalibrationCam(const char *name, const char *title)
    : fMeanNrPhotAvailable(kFALSE)
{
    fName  = name  ? name  : "MCalibrationCam";
    fTitle = title ? title : "Storage container for the Calibration Information in the camera";

    fPixels     = new TClonesArray("MCalibrationPix",0);
    fBlindPixel = new MCalibrationBlindPix();
    fPINDiode   = new MCalibrationPINDiode();
}

// --------------------------------------------------------------------------
//
// Delete the TClonesArray or MCalibrationPix's
// Delete the MCalibrationPINDiode and the MCalibrationBlindPix
//
MCalibrationCam::~MCalibrationCam()
{

  //
  // delete fPixels should delete all Objects stored inside
  // 
  delete fPixels;
  delete fBlindPixel;
  delete fPINDiode;
}

// --------------------------------------------------------------------------
//
// This function return the size of the FILLED MCalibrationCam
// It is NOT the size of the array fPixels !!!
// Note that with every call to AddPixels, GetSize() might change
//
Int_t MCalibrationCam::GetSize() const
{

  //
  // Here it is important to use GetEntriesFast, 
  // We get the array index of the last "filled" fPixel entry
  // (Not the number of "filled" fPixels!!)
  //
  return fPixels->GetEntriesFast();
}

// --------------------------------------------------------------------------
//
// Get i-th pixel (pixel number)
//
MCalibrationPix &MCalibrationCam::operator[](Int_t i)
{
    return *static_cast<MCalibrationPix*>(fPixels->UncheckedAt(i));
}

// --------------------------------------------------------------------------
//
// Get i-th pixel (pixel number)
//
MCalibrationPix &MCalibrationCam::operator[](Int_t i) const
{
    return *static_cast<MCalibrationPix*>(fPixels->UncheckedAt(i));
}

Bool_t MCalibrationCam::AddPixel(Int_t idx) 
{

  //
  // Check bounds first
  //
  if (!CheckBounds(idx))
    InitSize(idx);

  //
  // in case, new cannot allocate memory, 
  // return FALSE
  //
  return (new ((*fPixels)[idx]) MCalibrationPix(idx));

}


// --------------------------------------------------------------------------
//
// Return a pointer to the pixel with the requested idx. 
// NULL if it doesn't exist. 
//
MCalibrationPix *MCalibrationCam::GetCalibrationPix(Int_t idx) const
{
   if (idx<0)
     return NULL;

  if (!CheckBounds(idx))
    return NULL;
  
//    return (*this)[idx];
  return (MCalibrationPix*)fPixels->At(idx);
}



// --------------------------------------------------------------------------
//
// Check if position i is inside bounds
//
Bool_t MCalibrationCam::CheckBounds(Int_t i) const 
{
    return i < fPixels->GetEntriesFast();
} 

Bool_t MCalibrationCam::IsPixelUsed(Int_t idx) const 
{
  if (!CheckBounds(idx))
    return kFALSE;

  return (&(*this)[idx]);
}

Bool_t MCalibrationCam::IsPixelFitted(Int_t idx) const 
{
  if (!IsPixelUsed(idx))
    return kFALSE;

  return ((*this)[idx].GetRQ() > 0. && (*this)[idx].GetErrRQ() > 0.);
}



void MCalibrationCam::Clear(Option_t *o)
{
    fPixels->ForEach(TObject, Clear)();
}

//
// Perform the fits of the charges 
// If i=-1, then all pixels will be fitted
// Otherwise only the one with index i
// 
// The number of succesful fits is returned
//
UShort_t MCalibrationCam::FitQ(Int_t i)
{

  UShort_t nsuccess = 0;

  // invalidate i if it exceeds the number of entries in fPixels
  if (i > fPixels->GetEntriesFast())      
    {
      *fLog << warn << "Tried to fit pixel out of allowed range " << endl;
      return 0;
    }
  
  if (i == -1)    // loop over all events
    {

      Int_t id = 0;
      
      TIter Next(fPixels);
      MCalibrationPix *pix;
      while ((pix=(MCalibrationPix*)Next()))
        {
          if (!IsPixelUsed(id++))
            continue;
          if (pix->FitQ())
            nsuccess++;
        }
    }
  else                  // fit only the pixel with index i
    {
      if (!IsPixelUsed(i))
        return 0;
      if((*this)[i].FitQ())
        nsuccess++;
    }
  
  return nsuccess;
  
}

//
// Perform the fits of the charges of all pixels 
// plus the blind pixel and the PIN Diode
// 
// The number of succesful fits is returned
//
UShort_t MCalibrationCam::FitAllQ()
{

  // FIXME: Once the blind pixel is fully working, 
  //        there must be some penalty in case the
  //        fit does not succeed
  //
  UShort_t nsuccess = 0;

  Int_t id = 0;

  TIter Next(fPixels);
  MCalibrationPix *pix;
  while ((pix=(MCalibrationPix*)Next()))
    {
      if (!IsPixelUsed(id++))
        continue;
      if (pix->FitQ())
        nsuccess++;
    }
  
  if (fBlindPixel->FitQ())
        nsuccess++;

  if (fPINDiode->FitQ())
        nsuccess++;

  return nsuccess;

}



//
// Perform the fits of the arrival times
// If i=-1, then all pixels will be fitted
// Otherwise only the one with index i
// 
// The number of succesful fits is returned
//
UShort_t MCalibrationCam::FitT(Int_t i)
{

  UShort_t nsuccess = 0;

  // invalidate i if it exceeds the number of entries in fPixels
  if (i > fPixels->GetEntriesFast())      
    {
      *fLog << warn << "Tried to fit pixel out of allowed range " << endl;
      return 0;
    }
  
  if (i == -1)                // loop over all events
    {

      Int_t id = 0;
      
      TIter Next(fPixels);
      MCalibrationPix *pix;
      while ((pix=(MCalibrationPix*)Next()))
        {
          if (!IsPixelUsed(id++))
            continue;
          if (pix->FitT())
            nsuccess++;
        }
    }
  else                     // fit only the pixel with index i
    {
      if((*this)[i].FitT())
        nsuccess++;
    }
  
  return nsuccess;

}


//
// Perform the fits of the times of all pixels 
// plus the blind pixel and the PIN Diode
// 
// The number of succesful fits is returned
//
UShort_t MCalibrationCam::FitAllT()
{

  // FIXME: Once the blind pixel is fully working, 
  //        there must be some penalty in case the
  //        fit does not succeed
  //
  UShort_t nsuccess = 0;

  Int_t id = 0;

  TIter Next(fPixels);
  MCalibrationPix *pix;
  while ((pix=(MCalibrationPix*)Next()))
    {
      if (!IsPixelUsed(id++))
        continue;
      if (pix->FitT())
        nsuccess++;
    }
  
  if (fBlindPixel->FitT())
        nsuccess++;

  if (fPINDiode->FitT())
       nsuccess++;

  return nsuccess;

}



void MCalibrationCam::CutEdges()
{

  fBlindPixel->GetHist()->CutAllEdges();

  TIter Next(fPixels);
  MCalibrationPix *pix;
  while ((pix=(MCalibrationPix*)Next()))
    {
      pix->GetHist()->CutAllEdges();
    }

  return;
}

// ---------------------------------------------------------------------
//
// This function simply allocates memory via the ROOT command:
// (TObject**) TStorage::ReAlloc(fCont, newSize * sizeof(TObject*),
//                                       fSize * sizeof(TObject*));
// newSize corresponds to size in our case
// fSize is the old size (in most cases: 0)
//
void MCalibrationCam::InitSize(Int_t size)
{
  
  //
  // check if we have already initialized to size
  //
  if (CheckBounds(size))
    return;
  
  //
  // It is important to use Expand and NOT ExpandCreate, because 
  // we want to keep all pixel not used with a NULL pointer.
  //
  fPixels->Expand(size);
  //
  // Set all entries to the null pointer.  
  // Later we fill the array per pixId with the contruction: new (fPixels[i]) MCalibrationPix(pixid)
  // To check, if pixels is already filled, we test the NULL pointer (see IsPixelUsed)
  //
  for (Int_t i=0; i< size; i++)
    {
      MCalibrationPix *pix = &(*this)[i];
      pix = NULL;
    }
}
  
void MCalibrationCam::Print(Option_t *o) const
{
    *fLog << all << GetDescriptor() << ":" << endl;
    int id = 0;

    TIter Next(fPixels);
    MCalibrationPix *pix;
    while ((pix=(MCalibrationPix*)Next()))
    {
        if (!IsPixelUsed(id))
          continue;

        *fLog << id << ": ";
        *fLog << pix->GetQ() << " " << pix->GetRQ() << endl;

        id++;
    }
}


Bool_t MCalibrationCam::GetPixelContent(Double_t &val, Int_t idx, const MGeomCam &cam, Int_t type) const
{

  if (!IsPixelFitted(idx))
    return kFALSE;

  switch (type)
    {
    case 0:
       val = (*this)[idx].GetQ();
        break;
    case 1:
       val = (*this)[idx].GetErrQ();
        break;
    case 2:
       val = (*this)[idx].GetSigmaQ();
        break;
    case 3:
       val = (*this)[idx].GetErrSigmaQ();
        break;
    case 4:
       val = (*this)[idx].GetQProb();
        break;
    case 5:
       val = (*this)[idx].GetT();
        break;
    case 6:
       val = (*this)[idx].GetErrT();
        break;
    case 7:
       val = (*this)[idx].GetTProb();
        break;
    case 8:
       val = (*this)[idx].GetPed();
        break;
    case 9:
       val = (*this)[idx].GetPedRms();
        break;
    case 10:
      val = (*this)[idx].GetRQ();
        break;
    case 11:
      val = (*this)[idx].GetErrRQ();
      break;
    }
  return val>=0;
}

void MCalibrationCam::DrawPixelContent(Int_t num) const
{
    *fLog << warn << "MCalibrationCam::DrawPixelContent - not available." << endl;
}

Bool_t MCalibrationCam::CalcNrPhotInnerPixel()
{
  if (!fBlindPixel->IsValid())
    return kFALSE;
  
  const Float_t mean = fBlindPixel->GetLambda();
  //  const Float_t merr = fBlindPixel->GetErrLambda();
  
  switch (fColor)
    {
    case kECGreen:
      fMeanNrPhotInnerPix = mean / 
                           (gkCalibrationBlindPixelQEGreen
                            *TMath::Power(10,gkCalibrationBlindPixelAttGreen)
                            *gkCalibrationBlindPixelArea);
      break;
    case kECBlue:
      fMeanNrPhotInnerPix = mean / 
                           (gkCalibrationBlindPixelQEBlue
                            *TMath::Power(10,gkCalibrationBlindPixelAttBlue)
                            *gkCalibrationBlindPixelArea);
      break;
    case kECUV:
      fMeanNrPhotInnerPix = mean / 
                           (gkCalibrationBlindPixelQEUV
                            *TMath::Power(10,gkCalibrationBlindPixelAttUV)
                            *gkCalibrationBlindPixelArea);
      break;
    }

  fMeanNrPhotAvailable = kTRUE;
  return kTRUE;
}



Bool_t MCalibrationCam::GetConversionFactor(Int_t ipx, Float_t &mean, Float_t &err)
{
  
  if (ipx < 0 || !IsPixelFitted(ipx))
    return kFALSE;

  if (!fMeanNrPhotAvailable)
    if (!CalcNrPhotInnerPixel())
      return kFALSE;

  mean = fMeanNrPhotInnerPix / (*this)[ipx].GetRQ();

  //
  // Not yet ready , sorry 
  //
  err  = 100000000000.;

  return kTRUE;
}