#include <cmath>
#include "Stntuple/alg/TStnPhotonID.hh"
#include "Stntuple/obj/TStnPhoton.hh"

ClassImp(TStnPhotonID)

//_____________________________________________________________________________
TStnPhotonID::TStnPhotonID(const char* Name, const char* Title):TNamed(Name,Title) {
  Reset();
}


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

//_____________________________________________________________________________
int TStnPhotonID::IDWord(TStnPhoton* photon) {
  Int_t      id_word = 0;

  float et = photon->Etc();
  float e = photon->ECorr();
  float ces2 = (photon->CesStripE2()>photon->CesWireE2()? 
		photon->CesStripE2() : photon->CesWireE2() );
  ces2 = ces2*photon->SinTheta();
  float ces2cut = (et<18.0 ? 0.14*et : 2.4 + 0.01*et);

  if(photon->Detector() != 0) id_word |= fCentralBit;
  if(fPrintLevel>0) printf("PhotonID  Cent  =%9i  Cut=%7i  word=0x%08x\n",
			   photon->Detector(),0,id_word&fUseMask);

  if(photon->Detector() != 1) id_word |= fPlugBit;

  if(et<fMinEt) id_word |= fEtBit;
  if(fPrintLevel>0) printf("PhotonID  Et    =%9.3f  Cut=%7.3f  word=0x%08x\n",
			   et,fMinEt,id_word&fUseMask);

  if(fabs(photon->DetEta())>fMaxDetEta) id_word |= fDetEtaBit;
  if(fPrintLevel>0) printf("PhotonID  DetEta=%9.3f  Cut=%7.3f  word=0x%08x\n",
			   photon->DetEta(),fMaxDetEta,id_word&fUseMask);

  if(fabs(photon->Eta())>fMaxEta) id_word |= fEtaBit;
  if(fPrintLevel>0) printf("PhotonID     Eta=%9.3f  Cut=%7.3f  word=0x%08x\n",
			   photon->Eta(),fMaxEta,id_word&fUseMask);

  if(photon->HadEm()>(fMaxHadEm + fHadEmSlope*e)) 
    id_word |= fHadEmBit;
  if(fPrintLevel>0) printf("PhotonID  HadEm =%9.3f  Cut=%7.3f  word=0x%08x\n",
			   photon->HadEm(),(fMaxHadEm + fHadEmSlope*e),
			   id_word&fUseMask);

  bool qces = photon->Chi2()>0.0 && 
    ( (et<fMaxCesAvgChi2Break && photon->Chi2()<fMaxCesAvgChi2     ) ||
      (et>fMaxCesAvgChi2Break && photon->Chi2()<fMaxCesAvgChi2High )    );
  if(!qces) id_word |= fCesAvgChi2Bit;
  if(fPrintLevel>0) printf("PhotonID  Chi2  =%9.3f  Cut=%7.3f  word=0x%08x\n",
			   photon->Chi2(),fMaxCesAvgChi2,id_word&fUseMask);

  if(fabs(photon->XCes())>fMaxXCes) id_word |= fXCesBit;
  if(fPrintLevel>0) printf("PhotonID  CesX  =%9.3f  Cut=%7.3f  word=0x%08x\n",
			   photon->XCes(),fMaxXCes,id_word&fUseMask);

  if(fabs(photon->ZCes())<fMinZCes || fabs(photon->ZCes())>fMaxZCes) 
    id_word |= fZCesBit;
  if(fPrintLevel>0) printf("PhotonID  CesZ  =%9.3f  Cut=%7.3f  word=0x%08x\n",
			   photon->ZCes(),fMinZCes,id_word&fUseMask);


  if(photon->LShr()>fMaxLshr) id_word |= fLshrBit;
  if(fPrintLevel>0) printf("PhotonID  LShr  =%9.3f  Cut=%7.3f  word=0x%08x\n",
			   photon->LShr(),fMaxLshr,id_word&fUseMask);


  float isoCut = ( et<fMaxIso4Break ? 
		   fMaxIso4Low+fMaxIso4SlideLow*et : 
		   fMaxIso4+fMaxIso4Slide*(et-fMaxIso4Break) );
  if(photon->EIso4(2)>isoCut) id_word |= fIsoBit;
  if(fPrintLevel>0) printf("PhotonID  CalIso=%9.3f  Cut=%7.3f  word=0x%08x\n",
			   photon->EIso4(2),isoCut,id_word&fUseMask);

  if(photon->N3d()>fMaxN3D) id_word |= fN3DBit;
  if(fPrintLevel>0) printf("PhotonID  N3d   =%9i  Cut=%7i  word=0x%08x\n",
			   photon->N3d(),fMaxN3D,id_word&fUseMask);

  if(photon->N3d()>0 && photon->Pt()>(fMaxPt+fMaxPtSlide*et)) 
    id_word |= fPtBit;
  if(fPrintLevel>0) printf("PhotonID  Pt    =%9.3f  Cut=%7.3f  word=0x%08x\n",
			   photon->Pt(),(fMaxPt+fMaxPtSlide*et),
			   id_word&fUseMask);

  if(photon->SumPt4()>(fMaxSumPt4+fMaxSumPt4Slide*et)) id_word |= fTIsoBit;
  if(fPrintLevel>0) printf("PhotonID  SumPt =%9.3f  Cut=%7.3f  word=0x%08x\n",
			   photon->SumPt4(),(fMaxSumPt4+fMaxSumPt4Slide*et),
			   id_word&fUseMask);

  if(ces2>ces2cut) id_word |= fCes2Bit;
  if(fPrintLevel>0) printf("PhotonID  Ces2  =%9.3f  Cut=%7.3f  word=0x%08x\n",
			   ces2,ces2cut,id_word&fUseMask);

  if(photon->PesEnergy()<fMinPesE) id_word |= fPesEBit;
  if(fPrintLevel>0) printf("PhotonID  PesE  =%9.3f  Cut=%7.3f  word=0x%08x\n",
			   photon->PesEnergy(),fMinPesE,id_word&fUseMask);

  if(fabs(photon->DetEta())<fMinPesFid) id_word |= fPesFidBit;
  if(fPrintLevel>0) printf("PhotonID  DetEta=%9.3f  Cut=%7.3f  word=0x%08x\n",
			   photon->DetEta(),fMinPesFid,id_word&fUseMask);

  if(fabs(photon->DetEta())>fMaxPesFid) id_word |= fPesFidBit;
  if(fPrintLevel>0) printf("PhotonID  DetEta=%9.3f  Cut=%7.3f  word=0x%08x\n",
			   photon->DetEta(),fMaxPesFid,id_word&fUseMask);

  if(photon->Chi3x3()>fMaxPemChi2) id_word |= fPemChi2Bit;
  if(fPrintLevel>0) printf("PhotonID  Chi3x3=%9.3f  Cut=%7.3f  word=0x%08x\n",
			   photon->Chi3x3(),fMaxPemChi2,id_word&fUseMask);

  if(photon->Prof5by9U()<fMinPes5By9) id_word |= fPes5By9Bit;
  if(fPrintLevel>0) printf("PhotonID  5x9 U =%9.3f  Cut=%7.3f  word=0x%08x\n",
			   photon->Prof5by9U(),fMinPes5By9,id_word&fUseMask);

  if(photon->Prof5by9V()<fMinPes5By9) id_word |= fPes5By9Bit;
  if(fPrintLevel>0) printf("PhotonID  5x9 V =%9.3f  Cut=%7.3f  word=0x%08x\n",
			   photon->Prof5by9V(),fMinPes5By9,id_word&fUseMask);

  float pemHadEmCut = 0.05;
  if(e>100.0) pemHadEmCut += 0.026*log(e/100.0);
  if(photon->HadEm()>pemHadEmCut) id_word |= fPemHadEmBit;
  if(fPrintLevel>0) printf("PhotonID  PEMH/E=%9.3f  Cut=%7.3f  word=0x%08x\n",
			   photon->HadEm(),pemHadEmCut,id_word&fUseMask);

  int idret = (id_word & fUseMask);
  if(fPrintLevel>0) 
    printf("PhotonID word=0x%08x, mask=0x%08x, masked word=0x%08x \n",
	   id_word,fUseMask,idret);
  return idret;
}

