Merge pull request #120 from EmyrClement/hadronEfficiency

Hadron leg efficiency

interface/GlobalVariables.h

@@ -16,6 +16,7 @@
 #include "TH1D.h"
 #include "TH2F.h"
 #include "TH3F.h"
+#include "TEfficiency.h"
 #include <string>
 
 
@@ -51,6 +52,9 @@ struct Globals {
 
 	//jets
 	static JetAlgorithm::value jetAlgorithm;
+	static boost::shared_ptr<TEfficiency> hadronTriggerLegEfficiencyHistogram_nonIsoJets;
+	static boost::shared_ptr<TEfficiency> hadronTriggerLegEfficiencyHistogram_isoJets;
+	static boost::shared_ptr<TEfficiency> hadronTriggerLegEfficiencyHistogram_isoPFJets;
 	static boost::array<boost::shared_ptr<TF1>, 12> bL7Corrections;
 	static boost::array<boost::shared_ptr<TF1>, 12> lightL7Corrections;
 	static int JESsystematic;

interface/Python/ConfigFile.h

@@ -16,6 +16,7 @@
 #include "TH1D.h"
 #include "TH2F.h"
 #include "TH3F.h"
+#include "TEfficiency.h"
 #include <boost/shared_ptr.hpp>
 #include <boost/array.hpp>
 
@@ -33,6 +34,7 @@ class ConfigFile {
 	std::string MuonScaleFactorsFile() const;
     std::string bJetResoFile() const;
     std::string lightJetResoFile() const;
+    std::string hadronTriggerFile() const;
 	bool useHitFit() const;
 	bool fitterOutputFlag() const;
 	const std::vector<std::string>& inputFiles() const;
@@ -71,6 +73,8 @@ class ConfigFile {
 	std::string muonScaleFactorsFile_;
     std::string bJetResoFile_;
     std::string lightJetResoFile_;
+    bool getHadronTriggerFromFile_;
+    std::string hadronTriggerFile_;
 	bool useHitFit_;
 	bool fitterOutputFlag_;
 	std::vector<std::string> inputFiles_;
@@ -100,6 +104,7 @@ class ConfigFile {
 	boost::shared_ptr<TH2F> getMuonIdIsoScaleFactorsHistogram(std::string muonScaleFactorsFile);
 	boost::shared_ptr<TH3F> getMuonTriggerScaleFactorsHistogram(std::string muonScaleFactorsFile);
 	boost::array<boost::shared_ptr<TF1>, 12> getL7Correction(std::string correctionFile);
+	void getHadronTriggerLegHistogram(std::string hadronTriggerFile);
 };
 
 } /* namespace BAT */

interface/RecoObjects/Jet.h

@@ -9,6 +9,7 @@
 #define JET_H_
 #include "Particle.h"
 #include "TF1.h"
+#include "TEfficiency.h"
 #include <vector>
 #include <string>
 #include <boost/array.hpp>
@@ -135,6 +136,9 @@ class Jet: public Particle {
 
 	bool FirstDataLoosePFJetID() const;
 	bool FirstDataLooseCaloJetID() const;
+
+	double getEfficiencyCorrection( int scale_factor_systematic ) const;
+
 private:
 	JetAlgorithm::value usedAlgorithm;
 	double electromagneticFraction;
@@ -168,7 +172,6 @@ class Jet: public Particle {
 	ParticlePointer matchedGeneratedJet;
 	ParticlePointer unsmearedJet;
 	ParticlePointer smearedJet;
-
 };
 
 typedef boost::shared_ptr<Jet> JetPointer;

python/master_2011_53X_cfg.py

@@ -85,6 +85,7 @@
             'JetSmearingSystematic':0,
             'PUFile':'Data_PUDist_2011Full_central_68000mb_June2014.root',
             'MuonScaleFactors':'MuonEfficiencies_SF_2011_53X_DataMC.root',
+            'hadronTriggerFile':'hadronLegEfficiencies_electron.root',
             'BTagSystematic':0,
             'LightTagSystematic':0,
             'custom_file_suffix':'',
@@ -160,6 +161,8 @@ def getAnalysisSettings(analysisMode):
 PUFile = toolsFolder + "data/" + settings['PUFile']
 getMuonScaleFactorsFromFile = True
 MuonScaleFactorsFile = toolsFolder + "data/" + settings['MuonScaleFactors']
+getHadronTriggerFromFile = True
+hadronTriggerFile = toolsFolder + "data/" + settings['hadronTriggerFile']
 ElectronScaleFactorSystematic = settings['ElectronScaleFactorSystematic']
 MuonScaleFactorSystematic = settings['MuonScaleFactorSystematic']
 #JES Systematic, the +/- number of uncertainties to vary the jets with

python/master_2012_cfg.py

@@ -93,6 +93,7 @@
             'JetSmearingSystematic':0,
             'PUFile':'PileUp_2012_69300_truth_finebin.root',
             'MuonScaleFactors':'nofile.root',
+            'PUFile':'noFile.root',
             'BTagSystematic':0,
             'LightTagSystematic':0,
             'custom_file_suffix':'',
@@ -168,6 +169,8 @@ def getAnalysisSettings(analysisMode):
 PUFile = toolsFolder + "data/" + settings['PUFile']
 getMuonScaleFactorsFromFile = False
 MuonScaleFactorsFile = toolsFolder + "data/" + settings['MuonScaleFactors']
+getHadronTriggerFromFile = True
+hadronTriggerFile = toolsFolder + "data/" + settings['hadronTriggerFile']
 ElectronScaleFactorSystematic = settings['ElectronScaleFactorSystematic']
 MuonScaleFactorSystematic = settings['MuonScaleFactorSystematic']
 #JES Systematic, the +/- number of uncertainties to vary the jets with

src/Analysers/TTbar_plus_X_analyser.cpp

@@ -33,6 +33,7 @@ void TTbar_plus_X_analyser::ePlusJetsSignalAnalysis(const EventPtr event) {
 	if (topEplusJetsRefSelection_->passesFullSelectionExceptLastTwoSteps(event)) {
 		const JetCollection jets(topEplusJetsRefSelection_->cleanedJets(event));
 		const JetCollection bJets(topEplusJetsRefSelection_->cleanedBJets(event));
+		const JetPointer fourthJet = jets[3];
 		unsigned int numberOfBjets(bJets.size());
 		vector<double> bjetWeights;
 		if (event->isRealData()) {
@@ -51,7 +52,10 @@ void TTbar_plus_X_analyser::ePlusJetsSignalAnalysis(const EventPtr event) {
 
 		const LeptonPointer signalLepton = topEplusJetsRefSelection_->signalLepton(event);
 		const ElectronPointer signalElectron(boost::static_pointer_cast<Electron>(signalLepton));
-		double efficiencyCorrection = event->isRealData() ? 1. : signalElectron->getEfficiencyCorrection(false, Globals::ElectronScaleFactorSystematic, event->runnumber());
+
+		double hadronTriggerLegCorrection = event->isRealData() ? 1. : fourthJet->getEfficiencyCorrection( Globals::ElectronScaleFactorSystematic );
+
+		double efficiencyCorrection = event->isRealData() ? 1. : signalElectron->getEfficiencyCorrection(false, Globals::ElectronScaleFactorSystematic, event->runnumber()) * hadronTriggerLegCorrection;
 
 		for (unsigned int weightIndex = 0; weightIndex < bjetWeights.size(); ++weightIndex) {
 			double bjetWeight = bjetWeights.at(weightIndex);
@@ -146,6 +150,7 @@ void TTbar_plus_X_analyser::ePlusJetsQcdAnalysis(const EventPtr event) {
 			TTbarEPlusJetsReferenceSelection::AtLeastFourGoodJets)) {
 		const JetCollection jets(qcdNonIsoElectronSelection_->cleanedJets(event));
 		const JetCollection bJets(qcdNonIsoElectronSelection_->cleanedBJets(event));
+		const JetPointer fourthJet = jets[3];
 		unsigned int numberOfBjets(bJets.size());
 		vector<double> bjetWeights;
 		if (event->isRealData()) {
@@ -163,7 +168,8 @@ void TTbar_plus_X_analyser::ePlusJetsQcdAnalysis(const EventPtr event) {
 		unsigned int prescale(qcdNonIsoElectronSelection_->prescale(event));
 		const LeptonPointer signalLepton = qcdNonIsoElectronSelection_->signalLepton(event);
 		const ElectronPointer signalElectron(boost::static_pointer_cast<Electron>(signalLepton));
-		double efficiencyCorrection = event->isRealData() ? 1. : signalElectron->getEfficiencyCorrection(false, Globals::ElectronScaleFactorSystematic, event->runnumber()); // should really be  qcd=true, not false, and should have a loop in Electron.cpp to get scale factors for QCD selection, but it doesn't exist (yet)
+		double hadronTriggerLegCorrection = event->isRealData() ? 1. : fourthJet->getEfficiencyCorrection( Globals::ElectronScaleFactorSystematic );
+		double efficiencyCorrection = event->isRealData() ? 1. : signalElectron->getEfficiencyCorrection(false, Globals::ElectronScaleFactorSystematic, event->runnumber()) * hadronTriggerLegCorrection;
 
 		qcdNonIsoElectronAnalyser_->setPrescale(prescale);
 		metAnalyserqcdNonIsoElectronSelection_->setPrescale(prescale);
@@ -237,6 +243,7 @@ void TTbar_plus_X_analyser::ePlusJetsQcdAnalysis(const EventPtr event) {
 			TTbarEPlusJetsReferenceSelection::AtLeastFourGoodJets)) {
 		const JetCollection jets(qcdConversionSelection_->cleanedJets(event));
 		const JetCollection bJets(qcdConversionSelection_->cleanedBJets(event));
+		const JetPointer fourthJet = jets[3];
 		unsigned int numberOfBjets(bJets.size());
 		vector<double> bjetWeights;
 		if (event->isRealData()) {
@@ -253,7 +260,8 @@ void TTbar_plus_X_analyser::ePlusJetsQcdAnalysis(const EventPtr event) {
 		unsigned int prescale(qcdConversionSelection_->prescale(event));
 		const LeptonPointer signalLepton = qcdConversionSelection_->signalLepton(event);
 		const ElectronPointer signalElectron(boost::static_pointer_cast<Electron>(signalLepton));
-		double efficiencyCorrection = event->isRealData() ? 1. : signalElectron->getEfficiencyCorrection(false, Globals::ElectronScaleFactorSystematic, event->runnumber());
+		double hadronTriggerLegCorrection = event->isRealData() ? 1. : fourthJet->getEfficiencyCorrection( Globals::ElectronScaleFactorSystematic );
+		double efficiencyCorrection = event->isRealData() ? 1. : signalElectron->getEfficiencyCorrection(false, Globals::ElectronScaleFactorSystematic, event->runnumber()) * hadronTriggerLegCorrection;
 
 		qcdConversionsElectronAnalyser_->setPrescale(prescale);
 		metAnalyserqcdConversionSelection_->setPrescale(prescale);
@@ -343,7 +351,9 @@ void TTbar_plus_X_analyser::ePlusJetsQcdAnalysis(const EventPtr event) {
 		unsigned int prescale(qcdPFRelIsoEPlusJetsSelection_->prescale(event));
 		const LeptonPointer signalLepton = qcdPFRelIsoEPlusJetsSelection_->signalLepton(event);
 		const ElectronPointer signalElectron(boost::static_pointer_cast<Electron>(signalLepton));
-		double efficiencyCorrection = event->isRealData() ? 1. : signalElectron->getEfficiencyCorrection(false, Globals::ElectronScaleFactorSystematic, event->runnumber());
+		const JetPointer fourthJet = jets[3];
+		double hadronTriggerLegCorrection = event->isRealData() ? 1. : fourthJet->getEfficiencyCorrection( Globals::ElectronScaleFactorSystematic );
+		double efficiencyCorrection = event->isRealData() ? 1. : signalElectron->getEfficiencyCorrection(false, Globals::ElectronScaleFactorSystematic, event->runnumber()) * hadronTriggerLegCorrection;
 
 		qcdEPlusjetsPFRelIsoElectronAnalyser_->setPrescale(prescale);
 		for (unsigned int weightIndex = 0; weightIndex < bjetWeights.size(); ++weightIndex) {

src/GlobalsVariables.cpp

@@ -40,6 +40,9 @@ boost::shared_ptr<TH3F> Globals::muonTriggerScaleFactorsHistogram = boost::share
 
 //jets
 JetAlgorithm::value Globals::jetAlgorithm = JetAlgorithm::PF2PAT;
+boost::shared_ptr<TEfficiency> Globals::hadronTriggerLegEfficiencyHistogram_nonIsoJets = boost::shared_ptr<TEfficiency>(new TEfficiency());
+boost::shared_ptr<TEfficiency> Globals::hadronTriggerLegEfficiencyHistogram_isoJets = boost::shared_ptr<TEfficiency>(new TEfficiency());
+boost::shared_ptr<TEfficiency> Globals::hadronTriggerLegEfficiencyHistogram_isoPFJets = boost::shared_ptr<TEfficiency>(new TEfficiency());
 int Globals::JESsystematic = 0;
 int Globals::BJetSystematic = 0;
 int Globals::LightJetSystematic = 0;

src/Python/ConfigFile.cpp

@@ -33,6 +33,8 @@ ConfigFile::ConfigFile(int argc, char **argv) :
 		muonScaleFactorsFile_(PythonParser::getAttributeFromPyObject<string>(config, "MuonScaleFactorsFile")), //
 		bJetResoFile_(PythonParser::getAttributeFromPyObject<string>(config, "bJetResoFile")), //
 		lightJetResoFile_(PythonParser::getAttributeFromPyObject<string>(config, "lightJetResoFile")), //
+		getHadronTriggerFromFile_(PythonParser::getAttributeFromPyObject<bool>(config, "getHadronTriggerFromFile")), //
+		hadronTriggerFile_(PythonParser::getAttributeFromPyObject<string>(config, "hadronTriggerFile")), //
 		useHitFit_(PythonParser::getAttributeFromPyObject<bool>(config, "useHitFit")), //
 		fitterOutputFlag_(PythonParser::getAttributeFromPyObject<bool>(config, "produceFitterASCIIoutput")), //
 		inputFiles_(PythonParser::getVectorFromPythonObject(config, "inputFiles")), //
@@ -100,6 +102,7 @@ boost::program_options::variables_map ConfigFile::getParameters(int argc, char**
 	desc.add_options()("MuonScaleFactorsFile", value<std::string>(), "set input file for muon scale factors");
 	desc.add_options()("bJetResoFile", value<std::string>(), "set input root file for b-jet L7 resolutions");
 	desc.add_options()("lightJetResoFile", value<std::string>(), "set input root file for light jet L7 resolutions");
+	desc.add_options()("getHadronTriggerFromFile", value<bool>(), "state whether we are getting the electron trigger hadron leg efficiencies from a file or not");
 	desc.add_options()("fitter", value<bool>(), "turn on the fitter (HitFit)");
 	desc.add_options()("fitterASCIIoutput", value<bool>(), "produce full kinematic fit output in ASCII format");
 	desc.add_options()("TQAFPath", value<std::string>(),
@@ -205,6 +208,13 @@ string ConfigFile::MuonScaleFactorsFile() const {
 		return muonScaleFactorsFile_;
 }
 
+string ConfigFile::hadronTriggerFile() const {
+	if (programOptions.count("hadronTriggerFile"))
+		return programOptions["hadronTriggerFile"].as<std::string>();
+	else
+		return hadronTriggerFile_;
+}
+
 string ConfigFile::datasetInfoFile() const {
 	if (programOptions.count("datasetInfoFile"))
 		return programOptions["datasetInfoFile"].as<std::string>();
@@ -417,6 +427,11 @@ void ConfigFile::loadIntoMemory() {
 	std::cout << "ConfigFile.cpp: Globals::ElectronScaleFactorSystematic = " << Globals::ElectronScaleFactorSystematic << std::endl;
 	std::cout << "ConfigFile.cpp: Globals::MuonScaleFactorSystematic = " << Globals::MuonScaleFactorSystematic << std::endl;
 
+	if ( Globals::energyInTeV == 7 && getHadronTriggerFromFile_ ) {
+		std::cout << "Getting electron trigger hadron leg efficiencies from file " << hadronTriggerFile() << "." << std::endl;
+		getHadronTriggerLegHistogram(hadronTriggerFile());
+	}
+
 	//JES systematic
 	Globals::JESsystematic = jesSystematic();
 
@@ -491,6 +506,21 @@ boost::shared_ptr<TH3F> ConfigFile::getMuonTriggerScaleFactorsHistogram(std::str
 	return triggerHistogram;
 }
 
+void ConfigFile::getHadronTriggerLegHistogram(std::string hadronTriggerFile) {
+	using namespace std;
+
+	if (!boost::filesystem::exists(hadronTriggerFile)) {
+		cerr << "ConfigFile::getHadronTriggerLegHistogram(" << hadronTriggerFile << "): could not find file" << endl;
+		throw "Could not find hadron leg efficiency histogram file in " + hadronTriggerFile;
+	}
+
+	boost::scoped_ptr<TFile> file(TFile::Open(hadronTriggerFile.c_str()));
+	Globals::hadronTriggerLegEfficiencyHistogram_nonIsoJets = (boost::shared_ptr<TEfficiency>) ((TEfficiency*) file->Get("data_1")->Clone("data_1"));
+	Globals::hadronTriggerLegEfficiencyHistogram_isoJets = (boost::shared_ptr<TEfficiency>) ((TEfficiency*) file->Get("data_2")->Clone("data_2"));
+	Globals::hadronTriggerLegEfficiencyHistogram_isoPFJets = (boost::shared_ptr<TEfficiency>) ((TEfficiency*) file->Get("data_3")->Clone("data_3"));
+	file->Close();
+}
+
 unsigned int ConfigFile::nTupleVersion() const {
 	return nTupleVersion_;
 }

src/RecoObjects/Jet.cpp

@@ -525,5 +525,50 @@ double Jet::getBTagDiscriminator(BtagAlgorithm::value type) const {
 	return btag_discriminators.at((unsigned int) type);
 }
 
+double Jet::getEfficiencyCorrection( int scale_factor_systematic ) const {
+
+	if (Globals::energyInTeV == 7) {
+		double correction(1.);
+		double jetPt = pt();
+
+		if ( jetPt >= 100 ) return 1.0;
+
+		unsigned int binNumber = Globals::hadronTriggerLegEfficiencyHistogram_nonIsoJets->FindFixBin( jetPt );
+
+		double correction_nonIsoJets = 1.;
+		double lumi_nonIsoJets = Globals::luminosity * 0.076;
+		double correction_isoJets = 1.;
+		double lumi_isoJets = Globals::luminosity * 0.764;
+		double correction_isoPFJets = 1.;
+		double lumi_isoPFJets = Globals::luminosity * 0.16;
+
+
+		switch (scale_factor_systematic) {
+			case -1:
+				correction_nonIsoJets = Globals::hadronTriggerLegEfficiencyHistogram_nonIsoJets->GetEfficiency( binNumber ) - Globals::hadronTriggerLegEfficiencyHistogram_nonIsoJets->GetEfficiencyErrorLow( binNumber );
+				correction_isoJets = Globals::hadronTriggerLegEfficiencyHistogram_isoJets->GetEfficiency( binNumber ) - Globals::hadronTriggerLegEfficiencyHistogram_isoJets->GetEfficiencyErrorLow( binNumber );
+				correction_isoPFJets = Globals::hadronTriggerLegEfficiencyHistogram_isoPFJets->GetEfficiency( binNumber ) - Globals::hadronTriggerLegEfficiencyHistogram_isoPFJets->GetEfficiencyErrorLow( binNumber );
+				break;
+			case 1:
+				correction_nonIsoJets = Globals::hadronTriggerLegEfficiencyHistogram_nonIsoJets->GetEfficiency( binNumber ) + Globals::hadronTriggerLegEfficiencyHistogram_nonIsoJets->GetEfficiencyErrorUp( binNumber );
+				correction_isoJets = Globals::hadronTriggerLegEfficiencyHistogram_isoJets->GetEfficiency( binNumber ) + Globals::hadronTriggerLegEfficiencyHistogram_isoJets->GetEfficiencyErrorUp( binNumber );
+				correction_isoPFJets = Globals::hadronTriggerLegEfficiencyHistogram_isoPFJets->GetEfficiency( binNumber ) + Globals::hadronTriggerLegEfficiencyHistogram_isoPFJets->GetEfficiencyErrorUp( binNumber );
+				break;
+			default:
+				correction_nonIsoJets = Globals::hadronTriggerLegEfficiencyHistogram_nonIsoJets->GetEfficiency( binNumber );
+				correction_isoJets = Globals::hadronTriggerLegEfficiencyHistogram_isoJets->GetEfficiency( binNumber );
+				correction_isoPFJets = Globals::hadronTriggerLegEfficiencyHistogram_isoPFJets->GetEfficiency( binNumber );				
+		}
+
+		correction = ((correction_nonIsoJets * lumi_nonIsoJets) + (correction_isoJets * lumi_isoJets) + (correction_isoPFJets * lumi_isoPFJets)) / (Globals::luminosity);
+
+		return correction;
+	}
+	else {
+		return 1.;
+	}
+}
+
 }
 
+