diff --git a/PhysicsAnalysis/ElectronPhotonID/ElectronEfficiencyCorrection/CMakeLists.txt b/PhysicsAnalysis/ElectronPhotonID/ElectronEfficiencyCorrection/CMakeLists.txt
index bf066303bd7629ca4369ead0ca5780d2157716a4..9aa502462ad1997d59a2812085d9f86536225cb3 100644
--- a/PhysicsAnalysis/ElectronPhotonID/ElectronEfficiencyCorrection/CMakeLists.txt
+++ b/PhysicsAnalysis/ElectronPhotonID/ElectronEfficiencyCorrection/CMakeLists.txt
@@ -1,3 +1,4 @@
+# $Id: CMakeLists.txt 794295 2017-01-28 13:21:22Z kkoeneke $
################################################################################
# Package: ElectronEfficiencyCorrection
################################################################################
@@ -5,50 +6,91 @@
# Declare the package name:
atlas_subdir( ElectronEfficiencyCorrection )
+# Extra dependencies when not in AnalysisBase:
+set( extra_deps )
+if( NOT XAOD_STANDALONE )
+ set( extra_deps Control/AthenaBaseComps GaudiKernel )
+else()
+ set( extra_deps Control/xAODRootAccess )
+endif()
+
# Declare the package's dependencies:
-atlas_depends_on_subdirs( PUBLIC
- Control/AthContainers
- Control/AthToolSupport/AsgTools
- Event/xAOD/xAODEgamma
- PhysicsAnalysis/AnalysisCommon/PATCore
- PhysicsAnalysis/AnalysisCommon/PATInterfaces
- PRIVATE
- Tools/PathResolver
- Control/AthenaBaseComps
- Control/CxxUtils
- Event/xAOD/xAODCaloEvent
- Event/xAOD/xAODCore
- Event/xAOD/xAODEventInfo
- Event/xAOD/xAODTracking
- GaudiKernel )
+atlas_depends_on_subdirs(
+ PUBLIC
+ Control/AthContainers
+ Control/AthToolSupport/AsgTools
+ Event/xAOD/xAODEgamma
+ PhysicsAnalysis/AnalysisCommon/PATCore
+ PhysicsAnalysis/AnalysisCommon/PATInterfaces
+ PRIVATE
+ PhysicsAnalysis/Interfaces/AsgAnalysisInterfaces
+ Control/AthAnalysisBaseComps
+ Tools/PathResolver
+ Control/CxxUtils
+ Event/xAOD/xAODCaloEvent
+ Event/xAOD/xAODCore
+ Event/xAOD/xAODEventInfo
+ Event/xAOD/xAODTracking
+ Event/xAOD/xAODMetaData
-# External dependencies:
-find_package( Boost COMPONENTS filesystem thread system )
-find_package( ROOT COMPONENTS Core Tree MathCore Hist RIO pthread MathMore Minuit Minuit2 Matrix Physics HistPainter Rint PyROOT )
+ ${extra_deps} )
-# tag ROOTBasicLibs was not recognized in automatic conversion in cmt2cmake
+# External dependencies:
+find_package( Boost )
+find_package( ROOT COMPONENTS Core MathCore Hist RIO MathMore )
# Component(s) in the package:
atlas_add_library( ElectronEfficiencyCorrectionLib
- src/*.cxx
- Root/*.cxx
- PUBLIC_HEADERS ElectronEfficiencyCorrection
- INCLUDE_DIRS ${Boost_INCLUDE_DIRS} ${ROOT_INCLUDE_DIRS}
- LINK_LIBRARIES ${Boost_LIBRARIES} ${ROOT_LIBRARIES} AthContainers AsgTools xAODEgamma PATInterfaces PATCoreLib
- PRIVATE_LINK_LIBRARIES AthenaBaseComps xAODCaloEvent xAODCore xAODEventInfo xAODTracking GaudiKernel PathResolver )
-
-atlas_add_component( ElectronEfficiencyCorrection
- src/components/*.cxx
- INCLUDE_DIRS ${Boost_INCLUDE_DIRS} ${ROOT_INCLUDE_DIRS}
- LINK_LIBRARIES ${Boost_LIBRARIES} ${ROOT_LIBRARIES} AthContainers AsgTools xAODEgamma PATCoreLib PATInterfaces AthenaBaseComps xAODCaloEvent xAODCore xAODEventInfo xAODTracking GaudiKernel PathResolver ElectronEfficiencyCorrectionLib )
+ ElectronEfficiencyCorrection/*.h Root/*.cxx
+ PUBLIC_HEADERS ElectronEfficiencyCorrection
+ INCLUDE_DIRS ${Boost_INCLUDE_DIRS} ${ROOT_INCLUDE_DIRS}
+ LINK_LIBRARIES ${Boost_LIBRARIES} ${ROOT_LIBRARIES} AthContainers AsgTools
+ xAODEgamma PATInterfaces PATCoreLib
+ PRIVATE_LINK_LIBRARIES xAODCaloEvent xAODCore xAODEventInfo xAODTracking
+ xAODMetaData PathResolver )
+
+if( NOT XAOD_STANDALONE )
+ atlas_add_component( ElectronEfficiencyCorrection
+ src/*.h src/*.cxx src/components/*.cxx
+ LINK_LIBRARIES xAODEgamma PATInterfaces AthenaBaseComps xAODCore
+ xAODEventInfo GaudiKernel ElectronEfficiencyCorrectionLib )
+endif()
atlas_add_dictionary( ElectronEfficiencyCorrectionDict
- ElectronEfficiencyCorrection/ElectronEfficiencyCorrectionDict.h
- ElectronEfficiencyCorrection/selection.xml
- INCLUDE_DIRS ${Boost_INCLUDE_DIRS} ${ROOT_INCLUDE_DIRS}
- LINK_LIBRARIES ${Boost_LIBRARIES} ${ROOT_LIBRARIES} AthContainers AsgTools xAODEgamma PATCoreLib PATInterfaces AthenaBaseComps CxxUtils xAODCaloEvent xAODCore xAODEventInfo xAODTracking GaudiKernel PathResolver ElectronEfficiencyCorrectionLib )
+ ElectronEfficiencyCorrection/ElectronEfficiencyCorrectionDict.h
+ ElectronEfficiencyCorrection/selection.xml
+ LINK_LIBRARIES ElectronEfficiencyCorrectionLib )
+
+# Utilities provided by the package:
+atlas_add_executable( EgEfficiencyCorr_mem_check
+ util/EgEfficiencyCorr_mem_check.cxx
+ LINK_LIBRARIES AsgTools ElectronEfficiencyCorrectionLib )
+
+# AnalysisBase-only utilities:
+if( XAOD_STANDALONE )
+
+ atlas_add_executable( EgEfficiencyCorr_testCorrelationModels
+ util/testCorrelationModels.cxx
+ INCLUDE_DIRS ${ROOT_INCLUDE_DIRS}
+ LINK_LIBRARIES ${ROOT_LIBRARIES} xAODRootAccess xAODEventInfo xAODEgamma
+ xAODCore xAODCaloEvent xAODTracking AsgTools PATInterfaces
+ ElectronEfficiencyCorrectionLib )
+
+ atlas_add_executable( EgEfficiencyCorr_testEgEfficiencyCorr
+ util/testEgEfficiencyCorr.cxx
+ INCLUDE_DIRS ${ROOT_INCLUDE_DIRS}
+ LINK_LIBRARIES ${ROOT_LIBRARIES} xAODRootAccess xAODEventInfo xAODEgamma
+ xAODCore ElectronEfficiencyCorrectionLib )
+
+ atlas_add_executable( EgEfficiencyCorr_testEgEfficiencyCorrWithoutFile
+ util/testEgEfficiencyCorrWithoutFile.cxx
+ INCLUDE_DIRS ${ROOT_INCLUDE_DIRS}
+ LINK_LIBRARIES ${ROOT_LIBRARIES} xAODRootAccess xAODEventInfo xAODEgamma
+ xAODCore xAODCaloEvent xAODTracking AsgTools PATInterfaces
+ ElectronEfficiencyCorrectionLib )
+
+endif()
# Install files from the package:
atlas_install_joboptions( share/*.py )
atlas_install_data( data/*.root data/*.txt )
-
diff --git a/PhysicsAnalysis/ElectronPhotonID/ElectronEfficiencyCorrection/ElectronEfficiencyCorrection/AsgElectronEfficiencyCorrectionTool.h b/PhysicsAnalysis/ElectronPhotonID/ElectronEfficiencyCorrection/ElectronEfficiencyCorrection/AsgElectronEfficiencyCorrectionTool.h
index 271548dcf4a5e13bc3b412e153aed7cc322d6ba1..262346abf40341fb684d5b5442371f698a427560 100644
--- a/PhysicsAnalysis/ElectronPhotonID/ElectronEfficiencyCorrection/ElectronEfficiencyCorrection/AsgElectronEfficiencyCorrectionTool.h
+++ b/PhysicsAnalysis/ElectronPhotonID/ElectronEfficiencyCorrection/ElectronEfficiencyCorrection/AsgElectronEfficiencyCorrectionTool.h
@@ -27,6 +27,8 @@
#include "AthContainers/AuxElement.h"
//xAOD includes
#include "AsgTools/AsgTool.h"
+#include "AsgTools/AsgMetadataTool.h"
+
#include "PATInterfaces/ISystematicsTool.h"
#include "PATInterfaces/SystematicRegistry.h"
#include "PATInterfaces/CorrectionCode.h"
@@ -35,8 +37,7 @@
#include "xAODEgamma/ElectronFwd.h"
class AsgElectronEfficiencyCorrectionTool
- : virtual public IAsgElectronEfficiencyCorrectionTool,
- public asg::AsgTool
+ : virtual public IAsgElectronEfficiencyCorrectionTool, public asg::AsgMetadataTool
{
ASG_TOOL_CLASS(AsgElectronEfficiencyCorrectionTool, IAsgElectronEfficiencyCorrectionTool)
@@ -55,6 +56,11 @@ public:
/// Gaudi Service Interface method implementations
virtual StatusCode finalize();
+ // Introducing to check if METADATA working
+ virtual StatusCode beginInputFile();
+ virtual StatusCode beginEvent();
+ virtual StatusCode endInputFile();
+
// Main methods from IUserDataCalcTool
public:
@@ -96,6 +102,12 @@ public:
// Private member variables
private:
+ // To check if the metadat can be retrieved
+ bool m_metadata_retrieved = false;
+
+ // Get the simulation type from metadata
+ StatusCode get_simType_from_metadata(PATCore::ParticleDataType::DataType& result) const;
+
/// The main calculate method: the actual correction factors are determined here
const Root::TResult& calculate( const xAOD::Electron& egam, const unsigned int runnumber, int ¤tElectronSimplifiedUncorrSystRegion, int& currentElectronUncorrSystRegion ) const ;
diff --git a/PhysicsAnalysis/ElectronPhotonID/ElectronEfficiencyCorrection/ElectronEfficiencyCorrection/TElectronEfficiencyCorrectionTool.h b/PhysicsAnalysis/ElectronPhotonID/ElectronEfficiencyCorrection/ElectronEfficiencyCorrection/TElectronEfficiencyCorrectionTool.h
index 31cfcd935979be6326ff92d3527bac15fae13aaf..a67dce1cf35074ec5e06616700452c6159c2db19 100644
--- a/PhysicsAnalysis/ElectronPhotonID/ElectronEfficiencyCorrection/ElectronEfficiencyCorrection/TElectronEfficiencyCorrectionTool.h
+++ b/PhysicsAnalysis/ElectronPhotonID/ElectronEfficiencyCorrection/ElectronEfficiencyCorrection/TElectronEfficiencyCorrectionTool.h
@@ -111,76 +111,45 @@ namespace Root {
/// Load all histograms from the input file(s)
int getHistograms();
int getHistogramInDirectory( TKey *key );
- int setupHistogramsInFolder( TObjArray* dirNameArray, int lastIdx );
+ int setupHistogramsInFolder( const TObjArray& dirNameArray, int lastIdx );
void calcDetailLevels(TH1D *eig) ;
- std::vector<TObjArray*> *buildToyMCTable(TObjArray *sf, TObjArray *eig, TObjArray* stat, TObjArray* uncorr, std::vector<TObjArray*> *corr);
- std::vector<TH2D*> *buildSingleToyMC(TH2D *sf, TH2D* stat, TH2D* uncorr, TObjArray *corr);
- TH2D *buildSingleCombToyMC(TH2D *sf, TH2D* stat, TH2D* uncorr, TObjArray *corr);
+ std::vector<TObjArray> buildToyMCTable (const TObjArray &sf, const TObjArray &eig,
+ const TObjArray& stat, const TObjArray& uncorr, const std::vector<TObjArray> &corr);
+
+ std::vector<TH2D*> buildSingleToyMC(TH2D *sf, TH2D* stat, TH2D* uncorr, const TObjArray& corr);
+
+ TH2D* buildSingleCombToyMC(TH2D *sf, TH2D* stat, TH2D* uncorr, const TObjArray& corr);
/// Fill and interpret the setup, depending on which histograms are found in the input file(s)
- int setup( TObjArray* hist,
- std::vector< TObjArray* >& histList,
+ int setup( const TObjArray& hist,
+ std::vector< TObjArray >& histList,
std::vector< unsigned int >& beginRunNumberList,
std::vector< unsigned int >& endRunNumberList );
- int setupSys( std::vector<TObjArray*> *hist,
- std::vector< std::vector< TObjArray* > *> & histList);
+ int setupSys( std::vector<TObjArray> & hist,
+ std::vector< std::vector< TObjArray>> & histList);
void printDefaultReturnMessage(TString reason, int line);
- /// A debug flag: if true, print out more statements
- TRandom3 *m_Rndm;
+ TRandom3 m_Rndm;
int m_randomCounter;
-
bool m_isInitialized;
-
/// The detail level
int m_detailLevel;
-
int m_toyMCSF;
-
+ //
/// The seed
int m_seed;
bool m_doToyMC;
bool m_doCombToyMC;
int m_nToyMC;
- //These are pointers to vectors of vector pointers to TobjectArray pointers ...
- //Nighmare to delete
- std::vector< std::vector<TObjArray*> *> *m_uncorrToyMCSystFull, *m_uncorrToyMCSystFast;
-
int m_sLevel[3];
int m_nSys;
int m_nSysMax;
-
-
- /// The list of file name(s)
- std::vector< std::string > m_corrFileNameList;
-
- /// List of run numbers where histgrams become valid for full simulation
- std::vector< unsigned int > m_begRunNumberList;
-
- /// List of run numbers where histgrams stop being valid for full simulation
- std::vector< unsigned int > m_endRunNumberList;
-
-
- int m_runNumBegin, m_runNumEnd;
- /// List of run numbers where histgrams become valid for fast simulation
- std::vector< unsigned int > m_begRunNumberListFastSim;
-
- /// List of run numbers where histgrams stop being valid for fast simulation
- std::vector< unsigned int > m_endRunNumberListFastSim;
-
-
- /// List of histograms for full Geant4 simulation
- std::map<int, std::vector< TObjArray* > > m_histList;
- std::vector< std::vector< TObjArray* > > m_sysList;
-
- /// List of histograms for fast simulation
- std::map<int, std::vector< TObjArray* > > m_fastHistList;
- std::vector< std::vector< TObjArray* > > m_fastSysList;
-
+ int m_runNumBegin;
+ int m_runNumEnd;
/// The prefix string for the result
std::string m_resultPrefix;
@@ -208,15 +177,34 @@ namespace Root {
int m_nSimpleUncorrSyst;
+ /// The position of the efficiency scale factor uncorrelated systematic uncertainty in the result
+ int m_position_globalBinNumber;
+
+
+ //description here ?
+ std::vector< std::vector<TObjArray>> m_uncorrToyMCSystFull;
+ std::vector< std::vector<TObjArray>> m_uncorrToyMCSystFast;
+ /// The list of file name(s)
+ std::vector< std::string > m_corrFileNameList;
+ /// List of run numbers where histgrams become valid for full simulation
+ std::vector< unsigned int > m_begRunNumberList;
+ /// List of run numbers where histgrams stop being valid for full simulation
+ std::vector< unsigned int > m_endRunNumberList;
+ /// List of run numbers where histgrams become valid for fast simulation
+ std::vector< unsigned int > m_begRunNumberListFastSim;
+ /// List of run numbers where histgrams stop being valid for fast simulation
+ std::vector< unsigned int > m_endRunNumberListFastSim;
+ /// List of histograms for full Geant4 simulation
+ std::map<int, std::vector< TObjArray > > m_histList;
+ std::vector< std::vector< TObjArray > > m_sysList;
+ /// List of histograms for fast simulation
+ std::map<int, std::vector< TObjArray > > m_fastHistList;
+ std::vector< std::vector< TObjArray > > m_fastSysList;
/// The positions of the efficiency scale factor correlated sustematic uncertainties in the result
std::vector<int> m_position_corrSys;
-
/// The positions of the toy MC scale factors
std::vector<int> m_position_uncorrToyMCSF;
-
- /// The position of the efficiency scale factor uncorrelated systematic uncertainty in the result
- int m_position_globalBinNumber;
-
+ ///The vector holding the keys
std::vector<int> m_keys;
}; // End: class definition
diff --git a/PhysicsAnalysis/ElectronPhotonID/ElectronEfficiencyCorrection/Root/AsgElectronEfficiencyCorrectionTool.cxx b/PhysicsAnalysis/ElectronPhotonID/ElectronEfficiencyCorrection/Root/AsgElectronEfficiencyCorrectionTool.cxx
index 4d75e19a9d042d824219c808953b90409dd18da3..dbbea7bb499e08a0628299b0dfebc15f57a5de82 100644
--- a/PhysicsAnalysis/ElectronPhotonID/ElectronEfficiencyCorrection/Root/AsgElectronEfficiencyCorrectionTool.cxx
+++ b/PhysicsAnalysis/ElectronPhotonID/ElectronEfficiencyCorrection/Root/AsgElectronEfficiencyCorrectionTool.cxx
@@ -32,6 +32,11 @@
// xAOD includes
#include "xAODEgamma/Electron.h"
#include "xAODEventInfo/EventInfo.h"
+#ifndef ROOTCORE
+#include <boost/algorithm/string.hpp>
+#include "AthAnalysisBaseComps/AthAnalysisHelper.h"
+#endif
+#include "xAODMetaData/FileMetaData.h"
#include "PathResolver/PathResolver.h"
#include "ElectronEfficiencyCorrection/TElectronEfficiencyCorrectionTool.h"
@@ -41,11 +46,11 @@
namespace correlationModel{
enum model { COMBMCTOYS =0,
- MCTOYS=1,
- FULL=2,
- SIMPLIFIED=3,
- TOTAL=4,
- SYST=5
+ MCTOYS=1,
+ FULL=2,
+ SIMPLIFIED=3,
+ TOTAL=4,
+ SYST=5
};
}
@@ -53,7 +58,7 @@ namespace correlationModel{
// Standard constructor
// =============================================================================
AsgElectronEfficiencyCorrectionTool::AsgElectronEfficiencyCorrectionTool(std::string myname) :
- AsgTool(myname),
+ asg::AsgMetadataTool(myname),
m_rootTool(0),
m_affectedSys(),
m_appliedSystematics(0),
@@ -87,7 +92,7 @@ AsgElectronEfficiencyCorrectionTool::AsgElectronEfficiencyCorrectionTool(std::st
declareProperty("ResultName", m_resultName = "", "The string for the result");
declareProperty("CorrelationModel", m_correlation_model_name = "SIMPLIFIED",
- "Uncertainty correlation model. At the moment TOTAL, FULL, SIMPLIFIED, SYST, MCTOYS and COMBMCTOYS are implemented. SIMPLIFIED is the default option.");
+ "Uncertainty correlation model. At the moment TOTAL, FULL, SIMPLIFIED, SYST, MCTOYS and COMBMCTOYS are implemented. SIMPLIFIED is the default option.");
declareProperty("NumberOfToys", m_number_of_toys = 100,"Number of ToyMC replica, affecting MCTOYS and COMBMCTOYS correlation models only.");
declareProperty("MCToySeed", m_seed_toys = 0,"Seed for ToyMC replica, affecting MCTOYS and COMBMCTOYS correlation models only." );
declareProperty("MCToyScale", m_scale_toys = 1,"Scales Toy systematics up by this factor, affecting MCTOYS and COMBMCTOYS correlation models only." );
@@ -125,7 +130,7 @@ AsgElectronEfficiencyCorrectionTool::initialize() {
// The user wants to overwrite the m_dataType
// Otherwise it will be full
// When we will use metadata m_dataType will be whatever the metadata says i.e Full or Fast
- // We will need to check for the existence of overwrite
+ // We will need to check for the existence of overwrite
// and if overwrite is not set we will set the m_dataType based on metadata
if (m_dataTypeOverwrite != -1) {
if (m_dataTypeOverwrite != static_cast<int> (PATCore::ParticleDataType::Full)
@@ -136,7 +141,7 @@ AsgElectronEfficiencyCorrectionTool::initialize() {
m_dataType = static_cast<PATCore::ParticleDataType::DataType>(m_dataTypeOverwrite);
}
- //Find the relevant input files
+ //Find the relevant input files
//Fill the vector with filename using keys if the user
//has not passed the full filename as a property
if (m_corrFileNameList.size() == 0) {
@@ -148,7 +153,7 @@ AsgElectronEfficiencyCorrectionTool::initialize() {
// Resolve the paths to the input files for the full Geant4 simualtion corrections
for (size_t i = 0; i < m_corrFileNameList.size(); ++i) {
- std::string filename = PathResolverFindCalibFile(m_corrFileNameList.at(i));
+ std::string filename = PathResolverFindCalibFile(m_corrFileNameList.at(i));
if (filename.empty()) {
ATH_MSG_ERROR("Could NOT resolve file name " << m_corrFileNameList.at(i));
return StatusCode::FAILURE;
@@ -176,7 +181,7 @@ AsgElectronEfficiencyCorrectionTool::initialize() {
ATH_MSG_ERROR("Could NOT find systematics Substring file name " << m_sysSubstring);
return StatusCode::FAILURE;
}
- }
+ }
//
m_rootTool->setResultPrefix(m_resultPrefix);
m_rootTool->setResultName(m_resultName);
@@ -287,7 +292,7 @@ AsgElectronEfficiencyCorrectionTool::initialize() {
// =============================================================================
// Athena finalize method
// =============================================================================
-StatusCode
+StatusCode
AsgElectronEfficiencyCorrectionTool::finalize() {
if (!(m_rootTool->finalize())) {
ATH_MSG_ERROR("Something went wrong at finalize!");
@@ -316,7 +321,7 @@ AsgElectronEfficiencyCorrectionTool::getEfficiencyScaleFactor(const xAOD::Electr
if (!randomrunnumber.isAvailable(*eventInfo)) {
efficiencyScaleFactor = 1;
ATH_MSG_WARNING(
- "Pileup tool not run before using ElectronEfficiencyTool! SFs do not reflect PU distribution in data");
+ "Pileup tool not run before using ElectronEfficiencyTool! SFs do not reflect PU distribution in data");
return CP::CorrectionCode::Error;
}
runnumber = randomrunnumber(*(eventInfo));
@@ -324,7 +329,7 @@ AsgElectronEfficiencyCorrectionTool::getEfficiencyScaleFactor(const xAOD::Electr
//
//Get the result
const Root::TResult& result = calculate(inputObject, runnumber, currentSimplifiedUncorrSystRegion,
- currentUncorrSystRegion);
+ currentUncorrSystRegion);
efficiencyScaleFactor = result.getScaleFactor();
//
// The default of the underlying tool is -999 , if we are in a valid range
@@ -337,7 +342,7 @@ AsgElectronEfficiencyCorrectionTool::getEfficiencyScaleFactor(const xAOD::Electr
return CP::CorrectionCode::Ok;
}
// ==============================================================================//
- //Systematic Variations
+ //Systematic Variations
//We pass only one variation per time
// The applied systemetic is always one
// Either is relevant and acquires a values
@@ -374,15 +379,15 @@ AsgElectronEfficiencyCorrectionTool::getEfficiencyScaleFactor(const xAOD::Electr
if (m_nCorrSyst == 0) {
if (appliedSystematics().matchSystematic(CP::SystematicVariation("EL_EFF_" + m_sysSubstring + "CorrUncertainty",1))) {
- sys = sqrt(result.getTotalUncertainty() * result.getTotalUncertainty()
- - result.getResult(4) * result.getResult(4)); // total -stat
- func(efficiencyScaleFactor, sys);
+ sys = sqrt(result.getTotalUncertainty() * result.getTotalUncertainty()
+ - result.getResult(4) * result.getResult(4)); // total -stat
+ func(efficiencyScaleFactor, sys);
}
if (appliedSystematics().matchSystematic(CP::SystematicVariation("EL_EFF_" + m_sysSubstring + "CorrUncertainty" ,-1))) {
-
- sys = -1* sqrt(result.getTotalUncertainty() * result.getTotalUncertainty()
- - result.getResult(4) * result.getResult(4)); // total -stat
- func(efficiencyScaleFactor, sys);
+
+ sys = -1* sqrt(result.getTotalUncertainty() * result.getTotalUncertainty()
+ - result.getResult(4) * result.getResult(4)); // total -stat
+ func(efficiencyScaleFactor, sys);
}
}else if (m_correlation_model == correlationModel::TOTAL) { // one "TOTAL" uncertainty
if (appliedSystematics().matchSystematic(CP::SystematicVariation("EL_EFF_" + m_sysSubstring +
@@ -396,8 +401,8 @@ AsgElectronEfficiencyCorrectionTool::getEfficiencyScaleFactor(const xAOD::Electr
"1NPCOR_PLUS_UNCOR" ,-1))) {
sys = -1*result.getTotalUncertainty();
func(efficiencyScaleFactor, sys);
- }
- }
+ }
+ }
// =======================================================================
// Then the uncorrelated, we just need to see if the applied matches the current electron pt and eta
if (m_correlation_model == correlationModel::FULL) {// The Full Model
@@ -433,16 +438,16 @@ AsgElectronEfficiencyCorrectionTool::getEfficiencyScaleFactor(const xAOD::Electr
func(efficiencyScaleFactor, sys);
}
}
-
+
//If it has not returned so far , it means we wants to do the correlated for the full models
- for (int i = 0; i < m_nCorrSyst; ++i) {/// number of correlated sources
+ for (int i = 0; i < m_nCorrSyst; ++i) {/// number of correlated sources
if (appliedSystematics().matchSystematic(CP::SystematicVariation("EL_EFF_" + m_sysSubstring +
- Form("CorrUncertaintyNP%d", i),1))) {
+ Form("CorrUncertaintyNP%d", i),1))) {
sys = result.getResult(5 + i);
func(efficiencyScaleFactor, sys);
}
if (appliedSystematics().matchSystematic(CP::SystematicVariation("EL_EFF_" + m_sysSubstring +
- Form("CorrUncertaintyNP%d", i),-1))) {
+ Form("CorrUncertaintyNP%d", i),-1))) {
sys = -1* result.getResult(5 + i);
func(efficiencyScaleFactor, sys);
}
@@ -472,10 +477,10 @@ AsgElectronEfficiencyCorrectionTool::isAffectedBySystematic(const CP::Systematic
}
CP::SystematicSet sys = affectingSystematics();
-
- if (m_correlation_model == correlationModel::MCTOYS
+
+ if (m_correlation_model == correlationModel::MCTOYS
|| m_correlation_model == correlationModel::COMBMCTOYS ){
- return(sys.begin()->ensembleContains(systematic)) ;
+ return(sys.begin()->ensembleContains(systematic)) ;
}
else{
return (sys.find(systematic) != sys.end());
@@ -489,7 +494,7 @@ AsgElectronEfficiencyCorrectionTool::affectingSystematics() const {
}
// Register the systematics with the registry and add them to the recommended list
-CP::SystematicCode
+CP::SystematicCode
AsgElectronEfficiencyCorrectionTool::registerSystematics() {
CP::SystematicRegistry ®istry = CP::SystematicRegistry::getInstance();
@@ -554,8 +559,8 @@ AsgElectronEfficiencyCorrectionTool::applySystematicVariation(const CP::Systemat
// TRANSFER RESULT OF UNDERLYING TOOL TO xAOD TOOL
// =============================================================================
const Root::TResult& AsgElectronEfficiencyCorrectionTool::calculate(const xAOD::Electron &egam, const unsigned int runnumber,
- int ¤tSimplifiedUncorrSystRegion, int ¤tUncorrSystRegion
- ) const {
+ int ¤tSimplifiedUncorrSystRegion, int ¤tUncorrSystRegion
+ ) const {
double cluster_eta(-9999.9);
double et(0.0);
@@ -667,17 +672,117 @@ AsgElectronEfficiencyCorrectionTool::InitSystematics() {
return CP::SystematicCode::Ok;
}
+//===============================================================================
+// begin input file
+//===============================================================================
+StatusCode AsgElectronEfficiencyCorrectionTool::beginInputFile(){
+
+ // User preference of dataType, already done in initialize
+ if (m_dataTypeOverwrite != -1) return StatusCode::SUCCESS;
+
+ PATCore::ParticleDataType::DataType dataType_metadata;
+ const StatusCode status = get_simType_from_metadata(dataType_metadata);
+
+ if (status == StatusCode::SUCCESS) {
+ //m_metadata_retrieved isn't useful (might remove it later)
+ m_metadata_retrieved = true;
+ ATH_MSG_DEBUG("metadata from new file: " << (dataType_metadata == PATCore::ParticleDataType::Data ? "data" : (dataType_metadata == PATCore::ParticleDataType::Full ? "full simulation" : "fast simulation")));
+
+ if (dataType_metadata != PATCore::ParticleDataType::Data) {
+
+ if (m_dataTypeOverwrite == -1) { m_dataType = dataType_metadata; }
+ else {ATH_MSG_DEBUG("Use should set the dataType, otherwise it will take FullSim Type");}
+ }
+ }
+
+ else { // not able to retrieve metadata
+ m_metadata_retrieved = false;
+ ATH_MSG_DEBUG("not able to retrieve metadata, please set the dataType");
+ }
+
+return StatusCode::SUCCESS;
+}
+//===============================================================================
+// end input file
+//===============================================================================
+StatusCode AsgElectronEfficiencyCorrectionTool::endInputFile(){
+ m_metadata_retrieved = false;
+ return StatusCode::SUCCESS;
+
+}
+//===============================================================================
+// end input file
+//===============================================================================
+StatusCode AsgElectronEfficiencyCorrectionTool::beginEvent(){
+
+return StatusCode::SUCCESS;
+}
+
+//===============================================================================
+// Get Simulation flavor (FastSim or FullSim) from METADATA
+//===============================================================================
+StatusCode
+AsgElectronEfficiencyCorrectionTool::get_simType_from_metadata(PATCore::ParticleDataType::DataType& result) const
+{
+ // adapted from https://svnweb.cern.ch/trac/atlasoff/browser/PhysicsAnalysis/AnalysisCommon/CPAnalysisExamples/trunk/Root/MetadataToolExample.cxx
+#ifndef ROOTCORE
+ //Determine MC/Data
+ std::string dataType("");
+ if ( (AthAnalysisHelper::retrieveMetadata("/TagInfo", "project_name", dataType, inputMetaStore())).isFailure() ){
+ ATH_MSG_DEBUG("Failure to retrieve the project_name, Either running on data or something is wrong?");
+ }
+ if (!(dataType == "IS_SIMULATION")) {
+ result = PATCore::ParticleDataType::Data;
+ ATH_MSG_DEBUG("Running on simulation");
+ return StatusCode::SUCCESS;
+ }
+
+ // Determine Fast/FullSim
+ if (dataType == "IS_SIMULATION") {
+ std::string simType("");
+ ATH_CHECK(AthAnalysisHelper::retrieveMetadata("/Simulation/Parameters", "SimulationFlavour", simType, inputMetaStore()));
+ boost::to_upper(simType);
+ result = (simType.find("ATLFASTII")==std::string::npos) ? PATCore::ParticleDataType::Full : PATCore::ParticleDataType::Fast;
+ return StatusCode::SUCCESS;
+ }
+#endif
+
+ //here's how things will work dual use, when file metadata is available in files
+ if (inputMetaStore()->contains<xAOD::FileMetaData>("FileMetaData")) {
+ const xAOD::FileMetaData* fmd = 0;
+ ATH_CHECK(inputMetaStore()->retrieve(fmd, "FileMetaData"));
+
+ std::string simType("");
+ const bool s = fmd->value(xAOD::FileMetaData::simFlavour, simType);
+
+ if (!s) {
+ ATH_MSG_DEBUG("no sim flavour from metadata: must be data");
+ result = PATCore::ParticleDataType::Data;
+ return StatusCode::SUCCESS;
+ }
+ else {
+ ATH_MSG_DEBUG("sim type = " + simType);
+ result = simType == "FullSim" ? PATCore::ParticleDataType::Full : PATCore::ParticleDataType::Fast;
+ return StatusCode::SUCCESS;
+ }
+ }
+ else { // no metadata in the file
+ ATH_MSG_DEBUG("no metadata found in the file");
+ return StatusCode::FAILURE;
+ }
+}
+
//===============================================================================
// Map Key Feature
//===============================================================================
-// Gets the correction filename from map
+// Gets the correction filename from map
StatusCode
AsgElectronEfficiencyCorrectionTool::getFile(const std::string& recokey, const std::string& idkey, const std::string& isokey, const std::string& trigkey) {
std::string key = convertToOneKey(recokey, idkey, isokey, trigkey);
std::string mapFileName = PathResolverFindCalibFile(m_mapFile);
std::string value = getValueByKey(mapFileName, key);
-
+
if (value != "") {
m_corrFileNameList.push_back(value);
} else {
@@ -689,51 +794,51 @@ AsgElectronEfficiencyCorrectionTool::getFile(const std::string& recokey, const s
}
return StatusCode::FAILURE;
}
-
+
ATH_MSG_DEBUG("Full File Name is " + value);
return StatusCode::SUCCESS;
}
-// Convert reco, ID, iso and trigger key values into a
+// Convert reco, ID, iso and trigger key values into a
// single key according to the map file key format
-std::string
-AsgElectronEfficiencyCorrectionTool::convertToOneKey(const std::string& recokey, const std::string& idkey, const std::string& isokey, const std::string& trigkey) const {
+std::string
+AsgElectronEfficiencyCorrectionTool::convertToOneKey(const std::string& recokey, const std::string& idkey, const std::string& isokey, const std::string& trigkey) const {
+
+ std::string key;
- std::string key;
-
// Reconstruction Key
if (recokey != ""){ key = recokey; }
-
+
// Identification Key
- if (idkey != "" && (recokey == "" && isokey == "" && trigkey == "")){ key = idkey; }
-
+ if (idkey != "" && (recokey == "" && isokey == "" && trigkey == "")){ key = idkey; }
+
// Isolation Key
- if ((idkey != "" && isokey != "") && recokey == "" && trigkey == ""){ key = std::string(idkey + "_" + isokey); }
-
+ if ((idkey != "" && isokey != "") && recokey == "" && trigkey == ""){ key = std::string(idkey + "_" + isokey); }
+
// Trigger Key
if (trigkey != "" && idkey != "") {
-
- // Trigger SF file with isolation
+
+ // Trigger SF file with isolation
if (isokey != "") {
- key = std::string (trigkey + "_" + idkey + "_" + isokey);
+ key = std::string (trigkey + "_" + idkey + "_" + isokey);
} else {
- // Trigger SF file without isolation
- key = std::string(trigkey + "_" + idkey);
+ // Trigger SF file without isolation
+ key = std::string(trigkey + "_" + idkey);
}
}
ATH_MSG_DEBUG("Full Key is " + key);
return key;
}
-// Retrieves the value from the provided map file as
+// Retrieves the value from the provided map file as
// associated with the provided key
-std::string
+std::string
AsgElectronEfficiencyCorrectionTool::getValueByKey(const std::string& mapFile, const std::string& key) {
std::string value;
if (read(mapFile).isFailure()) {
ATH_MSG_ERROR("Couldn't read Map File" + mapFile);
- return "" ;
+ return "" ;
}
if (getValue(key, value) == "") {
ATH_MSG_DEBUG("Error(" + key + ") not found ");
@@ -758,7 +863,7 @@ AsgElectronEfficiencyCorrectionTool::read(const std::string& strFile) {
getline(is,strLine);
int nPos = strLine.find('=');
-
+
if ((signed int)std::string::npos == nPos) continue; // no '=', invalid line;
std::string strKey = strLine.substr(0,nPos);
std::string strVal = strLine.substr(nPos + 1, strLine.length() - nPos + 1);
@@ -766,11 +871,11 @@ AsgElectronEfficiencyCorrectionTool::read(const std::string& strFile) {
}
return StatusCode::SUCCESS;
}
-// Retrieves the value from the map file if
-// the provided key is found. If the key has an
+// Retrieves the value from the map file if
+// the provided key is found. If the key has an
// association then, the actual retrieved value would
-// be assigned to the 2nd argument of this method
-std::string
+// be assigned to the 2nd argument of this method
+std::string
AsgElectronEfficiencyCorrectionTool::getValue(const std::string& strKey, std::string& strValue) {
std::map<std::string,std::string>::const_iterator i;
@@ -778,8 +883,7 @@ AsgElectronEfficiencyCorrectionTool::getValue(const std::string& strKey, std::st
if (i != m_map.end()) {
strValue = i->second;
- return strValue;
+ return strValue;
}
return "";
}
-
diff --git a/PhysicsAnalysis/ElectronPhotonID/ElectronEfficiencyCorrection/Root/ElectronChargeEfficiencyCorrectionTool.cxx b/PhysicsAnalysis/ElectronPhotonID/ElectronEfficiencyCorrection/Root/ElectronChargeEfficiencyCorrectionTool.cxx
new file mode 100644
index 0000000000000000000000000000000000000000..1fcbfc8277dc4668d9aec3074032c0d145687b2d
--- /dev/null
+++ b/PhysicsAnalysis/ElectronPhotonID/ElectronEfficiencyCorrection/Root/ElectronChargeEfficiencyCorrectionTool.cxx
@@ -0,0 +1,725 @@
+/*
+ Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration
+*/
+
+/**
+ @class ElectronChargeEfficiencyCorrectionTool
+ @brief Apply the correction for the electrons charge mis-ID different rates in MC/data
+
+ @author Giulia Gonella <giulia.gonella@cern.ch>
+ @date September 2015
+*/
+
+//#define XXXMSG(MSG) std::cout << __FILE__ << ":" << __LINE__ << ": " << MSG << std::endl;
+
+// Include this class's header
+// #include "ElectronChargeEfficiencyCorrectionTool/ElectronChargeEfficiencyCorrectionTool.h"
+#include "ElectronChargeEfficiencyCorrectionTool.h"
+
+// xAOD includes
+#include "PathResolver/PathResolver.h"
+// #include "xAODTruth/TruthEventContainer.h"
+// #include "xAODEgamma/EgammaTruthxAODHelpers.h"
+// #include "MCTruthClassifier/IMCTruthClassifier.h"
+#include "xAODEventInfo/EventInfo.h"
+#include "xAODEgamma/Electron.h"
+
+// ROOT includes
+#include "TFile.h"
+// #include "TRandom.h"
+
+// STL includes
+#include <stdlib.h> /* atoi */
+
+
+// =============================================================================
+// Standard constructor
+// =============================================================================
+CP::ElectronChargeEfficiencyCorrectionTool::ElectronChargeEfficiencyCorrectionTool(std::string name) :
+ AsgTool(name),
+ m_dataTypeOverwrite(-1),
+ m_eventInfoCollectionName("EventInfo"),
+ m_SF_SS(),
+ m_SF_OS(),
+ m_RunNumbers(),
+ m_useRandomRunNumber(true),
+ m_defaultRandomRunNumber(999999),
+ m_filename(""),
+ m_workingPoint(""),
+ m_eta_lowlimit(0.0),
+ m_eta_uplimit(0.0),
+ m_pt_lowlimit(0.0),
+ m_pt_uplimit(0.0),
+ m_gevmev(0.0),
+ m_truthCharge(0),
+ m_filtered_sys_sets(),
+ m_mySysConf(),
+ m_affectingSys(),
+ m_sf_decoration_name("chargeIDEffiSF"),
+ m_sfDec(0)
+{
+ // Declare the needed properties
+ declareProperty("CorrectionFileName", m_filename, "Name of the file with charge flipping rates");
+ declareProperty("WorkingPoint", m_workingPoint, "Name of working point folder in the file");
+ declareProperty("ScaleFactorDecorationName", m_sf_decoration_name);
+ declareProperty("ForceDataType", m_dataTypeOverwrite,
+ "Force the DataType of the electron to specified value (to circumvent problem of incorrect DataType for forward electrons in some old releases)");
+ declareProperty("EventInfoCollectionName", m_eventInfoCollectionName, "The EventInfo Collection Name");
+ declareProperty("UseRandomRunNumber", m_useRandomRunNumber);
+ declareProperty("DefaultRandomRunNumber", m_defaultRandomRunNumber);
+
+}
+
+// =============================================================================
+// Standard destructor
+// =============================================================================
+CP::ElectronChargeEfficiencyCorrectionTool::~ElectronChargeEfficiencyCorrectionTool()
+{
+ if(m_sfDec) delete m_sfDec;
+}
+
+
+// =============================================================================
+// Athena initialize method
+// =============================================================================
+StatusCode CP::ElectronChargeEfficiencyCorrectionTool::initialize()
+{
+ ATH_MSG_DEBUG("initializing");
+
+ // initialize the random number generator (used in case of charge flip approach)
+ //m_Rndm = new TRandom3(1);
+
+ if(m_sfDec) delete m_sfDec;
+ m_sfDec = new SG::AuxElement::Decorator< float>(m_sf_decoration_name);//xxxx
+
+ // Do some consistency checks //xxx
+ bool allOK(true);
+
+ // Stop here if the user configuration is wrong
+ if ( !allOK ) return StatusCode::FAILURE;
+
+ //Resolve the path to the input file for the charge flip rates
+ const std::string rootfilename = PathResolverFindCalibFile(m_filename);
+ if (m_filename.empty()) {
+ ATH_MSG_ERROR ( " PathResolver was not able to find the file ... aborting" );
+ return StatusCode::FAILURE;
+ }
+
+ // Getting the root file and histograms
+ TFile* rootFile = TFile::Open( rootfilename.c_str() );
+
+ // protection against bad file
+ if ( rootFile==0 ) {
+ ATH_MSG_ERROR ( " Was not able to open file: " << rootfilename << " ...... aborting" );
+ return StatusCode::FAILURE;
+ }
+
+ /////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
+ //
+ // explanation: attempt to loop generally over a file
+ // -- if certain SINGALWORD is present -- then this is taken as a signal,
+ // that this is another dimension... can be dynamically added.
+ // e.g.
+ // SFSyst<number>_RunNumber<minRN>-<maxRN>_Nvtx<minNvtx>-<maxNvtx>
+ // SFStat_RunNumber<minRN>-<maxRN>_Nvtx<minNvtx>-<maxNvtx>
+ // SFCentral_RunNumber<minRN>-<maxRN>_Nvtx<minNvtx>-<maxNvtx>
+
+ // Then can create a key that will dynamically give us access to a map:
+ // std::map<std::string key, std::vector<TH2D *>> m_SF_SS; // keys (e.g. RunNumber223333_319200_Nvtx0_10_Phi1.5_1.6) mapping to vector of SF histograms --> vector m_SF: 0=nominal, 1=stat, 2,3,4...n=syst
+ // std::map<std::string key, std::vector<TH2D *>> m_SF_OS; // keys (e.g. RunNumber223333_319200_Nvtx0_10_Phi1.5_1.6) mapping to vector of SF histograms --> vector m_SF: 0=nominal, 1=stat, 2,3,4...n=syst
+ // TFile* data/ChMisIDSF_TightLL_FixedCutTight.root
+ // KEY: TH2F SFCentral_RunNumber296939_311481_SS;1 SFCentral_RunNumber296939_311481_SS
+ // KEY: TH2F SFCentral_RunNumber296939_311481_OS;1 SFCentral_RunNumber296939_311481_OS
+ // KEY: TH2F STAT_RunNumber296939_311481_SS;1 STAT_RunNumber296939_311481_SS
+ // KEY: TH2F STAT_RunNumber296939_311481_OS;1 STAT_RunNumber296939_311481_OS
+ // KEY: TH2F SYST_RunNumber296939_311481_total_SS;1 SYST_RunNumber296939_311481_SS: total
+ // KEY: TH2F SYST_RunNumber296939_311481_total_OS;1 SYST_RunNumber296939_311481_OS: total
+
+
+ m_SF_SS.clear();
+ m_SF_OS.clear();
+ TList* keyListfolder = rootFile->GetListOfKeys();
+
+ for ( int j=0; j<keyListfolder->GetEntries(); j++ ){
+
+ std::string name = ( keyListfolder->At(j)->GetName() );
+ ATH_MSG_DEBUG("Got ROOT object with name: " << name);
+ if ( name.find(Form("SFCentral_") ) != std::string::npos){
+ ATH_MSG_VERBOSE("Found name 'SFCentral_' in ROOT object name");
+ // Check for opposite-sign (=opposite-charge)
+ bool isOS = false;
+ if ( name.find(Form("_OS") ) != std::string::npos ){
+ isOS = true;
+ ATH_MSG_VERBOSE("Found name '_OS' in ROOT object name");
+ }
+ if (isOS){
+ std::string histid = ( keyListfolder->At(j)->GetName() );
+ histid.erase(0,10);
+ histid.erase(histid.size()-3,3);// remove _SS, _OS
+ ATH_MSG_VERBOSE("Using histid: " << histid);
+
+ if ( histid.find("RunNumber") != std::string::npos ){
+ ATH_MSG_VERBOSE("Found name 'RunNumber' in histid");
+ std::string runlow = histid;
+ runlow.erase(histid.find(Form("RunNumber") ),9);
+ runlow.erase(runlow.find("_"),runlow.size() );
+ m_RunNumbers.push_back( static_cast<unsigned int>(atoi(runlow.c_str())) );
+ std::string runhigh = histid;
+ runhigh.erase(histid.find(Form("RunNumber") ),9);
+ runhigh.erase(0,runhigh.find("_")+1);
+ m_RunNumbers.push_back( static_cast<unsigned int>(atoi(runhigh.c_str())) );
+ }
+ ATH_MSG_VERBOSE("Using histid (OS hid): " << histid);
+ m_SF_OS[histid].push_back( (TH2D*)rootFile->Get( keyListfolder->At(j)->GetName() ));
+ }
+ else {
+ std::string histid = ( keyListfolder->At(j)->GetName() );
+ histid.erase(0,10);
+ histid.erase(histid.size()-3,3);// remove _SS, _OS
+ ATH_MSG_VERBOSE("Using histid (do we this in ? SS): " << histid);
+ m_SF_SS[histid].push_back( (TH2D*)rootFile->Get( keyListfolder->At(j)->GetName() ));
+ }
+ }///// if ( name.find(Form("SFCentral_") ) != std::string::npos)
+
+
+ /// STAT ERROR
+ if ( name.find(Form("STAT_") ) != std::string::npos ){
+ ATH_MSG_VERBOSE("Found name 'STAT_' in ROOT object name");
+ bool isOS = false;
+ if ( name.find(Form("_OS") ) != std::string::npos ){
+ isOS = true;
+ ATH_MSG_VERBOSE("Found name '_OS' in ROOT object name");
+ }
+ if ( isOS ){
+ std::string histid = ( keyListfolder->At(j)->GetName() );
+ histid.erase(0,5);
+ histid.erase(histid.size()-3,3);// remove _SS, _OS
+ ATH_MSG_VERBOSE("Using histid: " << histid);
+
+ if ( histid.find("RunNumber") != std::string::npos ){
+ ATH_MSG_VERBOSE("Found name 'RunNumber' in histid");
+ std::string runlow = histid;
+ runlow.erase(histid.find(Form("RunNumber") ),9);
+ runlow.erase(runlow.find("_"),runlow.size() );
+ m_RunNumbers.push_back( static_cast<unsigned int>(atoi(runlow.c_str())) );
+ std::string runhigh = histid;
+ runhigh.erase(histid.find(Form("RunNumber") ),9);
+ runhigh.erase(0,runhigh.find("_")+1);
+ m_RunNumbers.push_back( static_cast<unsigned int>(atoi(runhigh.c_str())) );
+ }
+ ATH_MSG_VERBOSE("Using histid (OS hid): " << histid);
+ m_SF_OS[histid].push_back( (TH2D*)rootFile->Get( keyListfolder->At(j)->GetName() ));
+ }
+ else {
+ std::string histid = ( keyListfolder->At(j)->GetName() );
+ ATH_MSG_VERBOSE("Found histid: " << histid);
+ histid.erase(0,5);
+ histid.erase(histid.size()-3,3);// remove _SS, _OS
+ ATH_MSG_VERBOSE("Using histid (do we this in ? SS): " << histid);
+ m_SF_SS[histid].push_back( (TH2D*)rootFile->Get( keyListfolder->At(j)->GetName() ));
+ }
+
+ }///// if ( name.find(Form("SYST") ) != std::string::npos)
+
+
+ /// STAT ERROR
+ if ( name.find(Form("SYST") ) != std::string::npos ){
+ ATH_MSG_VERBOSE("Found name 'SYST' in ROOT object name");
+ bool isOS = false;
+ if ( name.find(Form("_OS") ) != std::string::npos ){
+ isOS = true;
+ ATH_MSG_VERBOSE("Found name '_OS' in ROOT object name");
+ }
+ if ( isOS ){
+ std::string histid = ( keyListfolder->At(j)->GetName() );
+ histid.erase(0,4);
+ histid.erase(histid.size()-3,3);// remove _SS, _OS
+
+ std::string sysname = histid;
+ sysname.erase(sysname.find("_"),sysname.size());
+ m_systematics.push_back(sysname);
+
+ histid.erase(0,histid.find("_")+1);// remove _SS, _OS
+ ATH_MSG_VERBOSE("Using syst histid: " << histid);
+
+ if (histid.find("RunNumber") != std::string::npos ) {
+ std::string runlow = histid;
+ runlow.erase(histid.find(Form("RunNumber") ),9);
+ runlow.erase(runlow.find("_"),runlow.size() );
+ m_RunNumbers.push_back( static_cast<unsigned int>(atoi(runlow.c_str())) );
+ std::string runhigh = histid;
+ runhigh.erase(histid.find(Form("RunNumber") ),9);
+ runhigh.erase(0,runhigh.find("_")+1);
+ m_RunNumbers.push_back( static_cast<unsigned int>(atoi(runhigh.c_str())) );
+ }
+ ATH_MSG_VERBOSE("Using histid (OS hid): " << histid);
+ m_SF_OS[histid].push_back( (TH2D*)rootFile->Get( keyListfolder->At(j)->GetName() ));
+ }
+ else {
+ std::string histid = ( keyListfolder->At(j)->GetName() );
+ histid.erase(0,4);
+ histid.erase(histid.size()-3,3);// remove _SS, _OS
+ histid.erase(0,histid.find("_")+1);// remove _SS, _OS
+ ATH_MSG_VERBOSE("Using histid (sys ? SS): " << histid);
+ m_SF_SS[histid].push_back( (TH2D*)rootFile->Get( keyListfolder->At(j)->GetName() ));
+ }
+
+ }///end // if ( name.find(Form("SYST") ) != std::string::npos)
+
+ }
+
+ /////////// checks ... --> same vector length... all files there?
+
+ if ( m_SF_OS.size() <1 || m_SF_SS.size() <1 || m_SF_SS.size()!=m_SF_OS.size() ) {
+ ATH_MSG_ERROR("OS/SS SF vectors not filled or of different size. -- Problem with files. -- Report to <hn-atlas-EGammaWG@cern.ch>");
+ return CP::CorrectionCode::Error;
+}
+
+
+ /////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
+
+ // rootFile->Close();
+
+ // Determine the limits of validity
+
+ /// here: need to use iterator over map!!!
+ ATH_MSG_DEBUG("Having m_SF_OS.size() = " << m_SF_OS.size());
+ std::map<std::string,std::vector<TH2D*> >::iterator it = m_SF_OS.begin();
+
+ // Get the kinematic limits
+ m_eta_lowlimit = (*it).second.at(0)->GetYaxis()->GetXmin();
+ m_eta_uplimit = (*it).second.at(0)->GetYaxis()->GetXmax();
+ ATH_MSG_VERBOSE("|eta| limits " << m_eta_lowlimit << ", " << m_eta_uplimit);
+
+ m_pt_lowlimit = (*it).second.at(0)->GetXaxis()->GetXmin();
+ m_pt_uplimit = (*it).second.at(0)->GetXaxis()->GetXmax();
+ ATH_MSG_VERBOSE("pt limits " << m_pt_lowlimit << ", " << m_pt_uplimit);
+
+ // Check if the input file is in GeV or MeV
+ if ( m_pt_uplimit > 1500 ) {
+ ATH_MSG_VERBOSE("Rates in input file are in MeV");
+ m_gevmev = 1.;
+ }
+ else {
+ ATH_MSG_VERBOSE("Rates in input file are in GeV");
+ m_gevmev = 0.001;
+ }
+
+ // Systematics // dynamic too?
+ m_affectingSys = affectingSystematics();
+
+
+ return StatusCode::SUCCESS;
+}
+
+
+// =============================================================================
+// Athena finalize method
+// =============================================================================
+StatusCode CP::ElectronChargeEfficiencyCorrectionTool::finalize()
+{
+ return StatusCode::SUCCESS ;
+}
+
+
+
+//---------------------------------------------------------------------------------------
+// Get the scale factor for the electron
+//---------------------------------------------------------------------------------------
+
+//
+
+CP::CorrectionCode
+CP::ElectronChargeEfficiencyCorrectionTool::getEfficiencyScaleFactor(const xAOD::IParticle& part, double& sf) const {
+ ATH_MSG_DEBUG("In CP::ElectronChargeEfficiencyCorrectionTool::getEfficiencyScaleFactor(const xAOD::IParticle& part, double& sf) const");
+ if ( part.type() != xAOD::Type::Electron ){
+ ATH_MSG_ERROR("This function requires an electron to be passed. Failing!");
+ return CP::CorrectionCode::Error;
+ }
+ const xAOD::Electron& ele = static_cast<const xAOD::Electron&>(part);
+
+ // initialize the SF at 1
+ sf = 1.0;
+
+ // checking on the truth electron: up to now if this is not a good ele it's returning
+ bool goodEle = false;
+ CP::CorrectionCode goodEle_result = this->isGoodEle( ele, goodEle); //// ## Giulia: This is to change
+ if ( goodEle_result != CP::CorrectionCode::Ok ) {
+ sf = -999.0;
+ ATH_MSG_DEBUG("This is the check of goodeleCC in getscalefactor. Scale factor set to -999");
+ return goodEle_result;
+ }
+
+ if ( goodEle == false ) {
+ // electron is background electron and should not be corrected
+ return CP::CorrectionCode::Ok;
+ ATH_MSG_DEBUG("Here goodele is false but CC ok");
+ }
+
+ // taking reconstructed variables
+ int reco_ele_charge = ele.charge();
+ const double ele_pt = ele.pt()*m_gevmev;
+ const double ele_eta = std::abs(ele.eta());
+
+ // getting the truth charge
+ int truth_ele_charge = 9999;
+ CP::CorrectionCode charge_result = this->getEleTruthCharge( ele, truth_ele_charge); //// ## Giulia: second function to change
+ if ( charge_result != CP::CorrectionCode::Ok ) {
+ sf = -9999.0;
+ ATH_MSG_VERBOSE("This is check of geteletruthchargeCC in getscalefactor. Scale factor set to -9999");
+ return charge_result;
+ }
+
+ if ( truth_ele_charge == 0 ) {
+ ATH_MSG_DEBUG("Here truth charge is =0!!");
+ return CP::CorrectionCode::Ok;
+ }
+
+ ATH_MSG_DEBUG("Reco charge = " << reco_ele_charge << "; Truth charge = " << truth_ele_charge);
+
+ // getting the rates from file....
+ float retVal(0.0);
+
+ //////////////////////////////////////
+ // here determine, WHICH of the [histid] to choose (after cuuts on runnumber etc....)
+ std::string cutRunNumber = "all";
+
+ if (m_RunNumbers.size()>0) {
+ unsigned int runnumber = m_defaultRandomRunNumber;
+ ATH_MSG_DEBUG("RandomRunNumber: " << runnumber << " " << m_useRandomRunNumber);
+ if (m_useRandomRunNumber) {
+ const xAOD::EventInfo *eventInfo = evtStore()->retrieve< const xAOD::EventInfo> (m_eventInfoCollectionName);
+ if (!eventInfo) {
+ ATH_MSG_ERROR("Could not retrieve EventInfo object!");
+ sf = 1.0;
+ return CP::CorrectionCode::Error;
+ }
+ static const SG::AuxElement::Accessor<unsigned int> randomrunnumber("RandomRunNumber");
+ if (!randomrunnumber.isAvailable(*eventInfo)) {
+ sf = 1.0;
+ ATH_MSG_WARNING("Pileup tool not run before using ElectronEfficiencyTool! SFs do not reflect PU distribution in data");
+ return CP::CorrectionCode::Error;
+ }
+ runnumber = randomrunnumber(*(eventInfo));
+ }
+ ATH_MSG_DEBUG("Number of RunNumbers: " << m_RunNumbers.size() );
+
+ for ( std::size_t r=0; r<m_RunNumbers.size()-1; r++ ){
+ ATH_MSG_VERBOSE(m_RunNumbers.size()-1 << " " << m_RunNumbers.at(r) << " " << m_RunNumbers.at(r+1) << " " << runnumber);
+ if ( runnumber > (unsigned int)m_RunNumbers.at(r) && runnumber < (unsigned int)m_RunNumbers.at(r+1) ) {
+ cutRunNumber.clear();
+ cutRunNumber = Form("RunNumber%d_%d",m_RunNumbers.at(r) ,m_RunNumbers.at(r+1));
+ ATH_MSG_DEBUG(m_RunNumbers.at(r));
+ }
+ }
+
+ if (runnumber<m_RunNumbers.at(0) || (runnumber>m_RunNumbers.at(m_RunNumbers.size()-1) )) {
+ ATH_MSG_DEBUG("RunNumber not in valid RunNumber Range ");
+ sf = -9999.0;
+ return CP::CorrectionCode::OutOfValidityRange;
+ }
+
+ }
+
+ // Determine WHICH histograms to use here
+ const std::vector<TH2D*>& SShistograms = m_SF_SS.at(cutRunNumber.c_str());
+ const std::vector<TH2D*>& OShistograms = m_SF_OS.at(cutRunNumber.c_str());
+
+
+ // here check OS or SS
+ bool isOS=false;
+
+ if (truth_ele_charge * reco_ele_charge > 0) isOS=true;
+
+ if (isOS) {
+ retVal = this->getChargeFlipRate( ele_eta, ele_pt,OShistograms.at(0), sf);
+ if ( retVal != 0 ) {
+ sf = -9999.0;
+ return CP::CorrectionCode::OutOfValidityRange;
+ }
+ }
+ else {
+ ATH_MSG_DEBUG("Get SS his");
+ retVal = this->getChargeFlipRate( ele_eta, ele_pt, SShistograms.at(0), sf);
+ if ( retVal != 0 ) {
+ sf = -9999.0;
+ return CP::CorrectionCode::OutOfValidityRange;
+ }
+ }
+
+ ATH_MSG_VERBOSE("SF Rates---- . SF: " << sf );
+
+ // Systematics ------------------------------------------------------------------------------------------------------
+ double val_stat;
+
+ /// STAT
+ if (isOS) {
+ retVal = this->getChargeFlipRate( ele_eta, ele_pt,OShistograms.at(1), val_stat);
+ if ( retVal != 0 ) {
+ sf = -9999.0;
+ return CP::CorrectionCode::OutOfValidityRange;
+ }
+ }
+ else {
+ ATH_MSG_DEBUG("Get SS his");
+ retVal = this->getChargeFlipRate( ele_eta, ele_pt, SShistograms.at(1), val_stat);
+ if ( retVal != 0 ) {
+ sf = -9999.0;
+ return CP::CorrectionCode::OutOfValidityRange;
+ }
+ }
+
+ std::vector<float> systs;
+ double val_sys;
+ /// STAT
+ for (unsigned int s=2;s<OShistograms.size();s++){
+ if (isOS) {
+ retVal = this->getChargeFlipRate( ele_eta, ele_pt,OShistograms.at(s), val_sys);
+ if ( retVal != 0 ) {
+ sf = -9999.0;
+ return CP::CorrectionCode::OutOfValidityRange;
+ }
+ }
+ else {
+ ATH_MSG_DEBUG("Get SS his");
+ retVal = this->getChargeFlipRate( ele_eta, ele_pt, SShistograms.at(s), val_sys);
+ if ( retVal != 0 ) {
+ sf = -9999.0;
+ return CP::CorrectionCode::OutOfValidityRange;
+ }
+ }
+ systs.push_back(static_cast<float>(val_sys));
+ }
+
+
+ if ( m_mySysConf.size()==0 ) {
+ ATH_MSG_DEBUG(" ... nominal SF: " << sf);
+ }
+ else if (*(m_mySysConf.begin()) == SystematicVariation ("EL_CHARGEID_STAT", 1)) { sf=(sf+(val_stat)); ATH_MSG_DEBUG("SF after STATup = " << sf); }
+ else if (*(m_mySysConf.begin()) == SystematicVariation ("EL_CHARGEID_STAT", -1)) { sf=(sf-(val_stat)); ATH_MSG_DEBUG("SF after STATdown = " << sf); }
+ else {
+
+ for (unsigned int i=0;i<m_systematics.size();i++){
+ if (*(m_mySysConf.begin()) == SystematicVariation (Form("EL_CHARGEID_SYS%s",m_systematics.at(i).c_str()), 1)) { sf=(sf+(val_sys)); ATH_MSG_DEBUG("SF after SYSup = " << sf); }
+
+ if (*(m_mySysConf.begin()) == SystematicVariation (Form("EL_CHARGEID_SYS%s",m_systematics.at(i).c_str()), -1)) { sf=(sf-(val_sys)); ATH_MSG_DEBUG("SF after SYSup = " << sf); }
+
+ }
+
+ // else if (*(m_mySysConf.begin()) == SystematicVariation ("EL_CHARGEID_SYS" , -1)) { val_sf*=(1-(val_sf_sys*0.01)); ATH_MSG_DEBUG("SF after SYSdown = " << val_sf); }
+ // else ATH_MSG_ERROR("No systematic string found");
+ }
+
+ return CP::CorrectionCode::Ok;
+}
+
+
+
+//---------------------------------------------------------------------------------------
+// Decorate the electron with the scale factor
+//---------------------------------------------------------------------------------------
+
+
+CP::CorrectionCode
+CP::ElectronChargeEfficiencyCorrectionTool::applyEfficiencyScaleFactor(const xAOD::IParticle& part) const {
+ ATH_MSG_DEBUG("In CP::ElectronChargeEfficiencyCorrectionTool::applyEfficiencyScaleFactor(const xAOD::IParticle& part) const");
+ double sf = 0.0;
+ CP::CorrectionCode result = this->getEfficiencyScaleFactor(part,sf);
+ //Decorate the electron
+ (*m_sfDec)(part) = static_cast<float>(sf);
+ return result;
+
+}
+
+
+//// Giulia: This one won't exist!! ################################################ Kristin: But ok to just comment out for the moment
+
+// //---------------------------------------------------------------------------------------
+// // Get the efficiency for data only
+// //---------------------------------------------------------------------------------------
+
+// CP::CorrectionCode CP::ElectronChargeEfficiencyCorrectionTool::getDataEfficiency(const xAOD::Electron& ele, float& rate) {
+
+// const double ele_pt = ele.pt()*m_gevmev;
+// const double ele_eta = fabs(ele.eta());
+// rate = 999.;
+
+// float syst_err(-1.0);
+// float stat_err(-1.0);
+// float retVal(0.0);
+
+// retVal= this->getChargeFlipRate(ele_eta, ele_pt, m_correctionRates_data[0], rate);
+// if ( retVal != 0 ) return CP::CorrectionCode::OutOfValidityRange;
+
+// //.......sys.......
+// retVal = this->getChargeFlipRate( ele_eta, ele_pt, m_correctionRates_data[1], syst_err);
+// if ( retVal != 0 ) return CP::CorrectionCode::OutOfValidityRange;
+
+// //.......stat.......
+// retVal = this->getChargeFlipRate( ele_eta, ele_pt, m_correctionRates_data[2], stat_err);
+// if ( retVal != 0 ) return CP::CorrectionCode::OutOfValidityRange;
+
+// float total=sqrt( syst_err*syst_err + stat_err*stat_err );
+
+// ATH_MSG_VERBOSE("Rates---- data: " << rate);
+
+
+// // Systematics ------------------------------------------------------------------------------------------------------
+// if ( m_mySysConf.size() >1 ) ATH_MSG_ERROR("m_mySysConf.size() >1 !! By now just one systematic implemented!!");
+
+// if ( m_mySysConf.empty() ) ATH_MSG_VERBOSE("mySysConf is empty. NOMINAL value for rate!!!--->" << rate);
+
+// else if (*(m_mySysConf.begin()) == SystematicVariation ("ELE_ChargeMisID_STAT", 1)) {rate*=(1+(stat_err*0.01)); ATH_MSG_VERBOSE("Rate data after STATup = " << rate); }
+// else if (*(m_mySysConf.begin()) == SystematicVariation ("ELE_ChargeMisID_STAT", -1)) {rate*=(1-(stat_err*0.01)); ATH_MSG_VERBOSE("Rate data after STATdown = " << rate); }
+
+// else if (*(m_mySysConf.begin()) == SystematicVariation ("ELE_ChargeMisID_SYS" , 1)) {rate*=(1+(syst_err*0.01)); ATH_MSG_VERBOSE("Rate data after SYSup = " << rate); }
+// else if (*(m_mySysConf.begin()) == SystematicVariation ("ELE_ChargeMisID_SYS" , -1)) {rate*=(1-(syst_err*0.01)); ATH_MSG_VERBOSE("Rate data after SYSdown = " << rate); }
+
+// else if (*(m_mySysConf.begin()) == SystematicVariation ("ELE_ChargeMisID_TOT" , 1)) {rate*=(1+(total*0.01)); ATH_MSG_VERBOSE("Rate data after TOTup = " << rate); }
+// else if (*(m_mySysConf.begin()) == SystematicVariation ("ELE_ChargeMisID_TOT" , -1)) {rate*=(1-(total*0.01)); ATH_MSG_VERBOSE("Rate data after TOTdown = " << rate); }
+
+// else ATH_MSG_ERROR("No systematic string found");
+
+
+// return CP::CorrectionCode::Ok;
+// }
+
+
+
+
+/////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
+// Get the charge of the original electron
+// Giulia: not needed anymore!
+
+CP::CorrectionCode CP::ElectronChargeEfficiencyCorrectionTool::getEleTruthCharge( const xAOD::Electron& ele, int& truthcharge) const
+{
+ // if ( !(ele.isAvailable<int>("firstEgMotherPdgId")) || !(ele.isAvailable<int>("truthPdgId")) ) {
+ if ( !(ele.isAvailable<int>("firstEgMotherPdgId")) ) {
+ ATH_MSG_ERROR("Link not available for firstEgMotherPdgId...BAD!!!");
+ ATH_MSG_ERROR("Need to have present: ( !(ele.isAvailable<int>('firstEgMotherPdgId')) )" );
+ return CP::CorrectionCode::OutOfValidityRange;
+ }
+
+ //don't you need a - sign? like electron pdg is +11? or am I wrong?
+ truthcharge = (-1)*ele.auxdata<int>("firstEgMotherPdgId");
+ // Make truthcharge -1, 0, +1
+ truthcharge = (0 < truthcharge) - (truthcharge < 0);
+
+ return CP::CorrectionCode::Ok;
+
+}
+
+
+
+//Giulia: already changed -> need to be checked --> Kristin: Looks good!
+
+CP::CorrectionCode CP::ElectronChargeEfficiencyCorrectionTool::isGoodEle( const xAOD::Electron& ele, bool& goodele) const
+{
+
+ // good ele => (firstEgMotherPdgId) == 11 ## valid for both iso and conversion ele
+
+ goodele = false;
+ int firstEgPdgId = -9999;
+
+ if ( !(ele.isAvailable<int>("firstEgMotherPdgId")) ) {
+ ATH_MSG_ERROR( "firstEgMotherPdgId IS NOT AVAILABLE!!" );
+ return CP::CorrectionCode::OutOfValidityRange;
+ }
+ else {
+
+ firstEgPdgId = ele.auxdata<int>("firstEgMotherPdgId");
+
+ if ( std::abs(firstEgPdgId) != 11) {
+
+ goodele = false;
+ ATH_MSG_VERBOSE( " electron is not GOOD .... returning ...." );
+ return CP::CorrectionCode::Ok;
+
+ }
+ else {
+
+ goodele = true;
+ return CP::CorrectionCode::Ok;
+
+ }
+ }
+
+ return CP::CorrectionCode::OutOfValidityRange;
+}
+
+
+// Get the correction rate given pt (E), eta, histogram
+float CP::ElectronChargeEfficiencyCorrectionTool::getChargeFlipRate( double eta, double pt, TH2D *hrates, double& flipRate) const
+{
+ ATH_MSG_VERBOSE(" -> in: getChargeFlipRate(" << pt <<", " << eta << " TH2D, double&)");
+
+ if ( eta < m_eta_lowlimit || eta > m_eta_uplimit ) {
+
+ ATH_MSG_ERROR("Got an electron outside of the range of ETA validity " << eta);
+ return 1;
+ }
+
+ if ( pt < m_pt_lowlimit ) {
+
+ ATH_MSG_ERROR("Got an electron outside of the range of pt validity: pt lower than lower limit");
+ return 2;
+ }
+
+ if ( pt > m_pt_uplimit ) pt=m_pt_uplimit*0.999;
+
+ int bin2D = hrates->FindBin(pt, eta);
+ flipRate = hrates->GetBinContent(bin2D);
+
+ ATH_MSG_VERBOSE(" -> flipRate is " << flipRate <<", for histogram " << hrates->GetName() );
+
+
+ return 0;
+}
+
+/////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
+// returns whether this tool is affected by the given systematics
+bool CP::ElectronChargeEfficiencyCorrectionTool::isAffectedBySystematic( const SystematicVariation& systematic ) const {
+
+ CP::SystematicSet sys = affectingSystematics();
+ return sys.find (systematic) != sys.end ();
+}
+
+/////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
+// returns the list of all systematics this tool can be affected by
+
+CP::SystematicSet CP::ElectronChargeEfficiencyCorrectionTool::affectingSystematics() const{
+ CP::SystematicSet result;
+ result.insert (SystematicVariation ("EL_CHARGEID_STAT" , 1));
+ result.insert (SystematicVariation ("EL_CHARGEID_STAT" , -1));
+
+ for (unsigned int i=0;i< m_systematics.size();i++){
+ result.insert (SystematicVariation (Form("EL_CHARGEID_SYS%s",m_systematics.at(i).c_str()), 1));
+ result.insert (SystematicVariation (Form("EL_CHARGEID_SYS%s",m_systematics.at(i).c_str()), -1));
+ }
+ return result;
+}
+
+
+/////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
+// returns the list of all systematics this tool recommends to use
+CP::SystematicSet CP::ElectronChargeEfficiencyCorrectionTool::recommendedSystematics() const {
+
+ return affectingSystematics();
+
+}
+
+/////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
+// Gets a SystematicSet and filters it
+CP::SystematicCode CP::ElectronChargeEfficiencyCorrectionTool::applySystematicVariation( const SystematicSet& systConfig ) {
+
+ if (!SystematicSet::filterForAffectingSystematics(systConfig, m_affectingSys, m_mySysConf)){
+ ATH_MSG_ERROR("Unsupported combination of systematics passed to the tool! ");
+ return SystematicCode::Unsupported;
+
+ }
+
+ return SystematicCode::Ok;
+
+}
+
+/////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
diff --git a/PhysicsAnalysis/ElectronPhotonID/ElectronEfficiencyCorrection/Root/ElectronChargeEfficiencyCorrectionTool.h b/PhysicsAnalysis/ElectronPhotonID/ElectronEfficiencyCorrection/Root/ElectronChargeEfficiencyCorrectionTool.h
new file mode 100644
index 0000000000000000000000000000000000000000..68590025904bedc4f98f6c433dc683599593c41d
--- /dev/null
+++ b/PhysicsAnalysis/ElectronPhotonID/ElectronEfficiencyCorrection/Root/ElectronChargeEfficiencyCorrectionTool.h
@@ -0,0 +1,162 @@
+/*
+ Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration
+*/
+
+#ifndef ELECTRONCHARGECORRECTION__ELECTRONCHARGECORRECTIONTOOL__H
+#define ELECTRONCHARGECORRECTION__ELECTRONCHARGECORRECTIONTOOL__H
+
+#include "AsgTools/AsgTool.h"
+//#include "ElectronChargeEfficiencyCorrectionTool/IElectronChargeEfficiencyCorrectionTool.h"
+#include "AsgAnalysisInterfaces/IEfficiencyScaleFactorTool.h"
+#include "PATInterfaces/ISystematicsTool.h"
+#include "TH2.h"
+// #include "xAODTruth/TruthParticle.h"
+#include "AthContainers/AuxElement.h"
+#include "xAODEgamma/ElectronFwd.h"
+
+#include <string>
+#include <vector>
+#include <map>
+
+// Forward declaration(s):
+namespace xAOD {
+ class IParticle;
+}
+//class TH2D;
+
+
+namespace CP {
+
+ class ElectronChargeEfficiencyCorrectionTool
+ : virtual public CP::IEfficiencyScaleFactorTool,
+ public asg::AsgTool
+ {
+ ASG_TOOL_CLASS2(ElectronChargeEfficiencyCorrectionTool, CP::IEfficiencyScaleFactorTool, CP:: ISystematicsTool)
+
+ public:
+
+ /// Standard constructor
+ ElectronChargeEfficiencyCorrectionTool(const std::string name);
+
+
+ /// Standard destructor
+ virtual ~ElectronChargeEfficiencyCorrectionTool();
+
+ // METADATA // TO DO
+ // virtual StatusCode beginInputFile();
+ // virtual StatusCode beginEvent();
+ // virtual StatusCode endInputFile();
+
+
+
+ public:
+ /// Gaudi Service Interface method implementations
+ virtual StatusCode initialize();
+
+ /// Gaudi Service Interface method implementations
+ virtual StatusCode finalize();
+
+ /// Retrieve the Scale factor
+ virtual CP::CorrectionCode getEfficiencyScaleFactor(const xAOD::IParticle& part, double& sf) const;
+
+ /// Decorate the electron
+ virtual CP::CorrectionCode applyEfficiencyScaleFactor(const xAOD::IParticle& part) const;
+
+ /// Systematics
+ // void applySysVariation();
+
+ /// Returns whether this tool is affected by the given systematics
+ virtual bool isAffectedBySystematic( const SystematicVariation& systematic ) const;
+
+ /// Returns the list of all systematics this tool can be affected by
+ virtual SystematicSet affectingSystematics() const;
+
+ /// Returns the list of all systematics this tool recommends to use
+ virtual CP::SystematicSet recommendedSystematics() const;
+
+ virtual CP::SystematicCode applySystematicVariation( const SystematicSet& systConfig );
+
+ CP::SystematicCode registerSystematics();
+
+ private:
+
+ /// Get the charge flip rate rate given pt, eta, histogram
+ float getChargeFlipRate( double eta, double pt, TH2D *hrates, double& flipRate) const;
+
+ /// Get the charge of the original electron
+ CP::CorrectionCode getEleTruthCharge( const xAOD::Electron& ele, int& truthcharge) const;
+
+ /// Return true if it's prompt ele
+ CP::CorrectionCode isGoodEle( const xAOD::Electron& ele, bool& goodele) const;
+
+ /// Force the data type to a given value
+ int m_dataTypeOverwrite;
+ //
+ //
+ /// The Event info collection name
+ std::string m_eventInfoCollectionName;
+
+ /// Histogram that holds the correction rates for Monte Carlo
+ std::map<std::string, std::vector<TH2D*> > m_SF_SS; // keys (e.g. RunNumber223333_319200_Nvtx0_10_Phi1.5_1.6) mapping to vector of SF histograms --> vector m_SF: 0=nominal, 1=stat, 2,3,4...n=syst
+ std::map<std::string, std::vector<TH2D*> > m_SF_OS; // keys (e.g. RunNumber223333_319200_Nvtx0_10_Phi1.5_1.6) mapping to vector of SF histograms --> vector m_SF: 0=nominal, 1=stat, 2,3,4...n=syst
+
+ //cuts // further variables to bin in
+ std::vector<unsigned int> m_RunNumbers;
+ bool m_useRandomRunNumber;
+ unsigned int m_defaultRandomRunNumber;
+
+ /// The name of the input file that contains the histograms
+ std::string m_filename;
+
+ /// The name of the folder defining the LH/iso working point
+ std::string m_workingPoint;
+
+ /// Lower limit of eta range where corrections are available; taken from histogram
+ double m_eta_lowlimit;
+
+ /// Upper limit of eta range where corrections are available; taken from histogram
+ double m_eta_uplimit;
+
+ /// Lower limit of pt range where corrections are available; taken from histogram
+ double m_pt_lowlimit;
+
+ /// Upper limit of pt range where corrections are available; taken from histogram
+ double m_pt_uplimit;
+
+ /// Factor for GeV <-> MeV switching
+ float m_gevmev;
+
+
+ /// Truth charge
+ mutable int m_truthCharge;
+
+ //const xAOD::TruthParticle *m_truthparticle;
+
+ /// Random number generator for throwing the dice
+ // TRandom3 *m_Rndm;
+
+ // Systematics
+ std::vector<std::string> m_systematics;
+
+ std::map<CP::SystematicSet, CP::SystematicSet> m_filtered_sys_sets;
+ //boost::unordered_map<SystematicSet,EffiCollection*> m_sf_sets;
+ CP::SystematicSet m_mySysConf;
+ CP::SystematicSet m_affectingSys;
+
+ /// Decorator
+ std::string m_sf_decoration_name;
+ SG::AuxElement::Decorator< float >* m_sfDec;
+
+
+ };
+
+} // End namespace CP
+
+#endif
+
+/// Apply the correction on a modifyable object
+//......virtual CP::CorrectionCode applyCorrection( xAOD::Electron& ele );
+
+/// Create a corrected copy from a constant muon
+//...virtual CP::CorrectionCode correctedCopy( const xAOD::Electron& input,
+//... xAOD::Electron*& output );
diff --git a/PhysicsAnalysis/ElectronPhotonID/ElectronEfficiencyCorrection/Root/TElectronEfficiencyCorrectionTool.cxx b/PhysicsAnalysis/ElectronPhotonID/ElectronEfficiencyCorrection/Root/TElectronEfficiencyCorrectionTool.cxx
index 22a734b19cf478751cd73d3448b74eeedf09a383..78485f96373222ae52e517b191e67e215c97e82e 100644
--- a/PhysicsAnalysis/ElectronPhotonID/ElectronEfficiencyCorrection/Root/TElectronEfficiencyCorrectionTool.cxx
+++ b/PhysicsAnalysis/ElectronPhotonID/ElectronEfficiencyCorrection/Root/TElectronEfficiencyCorrectionTool.cxx
@@ -71,7 +71,7 @@ namespace mapkey{
Root::TElectronEfficiencyCorrectionTool::TElectronEfficiencyCorrectionTool(const char *name) :
Root::TCalculatorToolBase(name),
asg::AsgMessaging(std::string(name)),
- m_Rndm(0),
+ m_Rndm(),
m_randomCounter(0),
m_isInitialized(false),
m_detailLevel(2),
@@ -80,8 +80,6 @@ Root::TElectronEfficiencyCorrectionTool::TElectronEfficiencyCorrectionTool(const
m_doToyMC(false),
m_doCombToyMC(false),
m_nToyMC(0),
- m_uncorrToyMCSystFull(0),
- m_uncorrToyMCSystFast(0),
m_nSys(0),
m_nSysMax(0),
m_runNumBegin(0),
@@ -97,6 +95,12 @@ Root::TElectronEfficiencyCorrectionTool::TElectronEfficiencyCorrectionTool(const
m_nSimpleUncorrSyst(0),
m_position_globalBinNumber(0)
{
+ //Setup the keys
+ m_keys.push_back(mapkey::sf);
+ m_keys.push_back(mapkey::stat);
+ m_keys.push_back(mapkey::eig);
+ m_keys.push_back(mapkey::uncorr);
+ //
}
// =============================================================================
@@ -105,95 +109,63 @@ Root::TElectronEfficiencyCorrectionTool::TElectronEfficiencyCorrectionTool(const
Root::TElectronEfficiencyCorrectionTool::~TElectronEfficiencyCorrectionTool() {
- for (auto const &tempit : m_histList) {
+ //Need some gymnastic to make sure that the TObjArray elements are owned
+ //and deleted ...
+ for (auto &tempit : m_histList) {
for (unsigned int i = 0; i < tempit.second.size(); ++i) {
- if (tempit.second.at(i)) {
- tempit.second.at(i)->SetOwner(kTRUE);
- delete tempit.second.at(i);
- }
+ tempit.second.at(i).SetOwner(kTRUE);
}
}
-
- for (auto const &tempit : m_fastHistList) {
+ for (auto &tempit : m_fastHistList) {
for (unsigned int i = 0; i < tempit.second.size(); ++i) {
- if (tempit.second.at(i)) {
- tempit.second.at(i)->SetOwner(kTRUE);
- delete tempit.second.at(i);
- }
+ tempit.second.at(i).SetOwner(kTRUE);
}
}
-
- for (auto const &tempit : m_sysList) {
- for (auto const &i : tempit) {
- i->SetOwner(kTRUE);
- delete i;
+ for (auto &tempit : m_sysList) {
+ for (auto &i : tempit) {
+ i.SetOwner(kTRUE);
}
}
- for (auto const &tempit : m_fastSysList) {
- for (auto const &i : tempit) {
- i->SetOwner(kTRUE);
- delete i;
+ for (auto &tempit : m_fastSysList) {
+ for (auto &i : tempit) {
+ i.SetOwner(kTRUE);
}
}
-
- if (m_Rndm) {
- delete m_Rndm;
- }
- //These are pointers to vectors of vector pointers to TobjectArray pointers ...
- //Nighmare to delete. They are not always allocated only under certain conditions
- if(m_uncorrToyMCSystFull){
- for (auto const &tempit : *m_uncorrToyMCSystFull) {
-
- for (auto const& i : *tempit) {
- i->SetOwner(kTRUE);
- delete i;
- }
- delete tempit;
+ for (auto &tempit : m_uncorrToyMCSystFull) {
+ for (auto & i : tempit) {
+ i.SetOwner(kTRUE);
}
- delete m_uncorrToyMCSystFull;
- }
-
- if(m_uncorrToyMCSystFast){
- for (auto const &tempit : *m_uncorrToyMCSystFast) {
- for (auto const& i : *tempit) {
- i->SetOwner(kTRUE);
- delete i;
- }
- delete tempit;
+ }
+ for (auto &tempit : m_uncorrToyMCSystFast) {
+ for (auto & i : tempit) {
+ i.SetOwner(kTRUE);
}
- delete m_uncorrToyMCSystFast;
}
}
-
-
// =============================================================================
// Initialize this correction tool
// =============================================================================
-int
-Root::TElectronEfficiencyCorrectionTool::initialize() {
+int Root::TElectronEfficiencyCorrectionTool::initialize() {
// use an int as a StatusCode
int sc(1);
-
ATH_MSG_DEBUG(" (file: " << __FILE__ << ", line: "
<< __LINE__ << ") " << "Debug flag set. Printing verbose output!");
+ //Avoid double init
if (m_isInitialized) {
ATH_MSG_ERROR(" (file: " << __FILE__ << ", line: " << __LINE__ << ") " << "tool initialized twice!");
return 0;
}
-
+ //Check if files are present
if (m_corrFileNameList.size() == 0) {
ATH_MSG_ERROR(" (file: " << __FILE__ << ", line: " << __LINE__ << ") " << " No file added!");
return 0;
}
-
ATH_MSG_INFO("Initializing tool with " << m_corrFileNameList.size() << " configuration file(s)");
-
- ///////////////////////////////////////////////////
+ //
// Check if the first file can be opened (needed for auto-setting of the seed based on the md5-sum of the file)
const std::unique_ptr<char> fname(gSystem->ExpandPathName(m_corrFileNameList[0].c_str()));
-
std::unique_ptr<TFile> rootFile_tmp = CxxUtils::make_unique<TFile> (fname.get(), "READ");
if (!rootFile_tmp) {
ATH_MSG_ERROR(
@@ -202,35 +174,29 @@ Root::TElectronEfficiencyCorrectionTool::initialize() {
}
rootFile_tmp->Close();
// close it back again
- ///////////////////////////////////////////////////
-
+ //
+ //invalid input requested
if (m_doToyMC && m_doCombToyMC) {
ATH_MSG_ERROR(" (file: " << __FILE__ << ", line: " << __LINE__ << ") " << " Both regular and combined toy MCs booked!"
<< " Only use one!");
return 0;
}
- m_keys.push_back(mapkey::sf);
- m_keys.push_back(mapkey::stat);
- m_keys.push_back(mapkey::eig);
- m_keys.push_back(mapkey::uncorr);
- //
- // initialize the random number generated if toyMC propagation booked
+
+ // initialize the random number generator if toyMC propagation booked
if (m_doToyMC || m_doCombToyMC) {
if (m_seed == 0) {
std::unique_ptr<TMD5> tmd=CxxUtils::make_unique<TMD5>();
- m_seed = *reinterpret_cast<const int *>(tmd->FileChecksum(fname.get())->AsString());
+ const char* tmd_as_string=tmd->FileChecksum(fname.get())->AsString();
+ m_seed = *reinterpret_cast<const int*>(tmd_as_string);
ATH_MSG_INFO("Seed (automatically) set to " << m_seed);
}else {
ATH_MSG_INFO("Seed set to " << m_seed);
}
- m_Rndm = new TRandom3(m_seed);
-
- m_uncorrToyMCSystFull = new std::vector< std::vector<TObjArray *> * >;
- m_uncorrToyMCSystFast = new std::vector< std::vector<TObjArray *> *>;
+ m_Rndm= TRandom3(m_seed);
}
//
// Load the needed histograms
- if (0 == this->getHistograms()) {
+ if (0 == getHistograms()) {
ATH_MSG_ERROR(" (file: " << __FILE__ << ", line: " << __LINE__ << ") " << "! Problem when calling getHistograms()");
return 0;
}
@@ -244,10 +210,10 @@ Root::TElectronEfficiencyCorrectionTool::initialize() {
ATH_MSG_DEBUG(" (file: " << __FILE__ << ", line: " << __LINE__ << ") "
<< "Found " << nRunNumbersFull << " run number ranges for full sim with a total of " <<
m_histList[mapkey::sf].size() << " scale factor histograms.");
-
+ //
// --------------------------------------------------------------------------
- // Register the cuts and check that the registration worked:
- // NOTE: THE ORDER IS IMPORTANT!!! Cut0 corresponds to bit 0, Cut1 to bit 1,...
+ // Register the output / results
+ // NOTE: THE ORDER IS IMPORTANT!!!
m_position_eff = m_result.addResult((m_resultPrefix + m_resultName).c_str(), "efficiency scale factor");
if (m_position_eff < 0) {
sc = 0; // Exceeded the number of allowed results
@@ -305,7 +271,7 @@ Root::TElectronEfficiencyCorrectionTool::initialize() {
m_result.setResult(m_position_nSys, 0);
m_result.setResult(m_position_uncorrSys, 0);
m_isInitialized = kTRUE;
-
+ // --------------------------------------------------------------------------
ATH_MSG_INFO("Tool succesfully initialized!");
return sc;
}
@@ -334,8 +300,7 @@ Root::TElectronEfficiencyCorrectionTool::calculate(const PATCore::ParticleDataTy
" Tool not correctly initialized or no scale factor file given. ");
return m_result;
}
- ATH_MSG_DEBUG(
- " (file: " << __FILE__ << ", line: " << __LINE__ << ") " << this->getName() << " n Systematics: " << m_nSysMax);
+ ATH_MSG_DEBUG("(file: " << __FILE__ << ", line: " << __LINE__ << ") " << this->getName() << " n Systematics: " << m_nSysMax);
// Reset the results to default values
m_result.setResult(m_position_eff, -999.0);
m_result.setResult(m_position_err, 1.0);
@@ -388,8 +353,8 @@ Root::TElectronEfficiencyCorrectionTool::calculate(const PATCore::ParticleDataTy
//The vector has as many entries as supported run period
//The TobjArray has 2D histos for high, low et, or forward
//The 2D Histo has the number we want.
- const std::map<int, std::vector< TObjArray* > >& currentmap = (isFastSim)? m_fastHistList : m_histList;
- std::map<int, std::vector< TObjArray* > >::const_iterator currentVector_itr = currentmap.find(mapkey::sf); //find the vector
+ const std::map<int, std::vector< TObjArray > >& currentmap = (isFastSim)? m_fastHistList : m_histList;
+ std::map<int, std::vector< TObjArray > >::const_iterator currentVector_itr = currentmap.find(mapkey::sf); //find the vector
//See if we can find a SF vector in the map and the corresponding TobjArray for this run period
if (currentVector_itr == currentmap.end()) {
if (this->msgLvl(MSG::DEBUG)){
@@ -401,7 +366,7 @@ Root::TElectronEfficiencyCorrectionTool::calculate(const PATCore::ParticleDataTy
return m_result;
}
//Get a reference (synonym) to this vector
- const std::vector<TObjArray*>& currentVector=currentVector_itr->second;
+ const std::vector<TObjArray>& currentVector=currentVector_itr->second;
if (currentVector.size()<=0 || runnumberIndex>= static_cast <signed> (currentVector.size())) {
if (this->msgLvl(MSG::DEBUG)){
printDefaultReturnMessage(TString::Format(
@@ -411,15 +376,13 @@ Root::TElectronEfficiencyCorrectionTool::calculate(const PATCore::ParticleDataTy
}
return m_result;
}
-
//
//
//Now we can get the corresponding Object array
- const TObjArray* currentObjectArray = currentVector.at(runnumberIndex);
- //
+ const TObjArray& currentObjectArray = currentVector.at(runnumberIndex);
//
//This is the number of entries in the TObjArray
- const int entries = currentObjectArray->GetEntries();
+ const int entries = currentObjectArray.GetEntries();
//Find the correct histogram in the TObjArray (Low, high Et or forward)
double xValue, yValue;
xValue = et;
@@ -431,7 +394,7 @@ Root::TElectronEfficiencyCorrectionTool::calculate(const PATCore::ParticleDataTy
TH2 *tmpHist(0);
for (int i = 0; i < entries ; ++i) {
block = kFALSE;
- tmpHist = (TH2 *) (currentObjectArray->At(i));
+ tmpHist = (TH2 *) (currentObjectArray.At(i));
//block if we are below minimum et
if (et < tmpHist->GetXaxis()->GetXmin()) {
@@ -489,7 +452,7 @@ Root::TElectronEfficiencyCorrectionTool::calculate(const PATCore::ParticleDataTy
//Now we have the index of the histogram for this region in the TObjectarray
TH2* currentHist(0);
if (index >= 0) {
- currentHist = static_cast<TH2*> (currentObjectArray->At(index));
+ currentHist = static_cast<TH2*> (currentObjectArray.At(index));
}
else {
if (this->msgLvl(MSG::DEBUG)){
@@ -526,7 +489,7 @@ Root::TElectronEfficiencyCorrectionTool::calculate(const PATCore::ParticleDataTy
if (currentVector_itr != currentmap.end()) {
//itr at the location of the vector, .second get the vector, at(runnumberIndex is the TObjectArray
// for the period , finaly get the hist at index (from above).
- TH1 *stat = static_cast<TH1*>(currentVector_itr->second.at(runnumberIndex)->At(index));
+ TH1 *stat = static_cast<TH1*>(currentVector_itr->second.at(runnumberIndex).At(index));
statErr = stat->GetBinContent(globalBinNumber);
m_result.setResult(m_position_statErr, statErr);
}
@@ -536,10 +499,10 @@ Root::TElectronEfficiencyCorrectionTool::calculate(const PATCore::ParticleDataTy
currentVector_itr = currentmap.find(mapkey::eig); //find the vector
if (currentVector_itr != currentmap.end()) {
//Check on the ObjArray
- if (currentVector_itr->second.at(runnumberIndex) != 0) {
+ if (currentVector_itr->second.at(runnumberIndex).GetEntriesFast()> 0) {
//itr at the location of the vector, .second get the vector, at(runnumberIndex is the TObjectArray
// for the period , finaly get the hist at index (from above).
- TH1 *eig = static_cast<TH1*>(currentVector_itr->second.at(runnumberIndex)->At(index));
+ TH1 *eig = static_cast<TH1*>(currentVector_itr->second.at(runnumberIndex).At(index));
m_sLevel[Root::TElectronEfficiencyCorrectionTool::simple] = 0;
m_sLevel[Root::TElectronEfficiencyCorrectionTool::medium] = 0;
m_sLevel[Root::TElectronEfficiencyCorrectionTool::detailed] = 0;
@@ -568,29 +531,25 @@ Root::TElectronEfficiencyCorrectionTool::calculate(const PATCore::ParticleDataTy
//And the obj array for high low etc
//We need to see if we can do something better here ...
//
- std::vector< std::vector< TObjArray * > > &sysList = (isFastSim) ? m_fastSysList : m_sysList;
+ std::vector< std::vector< TObjArray > > &sysList = (isFastSim) ? m_fastSysList : m_sysList;
static std::vector<double> corrSys(10); // Avoid some re-allocations of memory
corrSys.clear();
if (sysList.size() > (unsigned int) index) {
if (sysList.at(index).size() > (unsigned int) runnumberIndex) {
- const int sys_entries = sysList.at(index).at( runnumberIndex)->GetEntries();
+ const int sys_entries = sysList.at(index).at( runnumberIndex).GetEntries();
for (int sys = 0; sys < sys_entries; ++sys) {
- tmpHist = (TH2 *) sysList.at(index).at(runnumberIndex)->At(sys_entries - 1 - sys);
+ tmpHist = (TH2 *) sysList.at(index).at(runnumberIndex).At(sys_entries - 1 - sys);
corrSys.push_back(tmpHist->GetBinContent(globalBinNumber));
m_result.setResult(m_position_corrSys[(sys_entries - 1 - sys)], corrSys[sys]);
}
-
if (m_position_corrSys.size() > 0 && sys_entries<=1) {
if (m_result.getResult(m_position_corrSys[0]) == 0) {
m_result.setResult(m_position_corrSys[0], scaleFactorErr);
}
}
-
}
}
-
- //
//
//Do the uncorr error using the available info from above i.e index etc
//////////////////////////////////////////////////////////////////
@@ -599,8 +558,8 @@ Root::TElectronEfficiencyCorrectionTool::calculate(const PATCore::ParticleDataTy
currentVector_itr = currentmap.find(mapkey::uncorr); //find the vector
if (currentVector_itr != currentmap.end()) {
//Check on the ObjArray
- if (currentVector_itr->second.at(runnumberIndex) != 0) {
- TH1 *uncorr = static_cast<TH1*>(currentVector_itr->second.at(runnumberIndex)->At(index));
+ if (currentVector_itr->second.at(runnumberIndex).GetEntriesFast()>0) {
+ TH1 *uncorr = static_cast<TH1*>(currentVector_itr->second.at(runnumberIndex).At(index));
double valAdd = uncorr->GetBinContent(globalBinNumber);
val = sqrt(val * val + valAdd * valAdd);
for (int i = 0; i < m_sLevel[m_detailLevel]; ++i) {
@@ -617,23 +576,15 @@ Root::TElectronEfficiencyCorrectionTool::calculate(const PATCore::ParticleDataTy
//Here are the toys
//////////////////////////
if (m_doToyMC || m_doCombToyMC) {
- std::vector< std::vector<TObjArray *> *> *toyMCList;
- if (isFastSim) {
- toyMCList = m_uncorrToyMCSystFast;
- } else {
- toyMCList = m_uncorrToyMCSystFull;
- }
- if (toyMCList != 0) {
- if (toyMCList->size() > (unsigned int) runnumberIndex) {
+ std::vector<std::vector<TObjArray > >& toyMCList = ((isFastSim) ? m_uncorrToyMCSystFast : m_uncorrToyMCSystFull);
+ if (toyMCList.size() > (unsigned int) runnumberIndex) {
for (int toy = 0; toy < m_nToyMC; toy++) {
- if (toyMCList->at(runnumberIndex)->at(toy)->GetLast() >= index) {
+ if (toyMCList.at(runnumberIndex).at(toy).GetLast() >= index) {
m_result.setResult(m_position_uncorrToyMCSF.at(toy),
- ((TH2 *) toyMCList->at(runnumberIndex)->at(toy)->At(index))->GetBinContent(
- globalBinNumber));
+ ((TH2 *) toyMCList.at(runnumberIndex).at(toy).At(index))->GetBinContent(globalBinNumber));
}
}
}
- }
}
m_result.setResult(m_position_globalBinNumber, globalBinNumber);
return m_result;
@@ -644,9 +595,11 @@ Root::TElectronEfficiencyCorrectionTool::calculate(const PATCore::ParticleDataTy
// =============================================================================
void
Root::TElectronEfficiencyCorrectionTool::calcDetailLevels(TH1D *eig) {
+
m_sLevel[Root::TElectronEfficiencyCorrectionTool::simple] = 0;
m_sLevel[Root::TElectronEfficiencyCorrectionTool::medium] = 0;
m_sLevel[Root::TElectronEfficiencyCorrectionTool::detailed] = 0;
+
int nSys = eig->GetNbinsX() - 1;
double sign = 0;
// Calculate detail level
@@ -665,14 +618,14 @@ Root::TElectronEfficiencyCorrectionTool::calcDetailLevels(TH1D *eig) {
// =============================================================================
// Build the toyMC tables from inputs
// =============================================================================
-std::vector<TH2D *> *
-Root::TElectronEfficiencyCorrectionTool::buildSingleToyMC(TH2D *sf, TH2D *stat, TH2D *uncorr, TObjArray *corr) {
+std::vector<TH2D *>
+Root::TElectronEfficiencyCorrectionTool::buildSingleToyMC(TH2D *sf, TH2D *stat, TH2D *uncorr, const TObjArray& corr) {
ATH_MSG_DEBUG(" (file: " << __FILE__ << ", line: " << __LINE__ << ")! " << "entering function buildSingleToyMC");
- std::vector<TH2D*>* tmpHists = new std::vector<TH2D*>;
+ std::vector<TH2D*> tmpHists;
int nBins = (stat->GetNbinsX() + 2) * (stat->GetNbinsY() + 2);
for (int toy = 0; toy < m_nToyMC; toy++) {
- tmpHists->push_back((TH2D *) corr->At(0)->Clone());
+ tmpHists.push_back((TH2D *) corr.At(0)->Clone());
}
// Loop over all bins
for (int bin = 0; bin < nBins; bin++) {
@@ -686,15 +639,15 @@ Root::TElectronEfficiencyCorrectionTool::buildSingleToyMC(TH2D *sf, TH2D *stat,
// Add smaller correlated systematics as uncorrelated
for (int i = 0; i < m_sLevel[m_detailLevel]; ++i) {
- if (corr->At(i) != 0) {
- double valAdd = ((TH2D *) corr->At(i))->GetBinContent(bin);
+ if (corr.At(i) != 0) {
+ double valAdd = ((TH2D *) corr.At(i))->GetBinContent(bin);
val = sqrt(val * val + valAdd * valAdd);
}
}
for (int toy = 0; toy < m_nToyMC; toy++) {
- tmpHists->at(toy)->SetBinContent(bin, (val * m_Rndm->Gaus(0, 1)) + sf->GetBinContent(bin));
+ tmpHists.at(toy)->SetBinContent(bin, (val * m_Rndm.Gaus(0, 1)) + sf->GetBinContent(bin));
m_randomCounter++;
- tmpHists->at(toy)->SetDirectory(0);
+ tmpHists.at(toy)->SetDirectory(0);
}
}
return tmpHists;
@@ -703,18 +656,18 @@ Root::TElectronEfficiencyCorrectionTool::buildSingleToyMC(TH2D *sf, TH2D *stat,
// =============================================================================
// Build the combined toyMC tables from inputs
// =============================================================================
-TH2D *
-Root::TElectronEfficiencyCorrectionTool::buildSingleCombToyMC(TH2D *sf, TH2D *stat, TH2D *uncorr, TObjArray *corr) {
+TH2D*
+Root::TElectronEfficiencyCorrectionTool::buildSingleCombToyMC(TH2D *sf, TH2D *stat, TH2D *uncorr, const TObjArray& corr) {
ATH_MSG_DEBUG(" (file: " << __FILE__ << ", line: " << __LINE__ << ") " << "entering function buildSingleCombToyMC");
TH2D *tmpHist;
int nBins = (stat->GetNbinsX() + 2) * (stat->GetNbinsY() + 2);
- tmpHist = (TH2D *) corr->At(0)->Clone();
+ tmpHist = (TH2D *) corr.At(0)->Clone();
// Create random numbers for the corr. uncertainties
std::vector<double> rnd (m_nSys,0);
for (int s = 0; s < m_nSys; s++) {
- rnd[s] = m_Rndm->Gaus(0, 1);
+ rnd[s] = m_Rndm.Gaus(0, 1);
m_randomCounter++;
}
@@ -730,19 +683,19 @@ Root::TElectronEfficiencyCorrectionTool::buildSingleCombToyMC(TH2D *sf, TH2D *st
// Add smaller correlated systematics as uncorrelated
for (int s = 0; s < m_sLevel[m_detailLevel]; s++) {
- if (corr->At(s) != 0) {
- double valAdd = ((TH2D *) corr->At(s))->GetBinContent(bin);
+ if (corr.At(s) != 0) {
+ double valAdd = ((TH2D *) corr.At(s))->GetBinContent(bin);
val = sqrt(val * val + valAdd * valAdd);
}
}
- val = val * m_Rndm->Gaus(0, 1);
+ val = val * m_Rndm.Gaus(0, 1);
m_randomCounter++;
// Add larger correlated systematics
for (int s = m_sLevel[m_detailLevel]; s < m_nSys; s++) {
- if (corr->At(s) != 0) {
- val += ((TH2D *) corr->At(s))->GetBinContent(bin) * rnd[s];
+ if (corr.At(s) != 0) {
+ val += ((TH2D *) corr.At(s))->GetBinContent(bin) * rnd[s];
}
}
tmpHist->SetBinContent(bin, val + sf->GetBinContent(bin));
@@ -755,48 +708,44 @@ Root::TElectronEfficiencyCorrectionTool::buildSingleCombToyMC(TH2D *sf, TH2D *st
// =============================================================================
// Build the toyMC tables from inputs
// =============================================================================
-std::vector<TObjArray *> *
-Root::TElectronEfficiencyCorrectionTool::buildToyMCTable(TObjArray *sf, TObjArray *eig, TObjArray *stat,
- TObjArray *uncorr, std::vector<TObjArray *> *corr) {
+std::vector<TObjArray>
+Root::TElectronEfficiencyCorrectionTool::buildToyMCTable(const TObjArray& sf, const TObjArray& eig, const TObjArray& stat,
+ const TObjArray& uncorr, const std::vector<TObjArray>& corr) {
+
ATH_MSG_DEBUG(" (file: " << __FILE__ << ", line: " << __LINE__ << ") " << "entering function buildToyMCTable");
+ std::vector<TObjArray> tmpVec;
+ const int stat_entries = stat.GetEntries();
- const int stat_entries = stat->GetEntries();
- std::vector<TObjArray *>* tmpVec= new std::vector<TObjArray*>;
- TObjArray *tmpArray;
if (m_doCombToyMC) {
for (int toyMC = 0; toyMC < m_nToyMC; toyMC++) {
- tmpArray = new TObjArray();
+ TObjArray tmpArray;
for (int i = 0; i < stat_entries; ++i) {
- if (eig != 0 && uncorr != 0) {
- calcDetailLevels((TH1D *) eig->At(i));
- tmpArray->Add(buildSingleCombToyMC((TH2D *) sf->At(i), (TH2D *) stat->At(i), (TH2D *) uncorr->At(i),
- corr->at(i)));
+ if (eig.GetEntriesFast()>0 && uncorr.GetEntriesFast()>0) {
+ calcDetailLevels((TH1D *) eig.At(i));
+ tmpArray.Add(buildSingleCombToyMC((TH2D *) sf.At(i), (TH2D *) stat.At(i), (TH2D *) uncorr.At(i), corr.at(i)));
}else {
- tmpArray->Add(buildSingleCombToyMC((TH2D *) sf->At(i), (TH2D *) stat->At(i), 0, corr->at(i)));
+ tmpArray.Add(buildSingleCombToyMC((TH2D *) sf.At(i), (TH2D *) stat.At(i), 0, corr.at(i)));
}
}
- tmpVec->push_back(tmpArray);
+ tmpVec.push_back(tmpArray);
}
}else {
- std::vector< std::vector<TH2D *> *> *tmpVec2 = new std::vector< std::vector<TH2D *> *>;
+ std::vector< std::vector<TH2D*> > tmpVec2 ;
for (int i = 0; i < stat_entries; ++i) {
- calcDetailLevels((TH1D *) eig->At(i));
- tmpVec2->push_back(buildSingleToyMC((TH2D *) sf->At(i), (TH2D *) stat->At(i), (TH2D *) uncorr->At(i),
- corr->at(i)));
+ calcDetailLevels((TH1D *) eig.At(i));
+ tmpVec2.push_back(buildSingleToyMC((TH2D *) sf.At(i), (TH2D *) stat.At(i),
+ (TH2D *) uncorr.At(i),corr.at(i)));
}
for (int toy = 0; toy < m_nToyMC; toy++) {
- tmpArray = new TObjArray();
- for (unsigned int i = 0; i < tmpVec2->size(); ++i) {
- tmpArray->Add(tmpVec2->at(i)->at(toy));
+ TObjArray tmpArray;
+ for (unsigned int i = 0; i < tmpVec2.size(); ++i) {
+ tmpArray.Add(tmpVec2.at(i).at(toy));
}
- tmpVec->push_back(tmpArray);
+ tmpVec.push_back(tmpArray);
}
- delete tmpVec2;
}
-
ATH_MSG_DEBUG(" (file: " << __FILE__ << ", line: " << __LINE__ << ") " << "m_Rndm->Uniform(0, 1) after throwing " << m_randomCounter
- << " random numbers: " << m_Rndm->Uniform(0,
- 1));
+ << " random numbers: " << m_Rndm.Uniform(0,1));
return tmpVec;
}
@@ -805,8 +754,7 @@ Root::TElectronEfficiencyCorrectionTool::buildToyMCTable(TObjArray *sf, TObjArra
// =============================================================================
int
Root::TElectronEfficiencyCorrectionTool::getNbins(std::map<float, std::vector<float> > &pt_eta1) {
- std::vector<TObjArray *> tmpVec;
-
+ std::vector<TObjArray > tmpVec;
tmpVec = m_histList[mapkey::sf];
int nbinsTotal = 0;
pt_eta1.clear();
@@ -814,45 +762,41 @@ Root::TElectronEfficiencyCorrectionTool::getNbins(std::map<float, std::vector<fl
eta1.clear();
for (unsigned int ikey = 0; ikey < tmpVec.size(); ++ikey) {
- for (int entries = 0; entries < (tmpVec.at(ikey))->GetEntries(); entries++) {
+ for (int entries = 0; entries < (tmpVec.at(ikey)).GetEntries(); entries++) {
eta1.clear();
- TH2D *h_tmp = ((TH2D * ) (tmpVec.at(ikey))->At(entries));
+ TH2D *h_tmp = ((TH2D * ) (tmpVec.at(ikey)).At(entries));
int nbinsX = h_tmp->GetNbinsX();
int nbinsY = h_tmp->GetNbinsY();
for (int biny = 1; biny <= nbinsY + 1; biny++) {
eta1.push_back(h_tmp->GetYaxis()->GetBinLowEdge(biny));
- if (entries == (tmpVec.at(ikey))->GetEntries() - 1) {
+ if (entries == (tmpVec.at(ikey)).GetEntries() - 1) {
eta1.push_back(h_tmp->GetYaxis()->GetBinLowEdge(biny + 1));
}
}
for (int binx = 1; binx <= nbinsX; binx++) {
pt_eta1[h_tmp->GetXaxis()->GetBinLowEdge(binx)] = eta1;
- if (entries == (tmpVec.at(ikey))->GetEntries() - 1) {
+ if (entries == (tmpVec.at(ikey)).GetEntries() - 1) {
pt_eta1[h_tmp->GetXaxis()->GetBinLowEdge(binx + 1)] = eta1;
}
}
}
}
-
for (auto &i : pt_eta1) {
nbinsTotal += i.second.size();
}
-
return nbinsTotal;
}
// =============================================================================
// Helper function to retrieve the position of the first toy MC scale factor
// =============================================================================
-int
-Root::TElectronEfficiencyCorrectionTool::getFirstToyMCPosition() {
+int Root::TElectronEfficiencyCorrectionTool::getFirstToyMCPosition() {
if (!m_isInitialized) {
ATH_MSG_ERROR(" (file: " << __FILE__ << ", line: " << __LINE__ << ") " << "Tool not initialized.");
return -1;
}
-
if (m_nToyMC > 0) {
return m_position_uncorrToyMCSF.at(0);
}else {
@@ -862,13 +806,11 @@ Root::TElectronEfficiencyCorrectionTool::getFirstToyMCPosition() {
}
}
-int
-Root::TElectronEfficiencyCorrectionTool::getLastToyMCPosition() {
+int Root::TElectronEfficiencyCorrectionTool::getLastToyMCPosition() {
if (!m_isInitialized) {
ATH_MSG_ERROR(" (file: " << __FILE__ << ", line: " << __LINE__ << ") " << "Tool not initialized.");
return -1;
}
-
if (m_nToyMC > 0) {
return m_position_uncorrToyMCSF.at(m_nToyMC - 1);
} else {
@@ -877,14 +819,13 @@ Root::TElectronEfficiencyCorrectionTool::getLastToyMCPosition() {
return 0;
}
}
-
-int
-Root::TElectronEfficiencyCorrectionTool::getFirstCorrSysPosition() {
+//================================================================================
+/// Helpers
+int Root::TElectronEfficiencyCorrectionTool::getFirstCorrSysPosition() {
if (!m_isInitialized) {
ATH_MSG_ERROR(" (file: " << __FILE__ << ", line: " << __LINE__ << ") " << "Tool not initialized.");
return -1;
}
-
if (m_position_corrSys.size() > 0) {
return m_position_corrSys.at(0);
}else {
@@ -894,8 +835,7 @@ Root::TElectronEfficiencyCorrectionTool::getFirstCorrSysPosition() {
}
}
-int
-Root::TElectronEfficiencyCorrectionTool::getLastCorrSysPosition() {
+int Root::TElectronEfficiencyCorrectionTool::getLastCorrSysPosition() {
if (!m_isInitialized) {
ATH_MSG_ERROR(" (file: " << __FILE__ << ", line: " << __LINE__ << ") " << "Tool not initialized.");
return -1;
@@ -913,9 +853,7 @@ Root::TElectronEfficiencyCorrectionTool::getLastCorrSysPosition() {
return 0;
}
}
-
-int
-Root::TElectronEfficiencyCorrectionTool::getGlobalBinNumberPosition() {
+int Root::TElectronEfficiencyCorrectionTool::getGlobalBinNumberPosition() {
if (!m_isInitialized) {
ATH_MSG_ERROR(" (file: " << __FILE__ << ", line: " << __LINE__ << ") " << "Tool not initialized.");
return -1;
@@ -923,8 +861,7 @@ Root::TElectronEfficiencyCorrectionTool::getGlobalBinNumberPosition() {
return m_position_globalBinNumber;
}
-void
-Root::TElectronEfficiencyCorrectionTool::printResultMap() {
+void Root::TElectronEfficiencyCorrectionTool::printResultMap() {
if (!m_isInitialized) {
ATH_MSG_ERROR(" (file: " << __FILE__ << ", line: " << __LINE__ << ") " << "Tool not initialized.");
return;
@@ -934,17 +871,16 @@ Root::TElectronEfficiencyCorrectionTool::printResultMap() {
ATH_MSG_INFO(pos << " \t \t - \t " << m_result.getResultDescription(pos));
}
}
+//================================================================================
// =============================================================================
// Get the input histograms from the input files
// =============================================================================
-int
-Root::TElectronEfficiencyCorrectionTool::getHistograms() {
- ATH_MSG_DEBUG(" (file: " << __FILE__ << ", line: " << __LINE__ << ") " << "entering function getHistograms");
+int Root::TElectronEfficiencyCorrectionTool::getHistograms() {
+ ATH_MSG_DEBUG(" (file: " << __FILE__ << ", line: " << __LINE__ << ") " << "entering function getHistograms");
// Cache the current directory in the root file
TDirectory *origDir = gDirectory;
-
// -------------------------------------------------------
// Get the name of the first input ROOT file and
// interpret from that what we have:
@@ -974,11 +910,9 @@ Root::TElectronEfficiencyCorrectionTool::getHistograms() {
<< " and resultName: " << m_resultName);
}
-
// -------------------------------------------------------
// Get all ROOT files and histograms
// -------------------------------------------------------
-
for (unsigned int i = 0; i < m_corrFileNameList.size(); ++i) {
// Load the ROOT file
const std::unique_ptr<char> fname (gSystem->ExpandPathName(m_corrFileNameList[i].c_str()));
@@ -989,7 +923,6 @@ Root::TElectronEfficiencyCorrectionTool::getHistograms() {
m_corrFileNameList[i]);
return 0;
}
-
// Loop over all directories inside the root file (correspond to the run number ranges
TIter nextdir(rootFile->GetListOfKeys());
TKey *dir;
@@ -998,16 +931,17 @@ Root::TElectronEfficiencyCorrectionTool::getHistograms() {
obj = dir->ReadObj();
if (obj->IsA()->InheritsFrom("TDirectory")) {
// splits string by delimiter --> e.g RunNumber1_RunNumber2
- std::unique_ptr<TObjArray> dirNameArray(TString(obj->GetName()).Tokenize("_"));
+ //std::unique_ptr<TObjArray> dirNameArray(TString(obj->GetName()).Tokenize("_"));
+ TObjArray dirNameArray = *(TString(obj->GetName()).Tokenize("_"));
//// returns index of last string --> if one, the directory name does not contain any run numbers
- int lastIdx = dirNameArray->GetLast();
+ int lastIdx = dirNameArray.GetLast();
if (lastIdx != 1) {
ATH_MSG_ERROR(" (file: " << __FILE__ << ", line: " << __LINE__ << ") " << "The folder name seems to have the wrong format! Directory name:"
<< obj->GetName());
return 0;
}
rootFile->cd(obj->GetName());
- if (0 == this->setupHistogramsInFolder(dirNameArray.get(), lastIdx)) {
+ if (0 == this->setupHistogramsInFolder(dirNameArray, lastIdx)) {
ATH_MSG_ERROR(
" (file: " << __FILE__ << ", line: " << __LINE__ << ") " << "unable to setup the histograms in directory " << dir->GetName()
<< "in file " << m_corrFileNameList[i]);
@@ -1023,25 +957,23 @@ Root::TElectronEfficiencyCorrectionTool::getHistograms() {
gDirectory = origDir;
} // End: directory loop
} // End: file loop
-
return 1;
}
// =============================================================================
// Get the input histograms from a given folder/run number range
// =============================================================================
-int
-Root::TElectronEfficiencyCorrectionTool::setupHistogramsInFolder(TObjArray *dirNameArray, int lastIdx) {
+int Root::TElectronEfficiencyCorrectionTool::setupHistogramsInFolder(const TObjArray& dirNameArray, int lastIdx) {
ATH_MSG_DEBUG(" (file: " << __FILE__ << ", line: " << __LINE__ << ") " << "entering funtion setupHistogramsInFolder");
m_runNumBegin = -1;
- TString myBegRunNumString = ((TObjString *) dirNameArray->At(lastIdx - 1))->GetString();
+ TString myBegRunNumString = ((TObjString *) dirNameArray.At(lastIdx - 1))->GetString();
if (myBegRunNumString.IsDigit()) {
m_runNumBegin = myBegRunNumString.Atoi();
}
m_runNumEnd = -1;
- TString myEndRunNumString = ((TObjString *) dirNameArray->At(lastIdx))->GetString();
+ TString myEndRunNumString = ((TObjString *) dirNameArray.At(lastIdx))->GetString();
if (myEndRunNumString.IsDigit()) {
m_runNumEnd = myEndRunNumString.Atoi();
}
@@ -1052,26 +984,30 @@ Root::TElectronEfficiencyCorrectionTool::setupHistogramsInFolder(TObjArray *dirN
}
ATH_MSG_DEBUG(" (file: " << __FILE__ << ", line: " << __LINE__ << ") " << m_runNumBegin << " " << m_runNumEnd);
//
- std::unique_ptr<std::map<int, TObjArray *> > objsFull = CxxUtils::make_unique< std::map<int, TObjArray *> >();
- std::unique_ptr<std::map<int, TObjArray *> > objsFast = CxxUtils::make_unique< std::map<int, TObjArray *> >();;
/// setup pairs of obj arrays and keys --> e.g. "sf", new Array to take all SF Histos
+ std::map<int, TObjArray> objsFull;
+ std::map<int, TObjArray > objsFast;
for (unsigned int ikey = 0; ikey < m_keys.size(); ++ikey) {
- objsFull->insert(std::make_pair(m_keys.at(ikey), new TObjArray()));
- objsFast->insert(std::make_pair(m_keys.at(ikey), new TObjArray()));
- }
- objsFull->insert(std::make_pair(mapkey::sys, new TObjArray()));
- objsFast->insert(std::make_pair(mapkey::sys, new TObjArray()));
+ TObjArray dummyFull;
+ objsFull.insert(std::make_pair(m_keys.at(ikey),dummyFull));
+ TObjArray dummyFast;
+ objsFast.insert(std::make_pair(m_keys.at(ikey),dummyFast));
+ }
+ TObjArray dummyFull;
+ objsFull.insert(std::make_pair(mapkey::sys, dummyFull));
+ TObjArray dummyFast;
+ objsFast.insert(std::make_pair(mapkey::sys, dummyFast));
+ //---------------------------------------------------------
//
- std::unique_ptr<std::vector<TObjArray *> > sysObjsFull(new std::vector<TObjArray *>);
- std::unique_ptr<std::vector<TObjArray *> > sysObjsFast(new std::vector<TObjArray *>);
+ std::vector<TObjArray > sysObjsFull;
+ std::vector<TObjArray > sysObjsFast;
//
TIter nextkey(gDirectory->GetListOfKeys());
TKey *key;
TObject *obj(0);
TString tmpName = "";
- TObjArray* tmpArrayFull(0);
- TObjArray* tmpArrayFast(0);
int tmpCounter = 0;
+
//
//Loop of the keys
while ((key = (TKey *) nextkey())) {
@@ -1083,39 +1019,38 @@ Root::TElectronEfficiencyCorrectionTool::setupHistogramsInFolder(TObjArray *dirN
//Add all except the correlated
for (unsigned int ikey = 0; ikey < m_keys.size(); ++ikey) {
if (TString(obj->GetName()).EndsWith("_" + TString(mapkey::keytostring(m_keys.at(ikey))))) {
- objsFull->find(m_keys.at(ikey))->second->Add(obj);
+ objsFull.find(m_keys.at(ikey))->second.Add(obj);
}
}
//
tmpName = TString(obj->GetName());
- //Special treatment is for photons
+ //Special treatment , this is only for photons
if (tmpName.EndsWith("_sys")) {
- objsFull->find(mapkey::sys)->second->Add(obj);
+ objsFull.find(mapkey::sys)->second.Add(obj);
+ TObjArray tmpArrayFull;
+ tmpArrayFull.Add(obj);
+ sysObjsFull.push_back(tmpArrayFull);
tmpCounter++;
- if (tmpArrayFull == 0 && (objsFull->find(mapkey::sys)->second)) {
- tmpArrayFull = new TObjArray();
- tmpArrayFull->Add(obj);
- }else if (tmpArrayFull != 0) {
- tmpArrayFull->Add(obj);
- }
}
+ //
//See if we are dealing with correlated
if (tmpName.Contains("_corr")) {
- //The 1st correlated one corr0
- if (tmpName.EndsWith("corr0")) {
- tmpCounter = 0;
- if (tmpArrayFull != 0) {
- sysObjsFull->push_back(tmpArrayFull);
- }
- tmpArrayFull = new TObjArray();
- }//corr0
- if (tmpArrayFull != 0) {
- tmpArrayFull->Add(obj);
- }else {
- ATH_MSG_ERROR(" (file: " << __FILE__ << ", line: " << __LINE__ << ") " << "Tried to read systematic uncertainties for full sim scale factor,"
- << " but could NOT find a histogram with name ending in \"corr0\". Please check input file! ");
- return 0;
+ if (tmpName.EndsWith("corr0")) {
+ //This is the worse part ...
+ //corr0 triggers a few thing
+ //
+ //1st create a TObjectArray
+ //For high or low Pt (one for each ...)
+ TObjArray tmpArrayFull;
+ //Resgister it to the vector (so the vector has size 1 or 2)
+ sysObjsFull.push_back(tmpArrayFull);
+ //Reset the counter
+ tmpCounter=0;
}
+ //Add to the current last element of the sysObject full
+ //This will be Low Pt once and high Pt the other time
+ sysObjsFull.back().Add(obj);
+ //Now increase the counter
tmpCounter++;
}
if (tmpCounter > m_nSysMax) {
@@ -1125,27 +1060,36 @@ Root::TElectronEfficiencyCorrectionTool::setupHistogramsInFolder(TObjArray *dirN
else if (TString(obj->GetName()).Contains("AtlFast2")) {
for (unsigned int ikey = 0; ikey < m_keys.size(); ++ikey) {
if (TString(obj->GetName()).EndsWith("_" + TString(mapkey::keytostring(m_keys.at(ikey))))) {
- objsFast->find(m_keys.at(ikey))->second->Add(obj);
+ objsFast.find(m_keys.at(ikey))->second.Add(obj);
}
}
//See if we are dealing with correlated
tmpName = TString(obj->GetName());
+ //Special treatment , this is only for photons
+ if (tmpName.EndsWith("_sys")) {
+ objsFast.find(mapkey::sys)->second.Add(obj);
+ TObjArray tmpArrayFast;
+ tmpArrayFast.Add(obj);
+ sysObjsFast.push_back(tmpArrayFast);
+ tmpCounter++;
+ }
+ //
+ //See if we are dealing with correlated
if (tmpName.Contains("_corr")) {
- //The 1st correlated one corr0
- if (tmpName.EndsWith("corr0")) {
- tmpCounter = 0;
- if (tmpArrayFast != 0) {
- sysObjsFast->push_back(tmpArrayFast);
- }
- tmpArrayFast = new TObjArray();
- }//corr0
- if (tmpArrayFast != 0) {
- tmpArrayFast->Add(obj);
- }else {
- ATH_MSG_ERROR(" (file: " << __FILE__ << ", line: " << __LINE__ << ") " << "! Tried to read systematic uncertainties for fast sim scale factor, "
- << "but could NOT find a histogram with name ending in \"corr0\". Please check input file! ");
- return 0;
+ if (tmpName.EndsWith("corr0")) {
+ //This is the worse part ...
+ //corr0 triggers a few thing
+ //1st create a TObjectArray
+ TObjArray tmpArrayFast;
+ //Resgister it
+ sysObjsFast.push_back(tmpArrayFast);
+ //Reset the counter
+ tmpCounter=0;
}
+ //Add to the current last element of the sysObject full
+ //This will be Low Pt once and high Pt the other time
+ sysObjsFast.back().Add(obj);
+ //Now increase the counter
tmpCounter++;
}
if (tmpCounter > m_nSysMax) {
@@ -1159,53 +1103,43 @@ Root::TElectronEfficiencyCorrectionTool::setupHistogramsInFolder(TObjArray *dirN
}
}
}
- if (tmpArrayFull != 0) {
- sysObjsFull->push_back(tmpArrayFull);
- }
- if (tmpArrayFast != 0) {
- sysObjsFast->push_back(tmpArrayFast);
- }
- ATH_MSG_DEBUG(
- " (file: " << __FILE__ << ", line: " << __LINE__ << ") " << "setting up histograms for run ranges " <<
+ ATH_MSG_DEBUG(" (file: " << __FILE__ << ", line: " << __LINE__ << ") " << "setting up histograms for run ranges " <<
m_runNumEnd);
//
- //The setup here copies from the tmp in this function , to the class member variables.
+ //The setup here copies from the temporaties created in this function ,
+ //to the actual class member variables.
for (unsigned int ikey = 0; ikey < m_keys.size(); ++ikey) {
- if (objsFull->find(m_keys.at(ikey))->second->GetEntries() != 0) {
+ if (objsFull.find(m_keys.at(ikey))->second.GetEntries() != 0) {
if (0 ==
- this->setup(objsFull->find(m_keys.at(ikey))->second, m_histList[m_keys.at(ikey)], m_begRunNumberList,
- m_endRunNumberList)) {
- ATH_MSG_ERROR(" (file: " << __FILE__ << ", line: " << __LINE__ << ") " << "! Could NOT setup histogram " << m_keys.at(
- ikey)
- << " for full sim!");
+ setup(objsFull.find(m_keys.at(ikey))->second, m_histList[m_keys.at(ikey)],
+ m_begRunNumberList,m_endRunNumberList)) {
+ ATH_MSG_ERROR(" (file: " << __FILE__ << ", line: " << __LINE__ << ") " << "! Could NOT setup histogram "
+ << m_keys.at(ikey)<< " for full sim!");
return 0;
}
}
- if (objsFast->find(m_keys.at(ikey))->second->GetEntries() != 0) {
+ if (objsFast.find(m_keys.at(ikey))->second.GetEntries() != 0) {
if (0 ==
- this->setup(objsFast->find(m_keys.at(ikey))->second, m_fastHistList[m_keys.at(ikey)],
- m_begRunNumberListFastSim, m_endRunNumberListFastSim)) {
- ATH_MSG_ERROR(" (file: " << __FILE__ << ", line: " << __LINE__ << ") " << "! Could NOT setup histogram " << m_keys.at(
- ikey)
+ setup(objsFast.find(m_keys.at(ikey))->second, m_fastHistList[m_keys.at(ikey)],
+ m_begRunNumberListFastSim, m_endRunNumberListFastSim)) {
+ ATH_MSG_ERROR(" (file: " << __FILE__ << ", line: " << __LINE__ << ") " << "! Could NOT setup histogram " << m_keys.at(ikey)
<< " for fast sim");
return 0;
}
}
}
- for (unsigned int sys = 0; sys < sysObjsFast->size(); sys++) {
- m_fastSysList.resize(sysObjsFast->size());
- if (0 ==
- this->setup(sysObjsFast->at(sys), m_fastSysList[sys], m_begRunNumberListFastSim, m_endRunNumberListFastSim)) {
- ATH_MSG_ERROR(
- " (file: " << __FILE__ << ", line: " << __LINE__ << ") " <<
+ for (unsigned int sys = 0; sys < sysObjsFast.size(); sys++) {
+ m_fastSysList.resize(sysObjsFast.size());
+ if (0 == setup(sysObjsFast.at(sys), m_fastSysList[sys], m_begRunNumberListFastSim, m_endRunNumberListFastSim)) {
+ ATH_MSG_ERROR(" (file: " << __FILE__ << ", line: " << __LINE__ << ") " <<
"! Could NOT setup systematic histograms for fast sim");
return 0;
}
}
- for (unsigned int sys = 0; sys < sysObjsFull->size(); sys++) {
- m_sysList.resize(sysObjsFull->size());
- if (0 == this->setup(sysObjsFull->at(sys), m_sysList[sys], m_begRunNumberList, m_endRunNumberList)) {
+ for (unsigned int sys = 0; sys < sysObjsFull.size(); sys++) {
+ m_sysList.resize(sysObjsFull.size());
+ if (0 == setup(sysObjsFull.at(sys), m_sysList[sys], m_begRunNumberList, m_endRunNumberList)) {
ATH_MSG_ERROR(
" (file: " << __FILE__ << ", line: " << __LINE__ << ") " <<
"! Could NOT setup systematic histograms for fast sim");
@@ -1213,8 +1147,6 @@ Root::TElectronEfficiencyCorrectionTool::setupHistogramsInFolder(TObjArray *dirN
}
}
//
- //
- //
//Toys
ATH_MSG_DEBUG(" (file: " << __FILE__ << ", line: " << __LINE__ << ") " << this->getName()
<< "Checking for (m_doToyMC || m_doCombToyMC)");
@@ -1222,43 +1154,44 @@ Root::TElectronEfficiencyCorrectionTool::setupHistogramsInFolder(TObjArray *dirN
if (m_doToyMC || m_doCombToyMC) {
bool fullToysBooked = kFALSE;
bool fastToysBooked = kFALSE;
+ TObjArray dummy;
if (m_histList[mapkey::sf].size() > 0) {
- if (objsFull->find(mapkey::eig)->second->GetEntries() < 1 || objsFull->find(mapkey::stat)->second->GetEntries() < 1 ||
- objsFull->find(mapkey::uncorr)->second->GetEntries() < 1) {
- if (objsFull->find(mapkey::stat)->second->GetEntries() > 1 || objsFull->find(mapkey::sys)->second->GetEntries() > 1) {
- m_uncorrToyMCSystFull->push_back(buildToyMCTable(objsFull->find(mapkey::sf)->second, 0,
- objsFull->find(mapkey::stat)->second,
- 0, sysObjsFull.get()));
+ if (objsFull.find(mapkey::eig)->second.GetEntries() < 1 || objsFull.find(mapkey::stat)->second.GetEntries() < 1 ||
+ objsFull.find(mapkey::uncorr)->second.GetEntries() < 1) {
+ if (objsFull.find(mapkey::stat)->second.GetEntries() > 1 || objsFull.find(mapkey::sys)->second.GetEntries() > 1) {
+ m_uncorrToyMCSystFull.push_back(buildToyMCTable(objsFull.find(mapkey::sf)->second, dummy,
+ objsFull.find(mapkey::stat)->second,
+ dummy, sysObjsFull));
fullToysBooked = kTRUE;
}else {
ATH_MSG_INFO(
"! Toy MC error propagation booked, but not all needed histograms found in the output (For full sim). Skipping toy creation!");
}
}else {
- m_uncorrToyMCSystFull->push_back(buildToyMCTable(objsFull->find(mapkey::sf)->second, objsFull->find(mapkey::eig)->second,
- objsFull->find(mapkey::stat)->second,
- objsFull->find(mapkey::uncorr)->second, sysObjsFull.get()));
+ m_uncorrToyMCSystFull.push_back(buildToyMCTable(objsFull.find(mapkey::sf)->second, objsFull.find(mapkey::eig)->second,
+ objsFull.find(mapkey::stat)->second,
+ objsFull.find(mapkey::uncorr)->second, sysObjsFull));
fullToysBooked = kTRUE;
}
}
///// AF2
if (m_fastHistList[mapkey::sf].size() > 0) {
- if (objsFast->find(mapkey::eig)->second->GetEntries() < 1 || objsFast->find(mapkey::stat)->second->GetEntries() < 1 ||
- objsFast->find(mapkey::uncorr)->second->GetEntries() < 1) {
- if (objsFast->find(mapkey::stat)->second->GetEntries() > 1 || objsFast->find(mapkey::sys)->second->GetEntries() > 1) {
- m_uncorrToyMCSystFast->push_back(buildToyMCTable(objsFast->find(mapkey::sf)->second, 0,
- objsFast->find(mapkey::stat)->second,
- 0, sysObjsFast.get()));
+ if (objsFast.find(mapkey::eig)->second.GetEntries() < 1 || objsFast.find(mapkey::stat)->second.GetEntries() < 1 ||
+ objsFast.find(mapkey::uncorr)->second.GetEntries() < 1) {
+ if (objsFast.find(mapkey::stat)->second.GetEntries() > 1 || objsFast.find(mapkey::sys)->second.GetEntries() > 1) {
+ m_uncorrToyMCSystFast.push_back(buildToyMCTable(objsFast.find(mapkey::sf)->second, dummy,
+ objsFast.find(mapkey::stat)->second,
+ dummy, sysObjsFast));
fastToysBooked = kTRUE;
}else {
ATH_MSG_INFO(
"! Toy MC error propagation booked, but not all needed histograms found in the output (For fast sim). Skipping toy creation!");
}
}else {
- m_uncorrToyMCSystFast->push_back(buildToyMCTable(objsFast->find(mapkey::sf)->second, objsFast->find(mapkey::eig)->second,
- objsFast->find(mapkey::stat)->second,
- objsFast->find(mapkey::uncorr)->second, sysObjsFast.get()));
+ m_uncorrToyMCSystFast.push_back(buildToyMCTable(objsFast.find(mapkey::sf)->second, objsFast.find(mapkey::eig)->second,
+ objsFast.find(mapkey::stat)->second,
+ objsFast.find(mapkey::uncorr)->second, sysObjsFast));
fastToysBooked = kTRUE;
}
}
@@ -1279,18 +1212,18 @@ Root::TElectronEfficiencyCorrectionTool::setupHistogramsInFolder(TObjArray *dirN
// Fill and interpret the setup, depending on which histograms are found in the input file(s)
// =============================================================================
int
-Root::TElectronEfficiencyCorrectionTool::setup(TObjArray *hists,
- std::vector< TObjArray * > &histList,
+Root::TElectronEfficiencyCorrectionTool::setup(const TObjArray& hists,
+ std::vector< TObjArray> &histList,
std::vector< unsigned int > &beginRunNumberList,
std::vector< unsigned int > &endRunNumberList) {
- if (!hists) {
+ if (hists.GetEntriesFast()==0) {
ATH_MSG_ERROR(
"(file: " << __FILE__ << ", line: " << __LINE__ << ") " << "! Could NOT find histogram with name *_sf in folder");
return 0;
}
TH2 *tmpHist(0);
- for (int i = 0; i < hists->GetEntries(); ++i) {
- tmpHist = (TH2 *) hists->At(i);
+ for (int i = 0; i < hists.GetEntries(); ++i) {
+ tmpHist = (TH2 *) hists.At(i);
tmpHist->SetDirectory(0);
}
// Now, we have all the needed info. Fill the vectors accordingly
@@ -1309,10 +1242,8 @@ Root::TElectronEfficiencyCorrectionTool::setup(TObjArray *hists,
}else {
endRunNumberList.push_back(m_runNumEnd);
}
-
return 1;
}
-
// =============================================================================
// print a message that the default scale factor is returned
// =============================================================================
diff --git a/PhysicsAnalysis/ElectronPhotonID/ElectronEfficiencyCorrection/cmt/Makefile.RootCore b/PhysicsAnalysis/ElectronPhotonID/ElectronEfficiencyCorrection/cmt/Makefile.RootCore
index 25ebe814d0db79b4c88130a26d531592bbae11b4..f0313ffdfa888fb18d6c33ca65c5dc73e5db245f 100644
--- a/PhysicsAnalysis/ElectronPhotonID/ElectronEfficiencyCorrection/cmt/Makefile.RootCore
+++ b/PhysicsAnalysis/ElectronPhotonID/ElectronEfficiencyCorrection/cmt/Makefile.RootCore
@@ -3,7 +3,7 @@ PACKAGE = ElectronEfficiencyCorrection
PACKAGE_PRELOAD = Tree TreePlayer XMLParser XMLIO
PACKAGE_CXXFLAGS = -I/usr/include/libxml2
PACKAGE_BINFLAGS = -lxml2
-PACKAGE_DEP = AsgTools Asg_Boost xAODEgamma xAODEventInfo PATCore PATInterfaces PathResolver CxxUtils
+PACKAGE_DEP = ElectronPhotonSelectorTools AsgAnalysisInterfaces AsgTools Asg_Boost xAODEgamma xAODEventInfo PATCore PATInterfaces PathResolver CxxUtils xAODMetaData
PACKAGE_PEDANTIC = 1
PACKAGE_REFLEX = 1
include $(ROOTCOREDIR)/Makefile-common
diff --git a/PhysicsAnalysis/ElectronPhotonID/ElectronEfficiencyCorrection/cmt/requirements b/PhysicsAnalysis/ElectronPhotonID/ElectronEfficiencyCorrection/cmt/requirements
index 6bcc58f73d30517b8ee2e851544945dd7244e6a9..dee22d6c4815d89ba06126c10ea83630fa3e3174 100644
--- a/PhysicsAnalysis/ElectronPhotonID/ElectronEfficiencyCorrection/cmt/requirements
+++ b/PhysicsAnalysis/ElectronPhotonID/ElectronEfficiencyCorrection/cmt/requirements
@@ -11,25 +11,30 @@ use AtlasROOT AtlasROOT-* External
use PATCore PATCore-* PhysicsAnalysis/AnalysisCommon
use PATInterfaces PATInterfaces-* PhysicsAnalysis/AnalysisCommon
use AsgTools AsgTools-* Control/AthToolSupport
-use AthContainers AthContainers-* Control
+use AthContainers AthContainers-* Control
use xAODEgamma xAODEgamma-* Event/xAOD
use AtlasBoost AtlasBoost-* External
## put here your package dependencies...
private
+use AsgAnalysisInterfaces AsgAnalysisInterfaces-* PhysicsAnalysis/Interfaces
+use AthAnalysisBaseComps AthAnalysisBaseComps-* Control
+use AsgTools AsgTools-* Control/AthToolSupport
use GaudiInterface GaudiInterface-* External
use xAODEventInfo xAODEventInfo-* Event/xAOD
-use CxxUtils CxxUtils-* Control
+use CxxUtils CxxUtils-* Control
use PathResolver PathResolver-* Tools
use AthenaBaseComps AthenaBaseComps-* Control
use xAODCore xAODCore-* Event/xAOD
#next two are needed to get checkreq to shhh about the util stuff... this is annoying, but that's life! (will@cern)
use xAODCaloEvent xAODCaloEvent-* Event/xAOD
-use xAODTracking xAODTracking-* Event/xAOD
+use xAODTracking xAODTracking-* Event/xAOD
+use xAODMetaData xAODMetaData-* Event/xAOD
+use ElectronPhotonSelectorTools ElectronPhotonSelectorTools-* PhysicsAnalysis/ElectronPhotonID
end_private
##
-## macros
+## macros
apply_tag ROOTBasicLibs
apply_tag ROOTMathLibs
@@ -49,9 +54,8 @@ apply_pattern declare_python_modules files="*.py"
apply_pattern declare_joboptions files="*.py"
## install the ROOT files (in the data/ directory)
-apply_pattern declare_calib files="../data/*.root ../data/*.txt"
+apply_pattern declare_calib files="../data/*.root ../data/*.txt"
## For reflex dictionary generation
use AtlasReflex AtlasReflex-* External -no_auto_imports
apply_pattern lcgdict dict=ElectronEfficiencyCorrection selectionfile=selection.xml headerfiles="..\/ElectronEfficiencyCorrection/ElectronEfficiencyCorrectionDict.h"
-
diff --git a/PhysicsAnalysis/ElectronPhotonID/ElectronEfficiencyCorrection/src/EleChargeAlg.cxx b/PhysicsAnalysis/ElectronPhotonID/ElectronEfficiencyCorrection/src/EleChargeAlg.cxx
new file mode 100644
index 0000000000000000000000000000000000000000..f1143ee0b7247c38faa2fb06a02410b80452aa99
--- /dev/null
+++ b/PhysicsAnalysis/ElectronPhotonID/ElectronEfficiencyCorrection/src/EleChargeAlg.cxx
@@ -0,0 +1,138 @@
+/*
+ Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration
+*/
+
+// ElectronChargeCorrection includes
+#include "EleChargeAlg.h"
+#include "xAODEventInfo/EventInfo.h"
+#include "xAODEgamma/ElectronContainer.h"
+// #include "xAODTruth/TruthEventContainer.h"
+// #include "xAODTruth/TruthVertex.h"
+// #include "xAODEgamma/EgammaTruthxAODHelpers.h"
+#include "GaudiKernel/ITHistSvc.h"
+#include "PATInterfaces/ISystematicsTool.h"
+
+// #include "../ElectronChargeEfficiencyCorrectionTool/IElectronChargeEfficiencyCorrectionTool.h"
+#include "AsgAnalysisInterfaces/IEfficiencyScaleFactorTool.h"
+
+
+
+EleChargeAlg::EleChargeAlg( const std::string& name, ISvcLocator* pSvcLocator ) :
+ AthAlgorithm( name, pSvcLocator ),
+ m_syst(),
+ m_eccTool("CP::ElectronChargeEfficiencyCorrectionTool/ElectronChargeEfficiencyCorrectionTool", this),
+ m_eleContName("Electrons")
+{
+ declareProperty("ElectronChargeEfficiencyCorrectionTool", m_eccTool, "The private ElectronChargeEfficiencyCorrectionTool" );
+ declareProperty("ElectronContainerName", m_eleContName, "The name of the input electron container");
+}
+
+
+EleChargeAlg::~EleChargeAlg() {}
+
+
+StatusCode EleChargeAlg::initialize() {
+
+ ATH_MSG_INFO ("Initializing " << name() << "...");
+ ATH_CHECK( m_eccTool.retrieve() );
+ m_syst = m_eccTool->affectingSystematics();
+ for( const auto& variation : m_syst ) {
+ ATH_MSG_INFO(" Affecting systematics: " << variation.name());
+ }
+ return StatusCode::SUCCESS;
+}
+
+StatusCode EleChargeAlg::finalize() {
+ ATH_MSG_INFO ("Finalizing " << name() << "...");
+ ATH_CHECK( m_eccTool.release() );
+
+ return StatusCode::SUCCESS;
+}
+
+StatusCode EleChargeAlg::execute() {
+
+ //----------------------------
+ // Event information
+ //---------------------------
+
+ const xAOD::EventInfo* eventInfo = 0; //NOTE: Everything that comes from the storegate direct from the input files is const!
+
+ // ask the event store to retrieve the xAOD EventInfo container
+
+ CHECK( evtStore()->retrieve( eventInfo) );
+ // if there is only one container of that type in the xAOD (as with the EventInfo container), you do not need to pass
+ // the key name, the default will be taken as the only key name in the xAOD
+
+ // check if data or MC
+ if(!eventInfo->eventType(xAOD::EventInfo::IS_SIMULATION ) ){
+ ATH_MSG_DEBUG( "DATA. Will stop processing this algorithm for the current event.");
+ return StatusCode::SUCCESS; // stop this algorithms execute() for this event, here only interested in MC
+ }
+
+
+ //---------
+ // electrons
+ //---------
+
+ const xAOD::ElectronContainer* electrons = 0;
+ CHECK( evtStore()->retrieve( electrons, m_eleContName.value() ) );
+
+ ATH_MSG_VERBOSE ("Executing " << name() << "... in event with: " << electrons->size() << " electrons");
+
+ // Loop over all electrons in this container
+ for ( const auto* electron : *electrons ) {
+
+ if ( electron->pt() < 25000 ) continue;
+ if ( std::abs(electron->eta()) > 2.5 ) continue;
+
+ ATH_MSG_DEBUG("-------------------------------------------------------------------" );
+ ATH_MSG_DEBUG(" ... electron with pt = " << electron->pt() << " , eta: " << fabs(electron->eta()) );
+
+ double sf=-999;
+ //float rate=-999;
+ CP::SystematicSet syst_set;
+
+ CP::SystematicCode tmpSysResult = m_eccTool->applySystematicVariation(syst_set);
+ ATH_MSG_DEBUG("applySystematicVariation returned CP::SystematicCode = " << tmpSysResult);
+ if ( tmpSysResult != CP::SystematicCode::Ok ) {
+ ATH_MSG_ERROR("ElectronChargeEfficiencyCorrectionTool did NOT report a CP::SystematicCode::Ok when calling applySystematicVariation");
+ return StatusCode::FAILURE;
+ }
+
+ CP::CorrectionCode tmpResult = m_eccTool->getEfficiencyScaleFactor(*electron,sf);
+ ATH_MSG_DEBUG("Nominal value SF = " << sf << " and CP::CorrectionCode = " << tmpResult);
+ if ( tmpResult == CP::CorrectionCode::Error ) {
+ ATH_MSG_ERROR("ElectronChargeEfficiencyCorrectionTool reported a CP::CorrectionCode::Error");
+ return StatusCode::FAILURE;
+ }
+
+ for( const auto& variation : m_syst ) {
+ double tmpSF=0;
+ ///float tmprate=0;
+ CP::SystematicSet syst_set1;
+ syst_set1.insert( variation );
+
+ tmpSysResult = m_eccTool->applySystematicVariation(syst_set1);
+ ATH_MSG_DEBUG("applySystematicVariation (second time) returned CP::SystematicCode = " << tmpSysResult);
+ if ( tmpSysResult != CP::SystematicCode::Ok ) {
+ ATH_MSG_ERROR("ElectronChargeEfficiencyCorrectionTool did NOT report a CP::SystematicCode::Ok when calling applySystematicVariation");
+ return StatusCode::FAILURE;
+ }
+ tmpResult = m_eccTool->getEfficiencyScaleFactor(*electron, tmpSF);
+ ATH_MSG_DEBUG("getEfficiencyScaleFactor returned CP::CorrectionCode = " << tmpResult);
+ if ( tmpResult == CP::CorrectionCode::Error ) {
+ ATH_MSG_ERROR("ElectronChargeEfficiencyCorrectionTool reported a CP::CorrectionCode::Error when calling getEfficiencyScaleFactor");
+ return StatusCode::FAILURE;
+ }
+ // //m_eccTool->getDataEfficiency(*electron, tmprate);
+ ATH_MSG_DEBUG(" Systematic "<< variation.name() << ":" << " SF-> " << tmpSF);
+ // //ATH_MSG_DEBUG("Systematic "<< variation.name() << ":" << " data rate-> " << tmprate);
+ }
+
+
+ } // end for loop over electrons
+
+
+ return StatusCode::SUCCESS;
+
+} //-------end of execute
diff --git a/PhysicsAnalysis/ElectronPhotonID/ElectronEfficiencyCorrection/src/EleChargeAlg.h b/PhysicsAnalysis/ElectronPhotonID/ElectronEfficiencyCorrection/src/EleChargeAlg.h
new file mode 100644
index 0000000000000000000000000000000000000000..7fe1c95c3297d747189c7857ef12ab8e6c51a16b
--- /dev/null
+++ b/PhysicsAnalysis/ElectronPhotonID/ElectronEfficiencyCorrection/src/EleChargeAlg.h
@@ -0,0 +1,47 @@
+/*
+ Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration
+*/
+
+#ifndef ELECTRONCHARGECORRECTION_ELECHARGEALG_H
+#define ELECTRONCHARGECORRECTION_ELECHARGEALG_H 1
+
+#include "AthenaBaseComps/AthAlgorithm.h"
+#include "GaudiKernel/ToolHandle.h"
+// #include "xAODTruth/TruthParticle.h"
+
+#include "TTree.h"
+#include "TLorentzVector.h"
+
+#include <vector>
+
+// #include "MCTruthClassifier/IMCTruthClassifier.h"
+#include "PATInterfaces/SystematicSet.h"
+
+namespace CP{
+ // class IElectronChargeEfficiencyCorrectionTool;
+ class IEfficiencyScaleFactorTool;
+}
+
+class EleChargeAlg: public ::AthAlgorithm {
+ public:
+ EleChargeAlg( const std::string& name, ISvcLocator* pSvcLocator);
+
+ virtual ~EleChargeAlg();
+
+ virtual StatusCode initialize();
+ virtual StatusCode execute();
+ virtual StatusCode finalize();
+
+
+ private:
+ /// The systematics set that will contain all the systematic variations
+ CP::SystematicSet m_syst;
+
+ /// The handle to the charge-misid scale-factor tool
+ ToolHandle<CP::IEfficiencyScaleFactorTool> m_eccTool;
+
+ /// The name of the input electron container
+ StringProperty m_eleContName;
+};
+
+#endif //> !ELECTRONCHARGECORRECTION_ELECHARGEALG_H
diff --git a/PhysicsAnalysis/ElectronPhotonID/ElectronEfficiencyCorrection/util/testEgChargeCorr.cxx b/PhysicsAnalysis/ElectronPhotonID/ElectronEfficiencyCorrection/util/testEgChargeCorr.cxx
new file mode 100644
index 0000000000000000000000000000000000000000..96859ba18460c85fcff3ee7ebf329c16379956f3
--- /dev/null
+++ b/PhysicsAnalysis/ElectronPhotonID/ElectronEfficiencyCorrection/util/testEgChargeCorr.cxx
@@ -0,0 +1,218 @@
+/*
+ Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration
+*/
+
+// System include(s):
+#include <memory>
+#include <cstdlib>
+#include <string>
+
+// ROOT include(s):
+#include <TFile.h>
+#include <TError.h>
+#include <TString.h>
+
+// Infrastructure include(s):
+#ifdef ROOTCORE
+# include "xAODRootAccess/Init.h"
+# include "xAODRootAccess/TEvent.h"
+# include "xAODRootAccess/TStore.h"
+#endif // ROOTCORE
+
+// EDM include(s):
+#include "xAODEventInfo/EventInfo.h"
+#include "xAODEgamma/ElectronContainer.h"
+#include "xAODEgamma/Egamma.h"
+#include "../Root/ElectronChargeEfficiencyCorrectionTool.h"
+#include "xAODCore/ShallowCopy.h"
+
+//Asg includes
+#include "AsgTools/AsgMessaging.h"
+#include "PATInterfaces/SystematicsUtil.h"
+
+#include <ElectronPhotonSelectorTools/AsgElectronLikelihoodTool.h>
+
+
+#define CHECK( ARG ) \
+do { \
+ const bool result = ARG; \
+ if( ! result ) { \
+ ::Error( APP_NAME, "Failed to execute: \"%s\"", \
+#ARG ); \
+ return 1; \
+ } \
+ } while( false )
+
+
+int main( int argc, char* argv[] ) {
+
+ // The application's name:
+ const char* APP_NAME = argv[ 0 ];
+
+ // Check if we received a file name:
+ if( argc < 2 ) {
+ Error( APP_NAME, "No file name received!" );
+ Error( APP_NAME, " Usage: %s [xAOD file name]", APP_NAME );
+ return 1;
+ }
+
+ double SF_chargeID=0;
+ // double SF_chargeMisID=0;
+ int SF_nevents=0;
+
+ double n_chargeID=0;
+ double n_chargeMisID=0;
+
+ // Initialise the application:
+ CHECK( xAOD::Init( APP_NAME ) );
+
+ // Set message level
+ MSG::Level mylevel=MSG::VERBOSE;
+ //this.msg().setLevel(mylevel);
+
+ // Open the input file:
+ const TString fileName = argv[ 1 ];
+ Info( APP_NAME, "Opening file: %s", fileName.Data() );
+ std::auto_ptr< TFile > ifile( TFile::Open( fileName, "READ" ) );
+ CHECK( ifile.get() );
+
+ // Create a TEvent object:
+ //xAOD::TEvent event( xAOD::TEvent::kBranchAccess );
+ xAOD::TEvent event( xAOD::TEvent::kClassAccess );
+ CHECK( event.readFrom( ifile.get() ) );
+ Info( APP_NAME, "Number of events in the file: %i",
+ static_cast< int >( event.getEntries() ) );
+
+ std::cout << "=="<<std::endl;
+
+
+ // Decide how many events to run over:
+ Long64_t entries = event.getEntries();
+ if( argc > 2 ) {
+ const Long64_t e = atoll( argv[ 2 ] );
+ if( e < entries ) {
+ entries = e;
+ }
+ }
+
+ //Likelihood
+ CP::ElectronChargeEfficiencyCorrectionTool myEgCorrections ("myEgCorrections");
+ CHECK( myEgCorrections.setProperty("OutputLevel", mylevel) );
+
+ std::string inputFile = "data/ChMisIDSF_TightLL_FixedCutTight.root" ;
+ CHECK( myEgCorrections.setProperty("CorrectionFileName",inputFile) );
+ CHECK( myEgCorrections.setProperty("ForceDataType",1) );
+
+ CHECK( myEgCorrections.setProperty("DefaultRandomRunNumber", (unsigned int)100000 ) );
+ CHECK( myEgCorrections.setProperty("UseRandomRunNumber",false) );
+ myEgCorrections.initialize();
+
+ AsgElectronLikelihoodTool * m_LHToolTight = new AsgElectronLikelihoodTool("m_LHToolTight");
+ m_LHToolTight->setProperty("primaryVertexContainer","PrimaryVertices").ignore();
+ m_LHToolTight->setProperty("ConfigFile","ElectronPhotonSelectorTools/offline/mc15_20160512/ElectronLikelihoodTightOfflineConfig2016_Smooth.conf").ignore();
+ m_LHToolTight->initialize();
+
+
+ // Get a list of systematics
+ CP::SystematicSet recSysts = myEgCorrections.recommendedSystematics();
+ // Convert into a simple list
+ std::vector<CP::SystematicSet> sysList = CP::make_systematics_vector(recSysts);
+ std::cout << "=="<<std::endl;
+ // Loop over the events:
+ // entries = 1;
+
+ Long64_t entry = 0;
+ for( ; entry < entries; ++entry ) {
+
+ // Tell the object which entry to look at:
+ event.getEntry( entry );
+
+ if (entry<10) std::cout << "=================NEXT EVENT==========================" << std::endl;
+ Info (APP_NAME,"Electron" );
+
+
+ const xAOD::EventInfo* event_info = 0;
+ CHECK( event.retrieve( event_info, "EventInfo" ) );
+
+ const xAOD::ElectronContainer* electrons = 0;
+ CHECK( event.retrieve(electrons, "Electrons") );
+
+ // Loop over systematics
+ // for(const auto& sys : sysList){
+
+ // if (entry<10) Info(APP_NAME, "Processing syst: %s", sys.name().c_str());
+
+ // Configure the tool for this systematic
+ // CHECK( myEgCorrections.applySystematicVariation(sys) );
+ // Info(APP_NAME, "Applied syst: %s",
+ // myEgCorrections.affectingSystematics().name().c_str());
+
+ // Create shallow copy for this systematic
+ std::pair< xAOD::ElectronContainer*, xAOD::ShallowAuxContainer* > electrons_shallowCopy =
+ xAOD::shallowCopyContainer( *electrons );
+
+ //Iterate over the shallow copy
+ xAOD::ElectronContainer* elsCorr = electrons_shallowCopy.first;
+ xAOD::ElectronContainer::iterator el_it = elsCorr->begin();
+ xAOD::ElectronContainer::iterator el_it_last = elsCorr->end();
+
+ //double SF_chargeID=0;
+ //double SF_chargeMisID=0;
+ // int SF_nevents=0;
+
+ unsigned int i = 0;
+ double SF = 0;
+ for (; el_it != el_it_last; ++el_it, ++i) {
+
+ xAOD::Electron* el = *el_it;
+ if (el->pt() < 25000) continue;//skip electrons outside of recommendations
+ if (fabs(el->eta()) > 2.4) continue;//skip electrons outside of recommendations
+ if(!m_LHToolTight->accept(el)) continue;
+
+ SF_nevents++;
+
+ // std::cout << "Electron " << i << std::endl;
+ std::cout << "xAOD/raw pt = " << el->pt() << ", eta: "
+ <<abs(el->eta()) <<" , " << el->caloCluster()->etaBE(2) << std::endl;
+
+ Info (APP_NAME,"Electron #%d", i);
+
+ if(myEgCorrections.getEfficiencyScaleFactor(*el,SF) != CP::CorrectionCode::Ok){
+ Error( APP_NAME, "Problem in getEfficiencyScaleFactor");
+ return EXIT_FAILURE;
+ }
+
+ if(myEgCorrections.applyEfficiencyScaleFactor(*el) != CP::CorrectionCode::Ok){
+ Error( APP_NAME, "Problem in applyEfficiencyScaleFactor");
+ return EXIT_FAILURE;
+ }
+
+ Info( APP_NAME, "===>>> Resulting SF (from get function) %f, (from apply function) %f",
+ SF, el->auxdata< float >("SF"));
+
+ SF_chargeID=SF_chargeID+SF;
+
+ int truthcharge = (-1)*el->auxdata<int>("firstEgMotherPdgId");
+ if (el->auxdata<int>("truthType")) { std::cout << el->charge() << " " << el->auxdata<int>("firstEgMotherPdgId") << std::endl; }
+ if ( el->charge() * truthcharge < 0 ) {
+ Info( APP_NAME, "===>>> MISID %f ",SF);
+ n_chargeMisID++;
+ }
+ else n_chargeID++;
+ }
+
+ // }
+
+ }
+ Info( APP_NAME, "===>>> done processing event #%lld ",entry);
+ Info( APP_NAME, "===>>> processed #%d electrons",SF_nevents);
+ Info( APP_NAME, "===>>> compared to #%f (from Charge MisId SF)",SF_chargeID);
+ Info( APP_NAME, "===>>> compared to #%f and #%f ",n_chargeID,n_chargeMisID);
+
+
+
+ CHECK( myEgCorrections.finalize() );
+
+ // Return gracefully:
+ return 0;
+}
diff --git a/PhysicsAnalysis/ElectronPhotonID/ElectronEfficiencyCorrection/util/testEgEfficiencyCorr.cxx b/PhysicsAnalysis/ElectronPhotonID/ElectronEfficiencyCorrection/util/testEgEfficiencyCorr.cxx
index 8b4c28402854711ce5fa798c7e3f6a54825ef964..0273f4c1af79f0e33aa19ceb275035847a56b3d9 100644
--- a/PhysicsAnalysis/ElectronPhotonID/ElectronEfficiencyCorrection/util/testEgEfficiencyCorr.cxx
+++ b/PhysicsAnalysis/ElectronPhotonID/ElectronEfficiencyCorrection/util/testEgEfficiencyCorr.cxx
@@ -83,10 +83,13 @@ int main( int argc, char* argv[] ) {
//Likelihood
AsgElectronEfficiencyCorrectionTool myEgCorrections ("myEgCorrections");
- std::vector<std::string> inputFiles{"ElectronEfficiencyCorrection/efficiencySF.offline.Loose.2012.8TeV.rel17p2.v07.root"} ;
+ //std::vector<std::string> inputFiles{"ElectronEfficiencyCorrection/efficiencySF.offline.Loose.2012.8TeV.rel17p2.v07.root"} ;
+ std::vector<std::string> inputFiles{"ElectronEfficiencyCorrection/efficiencySF.offline.MediumLLH_d0z0_v11.root"} ;
CHECK( myEgCorrections.setProperty("CorrectionFileNameList",inputFiles) );
- CHECK( myEgCorrections.setProperty("ForceDataType",1) );
- CHECK( myEgCorrections.setProperty("CorrelationModel", "MCTOYS100" ));
+
+ CHECK( myEgCorrections.setProperty("ForceDataType",3) );
+ CHECK( myEgCorrections.setProperty("CorrelationModel", "SIMPLIFIED" ));
+ CHECK( myEgCorrections.setProperty("UseRandomRunNumber", false ));
CHECK( myEgCorrections.initialize() );
@@ -110,7 +113,8 @@ int main( int argc, char* argv[] ) {
CHECK( event.retrieve( event_info, "EventInfo" ) );
const xAOD::ElectronContainer* electrons = 0;
- CHECK( event.retrieve(electrons, "ElectronCollection") );
+ //CHECK( event.retrieve(electrons, "ElectronCollection") ); //For DAOD
+ CHECK( event.retrieve(electrons, "Electrons") );
// Loop over systematics
for(const auto& sys : sysList){