//_____________________________________________________________________________
Int_t TStnPhotonID::StdCentralIDWord() {
  Reset();
  int i=0;
  i |= fCentralBit;
  i |= fEtBit;
  i |= fHadEmBit;
  i |= fCesAvgChi2Bit;
  i |= fXCesBit;
  i |= fZCesBit;
  i |= fIsoBit;
  i |= fN3DBit;
  i |= fPtBit;
  i |= fTIsoBit;
  i |= fCes2Bit;
  SetUseMask(i);
  return i;
}

//_____________________________________________________________________________
Int_t TStnPhotonID::StdCentralLooseIDWord() {
  Reset();
  int i=0;
  i |= fCentralBit;
  i |= fEtBit;
  i |= fHadEmBit;
  SetMaxHadEm(0.125);
  SetHadEmSlope(0.0);
  i |= fXCesBit;
  i |= fZCesBit;
  i |= fIsoBit;
  SetMaxIso4(3.0);
  SetMaxIso4SlideLow(0.15);
  i |= fPtBit;
  SetMaxPt(0.0);
  SetMaxPtSlide(0.25);
  i |= fTIsoBit;
  SetMaxSumPt4(5.0);
  SetMaxSumPt4Slide(0.0);

  SetUseMask(i);
  return i;
}

//_____________________________________________________________________________
Int_t TStnPhotonID::StdPlugIDWord() {
  Reset();
  int i=0;
  i |= fPlugBit;
  i |= fEtBit;
  i |= fIsoBit;
  i |= fTIsoBit;
  i |= fPesFidBit;
  i |= fPemChi2Bit;
  i |= fPes5By9Bit;
  i |= fPemHadEmBit;
  SetUseMask(i);
  return i;
}

