VHbbDataFormats/interface/TriggerWeight.h

#ifndef TRIGGERWEIGHT_H
#define TRIGGERWEIGHT_H

#include "FWCore/ParameterSet/interface/ProcessDesc.h"
#include "FWCore/PythonParameterSet/interface/PythonProcessDesc.h"
//#include "VHbbAnalysis/VHbbDataFormats/interface/TriggerZnunuCurve.h"
#include "VHbbAnalysis/VHbbDataFormats/interface/MultiThresholdEfficiency.h"
#include <TH1F.h>
#include <TF1.h>
#include <TFile.h>
#include <TTree.h>
#include <iostream>

class TriggerWeight
{
public:  
  TriggerWeight(const edm::ParameterSet& ana) : combiner2Thr(2), combiner1Thr(1)
  {
   tscaleHLTmu=openFile(ana,"hltMuFileName");
   tscaleIDmu=openFile(ana,"idMuFileName");
   tscaleHLTele1=openFile(ana,"hltEle1FileName");
   tscaleHLTele2=openFile(ana,"hltEle2FileName");
   tscaleHLTele1Aug=openFile(ana,"hltEle1AugFileName");
   tscaleHLTele2Aug=openFile(ana,"hltEle2AugFileName");
   tscaleID80Ele=openFile(ana,"idEle80FileName");
   tscaleID95Ele=openFile(ana,"idEle95FileName");
   tscaleHLTeleJet1=openFile(ana,"hltJetEle1FileName");
   tscaleHLTeleJet2=openFile(ana,"hltJetEle2FileName");
   tscaleRecoEle=openFile(ana,"recoEleFileName");
//   tscalePFMHTele=openFile(ana,"hltPFMHTEleFileName");
   tscaleSingleEleMay=openFile(ana,"hltSingleEleMayFileName");
   tscaleSingleEleV4=openFile(ana,"hltSingleEleV4FileName");
   tscaleHLTmuOr30=openFile(ana,"hltMuOr30FileName");
   tscaleHLTmuonCrossTrig = openFile(ana,"hltMuCrossTrig");
   tscaleHLTelectronCrossTrig = openFile(ana,"hltEleCrossTrig");

   tscaleSingleEle2012Awp95=openFile(ana,"hltSingleEle2012Awp95");
   tscaleSingleEle2012Awp80=openFile(ana,"hltSingleEle2012Awp80");
   tscaleSingleMuon2012A=openFile(ana,"hltSingleMuon2012A");
   tscaleDoubleEle2012A_leg8=openFile(ana,"hltDoubleEle2012A_leg8");
   tscaleDoubleEle2012A_leg17=openFile(ana,"hltDoubleEle2012A_leg17");
   tscaleDoubleMuon2012A_leg8=openFile(ana,"hltDoubleMuon2012A_leg8");
   tscaleDoubleMuon2012A_leg17=openFile(ana,"hltDoubleMuon2012A_leg17");

   tscaleDoubleMuon2012A_dZ=openFile(ana,"hltDoubleMuon2012A_dZ");
   tscaleDoubleEle2012A_dZ=openFile(ana,"hltDoubleEle2012A_dZ");

   tscaleMuPlusWCandPt2012A_legMu=openFile(ana,"hltMuPlusWCandPt2012A_legMu");  
   tscaleMuPlusWCandPt2012A_legW=openFile(ana,"hltMuPlusWCandPt2012A_legW");  
   
   tscaleMuID2012A=openFile(ana,"idMu2012A");  
   tscaleEleID2012A=openFile(ana,"idEle2012A");  
   tscaleEleID2012Awp80=openFile(ana,"idEle2012Awp80");  

   if(tscaleHLTmu == 0 || tscaleIDmu == 0) 
    {
      std::cout << "ERROR: cannot load Muon Trigger efficiencies" << std::endl;
    } 
   

  }
 
static TTree * openFile(const edm::ParameterSet& ana, const char * name)
  {
   TFile *hltMuFile = new TFile (ana.getParameter<std::string> (name).c_str(),"read");
   if(hltMuFile)   return (TTree*) hltMuFile->Get("tree");
   else return 0;
  }

static  std::pair<float,float> efficiencyFromPtEta(float pt1, float eta1, TTree *t)
  {
//     std::cout << "here " << t << " pt1 "  << pt1 << " eta1 " << eta1 <<  std::endl;
     float s1 = 1.,err=1.;
     std::pair<float,float> r(s1,err);
     if(!t) return r;
     float ptMin,ptMax,etaMin,etaMax,scale,error;
     int count = 0;
     t->SetBranchAddress("ptMin",&ptMin);
     t->SetBranchAddress("ptMax",&ptMax);
     t->SetBranchAddress("etaMin",&etaMin);
     t->SetBranchAddress("etaMax",&etaMax);
     t->SetBranchAddress("scale",&scale);
     t->SetBranchAddress("error",&error);
    float lastPtBin = 200;
/*    for(int jentry = 0; jentry < t->GetEntries(); jentry++)
      {
        t->GetEntry(jentry);
         if(ptMax >= lastPtBin) lastPtBin =ptMax;
      }*/
    for(int jentry = 0; jentry < t->GetEntries(); jentry++)
      {
        t->GetEntry(jentry);
        if(ptMax==lastPtBin) ptMax=1e99;
        if((pt1 > ptMin) && (pt1 < ptMax) && (eta1 > etaMin) && (eta1 < etaMax))
          {
            s1 = scale;
            err=error;
            count++;
          }
      }

    if(count == 0)
      {
        return r;
      }

     r.first=s1;
     r.second = err;
    return (r);
  }
 
