///////////////////////////////////////////////////////////////////////////
//  Code for applying Trigger Scale Factor (SF) in electroweak analyses  //
//  Marilyn Marx (marx@cern.ch)                                          //
//  Junjie Zhu (junjie@umich.edu)                                        //
//  (August 17, 2011)                                                    //
//  (November 12, 2011) - updated by marx                                //
///////////////////////////////////////////////////////////////////////////
//
// 1) For single boson, WW, ZZ or other analyses:                
// If your analysis ONLY applies COMBINED muon efficiency/SF numbers, call SF function as on line 163.       
//
///////////////////////////////////////////////////////////////////////////
//
// 2) For WZ analysis:
// If your analysis can apply TIGHT AND COMBINED muon numbers, call SF function as on line 215,
// where runnumber is an int associating to a real data runnumber in mc, 
// useGeV is a bool you set depending on whether quantities are in MeV (useGeV=false) or GeV (useGeV=true) in your code,
// channel is a string that can take values like "TriMu", "TriEm", "DiMuEm", "DiEmMu"
// and zl1 and zl2 are ints determining the position of the Z leptons in the lepton vector (in ambiguous cases like eee and mmm channels).
// If the useTight bool is set to true, tight numbers are applied to Z leptons and combined numbers to W leptons.
// If the useTight bool is set to false, combined numbers are applied to Z and W leptons.                         
// (N.B. The latter is not the same as calling GetTriggerSF() as in 1) because the WZ analysis does not require exactly 3 leptons.)
//
///////////////////////////////////////////////////////////////////////////

#include "LeptonTriggerSF.h"
#include "electron_SF_VV.h"

#include "TFile.h"
#include "TROOT.h"
#include <iostream>
#include <vector>
#include <stdio.h>
#include <math.h>

using namespace std;
using namespace TrigMuonEff;
using TrigMuonEff::SFDataPeriod;
const double pi = acos(-1.);

void LeptonTriggerSF::closefile(){
  
  file_Muon_Trig_Eff->Close();
  delete file_Muon_Trig_Eff;
}

