#include <cstdio>

#include "Stntuple/photon/TDiCesCpr.hh"

 ClassImp(TDiCesCpr)

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

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

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


//_____________________________________________________________________________
Float_t TDiCesCpr::AddPhoton(TStnEvent* event, TStnPhoton* pho0, 
			   TStnPhoton* pho1) {
  SetEvent(event);

  fCesFlag[0] = OkForCes(pho0);
  fCesFlag[1] = OkForCes(pho1);
  fCesWt[0] = CesWeight(pho0, fCesProbS[0], fCesProbB[0]);
  fCesWt[1] = CesWeight(pho1, fCesProbS[1], fCesProbB[1]);
  fCesHit[0] = (fCesWt[0]>0.0);
  fCesHit[1] = (fCesWt[1]>0.0);
  fCprFlag[0] = OkForCpr(pho0);
  fCprFlag[1] = OkForCpr(pho1);
  fCprWt[0] = CprWeight(pho0, fCprProbS[0], fCprProbB[0]);
  fCprWt[1] = CprWeight(pho1, fCprProbS[1], fCprProbB[1]);
  fCprHit[0] = (fCprWt[0]<0.0);
  fCprHit[1] = (fCprWt[1]<0.0);

  float ps0,ps1,pb0,pb1;
  bool q0,q1;
  bool OK0 = false, OK1=false;

  if(fCesFlag[0]==0) {
    ps0 = fCesProbS[0];
    pb0 = fCesProbB[0];
    q0  = fCesHit[0];
    OK0 = true;
  } else if(fCprFlag[0]==0) {
    ps0 = fCprProbS[0];
    pb0 = fCprProbB[0];
    q0  = fCprHit[0];
    OK0 = true;
  }

  if(fCesFlag[1]==0) {
    ps1 = fCesProbS[1];
    pb1 = fCesProbB[1];
    q1  = fCesHit[1];
    OK1 = true;
  } else if(fCprFlag[1]==0) {
    ps1 = fCprProbS[1];
    pb1 = fCprProbB[1];
    q1  = fCprHit[1];
    OK1 = true;
  }

  if(OK0 && OK1) {
    fWeight = fMat.AddEvent(ps0, pb0, ps1, pb1, q0, q1);
    fNUsed++;
  } else {
    fWeight = 0.0;
    fNCut++;
    fMat.AddNonFiducialEvent();
  }
  fNEvent++;

  return 0;
}


//_____________________________________________________________________________
Float_t TDiCesCpr::Fraction(Int_t ind) { 
  return fMat.GetEbe(ind)/fMat.GetEbeTotal();
}
Float_t TDiCesCpr::FractionErr(Int_t ind){ 
  return fMat.GetEbeSig(ind)/fMat.GetEbeTotal();
}
Float_t TDiCesCpr::Events(Int_t ind, Bool_t scale) { 
  return fMat.GetEbe(ind)*(scale?float(NEvent())/float(NUsed()):1.0);
}
Float_t TDiCesCpr::EventsErr(Int_t ind, Bool_t scale){ 
  return fMat.GetEbeSig(ind)*(scale?float(NEvent())/float(NUsed()):1.0);
}

//_____________________________________________________________________________
void TDiCesCpr::Print() {
  printf("Results for diphoton Ces/CPR method, %s\n",GetName());
  TCesCprBase::Print();

  printf(" Observed   Signal                          Events\n");
  printf(" Events     Fraction        Events          Scaled\n");
  for(int i=0; i<4; i++) {
    printf("%4d   %6.3f +-%6.3f    %8.3f +-%6.3f   %8.3f +-%6.3f\n",
	   fMat.GetObs(i),
	   Fraction(i),FractionErr(i),
	   Events(i,false),EventsErr(i,false),
	   Events(i),EventsErr(i));
  }
}

//_____________________________________________________________________________
void TDiCesCpr::Clear() {
  TCesCprBase::Clear();
  fMat.Clear();
}