  float scaleMuIsoHLT(float pt1, float eta1)
  {
    return efficiencyFromPtEta(pt1,eta1,tscaleHLTmu).first;
  }
  
  
  float scaleMuID(float pt1, float eta1)
  {
    return efficiencyFromPtEta(pt1,eta1,tscaleIDmu).first;
  }

        float scaleMuCrossTrig(float pt1, float eta1)
        {
                return efficiencyFromPtEta(pt1,eta1,tscaleHLTmuonCrossTrig).first;
        }
        
        float scaleEleCrossTrig(float pt1, float eta1)
        {
                return efficiencyFromPtEta(pt1,eta1,tscaleHLTelectronCrossTrig).first;
        }
        

double scaleDoubleEle17Ele8Aug( std::vector<float> pt, std::vector<float> eta )
{
   std::vector< std::vector<float> > allEleWithEffs;
for(unsigned int j=0; j< pt.size() && j < 10 ; j++)
 {
  std::vector<float> thisEleEffs;
  thisEleEffs.push_back(efficiencyFromPtEta(pt[j],eta[j],tscaleHLTele1Aug).first);
  thisEleEffs.push_back(efficiencyFromPtEta(pt[j],eta[j],tscaleHLTele2Aug).first);

  allEleWithEffs.push_back(thisEleEffs);
 }

  return   combiner2Thr.weight<Trigger1High2Loose>(allEleWithEffs);

}


double scaleDoubleEle17Ele8( std::vector<float> pt, std::vector<float> eta )
{
   std::vector< std::vector<float> > allEleWithEffs;
for(unsigned int j=0; j< pt.size() && j<10; j++)
 {
  std::vector<float> thisEleEffs;
  thisEleEffs.push_back(efficiencyFromPtEta(pt[j],eta[j],tscaleHLTele1).first);
  thisEleEffs.push_back(efficiencyFromPtEta(pt[j],eta[j],tscaleHLTele2).first);

  allEleWithEffs.push_back(thisEleEffs);
 }

  return   combiner2Thr.weight<Trigger1High2Loose>(allEleWithEffs);

}
double muId2012A( float pt, float eta){    return efficiencyFromPtEta(pt,eta,tscaleMuID2012A).first;}
double eleId2012A( float pt, float eta){    return efficiencyFromPtEta(pt,eta,tscaleEleID2012A).first;}
double eleId2012Awp80( float pt, float eta){    return efficiencyFromPtEta(pt,eta,tscaleEleID2012Awp80).first;}

double scaleSingleEleMay( float pt, float eta){    return efficiencyFromPtEta(pt,eta,tscaleSingleEleMay).first;}
double scaleSingleEleV4( float pt, float eta){   return efficiencyFromPtEta(pt,eta,tscaleSingleEleV4).first; }
double scaleID80Ele( float pt, float eta) {      return efficiencyFromPtEta(pt,eta,tscaleID80Ele).first; }
double scaleID95Ele( float pt, float eta) {      return efficiencyFromPtEta(pt,eta,tscaleID95Ele).first; }
double scaleRecoEle( float pt, float eta){     return efficiencyFromPtEta(pt,eta,tscaleRecoEle).first; }
double scalePFMHTEle( float MetPFPt){
    double weightPFMHTrigger=0.;

    //FIXME: read from file
    if(MetPFPt>0. && MetPFPt<5.) weightPFMHTrigger=0.305;
    if(MetPFPt>5. && MetPFPt<10.) weightPFMHTrigger=0.351;
    if(MetPFPt>10. && MetPFPt<15.) weightPFMHTrigger=0.461;
    if(MetPFPt>15. && MetPFPt<20.) weightPFMHTrigger=0.572;
    if(MetPFPt>20. && MetPFPt<25.) weightPFMHTrigger=0.713;
    if(MetPFPt>25. && MetPFPt<30.) weightPFMHTrigger=0.844;
    if(MetPFPt>30. && MetPFPt<35.) weightPFMHTrigger=0.914;
    if(MetPFPt>35. && MetPFPt<40.) weightPFMHTrigger=0.939;
    if(MetPFPt>40. && MetPFPt<45.) weightPFMHTrigger=0.981;
    if(MetPFPt>45. && MetPFPt<50.) weightPFMHTrigger=0.982;
    if(MetPFPt>50. && MetPFPt<60.) weightPFMHTrigger=0.993;
    if(MetPFPt>60. && MetPFPt<70.) weightPFMHTrigger=0.995;
    if(MetPFPt>70. && MetPFPt<100.) weightPFMHTrigger=0.995;
    if(MetPFPt>100.) weightPFMHTrigger=1.;
    return weightPFMHTrigger;
}


double scaleJet30Jet25( std::vector<float> pt, std::vector<float> eta)
{

  std::vector< std::vector<float> > allJetsWithEffs;
for(unsigned int j=0; j< pt.size() && j < 10 ; j++)
 {
  std::vector<float> thisJetEffs;
  thisJetEffs.push_back(efficiencyFromPtEta(pt[j],eta[j],tscaleHLTeleJet1).first);
  thisJetEffs.push_back(efficiencyFromPtEta(pt[j],eta[j],tscaleHLTeleJet2).first);
//  std::cout << " jet pt " << pt[j] << " eta " << eta[j] << " eff1 "  <<  thisJetEffs[0] << " eff2 " << thisJetEffs[1] << std::endl;  
  allJetsWithEffs.push_back(thisJetEffs);
  
 }
 float res = combiner2Thr.weight<Trigger1High2Loose>(allJetsWithEffs);
// std::cout << "Result is " << res << std::endl;
 return res;
//  return   combiner2Thr.weight<Trigger1High2Loose>(allJetsWithEffs);

}
/*
TF1 fpt("f","1-exp(-0.157*(x-19.3))", 0., 9999999.);

MET80:
TF1 fmet80("f","1/ (1 + exp( -0.0709 * (x - 100.7)))", 0., 9999999.);

MET100:
TF1 fmet100("f","1/ (1 + exp( -0.0679 * (x - 128.8)))", 0., 9999999.);
*/

//LP curve used for MET
double  scaleMetHLT( double met){
    return 1. / (1. +  ( exp( 0.059486 * ( 123.27 - met ))));
}

//MET80 component of the factorized JET+MET trigger
double scaleMET80(double et)
{
 return 1. / (1. +  exp( -0.0709 * (et - 100.7)));
}

//MET100 component
double scaleMET100(double et)
{
 return 1. / (1. +  exp( -0.0679 * (et - 128.8)));
}

//Single jet20 efficiency for MET+2CJet20
double jet20efficiency( double pt)
{
 if(pt < 10 ) return 0;
 return 1. - exp(-0.157*(pt-19.3));
}

//combined 2 jets efficiency out of N jets, using jet20 efficiency curve
double scale2CentralJet( std::vector<float> pt, std::vector<float> eta)
{

  std::vector< std::vector<float> > allJetsWithEffs;
for(unsigned int j=0; j< pt.size() && j < 10; j++)
 {
  if(fabs(eta[j]) < 2.5)
  {
   std::vector<float> thisJetEffs;
   thisJetEffs.push_back(jet20efficiency(pt[j]));
   allJetsWithEffs.push_back(thisJetEffs);
  }

 }

 return   combiner1Thr.weight<Trigger2SingleThr>(allJetsWithEffs);
}

//New MET 150 
double scaleMET150(double et)
{ 
  return 1./ (1. +  exp( -0.129226 * (et - 156.699)));
}

float scaleMuOr30IsoHLT(float pt1, float eta1)
{
   return efficiencyFromPtEta(pt1,eta1,tscaleHLTmuOr30).first;
}


//For 2012A HLT_DiCentralPFJet30_PFMHT80, valid for pfMET > 100 GeV:
double scaleDiJet30MHT80_2012A(double x)
{
  if(x<100) return 0;
  return   (1e0 - exp(-0.04197*(x-75.73))) * 0.9721 ;
}
//For 2012B HLT_DiCentralJetSumpT100_dPhi05_DiCentralPFJet60_25_PFMET100_HBHENoiseCleaned, valid for pfMET > 100 GeV:
double scaleSumpT100MET100_2012B(double x)
{
  if(x<100) return 0;
  return   (1e0 - exp(-0.06704*(x-96.84))) * 0.9199 ;
}
//For 2012A+B HLT_PFMET150, valid for pfMET > 150 GeV:
double scalePFMET150_2012AB(double x)
{
  if(x<150) return 0;
  return  (1e0 - exp(-0.07135*(x-147.4))) * 0.9707;
}


//For 2012A HLT_PFMET150 OR HLT_DiCentralPFJet30_PFMHT80, valid for pfMET > 100 GeV:
double scalePFMET150orDiJetMET_2012A(double x)
{
  if(x<100) return 0;
  return  (1e0 - exp(-0.0412*(x-75.52))) * 0.9772;
}

//For 2012B HLT_PFMET150 OR HLT_DiCentralJetSumpT100_dPhi05_DiCentralPFJet60_25_PFMET100_HBHENoiseCleaned, valid for pfMET > 100 GeV:
double scalePFMET150orDiJetMET_2012B(double x)
{
  if(x<100) return 0;
  return  (1e0 - exp(-0.05482*(x-95.59))) * 0.9702;
}


//For 2012C HLT_PFMET150 OR HLT_DiCentralJetSumpT100_dPhi05_DiCentralPFJet60_25_PFMET100_HBHENoiseCleaned, valid for pfMET > 100 GeV:
double scalePFMET150orDiJetMET_2012C(double x)
{
  if(x<100) return 0;
  return  (1e0 - exp(-0.05627*(x-95.15))) * 0.9659;
}


float doubleEle2012A( float pt1, float eta1, float pt2, float eta2)
{
// std::cout << "di ele" << std::endl;
 float eff1_17 = efficiencyFromPtEta(pt1,eta1,tscaleDoubleEle2012A_leg17).first;
 float eff2_17 = efficiencyFromPtEta(pt2,eta2,tscaleDoubleEle2012A_leg17).first;
 float eff1_8 = efficiencyFromPtEta(pt1,eta1,tscaleDoubleEle2012A_leg8).first;
 float eff2_8 = efficiencyFromPtEta(pt2,eta2,tscaleDoubleEle2012A_leg8).first;
// std::cout << tscaleDoubleEle2012A_dZ << std::endl;
 float eff_dz = efficiencyFromPtEta(eta1,eta2,tscaleDoubleEle2012A_dZ).first;  // despite the name pt,eta is actually eta1,eta2
 
 return (eff1_17 * eff2_8 + eff2_17 * eff1_8 - eff1_17*eff2_17)*eff_dz;

}
float doubleMuon2012A( float pt1, float eta1, float pt2, float eta2)
{
// std::cout << "di mu" << std::endl;
 float eff1_17 = efficiencyFromPtEta(pt1,eta1,tscaleDoubleMuon2012A_leg17).first;
 float eff2_17 = efficiencyFromPtEta(pt2,eta2,tscaleDoubleMuon2012A_leg17).first;
 float eff1_8 = efficiencyFromPtEta(pt1,eta1,tscaleDoubleMuon2012A_leg8).first;
 float eff2_8 = efficiencyFromPtEta(pt2,eta2,tscaleDoubleMuon2012A_leg8).first;
// std::cout << tscaleDoubleMuon2012A_dZ << std::endl;
 float eff_dz = efficiencyFromPtEta(eta1,eta2,tscaleDoubleMuon2012A_dZ).first;  // despite the name pt,eta is actually eta1,eta2

 return (eff1_17 * eff2_8 + eff2_17 * eff1_8 - eff1_17*eff2_17)*eff_dz;

}


float muPlusWCandPt2012A_legW( float pt, float eta){    return efficiencyFromPtEta(pt,eta,tscaleMuPlusWCandPt2012A_legW).first;}
float muPlusWCandPt2012A_legMu( float pt, float eta){    return efficiencyFromPtEta(pt,eta,tscaleMuPlusWCandPt2012A_legMu).first;}
float singleEle2012Awp80( float pt, float eta){    return efficiencyFromPtEta(pt,eta,tscaleSingleEle2012Awp80).first;}
float singleEle2012Awp95( float pt, float eta){    return efficiencyFromPtEta(pt,eta,tscaleSingleEle2012Awp95).first;}
float singleMuon2012A( float pt, float eta){    return efficiencyFromPtEta(pt,eta,tscaleSingleMuon2012A).first;}


private:
  TTree * tscaleHLTele1;
  TTree * tscaleHLTele2;
  TTree * tscaleHLTeleJet1;
  TTree * tscaleHLTeleJet2;
  TTree * tscaleID80Ele;
  TTree * tscaleID95Ele;
  TTree * tscaleRecoEle;
  TTree * tscaleHLTmuOr30;