void LeptonTriggerSF::initialize() {
  // open root file
  cout << "LeptonTriggerSF : Initializing" << endl;
  file_Muon_Trig_Eff = new TFile((dir+"muon_trigger_sf.root").c_str());
  if(file_Muon_Trig_Eff) cout << "LeptonTriggerSF : File " << file_Muon_Trig_Eff->GetName() << " was found." << endl;
  else cout << "LeptonTriggerSF : ERROR Rootfile failed to open" << endl;

  // read histograms
  _MuonTrigEffB2toI_EtaPhi_Barrel_Data = (TH2D *)(file_Muon_Trig_Eff->Get("mu18_MG/etaphi_barrel_eff_data_periodB2toI_EF_mu18_MG"));
  _MuonTrigEffB2toI_EtaPhi_Barrel_MC = (TH2D *)(file_Muon_Trig_Eff->Get("mu18_MG/etaphi_barrel_eff_mc_periodB2toI_EF_mu18_MG"));
  _MuonTrigEffB2toI_EtaPhi_Endcap_Data = (TH2D *)(file_Muon_Trig_Eff->Get("mu18_MG/etaphi_endcap_eff_data_periodB2toI_EF_mu18_MG"));
  _MuonTrigEffB2toI_EtaPhi_Endcap_MC = (TH2D *)(file_Muon_Trig_Eff->Get("mu18_MG/etaphi_endcap_eff_mc_periodB2toI_EF_mu18_MG"));
  _MuonTrigEffJtoM_EtaPhi_Barrel_Data = (TH2D *)(file_Muon_Trig_Eff->Get("mu18_MG_medium/etaphi_barrel_eff_data_periodJtoM_EF_mu18_MG_medium"));
  _MuonTrigEffJtoM_EtaPhi_Barrel_MC = (TH2D *)(file_Muon_Trig_Eff->Get("mu18_MG_medium/etaphi_barrel_eff_mc_periodJtoM_EF_mu18_MG_medium"));
  _MuonTrigEffJtoM_EtaPhi_Endcap_Data = (TH2D *)(file_Muon_Trig_Eff->Get("mu18_MG_medium/etaphi_endcap_eff_data_periodJtoM_EF_mu18_MG_medium"));
  _MuonTrigEffJtoM_EtaPhi_Endcap_MC = (TH2D *)(file_Muon_Trig_Eff->Get("mu18_MG_medium/etaphi_endcap_eff_mc_periodJtoM_EF_mu18_MG_medium"));

  _MuonTrigEffJ_EtaPhi_Barrel_Data = (TH2D *)(file_Muon_Trig_Eff->Get("mu18_MG_medium_J/etaphi_barrel_eff_data_periodJ_EF_mu18_MG_medium"));
  _MuonTrigEffJ_EtaPhi_Barrel_MC = (TH2D *)(file_Muon_Trig_Eff->Get("mu18_MG_medium_J/etaphi_barrel_eff_mc_periodJ_EF_mu18_MG_medium"));
  _MuonTrigEffJ_EtaPhi_Endcap_Data = (TH2D *)(file_Muon_Trig_Eff->Get("mu18_MG_medium_J/etaphi_endcap_eff_data_periodJ_EF_mu18_MG_medium"));
  _MuonTrigEffJ_EtaPhi_Endcap_MC = (TH2D *)(file_Muon_Trig_Eff->Get("mu18_MG_medium_J/etaphi_endcap_eff_mc_periodJ_EF_mu18_MG_medium"));
  _MuonTrigEffK_EtaPhi_Barrel_Data = (TH2D *)(file_Muon_Trig_Eff->Get("mu18_MG_medium_K/etaphi_barrel_eff_data_periodK_EF_mu18_MG_medium"));
  _MuonTrigEffK_EtaPhi_Barrel_MC = (TH2D *)(file_Muon_Trig_Eff->Get("mu18_MG_medium_K/etaphi_barrel_eff_mc_periodK_EF_mu18_MG_medium"));
  _MuonTrigEffK_EtaPhi_Endcap_Data = (TH2D *)(file_Muon_Trig_Eff->Get("mu18_MG_medium_K/etaphi_endcap_eff_data_periodK_EF_mu18_MG_medium"));
  _MuonTrigEffK_EtaPhi_Endcap_MC = (TH2D *)(file_Muon_Trig_Eff->Get("mu18_MG_medium_K/etaphi_endcap_eff_mc_periodK_EF_mu18_MG_medium"));
  _MuonTrigEffK1_EtaPhi_Barrel_Data = (TH2D *)(file_Muon_Trig_Eff->Get("mu18_MG_medium_K1/etaphi_barrel_eff_data_periodK1_EF_mu18_MG_medium"));
  _MuonTrigEffK1_EtaPhi_Barrel_MC = (TH2D *)(file_Muon_Trig_Eff->Get("mu18_MG_medium_K1/etaphi_barrel_eff_mc_periodK1_EF_mu18_MG_medium"));
  _MuonTrigEffK1_EtaPhi_Endcap_Data = (TH2D *)(file_Muon_Trig_Eff->Get("mu18_MG_medium_K1/etaphi_endcap_eff_data_periodK1_EF_mu18_MG_medium"));
  _MuonTrigEffK1_EtaPhi_Endcap_MC = (TH2D *)(file_Muon_Trig_Eff->Get("mu18_MG_medium_K1/etaphi_endcap_eff_mc_periodK1_EF_mu18_MG_medium"));
  _MuonTrigEffK2_K6_EtaPhi_Barrel_Data = (TH2D *)(file_Muon_Trig_Eff->Get("mu18_MG_medium_K2toK6/etaphi_barrel_eff_data_periodK2toK6_EF_mu18_MG_medium"));
  _MuonTrigEffK2_K6_EtaPhi_Barrel_MC = (TH2D *)(file_Muon_Trig_Eff->Get("mu18_MG_medium_K2toK6/etaphi_barrel_eff_mc_periodK2toK6_EF_mu18_MG_medium"));
  _MuonTrigEffK2_K6_EtaPhi_Endcap_Data = (TH2D *)(file_Muon_Trig_Eff->Get("mu18_MG_medium_K2toK6/etaphi_endcap_eff_data_periodK2toK6_EF_mu18_MG_medium"));
  _MuonTrigEffK2_K6_EtaPhi_Endcap_MC = (TH2D *)(file_Muon_Trig_Eff->Get("mu18_MG_medium_K2toK6/etaphi_endcap_eff_mc_periodK2toK6_EF_mu18_MG_medium"));
  _MuonTrigEffL_EtaPhi_Barrel_Data = (TH2D *)(file_Muon_Trig_Eff->Get("mu18_MG_medium_L/etaphi_barrel_eff_data_periodL_EF_mu18_MG_medium"));
  _MuonTrigEffL_EtaPhi_Barrel_MC = (TH2D *)(file_Muon_Trig_Eff->Get("mu18_MG_medium_L/etaphi_barrel_eff_mc_periodL_EF_mu18_MG_medium"));
  _MuonTrigEffL_EtaPhi_Endcap_Data = (TH2D *)(file_Muon_Trig_Eff->Get("mu18_MG_medium_L/etaphi_endcap_eff_data_periodL_EF_mu18_MG_medium"));
  _MuonTrigEffL_EtaPhi_Endcap_MC = (TH2D *)(file_Muon_Trig_Eff->Get("mu18_MG_medium_L/etaphi_endcap_eff_mc_periodL_EF_mu18_MG_medium"));
  _MuonTrigEffL3_L4_EtaPhi_Barrel_Data = (TH2D *)(file_Muon_Trig_Eff->Get("mu18_MG_medium_L3toL4/etaphi_barrel_eff_data_periodL3toL4_EF_mu18_MG_medium"));
  _MuonTrigEffL3_L4_EtaPhi_Barrel_MC = (TH2D *)(file_Muon_Trig_Eff->Get("mu18_MG_medium_L3toL4/etaphi_barrel_eff_mc_periodL3toL4_EF_mu18_MG_medium"));
  _MuonTrigEffL3_L4_EtaPhi_Endcap_Data = (TH2D *)(file_Muon_Trig_Eff->Get("mu18_MG_medium_L3toL4/etaphi_endcap_eff_data_periodL3toL4_EF_mu18_MG_medium"));
  _MuonTrigEffL3_L4_EtaPhi_Endcap_MC = (TH2D *)(file_Muon_Trig_Eff->Get("mu18_MG_medium_L3toL4/etaphi_endcap_eff_mc_periodL3toL4_EF_mu18_MG_medium"));
  _MuonTrigEffLwoL3_L4_EtaPhi_Barrel_Data = (TH2D *)(file_Muon_Trig_Eff->Get("mu18_MG_medium_LwoL3L4/etaphi_barrel_eff_data_periodLwoL3L4_EF_mu18_MG_medium"));
  _MuonTrigEffLwoL3_L4_EtaPhi_Barrel_MC = (TH2D *)(file_Muon_Trig_Eff->Get("mu18_MG_medium_LwoL3L4/etaphi_barrel_eff_mc_periodLwoL3L4_EF_mu18_MG_medium"));
  _MuonTrigEffLwoL3_L4_EtaPhi_Endcap_Data = (TH2D *)(file_Muon_Trig_Eff->Get("mu18_MG_medium_LwoL3L4/etaphi_endcap_eff_data_periodLwoL3L4_EF_mu18_MG_medium"));
  _MuonTrigEffLwoL3_L4_EtaPhi_Endcap_MC = (TH2D *)(file_Muon_Trig_Eff->Get("mu18_MG_medium_LwoL3L4/etaphi_endcap_eff_mc_periodLwoL3L4_EF_mu18_MG_medium"));
  _MuonTrigEffM_EtaPhi_Barrel_Data = (TH2D *)(file_Muon_Trig_Eff->Get("mu18_MG_medium_M/etaphi_barrel_eff_data_periodM_EF_mu18_MG_medium"));
  _MuonTrigEffM_EtaPhi_Barrel_MC = (TH2D *)(file_Muon_Trig_Eff->Get("mu18_MG_medium_M/etaphi_barrel_eff_mc_periodM_EF_mu18_MG_medium"));
  _MuonTrigEffM_EtaPhi_Endcap_Data = (TH2D *)(file_Muon_Trig_Eff->Get("mu18_MG_medium_M/etaphi_endcap_eff_data_periodM_EF_mu18_MG_medium"));
  _MuonTrigEffM_EtaPhi_Endcap_MC = (TH2D *)(file_Muon_Trig_Eff->Get("mu18_MG_medium_M/etaphi_endcap_eff_mc_periodM_EF_mu18_MG_medium"));
  if(!(_MuonTrigEffB2toI_EtaPhi_Endcap_MC&&_MuonTrigEffB2toI_EtaPhi_Endcap_Data&&_MuonTrigEffB2toI_EtaPhi_Barrel_MC&&_MuonTrigEffB2toI_EtaPhi_Barrel_Data&&_MuonTrigEffJtoM_EtaPhi_Endcap_MC&&_MuonTrigEffJtoM_EtaPhi_Endcap_Data&&_MuonTrigEffJtoM_EtaPhi_Barrel_MC&&_MuonTrigEffJtoM_EtaPhi_Barrel_Data)){    
    cout << "Not all histograms could be initialized, exiting." << endl;
    exit (1);
  }
  if(!(_MuonTrigEffJ_EtaPhi_Barrel_Data&&_MuonTrigEffJ_EtaPhi_Barrel_MC&&_MuonTrigEffJ_EtaPhi_Endcap_Data&&_MuonTrigEffJ_EtaPhi_Endcap_MC&&_MuonTrigEffK_EtaPhi_Barrel_Data&&_MuonTrigEffK_EtaPhi_Barrel_MC&&_MuonTrigEffK_EtaPhi_Endcap_Data&&_MuonTrigEffK_EtaPhi_Endcap_MC&&_MuonTrigEffK1_EtaPhi_Barrel_Data&&_MuonTrigEffK1_EtaPhi_Barrel_MC&&_MuonTrigEffK1_EtaPhi_Endcap_Data&&_MuonTrigEffK1_EtaPhi_Endcap_MC&&_MuonTrigEffK2_K6_EtaPhi_Barrel_Data&&_MuonTrigEffK2_K6_EtaPhi_Barrel_MC&&_MuonTrigEffK2_K6_EtaPhi_Endcap_Data&&_MuonTrigEffK2_K6_EtaPhi_Endcap_MC&&_MuonTrigEffL_EtaPhi_Barrel_Data&&_MuonTrigEffL_EtaPhi_Barrel_MC&&_MuonTrigEffL_EtaPhi_Endcap_Data&&_MuonTrigEffL_EtaPhi_Endcap_MC&&_MuonTrigEffL3_L4_EtaPhi_Barrel_Data&&_MuonTrigEffL3_L4_EtaPhi_Barrel_MC&&_MuonTrigEffL3_L4_EtaPhi_Endcap_Data&&_MuonTrigEffL3_L4_EtaPhi_Endcap_MC&&_MuonTrigEffLwoL3_L4_EtaPhi_Barrel_Data&&_MuonTrigEffLwoL3_L4_EtaPhi_Barrel_MC&&_MuonTrigEffLwoL3_L4_EtaPhi_Endcap_Data&&_MuonTrigEffLwoL3_L4_EtaPhi_Endcap_MC&&_MuonTrigEffM_EtaPhi_Barrel_Data&&_MuonTrigEffM_EtaPhi_Barrel_MC&&_MuonTrigEffM_EtaPhi_Endcap_Data&&_MuonTrigEffM_EtaPhi_Endcap_MC)){    
    cout << "Not all histograms could be initialized, exiting." << endl;
    exit (1);
  }
  else cout << "LeptonTriggerSF : Initialization of LeptonTriggerSF successful " << endl;

  m_minimum_phi = _MuonTrigEffJtoM_EtaPhi_Barrel_Data->GetYaxis()->GetXmin();
  return;
}

