From dcf27b2b17483bcd8c3dd9a77922a8d59935cddc Mon Sep 17 00:00:00 2001
From: Graham Richard Lee <graham.richard.lee@cern.ch>
Date: Wed, 5 Dec 2018 13:53:05 +0000
Subject: [PATCH] InDetTruthVertexValidation merge request
Former-commit-id: 59dcc6cd30155480e0f05ad6007d8e3d9efbeb77
---
.../InDetTruthVertexValidation/CMakeLists.txt | 62 +++
.../IInDetVertexTruthMatchTool.h | 32 ++
.../InDetVertexTruthMatchTool.h | 54 +++
.../InDetVertexTruthMatchUtils.h | 53 +++
.../Root/InDetVertexTruthMatchTool.cxx | 422 ++++++++++++++++++
.../Root/InDetVertexTruthMatchUtils.cxx | 128 ++++++
.../InDetTruthVertexValidation/Root/LinkDef.h | 23 +
.../share/VertexTruthMatch_jobOptions.py | 27 ++
.../src/InDetVertexTruthMatchAlgorithm.cxx | 43 ++
.../src/InDetVertexTruthMatchAlgorithm.h | 34 ++
.../InDetTruthVertexValidation_entries.cxx | 13 +
.../InDetTruthVertexValidation_load.cxx | 7 +
.../util/VertexTruthMatchTest.cxx | 145 ++++++
13 files changed, 1043 insertions(+)
create mode 100644 InnerDetector/InDetValidation/InDetTruthVertexValidation/CMakeLists.txt
create mode 100644 InnerDetector/InDetValidation/InDetTruthVertexValidation/InDetTruthVertexValidation/IInDetVertexTruthMatchTool.h
create mode 100644 InnerDetector/InDetValidation/InDetTruthVertexValidation/InDetTruthVertexValidation/InDetVertexTruthMatchTool.h
create mode 100644 InnerDetector/InDetValidation/InDetTruthVertexValidation/InDetTruthVertexValidation/InDetVertexTruthMatchUtils.h
create mode 100644 InnerDetector/InDetValidation/InDetTruthVertexValidation/Root/InDetVertexTruthMatchTool.cxx
create mode 100644 InnerDetector/InDetValidation/InDetTruthVertexValidation/Root/InDetVertexTruthMatchUtils.cxx
create mode 100644 InnerDetector/InDetValidation/InDetTruthVertexValidation/Root/LinkDef.h
create mode 100644 InnerDetector/InDetValidation/InDetTruthVertexValidation/share/VertexTruthMatch_jobOptions.py
create mode 100644 InnerDetector/InDetValidation/InDetTruthVertexValidation/src/InDetVertexTruthMatchAlgorithm.cxx
create mode 100644 InnerDetector/InDetValidation/InDetTruthVertexValidation/src/InDetVertexTruthMatchAlgorithm.h
create mode 100644 InnerDetector/InDetValidation/InDetTruthVertexValidation/src/components/InDetTruthVertexValidation_entries.cxx
create mode 100644 InnerDetector/InDetValidation/InDetTruthVertexValidation/src/components/InDetTruthVertexValidation_load.cxx
create mode 100644 InnerDetector/InDetValidation/InDetTruthVertexValidation/util/VertexTruthMatchTest.cxx
diff --git a/InnerDetector/InDetValidation/InDetTruthVertexValidation/CMakeLists.txt b/InnerDetector/InDetValidation/InDetTruthVertexValidation/CMakeLists.txt
new file mode 100644
index 00000000000..8d9c36b25eb
--- /dev/null
+++ b/InnerDetector/InDetValidation/InDetTruthVertexValidation/CMakeLists.txt
@@ -0,0 +1,62 @@
+################################################################################
+# Package: InDetTruthVertexValidation
+################################################################################
+
+# Declare the package name:
+atlas_subdir( InDetTruthVertexValidation )
+
+# Extra dependencies, based on the build environment:
+set( extra_deps )
+set( extra_libs )
+if( XAOD_STANDALONE )
+ set( extra_deps Control/xAODRootAccess
+ PhysicsAnalysis/D3PDTools/EventLoop )
+ set( extra_libs xAODRootAccess EventLoop )
+else()
+ set( extra_deps PRIVATE Control/AthenaBaseComps PhysicsAnalysis/POOLRootAccess GaudiKernel )
+ set( extra_libs AthAnalysisBaseCompsLib )
+endif()
+
+
+# Declare the package's dependencies:
+atlas_depends_on_subdirs( PUBLIC
+ Control/AthToolSupport/AsgTools
+ Event/xAOD/xAODTracking
+ Event/xAOD/xAODTruth
+ Event/EventPrimitives
+ Event/xAOD/xAODEventInfo
+ ${extra_deps} )
+
+# External dependencies:
+find_package( Eigen )
+find_package( ROOT COMPONENTS Core Geometry Tree MathCore Hist RIO pthread TBB )
+
+# Generate a CINT dictionary source file:
+atlas_add_root_dictionary( InDetTruthVertexValidationLib _cintDictSource
+ ROOT_HEADERS Root/LinkDef.h
+ EXTERNAL_PACKAGES ROOT )
+
+# Component(s) in the package:
+atlas_add_library( InDetTruthVertexValidationLib
+ Root/*.cxx ${_cintDictSource}
+ PUBLIC_HEADERS InDetTruthVertexValidation
+ PRIVATE_INCLUDE_DIRS ${ROOT_INCLUDE_DIRS} ${EIGEN_INCLUDE_DIRS} ${TBB_INCLUDE_DIRS}
+ PRIVATE_LINK_LIBRARIES ${ROOT_LIBRARIES} ${EIGEN_LIBRARIES} EventPrimitives xAODTracking AsgTools xAODEventInfo xAODTruth ${extra_libs})
+
+if( NOT XAOD_STANDALONE )
+ atlas_add_component( InDetTruthVertexValidation
+ src/*.cxx src/*.h
+ src/components/*.cxx
+ INCLUDE_DIRS ${ROOT_INCLUDE_DIRS} ${TBB_INCLUDE_DIRS} ${EIGEN_INCLUDE_DIRS}
+ LINK_LIBRARIES ${ROOT_LIBRARIES} ${EIGEN_LIBRARIES} ${extra_libs} EventPrimitives xAODTracking AthenaBaseComps AsgTools xAODEventInfo xAODTruth InDetTruthVertexValidationLib )
+endif()
+
+# Install files from the package:
+atlas_install_joboptions( share/*.py )
+
+if( XAOD_STANDALONE )
+ atlas_add_executable( VertexTruthMatchTest
+ util/VertexTruthMatchTest.cxx
+ LINK_LIBRARIES InDetTruthVertexValidationLib xAODTracking xAODEventInfo xAODRootAccess ${ROOT_LIBRARIES} ${EIGEN_LIBRARIES} ${extra_libs})
+endif()
+
diff --git a/InnerDetector/InDetValidation/InDetTruthVertexValidation/InDetTruthVertexValidation/IInDetVertexTruthMatchTool.h b/InnerDetector/InDetValidation/InDetTruthVertexValidation/InDetTruthVertexValidation/IInDetVertexTruthMatchTool.h
new file mode 100644
index 00000000000..3302f3c0150
--- /dev/null
+++ b/InnerDetector/InDetValidation/InDetTruthVertexValidation/InDetTruthVertexValidation/IInDetVertexTruthMatchTool.h
@@ -0,0 +1,32 @@
+/*
+ Copyright (C) 2002-2018 CERN for the benefit of the ATLAS collaboration
+*/
+
+#ifndef IInDetVertexTruthMatchTool_h
+#define IInDetVertexTruthMatchTool_h
+
+// Framework include(s):
+#include "AsgTools/IAsgTool.h"
+
+// EDM include(s):
+#include "xAODTracking/VertexContainer.h"
+
+
+/** Class for vertex truth matching.
+ * Match reconstructed vertices to truth level interactions vertices
+ * through the chain: track -> particle -> genEvent -> genVertex
+ * Categorize reconstructed vertices depending on their composition.
+ */
+
+class IInDetVertexTruthMatchTool : public virtual asg::IAsgTool {
+
+ASG_TOOL_INTERFACE( IInDetVertexTruthMatchTool )
+
+public:
+
+//take const collection of vertices, match them, and decorate with matching info
+ virtual StatusCode matchVertices( const xAOD::VertexContainer & vxContainer ) = 0;
+
+};
+
+#endif
diff --git a/InnerDetector/InDetValidation/InDetTruthVertexValidation/InDetTruthVertexValidation/InDetVertexTruthMatchTool.h b/InnerDetector/InDetValidation/InDetTruthVertexValidation/InDetTruthVertexValidation/InDetVertexTruthMatchTool.h
new file mode 100644
index 00000000000..6f11c026103
--- /dev/null
+++ b/InnerDetector/InDetValidation/InDetTruthVertexValidation/InDetTruthVertexValidation/InDetVertexTruthMatchTool.h
@@ -0,0 +1,54 @@
+/*
+ Copyright (C) 2002-2018 CERN for the benefit of the ATLAS collaboration
+*/
+
+#ifndef InDetVertexTruthMatchTool_h
+#define InDetVertexTruthMatchTool_h
+
+// Framework include(s):
+#include "AsgTools/AsgTool.h"
+
+// EDM include(s):
+#include "xAODTracking/VertexContainerFwd.h"
+#include "xAODTruth/TruthParticleFwd.h"
+#include "InDetTruthVertexValidation/IInDetVertexTruthMatchTool.h"
+
+/** Class for vertex truth matching.
+ * Match reconstructed vertices to truth level interactions vertices
+ * through the chain: track -> particle -> genEvent -> genVertex
+ * Categorize reconstructed vertices depending on their composition.
+ */
+
+
+class InDetVertexTruthMatchTool : public virtual IInDetVertexTruthMatchTool,
+ public asg::AsgTool {
+
+ ASG_TOOL_CLASS( InDetVertexTruthMatchTool, IInDetVertexTruthMatchTool )
+
+ public:
+
+ InDetVertexTruthMatchTool( const std::string & name );
+
+ virtual StatusCode initialize() override final;
+
+ //take const collection of vertices, match them, and decorate with matching info
+ virtual StatusCode matchVertices( const xAOD::VertexContainer & vxContainer ) override;
+
+ private:
+
+ const xAOD::TrackParticleContainer* findTrackParticleContainer( const xAOD::VertexContainer& vxContainer );
+
+ //required MC match probability to consider track a good match
+ float m_trkMatchProb;
+ //relative weight threshold to consider vertex matched
+ float m_vxMatchWeight;
+ //pt cut to use on tracks
+ float m_trkPtCut;
+
+ //private methods to check if particles are good to use
+ bool pass( const xAOD::TruthParticle & truthPart ) const;
+ bool pass( const xAOD::TrackParticle & trackPart ) const;
+
+};
+
+#endif
diff --git a/InnerDetector/InDetValidation/InDetTruthVertexValidation/InDetTruthVertexValidation/InDetVertexTruthMatchUtils.h b/InnerDetector/InDetValidation/InDetTruthVertexValidation/InDetTruthVertexValidation/InDetVertexTruthMatchUtils.h
new file mode 100644
index 00000000000..5871724d924
--- /dev/null
+++ b/InnerDetector/InDetValidation/InDetTruthVertexValidation/InDetTruthVertexValidation/InDetVertexTruthMatchUtils.h
@@ -0,0 +1,53 @@
+/*
+ Copyright (C) 2002-2018 CERN for the benefit of the ATLAS collaboration
+*/
+
+#ifndef InDetVertexTruthMatchUtils_h
+#define InDetVertexTruthMatchUtils_h
+
+#include "xAODTracking/VertexContainerFwd.h"
+#include "xAODTruth/TruthEventBaseContainer.h"
+
+namespace InDetVertexTruthMatchUtils {
+
+ //Namespace for useful analysis things on the truth matching decorations applied to the VertexContainer
+ //Can be called by algorithms inside or outside Athena to do their analysis
+
+ //How the matching info is stored; link to the truth vertex, a float with its contribution to the relative track weight, and a float with its contribution to the track sumpt2 of the truth vertex
+ typedef std::tuple<ElementLink<xAOD::TruthEventBaseContainer>, float, float> VertexTruthMatchInfo;
+
+ //type codes for vertex matching on all vertices
+ enum VertexMatchType {
+ MATCHED, // > threshold (default 70%) from one truth interaction
+ MERGED, // not matched
+ SPLIT, // highest weight truth interaction contributes to >1 vtx (vtx with highest fraction of sumpT2 remains matched/merged)
+ FAKE, // highest contribution is fake (if pile-up MC info not present those tracks end up as "fakes")
+ DUMMY, // is the dummy vertex
+ NTYPES // access to number of types
+};
+
+ //Classification for true hard-scatter interaction
+ enum HardScatterType {
+ CLEAN, // matched
+ LOWPU, // merged, but dominant contributor
+ HIGHPU, // merged and dominated by pile-up or fakes
+ HSSPLIT, // split
+ NONE, // not found
+ NHSTYPES
+ };
+
+ //Find the best reco vertex matching the hard scatter event
+ //returns 0 in case can't find one
+ const xAOD::Vertex * bestHardScatterMatch( const xAOD::VertexContainer & vxContainer );
+
+ //pointers to all reco vertices for which the hard scatter has some contribution
+ //includes the index in the matching info where to find it
+ const std::vector<std::pair<const xAOD::Vertex*, size_t> > hardScatterMatches( const xAOD::VertexContainer & vxContainer );
+
+ //find the hard scatter type
+ HardScatterType classifyHardScatter( const xAOD::VertexContainer & vxContainer );
+
+}
+
+
+#endif
diff --git a/InnerDetector/InDetValidation/InDetTruthVertexValidation/Root/InDetVertexTruthMatchTool.cxx b/InnerDetector/InDetValidation/InDetTruthVertexValidation/Root/InDetVertexTruthMatchTool.cxx
new file mode 100644
index 00000000000..2869b9d284c
--- /dev/null
+++ b/InnerDetector/InDetValidation/InDetTruthVertexValidation/Root/InDetVertexTruthMatchTool.cxx
@@ -0,0 +1,422 @@
+/*
+ Copyright (C) 2002-2018 CERN for the benefit of the ATLAS collaboration
+*/
+
+#include "InDetTruthVertexValidation/InDetVertexTruthMatchTool.h"
+
+#include "xAODTracking/TrackParticleContainer.h"
+#include "xAODTruth/TruthEventContainer.h"
+
+#include "InDetTruthVertexValidation/InDetVertexTruthMatchUtils.h"
+using namespace InDetVertexTruthMatchUtils;
+
+InDetVertexTruthMatchTool::InDetVertexTruthMatchTool( const std::string & name ) : asg::AsgTool(name) {
+ declareProperty("trackMatchProb", m_trkMatchProb = 0.7 );
+ declareProperty("vertexMatchWeight", m_vxMatchWeight = 0.7 );
+ declareProperty("trackPtCut", m_trkPtCut = 100. );
+}
+
+StatusCode InDetVertexTruthMatchTool::initialize() {
+ ATH_MSG_INFO("Initializing");
+
+ return StatusCode::SUCCESS;
+}
+
+
+namespace {
+//Helper methods for this file only
+
+//output typedef
+// this is defined like this in InDetTruthMatchUtils:
+// typedef std::pair<ElementLink<xAOD::TruthEventBaseContainer>, float> VertexTruthMatchInfo;
+// std::pair<ElementLink<>, T> is special templated pair in ElementLink.h even
+// pair of link to a truth event with relative weight of matched tracks
+
+
+// Create a truth map by decoarting truth particles with a back link to the truth event they live in
+// Needed because the track->truth assoc gives us the particles but they don't store event normally
+// Add as decoration to avoid full loop for every track ( this time only once per event )
+// Use a vector so any number of truth event collections can be used at once -- but the pointers need to be valid
+void createTruthMap(std::vector<const xAOD::TruthEventBaseContainer *> truthEventContainers ) {
+
+ xAOD::TruthParticle::Decorator<ElementLink<xAOD::TruthEventBaseContainer> > backLinkDecor("TruthEventLink");
+
+ for ( auto cit : truthEventContainers ) {
+
+ const xAOD::TruthEventBaseContainer & truthEvents = *cit;
+
+ for ( size_t i = 0; i < truthEvents.size(); ++i) {
+
+ for ( auto & tkit : truthEvents[i]->truthParticleLinks() ) { //std::vector<ElementLink...
+
+ const ElementLink<xAOD::TruthEventBaseContainer> elLink = ElementLink<xAOD::TruthEventBaseContainer>( truthEvents, i );
+
+ if (elLink.isValid() && tkit.isValid()) {
+ backLinkDecor(**tkit) = elLink;
+ }
+
+ }
+
+ }
+
+ }
+
+}
+
+void createTrackTruthMap(std::vector<const xAOD::TruthEventBaseContainer *> truthEventContainers,
+ const xAOD::TrackParticleContainer & trackParticleContainer,
+ float matchCut)
+{
+
+ createTruthMap(truthEventContainers);
+
+ xAOD::TrackParticle::ConstAccessor<ElementLink<xAOD::TruthParticleContainer> > trk_truthPartAcc("truthParticleLink");
+ xAOD::TrackParticle::ConstAccessor<float> trk_truthProbAcc("truthMatchProbability");
+ xAOD::TruthParticle::Decorator<ElementLink<xAOD::TruthEventBaseContainer> > backLinkDecor("TruthEventLink");
+ xAOD::TrackParticle::Decorator<ElementLink<xAOD::TruthEventBaseContainer> > trackLinkDecor("TrackEventLink");
+
+ int nGood = 0;
+ int nMatch = 0;
+ int nLink = 0;
+ for (auto trk : trackParticleContainer)
+ {
+ {
+ nGood++;
+ if (trk_truthPartAcc.isAvailable(*trk) && trk_truthProbAcc.isAvailable(*trk)
+ && trk_truthPartAcc(*trk).isValid() && trk_truthProbAcc(*trk) >= matchCut)
+ {
+ nMatch++;
+ auto truthParticle = trk_truthPartAcc(*trk);
+ if (backLinkDecor.isAvailable(**truthParticle) && backLinkDecor(**truthParticle).isValid())
+ {
+ nLink++;
+ trackLinkDecor(*trk) = backLinkDecor(**truthParticle);
+ }
+ }
+ }
+ }
+ // won't compile, no idea why
+ //ATH_MSG_DEBUG("Linked/Matched/Good/All: " << nLink << " / " << nMatch << " / " << nGood << " / " << trackParticleContainer.size());
+}
+
+//In the vector of match info, find the element corresponding to link and return its index; create a new one if necessary
+size_t indexOfMatchInfo( std::vector<VertexTruthMatchInfo> & matches, ElementLink<xAOD::TruthEventBaseContainer> & link ) {
+ for ( size_t i = 0; i < matches.size(); ++i ) {
+ if ( link.key() == std::get<0>(matches[i]).key() && link.index() == std::get<0>(matches[i]).index() )
+ return i;
+ }
+ // This is the first time we've seen this truth vertex, so make a new entry
+ matches.push_back( std::make_tuple( link, 0., 0. ) );
+ return matches.size() - 1;
+}
+
+
+//for sorting the container -> highest relative match weight first
+bool compareMatchPair(const VertexTruthMatchInfo a, const VertexTruthMatchInfo b ) { return std::get<1>(a) > std::get<1>(b); }
+
+}
+
+
+const xAOD::TrackParticleContainer*
+InDetVertexTruthMatchTool::findTrackParticleContainer( const xAOD::VertexContainer& vxContainer )
+{
+ for (auto vtx : vxContainer)
+ {
+ for (const ElementLink<xAOD::TrackParticleContainer>& tpLink : vtx->trackParticleLinks())
+ {
+ if (tpLink.isValid())
+ {
+ return tpLink.getStorableObjectPointer();
+ }
+ }
+ }
+ return 0;
+}
+
+StatusCode InDetVertexTruthMatchTool::matchVertices( const xAOD::VertexContainer & vxContainer ) {
+
+ ATH_MSG_DEBUG("Start vertex matching");
+
+ // Identify MC vertices to match to -- this is the collection for hard scatter
+ const xAOD::TruthEventBaseContainer * truthEvents = nullptr;
+ if ( evtStore()->contains<xAOD::TruthEventBaseContainer>( "TruthEvents" ) )
+ ATH_CHECK( evtStore()->retrieve( truthEvents, "TruthEvents" ) );
+ else
+ ATH_CHECK( evtStore()->retrieve( truthEvents, "TruthEvent" ) );
+
+ std::vector<const xAOD::TruthEventBaseContainer *> truthContainers;
+ truthContainers.push_back( truthEvents );
+
+ ATH_MSG_DEBUG("Found Hard Scatter collection");
+
+ // These are the pile-up truth -- don't want to fail if they don't exist
+ const xAOD::TruthEventBaseContainer * truthPileup = nullptr;
+ if ( evtStore()->contains<xAOD::TruthEventBaseContainer>( "TruthPileupEvents" ) )
+ ATH_CHECK( evtStore()->retrieve( truthPileup, "TruthPileupEvents" ) );
+ if (truthPileup)
+ truthContainers.push_back( truthPileup );
+
+ ATH_MSG_DEBUG("Found Pileup collection");
+
+ // Find the trackParticle container associated with our reconstructed vertices
+ // We could pass this, but it would break the original interface...
+ const xAOD::TrackParticleContainer* tkContainer = findTrackParticleContainer(vxContainer);
+ if (!tkContainer)
+ {
+ ATH_MSG_WARNING("Vertex container has no vertices with valid TrackParticle links");
+ return StatusCode::SUCCESS;
+ }
+
+ ATH_MSG_DEBUG("Found track collection");
+
+ // create the particle links to events to avoid excessive looping
+ // also decorate reconstructed tracks passing selection with truthEvent links
+ createTrackTruthMap( truthContainers, *tkContainer, m_trkMatchProb );
+
+ // Accessor for the links we just created
+ xAOD::TruthParticle::Decorator<ElementLink<xAOD::TruthEventBaseContainer> > backLinkDecor("TruthEventLink");
+
+ //setup decorators for truth matching info
+ xAOD::Vertex::Decorator<std::vector<VertexTruthMatchInfo> > matchInfoDecor("TruthEventMatchingInfos");
+ xAOD::Vertex::Decorator<std::vector<VertexTruthMatchInfo> > rawMatchInfoDecor("TruthEventRawMatchingInfos");
+ xAOD::Vertex::Decorator<VertexMatchType> matchTypeDecor("VertexMatchType");
+ xAOD::Vertex::Decorator<std::vector<ElementLink<xAOD::VertexContainer> > > splitPartnerDecor("SplitPartners");
+
+ //setup accessors
+ // can switch to built in method in xAOD::Vertex once don't have to deal with changing names anymore
+ xAOD::Vertex::ConstAccessor<xAOD::Vertex::TrackParticleLinks_t> trkAcc("trackParticleLinks");
+ xAOD::Vertex::ConstAccessor<std::vector<float> > weightAcc("trackWeights");
+
+ xAOD::TrackParticle::ConstAccessor<ElementLink<xAOD::TruthParticleContainer> > trk_truthPartAcc("truthParticleLink");
+ xAOD::TrackParticle::ConstAccessor<float> trk_truthProbAcc("truthMatchProbability");
+
+ xAOD::TrackParticle::Decorator<ElementLink<xAOD::VertexContainer> > trk_recoVtx("RecoVertex");
+ xAOD::TrackParticle::Decorator<float> trk_wtVtx("WeightVertex");
+
+ //some variables to store
+ size_t ntracks;
+ xAOD::VxType::VertexType vxType;
+
+ ATH_MSG_DEBUG("Starting Loop on Vertices");
+
+ //=============================================================================
+ //First loop over vertices: get tracks, then TruthParticles, and store relative
+ //weights of contribution from each TruthEvent
+ //=============================================================================
+ size_t vxEntry = 0;
+
+ for ( auto vxit : vxContainer.stdcont() ) {
+ vxEntry++;
+ vxType = static_cast<xAOD::VxType::VertexType>( vxit->vertexType() );
+ if (vxType == xAOD::VxType::NoVtx) {
+ //skip dummy vertices -> match info will be empty vector if someone tries to access later
+ //type will be set to dummy
+ continue;
+ }
+
+
+ //create the vector we will add as matching info decoration later
+ std::vector<VertexTruthMatchInfo> matchinfo;
+ std::vector<VertexTruthMatchInfo> rawMatchinfo; //not normalized to one for each reco vertex
+
+ //if don't have track particles
+ if (!trkAcc.isAvailable(*vxit) || !weightAcc.isAvailable(*vxit) ) {
+ ATH_MSG_DEBUG("trackParticles or trackWeights not available, setting fake");
+ // Add invalid link for fakes
+ matchinfo.push_back( std::make_tuple( ElementLink<xAOD::TruthEventBaseContainer>(), 1., 0. ) );
+ matchInfoDecor( *vxit ) = matchinfo;
+ rawMatchinfo.push_back( std::make_tuple( ElementLink<xAOD::TruthEventBaseContainer>(), 1., 0. ) );
+ rawMatchInfoDecor( *vxit ) = rawMatchinfo;
+ continue;
+ }
+
+ //things we need to do the matching
+ const xAOD::Vertex::TrackParticleLinks_t & trkParts = trkAcc( *vxit );
+ ntracks = trkParts.size();
+ const std::vector<float> & trkWeights = weightAcc( *vxit );
+
+ //double check
+ if ( trkWeights.size() != ntracks ) {
+ ATH_MSG_DEBUG("Vertex without same number of tracks and trackWeights, setting fake");
+ matchinfo.push_back( std::make_tuple( ElementLink<xAOD::TruthEventBaseContainer>(), 1., 0. ) );
+ matchInfoDecor( *vxit ) = matchinfo;
+ rawMatchinfo.push_back( std::make_tuple( ElementLink<xAOD::TruthEventBaseContainer>(), 1., 0. ) );
+ rawMatchInfoDecor( *vxit ) = rawMatchinfo;
+ continue;
+ }
+
+ ATH_MSG_DEBUG("Matching new vertex at (" << vxit->x() << ", " << vxit->y() << ", " << vxit->z() << ")" << " with " << ntracks << " tracks:");
+
+ float totalWeight = 0.;
+ float totalFake = 0.;
+
+ //loop element link to track particle
+ for ( size_t t = 0; t < ntracks; ++t ) {
+ const xAOD::TrackParticle & trk = **trkParts[t];
+
+ totalWeight += trkWeights[t];
+ trk_recoVtx(trk) = ElementLink<xAOD::VertexContainer>(vxContainer, vxEntry - 1);
+ trk_wtVtx(trk) = trkWeights[t];
+
+ const ElementLink<xAOD::TruthParticleContainer> & truthPartLink = trk_truthPartAcc( trk );
+ float prob = trk_truthProbAcc( trk );
+
+ if (truthPartLink.isValid() && prob > m_trkMatchProb) {
+ const xAOD::TruthParticle & truthPart = **truthPartLink;
+ //check if the truth particle is "good"
+ if ( pass( truthPart) ) {
+ ElementLink<xAOD::TruthEventBaseContainer> match = backLinkDecor( truthPart );
+ //check we have an actual link
+ if ( match.isValid() ) {
+ size_t matchIdx = indexOfMatchInfo( matchinfo, match );
+ std::get<1>(matchinfo[matchIdx]) += trkWeights[t];
+ std::get<2>(matchinfo[matchIdx]) += (trk.pt()/1000.) * (trk.pt()/1000.) * trkWeights[t];
+ matchIdx = indexOfMatchInfo( rawMatchinfo, match );
+ std::get<1>(rawMatchinfo[matchIdx]) += trkWeights[t];
+ std::get<2>(rawMatchinfo[matchIdx]) += (trk.pt()/1000.) * (trk.pt()/1000.) * trkWeights[t];
+ } else {
+ totalFake += trkWeights[t];
+ }
+
+ } else {
+ //truth particle failed cuts -> add to fakes
+ totalFake += trkWeights[t];
+ }
+ } else {
+ //not valid or low matching probability -> add to fakes
+ totalFake += trkWeights[t];
+ }
+ }//end loop over tracks in vertex
+
+ //finalize the match info vector
+ if ( totalWeight < 1e-12 ) { // in case we don't have any passing tracks we want to make sure labeled fake
+ ATH_MSG_DEBUG(" Declaring vertex fully fake (no passing tracks included)");
+ totalWeight = 1.;
+ totalFake = 1.;
+ }
+ if ( totalFake > 0. )
+ {
+ matchinfo.push_back( std::make_tuple( ElementLink<xAOD::TruthEventBaseContainer>(), totalFake, 0. ) );
+ rawMatchinfo.push_back( std::make_tuple( ElementLink<xAOD::TruthEventBaseContainer>(), totalFake, 0. ) );
+ }
+
+ for ( auto & mit : matchinfo ) {
+ std::get<1>(mit) /= totalWeight;
+ }
+ std::sort( matchinfo.begin(), matchinfo.end(), compareMatchPair );
+ std::sort( rawMatchinfo.begin(), rawMatchinfo.end(), compareMatchPair );
+ matchInfoDecor( *vxit ) = matchinfo;
+ rawMatchInfoDecor( *vxit ) = rawMatchinfo;
+ }
+
+ //After first loop, all vertices have been decorated with their vector of match info (link to TruthEvent paired with weight)
+ //now we want to use that information from the whole collection to assign types
+
+ //keep track of whether a type is assigned
+ //useful since looking for splits involves a double loop, and then setting types ahead in the collection
+ std::vector<bool> assignedType( vxContainer.size(), false );
+
+ for ( size_t i = 0; i < vxContainer.size(); ++i ) {
+
+ if ( assignedType[i] ) continue; // make sure we don't reclassify vertices already found in the split loop below
+
+ std::vector<VertexTruthMatchInfo> & info = matchInfoDecor( *vxContainer[i] );
+ if (info.size() == 0) {
+ matchTypeDecor( *vxContainer[i] ) = DUMMY;
+ } else if ( !std::get<0>(info[0]).isValid() ) {
+ matchTypeDecor( *vxContainer[i] ) = FAKE;
+ } else if ( std::get<1>(info[0]) > m_vxMatchWeight ) {
+ matchTypeDecor( *vxContainer[i] ) = MATCHED;
+ } else {
+ matchTypeDecor( *vxContainer[i] ) = MERGED;
+ }
+
+ //check for splitting
+ if ( matchTypeDecor( *vxContainer[i] ) == MATCHED || matchTypeDecor( *vxContainer[i] ) == MERGED ) {
+ std::vector<size_t> foundSplits;
+ for ( size_t j = i + 1; j < vxContainer.size(); ++j ) {
+ std::vector<VertexTruthMatchInfo> & info2 = matchInfoDecor( *vxContainer[j] );
+ //check second vertex is not dummy or fake, and that it has same elementlink as first vertex
+ //equality test is in code but doesnt seem to work for ElementLinks that I have?
+ //so i am just checking that the contianer key hash and the index are the same
+ if (matchTypeDecor( *vxContainer[j] ) == FAKE || matchTypeDecor( *vxContainer[j] ) == DUMMY) continue;
+ if (info2.size() > 0 && std::get<0>(info2[0]).isValid() && std::get<0>(info[0]).key() == std::get<0>(info2[0]).key() && std::get<0>(info[0]).index() == std::get<0>(info2[0]).index() ) {
+ //add split links; first between first one found and newest one
+ splitPartnerDecor( *vxContainer[i] ).push_back( ElementLink<xAOD::VertexContainer>( vxContainer, j ) );
+ splitPartnerDecor( *vxContainer[j] ).push_back( ElementLink<xAOD::VertexContainer>( vxContainer, i ) );
+ //then between any others we found along the way
+ for ( auto k : foundSplits ) { //k is a size_t in the vector of splits
+ splitPartnerDecor( *vxContainer[k] ).push_back( ElementLink<xAOD::VertexContainer>( vxContainer, j ) );
+ splitPartnerDecor( *vxContainer[j] ).push_back( ElementLink<xAOD::VertexContainer>( vxContainer, k ) );
+ }
+ //then keep track that we found this one
+ foundSplits.push_back(j);
+ } //if the two vertices match to same TruthEvent
+ }//inner loop over vertices
+
+ // Correct labelling of split vertices - keep highest sumpt2 vertex labelled as matched/merged
+ float maxSumpT2 = std::get<2>( matchInfoDecor( *vxContainer[i] )[0] );
+ size_t indexOfMax = i;
+ for ( auto l : foundSplits ) {
+ if ( std::get<2>( matchInfoDecor( *vxContainer[l] )[0] ) > maxSumpT2 ){
+ maxSumpT2 = std::get<2>( matchInfoDecor( *vxContainer[l] )[0] );
+ indexOfMax = l;
+ } else {
+ matchTypeDecor( *vxContainer[l] ) = SPLIT;
+ assignedType[l] = true;
+ }
+ }
+ if ( indexOfMax!=i ) matchTypeDecor( *vxContainer[i] ) = SPLIT;
+ } //if matched or merged
+ } //outer loop
+
+ //DEBUG MATCHING
+ if (msgLvl(MSG::DEBUG)) {
+ for (const auto &vxit : vxContainer.stdcont() ) {
+ ATH_MSG_DEBUG("Matched vertex (index " << (*vxit).index() << ") to type " << matchTypeDecor(*vxit) << " with following info of size " << matchInfoDecor(*vxit).size() << ":");
+ for (const auto &vit : matchInfoDecor(*vxit) ) {
+ if ( std::get<0>(vit).isValid() ) {
+ ATH_MSG_DEBUG(" GenEvent type " << (* std::get<0>(vit))->type() << ", index " << std::get<0>(vit).index() << " with relative weight " << std::get<1>(vit) );
+ } else {
+ ATH_MSG_DEBUG(" Fakes with relative weight " << std::get<1>(vit) );
+ }
+ }
+ if (matchTypeDecor(*vxit) == SPLIT) {
+ ATH_MSG_DEBUG(" Split partners are:");
+ for (const auto &split : splitPartnerDecor( *vxit ) ) {
+ if ( split.isValid() )
+ ATH_MSG_DEBUG(" Vertex " << split.index());
+ else
+ ATH_MSG_DEBUG(" ERROR");
+ }
+ }
+ }
+ }
+
+ return StatusCode::SUCCESS;
+
+}
+
+
+//Set up any cuts on either the tracks or truth particles to allow here
+//A failing track is removed from consideration entirely
+//If a passing track matches to a failing truth particle it will be considered "fake"
+
+bool InDetVertexTruthMatchTool::pass( const xAOD::TruthParticle & truthPart ) const {
+
+ //remove the registered secondaries
+ if ( truthPart.barcode() > 200000 ) return false;
+
+ return true;
+
+}
+
+/*
+bool InDetVertexTruthMatchTool::pass( const xAOD::TrackParticle & trackPart ) {
+
+ if( trackPart.pt() < m_trkPtCut ) return false;
+
+ return true;
+
+}
+*/
diff --git a/InnerDetector/InDetValidation/InDetTruthVertexValidation/Root/InDetVertexTruthMatchUtils.cxx b/InnerDetector/InDetValidation/InDetTruthVertexValidation/Root/InDetVertexTruthMatchUtils.cxx
new file mode 100644
index 00000000000..91f5fc1bfe1
--- /dev/null
+++ b/InnerDetector/InDetValidation/InDetTruthVertexValidation/Root/InDetVertexTruthMatchUtils.cxx
@@ -0,0 +1,128 @@
+/*
+ Copyright (C) 2002-2018 CERN for the benefit of the ATLAS collaboration
+*/
+
+#include "InDetTruthVertexValidation/InDetVertexTruthMatchUtils.h"
+#include "xAODTracking/VertexContainer.h"
+
+namespace InDetVertexTruthMatchUtils {
+
+namespace {
+
+bool isHardScatterEvent( const ElementLink<xAOD::TruthEventBaseContainer> evlink ) {
+ //good link of type "TruthEvent" (not truthpileupevent) and the first one in the collection
+ return ( evlink.isValid() && (*evlink)->type() == xAOD::Type::TruthEvent && evlink.index() == 0 );
+}
+
+}
+
+//find highest weight match to the hard scatter interaction
+//in the case of splitting user can take care to check this is really "best"
+// --> example:
+// vertex 0 is 90% hard scatter with many tracks, and 10% another vertex
+// vertex 1 is 100% hard scatter with few tracks.
+// this algorithm will pick vertex 1, but user may want to use vertex 0 (accessible from splitting link decoration)
+const xAOD::Vertex * bestHardScatterMatch( const xAOD::VertexContainer & vxContainer ) {
+ //accessors for decorations
+ xAOD::Vertex::Decorator<std::vector<VertexTruthMatchInfo> > matchInfoDecor("TruthEventMatchingInfos");
+
+ std::vector<std::pair<const xAOD::Vertex*, size_t> > allMatches = hardScatterMatches( vxContainer );
+
+ //look for the highest weight
+ float bestw = 0.;
+ const xAOD::Vertex * best = nullptr;
+ for ( auto it : allMatches ) {
+ std::vector<VertexTruthMatchInfo> & info = matchInfoDecor( *(it.first) );
+ //go to the index of the HS match in the info vector, and check the weight
+ if ( std::get<1>(info[ it.second ]) > bestw ) {
+ bestw = std::get<1>(info[ it.second ]);
+ best = it.first;
+ }
+ }
+
+ //if didn't find any matches will return 0
+ return best;
+}
+
+//Find all hard scatter matches
+const std::vector<std::pair<const xAOD::Vertex*, size_t> > hardScatterMatches( const xAOD::VertexContainer & vxContainer ) {
+ //accessors for decorations
+ xAOD::Vertex::Decorator<std::vector<VertexTruthMatchInfo> > matchInfoDecor("TruthEventMatchingInfos");
+
+ //return vector
+ std::vector<std::pair<const xAOD::Vertex*, size_t> > result;
+
+ //loop and look
+ for ( auto vxit : vxContainer ) {
+ const std::vector<VertexTruthMatchInfo> & info = matchInfoDecor( *vxit );
+ for ( size_t i = 0; i < info.size(); ++i ) {
+ if ( isHardScatterEvent( std::get<0>(info[i]) ) ) {
+ result.push_back( std::make_pair(vxit, i) );
+ break;
+ }
+ }
+
+ }
+ return result;
+}
+
+HardScatterType classifyHardScatter( const xAOD::VertexContainer & vxContainer ) {
+
+ // return a vector of reco vertices with contributions from the hard-scatter
+ // if the second element in the pair is zero, it means tracks from the hard-scatter
+ // contribute more weight to that reco vertex than any other interaction
+ std::vector<std::pair<const xAOD::Vertex*, size_t> > matches = hardScatterMatches( vxContainer );
+
+ // information indicating whether a given reco vertex is:
+ // MATCHED, MERGED, SPLIT or FAKE
+ xAOD::Vertex::Decorator<VertexMatchType> matchTypeDecor("VertexMatchType");
+
+ if ( matches.size() == 0 ) {
+ // No reco vertices with tracks from hard-scatter
+ return NONE;
+ } else if ( matches.size() == 1 ) {
+ // A single reco vertex has tracks from hard-scatter
+ const VertexMatchType & type = matchTypeDecor( *(matches[0].first) );
+ //check if the index in the matching info for HS is first, then can assign clean/lowpu based on match/merge
+ if ( matches[0].second == 0 && type == MATCHED ) {
+ // HS made largest contribution and vertex matched
+ return CLEAN;
+ } else if ( matches[0].second == 0 && type == MERGED ) {
+ // HS made largest contribution and vertex merged
+ return LOWPU;
+ } else {
+ // HS did not make largest contribution or vertex is fake or vertex is part of non-HS split pair
+ return HIGHPU;
+ }
+
+ } else {
+ // multiple reco vertices have tracks from hard-scatter
+ // count how many have hard-scatter tracks as largest contribution
+ int num_main = 0;
+ const xAOD::Vertex* mainVtx = nullptr;
+ for ( auto it : matches ) {
+ if ( it.second == 0 ) // is hard-scatter the largest contribution?
+ {
+ num_main++;
+ mainVtx = it.first;
+ }
+ }
+ if (num_main == 0 ) {
+ return HIGHPU;
+ } else if (num_main == 1 ) {
+ const VertexMatchType& type = matchTypeDecor(*mainVtx);
+ if (type == MATCHED)
+ return CLEAN;
+ else if (type == MERGED)
+ return LOWPU;
+ else
+ return HIGHPU;
+ } else {
+ return HSSPLIT;
+ }
+ }
+
+ return NONE;
+}
+
+}
diff --git a/InnerDetector/InDetValidation/InDetTruthVertexValidation/Root/LinkDef.h b/InnerDetector/InDetValidation/InDetTruthVertexValidation/Root/LinkDef.h
new file mode 100644
index 00000000000..6477846377b
--- /dev/null
+++ b/InnerDetector/InDetValidation/InDetTruthVertexValidation/Root/LinkDef.h
@@ -0,0 +1,23 @@
+/*
+ Copyright (C) 2002-2018 CERN for the benefit of the ATLAS collaboration
+*/
+
+#ifndef INDETRUTHVERTEXVALIDATION_LINKDEF_H
+#define INDETRUTHVERTEXVALIDATION_LINKDEF_H
+
+#include "InDetTruthVertexValidation/InDetVertexTruthMatchUtils.h"
+#include "InDetTruthVertexValidation/InDetVertexTruthMatchTool.h"
+
+#ifdef __CINT__
+
+#pragma link off all globals;
+#pragma link off all classes;
+#pragma link off all functions;
+#pragma link C++ nestedclass;
+
+#pragma link C++ class InDetVertexTruthMatchTool;
+#pragma link C++ class InDetVertexTruthMatchUtils;
+
+#endif
+
+#endif
\ No newline at end of file
diff --git a/InnerDetector/InDetValidation/InDetTruthVertexValidation/share/VertexTruthMatch_jobOptions.py b/InnerDetector/InDetValidation/InDetTruthVertexValidation/share/VertexTruthMatch_jobOptions.py
new file mode 100644
index 00000000000..4e54cdd4ea3
--- /dev/null
+++ b/InnerDetector/InDetValidation/InDetTruthVertexValidation/share/VertexTruthMatch_jobOptions.py
@@ -0,0 +1,27 @@
+# Copyright (C) 2002-2018 CERN for the benefit of the ATLAS collaboration
+
+# Set up the file reading:
+FNAME = "/afs/cern.ch/atlas/project/PAT/xAODs/r5591/mc14_8TeV.117050.PowhegPythia_P2011C_ttbar.recon.AOD.e1727_s1933_s1911_r5591/AOD.01494882._111853.pool.root.1"
+
+import AthenaPoolCnvSvc.ReadAthenaPool
+ServiceMgr.EventSelector.InputCollections = [ FNAME ]
+
+# Access the algorithm sequence:
+from AthenaCommon.AlgSequence import AlgSequence
+theJob = AlgSequence()
+
+from InDetTruthVertexValidation.InDetTruthVertexValidationConf import InDetVertexTruthMatchTool
+tool = InDetVertexTruthMatchTool()
+tool.OutputLevel=DEBUG
+ToolSvc += tool
+
+# Add the test algorithm:
+from InDetTruthVertexValidation.InDetTruthVertexValidationConf import InDetVertexTruthMatchAlgorithm
+alg = InDetVertexTruthMatchAlgorithm()
+theJob += alg
+
+# Do some additional tweaking:
+from AthenaCommon.AppMgr import theApp
+theApp.EvtMax = 100
+ServiceMgr.MessageSvc.OutputLevel = INFO
+ServiceMgr.MessageSvc.defaultLimit = 1000000
diff --git a/InnerDetector/InDetValidation/InDetTruthVertexValidation/src/InDetVertexTruthMatchAlgorithm.cxx b/InnerDetector/InDetValidation/InDetTruthVertexValidation/src/InDetVertexTruthMatchAlgorithm.cxx
new file mode 100644
index 00000000000..1784403e25d
--- /dev/null
+++ b/InnerDetector/InDetValidation/InDetTruthVertexValidation/src/InDetVertexTruthMatchAlgorithm.cxx
@@ -0,0 +1,43 @@
+/*
+ Copyright (C) 2002-2018 CERN for the benefit of the ATLAS collaboration
+*/
+
+#include "InDetVertexTruthMatchAlgorithm.h"
+
+#include "xAODTracking/VertexContainer.h"
+
+
+InDetVertexTruthMatchAlgorithm::InDetVertexTruthMatchAlgorithm( const std::string& name, ISvcLocator* svcLoc )
+ : AthAlgorithm( name, svcLoc ),
+ m_matchTool( "InDetVertexTruthMatchTool/InDetVertexTruthMatchTool", this ) {
+
+ declareProperty( "SGKey", m_sgKey = "PrimaryVertices" );
+
+ declareProperty( "VertexTruthMatchTool", m_matchTool );
+}
+
+StatusCode InDetVertexTruthMatchAlgorithm::initialize() {
+
+ // Greet the user:
+ ATH_MSG_INFO( "Initialising - Package version: " << PACKAGE_VERSION );
+ ATH_MSG_DEBUG( "SGKey = " << m_sgKey );
+
+ // Retrieve the tools:
+ ATH_CHECK( m_matchTool.retrieve() );
+
+ // Return gracefully:
+ return StatusCode::SUCCESS;
+}
+
+StatusCode InDetVertexTruthMatchAlgorithm::execute() {
+
+ //Retrieve the vertices:
+ const xAOD::VertexContainer * vxContainer = 0;
+ ATH_CHECK( evtStore()->retrieve( vxContainer, m_sgKey ) );
+
+ //pass to the tool for decoration:
+ ATH_CHECK( m_matchTool->matchVertices( *vxContainer ) );
+
+ return StatusCode::SUCCESS;
+
+}
diff --git a/InnerDetector/InDetValidation/InDetTruthVertexValidation/src/InDetVertexTruthMatchAlgorithm.h b/InnerDetector/InDetValidation/InDetTruthVertexValidation/src/InDetVertexTruthMatchAlgorithm.h
new file mode 100644
index 00000000000..2d0e7e98a9e
--- /dev/null
+++ b/InnerDetector/InDetValidation/InDetTruthVertexValidation/src/InDetVertexTruthMatchAlgorithm.h
@@ -0,0 +1,34 @@
+/*
+ Copyright (C) 2002-2018 CERN for the benefit of the ATLAS collaboration
+*/
+
+#ifndef INDETVERTEXTRUTHMATCHALGORITHM_H
+#define INDETVERTEXTRUTHMATCHALGORITHM_H
+
+// Gaudi/Athena include(s):
+#include "AthenaBaseComps/AthAlgorithm.h"
+#include "AsgTools/ToolHandle.h"
+
+// Local include(s):
+#include "InDetTruthVertexValidation/IInDetVertexTruthMatchTool.h"
+
+class InDetVertexTruthMatchAlgorithm : public AthAlgorithm {
+
+ public:
+ /// Regular Algorithm constructor
+ InDetVertexTruthMatchAlgorithm( const std::string& name, ISvcLocator* svcLoc );
+
+ /// Function initialising the algorithm
+ virtual StatusCode initialize() override final;
+ /// Function executing the algorithm
+ virtual StatusCode execute() override final;
+
+ private:
+ /// StoreGate key for the muon container to investigate
+ std::string m_sgKey;
+
+ ToolHandle<IInDetVertexTruthMatchTool> m_matchTool;
+
+};
+
+#endif
diff --git a/InnerDetector/InDetValidation/InDetTruthVertexValidation/src/components/InDetTruthVertexValidation_entries.cxx b/InnerDetector/InDetValidation/InDetTruthVertexValidation/src/components/InDetTruthVertexValidation_entries.cxx
new file mode 100644
index 00000000000..140e6664726
--- /dev/null
+++ b/InnerDetector/InDetValidation/InDetTruthVertexValidation/src/components/InDetTruthVertexValidation_entries.cxx
@@ -0,0 +1,13 @@
+/*
+ Copyright (C) 2002-2018 CERN for the benefit of the ATLAS collaboration
+*/
+
+// Gaudi/Athena include(s):
+//#include "GaudiKernel/DeclareFactoryEntries.h"
+
+// Local include(s):
+#include "InDetTruthVertexValidation/InDetVertexTruthMatchTool.h"
+#include "../InDetVertexTruthMatchAlgorithm.h"
+
+DECLARE_COMPONENT( InDetVertexTruthMatchTool )
+DECLARE_COMPONENT( InDetVertexTruthMatchAlgorithm )
diff --git a/InnerDetector/InDetValidation/InDetTruthVertexValidation/src/components/InDetTruthVertexValidation_load.cxx b/InnerDetector/InDetValidation/InDetTruthVertexValidation/src/components/InDetTruthVertexValidation_load.cxx
new file mode 100644
index 00000000000..b1a3b84230d
--- /dev/null
+++ b/InnerDetector/InDetValidation/InDetTruthVertexValidation/src/components/InDetTruthVertexValidation_load.cxx
@@ -0,0 +1,7 @@
+/*
+ Copyright (C) 2002-2018 CERN for the benefit of the ATLAS collaboration
+*/
+
+#include "GaudiKernel/LoadFactoryEntries.h"
+
+LOAD_FACTORY_ENTRIES(InDetTruthVertexValidation)
diff --git a/InnerDetector/InDetValidation/InDetTruthVertexValidation/util/VertexTruthMatchTest.cxx b/InnerDetector/InDetValidation/InDetTruthVertexValidation/util/VertexTruthMatchTest.cxx
new file mode 100644
index 00000000000..0daa864fa20
--- /dev/null
+++ b/InnerDetector/InDetValidation/InDetTruthVertexValidation/util/VertexTruthMatchTest.cxx
@@ -0,0 +1,145 @@
+/*
+ Copyright (C) 2002-2018 CERN for the benefit of the ATLAS collaboration
+*/
+
+//ROOT
+#include <TFile.h>
+#include <TH1.h>
+
+#ifdef ROOTCORE
+#include "xAODRootAccess/Init.h"
+#include "xAODRootAccess/TEvent.h"
+#include "xAODRootAccess/TStore.h"
+#endif
+
+#include "xAODEventInfo/EventInfo.h"
+#include "xAODTracking/VertexContainer.h"
+#include "xAODTracking/TrackParticleContainer.h"
+
+#include "InDetTruthVertexValidation/InDetVertexTruthMatchTool.h"
+#include "InDetTruthVertexValidation/InDetVertexTruthMatchUtils.h"
+
+int main( int argc, char* argv[] ) {
+
+ const char * APP_NAME = argv[0];
+
+ //check for filename
+ if( argc < 2 ) {
+ Error( APP_NAME, "No file name argument" );
+ Error( APP_NAME, " Usage: %s xAOD_file_name [nEntries]",APP_NAME );
+ return 1;
+ }
+
+ //init app
+ bool result = xAOD::Init( APP_NAME );
+ if(!result) {
+ Error( APP_NAME, "Failed to init");
+ return 1;
+ }
+
+ //open file
+ const char * fileName = argv[1];
+ Info(APP_NAME, "Opening file: %s", fileName );
+ TFile * file = new TFile(fileName,"READ");
+
+ xAOD::TEvent event( xAOD::TEvent::kClassAccess );
+ result = event.readFrom( file );
+ if(!result) {
+ Error( APP_NAME, "Failed to have TEvent read file");
+ return 1;
+ }
+
+ //create transient object store
+ xAOD::TStore store;
+
+ Long64_t entries = event.getEntries();
+ if( argc > 2 ) {
+ const Long64_t e = atoll( argv[2] );
+ if( e < entries )
+ entries = e;
+ }
+
+ //Tool to test
+ InDetVertexTruthMatchTool matchTool("VertexTruthMatchTool");
+
+ //config
+ matchTool.msg().setLevel( MSG::DEBUG );
+ matchTool.setProperty( "trackMatchProb", 0.7 );
+ matchTool.setProperty( "vertexMatchWeight", 0.7 );
+ matchTool.initialize();
+
+ TFile fout("output.root","RECREATE");
+ fout.cd();
+
+ //some inclusive test plots
+ TH1F * h_cat = new TH1F("vertexMatchTypes","Number of vertices of each type",4,-0.5,3.5);
+ h_cat->GetXaxis()->SetBinLabel(1, "Matched");
+ h_cat->GetXaxis()->SetBinLabel(2, "Merged");
+ h_cat->GetXaxis()->SetBinLabel(3, "Split");
+ h_cat->GetXaxis()->SetBinLabel(4, "Fake");
+
+
+ //plot some things for the best reco match to hard-scatter event
+ TH1F * h_cat_hs = new TH1F("vertexTypes_hardScatter","Number of events of each type",5,-0.5,4.5);
+ h_cat_hs->GetXaxis()->SetBinLabel(1, "Clean");
+ h_cat_hs->GetXaxis()->SetBinLabel(2, "Low PU");
+ h_cat_hs->GetXaxis()->SetBinLabel(3, "High PU");
+ h_cat_hs->GetXaxis()->SetBinLabel(4, "Split");
+ h_cat_hs->GetXaxis()->SetBinLabel(5, "No match");
+
+ //Loop
+ for( Long64_t entry = 0; entry < entries; ++entry ) {
+
+ event.getEntry( entry );
+
+ const xAOD::EventInfo * ei = nullptr;
+ result = event.retrieve( ei, "EventInfo" );
+ if(!result) {
+ Error( APP_NAME, "Failed to retrieve EventInfo on entry %i", static_cast<int>(entry));
+ return 1;
+ }
+
+ if( entry % 100 == 0 ) {
+ Info( APP_NAME,
+ "===Start event %i, run %i; %i events processed so far===",
+ static_cast<int>( ei->eventNumber() ),
+ static_cast<int>( ei->runNumber() ),
+ static_cast<int>( entry ) );
+ }
+
+ //get vertices
+ const xAOD::VertexContainer * vxContainer = nullptr;
+ result = event.retrieve( vxContainer, "PrimaryVertices" );
+ if(!result) {
+ Error( APP_NAME, "Failed to retrieve PrimaryVertices on entry %i", static_cast<int>(entry));
+ return 1;
+ }
+
+ int status = matchTool.matchVertices( *vxContainer );
+ if(!status) {
+ Error( APP_NAME, "Bad status from matching tool on entry %i", static_cast<int>(entry));
+ continue;
+ }
+
+ //accessor for the matching type
+ xAOD::Vertex::Decorator<InDetVertexTruthMatchUtils::VertexMatchType> matchTypeDecor("VertexMatchType");
+
+ //fill the category histogram for all vertices
+ for( auto vxit : *vxContainer ) {
+ h_cat->Fill( matchTypeDecor( *vxit ) );
+ }
+
+ //get the matching type for hard scatter and fill
+ h_cat_hs->Fill( InDetVertexTruthMatchUtils::classifyHardScatter( *vxContainer ) );
+
+ //get the best match for the hard scatter
+ //const xAOD::Vertex * hsMatch = InDetVertexTruthMatchUtils::bestHardScatterMatch( *vxContainer );
+
+ }
+
+ h_cat->Write();
+ h_cat_hs->Write();
+ fout.Close();
+
+ return 0;
+}
--
GitLab