  TTree * tscaleSingleEle2012Awp95;
  TTree * tscaleSingleEle2012Awp80;
  TTree * tscaleSingleMuon2012A;
  TTree * tscaleDoubleEle2012A_leg8;
  TTree * tscaleDoubleEle2012A_leg17;
  TTree * tscaleDoubleMuon2012A_leg8;
  TTree * tscaleDoubleMuon2012A_leg17;
  TTree * tscaleMuPlusWCandPt2012A_legMu;
  TTree * tscaleMuPlusWCandPt2012A_legW;
  TTree * tscaleDoubleEle2012A_dZ;
  TTree * tscaleDoubleMuon2012A_dZ;
//  TTree * tscalePFMHTele;
  TTree * tscaleSingleEleMay;
  TTree * tscaleSingleEleV4;

  TTree * tscaleHLTele1Aug;
  TTree * tscaleHLTele2Aug;

  TTree * tscaleMuID2012A;
  TTree * tscaleEleID2012A;
  TTree * tscaleEleID2012Awp80;

  TTree * tscaleHLTmu;
        TTree * tscaleIDmu;
        TTree * tscaleHLTmuonCrossTrig;
        TTree * tscaleHLTelectronCrossTrig;
  
  MultiThresholdEfficiency combiner2Thr; 
  MultiThresholdEfficiency combiner1Thr; 
};