/// Constructor
LeptonTriggerSF::LeptonTriggerSF(string path){
  dir=path;
  inGeV=false;
  muThreshold=20000.;
  emThreshold=25000.;
  initialize();
}

/// Destructor 
LeptonTriggerSF::~LeptonTriggerSF() {
  delete _MuonTrigEffB2toI_EtaPhi_Barrel_Data;
  delete _MuonTrigEffB2toI_EtaPhi_Barrel_MC;
  delete _MuonTrigEffB2toI_EtaPhi_Endcap_Data;
  delete _MuonTrigEffB2toI_EtaPhi_Endcap_MC;
  delete _MuonTrigEffJtoM_EtaPhi_Barrel_Data;
  delete _MuonTrigEffJtoM_EtaPhi_Barrel_MC;
  delete _MuonTrigEffJtoM_EtaPhi_Endcap_Data;
  delete _MuonTrigEffJtoM_EtaPhi_Endcap_MC;

  delete _MuonTrigEffJ_EtaPhi_Barrel_Data;
  delete _MuonTrigEffJ_EtaPhi_Barrel_MC;
  delete _MuonTrigEffJ_EtaPhi_Endcap_Data;
  delete _MuonTrigEffJ_EtaPhi_Endcap_MC;
  delete _MuonTrigEffK_EtaPhi_Barrel_Data;
  delete _MuonTrigEffK_EtaPhi_Barrel_MC;
  delete _MuonTrigEffK_EtaPhi_Endcap_Data;
  delete _MuonTrigEffK_EtaPhi_Endcap_MC;
  delete _MuonTrigEffK1_EtaPhi_Barrel_Data;
  delete _MuonTrigEffK1_EtaPhi_Barrel_MC;
  delete _MuonTrigEffK1_EtaPhi_Endcap_Data;
  delete _MuonTrigEffK1_EtaPhi_Endcap_MC;
  delete _MuonTrigEffK2_K6_EtaPhi_Barrel_Data;
  delete _MuonTrigEffK2_K6_EtaPhi_Barrel_MC;
  delete _MuonTrigEffK2_K6_EtaPhi_Endcap_Data;
  delete _MuonTrigEffK2_K6_EtaPhi_Endcap_MC;
  delete _MuonTrigEffL_EtaPhi_Barrel_Data;
  delete _MuonTrigEffL_EtaPhi_Barrel_MC;
  delete _MuonTrigEffL_EtaPhi_Endcap_Data;
  delete _MuonTrigEffL_EtaPhi_Endcap_MC;
  delete _MuonTrigEffL3_L4_EtaPhi_Barrel_Data;
  delete _MuonTrigEffL3_L4_EtaPhi_Barrel_MC;
  delete _MuonTrigEffL3_L4_EtaPhi_Endcap_Data;
  delete _MuonTrigEffL3_L4_EtaPhi_Endcap_MC;
  delete _MuonTrigEffLwoL3_L4_EtaPhi_Barrel_Data;
  delete _MuonTrigEffLwoL3_L4_EtaPhi_Barrel_MC;
  delete _MuonTrigEffLwoL3_L4_EtaPhi_Endcap_Data;
  delete _MuonTrigEffLwoL3_L4_EtaPhi_Endcap_MC;
  delete _MuonTrigEffM_EtaPhi_Barrel_Data;
  delete _MuonTrigEffM_EtaPhi_Barrel_MC;
  delete _MuonTrigEffM_EtaPhi_Endcap_Data;
  delete _MuonTrigEffM_EtaPhi_Endcap_MC;

  closefile(); 

}


double LeptonTriggerSF::GetTriggerSF(int runnumber,bool useGeV, std::vector<TLorentzVector> muons, std::vector<TLorentzVector> electrons){
  
  SFDataPeriod period = getDataPeriod(runnumber);
  if(period==perUnDefined) {
    cout << "RunNumber is not in 2011 dataset. No scale factors calculated. Please use Runnumber between 178044-191933" << endl;
    return 0.;
  }

  setThresholds(useGeV);

  double rate_not_fired_data = 1.;
  double rate_not_fired_mc = 1.;
  
  for(unsigned int ielec = 0; ielec < electrons.size(); ielec ++) {
    double eff_data = 0., eff_mc = 0.;    

    // get efficiency from data
    eff_data= EmEff_Data(electrons[ielec]);
    // get efficiency from MC
    eff_mc=  EmEff_MC(electrons[ielec]);
    
    rate_not_fired_data *= (1-eff_data);
    rate_not_fired_mc *= (1-eff_mc);
  }
  
  for(unsigned int imuon = 0; imuon < muons.size(); imuon ++) {
    double eff_data = 0., eff_mc = 0.;
    
    if(muons[imuon].Pt() < muThreshold){
      eff_data=0;
      eff_mc=0;
    }

    else{
      // get efficiency from data
      eff_data= MuEff_CB(period,true, muons[imuon]);    
      // get efficiency from MC
      eff_mc=  MuEff_CB(period,false, muons[imuon]);
    }

    rate_not_fired_data *= (1-eff_data);
    rate_not_fired_mc *= (1-eff_mc);
  }

  double event_SF = 1.;
  
  // prevent events with no triggered electrons or muons
  if ((electrons.size() != 0 || muons.size() != 0) && ((1-rate_not_fired_mc) != 0))
    event_SF = (1-rate_not_fired_data)/(1-rate_not_fired_mc);
  
  return event_SF; 
}
  

double LeptonTriggerSF::GetTriggerSF(int runnumber,bool useGeV, std::vector<TLorentzVector> muons,std::vector<TLorentzVector> electrons, TString channel , int zl1, int zl2, bool useTight){

  SFDataPeriod period = getDataPeriod(runnumber);
  if(period==perUnDefined) {
    cout << "RunNumber is not in 2011 dataset. No scale factors calculated. Please use Runnumber between 178044-191933" << endl;
    return 0.;
  }

  setThresholds(useGeV);

  int ichannel=-1;
  if(channel.Contains("TriMu")) ichannel=1;
  if(channel.Contains("TriEm")) ichannel=2;
  if(channel.Contains("DiMuEm")) ichannel=3;
  if(channel.Contains("DiEmMu")) ichannel=4;
  
  if(ichannel==1 && useTight) return TriMuSFTight(period,muons, zl1, zl2);
  if(ichannel==3 && useTight) return DiMuEmSFTight(period,muons,electrons.at(0),zl1,zl2);
  if(ichannel==1 && !useTight) return TriMuSF(period,muons, zl1, zl2);	
  if(ichannel==3 && !useTight) return DiMuEmSF(period,muons,electrons.at(0),zl1,zl2);

  if(ichannel==2) return TriEmSF(electrons,zl1,zl2 );
  if(ichannel==4) return DiEmMuSF(period,electrons,muons.at(0),zl1,zl2);

  return 0.;
}


double LeptonTriggerSF::TriMuSF(SFDataPeriod period,std::vector<TLorentzVector> muons,int zl1, int zl2) const{
  
  double mu1Eff_data= MuEff_CB(period,true, muons.at(zl1));   
  double mu2Eff_data= MuEff_CB(period,true, muons.at(zl2));
  
  double mu1Eff_mc= MuEff_CB(period,false, muons.at(zl1));   
  double mu2Eff_mc= MuEff_CB(period,false, muons.at(zl2));
  
  if(muons.at(zl1).Pt() < muThreshold){
    mu1Eff_data=0;
    mu1Eff_mc=0;
  }
  
  if(muons.at(zl2).Pt() < muThreshold){
    mu2Eff_data = 0;
    mu2Eff_mc = 0;
  }
  
  int lastmu=0;
  
  for(int i=0; i<10; ++i){
    if(i==zl1 || i==zl2)continue;
    lastmu=i;
    break;
  }  
  
  double mu3Eff_mc= MuEff_CB(period,false, muons.at(lastmu));
  double mu3Eff_data= MuEff_CB(period,true, muons.at(lastmu));
  
  if(muons.at(lastmu).Pt() < muThreshold){
    mu3Eff_data = 0; 
    mu3Eff_mc = 0; 
  }
  
  double Eff_data =  1. - (1.- mu1Eff_data) * (1.- mu2Eff_data) * (1.- mu3Eff_data);    
  double Eff_mc =  1. - (1.- mu1Eff_mc) * (1.-mu2Eff_mc) * (1.-mu3Eff_mc);
 
  double TriMuEventSF = 0;
  if(Eff_data == 0 || Eff_mc == 0) TriMuEventSF = 1.;
  else TriMuEventSF = (Eff_data / Eff_mc);  
  return TriMuEventSF;
}