//_____________________________________________________________________________
Int_t TStnPhotonID::StdPlugLooseIDWord() {
  Reset();
  int i=0;
  i |= fPlugBit;
  i |= fEtBit;
  i |= fHadEmBit;
  SetMaxHadEm(0.125);
  SetHadEmSlope(0.0);
  i |= fIsoBit;
  SetMaxIso4(3.0);
  SetMaxIso4SlideLow(0.15);
  i |= fTIsoBit;
  SetMaxSumPt4(5.0);
  SetMaxSumPt4Slide(0.0);
  i |= fPesFidBit;

  SetUseMask(i);
  return i;
}


//_____________________________________________________________________________
void TStnPhotonID::Reset() {
  fMinEt = 25.0;
  fMaxDetEta = 999.0;
  fMaxEta = 999.0;
  fMaxHadEm = 0.055;
  fHadEmSlope = 0.00045;
  fMaxCesAvgChi2 = 20.0;
  fMaxCesAvgChi2Break = 9999.0;
  fMaxCesAvgChi2High = 1.0e30;
  fMaxXCes = 21.0;
  fMinZCes = 9.0;
  fMaxZCes = 230.0;
  fMaxLshr = 0.2;
  fMaxIso4Low = 0.0;
  fMaxIso4SlideLow = 0.1;
  fMaxIso4Break = 20.0;
  fMaxIso4 = 2.0;
  fMaxIso4Slide = 0.02;
  fMaxN3D = 1;
  fMaxPt = 1.0;
  fMaxPtSlide = 0.005;
  fMaxSumPt4 = 2.0;
  fMaxSumPt4Slide = 0.005;
  fMinPesE = 0.001;
  fMinPesFid = 1.2;
  fMaxPesFid = 2.8;
  fMaxPemChi2 = 10.0;
  fMinPes5By9 = 0.65;

  fUseMask = 0xFFFFFFFF;
  fPrintLevel=0;
}


