#include <cstdio>

#include "Stntuple/photon/TCesCpr.hh"

 ClassImp(TCesCpr)

//_____________________________________________________________________________
TCesCpr::TCesCpr()
{
  Clear();
}

//_____________________________________________________________________________
TCesCpr::TCesCpr(const char* name)
{
  Clear();
  SetName(name);
}

//_____________________________________________________________________________
TCesCpr::~TCesCpr() {
}


//_____________________________________________________________________________
Float_t TCesCpr::AddPhoton(TStnEvent* event, TStnPhoton* pho) {

  SetEvent(event);

  fCesFlag = OkForCes(pho);
  fCprFlag = OkForCpr(pho);
  fCesWt = CesWeight(pho, fCesProbS, fCesProbB);
  fCprWt = CprWeight(pho, fCprProbS, fCprProbB);
  fCesHit = (fCesWt>0.0);
  fCprHit = (fCprWt<0.0);

  float fWeight = 0.0;
  if(fCesFlag==0) {
    fWeight = fCesWt;
  } else if (fCprFlag==0) {
    fWeight = fCprWt;
  }

  if(fabs(fWeight)>1e-6){
    fSum += fWeight;
    fSum2 += pow(fWeight,2);
    fSum2Bg += pow(1.0-fWeight,2);
    fNUsed++;
  } else {
    fNCut++;
  }
  fNEvent++;
  return fWeight;
}

//_____________________________________________________________________________
Float_t TCesCpr::SignalFraction() { return fSum/NUsed();}
Float_t TCesCpr::SignalFractionErr(){ return sqrt(fSum2)/NUsed();}
Float_t TCesCpr::BGFraction(){ return 1.0-SignalFraction();}
Float_t TCesCpr::BGFractionErr(){ return sqrt(fSum2Bg)/NUsed();}
Float_t TCesCpr::Signal(Bool_t scale){
  return fSum*(scale?float(NEvent())/float(NUsed()):1.0);
}
Float_t TCesCpr::SignalErr(Bool_t scale){
  return sqrt(fSum2)*(scale?float(NEvent())/float(NUsed()):1.0);
}
Float_t TCesCpr::BG(Bool_t scale){
  return (fNUsed - fSum)*(scale?float(NEvent())/float(NUsed()):1.0);
}
Float_t TCesCpr::BGErr(Bool_t scale){
  return sqrt(fSum2Bg)*(scale?float(NEvent())/float(NUsed()):1.0);
}

//_____________________________________________________________________________
void TCesCpr::Print() {
  printf("Results for Ces/CPR method, %s\n",GetName());
  TCesCprBase::Print();
  printf("Signal fraction = %6.3f +-%6.3f \n",SignalFraction(),
	 SignalFractionErr());
  printf("BG fraction     = %6.3f +-%6.3f \n",BGFraction(),
	 BGFractionErr());
  printf("Result for useable photons:\n");
  printf("Signal     = %7.1f +-%7.1f \n",Signal(false),
	 SignalErr(false));
  printf("Background = %7.1f +-%7.1f \n",BG(false),
	 BGErr(false));
  printf("Scaled to total numberof phtons:\n");
  printf("Signal     = %7.1f +-%7.1f \n",Signal(),
	 SignalErr());
  printf("Background = %7.1f +-%7.1f \n",BG(),
	 BGErr());
}

//_____________________________________________________________________________
void TCesCpr::Clear() {
  TCesCprBase::Clear();
  fSum = 0.0;
  fSum2 = 0.0;
  fSum2Bg = 0.0;
}