double LeptonTriggerSF::TriMuSFTight(SFDataPeriod period,std::vector<TLorentzVector> muons,int zl1, int zl2) const{
  
  double mu1Eff_data= MuEff_Tight(period,true,muons.at(zl1));   
  double mu2Eff_data= MuEff_Tight(period,true, muons.at(zl2));
  
  double mu1Eff_mc= MuEff_Tight(period,false, muons.at(zl1));   
  double mu2Eff_mc= MuEff_Tight(period,false, muons.at(zl2));
  
  if(muons.at(zl1).Pt() < muThreshold){
    mu1Eff_data=0;
    mu1Eff_mc=0;
  }
  
  if(muons.at(zl2).Pt() < muThreshold){
    mu2Eff_data = 0;
    mu2Eff_mc = 0;
  }
  
  int lastmu=0;
  
  for(int i=0; i<10; ++i){
    if(i==zl1 || i==zl2)continue;
    lastmu=i;
    break;
  }  
  
  double mu3Eff_mc= MuEff_CB(period,false, muons.at(lastmu));
  double mu3Eff_data= MuEff_CB(period,true, muons.at(lastmu));
  
  if(muons.at(lastmu).Pt() < muThreshold){
    mu3Eff_data = 0; 
    mu3Eff_mc = 0; 
  }
  
  double Eff_data =  1. - (1.- mu1Eff_data) * (1.- mu2Eff_data) * (1.- mu3Eff_data);    
  double Eff_mc =  1. - (1.- mu1Eff_mc) * (1.-mu2Eff_mc) * (1.-mu3Eff_mc);    
  
  double TriMuEventSF = 0;
  if(Eff_data == 0 || Eff_mc == 0) TriMuEventSF = 1.;
  else TriMuEventSF = (Eff_data / Eff_mc);  
  return TriMuEventSF;
}

double LeptonTriggerSF::TriEmSF(std::vector<TLorentzVector> electrons,int zl1, int zl2) const{
  
  double em1Eff_data= EmEff_Data(electrons.at(zl1));   
  double em2Eff_data= EmEff_Data(electrons.at(zl2));
  
  double em1Eff_mc= EmEff_MC(electrons.at(zl1));   
  double em2Eff_mc= EmEff_MC(electrons.at(zl2));

  if(electrons.at(zl1).Pt() < emThreshold){
    em1Eff_data=0;   
    em1Eff_mc=0;   
  }
  if(electrons.at(zl2).Pt() < emThreshold) {
    em2Eff_data=0;   
    em2Eff_mc=0;  
  }
  
  int lastem=0;
  
  for(int i=0; i<10; i++){
    if(i==zl1 || i==zl2)continue;
    lastem=i;
    break;
  }   
  
  double em3Eff_mc= EmEff_MC(electrons.at(lastem));   
  double em3Eff_data= EmEff_Data(electrons.at(lastem));   

  if(electrons.at(lastem).Pt() < emThreshold){
    em3Eff_data=0;   
    em3Eff_mc=0;  
  }
  
  double Eff_data =  1. - (1.- em1Eff_data) * (1.- em2Eff_data) * (1.- em3Eff_data);    
  double Eff_mc =  1. - (1.- em1Eff_mc) * (1.-em2Eff_mc) * (1.-em3Eff_mc);   
  
  double TriEmEventSF = 0;
  if(Eff_data == 0 || Eff_mc == 0) TriEmEventSF = 1.;
  else TriEmEventSF = (Eff_data / Eff_mc);
  return TriEmEventSF;
}

double LeptonTriggerSF::DiEmMuSF(SFDataPeriod period,std::vector<TLorentzVector> electrons, TLorentzVector muon,  int zl1, int zl2) const{

  double em1Eff_data= EmEff_Data(electrons.at(zl1));   
  double em2Eff_data= EmEff_Data(electrons.at(zl2));
  double em1Eff_mc= EmEff_MC(electrons.at(zl1));   
  double em2Eff_mc= EmEff_MC(electrons.at(zl2));

  if(electrons.at(zl1).Pt() < emThreshold){
    em1Eff_data=0;   
    em1Eff_mc=0;     
  }
  if( electrons.at(zl2).Pt() < emThreshold){
    em2Eff_data=0;   
    em2Eff_mc=0; 
  }

  double mu1Eff_mc= MuEff_CB(period,false, muon);
  double mu1Eff_data= MuEff_CB(period,true, muon);

  if( muon.Pt() < muThreshold) {
    mu1Eff_data=0;   
    mu1Eff_mc=0;
  }

  double Eff_data =  1. - (1.- em1Eff_data) * (1.- em2Eff_data) * (1.- mu1Eff_data);    
  double Eff_mc =  1. - (1.- em1Eff_mc) * (1.-em2Eff_mc) * (1.-mu1Eff_mc);    
  
  double DiEmMuEventSF = 0;
  if(Eff_data == 0 || Eff_mc == 0) DiEmMuEventSF = 1.;
  else DiEmMuEventSF = (Eff_data / Eff_mc);
  return DiEmMuEventSF;
}

double LeptonTriggerSF::DiMuEmSF(SFDataPeriod period,std::vector<TLorentzVector> muons, TLorentzVector electron,  int zl1, int zl2) const{

  double mu1Eff_data= MuEff_CB(period,true, muons.at(zl1));   
  double mu2Eff_data= MuEff_CB(period,true, muons.at(zl2));
  double mu1Eff_mc= MuEff_CB(period,false, muons.at(zl1));   
  double mu2Eff_mc= MuEff_CB(period,false, muons.at(zl2));
  
  if(muons.at(zl1).Pt() < muThreshold) {
    mu1Eff_data=0;   
    mu1Eff_mc=0;
  }
  if(muons.at(zl2).Pt() < muThreshold) {
    mu2Eff_data=0;   
    mu2Eff_mc=0;
  }
  
  double em1Eff_data= EmEff_Data(electron);
  double em1Eff_mc= EmEff_MC(electron);  
  
  if( electron.Pt() < emThreshold) {
    em1Eff_data=0;   
    em1Eff_mc=0;    
  }

  double Eff_data =  1. - (1.- em1Eff_data) * (1.- mu1Eff_data) * (1.- mu2Eff_data);    
  double Eff_mc =  1. - (1.- em1Eff_mc) * (1.-mu1Eff_mc) * (1.-mu2Eff_mc);    
  
  double DiMuEmEventSF = 0;
  if(Eff_data == 0 || Eff_mc == 0) DiMuEmEventSF = 1.;
  else DiMuEmEventSF = (Eff_data / Eff_mc); 
  return DiMuEmEventSF;
}

double LeptonTriggerSF::DiMuEmSFTight(SFDataPeriod period,std::vector<TLorentzVector> muons, TLorentzVector electron,  int zl1, int zl2) const{

  double mu1Eff_data= MuEff_Tight(period,true, muons.at(zl1));   
  double mu2Eff_data= MuEff_Tight(period,true, muons.at(zl2));
  double mu1Eff_mc= MuEff_Tight(period,false, muons.at(zl1));   
  double mu2Eff_mc= MuEff_Tight(period,false, muons.at(zl2));
  
  if(muons.at(zl1).Pt() < muThreshold) {
    mu1Eff_data=0;   
    mu1Eff_mc=0;
  }
  if(muons.at(zl2).Pt() < muThreshold) {
    mu2Eff_data=0;   
    mu2Eff_mc=0;
  }
  
  double em1Eff_data= EmEff_Data(electron);
  double em1Eff_mc= EmEff_MC(electron);  
  
  if( electron.Pt() < emThreshold) {
    em1Eff_data=0;   
    em1Eff_mc=0;    
  }

  double Eff_data =  1. - (1.- em1Eff_data) * (1.- mu1Eff_data) * (1.- mu2Eff_data);    
  double Eff_mc =  1. - (1.- em1Eff_mc) * (1.-mu1Eff_mc) * (1.-mu2Eff_mc);    
  
  double DiMuEmEventSF = 0;
  if(Eff_data == 0 || Eff_mc == 0) DiMuEmEventSF = 1.;
  else DiMuEmEventSF = (Eff_data / Eff_mc); 
  return DiMuEmEventSF;
}