#endif
Revision:	1.22.4.1
Committed:	Mon Feb 11 19:43:16 2013 UTC (12 years, 2 months ago) by wilken
Content type:	text/plain
Branch:	V42TauCandidate
Changes since 1.22:	+17 -1 lines
Log Message:	cross trigger
#	User	Rev	Content
1	arizzi	1.1	#ifndef TRIGGERWEIGHT_H
2			#define TRIGGERWEIGHT_H
3
4			#include "FWCore/ParameterSet/interface/ProcessDesc.h"
5			#include "FWCore/PythonParameterSet/interface/PythonProcessDesc.h"
6	arizzi	1.11	//#include "VHbbAnalysis/VHbbDataFormats/interface/TriggerZnunuCurve.h"
7	arizzi	1.4	#include "VHbbAnalysis/VHbbDataFormats/interface/MultiThresholdEfficiency.h"
8	arizzi	1.1	#include <TH1F.h>
9	arizzi	1.3	#include <TF1.h>
10	arizzi	1.1	#include <TFile.h>
11			#include <TTree.h>
12			#include <iostream>
13	wilken	1.22.4.1
14	arizzi	1.1	class TriggerWeight
15			{
16			public:
17	arizzi	1.11	TriggerWeight(const edm::ParameterSet& ana) : combiner2Thr(2), combiner1Thr(1)
18	arizzi	1.1	{
19	arizzi	1.4	tscaleHLTmu=openFile(ana,"hltMuFileName");
20			tscaleIDmu=openFile(ana,"idMuFileName");
21			tscaleHLTele1=openFile(ana,"hltEle1FileName");
22			tscaleHLTele2=openFile(ana,"hltEle2FileName");
23	arizzi	1.9	tscaleHLTele1Aug=openFile(ana,"hltEle1AugFileName");
24			tscaleHLTele2Aug=openFile(ana,"hltEle2AugFileName");
25	arizzi	1.5	tscaleID80Ele=openFile(ana,"idEle80FileName");
26			tscaleID95Ele=openFile(ana,"idEle95FileName");
27	arizzi	1.4	tscaleHLTeleJet1=openFile(ana,"hltJetEle1FileName");
28			tscaleHLTeleJet2=openFile(ana,"hltJetEle2FileName");
29	arizzi	1.5	tscaleRecoEle=openFile(ana,"recoEleFileName");
30			// tscalePFMHTele=openFile(ana,"hltPFMHTEleFileName");
31			tscaleSingleEleMay=openFile(ana,"hltSingleEleMayFileName");
32			tscaleSingleEleV4=openFile(ana,"hltSingleEleV4FileName");
33	arizzi	1.6	tscaleHLTmuOr30=openFile(ana,"hltMuOr30FileName");
34	wilken	1.22.4.1	tscaleHLTmuonCrossTrig = openFile(ana,"hltMuCrossTrig");
35			tscaleHLTelectronCrossTrig = openFile(ana,"hltEleCrossTrig");
36	arizzi	1.4
37	arizzi	1.12	tscaleSingleEle2012Awp95=openFile(ana,"hltSingleEle2012Awp95");
38			tscaleSingleEle2012Awp80=openFile(ana,"hltSingleEle2012Awp80");
39			tscaleSingleMuon2012A=openFile(ana,"hltSingleMuon2012A");
40			tscaleDoubleEle2012A_leg8=openFile(ana,"hltDoubleEle2012A_leg8");
41			tscaleDoubleEle2012A_leg17=openFile(ana,"hltDoubleEle2012A_leg17");
42			tscaleDoubleMuon2012A_leg8=openFile(ana,"hltDoubleMuon2012A_leg8");
43	arizzi	1.15	tscaleDoubleMuon2012A_leg17=openFile(ana,"hltDoubleMuon2012A_leg17");
44	arizzi	1.13
45			tscaleDoubleMuon2012A_dZ=openFile(ana,"hltDoubleMuon2012A_dZ");
46			tscaleDoubleEle2012A_dZ=openFile(ana,"hltDoubleEle2012A_dZ");
47	arizzi	1.12
48			tscaleMuPlusWCandPt2012A_legMu=openFile(ana,"hltMuPlusWCandPt2012A_legMu");
49			tscaleMuPlusWCandPt2012A_legW=openFile(ana,"hltMuPlusWCandPt2012A_legW");
50	arizzi	1.16
51			tscaleMuID2012A=openFile(ana,"idMu2012A");
52			tscaleEleID2012A=openFile(ana,"idEle2012A");
53			tscaleEleID2012Awp80=openFile(ana,"idEle2012Awp80");
54	arizzi	1.12
55	arizzi	1.1	if(tscaleHLTmu == 0 \|\| tscaleIDmu == 0)
56			{
57			std::cout << "ERROR: cannot load Muon Trigger efficiencies" << std::endl;
58			}
59	arizzi	1.4
60
61	arizzi	1.1
62			}
63	arizzi	1.4
64	arizzi	1.10	static TTree * openFile(const edm::ParameterSet& ana, const char * name)
65	arizzi	1.4	{
66	arizzi	1.5	TFile *hltMuFile = new TFile (ana.getParameter<std::string> (name).c_str(),"read");
67	arizzi	1.4	if(hltMuFile) return (TTree*) hltMuFile->Get("tree");
68			else return 0;
69			}
70
71	arizzi	1.10	static std::pair<float,float> efficiencyFromPtEta(float pt1, float eta1, TTree *t)
72	arizzi	1.4	{
73	arizzi	1.15	// std::cout << "here " << t << " pt1 " << pt1 << " eta1 " << eta1 << std::endl;
74	arizzi	1.4	float s1 = 1.,err=1.;
75			std::pair<float,float> r(s1,err);
76			if(!t) return r;
77			float ptMin,ptMax,etaMin,etaMax,scale,error;
78			int count = 0;
79			t->SetBranchAddress("ptMin",&ptMin);
80			t->SetBranchAddress("ptMax",&ptMax);
81			t->SetBranchAddress("etaMin",&etaMin);
82			t->SetBranchAddress("etaMax",&etaMax);
83			t->SetBranchAddress("scale",&scale);
84			t->SetBranchAddress("error",&error);
85	arizzi	1.7	float lastPtBin = 200;
86			/* for(int jentry = 0; jentry < t->GetEntries(); jentry++)
87			{
88			t->GetEntry(jentry);
89			if(ptMax >= lastPtBin) lastPtBin =ptMax;
90			}*/
91	arizzi	1.4	for(int jentry = 0; jentry < t->GetEntries(); jentry++)
92			{
93			t->GetEntry(jentry);
94	arizzi	1.7	if(ptMax==lastPtBin) ptMax=1e99;
95	arizzi	1.4	if((pt1 > ptMin) && (pt1 < ptMax) && (eta1 > etaMin) && (eta1 < etaMax))
96			{
97			s1 = scale;
98			err=error;
99			count++;
100			}
101			}
102
103	sethzenz	1.18	if(count == 0)
104	arizzi	1.4	{
105			return r;
106			}
107
108			r.first=s1;
109			r.second = err;
110			return (r);
111			}
112	arizzi	1.1
113			float scaleMuIsoHLT(float pt1, float eta1)
114			{
115	arizzi	1.4	return efficiencyFromPtEta(pt1,eta1,tscaleHLTmu).first;
116	arizzi	1.1	}
117
118
119
120			float scaleMuID(float pt1, float eta1)
121			{
122	arizzi	1.4	return efficiencyFromPtEta(pt1,eta1,tscaleIDmu).first;
123	arizzi	1.1	}
124
125	wilken	1.22.4.1	float scaleMuCrossTrig(float pt1, float eta1)
126			{
127			return efficiencyFromPtEta(pt1,eta1,tscaleHLTmuonCrossTrig).first;
128			}
129
130			float scaleEleCrossTrig(float pt1, float eta1)
131			{
132			return efficiencyFromPtEta(pt1,eta1,tscaleHLTelectronCrossTrig).first;
133			}
134
135	arizzi	1.3
136	arizzi	1.9	double scaleDoubleEle17Ele8Aug( std::vector<float> pt, std::vector<float> eta )
137			{
138			std::vector< std::vector<float> > allEleWithEffs;
139	arizzi	1.21	for(unsigned int j=0; j< pt.size() && j < 10 ; j++)
140	arizzi	1.9	{
141			std::vector<float> thisEleEffs;
142			thisEleEffs.push_back(efficiencyFromPtEta(pt[j],eta[j],tscaleHLTele1Aug).first);
143			thisEleEffs.push_back(efficiencyFromPtEta(pt[j],eta[j],tscaleHLTele2Aug).first);
144
145			allEleWithEffs.push_back(thisEleEffs);
146			}
147
148			return combiner2Thr.weight<Trigger1High2Loose>(allEleWithEffs);
149
150			}
151
152
153	arizzi	1.5	double scaleDoubleEle17Ele8( std::vector<float> pt, std::vector<float> eta )
154			{
155			std::vector< std::vector<float> > allEleWithEffs;
156	arizzi	1.21	for(unsigned int j=0; j< pt.size() && j<10; j++)
157	arizzi	1.5	{
158			std::vector<float> thisEleEffs;
159			thisEleEffs.push_back(efficiencyFromPtEta(pt[j],eta[j],tscaleHLTele1).first);
160			thisEleEffs.push_back(efficiencyFromPtEta(pt[j],eta[j],tscaleHLTele2).first);
161	arizzi	1.9
162	arizzi	1.5	allEleWithEffs.push_back(thisEleEffs);
163			}
164
165			return combiner2Thr.weight<Trigger1High2Loose>(allEleWithEffs);
166
167			}
168	arizzi	1.16	double muId2012A( float pt, float eta){ return efficiencyFromPtEta(pt,eta,tscaleMuID2012A).first;}
169			double eleId2012A( float pt, float eta){ return efficiencyFromPtEta(pt,eta,tscaleEleID2012A).first;}
170			double eleId2012Awp80( float pt, float eta){ return efficiencyFromPtEta(pt,eta,tscaleEleID2012Awp80).first;}
171	arizzi	1.5
172			double scaleSingleEleMay( float pt, float eta){ return efficiencyFromPtEta(pt,eta,tscaleSingleEleMay).first;}
173			double scaleSingleEleV4( float pt, float eta){ return efficiencyFromPtEta(pt,eta,tscaleSingleEleV4).first; }
174			double scaleID80Ele( float pt, float eta) { return efficiencyFromPtEta(pt,eta,tscaleID80Ele).first; }
175			double scaleID95Ele( float pt, float eta) { return efficiencyFromPtEta(pt,eta,tscaleID95Ele).first; }
176			double scaleRecoEle( float pt, float eta){ return efficiencyFromPtEta(pt,eta,tscaleRecoEle).first; }
177			double scalePFMHTEle( float MetPFPt){
178			double weightPFMHTrigger=0.;
179
180			//FIXME: read from file
181	dlopes	1.8	if(MetPFPt>0. && MetPFPt<5.) weightPFMHTrigger=0.305;
182			if(MetPFPt>5. && MetPFPt<10.) weightPFMHTrigger=0.351;
183			if(MetPFPt>10. && MetPFPt<15.) weightPFMHTrigger=0.461;
184			if(MetPFPt>15. && MetPFPt<20.) weightPFMHTrigger=0.572;
185			if(MetPFPt>20. && MetPFPt<25.) weightPFMHTrigger=0.713;
186			if(MetPFPt>25. && MetPFPt<30.) weightPFMHTrigger=0.844;
187			if(MetPFPt>30. && MetPFPt<35.) weightPFMHTrigger=0.914;
188			if(MetPFPt>35. && MetPFPt<40.) weightPFMHTrigger=0.939;
189			if(MetPFPt>40. && MetPFPt<45.) weightPFMHTrigger=0.981;
190			if(MetPFPt>45. && MetPFPt<50.) weightPFMHTrigger=0.982;
191			if(MetPFPt>50. && MetPFPt<60.) weightPFMHTrigger=0.993;
192			if(MetPFPt>60. && MetPFPt<70.) weightPFMHTrigger=0.995;
193			if(MetPFPt>70. && MetPFPt<100.) weightPFMHTrigger=0.995;
194			if(MetPFPt>100.) weightPFMHTrigger=1.;
195	arizzi	1.5	return weightPFMHTrigger;
196			}
197
198	arizzi	1.11
199	arizzi	1.5	double scaleJet30Jet25( std::vector<float> pt, std::vector<float> eta)
200	arizzi	1.4	{
201	arizzi	1.5
202			std::vector< std::vector<float> > allJetsWithEffs;
203	arizzi	1.22	for(unsigned int j=0; j< pt.size() && j < 10 ; j++)
204	arizzi	1.5	{
205			std::vector<float> thisJetEffs;
206			thisJetEffs.push_back(efficiencyFromPtEta(pt[j],eta[j],tscaleHLTeleJet1).first);
207			thisJetEffs.push_back(efficiencyFromPtEta(pt[j],eta[j],tscaleHLTeleJet2).first);
208	arizzi	1.11	// std::cout << " jet pt " << pt[j] << " eta " << eta[j] << " eff1 " << thisJetEffs[0] << " eff2 " << thisJetEffs[1] << std::endl;
209	arizzi	1.5	allJetsWithEffs.push_back(thisJetEffs);
210	arizzi	1.11
211			}
212			float res = combiner2Thr.weight<Trigger1High2Loose>(allJetsWithEffs);
213			// std::cout << "Result is " << res << std::endl;
214			return res;
215			// return combiner2Thr.weight<Trigger1High2Loose>(allJetsWithEffs);
216
217			}
218			/*
219			TF1 fpt("f","1-exp(-0.157*(x-19.3))", 0., 9999999.);
220
221			MET80:
222			TF1 fmet80("f","1/ (1 + exp( -0.0709 * (x - 100.7)))", 0., 9999999.);
223
224			MET100:
225			TF1 fmet100("f","1/ (1 + exp( -0.0679 * (x - 128.8)))", 0., 9999999.);
226			*/
227
228			//LP curve used for MET
229			double scaleMetHLT( double met){
230			return 1. / (1. + ( exp( 0.059486 * ( 123.27 - met ))));
231			}
232
233			//MET80 component of the factorized JET+MET trigger
234			double scaleMET80(double et)
235			{
236			return 1. / (1. + exp( -0.0709 * (et - 100.7)));
237			}
238
239			//MET100 component
240			double scaleMET100(double et)
241			{
242			return 1. / (1. + exp( -0.0679 * (et - 128.8)));
243			}
244
245			//Single jet20 efficiency for MET+2CJet20
246			double jet20efficiency( double pt)
247			{
248			if(pt < 10 ) return 0;
249			return 1. - exp(-0.157*(pt-19.3));
250			}
251
252			//combined 2 jets efficiency out of N jets, using jet20 efficiency curve
253			double scale2CentralJet( std::vector<float> pt, std::vector<float> eta)
254			{
255
256			std::vector< std::vector<float> > allJetsWithEffs;
257	arizzi	1.22	for(unsigned int j=0; j< pt.size() && j < 10; j++)
258	arizzi	1.11	{
259			if(fabs(eta[j]) < 2.5)
260			{
261			std::vector<float> thisJetEffs;
262			thisJetEffs.push_back(jet20efficiency(pt[j]));
263			allJetsWithEffs.push_back(thisJetEffs);
264			}
265
266	arizzi	1.5	}
267
268	arizzi	1.11	return combiner1Thr.weight<Trigger2SingleThr>(allJetsWithEffs);
269			}
270
271			//New MET 150
272			double scaleMET150(double et)
273			{
274			return 1./ (1. + exp( -0.129226 * (et - 156.699)));
275	arizzi	1.4	}
276	arizzi	1.3
277	arizzi	1.6	float scaleMuOr30IsoHLT(float pt1, float eta1)
278			{
279			return efficiencyFromPtEta(pt1,eta1,tscaleHLTmuOr30).first;
280			}
281	arizzi	1.3
282
283	degrutto	1.19	//For 2012A HLT_DiCentralPFJet30_PFMHT80, valid for pfMET > 100 GeV:
284	arizzi	1.14	double scaleDiJet30MHT80_2012A(double x)
285			{
286	degrutto	1.19	if(x<100) return 0;
287	degrutto	1.20	return (1e0 - exp(-0.04197(x-75.73))) 0.9721 ;
288	arizzi	1.14	}
289	degrutto	1.19	//For 2012B HLT_DiCentralJetSumpT100_dPhi05_DiCentralPFJet60_25_PFMET100_HBHENoiseCleaned, valid for pfMET > 100 GeV:
290	arizzi	1.14	double scaleSumpT100MET100_2012B(double x)
291			{
292	degrutto	1.19	if(x<100) return 0;
293	degrutto	1.20	return (1e0 - exp(-0.06704(x-96.84))) 0.9199 ;
294	arizzi	1.14	}
295			//For 2012A+B HLT_PFMET150, valid for pfMET > 150 GeV:
296			double scalePFMET150_2012AB(double x)
297			{
298	degrutto	1.17	if(x<150) return 0;
299	degrutto	1.20	return (1e0 - exp(-0.07135(x-147.4))) 0.9707;
300	degrutto	1.19	}
301
302
303			//For 2012A HLT_PFMET150 OR HLT_DiCentralPFJet30_PFMHT80, valid for pfMET > 100 GeV:
304			double scalePFMET150orDiJetMET_2012A(double x)
305			{
306			if(x<100) return 0;
307	degrutto	1.20	return (1e0 - exp(-0.0412(x-75.52))) 0.9772;
308	degrutto	1.19	}
309
310			//For 2012B HLT_PFMET150 OR HLT_DiCentralJetSumpT100_dPhi05_DiCentralPFJet60_25_PFMET100_HBHENoiseCleaned, valid for pfMET > 100 GeV:
311			double scalePFMET150orDiJetMET_2012B(double x)
312			{
313			if(x<100) return 0;
314	degrutto	1.20	return (1e0 - exp(-0.05482(x-95.59))) 0.9702;
315	arizzi	1.14	}
316
317	degrutto	1.19
318			//For 2012C HLT_PFMET150 OR HLT_DiCentralJetSumpT100_dPhi05_DiCentralPFJet60_25_PFMET100_HBHENoiseCleaned, valid for pfMET > 100 GeV:
319			double scalePFMET150orDiJetMET_2012C(double x)
320			{
321			if(x<100) return 0;
322	degrutto	1.20	return (1e0 - exp(-0.05627(x-95.15))) 0.9659;
323	degrutto	1.19	}
324
325
326
327
328
329	arizzi	1.12	float doubleEle2012A( float pt1, float eta1, float pt2, float eta2)
330			{
331	arizzi	1.15	// std::cout << "di ele" << std::endl;
332	arizzi	1.12	float eff1_17 = efficiencyFromPtEta(pt1,eta1,tscaleDoubleEle2012A_leg17).first;
333			float eff2_17 = efficiencyFromPtEta(pt2,eta2,tscaleDoubleEle2012A_leg17).first;
334			float eff1_8 = efficiencyFromPtEta(pt1,eta1,tscaleDoubleEle2012A_leg8).first;
335			float eff2_8 = efficiencyFromPtEta(pt2,eta2,tscaleDoubleEle2012A_leg8).first;
336	arizzi	1.15	// std::cout << tscaleDoubleEle2012A_dZ << std::endl;
337	arizzi	1.13	float eff_dz = efficiencyFromPtEta(eta1,eta2,tscaleDoubleEle2012A_dZ).first; // despite the name pt,eta is actually eta1,eta2
338	arizzi	1.12
339	arizzi	1.13	return (eff1_17 * eff2_8 + eff2_17 * eff1_8 - eff1_17eff2_17)eff_dz;
340	arizzi	1.12
341			}
342			float doubleMuon2012A( float pt1, float eta1, float pt2, float eta2)
343			{
344	arizzi	1.15	// std::cout << "di mu" << std::endl;
345	arizzi	1.12	float eff1_17 = efficiencyFromPtEta(pt1,eta1,tscaleDoubleMuon2012A_leg17).first;
346			float eff2_17 = efficiencyFromPtEta(pt2,eta2,tscaleDoubleMuon2012A_leg17).first;
347			float eff1_8 = efficiencyFromPtEta(pt1,eta1,tscaleDoubleMuon2012A_leg8).first;
348			float eff2_8 = efficiencyFromPtEta(pt2,eta2,tscaleDoubleMuon2012A_leg8).first;
349	arizzi	1.15	// std::cout << tscaleDoubleMuon2012A_dZ << std::endl;
350	arizzi	1.13	float eff_dz = efficiencyFromPtEta(eta1,eta2,tscaleDoubleMuon2012A_dZ).first; // despite the name pt,eta is actually eta1,eta2
351	arizzi	1.12
352	arizzi	1.13	return (eff1_17 * eff2_8 + eff2_17 * eff1_8 - eff1_17eff2_17)eff_dz;
353	arizzi	1.12
354			}
355
356
357
358			float muPlusWCandPt2012A_legW( float pt, float eta){ return efficiencyFromPtEta(pt,eta,tscaleMuPlusWCandPt2012A_legW).first;}
359			float muPlusWCandPt2012A_legMu( float pt, float eta){ return efficiencyFromPtEta(pt,eta,tscaleMuPlusWCandPt2012A_legMu).first;}
360			float singleEle2012Awp80( float pt, float eta){ return efficiencyFromPtEta(pt,eta,tscaleSingleEle2012Awp80).first;}
361			float singleEle2012Awp95( float pt, float eta){ return efficiencyFromPtEta(pt,eta,tscaleSingleEle2012Awp95).first;}
362			float singleMuon2012A( float pt, float eta){ return efficiencyFromPtEta(pt,eta,tscaleSingleMuon2012A).first;}
363
364
365	arizzi	1.1	private:
366	arizzi	1.4	TTree * tscaleHLTele1;
367			TTree * tscaleHLTele2;
368			TTree * tscaleHLTeleJet1;
369			TTree * tscaleHLTeleJet2;
370	arizzi	1.5	TTree * tscaleID80Ele;
371			TTree * tscaleID95Ele;
372			TTree * tscaleRecoEle;
373	arizzi	1.6	TTree * tscaleHLTmuOr30;
374
375	arizzi	1.12	TTree * tscaleSingleEle2012Awp95;
376			TTree * tscaleSingleEle2012Awp80;
377			TTree * tscaleSingleMuon2012A;
378			TTree * tscaleDoubleEle2012A_leg8;
379			TTree * tscaleDoubleEle2012A_leg17;
380			TTree * tscaleDoubleMuon2012A_leg8;
381			TTree * tscaleDoubleMuon2012A_leg17;
382			TTree * tscaleMuPlusWCandPt2012A_legMu;
383			TTree * tscaleMuPlusWCandPt2012A_legW;
384	arizzi	1.13	TTree * tscaleDoubleEle2012A_dZ;
385			TTree * tscaleDoubleMuon2012A_dZ;
386	arizzi	1.5	// TTree * tscalePFMHTele;
387			TTree * tscaleSingleEleMay;
388			TTree * tscaleSingleEleV4;
389	arizzi	1.4
390	arizzi	1.9	TTree * tscaleHLTele1Aug;
391			TTree * tscaleHLTele2Aug;
392
393	arizzi	1.16	TTree * tscaleMuID2012A;
394			TTree * tscaleEleID2012A;
395			TTree * tscaleEleID2012Awp80;
396
397	arizzi	1.1	TTree * tscaleHLTmu;
398	wilken	1.22.4.1	TTree * tscaleIDmu;
399			TTree * tscaleHLTmuonCrossTrig;
400			TTree * tscaleHLTelectronCrossTrig;
401
402	arizzi	1.4	MultiThresholdEfficiency combiner2Thr;
403	arizzi	1.11	MultiThresholdEfficiency combiner1Thr;
404	arizzi	1.1	};
405
406			#endif