//_____________________________________________________________________________
void TStnPhotonID::Print(const char* Opt) const {

  printf(" bit  2: fMinEt             = %12.5f\n"  ,fMinEt          );
  printf(" bit  3: fMaxDetEta         = %12.5f\n"  ,fMaxDetEta      );
  printf(" bit  4: fMaxEta            = %12.5f\n"  ,fMaxEta         );
  printf(" bit  5: fMaxHadEm          = %12.5f\n"  ,fMaxHadEm       );
  printf(" bit  5: fHadEmSlope        = %12.5f\n"  ,fHadEmSlope     );
  printf(" bit  6: fMaxCesAvgChi2     = %12.5f\n"  ,fMaxCesAvgChi2  );
  printf(" bit  7: fMaxXCes           = %12.5f\n"  ,fMaxXCes        );
  printf(" bit  8: fMinZCes           = %12.5f\n"  ,fMinZCes        );
  printf(" bit  8: fMaxZCes           = %12.5f\n"  ,fMaxZCes        );
  printf(" bit  9: fMaxLshr           = %12.5f\n"  ,fMaxLshr        );
  printf(" bit 10: fMaxIso4           = %12.5f\n"  ,fMaxIso4        );
  printf(" bit 10: fMaxIso4Slide      = %12.5f\n"  ,fMaxIso4Slide   );
  printf(" bit 10: fMaxIso4Break      = %12.5f\n"  ,fMaxIso4Break   );
  printf(" bit 10: fMaxIso4Low        = %12.5f\n"  ,fMaxIso4Low     );
  printf(" bit 10: fMaxIso4SlideLow   = %12.5f\n"  ,fMaxIso4SlideLow);
  printf(" bit 11: fMaxN3D            = %12.5f\n"  ,fMaxN3D         );
  printf(" bit 12: fMaxPt             = %12.5f\n"  ,fMaxPt          );
  printf(" bit 12: fMaxPtSlide        = %12.5f\n"  ,fMaxPtSlide     );
  printf(" bit 13: fMaxSumPt4         = %12.5f\n"  ,fMaxSumPt4      );
  printf(" bit 13: fMaxSumPt4Slide    = %12.5f\n"  ,fMaxSumPt4Slide );
  printf(" bit 15: fMinPesE           = %12.5f\n"  ,fMinPesE        );
  printf(" bit 16: fMinPesFid         = %12.5f\n"  ,fMinPesFid      );
  printf(" bit 16: fMaxPesFid         = %12.5f\n"  ,fMaxPesFid      );
  printf(" bit 17: fMaxPemChi2        = %12.5f\n"  ,fMaxPemChi2     );
  printf(" bit 18: fMinPes5By9        = %12.5f\n"  ,fMinPes5By9     );

  printf(" bit  0: Central     %s\n",((fUseMask&(1<< 0))>0?"    used":"not used") );
  printf(" bit  1: Plug        %s\n",((fUseMask&(1<< 1))>0?"    used":"not used") );
  printf(" bit  2: Et          %s\n",((fUseMask&(1<< 2))>0?"    used":"not used") );
  printf(" bit  3: DetEta      %s\n",((fUseMask&(1<< 3))>0?"    used":"not used") );
  printf(" bit  4: Eta         %s\n",((fUseMask&(1<< 4))>0?"    used":"not used") );
  printf(" bit  5: HadEm       %s\n",((fUseMask&(1<< 5))>0?"    used":"not used") );
  printf(" bit  6: CesAvgChi2  %s\n",((fUseMask&(1<< 6))>0?"    used":"not used") );
  printf(" bit  7: XCes        %s\n",((fUseMask&(1<< 7))>0?"    used":"not used") );
  printf(" bit  8: ZCes        %s\n",((fUseMask&(1<< 8))>0?"    used":"not used") );
  printf(" bit  9: Lshr        %s\n",((fUseMask&(1<< 9))>0?"    used":"not used") );
  printf(" bit 10: Iso         %s\n",((fUseMask&(1<<10))>0?"    used":"not used") );
  printf(" bit 11: N3D         %s\n",((fUseMask&(1<<11))>0?"    used":"not used") );
  printf(" bit 12: Pt          %s\n",((fUseMask&(1<<12))>0?"    used":"not used") );
  printf(" bit 13: TIso        %s\n",((fUseMask&(1<<13))>0?"    used":"not used") );
  printf(" bit 14: Ces2        %s\n",((fUseMask&(1<<14))>0?"    used":"not used") );
  printf(" bit 15: PesE        %s\n",((fUseMask&(1<<15))>0?"    used":"not used") );
  printf(" bit 16: PesFid      %s\n",((fUseMask&(1<<16))>0?"    used":"not used") );
  printf(" bit 17: PemChi2     %s\n",((fUseMask&(1<<17))>0?"    used":"not used") );
  printf(" bit 18: Pes5By9     %s\n",((fUseMask&(1<<18))>0?"    used":"not used") );
  printf(" bit 19: PesHadEm    %s\n",((fUseMask&(1<<19))>0?"    used":"not used") );


}