double LeptonTriggerSF::MuEff_CB(SFDataPeriod period,bool isData, TLorentzVector muon) const{

  const double mu_CB_eta=muon.Eta();
  double mu_CB_phi= check_Phi_range(muon.Phi());  
  //double mu_pt=muon.Pt();
  //if(!inGeV) mu_pt=mu_pt/1000.;

  Int_t bin = -1;
  double eff = 0;

  if(mu_CB_phi < m_minimum_phi) mu_CB_phi+=2*pi; // fix phi range for unusual binning
  
  if(fabs(mu_CB_eta) < 1.05) {  //barrel
    if(isData){
      if(period==per2011B2_I)bin = _MuonTrigEffB2toI_EtaPhi_Barrel_Data->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011B2_I)eff = _MuonTrigEffB2toI_EtaPhi_Barrel_Data->GetBinContent(bin);
      if(period==per2011J_M)bin = _MuonTrigEffJtoM_EtaPhi_Barrel_Data->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011J_M)eff = _MuonTrigEffJtoM_EtaPhi_Barrel_Data->GetBinContent(bin);

      if(period==per2011J)bin = _MuonTrigEffJ_EtaPhi_Barrel_Data->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011J)eff = _MuonTrigEffJ_EtaPhi_Barrel_Data->GetBinContent(bin);
      if(period==per2011K)bin = _MuonTrigEffK_EtaPhi_Barrel_Data->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011K)eff = _MuonTrigEffK_EtaPhi_Barrel_Data->GetBinContent(bin);
      if(period==per2011K1)bin = _MuonTrigEffK1_EtaPhi_Barrel_Data->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011K1)eff = _MuonTrigEffK1_EtaPhi_Barrel_Data->GetBinContent(bin);
      if(period==per2011K2_K6)bin = _MuonTrigEffK2_K6_EtaPhi_Barrel_Data->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011K2_K6)eff = _MuonTrigEffK2_K6_EtaPhi_Barrel_Data->GetBinContent(bin);
      if(period==per2011L)bin = _MuonTrigEffL_EtaPhi_Barrel_Data->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011L)eff = _MuonTrigEffL_EtaPhi_Barrel_Data->GetBinContent(bin);
      if(period==per2011L3_L4)bin = _MuonTrigEffL3_L4_EtaPhi_Barrel_Data->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011L3_L4)eff = _MuonTrigEffL3_L4_EtaPhi_Barrel_Data->GetBinContent(bin);
      if(period==per2011LwoL3_L4)bin = _MuonTrigEffLwoL3_L4_EtaPhi_Barrel_Data->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011LwoL3_L4)eff = _MuonTrigEffLwoL3_L4_EtaPhi_Barrel_Data->GetBinContent(bin);
      if(period==per2011M)bin = _MuonTrigEffM_EtaPhi_Barrel_Data->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011M)eff = _MuonTrigEffM_EtaPhi_Barrel_Data->GetBinContent(bin);
    }
    else{
       if(period==per2011B2_I)bin = _MuonTrigEffB2toI_EtaPhi_Barrel_MC->FindBin(mu_CB_eta, mu_CB_phi);
       if(period==per2011B2_I)eff = _MuonTrigEffB2toI_EtaPhi_Barrel_MC->GetBinContent(bin);
       if(period==per2011J_M)bin = _MuonTrigEffJtoM_EtaPhi_Barrel_MC->FindBin(mu_CB_eta, mu_CB_phi);
       if(period==per2011J_M)eff = _MuonTrigEffJtoM_EtaPhi_Barrel_MC->GetBinContent(bin);

       if(period==per2011J)bin = _MuonTrigEffJ_EtaPhi_Barrel_MC->FindBin(mu_CB_eta, mu_CB_phi);
       if(period==per2011J)eff = _MuonTrigEffJ_EtaPhi_Barrel_MC->GetBinContent(bin);
       if(period==per2011K)bin = _MuonTrigEffK_EtaPhi_Barrel_MC->FindBin(mu_CB_eta, mu_CB_phi);
       if(period==per2011K)eff = _MuonTrigEffK_EtaPhi_Barrel_MC->GetBinContent(bin);
       if(period==per2011K1)bin = _MuonTrigEffK1_EtaPhi_Barrel_MC->FindBin(mu_CB_eta, mu_CB_phi);
       if(period==per2011K1)eff = _MuonTrigEffK1_EtaPhi_Barrel_MC->GetBinContent(bin);
       if(period==per2011K2_K6)bin = _MuonTrigEffK2_K6_EtaPhi_Barrel_MC->FindBin(mu_CB_eta, mu_CB_phi);
       if(period==per2011K2_K6)eff = _MuonTrigEffK2_K6_EtaPhi_Barrel_MC->GetBinContent(bin);
       if(period==per2011L)bin = _MuonTrigEffL_EtaPhi_Barrel_MC->FindBin(mu_CB_eta, mu_CB_phi);
       if(period==per2011L)eff = _MuonTrigEffL_EtaPhi_Barrel_MC->GetBinContent(bin);
       if(period==per2011L3_L4)bin = _MuonTrigEffL3_L4_EtaPhi_Barrel_MC->FindBin(mu_CB_eta, mu_CB_phi);
       if(period==per2011L3_L4)eff = _MuonTrigEffL3_L4_EtaPhi_Barrel_MC->GetBinContent(bin);
       if(period==per2011LwoL3_L4)bin = _MuonTrigEffLwoL3_L4_EtaPhi_Barrel_MC->FindBin(mu_CB_eta, mu_CB_phi);
       if(period==per2011LwoL3_L4)eff = _MuonTrigEffLwoL3_L4_EtaPhi_Barrel_MC->GetBinContent(bin);
       if(period==per2011M)bin = _MuonTrigEffM_EtaPhi_Barrel_MC->FindBin(mu_CB_eta, mu_CB_phi);
       if(period==per2011M)eff = _MuonTrigEffM_EtaPhi_Barrel_MC->GetBinContent(bin);
    }
  }
  else { //endcap
    if(isData){ 
      if(period==per2011B2_I)bin = _MuonTrigEffB2toI_EtaPhi_Endcap_Data->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011B2_I)eff = _MuonTrigEffB2toI_EtaPhi_Endcap_Data->GetBinContent(bin);
      if(period==per2011J_M)bin = _MuonTrigEffJtoM_EtaPhi_Endcap_Data->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011J_M)eff = _MuonTrigEffJtoM_EtaPhi_Endcap_Data->GetBinContent(bin);

      if(period==per2011J)bin = _MuonTrigEffJ_EtaPhi_Endcap_Data->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011J)eff = _MuonTrigEffJ_EtaPhi_Endcap_Data->GetBinContent(bin);
      if(period==per2011K)bin = _MuonTrigEffK_EtaPhi_Endcap_Data->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011K)eff = _MuonTrigEffK_EtaPhi_Endcap_Data->GetBinContent(bin);
      if(period==per2011K1)bin = _MuonTrigEffK1_EtaPhi_Endcap_Data->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011K1)eff = _MuonTrigEffK1_EtaPhi_Endcap_Data->GetBinContent(bin);
      if(period==per2011K2_K6)bin = _MuonTrigEffK2_K6_EtaPhi_Endcap_Data->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011K2_K6)eff = _MuonTrigEffK2_K6_EtaPhi_Endcap_Data->GetBinContent(bin);
      if(period==per2011L)bin = _MuonTrigEffL_EtaPhi_Endcap_Data->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011L)eff = _MuonTrigEffL_EtaPhi_Endcap_Data->GetBinContent(bin);
      if(period==per2011L3_L4)bin = _MuonTrigEffL3_L4_EtaPhi_Endcap_Data->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011L3_L4)eff = _MuonTrigEffL3_L4_EtaPhi_Endcap_Data->GetBinContent(bin);
      if(period==per2011LwoL3_L4)bin = _MuonTrigEffLwoL3_L4_EtaPhi_Endcap_Data->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011LwoL3_L4)eff = _MuonTrigEffLwoL3_L4_EtaPhi_Endcap_Data->GetBinContent(bin);
      if(period==per2011M)bin = _MuonTrigEffM_EtaPhi_Endcap_Data->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011M)eff = _MuonTrigEffM_EtaPhi_Endcap_Data->GetBinContent(bin);      
    }
    else{
      if(period==per2011B2_I)bin = _MuonTrigEffB2toI_EtaPhi_Endcap_MC->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011B2_I)eff = _MuonTrigEffB2toI_EtaPhi_Endcap_MC->GetBinContent(bin);
      if(period==per2011J_M)bin = _MuonTrigEffJtoM_EtaPhi_Endcap_MC->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011J_M)eff = _MuonTrigEffJtoM_EtaPhi_Endcap_MC->GetBinContent(bin);

      if(period==per2011J)bin = _MuonTrigEffJ_EtaPhi_Endcap_MC->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011J)eff = _MuonTrigEffJ_EtaPhi_Endcap_MC->GetBinContent(bin);
      if(period==per2011K)bin = _MuonTrigEffK_EtaPhi_Endcap_MC->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011K)eff = _MuonTrigEffK_EtaPhi_Endcap_MC->GetBinContent(bin);
      if(period==per2011K1)bin = _MuonTrigEffK1_EtaPhi_Endcap_MC->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011K1)eff = _MuonTrigEffK1_EtaPhi_Endcap_MC->GetBinContent(bin);
      if(period==per2011K2_K6)bin = _MuonTrigEffK2_K6_EtaPhi_Endcap_MC->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011K2_K6)eff = _MuonTrigEffK2_K6_EtaPhi_Endcap_MC->GetBinContent(bin);
      if(period==per2011L)bin = _MuonTrigEffL_EtaPhi_Endcap_MC->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011L)eff = _MuonTrigEffL_EtaPhi_Endcap_MC->GetBinContent(bin);
      if(period==per2011L3_L4)bin = _MuonTrigEffL3_L4_EtaPhi_Endcap_MC->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011L3_L4)eff = _MuonTrigEffL3_L4_EtaPhi_Endcap_MC->GetBinContent(bin);
      if(period==per2011LwoL3_L4)bin = _MuonTrigEffLwoL3_L4_EtaPhi_Endcap_MC->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011LwoL3_L4)eff = _MuonTrigEffLwoL3_L4_EtaPhi_Endcap_MC->GetBinContent(bin);
      if(period==per2011M)bin = _MuonTrigEffM_EtaPhi_Endcap_MC->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011M)eff = _MuonTrigEffM_EtaPhi_Endcap_MC->GetBinContent(bin);      
    }
  }

  return eff;
}

