Helper Functions¶
-
namespace
HelperFunctions
¶ Functions
-
MsgStream &
msg
(MSG::Level lvl = MSG::INFO)¶ Static object that provides athena-based message logging functionality
-
bool
passPrimaryVertexSelection
(const xAOD::VertexContainer *vertexContainer, int Ntracks = 2)¶
-
int
countPrimaryVertices
(const xAOD::VertexContainer *vertexContainer, int Ntracks = 2)¶
-
const xAOD::Vertex *
getPrimaryVertex
(const xAOD::VertexContainer *vertexContainer, MsgStream &msg)¶
-
const xAOD::Vertex *
getPrimaryVertex
(const xAOD::VertexContainer *vertexContainer)¶
-
float
getPrimaryVertexZ
(const xAOD::Vertex *pvx)¶
-
int
getPrimaryVertexLocation
(const xAOD::VertexContainer *vertexContainer, MsgStream &msg)¶
-
int
getPrimaryVertexLocation
(const xAOD::VertexContainer *vertexContainer)¶
-
bool
applyPrimaryVertexSelection
(const xAOD::JetContainer *jets, const xAOD::VertexContainer *vertices)¶
-
float
GetBTagMV2c20_Cut
(int efficiency)¶
-
std::string
GetBTagMV2c20_CutStr
(int efficiency)¶
-
std::string
replaceString
(std::string subjet, const std::string &search, const std::string &replace)¶
-
std::vector<TString>
SplitString
(TString &orig, const char separator)¶
-
float
dPhi
(float phi1, float phi2)¶
-
bool
has_exact
(const std::string input, const std::string flag)¶
-
std::size_t
string_pos
(const std::string &haystack, const std::string &needle, unsigned int N)¶ Function which returns the position of the n-th occurence of a character in a string searching backwards. Returns -1 if no occurencies are found.
Source: http://stackoverflow.com/questions/18972258/index-of-nth-occurrence-of-the-string
-
std::string
parse_wp
(const std::string &type, const std::string &config_name, MsgStream &msg)¶ Function which returns the WP for ISO/ID from a config file. Returns empty string if no WP is found.
-
std::string
parse_wp
(const std::string &type, const std::string &config_name)¶
-
StatusCode
isAvailableMetaData
(TTree *metaData)¶
-
bool
isFilePrimaryxAOD
(TFile *inputFile)¶
-
std::vector<TLorentzVector>
jetReclustering
(const xAOD::JetContainer *jets, double radius = 1.0, double fcut = 0.05, fastjet::JetAlgorithm rc_alg = fastjet::antikt_algorithm)¶
-
std::vector<TLorentzVector>
jetTrimming
(const xAOD::JetContainer *jets, double radius = 0.3, double fcut = 0.05, fastjet::JetAlgorithm s_alg = fastjet::kt_algorithm)¶
-
TLorentzVector
jetTrimming
(const xAOD::Jet *jet, double radius = 0.3, double fcut = 0.05, fastjet::JetAlgorithm s_alg = fastjet::kt_algorithm)¶
-
bool
sort_pt
(const xAOD::IParticle *partA, const xAOD::IParticle *partB)¶
-
std::vector<CP::SystematicSet>
getListofSystematics
(const CP::SystematicSet inSysts, std::string systNames, float systVal, MsgStream &msg)¶ Get a list of systematics.
- Parameters
inSysts
: systematics set retrieved from the toolsystNames
: comma separated list of wanted systematics names, use “Nominal” for nominal and “All” for all systematicssystVal
: continuous systematics sigma valuemsg
: the MsgStream object with appropriate level for debugging
-
void
writeSystematicsListHist
(const std::vector<CP::SystematicSet> &systs, std::string histName, TFile *file)¶
- template <typename T>
-
std::string
type_name
(bool useXAOD = true)¶
- template <typename T1, typename T2>
-
StatusCode
makeSubsetCont
(T1 *&intCont, T2 *&outCont, MsgStream &msg, const std::string &flagSelect = "", HelperClasses::ToolName tool_name = HelperClasses::ToolName::DEFAULT)¶ Function to copy a subset of a generic input xAOD container into a generic output xAOD container.
- Author
- Marco Milesi (marco.milesi@cern.ch) If the optional parameters aren’t specified, the function will just make a full copy of the input container into the output one.
- Parameters
intCont
: input containeroutCont
: output containerflagSelect
: (optional) the name of the decoration for objects passing a certain selection (e.g. “passSel”, “overlaps” ...). When explicitly specified, it must not be empty.tool_name
: (optional) an enum specifying the tool type which is calling this function (definition inHelperClasses::ToolName
)
- template <typename T1, typename T2>
-
StatusCode
makeSubsetCont
(T1 *&intCont, T2 *&outCont, const std::string &flagSelect = "", HelperClasses::ToolName tool_name = HelperClasses::ToolName::DEFAULT)¶
- template <typename T>
-
StatusCode
retrieve
(T *&cont, std::string name, xAOD::TEvent *event, xAOD::TStore *store, MsgStream &msg)¶ Retrieve an arbitrary object from TStore / TEvent.
This tries to make your life simple by providing a one-stop container retrieval shop for all types.
Example Usage:
const xAOD::JetContainer* jets(0); // look for "AntiKt10LCTopoJets" in both TEvent and TStore ANA_CHECK( HelperFunctions::retrieve(jets, "AntiKt10LCTopoJets", m_event, m_store) ); // look for "AntiKt10LCTopoJets" in only TStore ANA_CHECK( HelperFunctions::retrieve(jets, "AntiKt10LCTopoJets", 0, m_store) ); // look for "AntiKt10LCTopoJets" in only TEvent, enable verbose output ANA_CHECK( HelperFunctions::retrieve(jets, "AntiKt10LCTopoJets", m_event, 0, msg()) );
Checking Order:
- start by checking TStore
- check if store contains ‘xAOD::JetContainer’ named ‘name’
- attempt to retrieve from store
- return if failure
- check if store contains ‘xAOD::JetContainer’ named ‘name’
- next check TEvent
- check if event contains ‘xAOD::JetContainer’ named ‘name’
- attempt to retrieve from event
- return if failure
- return FAILURE
- check if event contains ‘xAOD::JetContainer’ named ‘name’
- return SUCCESS (should never reach this last line)
- Parameters
cont
: pass in a pointer to the object to store the retrieved container inname
: the name of the object to look upevent
: the TEvent, usually wk()->xaodEvent(). Set to 0 to not search TEvent.store
: the TStore, usually wk()->xaodStore(). Set to 0 to not search TStore.msg
: the MsgStream object with appropriate level for debugging
- start by checking TStore
- template <typename T>
-
StatusCode
retrieve
(T *&cont, std::string name, xAOD::TEvent *event, xAOD::TStore *store)¶
- template <typename T>
-
StatusCode HelperFunctions::__attribute__((deprecated("retrieve<T>(..., bool) is deprecated. See https://github.com/UCATLAS/xAODAnaHelpers/pull/882")))
- template <typename T>
-
bool
isAvailable
(std::string name, xAOD::TEvent *event, xAOD::TStore *store, MsgStream &msg)¶ Return true if an arbitrary object from TStore / TEvent is available.
This tries to make your life simple by providing a one-stop container check shop for all types
Example Usage:
const xAOD::JetContainer* jets(0); // look for "AntiKt10LCTopoJets" in both TEvent and TStore HelperFunctions::isAvailable<xAOD::JetContainer>("AntiKt10LCTopoJets", m_event, m_store) // look for "AntiKt10LCTopoJets" in only TStore HelperFunctions::isAvailable<xAOD::JetContainer>("AntiKt10LCTopoJets", 0, m_store) // look for "AntiKt10LCTopoJets" in only TEvent, enable verbose output HelperFunctions::isAvailable<xAOD::JetContainer>("AntiKt10LCTopoJets", m_event, 0, MSG::VERBOSE)
- Parameters
name
: the name of the object to look upevent
: the TEvent, usually wk()->xaodEvent(). Set to 0 to not search TEvent.store
: the TStore, usually wk()->xaodStore(). Set to 0 to not search TStore.msg
: the MsgStream object with appropriate level for debugging
- template <typename T>
-
bool
isAvailable
(std::string name, xAOD::TEvent *event, xAOD::TStore *store)¶
- template <class T>
-
const T *
getLink
(const xAOD::IParticle *particle, std::string name)¶ Access to element link to object of type T stored in auxdata.
- template <typename T>
-
T
sort_container_pt
(T *inCont)¶
- template <typename T>
-
const T
sort_container_pt
(const T *inCont)¶
-
bool
found_non_dummy_sys
(std::vector<std::string> *sys_list)¶
- template <typename T1, typename T2, typename T3>
-
StatusCode
makeDeepCopy
(xAOD::TStore *m_store, std::string containerName, const T1 *cont)¶ Make a deep copy of a container and put it in the TStore.
This is a very powerful templating function. The point is to remove the triviality of making deep copies by specifying all that is needed. The best way is to demonstrate via example:
const xAOD::JetContainer* selected_jets(nullptr); ANA_CHECK( m_event->retrieve( selected_jets, "SelectedJets" )); ANA_CHECK( (HelperFunctions::makeDeepCopy<xAOD::JetContainer, xAOD::JetAuxContainer, xAOD::Jet>(m_store, "BaselineJets", selected_jets)));
- Template Parameters
T1
: The type of the container you’re going to deep copy intoT2
: The type of the aux container you’re going to deep copy intoT3
: The type of the object inside the container you’re going to deep copy
- Parameters
m_store
: A pointer to the TStore objectcontainerName
: The name of the container to create as output in the TStorecont
: The container to deep copy, it should be a container of pointers (IParticleContainer or ConstDataVector)
- template <typename T1, typename T2>
-
StatusCode
recordOutput
(xAOD::TEvent *m_event, xAOD::TStore *m_store, std::string containerName)¶ Copy a container from the TStore to be recorded in the TEvent (eg: to an output)
If you have a container in the TStore, this function will record it into the output for you without an issue. As an example:
ANA_CHECK( HelperFunctions::recordOutput<xAOD::JetContainer, xAOD::JetAuxContainer>(m_event, m_store, "BaselineJets"));
where we build off the previous example of making a deep copy (see
HelperFunctions::makeDeepCopy()
).- Template Parameters
T1
: The type of the container you’re going to recordT2
: The type of the aux container you’re going to record
- Parameters
m_event
: A pointer to the TEvent objectm_store
: A pointer to the TStore objectcontainerName
: The name of the container in the TStore to record to TEvent
- template <typename T_BR>
-
void
connectBranch
(std::string name, TTree *tree, const std::string &branch, std::vector<T_BR> **variable)¶
- template <typename T>
-
void
remove_duplicates
(std::vector<T> &vec)¶
Variables
-
StatusCode std::string HelperFunctions::name
-
StatusCode std::string xAOD::TEvent* HelperFunctions::event
-
StatusCode std::string xAOD::TEvent xAOD::TStore* HelperFunctions::store
-
StatusCode std::string xAOD::TEvent xAOD::TStore bool HelperFunctions::debug
= { return retrieve<T>(cont, name, event, store, msg())
-
struct
pt_sort
¶ Public Functions
-
bool
operator()
(const TLorentzVector &lhs, const TLorentzVector &rhs)¶
-
bool
operator()
(const TLorentzVector *lhs, const TLorentzVector *rhs)¶
-
bool
operator()
(const xAOD::IParticle &lhs, const xAOD::IParticle &rhs)¶
-
bool
operator()
(const xAOD::IParticle *lhs, const xAOD::IParticle *rhs)¶
-
bool
-
MsgStream &