//______________________________________________________________________________
void TStnPhotonID::Streamer(TBuffer &R__b)
{
   if (R__b.IsReading()) {
      Version_t R__v = R__b.ReadVersion(); if (R__v) { }
      TNamed::Streamer(R__b);
      R__b >> fMinEt;	      
      R__b >> fMaxDetEta;      
      R__b >> fMaxEta;	      
      R__b >> fMaxHadEm;	      
      R__b >> fHadEmSlope;     
      R__b >> fMaxCesAvgChi2;  
      R__b >> fMaxXCes;	      
      R__b >> fMinZCes;	      
      R__b >> fMaxZCes;	      
      R__b >> fMaxLshr;	      
      R__b >> fMaxIso4;	      
      R__b >> fMaxIso4Slide;   
      R__b >> fMaxIso4Break;   
      R__b >> fMaxIso4Low;     
      R__b >> fMaxIso4SlideLow;
      R__b >> fMaxN3D;	      
      R__b >> fMaxPt;	      
      R__b >> fMaxPtSlide;     
      R__b >> fMaxSumPt4;      
      R__b >> fMaxSumPt4Slide; 
      R__b >> fMinPesE;	      
      R__b >> fMinPesFid;      
      R__b >> fMaxPesFid;      
      R__b >> fMaxPemChi2;     
      R__b >> fMinPes5By9;     
      R__b >> fUseMask;     
   } else {
      R__b.WriteVersion(TStnPhotonID::IsA());
      TNamed::Streamer(R__b);
      R__b << fMinEt;	      
      R__b << fMaxDetEta;      
      R__b << fMaxEta;	      
      R__b << fMaxHadEm;	      
      R__b << fHadEmSlope;     
      R__b << fMaxCesAvgChi2;  
      R__b << fMaxXCes;	      
      R__b << fMinZCes;	      
      R__b << fMaxZCes;	      
      R__b << fMaxLshr;	      
      R__b << fMaxIso4;	      
      R__b << fMaxIso4Slide;   
      R__b << fMaxIso4Break;   
      R__b << fMaxIso4Low;     
      R__b << fMaxIso4SlideLow;
      R__b << fMaxN3D;	      
      R__b << fMaxPt;	      
      R__b << fMaxPtSlide;     
      R__b << fMaxSumPt4;      
      R__b << fMaxSumPt4Slide; 
      R__b << fMinPesE;	      
      R__b << fMinPesFid;      
      R__b << fMaxPesFid;      
      R__b << fMaxPemChi2;     
      R__b << fMinPes5By9;     
      R__b << fUseMask;     
   }
}