std::pair<double, double>  LeptonTriggerSF::MuEffErr_CB(SFDataPeriod period,bool isData, TLorentzVector muon) const{

  const double mu_CB_eta=muon.Eta();
  double mu_CB_phi= check_Phi_range(muon.Phi());  
  //double mu_pt=muon.Pt();
  //if(!inGeV) mu_pt=mu_pt/1000.;

  Int_t bin = -1;
  double eff = 0;
  double stat = 0;

  if(mu_CB_phi < m_minimum_phi) mu_CB_phi+=2*pi; // fix phi range for unusual binning
  
  if(fabs(mu_CB_eta) < 1.05) {  //barrel
    if(isData){
      if(period==per2011B2_I)bin = _MuonTrigEffB2toI_EtaPhi_Barrel_Data->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011B2_I)eff = _MuonTrigEffB2toI_EtaPhi_Barrel_Data->GetBinContent(bin);
      if(period==per2011B2_I)stat = _MuonTrigEffB2toI_EtaPhi_Barrel_Data->GetBinError(bin);
      if(period==per2011J_M)bin = _MuonTrigEffJtoM_EtaPhi_Barrel_Data->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011J_M)eff = _MuonTrigEffJtoM_EtaPhi_Barrel_Data->GetBinContent(bin);
      if(period==per2011J_M)stat = _MuonTrigEffJtoM_EtaPhi_Barrel_Data->GetBinError(bin);

      if(period==per2011J)bin = _MuonTrigEffJ_EtaPhi_Barrel_Data->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011J)eff = _MuonTrigEffJ_EtaPhi_Barrel_Data->GetBinContent(bin);
      if(period==per2011J)stat = _MuonTrigEffJ_EtaPhi_Barrel_Data->GetBinError(bin);
      if(period==per2011K)bin = _MuonTrigEffK_EtaPhi_Barrel_Data->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011K)eff = _MuonTrigEffK_EtaPhi_Barrel_Data->GetBinContent(bin);
      if(period==per2011K)stat = _MuonTrigEffK_EtaPhi_Barrel_Data->GetBinError(bin);
      if(period==per2011K1)bin = _MuonTrigEffK1_EtaPhi_Barrel_Data->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011K1)eff = _MuonTrigEffK1_EtaPhi_Barrel_Data->GetBinContent(bin);
      if(period==per2011K1)stat = _MuonTrigEffK1_EtaPhi_Barrel_Data->GetBinError(bin);
      if(period==per2011K2_K6)bin = _MuonTrigEffK2_K6_EtaPhi_Barrel_Data->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011K2_K6)eff = _MuonTrigEffK2_K6_EtaPhi_Barrel_Data->GetBinContent(bin);
      if(period==per2011K2_K6)stat = _MuonTrigEffK2_K6_EtaPhi_Barrel_Data->GetBinError(bin);
      if(period==per2011L)bin = _MuonTrigEffL_EtaPhi_Barrel_Data->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011L)eff = _MuonTrigEffL_EtaPhi_Barrel_Data->GetBinContent(bin);
      if(period==per2011L)stat = _MuonTrigEffL_EtaPhi_Barrel_Data->GetBinError(bin);
      if(period==per2011L3_L4)bin = _MuonTrigEffL3_L4_EtaPhi_Barrel_Data->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011L3_L4)eff = _MuonTrigEffL3_L4_EtaPhi_Barrel_Data->GetBinContent(bin);
      if(period==per2011L3_L4)stat = _MuonTrigEffL3_L4_EtaPhi_Barrel_Data->GetBinError(bin);
      if(period==per2011LwoL3_L4)bin = _MuonTrigEffLwoL3_L4_EtaPhi_Barrel_Data->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011LwoL3_L4)eff = _MuonTrigEffLwoL3_L4_EtaPhi_Barrel_Data->GetBinContent(bin);
      if(period==per2011LwoL3_L4)stat = _MuonTrigEffLwoL3_L4_EtaPhi_Barrel_Data->GetBinError(bin);
      if(period==per2011M)bin = _MuonTrigEffM_EtaPhi_Barrel_Data->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011M)eff = _MuonTrigEffM_EtaPhi_Barrel_Data->GetBinContent(bin);
      if(period==per2011M)stat = _MuonTrigEffM_EtaPhi_Barrel_Data->GetBinError(bin);
    }
    else{
      if(period==per2011B2_I)bin = _MuonTrigEffB2toI_EtaPhi_Barrel_MC->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011B2_I)eff = _MuonTrigEffB2toI_EtaPhi_Barrel_MC->GetBinContent(bin);
      if(period==per2011B2_I)stat = _MuonTrigEffB2toI_EtaPhi_Barrel_MC->GetBinError(bin);
      if(period==per2011J_M)bin = _MuonTrigEffJtoM_EtaPhi_Barrel_MC->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011J_M)eff = _MuonTrigEffJtoM_EtaPhi_Barrel_MC->GetBinContent(bin);
      if(period==per2011J_M)stat = _MuonTrigEffJtoM_EtaPhi_Barrel_MC->GetBinError(bin);

      if(period==per2011J)bin = _MuonTrigEffJ_EtaPhi_Barrel_MC->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011J)eff = _MuonTrigEffJ_EtaPhi_Barrel_MC->GetBinContent(bin);
      if(period==per2011J)stat = _MuonTrigEffJ_EtaPhi_Barrel_MC->GetBinError(bin);
      if(period==per2011K)bin = _MuonTrigEffK_EtaPhi_Barrel_MC->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011K)eff = _MuonTrigEffK_EtaPhi_Barrel_MC->GetBinContent(bin);
      if(period==per2011K)stat = _MuonTrigEffK_EtaPhi_Barrel_MC->GetBinError(bin);
      if(period==per2011K1)bin = _MuonTrigEffK1_EtaPhi_Barrel_MC->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011K1)eff = _MuonTrigEffK1_EtaPhi_Barrel_MC->GetBinContent(bin);
      if(period==per2011K1)stat = _MuonTrigEffK1_EtaPhi_Barrel_MC->GetBinError(bin);
      if(period==per2011K2_K6)bin = _MuonTrigEffK2_K6_EtaPhi_Barrel_MC->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011K2_K6)eff = _MuonTrigEffK2_K6_EtaPhi_Barrel_MC->GetBinContent(bin);
      if(period==per2011K2_K6)stat = _MuonTrigEffK2_K6_EtaPhi_Barrel_MC->GetBinError(bin);
      if(period==per2011L)bin = _MuonTrigEffL_EtaPhi_Barrel_MC->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011L)eff = _MuonTrigEffL_EtaPhi_Barrel_MC->GetBinContent(bin);
      if(period==per2011L)stat = _MuonTrigEffL_EtaPhi_Barrel_MC->GetBinError(bin);
      if(period==per2011L3_L4)bin = _MuonTrigEffL3_L4_EtaPhi_Barrel_MC->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011L3_L4)eff = _MuonTrigEffL3_L4_EtaPhi_Barrel_MC->GetBinContent(bin);
      if(period==per2011L3_L4)stat = _MuonTrigEffL3_L4_EtaPhi_Barrel_MC->GetBinError(bin);
      if(period==per2011LwoL3_L4)bin = _MuonTrigEffLwoL3_L4_EtaPhi_Barrel_MC->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011LwoL3_L4)eff = _MuonTrigEffLwoL3_L4_EtaPhi_Barrel_MC->GetBinContent(bin);
      if(period==per2011LwoL3_L4)stat = _MuonTrigEffLwoL3_L4_EtaPhi_Barrel_MC->GetBinError(bin);
      if(period==per2011M)bin = _MuonTrigEffM_EtaPhi_Barrel_MC->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011M)eff = _MuonTrigEffM_EtaPhi_Barrel_MC->GetBinContent(bin);
      if(period==per2011M)stat = _MuonTrigEffM_EtaPhi_Barrel_MC->GetBinError(bin);
    }
  }
  else { //endcap
    if(isData){ 
      if(period==per2011B2_I)bin = _MuonTrigEffB2toI_EtaPhi_Endcap_Data->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011B2_I)eff = _MuonTrigEffB2toI_EtaPhi_Endcap_Data->GetBinContent(bin);
      if(period==per2011B2_I)stat = _MuonTrigEffB2toI_EtaPhi_Endcap_Data->GetBinError(bin);
      if(period==per2011J_M)bin = _MuonTrigEffJtoM_EtaPhi_Endcap_Data->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011J_M)eff = _MuonTrigEffJtoM_EtaPhi_Endcap_Data->GetBinContent(bin);
      if(period==per2011J_M)stat = _MuonTrigEffJtoM_EtaPhi_Endcap_Data->GetBinError(bin);

      if(period==per2011J)bin = _MuonTrigEffJ_EtaPhi_Endcap_Data->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011J)eff = _MuonTrigEffJ_EtaPhi_Endcap_Data->GetBinContent(bin);
      if(period==per2011J)stat = _MuonTrigEffJ_EtaPhi_Endcap_Data->GetBinError(bin);
      if(period==per2011K)bin = _MuonTrigEffK_EtaPhi_Endcap_Data->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011K)eff = _MuonTrigEffK_EtaPhi_Endcap_Data->GetBinContent(bin);
      if(period==per2011K)stat = _MuonTrigEffK_EtaPhi_Endcap_Data->GetBinError(bin);
      if(period==per2011K1)bin = _MuonTrigEffK1_EtaPhi_Endcap_Data->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011K1)eff = _MuonTrigEffK1_EtaPhi_Endcap_Data->GetBinContent(bin);
      if(period==per2011K1)stat = _MuonTrigEffK1_EtaPhi_Endcap_Data->GetBinError(bin);
      if(period==per2011K2_K6)bin = _MuonTrigEffK2_K6_EtaPhi_Endcap_Data->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011K2_K6)eff = _MuonTrigEffK2_K6_EtaPhi_Endcap_Data->GetBinContent(bin);
      if(period==per2011K2_K6)stat = _MuonTrigEffK2_K6_EtaPhi_Endcap_Data->GetBinError(bin);
      if(period==per2011L)bin = _MuonTrigEffL_EtaPhi_Endcap_Data->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011L)eff = _MuonTrigEffL_EtaPhi_Endcap_Data->GetBinContent(bin);
      if(period==per2011L)stat = _MuonTrigEffL_EtaPhi_Endcap_Data->GetBinError(bin);
      if(period==per2011L3_L4)bin = _MuonTrigEffL3_L4_EtaPhi_Endcap_Data->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011L3_L4)eff = _MuonTrigEffL3_L4_EtaPhi_Endcap_Data->GetBinContent(bin);
      if(period==per2011L3_L4)stat = _MuonTrigEffL3_L4_EtaPhi_Endcap_Data->GetBinError(bin);
      if(period==per2011LwoL3_L4)bin = _MuonTrigEffLwoL3_L4_EtaPhi_Endcap_Data->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011LwoL3_L4)eff = _MuonTrigEffLwoL3_L4_EtaPhi_Endcap_Data->GetBinContent(bin);
      if(period==per2011LwoL3_L4)stat = _MuonTrigEffLwoL3_L4_EtaPhi_Endcap_Data->GetBinError(bin);
      if(period==per2011M)bin = _MuonTrigEffM_EtaPhi_Endcap_Data->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011M)eff = _MuonTrigEffM_EtaPhi_Endcap_Data->GetBinContent(bin);
      if(period==per2011M)stat = _MuonTrigEffM_EtaPhi_Endcap_Data->GetBinError(bin);
    }
    else{
      if(period==per2011B2_I)bin = _MuonTrigEffB2toI_EtaPhi_Endcap_MC->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011B2_I)eff = _MuonTrigEffB2toI_EtaPhi_Endcap_MC->GetBinContent(bin);
      if(period==per2011B2_I)stat = _MuonTrigEffB2toI_EtaPhi_Endcap_MC->GetBinError(bin);
      if(period==per2011J_M)bin = _MuonTrigEffJtoM_EtaPhi_Endcap_MC->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011J_M)eff = _MuonTrigEffJtoM_EtaPhi_Endcap_MC->GetBinContent(bin);
      if(period==per2011J_M)stat = _MuonTrigEffJtoM_EtaPhi_Endcap_MC->GetBinError(bin);

      if(period==per2011J)bin = _MuonTrigEffJ_EtaPhi_Endcap_MC->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011J)eff = _MuonTrigEffJ_EtaPhi_Endcap_MC->GetBinContent(bin);
      if(period==per2011J)stat = _MuonTrigEffJ_EtaPhi_Endcap_MC->GetBinError(bin);
      if(period==per2011K)bin = _MuonTrigEffK_EtaPhi_Endcap_MC->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011K)eff = _MuonTrigEffK_EtaPhi_Endcap_MC->GetBinContent(bin);
      if(period==per2011K)stat = _MuonTrigEffK_EtaPhi_Endcap_MC->GetBinError(bin);
      if(period==per2011K1)bin = _MuonTrigEffK1_EtaPhi_Endcap_MC->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011K1)eff = _MuonTrigEffK1_EtaPhi_Endcap_MC->GetBinContent(bin);
      if(period==per2011K1)stat = _MuonTrigEffK1_EtaPhi_Endcap_MC->GetBinError(bin);
      if(period==per2011K2_K6)bin = _MuonTrigEffK2_K6_EtaPhi_Endcap_MC->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011K2_K6)eff = _MuonTrigEffK2_K6_EtaPhi_Endcap_MC->GetBinContent(bin);
      if(period==per2011K2_K6)stat = _MuonTrigEffK2_K6_EtaPhi_Endcap_MC->GetBinError(bin);
      if(period==per2011L)bin = _MuonTrigEffL_EtaPhi_Endcap_MC->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011L)eff = _MuonTrigEffL_EtaPhi_Endcap_MC->GetBinContent(bin);
      if(period==per2011L)stat = _MuonTrigEffL_EtaPhi_Endcap_MC->GetBinError(bin);
      if(period==per2011L3_L4)bin = _MuonTrigEffL3_L4_EtaPhi_Endcap_MC->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011L3_L4)eff = _MuonTrigEffL3_L4_EtaPhi_Endcap_MC->GetBinContent(bin);
      if(period==per2011L3_L4)stat = _MuonTrigEffL3_L4_EtaPhi_Endcap_MC->GetBinError(bin);
      if(period==per2011LwoL3_L4)bin = _MuonTrigEffLwoL3_L4_EtaPhi_Endcap_MC->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011LwoL3_L4)eff = _MuonTrigEffLwoL3_L4_EtaPhi_Endcap_MC->GetBinContent(bin);
      if(period==per2011LwoL3_L4)stat = _MuonTrigEffLwoL3_L4_EtaPhi_Endcap_MC->GetBinError(bin);
      if(period==per2011M)bin = _MuonTrigEffM_EtaPhi_Endcap_MC->FindBin(mu_CB_eta, mu_CB_phi);
      if(period==per2011M)eff = _MuonTrigEffM_EtaPhi_Endcap_MC->GetBinContent(bin);
      if(period==per2011M)stat = _MuonTrigEffM_EtaPhi_Endcap_MC->GetBinError(bin);
    }
  }

  double syst = commonSystMTSG*eff;
  double err  = sqrt(stat*stat + syst*syst);

  return std::pair<double, double> ( eff, err );

}

double LeptonTriggerSF::MuEff_Tight(SFDataPeriod period,bool isData, TLorentzVector muon) const{

  const double mu_Tight_eta=muon.Eta();
  double mu_Tight_phi=check_Phi_range(muon.Phi());  
  //double mu_pt=muon.Pt();
  //if(!inGeV) mu_pt=mu_pt/1000.;

  Int_t bin = -1;
  double eff = 0;

  if(mu_Tight_phi < m_minimum_phi) mu_Tight_phi+=2*pi; // fix phi range for unusual binning

  if(fabs(mu_Tight_eta) < 1.05) {  //barrel
    if(isData){
      if(period==per2011B2_I)bin = _MuonTrigEffB2toI_EtaPhi_Barrel_Data->FindBin(mu_Tight_eta, mu_Tight_phi);
      if(period==per2011B2_I)eff = _MuonTrigEffB2toI_EtaPhi_Barrel_Data->GetBinContent(bin);
      if(period==per2011J_M)bin = _MuonTrigEffJtoM_EtaPhi_Barrel_Data->FindBin(mu_Tight_eta, mu_Tight_phi);
      if(period==per2011J_M)eff = _MuonTrigEffJtoM_EtaPhi_Barrel_Data->GetBinContent(bin);
    }
    else{
       if(period==per2011B2_I)bin = _MuonTrigEffB2toI_EtaPhi_Barrel_MC->FindBin(mu_Tight_eta, mu_Tight_phi);
       if(period==per2011B2_I)eff = _MuonTrigEffB2toI_EtaPhi_Barrel_MC->GetBinContent(bin);
       if(period==per2011J_M)bin = _MuonTrigEffJtoM_EtaPhi_Barrel_MC->FindBin(mu_Tight_eta, mu_Tight_phi);
       if(period==per2011J_M)eff = _MuonTrigEffJtoM_EtaPhi_Barrel_MC->GetBinContent(bin);
    }
  }
  else { //endcap
    if(isData){ 
       if(period==per2011B2_I)bin = _MuonTrigEffB2toI_EtaPhi_Endcap_Data->FindBin(mu_Tight_eta, mu_Tight_phi);
       if(period==per2011B2_I)eff = _MuonTrigEffB2toI_EtaPhi_Endcap_Data->GetBinContent(bin);
       if(period==per2011J_M)bin = _MuonTrigEffJtoM_EtaPhi_Endcap_Data->FindBin(mu_Tight_eta, mu_Tight_phi);
       if(period==per2011J_M)eff = _MuonTrigEffJtoM_EtaPhi_Endcap_Data->GetBinContent(bin);
    }
    else{
      if(period==per2011B2_I)bin = _MuonTrigEffB2toI_EtaPhi_Endcap_MC->FindBin(mu_Tight_eta, mu_Tight_phi);
      if(period==per2011B2_I)eff = _MuonTrigEffB2toI_EtaPhi_Endcap_MC->GetBinContent(bin);
      if(period==per2011J_M)bin = _MuonTrigEffJtoM_EtaPhi_Endcap_MC->FindBin(mu_Tight_eta, mu_Tight_phi);
      if(period==per2011J_M)eff = _MuonTrigEffJtoM_EtaPhi_Endcap_MC->GetBinContent(bin);
    }
  }

  return eff;
}


double LeptonTriggerSF::EmEff_Data(TLorentzVector electron) const{

  const double em_eta=electron.Eta();
  
  double m_trigger_em_Eff_data = ele_WZ_trigger_Eff_data(em_eta);
  return m_trigger_em_Eff_data;
}

double LeptonTriggerSF::EmEff_MC(TLorentzVector electron) const{

  const double em_eta=electron.Eta();
 
  double m_trigger_em_Eff_mc = ele_WZ_trigger_Eff_MC(em_eta);
  return m_trigger_em_Eff_mc;
}

void LeptonTriggerSF::setThresholds(bool useGeV){

  if(useGeV){
    inGeV=true;
    muThreshold=20.;
    emThreshold=25.;
  }
  return;
}

SFDataPeriod LeptonTriggerSF::getDataPeriod(Int_t runNumber){

  if(runNumber<1000000) {
    if(runNumber<178044) return perUnDefined;
    if(runNumber<=186493) return per2011B2_I;
    if(runNumber<=191933) return per2011J_M;
  } else {
    if(runNumber == (1000000 + 1)) return per2011J;
    if(runNumber == (1000000 + 2)) return per2011K;
    if(runNumber == (1000000 + 3)) return per2011K1;
    if(runNumber == (1000000 + 4)) return per2011K2_K6;
    if(runNumber == (1000000 + 5)) return per2011L;
    if(runNumber == (1000000 + 6)) return per2011L3_L4;
    if(runNumber == (1000000 + 7)) return per2011LwoL3_L4;
    if(runNumber == (1000000 + 8)) return per2011M;
  }
  

  return perUnDefined;
}

// Helper function to deal with possible phi ambiguity
double LeptonTriggerSF::check_Phi_range(double phi) const  {

  double newphi = phi;
  const double pi = acos(-1.);
  
  if (newphi > pi) {
    printf("<muon_SF>: WARNING: Muon phi %4.2f > pi! Using (phi-2*pi) \n", phi);
    newphi -= 2*pi;
  }
  if (newphi < -pi) {
    printf("<muon_SF>: WARNING: Muon phi %4.2f < -pi! Using (phi+2*pi) \n", phi);
    newphi += 2*pi;
  }
  
  return newphi;
}