diff --git a/InnerDetector/InDetRecTools/InDetSecVtxTruthMatchTool/CMakeLists.txt b/InnerDetector/InDetRecTools/InDetSecVtxTruthMatchTool/CMakeLists.txt
new file mode 100644
index 0000000000000000000000000000000000000000..5e9de9c0e9090aee90b1aeefac69ce74419034eb
--- /dev/null
+++ b/InnerDetector/InDetRecTools/InDetSecVtxTruthMatchTool/CMakeLists.txt
@@ -0,0 +1,24 @@
+# Copyright (C) 2002-2024 CERN for the benefit of the ATLAS collaboration
+
+# Declare the package name:
+atlas_subdir( InDetSecVtxTruthMatchTool )
+
+# External(s):
+find_package( ROOT COMPONENTS Core RIO Hist MathCore )
+
+atlas_add_library( InDetSecVtxTruthMatchToolLib
+ InDetSecVtxTruthMatchTool/*.h Root/*.h Root/*.cxx
+ PUBLIC_HEADERS InDetSecVtxTruthMatchTool
+ LINK_LIBRARIES PATCoreLib AsgAnalysisInterfaces xAODTracking xAODTruth)
+
+if( NOT XAOD_STANDALONE )
+ atlas_add_component( InDetSecVtxTruthMatchTool
+ src/components/*.cxx
+ LINK_LIBRARIES GaudiKernel AthenaBaseComps xAODTracking xAODTruth
+ InDetSecVtxTruthMatchToolLib)
+endif()
+
+atlas_add_dictionary( InDetSecVtxTruthMatchToolDict
+ InDetSecVtxTruthMatchTool/InDetSecVtxTruthMatchToolDict.h
+ InDetSecVtxTruthMatchTool/selection.xml
+ LINK_LIBRARIES InDetSecVtxTruthMatchToolLib )
diff --git a/InnerDetector/InDetRecTools/InDetSecVtxTruthMatchTool/InDetSecVtxTruthMatchTool/ATLAS_CHECK_THREAD_SAFETY b/InnerDetector/InDetRecTools/InDetSecVtxTruthMatchTool/InDetSecVtxTruthMatchTool/ATLAS_CHECK_THREAD_SAFETY
new file mode 100644
index 0000000000000000000000000000000000000000..4633f3803ed1117421017158787e508751d3ba84
--- /dev/null
+++ b/InnerDetector/InDetRecTools/InDetSecVtxTruthMatchTool/InDetSecVtxTruthMatchTool/ATLAS_CHECK_THREAD_SAFETY
@@ -0,0 +1 @@
+InnerDetector/InDetRecTools/InDetSecVtxTruthMatchTool
\ No newline at end of file
diff --git a/InnerDetector/InDetValidation/InDetSecVertexValidation/InDetSecVertexValidation/IInDetSecVertexTruthMatchTool.h b/InnerDetector/InDetRecTools/InDetSecVtxTruthMatchTool/InDetSecVtxTruthMatchTool/IInDetSecVtxTruthMatchTool.h
similarity index 55%
rename from InnerDetector/InDetValidation/InDetSecVertexValidation/InDetSecVertexValidation/IInDetSecVertexTruthMatchTool.h
rename to InnerDetector/InDetRecTools/InDetSecVtxTruthMatchTool/InDetSecVtxTruthMatchTool/IInDetSecVtxTruthMatchTool.h
index c12297efffd74c3537989950a0d403a83b23a930..eb4c813c22a0f09cd4cea5ff715698626ed049b8 100644
--- a/InnerDetector/InDetValidation/InDetSecVertexValidation/InDetSecVertexValidation/IInDetSecVertexTruthMatchTool.h
+++ b/InnerDetector/InDetRecTools/InDetSecVtxTruthMatchTool/InDetSecVtxTruthMatchTool/IInDetSecVtxTruthMatchTool.h
@@ -2,8 +2,8 @@
Copyright (C) 2002-2022 CERN for the benefit of the ATLAS collaboration
*/
-#ifndef IInDetSecVertexTruthMatchTool_h
-#define IInDetSecVertexTruthMatchTool_h
+#ifndef IInDetSecVtxTruthMatchTool_h
+#define IInDetSecVtxTruthMatchTool_h
// Framework include(s):
#include "AsgTools/IAsgTool.h"
@@ -19,17 +19,16 @@
* Categorize reconstructed vertices depending on their composition.
*/
-class IInDetSecVertexTruthMatchTool : public virtual asg::IAsgTool {
+class IInDetSecVtxTruthMatchTool : public virtual asg::IAsgTool {
-ASG_TOOL_INTERFACE( IInDetSecVertexTruthMatchTool )
+ASG_TOOL_INTERFACE( IInDetSecVtxTruthMatchTool )
public:
//take const collection of vertices, match them, and decorate with matching info
- virtual StatusCode matchVertices( const xAOD::VertexContainer & vtxContainer,
- const xAOD::TruthVertexContainer & truthVtxContainer ) = 0;
-
- virtual StatusCode labelTruthVertices( const xAOD::TruthVertexContainer & truthVtxContainer ) = 0;
+ virtual StatusCode matchVertices( std::vector<const xAOD::Vertex*> recoVerticesToMatch,
+ std::vector<const xAOD::TruthVertex*> truthVerticesToMatch,
+ const xAOD::TrackParticleContainer* trackParticles ) = 0;
};
diff --git a/InnerDetector/InDetRecTools/InDetSecVtxTruthMatchTool/InDetSecVtxTruthMatchTool/InDetSecVtxTruthMatchTool.h b/InnerDetector/InDetRecTools/InDetSecVtxTruthMatchTool/InDetSecVtxTruthMatchTool/InDetSecVtxTruthMatchTool.h
new file mode 100644
index 0000000000000000000000000000000000000000..e4422b42ed44fb79941686e0bd9aad1dcd0db7e8
--- /dev/null
+++ b/InnerDetector/InDetRecTools/InDetSecVtxTruthMatchTool/InDetSecVtxTruthMatchTool/InDetSecVtxTruthMatchTool.h
@@ -0,0 +1,124 @@
+/*
+ Copyright (C) 2002-2024 CERN for the benefit of the ATLAS collaboration
+*/
+
+#ifndef InDetSecVtxTruthMatchTool_h
+#define InDetSecVtxTruthMatchTool_h
+
+// Framework include(s):
+#include "AsgTools/AsgTool.h"
+#include "AsgTools/PropertyWrapper.h"
+
+// EDM include(s):
+#include "xAODTracking/VertexContainerFwd.h"
+#include "xAODTruth/TruthParticleFwd.h"
+#include "xAODTruth/TruthVertexFwd.h"
+#include "xAODTruth/TruthVertexContainer.h"
+#include "InDetSecVtxTruthMatchTool/IInDetSecVtxTruthMatchTool.h"
+
+// standard includes
+#include<vector>
+
+namespace InDetSecVtxTruthMatchUtils {
+
+ //Namespace for useful analysis things on the truth matching decorations applied to the VertexContainer
+
+ //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::TruthVertexContainer>, float, float> VertexTruthMatchInfo;
+
+ //type codes for vertex matching on all vertices
+ enum VertexMatchType {
+ Matched=0, // > threshold (default 50%) 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")
+ Other
+ };
+
+ //type codes for truth vertices
+ //NOTE: types are subsets of subsequent types
+ enum TruthVertexMatchType {
+ Reconstructable=0, // fiducial cuts, >= 2 charged daughters
+ Accepted, // >= 2 reco tracks
+ Seeded, // tracks pass cuts
+ Reconstructed, // matched to reco vtx
+ ReconstructedSplit // matched to multiple vtx
+ };
+
+ inline bool isMatched(int matchInfo) {
+ if (matchInfo & (0x1 << Matched)) return true;
+ return false;
+ }
+ inline bool isMerged(int matchInfo) {
+ if (matchInfo & (0x1 << Merged)) return true;
+ return false;
+ }
+ inline bool isSplit(int matchInfo) {
+ if (matchInfo & (0x1 << Split)) return true;
+ return false;
+ }
+ inline bool isFake(int matchInfo) {
+ if (matchInfo & (0x1 << Fake)) return true;
+ return false;
+ }
+ inline bool isOther(int matchInfo) {
+ if (matchInfo & (0x1 << Other)) return true;
+ return false;
+ }
+
+ inline bool isReconstructable(int matchInfo) {
+ if (matchInfo & (0x1 << Reconstructable)) return true;
+ return false;
+ }
+ inline bool isAccepted(int matchInfo) {
+ if (matchInfo & (0x1 << Accepted)) return true;
+ return false;
+ }
+ inline bool isSeeded(int matchInfo) {
+ if (matchInfo & (0x1 << Seeded)) return true;
+ return false;
+ }
+ inline bool isReconstructed(int matchInfo) {
+ if (matchInfo & (0x1 << Reconstructed)) return true;
+ return false;
+ }
+ inline bool isReconstructedSplit(int matchInfo) {
+ if (matchInfo & (0x1 << ReconstructedSplit)) return true;
+ return false;
+ }
+}
+
+/** Class for vertex truth matching.
+ * Match reconstructed vertices to truth level decay vertices
+ * Categorize reconstructed vertices depending on their composition.
+ */
+class InDetSecVtxTruthMatchTool : public virtual IInDetSecVtxTruthMatchTool,
+ public asg::AsgTool {
+
+ ASG_TOOL_CLASS( InDetSecVtxTruthMatchTool, IInDetSecVtxTruthMatchTool )
+
+ public:
+
+ InDetSecVtxTruthMatchTool( const std::string & name );
+
+ virtual StatusCode initialize() override final;
+
+ // take collection of vertices, match them, and decorate with matching info
+ virtual StatusCode matchVertices( std::vector<const xAOD::Vertex*> recoVerticesToMatch, std::vector<const xAOD::TruthVertex*> truthVerticesToMatch, const xAOD::TrackParticleContainer* trackParticles ) override;
+
+ private:
+
+ Gaudi::Property<float> m_trkMatchProb{this, "trackMatchProb", 0.5, "Required MC match probability to consider track a good match" };
+ Gaudi::Property<float> m_vxMatchWeight{this, "vertexMatchWeight", 0.5, "Relative weight threshold to consider vertex matched"};
+ Gaudi::Property<float> m_trkPtCut{this, "trackPtCut", 1000., "pt cut to apply on tracks"};
+ Gaudi::Property<std::string> m_selectedTrackFlag{this, "selectedTrackFlag", "is_selected", "Aux decoration on tracks for seeding efficiencies"};
+
+ //private methods to check if particles are good to use
+ //returns barcode of LLP production truth vertex
+ int checkProduction( const xAOD::TruthParticle& truthPart, std::vector<const xAOD::TruthVertex*> truthVerticesToMatch ) const;
+ void countReconstructibleDescendentParticles(const xAOD::TruthVertex& signalTruthVertex,
+ std::vector<const xAOD::TruthParticle*>& set, int counter) const;
+ std::vector<int> checkParticle( const xAOD::TruthParticle& part, const xAOD::TrackParticleContainer* tkCont ) const;
+};
+
+#endif
diff --git a/InnerDetector/InDetRecTools/InDetSecVtxTruthMatchTool/InDetSecVtxTruthMatchTool/InDetSecVtxTruthMatchToolDict.h b/InnerDetector/InDetRecTools/InDetSecVtxTruthMatchTool/InDetSecVtxTruthMatchTool/InDetSecVtxTruthMatchToolDict.h
new file mode 100644
index 0000000000000000000000000000000000000000..69cbdf1a526275ae71d631b1fc84a7432caf420e
--- /dev/null
+++ b/InnerDetector/InDetRecTools/InDetSecVtxTruthMatchTool/InDetSecVtxTruthMatchTool/InDetSecVtxTruthMatchToolDict.h
@@ -0,0 +1,16 @@
+/*
+ Copyright (C) 2002-2024 CERN for the benefit of the ATLAS collaboration
+*/
+
+#ifndef INDETSECVTXTRUTHMATCHTOOL_INDETSECVTXTRUTHMATCHTOOLDICT_H
+#define INDETSECVTXTRUTHMATCHTOOL_INDETSECVTXTRUTHMATCHTOOLDICT_H
+
+#if defined(__GCCXML__) and not defined(EIGEN_DONT_VECTORIZE)
+#define EIGEN_DONT_VECTORIZE
+#endif // __GCCXML__
+
+#include "InDetSecVtxTruthMatchTool/IInDetSecVtxTruthMatchTool.h"
+#include "InDetSecVtxTruthMatchTool/InDetSecVtxTruthMatchTool.h"
+
+#endif // not INDETSECVTXTRUTHMATCHTOOL_INDETSECVTXTRUTHMATCHTOOLDICT_H
+
diff --git a/InnerDetector/InDetRecTools/InDetSecVtxTruthMatchTool/InDetSecVtxTruthMatchTool/selection.xml b/InnerDetector/InDetRecTools/InDetSecVtxTruthMatchTool/InDetSecVtxTruthMatchTool/selection.xml
new file mode 100644
index 0000000000000000000000000000000000000000..a78768e094a96dd2bf40707c7d5dc0888fa16822
--- /dev/null
+++ b/InnerDetector/InDetRecTools/InDetSecVtxTruthMatchTool/InDetSecVtxTruthMatchTool/selection.xml
@@ -0,0 +1,5 @@
+<lcgdict>
+ <class name="IInDetSecVtxTruthMatchTool" />
+ <class name="InDetSecVtxTruthMatchTool" />
+</lcgdict>
+
diff --git a/InnerDetector/InDetRecTools/InDetSecVtxTruthMatchTool/Root/InDetSecVtxTruthMatchTool.cxx b/InnerDetector/InDetRecTools/InDetSecVtxTruthMatchTool/Root/InDetSecVtxTruthMatchTool.cxx
new file mode 100644
index 0000000000000000000000000000000000000000..cac1fa704cfd568f4de3a74c0670a229082c3268
--- /dev/null
+++ b/InnerDetector/InDetRecTools/InDetSecVtxTruthMatchTool/Root/InDetSecVtxTruthMatchTool.cxx
@@ -0,0 +1,422 @@
+/*
+ Copyright (C) 2002-2024 CERN for the benefit of the ATLAS collaboration
+*/
+#include "InDetSecVtxTruthMatchTool/InDetSecVtxTruthMatchTool.h"
+
+#include "xAODTracking/TrackParticleContainer.h"
+#include "xAODTruth/TruthEventContainer.h"
+
+using namespace InDetSecVtxTruthMatchUtils;
+
+InDetSecVtxTruthMatchTool::InDetSecVtxTruthMatchTool( const std::string & name ) : asg::AsgTool(name) {}
+
+StatusCode InDetSecVtxTruthMatchTool::initialize() {
+ ATH_MSG_INFO("Initializing InDetSecVtxTruthMatchTool");
+
+
+ return StatusCode::SUCCESS;
+}
+
+namespace {
+//Helper methods for this file only
+
+//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, const ElementLink<xAOD::TruthVertexContainer> & 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.emplace_back( link, 0., 0. );
+ return matches.size() - 1;
+}
+
+}
+
+StatusCode InDetSecVtxTruthMatchTool::matchVertices( std::vector<const xAOD::Vertex*> recoVerticesToMatch,
+ std::vector<const xAOD::TruthVertex*> truthVerticesToMatch,
+ const xAOD::TrackParticleContainer* trackParticles) {
+
+ ATH_MSG_DEBUG("Start vertex matching");
+
+ //setup decorators for truth matching info
+ static const xAOD::Vertex::Decorator<std::vector<VertexTruthMatchInfo> > matchInfoDecor("truthVertexMatchingInfos");
+ static const xAOD::Vertex::Decorator<int> recoMatchTypeDecor("vertexMatchType");
+ static const xAOD::Vertex::Decorator<std::vector<ElementLink<xAOD::VertexContainer> > > splitPartnerDecor("splitPartners");
+
+ const xAOD::Vertex::Decorator<float> fakeScoreDecor("fakeScore");
+ const xAOD::Vertex::Decorator<float> otherScoreDecor("otherScore");
+
+ //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");
+
+ ATH_MSG_DEBUG("Starting Loop on Vertices");
+
+ //=============================================================================
+ //First loop over vertices: get tracks, then TruthParticles, and store relative
+ //weights from each TruthVertex
+ //=============================================================================
+ for ( const xAOD::Vertex* vtx : recoVerticesToMatch ) {
+
+ //create the vector we will add as matching info decoration later
+ std::vector<VertexTruthMatchInfo> matchinfo;
+
+ const xAOD::Vertex::TrackParticleLinks_t & trkParts = trkAcc( *vtx );
+ size_t ntracks = trkParts.size();
+ const std::vector<float> & trkWeights = weightAcc( *vtx );
+
+ //if don't have track particles
+ if (!trkAcc.isAvailable(*vtx) || !weightAcc.isAvailable(*vtx) ) {
+ ATH_MSG_WARNING("trackParticles or trackWeights not available, vertex is missing info");
+ continue;
+ }
+ if ( trkWeights.size() != ntracks ) {
+ ATH_MSG_WARNING("Vertex without same number of tracks and trackWeights, vertex is missing info");
+ continue;
+ }
+
+ ATH_MSG_DEBUG("Matching new vertex at (" << vtx->x() << ", " << vtx->y() << ", " << vtx->z() << ")" << " with " << ntracks << " tracks, at index: " << vtx->index());
+
+ float totalWeight = 0.;
+ float totalPt = 0;
+ float otherPt = 0;
+ float fakePt = 0;
+
+ //loop over the tracks in the vertex
+ for ( size_t t = 0; t < ntracks; ++t ) {
+
+ ATH_MSG_DEBUG("Checking track number " << t);
+
+ if (!trkParts[t].isValid()) {
+ ATH_MSG_DEBUG("Track " << t << " is bad!");
+ continue;
+ }
+ const xAOD::TrackParticle & trk = **trkParts[t];
+
+ // store the contribution to total weight and pT
+ totalWeight += trkWeights[t];
+ totalPt += trk.pt();
+
+ // get the linked truth particle
+ if (!trk_truthPartAcc.isAvailable(trk)) {
+ ATH_MSG_DEBUG("The truth particle link decoration isn't available.");
+ continue;
+ }
+ const ElementLink<xAOD::TruthParticleContainer> & truthPartLink = trk_truthPartAcc( trk );
+ float prob = trk_truthProbAcc( trk );
+ ATH_MSG_DEBUG("Truth prob: " << prob);
+
+ // check the truth particle origin
+ if (truthPartLink.isValid() && prob > m_trkMatchProb) {
+ const xAOD::TruthParticle & truthPart = **truthPartLink;
+
+ int barcode = checkProduction(truthPart, truthVerticesToMatch);
+
+ //check if the truth particle is "good"
+ if ( barcode != -1 ) {
+ //track in vertex is linked to LLP descendant
+ //create link to truth vertex and add to matchInfo
+ auto it = std::find_if(truthVerticesToMatch.begin(), truthVerticesToMatch.end(),
+ [&](const auto& ele){ return ele->barcode() == barcode;} );
+
+ if(it == truthVerticesToMatch.end()) {
+ ATH_MSG_WARNING("Truth vertex with barcode " << barcode << " not found!");
+ }
+ else {
+ ElementLink<xAOD::TruthVertexContainer> elLink;
+ elLink.setElement(*it);
+ elLink.setStorableObject( *dynamic_cast<const xAOD::TruthVertexContainer*>( (*it)->container() ) );
+ size_t matchIdx = indexOfMatchInfo( matchinfo, elLink );
+
+ std::get<1>(matchinfo[matchIdx]) += trkWeights[t];
+ std::get<2>(matchinfo[matchIdx]) += trk.pt();
+ }
+ } else {
+ //truth particle failed cuts
+ ATH_MSG_DEBUG("Truth particle not from LLP decay.");
+ otherPt += trk.pt();
+ }
+ } else {
+ //not valid or low matching probability
+ ATH_MSG_DEBUG("Invalid or low prob truth link!");
+ fakePt += trk.pt();
+ }
+ }//end loop over tracks in vertex
+
+ // normalize by total weight and pT
+ std::for_each( matchinfo.begin(), matchinfo.end(), [&](VertexTruthMatchInfo& link)
+ {
+ std::get<1>(link) /= totalWeight;
+ std::get<2>(link) /= totalPt;
+ });
+
+ float fakeScore = fakePt/totalPt;
+ float otherScore = otherPt/totalPt;
+
+ matchInfoDecor ( *vtx ) = matchinfo;
+ fakeScoreDecor ( *vtx ) = fakeScore;
+ otherScoreDecor( *vtx ) = otherScore;
+ }
+
+ //After first loop, all vertices have been decorated with their vector of match info (link to TruthVertex 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( recoVerticesToMatch.size(), false );
+ static const xAOD::TruthVertex::Decorator<bool> isMatched("matchedToRecoVertex");
+ static const xAOD::TruthVertex::Decorator<bool> isSplit("vertexSplit");
+
+ for ( size_t i = 0; i < recoVerticesToMatch.size(); ++i ) {
+
+ int recoVertexMatchType = 0;
+
+ if ( assignedType[i] ) {
+ ATH_MSG_DEBUG("Vertex already defined as split.");
+ continue; // make sure we don't reclassify vertices already found in the split loop below
+ }
+
+ std::vector<VertexTruthMatchInfo> & info = matchInfoDecor( *recoVerticesToMatch[i] );
+ float fakeScore = fakeScoreDecor( *recoVerticesToMatch[i] );
+
+ if(fakeScore > m_vxMatchWeight) {
+ ATH_MSG_DEBUG("Vertex is fake.");
+ recoVertexMatchType = recoVertexMatchType | (0x1 << InDetSecVtxTruthMatchUtils::Fake);
+ } else if (info.size() == 1) {
+ if(std::get<2>(info[0]) > m_vxMatchWeight ) { // one truth matched vertex, sufficient weight
+ ATH_MSG_DEBUG("One true decay vertices matched with sufficient weight. Vertex is matched.");
+ recoVertexMatchType = recoVertexMatchType | (0x1 << InDetSecVtxTruthMatchUtils::Matched);
+ isMatched(**std::get<0>(info[0])) = true;
+ }
+ else {
+ ATH_MSG_DEBUG("One true decay vertices matched with insufficient weight. Vertex is other.");
+ recoVertexMatchType = recoVertexMatchType | (0x1 << InDetSecVtxTruthMatchUtils::Other);
+ }
+ } else if (info.size() >= 2 ) { // more than one true deacy vertices matched
+ ATH_MSG_DEBUG("Multiple true decay vertices matched. Vertex is merged.");
+ recoVertexMatchType = recoVertexMatchType | (0x1 << InDetSecVtxTruthMatchUtils::Merged);
+ std::for_each( info.begin(), info.end(), [](VertexTruthMatchInfo& link)
+ {
+ isMatched(**std::get<0>(link)) = true;
+ });
+ } else { // zero truth matched vertices, but not fake
+ ATH_MSG_DEBUG("Vertex is neither fake nor LLP. Marking as OTHER.");
+ recoVertexMatchType = recoVertexMatchType | (0x1 << InDetSecVtxTruthMatchUtils::Other);
+ }
+
+ recoMatchTypeDecor(*recoVerticesToMatch[i]) = recoVertexMatchType;
+
+ //check for splitting
+ if ( InDetSecVtxTruthMatchUtils::isMatched(recoMatchTypeDecor(*recoVerticesToMatch[i])) ||
+ InDetSecVtxTruthMatchUtils::isMerged(recoMatchTypeDecor(*recoVerticesToMatch[i])) ) {
+ std::vector<size_t> foundSplits;
+ for ( size_t j = i + 1; j < recoVerticesToMatch.size(); ++j ) {
+ std::vector<VertexTruthMatchInfo> & info2 = matchInfoDecor( *recoVerticesToMatch[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 (recoMatchTypeDecor( *recoVerticesToMatch[j] ) & (0x1 << InDetSecVtxTruthMatchUtils::Fake)) continue;
+ if (!info2.empty() && 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
+ ElementLink<xAOD::VertexContainer> splitLink_ij;
+ splitLink_ij.setElement( recoVerticesToMatch[j] );
+ splitLink_ij.setStorableObject( *dynamic_cast<const xAOD::VertexContainer*>(recoVerticesToMatch[j]->container()));
+ splitPartnerDecor( *recoVerticesToMatch[i] ).emplace_back(splitLink_ij);
+
+ ElementLink<xAOD::VertexContainer> splitLink_ji;
+ splitLink_ji.setElement( recoVerticesToMatch[i] );
+ splitLink_ji.setStorableObject( *dynamic_cast<const xAOD::VertexContainer*>(recoVerticesToMatch[i]->container()));
+ splitPartnerDecor( *recoVerticesToMatch[j] ).emplace_back(splitLink_ji);
+
+ //then between any others we found along the way
+ for ( auto k : foundSplits ) { //k is a size_t in the vector of splits
+ ElementLink<xAOD::VertexContainer> splitLink_kj;
+ splitLink_kj.setElement( recoVerticesToMatch[j] );
+ splitLink_kj.setStorableObject( *dynamic_cast<const xAOD::VertexContainer*>(recoVerticesToMatch[j]->container()));
+ splitPartnerDecor( *recoVerticesToMatch[k] ).emplace_back(splitLink_kj);
+
+ ElementLink<xAOD::VertexContainer> splitLink_jk;
+ splitLink_jk.setElement( recoVerticesToMatch[k] );
+ splitLink_jk.setStorableObject( *dynamic_cast<const xAOD::VertexContainer*>(recoVerticesToMatch[k]->container()));
+ splitPartnerDecor( *recoVerticesToMatch[j] ).emplace_back(splitLink_jk);
+ }
+ //then keep track that we found this one
+ foundSplits.push_back(j);
+ recoMatchTypeDecor( *recoVerticesToMatch[i] ) = recoMatchTypeDecor( *recoVerticesToMatch[i] ) | (0x1 << InDetSecVtxTruthMatchUtils::Split);
+ recoMatchTypeDecor( *recoVerticesToMatch[j] ) = recoMatchTypeDecor( *recoVerticesToMatch[j] ) | (0x1 << InDetSecVtxTruthMatchUtils::Split);
+ isSplit(**std::get<0>(info[0])) = true;
+ assignedType[j] = true;
+ } //if the two vertices match to same TruthVertex
+ }//inner loop over vertices
+ } //if matched or merged
+
+ } //outer loop
+
+ // now label truth vertices
+
+ ATH_MSG_DEBUG("Labeling truth vertices");
+
+ static const xAOD::TruthVertex::Decorator<int> truthMatchTypeDecor("truthVertexMatchType");
+
+ for(const xAOD::TruthVertex* truthVtx : truthVerticesToMatch) {
+
+ std::vector<const xAOD::TruthParticle*> reconstructibleParticles;
+ int counter = 0;
+ countReconstructibleDescendentParticles( *truthVtx, reconstructibleParticles, counter );
+
+ // hacky solution for keeping track of particles in the vertex
+ std::vector<int> particleInfo = {0,0,0};
+ std::vector<int> vertexInfo = {0,0,0};
+
+ for(size_t n = 0; n < reconstructibleParticles.size(); n++){
+ ATH_MSG_DEBUG("Checking daughter no. " << n);
+ const xAOD::TruthParticle* outPart = reconstructibleParticles.at(n);
+
+ if (trackParticles){
+ particleInfo = checkParticle(*outPart, trackParticles);
+
+ for(size_t h = 0; h < particleInfo.size(); h++){
+ vertexInfo.at(h) += particleInfo.at(h);
+ }
+ }
+ }
+
+ int truthMatchType = 0;
+ if( vertexInfo.at(0) > 1 && truthVtx->perp() < 320 && abs(truthVtx->z()) < 1500){
+ ATH_MSG_DEBUG("Vertex is reconstructable and in Inner Det region");
+ truthMatchType = truthMatchType | (0x1 << InDetSecVtxTruthMatchUtils::Reconstructable);
+ }
+ if( InDetSecVtxTruthMatchUtils::isReconstructable(truthMatchType) and vertexInfo.at(1) > 1){
+ ATH_MSG_DEBUG("Vertex has at least two tracks");
+ truthMatchType = truthMatchType | (0x1 << InDetSecVtxTruthMatchUtils::Accepted);
+ }
+ if(InDetSecVtxTruthMatchUtils::isAccepted(truthMatchType) and vertexInfo.at(2) > 1){
+ ATH_MSG_DEBUG("Vertex is has at least two tracks passing track selection: " << vertexInfo.at(2));
+ truthMatchType = truthMatchType | (0x1 << InDetSecVtxTruthMatchUtils::Seeded);
+ }
+ if(InDetSecVtxTruthMatchUtils::isSeeded(truthMatchType) and isMatched(*truthVtx)){
+ ATH_MSG_DEBUG("Vertex is matched to a reconstructed secVtx");
+ truthMatchType = truthMatchType | (0x1 << InDetSecVtxTruthMatchUtils::Reconstructed);
+ }
+ if(InDetSecVtxTruthMatchUtils::isSeeded(truthMatchType) and isSplit(*truthVtx)){
+ ATH_MSG_DEBUG("Vertex is matched to multiple secVtx");
+ truthMatchType = truthMatchType | (0x1 << InDetSecVtxTruthMatchUtils::ReconstructedSplit);
+ }
+ truthMatchTypeDecor(*truthVtx) = truthMatchType;
+ }
+ ATH_MSG_DEBUG("Done labeling truth vertices");
+
+ return StatusCode::SUCCESS;
+
+}
+
+std::vector<int> InDetSecVtxTruthMatchTool::checkParticle(const xAOD::TruthParticle &truthPart, const xAOD::TrackParticleContainer* trkCont) const {
+
+ xAOD::TrackParticle::ConstAccessor<char> trackPass(m_selectedTrackFlag);
+ xAOD::TrackParticle::ConstAccessor<ElementLink<xAOD::TruthParticleContainer> > trk_truthPartAcc("truthParticleLink");
+ xAOD::TrackParticle::ConstAccessor<float> trk_truthProbAcc("truthMatchProbability");
+
+ if(truthPart.pt() < m_trkPtCut){
+ ATH_MSG_DEBUG("Insufficient pt to reconstruct the particle");
+ return {0,0,0};
+ }
+ else{
+
+ for(const xAOD::TrackParticle* trkPart : *trkCont){
+ const ElementLink<xAOD::TruthParticleContainer> & truthPartLink = trk_truthPartAcc( *trkPart );
+ float matchProb = trk_truthProbAcc( *trkPart );
+
+ if(truthPartLink.isValid() && matchProb > m_trkMatchProb) {
+ const xAOD::TruthParticle& tmpPart = **truthPartLink;
+ if(tmpPart.barcode() == truthPart.barcode()) {
+ if(trackPass.isAvailable( *trkPart )) {
+ if(trackPass( *trkPart )) {
+ ATH_MSG_DEBUG("Particle has a track that passes track selection.");
+ return {1,1,1};
+ } else {
+ ATH_MSG_DEBUG("Particle has a track, but did not pass track selection.");
+ return {1,1,0};
+ }
+ } else {
+ ATH_MSG_DEBUG("Track selection decoration not available, calling the track selected");
+ return {1,1,1};
+ }
+ }
+ }
+ }
+ ATH_MSG_DEBUG("Particle has enough pt.");
+ return {1,0,0};
+
+ }
+ return {0,0,0};
+}
+
+// check if truth particle originated from decay of particle in the pdgIdList
+int InDetSecVtxTruthMatchTool::checkProduction( const xAOD::TruthParticle & truthPart, std::vector<const xAOD::TruthVertex*> truthVerticesToMatch ) const {
+
+ if (truthPart.nParents() == 0){
+ ATH_MSG_DEBUG("Particle has no parents (end of loop)");
+ return -1;
+ }
+ else{
+ const xAOD::TruthParticle * parent = truthPart.parent(0);
+ if(not parent) {
+ ATH_MSG_DEBUG("Particle parent is null");
+ return -1;
+ }
+ ATH_MSG_DEBUG("Parent ID: " << parent->pdgId());
+
+ const xAOD::TruthVertex* parentVertex = parent->decayVtx();
+ if(std::find(truthVerticesToMatch.begin(), truthVerticesToMatch.end(), parentVertex) != truthVerticesToMatch.end()) {
+ ATH_MSG_DEBUG("Found LLP decay.");
+ return parentVertex->barcode();
+ }
+ // recurse on parent
+ return checkProduction(*parent, truthVerticesToMatch);
+ }
+ return -1;
+}
+
+void InDetSecVtxTruthMatchTool::countReconstructibleDescendentParticles(const xAOD::TruthVertex& signalTruthVertex,
+ std::vector<const xAOD::TruthParticle*>& set, int counter) const {
+
+ counter++;
+
+ for( size_t itrk = 0; itrk < signalTruthVertex.nOutgoingParticles(); itrk++) {
+ const auto* particle = signalTruthVertex.outgoingParticle( itrk );
+ if( !particle ) continue;
+ // Recursively add descendents
+ if( particle->hasDecayVtx() ) {
+
+ TVector3 decayPos( particle->decayVtx()->x(), particle->decayVtx()->y(), particle->decayVtx()->z() );
+ TVector3 prodPos ( particle->prodVtx()->x(), particle->prodVtx()->y(), particle->prodVtx()->z() );
+
+ auto isInside = []( TVector3& v ) { return ( v.Perp() < 300. && std::abs( v.z() ) < 1500. ); };
+ auto isOutside = []( TVector3& v ) { return ( v.Perp() > 563. || std::abs( v.z() ) > 2720. ); };
+
+ const auto distance = (decayPos - prodPos).Mag();
+
+ if (counter > 100) {
+ ATH_MSG_WARNING("Vetoing particle that may be added recursively infinitely (potential loop in generator record");
+ break;
+ }
+
+ // consider track reconstructible if it travels at least 10mm
+ if( distance < 10.0 ) {
+ countReconstructibleDescendentParticles( *particle->decayVtx(), set , counter);
+ } else if( isInside ( prodPos ) && isOutside( decayPos ) && particle->isCharged() ) {
+ set.push_back( particle );
+ } else if( particle->isElectron() || particle->isMuon() ) {
+ set.push_back( particle );
+ }
+ } else {
+ if( !(particle->isCharged()) ) continue;
+ set.push_back( particle );
+ }
+ }
+
+ }
diff --git a/InnerDetector/InDetValidation/InDetSecVertexValidation/Root/LinkDef.h b/InnerDetector/InDetRecTools/InDetSecVtxTruthMatchTool/Root/LinkDef.h
similarity index 52%
rename from InnerDetector/InDetValidation/InDetSecVertexValidation/Root/LinkDef.h
rename to InnerDetector/InDetRecTools/InDetSecVtxTruthMatchTool/Root/LinkDef.h
index 0c34221b64a4f21f9b86b8090546dfd31f69d83d..7aac8d43f6a4ce79011eb00ee5ef3ac4e0a1a695 100644
--- a/InnerDetector/InDetValidation/InDetSecVertexValidation/Root/LinkDef.h
+++ b/InnerDetector/InDetRecTools/InDetSecVtxTruthMatchTool/Root/LinkDef.h
@@ -2,10 +2,10 @@
Copyright (C) 2002-2021 CERN for the benefit of the ATLAS collaboration
*/
-#ifndef INDETRUTHSECVERTEXVALIDATION_LINKDEF_H
-#define INDETRUTHSECVERTEXVALIDATION_LINKDEF_H
+#ifndef INDETRUTHSECVTXTRUTHMATCHTOOL_LINKDEF_H
+#define INDETRUTHSECVTXTRUTHMATCHTOOL_LINKDEF_H
-#include "InDetSecVertexValidation/InDetSecVertexTruthMatchTool.h"
+#include "InDetSecVtxTruthMatchTool/InDetSecVtxTruthMatchTool.h"
#ifdef __CINT__
@@ -14,7 +14,7 @@
#pragma link off all functions;
#pragma link C++ nestedclass;
-#pragma link C++ class InDetSecVertexTruthMatchTool;
+#pragma link C++ class InDetSecVtxTruthMatchTool;
#endif
diff --git a/InnerDetector/InDetRecTools/InDetSecVtxTruthMatchTool/src/components/InDetSecVtxTruthMatchTool_entries.cxx b/InnerDetector/InDetRecTools/InDetSecVtxTruthMatchTool/src/components/InDetSecVtxTruthMatchTool_entries.cxx
new file mode 100644
index 0000000000000000000000000000000000000000..b7e8fa3992a177ce2e1156fb690562ef608b7a12
--- /dev/null
+++ b/InnerDetector/InDetRecTools/InDetSecVtxTruthMatchTool/src/components/InDetSecVtxTruthMatchTool_entries.cxx
@@ -0,0 +1,11 @@
+/*
+ Copyright (C) 2002-2024 CERN for the benefit of the ATLAS collaboration
+*/
+
+// Gaudi/Athena include(s):
+//#include "GaudiKernel/DeclareFactoryEntries.h"
+
+// Local include(s):
+#include "InDetSecVtxTruthMatchTool/InDetSecVtxTruthMatchTool.h"
+
+DECLARE_COMPONENT( InDetSecVtxTruthMatchTool )
\ No newline at end of file
diff --git a/InnerDetector/InDetValidation/InDetSecVertexValidation/CMakeLists.txt b/InnerDetector/InDetValidation/InDetSecVertexValidation/CMakeLists.txt
index 4c66ee1e5a94fe6ef304642b63def6db115629b8..a8745449ecb0564b9d4bc6ff5368bc8459530640 100644
--- a/InnerDetector/InDetValidation/InDetSecVertexValidation/CMakeLists.txt
+++ b/InnerDetector/InDetValidation/InDetSecVertexValidation/CMakeLists.txt
@@ -14,30 +14,11 @@ endif()
# External dependencies:
find_package( ROOT COMPONENTS Core Tree Hist RIO )
-# Generate a CINT dictionary source file:
-atlas_add_root_dictionary( InDetSecVertexValidationLib _cintDictSource
- ROOT_HEADERS InDetSecVertexValidation/InDetSecVertexTruthMatchTool.h Root/LinkDef.h
- EXTERNAL_PACKAGES ROOT )
-
-# Component(s) in the package:
-atlas_add_library( InDetSecVertexValidationLib
- Root/*.cxx ${_cintDictSource}
- PUBLIC_HEADERS InDetSecVertexValidation
- PRIVATE_LINK_LIBRARIES ${extra_libs}
- LINK_LIBRARIES AsgTools xAODTracking xAODTruth )
-
if( NOT XAOD_STANDALONE )
atlas_add_component( InDetSecVertexValidation
src/*.cxx src/*.h
src/components/*.cxx
- LINK_LIBRARIES ${extra_libs} GaudiKernel xAODTracking AthenaBaseComps AsgTools InDetSecVertexValidationLib AthenaMonitoringLib TrkValHistUtils )
-endif()
-
-if( XAOD_STANDALONE )
- atlas_add_executable( SecVertexTruthMatchTest
- util/SecVertexTruthMatchTest.cxx
- INCLUDE_DIRS ${ROOT_INCLUDE_DIRS}
- LINK_LIBRARIES InDetSecVertexValidationLib xAODTracking xAODEventInfo xAODRootAccess AthenaMonitoringLib TrkValHistUtils ${ROOT_LIBRARIES} ${extra_libs} )
+ LINK_LIBRARIES ${extra_libs} GaudiKernel xAODTracking AthenaBaseComps AsgTools AthenaMonitoringLib TrkValHistUtils )
endif()
atlas_install_python_modules(python/*.py POST_BUILD_CMD ${ATLAS_FLAKE8})
diff --git a/InnerDetector/InDetValidation/InDetSecVertexValidation/InDetSecVertexValidation/ATLAS_CHECK_THREAD_SAFETY b/InnerDetector/InDetValidation/InDetSecVertexValidation/InDetSecVertexValidation/ATLAS_CHECK_THREAD_SAFETY
deleted file mode 100644
index c4b6ab6470835418aab6f71fd8e2c2a842cef9e8..0000000000000000000000000000000000000000
--- a/InnerDetector/InDetValidation/InDetSecVertexValidation/InDetSecVertexValidation/ATLAS_CHECK_THREAD_SAFETY
+++ /dev/null
@@ -1 +0,0 @@
-InnerDetector/InDetValidation/InDetSecVertexValidation
diff --git a/InnerDetector/InDetValidation/InDetSecVertexValidation/InDetSecVertexValidation/InDetSecVertexTruthMatchTool.h b/InnerDetector/InDetValidation/InDetSecVertexValidation/InDetSecVertexValidation/InDetSecVertexTruthMatchTool.h
deleted file mode 100644
index b26bbce5ea206c57ccc3961bde8f3946eebd4c21..0000000000000000000000000000000000000000
--- a/InnerDetector/InDetValidation/InDetSecVertexValidation/InDetSecVertexValidation/InDetSecVertexTruthMatchTool.h
+++ /dev/null
@@ -1,200 +0,0 @@
-/*
- Copyright (C) 2002-2022 CERN for the benefit of the ATLAS collaboration
-*/
-
-#ifndef InDetSecVertexTruthMatchTool_h
-#define InDetSecVertexTruthMatchTool_h
-
-// Framework include(s):
-#include "AsgTools/AsgTool.h"
-#include "GaudiKernel/ITHistSvc.h"
-
-
-
-// EDM include(s):
-#include "xAODTracking/VertexContainerFwd.h"
-#include "xAODTruth/TruthParticleFwd.h"
-#include "xAODTruth/TruthVertexFwd.h"
-#include "xAODTruth/TruthVertexContainer.h"
-#include "InDetSecVertexValidation/IInDetSecVertexTruthMatchTool.h"
-
-// standard includes
-#include<vector>
-#include <TH1.h>
-#include <TH2.h>
-#include <TH3.h>
-
-
-
-namespace InDetSecVertexTruthMatchUtils {
-
- //Namespace for useful analysis things on the truth matching decorations applied to the VertexContainer
-
- //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::TruthVertexContainer>, 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")
- LLP,
- OTHER,
- ALL,
- NTYPES // access to number of types
- };
-
- constexpr std::initializer_list<VertexMatchType> allTypes = { MATCHED, MERGED, SPLIT, FAKE, LLP, OTHER, ALL };
-
-
- //type codes for truth vertices
- //NOTE: types are subsets of subsequent types
- enum TruthVertexMatchType {
- INCLUSIVE, // all vertices
- RECONSTRUCTABLE, // fiducial cuts, >= 2 charged daughters
- ACCEPTED, // >= 2 reco tracks
- SEEDED, // tracks pass cuts
- RECONSTRUCTED, // matched to reco vtx
- RECONSTRUCTEDSPLIT, // matched to multiple vtx
- };
-
-}
-
-/** Class for vertex truth matching.
- * Match reconstructed vertices to truth level decay vertices
- * Categorize reconstructed vertices depending on their composition.
- */
-class InDetSecVertexTruthMatchTool : public virtual IInDetSecVertexTruthMatchTool,
- public asg::AsgTool {
-
- ASG_TOOL_CLASS( InDetSecVertexTruthMatchTool, IInDetSecVertexTruthMatchTool )
-
- public:
-
- InDetSecVertexTruthMatchTool( const std::string & name );
-
- virtual StatusCode initialize() override final;
- virtual StatusCode finalize() override;
-
- //take const collection of vertices, match them, and decorate with matching info
- virtual StatusCode matchVertices( const xAOD::VertexContainer& vtxContainer, const xAOD::TruthVertexContainer& truthVtxContainer ) override;
- //take const collection of truth vertices and decorate with type info
- virtual StatusCode labelTruthVertices( const xAOD::TruthVertexContainer & truthVtxContainer ) override;
-
- private:
-
- //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;
- //comma separated list of pdgids to consider in the matchiing
- std::string m_pdgIds;
- //turn on/off histogram output
- bool m_fillHist;
- //Augmentation string to add to the end of patterns.
- std::string m_AugString;
-
- //private methods to check if particles are good to use
- //returns barcode of LLP production truth vertex
- int checkProduction( const xAOD::TruthParticle& truthPart ) const;
- void countReconstructibleDescendentParticles(const xAOD::TruthVertex& signalTruthVertex,
- std::vector<const xAOD::TruthParticle*>& set, int counter) const;
- std::vector<int> checkParticle( const xAOD::TruthParticle& part, const xAOD::TrackParticleContainer &tkCont ) const;
-
- // (optional) write out histograms
- StatusCode fillRecoPlots( const xAOD::Vertex& secVtx );
- StatusCode fillTruthPlots( const xAOD::TruthVertex& truthVtx );
-
- //private internal variables (not properties)
- std::vector<int> m_pdgIdList;
-
- //histograms
- ITHistSvc* m_thistSvc{nullptr};
- TH1F* m_matchType = nullptr;
- TH1F* m_truth_Ntrk = nullptr;
-
- TH1F* m_truthInclusive_r = nullptr;
- TH1F* m_truthReconstructable_r = nullptr;
- TH1F* m_truthAccepted_r = nullptr;
- TH1F* m_truthSeeded_r = nullptr;
- TH1F* m_truthReconstructed_r = nullptr;
- TH1F* m_truthSplit_r = nullptr;
-
- TH1F* m_truthReconstructable_trkSel = nullptr;
- TH1F* m_truthReconstructed_trkSel = nullptr;
-
- std::map<std::string,TH1F*> m_positionRes_R;
- std::map<std::string,TH1F*> m_positionRes_Z;
- std::map<std::string,TH1F*> m_recoX;
- std::map<std::string,TH1F*> m_recoY;
- std::map<std::string,TH1F*> m_recoZ;
- std::map<std::string,TH1F*> m_recoR;
- std::map<std::string,TH1F*> m_recodistFromPV;
- std::map<std::string,TH1F*> m_recoPt;
- std::map<std::string,TH1F*> m_recoEta;
- std::map<std::string,TH1F*> m_recoPhi;
- std::map<std::string,TH1F*> m_recoMass;
- std::map<std::string,TH1F*> m_recoMu;
- std::map<std::string,TH1F*> m_recoChi2;
- std::map<std::string,TH1F*> m_recoDir;
- std::map<std::string,TH1F*> m_recoCharge;
- std::map<std::string,TH1F*> m_recoH;
- std::map<std::string,TH1F*> m_recoHt;
- std::map<std::string,TH1F*> m_recoMinOpAng;
- std::map<std::string,TH1F*> m_recoMaxOpAng;
- std::map<std::string,TH1F*> m_recoMaxDR;
- std::map<std::string,TH1F*> m_recoMinD0;
- std::map<std::string,TH1F*> m_recoMaxD0;
- std::map<std::string,TH1F*> m_recoNtrk;
-
- std::map<std::string,TH1F*> m_recoTrk_qOverP;
- std::map<std::string,TH1F*> m_recoTrk_theta;
- std::map<std::string,TH1F*> m_recoTrk_E;
- std::map<std::string,TH1F*> m_recoTrk_M;
- std::map<std::string,TH1F*> m_recoTrk_Pt;
- std::map<std::string,TH1F*> m_recoTrk_Px;
- std::map<std::string,TH1F*> m_recoTrk_Py;
- std::map<std::string,TH1F*> m_recoTrk_Pz;
- std::map<std::string,TH1F*> m_recoTrk_Eta;
- std::map<std::string,TH1F*> m_recoTrk_Phi;
- std::map<std::string,TH1F*> m_recoTrk_D0;
- std::map<std::string,TH1F*> m_recoTrk_Z0;
- std::map<std::string,TH1F*> m_recoTrk_errD0;
- std::map<std::string,TH1F*> m_recoTrk_errZ0;
- std::map<std::string,TH1F*> m_recoTrk_Chi2;
- std::map<std::string,TH1F*> m_recoTrk_nDoF;
- std::map<std::string,TH1F*> m_recoTrk_charge;
-
- std::map<std::string,TH1F*> m_matchScore_weight;
- std::map<std::string,TH1F*> m_matchScore_pt;
- std::map<std::string,TH1F*> m_matchedTruthID;
-
- TH1F* m_truthX = nullptr;
- TH1F* m_truthY = nullptr;
- TH1F* m_truthZ = nullptr;
- TH1F* m_truthR = nullptr;
- TH1F* m_truthdistFromPV = nullptr;
- TH1F* m_truthEta = nullptr;
- TH1F* m_truthPhi = nullptr;
- TH1F* m_truthNtrk_out = nullptr;
- TH1F* m_truthParent_E = nullptr;
- TH1F* m_truthParent_M = nullptr;
- TH1F* m_truthParent_Pt = nullptr;
- TH1F* m_truthParent_Eta = nullptr;
- TH1F* m_truthParent_Phi = nullptr;
- TH1F* m_truthParent_charge = nullptr;
- TH1F* m_truthParentProdX = nullptr;
- TH1F* m_truthParentProdY = nullptr;
- TH1F* m_truthParentProdZ = nullptr;
- TH1F* m_truthParentProdR = nullptr;
- TH1F* m_truthParentProddistFromPV = nullptr;
- TH1F* m_truthParentProdEta = nullptr;
- TH1F* m_truthParentProdPhi = nullptr;
-
-
-};
-
-#endif
diff --git a/InnerDetector/InDetValidation/InDetSecVertexValidation/Root/InDetSecVertexTruthMatchTool.cxx b/InnerDetector/InDetValidation/InDetSecVertexValidation/Root/InDetSecVertexTruthMatchTool.cxx
deleted file mode 100644
index 576c83d92663f240d9402a5f9506bcf8cd32072c..0000000000000000000000000000000000000000
--- a/InnerDetector/InDetValidation/InDetSecVertexValidation/Root/InDetSecVertexTruthMatchTool.cxx
+++ /dev/null
@@ -1,964 +0,0 @@
-/*
- Copyright (C) 2002-2022 CERN for the benefit of the ATLAS collaboration
-*/
-#include "AthenaKernel/Units.h"
-
-#include "InDetSecVertexValidation/InDetSecVertexTruthMatchTool.h"
-
-#include "xAODTracking/TrackParticleContainer.h"
-#include "xAODTruth/TruthEventContainer.h"
-
-using namespace InDetSecVertexTruthMatchUtils;
-using Athena::Units::GeV;
-
-InDetSecVertexTruthMatchTool::InDetSecVertexTruthMatchTool( const std::string & name ) : asg::AsgTool(name) {
- declareProperty("trackMatchProb", m_trkMatchProb = 0.5 );
- declareProperty("vertexMatchWeight", m_vxMatchWeight = 0.5 );
- declareProperty("trackPtCut", m_trkPtCut = 1000. );
- declareProperty("pdgIds", m_pdgIds = "36" );
- declareProperty("fillHist", m_fillHist = true );
- declareProperty("AugString", m_AugString = "" );
-}
-
-StatusCode InDetSecVertexTruthMatchTool::initialize() {
- ATH_MSG_INFO("Initializing InDetSecVertexTruthMatchTool");
-
- // build the vector of pdgIds from the input string
- std::stringstream ss(m_pdgIds);
- for (int i; ss >> i;) {
- m_pdgIdList.push_back(i);
- if (ss.peek() == ',')
- ss.ignore();
- }
-
- // histograms
- ATH_CHECK( service("THistSvc",m_thistSvc) );
-
- ////////////////////////////////////////////
- ////// Seconvery Vertex Histograms /////////
- ////////////////////////////////////////////
- m_matchType = new TH1F("reco_matchType","Vertex Match Type",8,-0.5,7.5);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/SecondaryVertex/matchType", m_matchType));
-
- for(auto type : allTypes) {
- std::string vertexType = "";
- switch (type)
- {
- case MATCHED:
- vertexType = "Matched";
- break;
- case MERGED:
- vertexType = "Merged";
- break;
- case SPLIT:
- vertexType = "Split";
- break;
- case LLP:
- vertexType = "LLP";
- break;
- case OTHER:
- vertexType = "Other";
- break;
- case FAKE:
- vertexType = "Fake";
- break;
- case ALL:
- vertexType = "All";
- break;
- default:
- break;
- }
-
-
- m_recoX[vertexType] = new TH1F(("recoZ" + vertexType).c_str(),"Reco vertex x [mm]",500,0,500);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/SecondaryVertex/" + vertexType + "/X", m_recoX[vertexType]));
- m_recoY[vertexType] = new TH1F(("recoY" + vertexType).c_str(),"Reco vertex y [mm]",500,0,500);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/SecondaryVertex/" + vertexType + "/Y", m_recoY[vertexType]));
- m_recoZ[vertexType] = new TH1F(("recoZ" + vertexType).c_str(),"Reco vertex z [mm]",500,0,500);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/SecondaryVertex/" + vertexType + "/Z", m_recoZ[vertexType]));
- m_recoR[vertexType] = new TH1F(("recoR" + vertexType).c_str(),"Reco vertex r [mm]",500,0,500);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/SecondaryVertex/" + vertexType + "/R", m_recoR[vertexType]));
- m_recoPt[vertexType] = new TH1F(("recoPt" + vertexType).c_str(),"Reco vertex Pt [GeV]",100,0,100);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/SecondaryVertex/" + vertexType + "/Pt", m_recoPt[vertexType]));
- m_recoEta[vertexType] = new TH1F(("recoEta" + vertexType).c_str(),"Reco vertex Eta ",100,-5,5);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/SecondaryVertex/" + vertexType + "/Eta", m_recoEta[vertexType]));
- m_recoPhi[vertexType] = new TH1F(("recoPhi" + vertexType).c_str(),"Reco vertex Phi ",64,-3.2,3.2);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/SecondaryVertex/" + vertexType + "/Phi", m_recoPhi[vertexType]));
- m_recoMass[vertexType] = new TH1F(("recoMass" + vertexType).c_str(),"Reco vertex Mass [GeV]",500,0,100);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/SecondaryVertex/" + vertexType + "/Mass", m_recoMass[vertexType]));
- m_recoMu[vertexType] = new TH1F(("recoMu" + vertexType).c_str(),"Reco vertex Red. Mass [GeV]",500,0,100);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/SecondaryVertex/" + vertexType + "/Mu", m_recoMu[vertexType]));
- m_recoChi2[vertexType] = new TH1F(("recoChi2" + vertexType).c_str(),"Reco vertex recoChi2",100,0,10);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/SecondaryVertex/" + vertexType + "/Chi2", m_recoChi2[vertexType]));
- m_recoDir[vertexType] = new TH1F(("recoDir" + vertexType).c_str(),"Reco vertex recoDirection",100,-1,1);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/SecondaryVertex/" + vertexType + "/Dir", m_recoDir[vertexType]));
- m_recoCharge[vertexType] = new TH1F(("recoCharge" + vertexType).c_str(),"Reco vertex recoCharge",20,-10,10);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/SecondaryVertex/" + vertexType + "/Charge", m_recoCharge[vertexType]));
- m_recoH[vertexType] = new TH1F(("recoH" + vertexType).c_str(),"Reco vertex H [GeV]",100,0,100);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/SecondaryVertex/" + vertexType + "/H", m_recoH[vertexType]));
- m_recoHt[vertexType] = new TH1F(("recoHt" + vertexType).c_str(),"Reco vertex Mass [GeV]",100,0,100);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/SecondaryVertex/" + vertexType + "/Ht", m_recoHt[vertexType]));
- m_recoMinOpAng[vertexType] = new TH1F(("recoMinOpAng" + vertexType).c_str(),"Reco vertex minOpAng",100,-1,1);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/SecondaryVertex/" + vertexType + "/MinOpAng", m_recoMinOpAng[vertexType]));
- m_recoMaxOpAng[vertexType] = new TH1F(("recoMaxOpAng" + vertexType).c_str(),"Reco vertex MaxOpAng",100,-1,1);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/SecondaryVertex/" + vertexType + "/MaxOpAng", m_recoMaxOpAng[vertexType]));
- m_recoMaxDR[vertexType] = new TH1F(("recoMaxDR" + vertexType).c_str(),"Reco vertex maxDR",100,0,10);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/SecondaryVertex/" + vertexType + "/MaxDR", m_recoMaxDR[vertexType]));
- m_recoMinD0[vertexType] = new TH1F(("recoMinD0" + vertexType).c_str(),"Reco vertex min d0 [mm]",100,0,100);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/SecondaryVertex/" + vertexType + "/MinD0", m_recoMinD0[vertexType]));
- m_recoMaxD0[vertexType] = new TH1F(("recoMaxD0" + vertexType).c_str(),"Reco vertex max d0 [mm]",100,0,100);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/SecondaryVertex/" + vertexType + "/MaxD0", m_recoMaxD0[vertexType]));
- m_recoNtrk[vertexType] = new TH1F(("recoNtrk" + vertexType).c_str(),"Reco vertex n tracks",30,0,30);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/SecondaryVertex/" + vertexType + "/Ntrk", m_recoNtrk[vertexType]));
-
- m_positionRes_R[vertexType] = new TH1F(("positionRes_R" + vertexType).c_str(),"Position resolution for vertices matched to truth decays",400,-20,20);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/SecondaryVertex/" + vertexType + "/positionRes_R", m_positionRes_R[vertexType]));
- m_positionRes_Z[vertexType] = new TH1F(("positionRes_Z" + vertexType).c_str(),"Position resolution for vertices matched to truth decays",400,-20,20);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/SecondaryVertex/" + vertexType + "/positionRes_Z", m_positionRes_Z[vertexType]));
- m_matchScore_weight[vertexType] = new TH1F(("matchScore_weight" + vertexType).c_str(),"Vertex Match Score (weight)",101,0,1.01);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/SecondaryVertex/" + vertexType + "/matchScore_weight", m_matchScore_weight[vertexType]));
- m_matchScore_pt[vertexType] = new TH1F(("matchScore_pt" + vertexType).c_str(),"Vertex Match Score (pT)",101,0,1.01);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/SecondaryVertex/" + vertexType + "/matchScore_pt", m_matchScore_pt[vertexType]));
- m_matchedTruthID[vertexType] = new TH1F(("matchedTruthID" + vertexType).c_str(),"Vertex Truth Match ID",100,0,100);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/SecondaryVertex/" + vertexType + "/matchedTruthID", m_matchedTruthID[vertexType]));
-
- // tracks
- m_recoTrk_qOverP[vertexType] = new TH1F(("recoTrk_qOverP" + vertexType).c_str(),"Reco track qOverP ",100,0,.01);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/SecondaryVertex/" + vertexType + "/Trk_qOverP", m_recoTrk_qOverP[vertexType]));
- m_recoTrk_theta[vertexType] = new TH1F(("recoTrk_theta" + vertexType).c_str(),"Reco track theta ",64,0,3.2);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/SecondaryVertex/" + vertexType + "/Trk_theta", m_recoTrk_theta[vertexType]));
- m_recoTrk_E[vertexType] = new TH1F(("recoTrk_E" + vertexType).c_str(),"Reco track E ",100,0,100);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/SecondaryVertex/" + vertexType + "/Trk_E", m_recoTrk_E[vertexType]));
- m_recoTrk_M[vertexType] = new TH1F(("recoTrk_M" + vertexType).c_str(),"Reco track M ",100,0,10);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/SecondaryVertex/" + vertexType + "/Trk_M", m_recoTrk_M[vertexType]));
- m_recoTrk_Pt[vertexType] = new TH1F(("recoTrk_Pt" + vertexType).c_str(),"Reco track Pt ",100,0,100);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/SecondaryVertex/" + vertexType + "/Trk_Pt", m_recoTrk_Pt[vertexType]));
- m_recoTrk_Px[vertexType] = new TH1F(("recoTrk_Px" + vertexType).c_str(),"Reco track Px ",100,0,100);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/SecondaryVertex/" + vertexType + "/Trk_Px", m_recoTrk_Px[vertexType]));
- m_recoTrk_Py[vertexType] = new TH1F(("recoTrk_Py" + vertexType).c_str(),"Reco track Py ",100,0,100);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/SecondaryVertex/" + vertexType + "/Trk_Py", m_recoTrk_Py[vertexType]));
- m_recoTrk_Pz[vertexType] = new TH1F(("recoTrk_Pz" + vertexType).c_str(),"Reco track Pz ",100,0,100);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/SecondaryVertex/" + vertexType + "/Trk_Pz", m_recoTrk_Pz[vertexType]));
- m_recoTrk_Eta[vertexType] = new TH1F(("recoTrk_Eta" + vertexType).c_str(),"Reco track Eta ",100,-5,5);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/SecondaryVertex/" + vertexType + "/Trk_Eta", m_recoTrk_Eta[vertexType]));
- m_recoTrk_Phi[vertexType] = new TH1F(("recoTrk_Phi" + vertexType).c_str(),"Reco track Phi ",63,-3.2,3.2);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/SecondaryVertex/" + vertexType + "/Trk_Phi", m_recoTrk_Phi[vertexType]));
- m_recoTrk_D0[vertexType] = new TH1F(("recoTrk_D0" + vertexType).c_str(),"Reco track D0 ",300,-300,300);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/SecondaryVertex/" + vertexType + "/Trk_D0", m_recoTrk_D0[vertexType]));
- m_recoTrk_Z0[vertexType] = new TH1F(("recoTrk_Z0" + vertexType).c_str(),"Reco track Z0 ",500,-500,500);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/SecondaryVertex/" + vertexType + "/Trk_Z0", m_recoTrk_Z0[vertexType]));
- m_recoTrk_errD0[vertexType] = new TH1F(("recoTrk_errD0" + vertexType).c_str(),"Reco track errD0 ",300,0,30);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/SecondaryVertex/" + vertexType + "/Trk_errD0", m_recoTrk_errD0[vertexType]));
- m_recoTrk_errZ0[vertexType] = new TH1F(("recoTrk_errZ0" + vertexType).c_str(),"Reco track errZ0 ",500,0,50);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/SecondaryVertex/" + vertexType + "/Trk_errZ0", m_recoTrk_errZ0[vertexType]));
- m_recoTrk_Chi2[vertexType] = new TH1F(("recoTrk_Chi2" + vertexType).c_str(),"Reco track Chi2 ",100,0,10);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/SecondaryVertex/" + vertexType + "/Trk_Chi2", m_recoTrk_Chi2[vertexType]));
- m_recoTrk_nDoF[vertexType] = new TH1F(("recoTrk_nDoF" + vertexType).c_str(),"Reco track nDoF ",100,0,100);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/SecondaryVertex/" + vertexType + "/Trk_nDoF", m_recoTrk_nDoF[vertexType]));
- m_recoTrk_charge[vertexType] = new TH1F(("recoTrk_charge" + vertexType).c_str(),"Reco track charge ",3,-1.5,1.5);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/SecondaryVertex/" + vertexType + "/Trk_charge", m_recoTrk_charge[vertexType]));
-
- }
-
- ////////////////////////////////////////////
- //////// Truth Vertex Histograms ///////////
- ////////////////////////////////////////////
- m_truthX = new TH1F("truth_X","truth vertex x [mm]",500,-500,500);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/TruthVertex/truthX", m_truthX));
- m_truthY = new TH1F("truth_Y","truth vertex y [mm]",500,-500,500);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/TruthVertex/truthY", m_truthY));
- m_truthZ = new TH1F("truth_Z","truth vertex z [mm]",500,-500,500);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/TruthVertex/truthZ", m_truthZ));
- m_truthR = new TH1F("truth_R","truth vertex r [mm]",6000,0,600);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/TruthVertex/truthR", m_truthR));
- m_truthdistFromPV = new TH1F("truth_distFromPV","truth vertex distFromPV [mm]",500,0,500);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/TruthVertex/truthdistFromPV", m_truthdistFromPV));
- m_truthEta = new TH1F("truth_Eta","truth veEtatex Eta ",100,-5,5);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/TruthVertex/truthEta", m_truthEta));
- m_truthPhi = new TH1F("truth_Phi","truth vePhitex Phi ",64,-3.2,3.2);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/TruthVertex/truthPhi", m_truthPhi));
- m_truthNtrk_out = new TH1F("truth_Ntrk_out","truth vertex n outgoing tracks",100,0,100);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/TruthVertex/truthNtrk_out", m_truthNtrk_out));
- m_truthParent_E = new TH1F("truth_Parent_E","Reco track E ",100,0,100);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/TruthVertex/truthParent_E", m_truthParent_E));
- m_truthParent_M = new TH1F("truth_Parent_M","Reco track M ",500,0,500);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/TruthVertex/truthParent_M", m_truthParent_M));
- m_truthParent_Pt = new TH1F("truth_Parent_Pt","Reco track Pt ",100,0,100);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/TruthVertex/truthParent_Pt", m_truthParent_Pt));
- m_truthParent_Eta = new TH1F("truth_Parent_Eta","Reco track Eta ",100,-5,5);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/TruthVertex/truthParent_Eta", m_truthParent_Eta));
- m_truthParent_Phi = new TH1F("truth_Parent_Phi","Reco track Phi ",63,-3.2,3.2);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/TruthVertex/truthParent_Phi", m_truthParent_Phi));
- m_truthParent_charge = new TH1F("truth_Parent_charge","Reco track charge ",3,-1,1);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/TruthVertex/truthParent_charge", m_truthParent_charge));
- m_truthParentProdX = new TH1F("truth_ParentProdX","truthParentProd vertex x [mm]",500,-500,500);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/TruthVertex/truthParentProdX", m_truthParentProdX));
- m_truthParentProdY = new TH1F("truth_ParentProdY","truthParentProd vertex y [mm]",500,-500,500);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/TruthVertex/truthParentProdY", m_truthParentProdY));
- m_truthParentProdZ = new TH1F("truth_ParentProdZ","truthParentProd vertex z [mm]",500,-500,500);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/TruthVertex/truthParentProdZ", m_truthParentProdZ));
- m_truthParentProdR = new TH1F("truth_ParentProdR","truthParentProd vertex r [mm]",6000,0,600);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/TruthVertex/truthParentProdR", m_truthParentProdR));
- m_truthParentProddistFromPV = new TH1F("truth_ParentProddistFromPV","truthParentProd vertex distFromPV [mm]",500,0,500);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/TruthVertex/truthParentProddistFromPV", m_truthParentProddistFromPV));
- m_truthParentProdEta = new TH1F("truth_ParentProdEta","truthParentProd veEtatex Eta ",100,-5,5);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/TruthVertex/truthParentProdEta", m_truthParentProdEta));
- m_truthParentProdPhi = new TH1F("truth_ParentProdPhi","truthParentProd vePhitex Phi ",64,-3.2,3.2);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/TruthVertex/truthParentProdPhi", m_truthParentProdPhi));
-
- m_truthInclusive_r = new TH1F("truth_R_Inclusive","Reconstructable Truth Vertices",6000,0,600);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/TruthVertex/Inclusive/R", m_truthInclusive_r));
- m_truthReconstructable_r = new TH1F("truth_R_Reconstructable","Truth Vertices in detector acceptance",6000,0,600);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/TruthVertex/Reconstructable/R", m_truthReconstructable_r));
- m_truthAccepted_r = new TH1F("truth_R_Accepted","Truth Vertices in detector acceptance",6000,0,600);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/TruthVertex/Accepted/R", m_truthAccepted_r));
- m_truthSeeded_r = new TH1F("truth_R_Seeded","Seedable Truth Vertices",6000,0,600);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/TruthVertex/Seeded/R", m_truthSeeded_r));
- m_truthReconstructed_r = new TH1F("truth_R_Reconstructed","Vertex with Match Score > 0.5",6000,0,600);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/TruthVertex/Reconstructed/R", m_truthReconstructed_r));
- m_truthSplit_r = new TH1F("truth_R_Split","Split Vertex",6000,0,600);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/TruthVertex/Split/R", m_truthSplit_r));
-
- // TODO: implement these plots
- m_truth_Ntrk = new TH1F("truth_Ntrk","Truth vertex n track pass tracks",30,0,30);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/TruthVertex/Inclusive/Ntrk", m_truth_Ntrk));
- m_truthReconstructable_trkSel = new TH1F("truth_Ntrk_Seeded","Seedable Truth Vertices",30,0,30);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/TruthVertex/Reconstructable/Ntrk", m_truthReconstructable_trkSel));
- m_truthReconstructed_trkSel = new TH1F("truth_Ntrk_Reconstructed","Vertex with Match Score > 0.5",30,0,30);
- ATH_CHECK( m_thistSvc->regHist("/VTXPLOTS/" + name() + "/TruthVertex/Reconstructed/Ntrk", m_truthReconstructed_trkSel));
- return StatusCode::SUCCESS;
-}
-
-StatusCode InDetSecVertexTruthMatchTool::finalize()
-{
-
- ATH_MSG_INFO("Finalizing InDetSecVertexTruthMatchTool");
- return StatusCode::SUCCESS;
-}
-
-namespace {
-//Helper methods for this file only
-
-//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, const ElementLink<xAOD::TruthVertexContainer> & 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.emplace_back( link, 0., 0. );
- return matches.size() - 1;
-}
-
-}
-
-StatusCode InDetSecVertexTruthMatchTool::matchVertices( const xAOD::VertexContainer & vtxContainer,
- const xAOD::TruthVertexContainer & truthVtxContainer) {
-
- ATH_MSG_DEBUG("Start vertex matching");
-
- //setup decorators for truth matching info
- static const xAOD::Vertex::Decorator<std::vector<VertexTruthMatchInfo> > matchInfoDecor("TruthVertexMatchingInfos");
- static const xAOD::Vertex::Decorator<VertexMatchType> matchTypeDecor("VertexMatchType");
- static const xAOD::Vertex::Decorator<std::vector<ElementLink<xAOD::VertexContainer> > > splitPartnerDecor("SplitPartners");
- static const xAOD::Vertex::Decorator<ElementLink<xAOD::TruthVertexContainer> > backLinkDecor("RecoToTruthLink");
-
- const xAOD::Vertex::Decorator<float> fakeScoreDecor("FakeScore");
- const xAOD::Vertex::Decorator<float> otherScoreDecor("OtherScore");
-
- //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");
-
- //some variables to store
- size_t ntracks;
- xAOD::VxType::VertexType vtxType;
-
- ATH_MSG_DEBUG("Starting Loop on Vertices");
-
- //=============================================================================
- //First loop over vertices: get tracks, then TruthParticles, and store relative
- //weights from each TruthVertex
- //=============================================================================
- size_t vtxEntry = 0;
-
- for ( const xAOD::Vertex* vtx : vtxContainer ) {
-
- vtxEntry++;
- ATH_MSG_DEBUG("Matching vertex number: " << vtxEntry << ".");
-
- vtxType = static_cast<xAOD::VxType::VertexType>( vtx->vertexType() );
-
- if(vtxType != xAOD::VxType::SecVtx ){
- ATH_MSG_DEBUG("Vertex not labeled as secondary");
- continue;
- }
-
- //create the vector we will add as matching info decoration later
- std::vector<VertexTruthMatchInfo> matchinfo;
-
- const xAOD::Vertex::TrackParticleLinks_t & trkParts = trkAcc( *vtx );
- ntracks = trkParts.size();
- const std::vector<float> & trkWeights = weightAcc( *vtx );
-
- //if don't have track particles
- if (!trkAcc.isAvailable(*vtx) || !weightAcc.isAvailable(*vtx) ) {
- ATH_MSG_WARNING("trackParticles or trackWeights not available, vertex is missing info");
- continue;
- }
- if ( trkWeights.size() != ntracks ) {
- ATH_MSG_WARNING("Vertex without same number of tracks and trackWeights, vertex is missing info");
- continue;
- }
-
- ATH_MSG_DEBUG("Matching new vertex at (" << vtx->x() << ", " << vtx->y() << ", " << vtx->z() << ")" << " with " << ntracks << " tracks, at index: " << vtx->index());
-
- float totalWeight = 0.;
- float totalPt = 0;
- float otherPt = 0;
- float fakePt = 0;
-
- //loop over the tracks in the vertex
- for ( size_t t = 0; t < ntracks; ++t ) {
-
- ATH_MSG_DEBUG("Checking track number " << t);
-
- if (!trkParts[t].isValid()) {
- ATH_MSG_DEBUG("Track " << t << " is bad!");
- continue;
- }
- const xAOD::TrackParticle & trk = **trkParts[t];
-
- // store the contribution to total weight and pT
- totalWeight += trkWeights[t];
- totalPt += trk.pt();
-
- // get the linked truth particle
- if (!trk_truthPartAcc.isAvailable(trk)) {
- ATH_MSG_DEBUG("The truth particle link decoration isn't available.");
- continue;
- }
- const ElementLink<xAOD::TruthParticleContainer> & truthPartLink = trk_truthPartAcc( trk );
- float prob = trk_truthProbAcc( trk );
- ATH_MSG_DEBUG("Truth prob: " << prob);
-
- // check the truth particle origin
- if (truthPartLink.isValid() && prob > m_trkMatchProb) {
- const xAOD::TruthParticle & truthPart = **truthPartLink;
-
- int barcode = checkProduction(truthPart);
-
- //check if the truth particle is "good"
- if ( barcode != -1 ) {
- //track in vertex is linked to LLP descendant
- //create link to truth vertex and add to matchInfo
- auto it = std::find_if(truthVtxContainer.begin(), truthVtxContainer.end(),
- [&](const auto& ele){ return ele->barcode() == barcode;} );
-
- if(it == truthVtxContainer.end()) {
- ATH_MSG_DEBUG("Truth vertex with barcode " << barcode << " not found!");
- }
-
- // get index for the elementLink
- size_t index = it - truthVtxContainer.begin();
-
- const ElementLink<xAOD::TruthVertexContainer> elLink = ElementLink<xAOD::TruthVertexContainer>( truthVtxContainer, index );
-
- size_t matchIdx = indexOfMatchInfo( matchinfo, elLink );
-
- std::get<1>(matchinfo[matchIdx]) += trkWeights[t];
- std::get<2>(matchinfo[matchIdx]) += trk.pt();
- } else {
- //truth particle failed cuts
- ATH_MSG_DEBUG("Truth particle not from LLP decay.");
- otherPt += trk.pt();
- }
- } else {
- //not valid or low matching probability
- ATH_MSG_DEBUG("Invalid or low prob truth link!");
- fakePt += trk.pt();
- }
- }//end loop over tracks in vertex
-
-
-
- // normalize by total weight and pT
- std::for_each( matchinfo.begin(), matchinfo.end(), [&](VertexTruthMatchInfo& link)
- {
- std::get<1>(link) /= totalWeight;
- std::get<2>(link) /= totalPt;
- });
-
- float fakeScore = fakePt/totalPt;
- float otherScore = otherPt/totalPt;
-
- matchInfoDecor ( *vtx ) = matchinfo;
- fakeScoreDecor ( *vtx ) = fakeScore;
- otherScoreDecor( *vtx ) = otherScore;
- }
-
- //After first loop, all vertices have been decorated with their vector of match info (link to TruthVertex 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( vtxContainer.size(), false );
- static const xAOD::TruthVertex::Decorator<bool> isMatched("VertexMatchedToTruth");
- static const xAOD::TruthVertex::Decorator<bool> isSplit("VertexSplit");
-
- for ( size_t i = 0; i < vtxContainer.size(); ++i ) {
-
- if ( assignedType[i] ) {
- ATH_MSG_DEBUG("Vertex already defined as split.");
- if(m_fillHist) {
- ATH_CHECK( fillRecoPlots( *vtxContainer[i] ) );
- }
- continue; // make sure we don't reclassify vertices already found in the split loop below
- }
-
- std::vector<VertexTruthMatchInfo> & info = matchInfoDecor( *vtxContainer[i] );
- float fakeScore = fakeScoreDecor( *vtxContainer[i] );
-
- if(fakeScore > m_vxMatchWeight) {
- ATH_MSG_DEBUG("Vertex is fake.");
- matchTypeDecor( *vtxContainer[i] ) = FAKE;
- } else if (info.size() == 1) {
- if(std::get<2>(info[0]) > m_vxMatchWeight ) { // one truth matched vertex, sufficient weight
- ATH_MSG_DEBUG("One true decay vertices matched with sufficient weight. Vertex is matched.");
- matchTypeDecor( *vtxContainer[i] ) = MATCHED;
- isMatched(**std::get<0>(info[0])) = true;
- }
- else {
- ATH_MSG_DEBUG("One true decay vertices matched with insufficient weight. Vertex is other.");
- matchTypeDecor( *vtxContainer[i] ) = OTHER;
- }
- } else if (info.size() >= 2 ) { // more than one true deacy vertices matched
- ATH_MSG_DEBUG("Multiple true decay vertices matched. Vertex is merged.");
- matchTypeDecor( *vtxContainer[i] ) = MERGED;
- std::for_each( info.begin(), info.end(), [](VertexTruthMatchInfo& link)
- {
- isMatched(**std::get<0>(link)) = true;
- });
- } else { // zero truth matched vertices, but not fake
- ATH_MSG_DEBUG("Vertex is neither fake nor LLP. Marking as OTHER.");
- matchTypeDecor( *vtxContainer[i] ) = OTHER;
- }
-
- //check for splitting
- // TODO: decorate linked truth vertices with isSplit
- if ( matchTypeDecor( *vtxContainer[i] ) == MATCHED || matchTypeDecor( *vtxContainer[i] ) == MERGED ) {
- std::vector<size_t> foundSplits;
- for ( size_t j = i + 1; j < vtxContainer.size(); ++j ) {
- std::vector<VertexTruthMatchInfo> & info2 = matchInfoDecor( *vtxContainer[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( *vtxContainer[j] ) == FAKE) continue;
- if (!info2.empty() && 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( *vtxContainer[i] ).emplace_back( vtxContainer, j );
- splitPartnerDecor( *vtxContainer[j] ).emplace_back( vtxContainer, 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( *vtxContainer[k] ).emplace_back( vtxContainer, j );
- splitPartnerDecor( *vtxContainer[j] ).emplace_back( vtxContainer, k );
- }
- //then keep track that we found this one
- foundSplits.push_back(j);
- matchTypeDecor( *vtxContainer[i] ) = SPLIT;
- matchTypeDecor( *vtxContainer[j] ) = SPLIT;
- isSplit(**std::get<0>(info[0])) = true;
- assignedType[j] = true;
- } //if the two vertices match to same TruthVertex
- }//inner loop over vertices
- } //if matched or merged
-
- if(m_fillHist) {
- ATH_CHECK( fillRecoPlots( *vtxContainer[i] ) );
- }
- } //outer loop
-
- return StatusCode::SUCCESS;
-
-}
-
-StatusCode InDetSecVertexTruthMatchTool::labelTruthVertices( const xAOD::TruthVertexContainer & truthVtxContainer) {
-
- ATH_MSG_DEBUG("Labeling truth vertices");
-
- const xAOD::TrackParticleContainer * trackPartCont = nullptr;
- const xAOD::TrackParticleContainer * largeD0TrackPartCont = nullptr;
-
- ATH_CHECK( evtStore()->retrieve( trackPartCont, "InDetTrackParticles" ) );
- if ( evtStore()->contains<xAOD::TrackParticleContainer>( "InDetLargeD0TrackParticles" ) )
- ATH_CHECK( evtStore()->retrieve( largeD0TrackPartCont, "InDetLargeD0TrackParticles" ) );
- else
- ATH_MSG_WARNING("No LRT container in input! Using standard tracks only.");
-
-
- static const xAOD::Vertex::Decorator<TruthVertexMatchType> matchTypeDecor("TruthVertexMatchType");
- xAOD::TruthVertex::Decorator<bool> isMatched("VertexMatchedToTruth");
- xAOD::TruthVertex::Decorator<bool> isSplit("VertexSplit");
-
- for(const xAOD::TruthVertex* truthVtx : truthVtxContainer) {
-
- if(truthVtx->nIncomingParticles() != 1){continue;}
- const xAOD::TruthParticle* truthPart = truthVtx->incomingParticle(0);
- if(not truthPart) continue;
- if(std::find(m_pdgIdList.begin(), m_pdgIdList.end(), std::abs(truthPart->pdgId())) == m_pdgIdList.end()) continue;
- if(truthVtx->nOutgoingParticles()<2){continue;} //Skipping vertices with only 1 outgoing particle.
- ATH_MSG_DEBUG("Analysing Truth Vertex " << truthVtx );
- std::vector<const xAOD::TruthParticle*> reconstructibleParticles;
- int counter = 0;
- countReconstructibleDescendentParticles( *truthVtx, reconstructibleParticles, counter );
-
- // temporary solution for keeping track of particles in the vertex
- std::vector<int> particleInfo = {0,0,0};
- std::vector<int> vertexInfo = {0,0,0};
-
- for(size_t n = 0; n < reconstructibleParticles.size(); n++){
- ATH_MSG_DEBUG("Checking daughter no. " << n);
- const xAOD::TruthParticle* outPart = reconstructibleParticles.at(n);
-
- if (trackPartCont){
- particleInfo = checkParticle(*outPart, *trackPartCont);
- ATH_MSG_DEBUG(particleInfo);
-
- for(size_t h = 0; h < particleInfo.size(); h++){
- vertexInfo.at(h) += particleInfo.at(h);
- }
- }
- if (largeD0TrackPartCont){
- particleInfo = checkParticle(*outPart, *largeD0TrackPartCont);
- ATH_MSG_DEBUG(particleInfo);
-
- // skip first value in the tuple, we already counted it
- // in the first loop
- for(size_t h = 1; h < particleInfo.size(); h++){
- vertexInfo.at(h) += particleInfo.at(h);
- }
- }
- }
-
- ATH_MSG_DEBUG("Info for this LLP decay: " << vertexInfo);
-
- matchTypeDecor(*truthVtx) = INCLUSIVE;
- if( vertexInfo.at(0) > 1 && truthVtx->perp() < 320 && abs(truthVtx->z()) < 1500){
- ATH_MSG_DEBUG("Vertex is reconstructable and in Inner Det region");
- matchTypeDecor(*truthVtx) = RECONSTRUCTABLE;
- }
- if( matchTypeDecor(*truthVtx) == RECONSTRUCTABLE and vertexInfo.at(1) > 1){
- ATH_MSG_DEBUG("Vertex has at least two tracks");
- matchTypeDecor(*truthVtx) = ACCEPTED;
- }
- if(matchTypeDecor(*truthVtx) == ACCEPTED and vertexInfo.at(2) > 1){
- ATH_MSG_DEBUG("Vertex is has at least two tracks passing track selection: " << vertexInfo.at(2));
- matchTypeDecor(*truthVtx) = SEEDED;
- }
- if(matchTypeDecor(*truthVtx) == SEEDED and isMatched(*truthVtx)){
- ATH_MSG_DEBUG("Vertex is matched to a reconstructed secVtx");
- matchTypeDecor(*truthVtx) = RECONSTRUCTED;
- }
- if(matchTypeDecor(*truthVtx) == SEEDED and isSplit(*truthVtx)){
- ATH_MSG_DEBUG("Vertex is matched to multiple secVtx");
- matchTypeDecor(*truthVtx) = RECONSTRUCTEDSPLIT;
- }
-
- if(m_fillHist) {
- ATH_CHECK( fillTruthPlots(*truthVtx) );
- }
-
- }
- ATH_MSG_DEBUG("Done labeling truth vertices");
-
- return StatusCode::SUCCESS;
-}
-
-std::vector<int> InDetSecVertexTruthMatchTool::checkParticle(const xAOD::TruthParticle &truthPart, const xAOD::TrackParticleContainer &trkCont) const {
-
- xAOD::TrackParticle::ConstAccessor<char> trackPass("is_selected"+m_AugString);
- xAOD::TrackParticle::ConstAccessor<ElementLink<xAOD::TruthParticleContainer> > trk_truthPartAcc("truthParticleLink");
- xAOD::TrackParticle::ConstAccessor<float> trk_truthProbAcc("truthMatchProbability");
-
- if(truthPart.pt() < m_trkPtCut){
- ATH_MSG_DEBUG("Insufficient pt to reconstruct the particle");
- return {0,0,0};
- }
- else{
-
- for(const xAOD::TrackParticle* trkPart : trkCont){
- const ElementLink<xAOD::TruthParticleContainer> & truthPartLink = trk_truthPartAcc( *trkPart );
- float matchProb = trk_truthProbAcc( *trkPart );
-
- if(truthPartLink.isValid() && matchProb > m_trkMatchProb) {
- const xAOD::TruthParticle& tmpPart = **truthPartLink;
- if(tmpPart.barcode() == truthPart.barcode()) {
- if(trackPass.isAvailable( *trkPart ) and trackPass( *trkPart )) {
- ATH_MSG_DEBUG("Particle has a track that passes track selection.");
- return {1,1,1};
- }
- else {
- ATH_MSG_DEBUG("Particle has a track, but did not pass track selection.");
- return {1,1,0};
- }
- }
- }
- }
- ATH_MSG_DEBUG("Particle has enough pt.");
- return {1,0,0};
-
- }
- return {0,0,0};
-}
-
-StatusCode InDetSecVertexTruthMatchTool::fillRecoPlots( const xAOD::Vertex& secVtx ) {
-
- // set of accessors for tracks and weights
- xAOD::Vertex::ConstAccessor<xAOD::Vertex::TrackParticleLinks_t> trkAcc("trackParticleLinks");
- xAOD::Vertex::ConstAccessor<std::vector<float> > weightAcc("trackWeights");
-
- // set of decorators for truth matching info
- const xAOD::Vertex::Decorator<VertexMatchType> matchTypeDecor("VertexMatchType");
- const xAOD::Vertex::Decorator<std::vector<VertexTruthMatchInfo> > matchInfoDecor("TruthVertexMatchingInfos");
-
- std::string vertexType = "";
- switch (matchTypeDecor(secVtx))
- {
- case MATCHED:
- vertexType = "Matched";
- break;
- case MERGED:
- vertexType = "Merged";
- break;
- case SPLIT:
- vertexType = "Split";
- break;
- case OTHER:
- vertexType = "Other";
- break;
- case FAKE:
- vertexType = "Fake";
- break;
- default:
- break;
- }
-
- m_matchType->Fill(matchTypeDecor(secVtx));
-
- TVector3 reco_pos(secVtx.x(), secVtx.y(), secVtx.z());
- float reco_r = reco_pos.Perp();
-
-
- size_t ntracks;
- const xAOD::Vertex::TrackParticleLinks_t & trkParts = trkAcc( secVtx );
- ntracks = trkParts.size();
-
- TLorentzVector sumP4(0,0,0,0);
- double H = 0.0;
- double HT = 0.0;
- int charge = 0;
- double minOpAng = -1.0* 1.e10;
- double maxOpAng = 1.0* 1.e10;
- double minD0 = 1.0* 1.e10;
- double maxD0 = 0.0;
- double maxDR = 0.0;
-
- xAOD::TrackParticle::ConstAccessor< std::vector< float > > accCovMatrixDiag( "definingParametersCovMatrixDiag" );
-
- ATH_MSG_DEBUG("Loop over tracks");
- for(size_t t = 0; t < ntracks; t++){
- if(!trkParts[t].isValid()){
- ATH_MSG_DEBUG("Track " << t << " is bad!");
- continue;
- }
- const xAOD::TrackParticle & trk = **trkParts[t];
-
- double trk_d0 = std::abs(trk.definingParameters()[0]);
- double trk_z0 = std::abs(trk.definingParameters()[1]);
-
- if(trk_d0 < minD0){ minD0 = trk_d0; }
- if(trk_d0 > maxD0){ maxD0 = trk_d0; }
-
- TLorentzVector vv;
- // TODO: use values computed w.r.t SV
- vv.SetPtEtaPhiM(trk.pt(),trk.eta(), trk.phi0(), trk.m());
- sumP4 += vv;
- H += vv.Vect().Mag();
- HT += vv.Pt();
-
- TLorentzVector v_minus_iv(0,0,0,0);
- for(size_t j = 0; j < ntracks; j++){
- if (j == t){ continue; }
- if(!trkParts[j].isValid()){
- ATH_MSG_DEBUG("Track " << j << " is bad!");
- continue;
- }
-
- const xAOD::TrackParticle & trk_2 = **trkParts[j];
-
- TLorentzVector tmp;
- // TODO: use values computed w.r.t. SV
- tmp.SetPtEtaPhiM(trk_2.pt(),trk_2.eta(), trk_2.phi0(), trk_2.m());
- v_minus_iv += tmp;
-
- if( j > t ) {
- double tm = vv * tmp / ( vv.Mag() * tmp.Mag() );
- if( minOpAng < tm ) minOpAng = tm;
- if( maxOpAng > tm ) maxOpAng = tm;
- }
- }
-
- double DR = vv.DeltaR(v_minus_iv);
- if( DR > maxDR ){ maxDR = DR;}
-
- charge += trk.charge();
-
- xAOD::TrackParticle::ConstAccessor<float> Trk_Chi2("chiSquared");
- xAOD::TrackParticle::ConstAccessor<float> Trk_nDoF("numberDoF");
-
- if ( Trk_Chi2.isAvailable(trk) && Trk_Chi2(trk) && Trk_nDoF.isAvailable(trk) && Trk_nDoF(trk) ) {
- m_recoTrk_Chi2[vertexType]->Fill( Trk_Chi2(trk) / Trk_nDoF(trk));
- m_recoTrk_nDoF[vertexType]->Fill( Trk_nDoF(trk) );
- }
- m_recoTrk_charge[vertexType]->Fill(trk.charge());
- m_recoTrk_errD0[vertexType]->Fill(trk.definingParametersCovMatrix()(0,0));
- m_recoTrk_errZ0[vertexType]->Fill(trk.definingParametersCovMatrix()(1,1));
-
- m_recoTrk_theta[vertexType]->Fill(trk.definingParameters()[3]);
- m_recoTrk_qOverP[vertexType]->Fill(trk.definingParameters()[4]);
- m_recoTrk_E[vertexType]->Fill(trk.e()/GeV);
- m_recoTrk_M[vertexType]->Fill(trk.m()/GeV);
- m_recoTrk_Pt[vertexType]->Fill(trk.pt()/GeV);
-
- m_recoTrk_Px[vertexType]->Fill(trk.p4().Px()/GeV);
- m_recoTrk_Py[vertexType]->Fill(trk.p4().Py()/GeV);
- m_recoTrk_Pz[vertexType]->Fill(trk.p4().Pz()/GeV);
-
- m_recoTrk_Eta[vertexType]->Fill(trk.eta());
- m_recoTrk_Phi[vertexType]->Fill(trk.phi0());
-
- m_recoTrk_D0[vertexType]->Fill(trk_d0);
- m_recoTrk_Z0[vertexType]->Fill(trk_z0);
- } // end loop over tracks
-
- const double dir = sumP4.Vect().Dot( reco_pos ) / sumP4.Vect().Mag() / reco_pos.Mag();
-
- xAOD::Vertex::ConstAccessor<float> Chi2("chiSquared");
- xAOD::Vertex::ConstAccessor<float> nDoF("numberDoF");
- m_recoX[vertexType]->Fill(secVtx.x());
- m_recoY[vertexType]->Fill(secVtx.y());
- m_recoZ[vertexType]->Fill(secVtx.z());
- m_recoR[vertexType]->Fill(reco_r);
- m_recoNtrk[vertexType]->Fill(ntracks);
- m_recoPt[vertexType]->Fill(sumP4.Pt() / GeV);
- m_recoEta[vertexType]->Fill(sumP4.Eta());
- m_recoPhi[vertexType]->Fill(sumP4.Phi());
- m_recoMass[vertexType]->Fill(sumP4.M() / GeV);
- m_recoMu[vertexType]->Fill(sumP4.M()/maxDR / GeV);
- m_recoChi2[vertexType]->Fill(Chi2(secVtx)/nDoF(secVtx));
- m_recoDir[vertexType]->Fill(dir);
- m_recoCharge[vertexType]->Fill(charge);
- m_recoH[vertexType]->Fill(H / GeV);
- m_recoHt[vertexType]->Fill(HT / GeV);
- m_recoMinOpAng[vertexType]->Fill(minOpAng);
- m_recoMaxOpAng[vertexType]->Fill(maxOpAng);
- m_recoMinD0[vertexType]->Fill(minD0);
- m_recoMaxD0[vertexType]->Fill(maxD0);
- m_recoMaxDR[vertexType]->Fill(maxDR);
-
- m_recoX["All"]->Fill(secVtx.x());
- m_recoY["All"]->Fill(secVtx.y());
- m_recoZ["All"]->Fill(secVtx.z());
- m_recoR["All"]->Fill(reco_r);
- m_recoNtrk["All"]->Fill(ntracks);
- m_recoPt["All"]->Fill(sumP4.Pt() / GeV);
- m_recoEta["All"]->Fill(sumP4.Eta());
- m_recoPhi["All"]->Fill(sumP4.Phi());
- m_recoMass["All"]->Fill(sumP4.M() / GeV);
- m_recoMu["All"]->Fill(sumP4.M()/maxDR / GeV);
- m_recoChi2["All"]->Fill(Chi2(secVtx)/nDoF(secVtx));
- m_recoDir["All"]->Fill(dir);
- m_recoCharge["All"]->Fill(charge);
- m_recoH["All"]->Fill(H / GeV);
- m_recoHt["All"]->Fill(HT / GeV);
- m_recoMinOpAng["All"]->Fill(minOpAng);
- m_recoMaxOpAng["All"]->Fill(maxOpAng);
- m_recoMinD0["All"]->Fill(minD0);
- m_recoMaxD0["All"]->Fill(maxD0);
- m_recoMaxDR["All"]->Fill(maxDR);
-
- std::vector<VertexTruthMatchInfo> truthmatchinfo;
- truthmatchinfo = matchInfoDecor(secVtx);
-
- // This includes all matched vertices, including splits
- if(not truthmatchinfo.empty()){
- float matchScore_weight = std::get<1>(truthmatchinfo.at(0));
- float matchScore_pt = std::get<2>(truthmatchinfo.at(0));
-
- m_matchScore_weight[vertexType]->Fill(matchScore_weight);
- m_matchScore_pt[vertexType]->Fill(matchScore_pt);
-
- ATH_MSG_DEBUG("Match Score and probability: " << matchScore_weight << " " << matchScore_pt/0.01);
-
- const ElementLink<xAOD::TruthVertexContainer>& truthVertexLink = std::get<0>(truthmatchinfo.at(0));
- const xAOD::TruthVertex& truthVtx = **truthVertexLink ;
-
- m_positionRes_R[vertexType]->Fill(reco_r - truthVtx.perp());
- m_positionRes_Z[vertexType]->Fill(secVtx.z() - truthVtx.z());
-
- m_recoX["LLP"]->Fill(secVtx.x());
- m_recoY["LLP"]->Fill(secVtx.y());
- m_recoZ["LLP"]->Fill(secVtx.z());
- m_recoR["LLP"]->Fill(reco_r);
- m_recoNtrk["LLP"]->Fill(ntracks);
- m_recoPt["LLP"]->Fill(sumP4.Pt() / GeV);
- m_recoEta["LLP"]->Fill(sumP4.Eta());
- m_recoPhi["LLP"]->Fill(sumP4.Phi());
- m_recoMass["LLP"]->Fill(sumP4.M() / GeV);
- m_recoMu["LLP"]->Fill(sumP4.M()/maxDR / GeV);
- m_recoChi2["LLP"]->Fill(Chi2(secVtx)/nDoF(secVtx));
- m_recoDir["LLP"]->Fill(dir);
- m_recoCharge["LLP"]->Fill(charge);
- m_recoH["LLP"]->Fill(H / GeV);
- m_recoHt["LLP"]->Fill(HT / GeV);
- m_recoMinOpAng["LLP"]->Fill(minOpAng);
- m_recoMaxOpAng["LLP"]->Fill(maxOpAng);
- m_recoMinD0["LLP"]->Fill(minD0);
- m_recoMaxD0["LLP"]->Fill(maxD0);
- m_recoMaxDR["LLP"]->Fill(maxDR);
-
- m_positionRes_R["LLP"]->Fill(reco_r - truthVtx.perp());
- m_positionRes_Z["LLP"]->Fill(secVtx.z() - truthVtx.z());
- m_matchScore_weight["LLP"]->Fill(matchScore_weight);
- m_matchScore_pt["LLP"]->Fill(matchScore_pt);
- }
-
- return StatusCode::SUCCESS;
-}
-
-StatusCode InDetSecVertexTruthMatchTool::fillTruthPlots( const xAOD::TruthVertex& truthVtx) {
-
- ATH_MSG_DEBUG("Plotting truth vertex");
-
- m_truthX->Fill(truthVtx.x());
- m_truthY->Fill(truthVtx.y());
- m_truthZ->Fill(truthVtx.z());
- m_truthR->Fill(truthVtx.perp());
- m_truthEta->Fill(truthVtx.eta());
- m_truthPhi->Fill(truthVtx.phi());
- m_truthNtrk_out->Fill(truthVtx.nOutgoingParticles());
-
- ATH_MSG_DEBUG("Plotting truth parent");
- const xAOD::TruthParticle& truthPart = *truthVtx.incomingParticle(0);
-
- m_truthParent_E->Fill(truthPart.e() / GeV);
- m_truthParent_M->Fill(truthPart.m() / GeV);
- m_truthParent_Pt->Fill(truthPart.pt() / GeV);
- m_truthParent_Phi->Fill(truthPart.phi());
- m_truthParent_Eta->Fill(truthPart.eta());
- m_truthParent_charge->Fill(truthPart.charge());
-
- ATH_MSG_DEBUG("Plotting truth prod vtx");
- if(truthPart.hasProdVtx()){
- const xAOD::TruthVertex & vertex = *truthPart.prodVtx();
-
- m_truthParentProdX->Fill(vertex.x());
- m_truthParentProdY->Fill(vertex.y());
- m_truthParentProdZ->Fill(vertex.z());
- m_truthParentProdR->Fill(vertex.perp());
- m_truthParentProdEta->Fill(vertex.eta());
- m_truthParentProdPhi->Fill(vertex.phi());
- }
-
- ATH_MSG_DEBUG("Plotting match types");
- static const xAOD::Vertex::Decorator<TruthVertexMatchType> matchTypeDecor("TruthVertexMatchType");
-
- m_truthInclusive_r->Fill(truthVtx.perp());
-
- if(matchTypeDecor(truthVtx) >= RECONSTRUCTABLE){
- m_truthReconstructable_r->Fill(truthVtx.perp());
- }
- if(matchTypeDecor(truthVtx) >= ACCEPTED){
- m_truthAccepted_r->Fill(truthVtx.perp());
- }
- if(matchTypeDecor(truthVtx) >= SEEDED){
- m_truthSeeded_r->Fill(truthVtx.perp());
- }
- if(matchTypeDecor(truthVtx) >= RECONSTRUCTED){
- m_truthReconstructed_r->Fill(truthVtx.perp());
- }
- if(matchTypeDecor(truthVtx) >= RECONSTRUCTEDSPLIT){
- m_truthSplit_r->Fill(truthVtx.perp());
- }
- return StatusCode::SUCCESS;
-}
-
-
-// check if truth particle originated from decay of particle in the pdgIdList
-int InDetSecVertexTruthMatchTool::checkProduction( const xAOD::TruthParticle & truthPart ) const {
-
- if (truthPart.nParents() == 0){
- ATH_MSG_DEBUG("Particle has no parents (end of loop)");
- return -1;
- }
- else{
- const xAOD::TruthParticle * parent = truthPart.parent(0);
- if(not parent) {
- ATH_MSG_DEBUG("Particle parent is null");
- return -1;
- }
- ATH_MSG_DEBUG("Parent ID: " << parent->pdgId());
-
- if(std::find(m_pdgIdList.begin(), m_pdgIdList.end(), std::abs(parent->pdgId())) != m_pdgIdList.end()) {
- ATH_MSG_DEBUG("Found LLP decay.");
- const xAOD::TruthVertex* vertex = parent->decayVtx();
- return vertex->barcode();
- }
- // recurse on parent
- return checkProduction(*parent);
- }
- return -1;
-}
-
-void InDetSecVertexTruthMatchTool::countReconstructibleDescendentParticles(const xAOD::TruthVertex& signalTruthVertex,
- std::vector<const xAOD::TruthParticle*>& set, int counter) const {
-
- counter++;
-
- for( size_t itrk = 0; itrk < signalTruthVertex.nOutgoingParticles(); itrk++) {
- const auto* particle = signalTruthVertex.outgoingParticle( itrk );
- if( !particle ) continue;
- // Recursively add descendents
- if( particle->hasDecayVtx() ) {
-
- TVector3 decayPos( particle->decayVtx()->x(), particle->decayVtx()->y(), particle->decayVtx()->z() );
- TVector3 prodPos ( particle->prodVtx()->x(), particle->prodVtx()->y(), particle->prodVtx()->z() );
-
- auto isInside = []( TVector3& v ) { return ( v.Perp() < 300. && std::abs( v.z() ) < 1500. ); };
- auto isOutside = []( TVector3& v ) { return ( v.Perp() > 563. || std::abs( v.z() ) > 2720. ); };
-
- const auto distance = (decayPos - prodPos).Mag();
-
- if (counter > 100) {
- ATH_MSG_WARNING("Vetoing particle that may be added recursively infinitely (potential loop in generator record");
- break;
- }
-
- // consider track reconstructible if it travels at least 10mm
- if( distance < 10.0 ) {
- countReconstructibleDescendentParticles( *particle->decayVtx(), set , counter);
- } else if( isInside ( prodPos ) && isOutside( decayPos ) && particle->isCharged() ) {
- set.push_back( particle );
- } else if( particle->isElectron() || particle->isMuon() ) {
- set.push_back( particle );
- }
- } else {
- if( !(particle->isCharged()) ) continue;
- set.push_back( particle );
- }
- }
-
- }
diff --git a/InnerDetector/InDetValidation/InDetSecVertexValidation/src/InDetSecVertexTruthMatchAlgorithm.cxx b/InnerDetector/InDetValidation/InDetSecVertexValidation/src/InDetSecVertexTruthMatchAlgorithm.cxx
deleted file mode 100644
index 157c24ef731fb6bcb4f5d132200430743643afc1..0000000000000000000000000000000000000000
--- a/InnerDetector/InDetValidation/InDetSecVertexValidation/src/InDetSecVertexTruthMatchAlgorithm.cxx
+++ /dev/null
@@ -1,50 +0,0 @@
-/*
- Copyright (C) 2002-2022 CERN for the benefit of the ATLAS collaboration
-*/
-
-#include "InDetSecVertexTruthMatchAlgorithm.h"
-
-#include "xAODTracking/VertexContainer.h"
-#include "xAODTruth/TruthVertexContainer.h"
-
-
-InDetSecVertexTruthMatchAlgorithm::InDetSecVertexTruthMatchAlgorithm( const std::string& name, ISvcLocator* svcLoc )
- : AthAlgorithm( name, svcLoc ),
- m_matchTool( "InDetSecVertexTruthMatchTool/InDetSecVertexTruthMatchTool", this ) {
-
- declareProperty( "SecVertexSGKey", m_secVtxSGKey = "VrtSecInclusive_SecondaryVertices" );
- declareProperty( "TruthVertexSGKey", m_truthVtxSGKey = "TruthBSMWithDecayVertices" );
- declareProperty( "VertexTruthMatchTool", m_matchTool );
-}
-
-StatusCode InDetSecVertexTruthMatchAlgorithm::initialize() {
-
- // Greet the user:
- ATH_MSG_INFO( "Initialising" );
- ATH_MSG_DEBUG( "SecVertexSGKey = " << m_secVtxSGKey );
- ATH_MSG_DEBUG( "TruthVertexSGKey = " << m_truthVtxSGKey );
-
- // Retrieve the tools:
- ATH_CHECK( m_matchTool.retrieve() );
-
- // Return gracefully:
- return StatusCode::SUCCESS;
-}
-
-StatusCode InDetSecVertexTruthMatchAlgorithm::execute() {
-
- //Retrieve the vertices:
- const xAOD::VertexContainer * inSecVert = nullptr;
- ATH_CHECK( evtStore()->retrieve( inSecVert, m_secVtxSGKey ) );
-
- //Retrieve truth vertices:
- const xAOD::TruthVertexContainer * inTruthVert = nullptr;
- ATH_CHECK( evtStore()->retrieve( inTruthVert, m_truthVtxSGKey ) );
-
- //pass to the tool for decoration:
- ATH_CHECK( m_matchTool->matchVertices( *inSecVert, *inTruthVert ) );
- ATH_CHECK( m_matchTool->labelTruthVertices( *inTruthVert ) );
-
- return StatusCode::SUCCESS;
-
-}
diff --git a/InnerDetector/InDetValidation/InDetSecVertexValidation/src/InDetSecVertexTruthMatchAlgorithm.h b/InnerDetector/InDetValidation/InDetSecVertexValidation/src/InDetSecVertexTruthMatchAlgorithm.h
deleted file mode 100644
index 664b5e56d6c84961f4d6893fee0714bb35791b78..0000000000000000000000000000000000000000
--- a/InnerDetector/InDetValidation/InDetSecVertexValidation/src/InDetSecVertexTruthMatchAlgorithm.h
+++ /dev/null
@@ -1,35 +0,0 @@
-/*
- Copyright (C) 2002-2021 CERN for the benefit of the ATLAS collaboration
-*/
-
-#ifndef INDETSECVERTEXTRUTHMATCHALGORITHM_H
-#define INDETSECVERTEXTRUTHMATCHALGORITHM_H
-
-// Gaudi/Athena include(s):
-#include "AthenaBaseComps/AthAlgorithm.h"
-#include "AsgTools/ToolHandle.h"
-
-// Local include(s):
-#include "InDetSecVertexValidation/IInDetSecVertexTruthMatchTool.h"
-
-class InDetSecVertexTruthMatchAlgorithm : public AthAlgorithm {
-
- public:
- /// Regular Algorithm constructor
- InDetSecVertexTruthMatchAlgorithm( 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 vertex container to investigate
- std::string m_secVtxSGKey;
- std::string m_truthVtxSGKey;
-
- ToolHandle<IInDetSecVertexTruthMatchTool> m_matchTool;
-
-};
-
-#endif
diff --git a/InnerDetector/InDetValidation/InDetSecVertexValidation/src/components/InDetSecVertexValidation_entries.cxx b/InnerDetector/InDetValidation/InDetSecVertexValidation/src/components/InDetSecVertexValidation_entries.cxx
index 7aa7e5189c52599b8bcc3d83de11cb59a2ac3a58..fbbbc38f6c2b3a9ca28723b2dded8092e8d6acca 100644
--- a/InnerDetector/InDetValidation/InDetSecVertexValidation/src/components/InDetSecVertexValidation_entries.cxx
+++ b/InnerDetector/InDetValidation/InDetSecVertexValidation/src/components/InDetSecVertexValidation_entries.cxx
@@ -1,15 +1,11 @@
/*
- Copyright (C) 2002-2021 CERN for the benefit of the ATLAS collaboration
+ Copyright (C) 2002-2024 CERN for the benefit of the ATLAS collaboration
*/
// Gaudi/Athena include(s):
//#include "GaudiKernel/DeclareFactoryEntries.h"
// Local include(s):
-#include "InDetSecVertexValidation/InDetSecVertexTruthMatchTool.h"
-#include "../InDetSecVertexTruthMatchAlgorithm.h"
#include "../PhysValSecVtx.h"
-DECLARE_COMPONENT( InDetSecVertexTruthMatchTool )
-DECLARE_COMPONENT( InDetSecVertexTruthMatchAlgorithm )
DECLARE_COMPONENT( PhysValSecVtx )
\ No newline at end of file
diff --git a/InnerDetector/InDetValidation/InDetSecVertexValidation/util/VertexTruthMatchTest.cxx b/InnerDetector/InDetValidation/InDetSecVertexValidation/util/VertexTruthMatchTest.cxx
deleted file mode 100644
index 5874e20cc7e154a6afb0dbc16b813e72d892fe53..0000000000000000000000000000000000000000
--- a/InnerDetector/InDetValidation/InDetSecVertexValidation/util/VertexTruthMatchTest.cxx
+++ /dev/null
@@ -1,143 +0,0 @@
-/*
- Copyright (C) 2002-2021 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 "InDetSecVertexValidation/InDetSecVertexTruthMatchTool.h"
-#include "InDetSecVertexValidation/InDetSecVertexTruthMatchUtils.h"
-
-#include "AsgMessaging/MessageCheck.h"
-
-int main( int argc, char* argv[] ) {
-
- using namespace asg::msgUserCode;
-
- 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
- InDetSecVertexTruthMatchTool matchTool("VertexTruthMatchTool");
-
- //config
- matchTool.msg().setLevel( MSG::DEBUG );
- ATH_CHECK( matchTool.setProperty( "trackMatchProb", 0.7 ) );
- ATH_CHECK( matchTool.setProperty( "vertexMatchWeight", 0.5 ) );
- ATH_CHECK( matchTool.setProperty( "pdgIds", "36" ) );
- ATH_CHECK( 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");
-
-
- //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, "VrtSecInclusive_SecondaryVertices" );
- if(!result) {
- Error( APP_NAME, "Failed to retrieve VrtSecInclusive_SecondaryVertices on entry %i", static_cast<int>(entry));
- return 1;
- }
- const xAOD::TruthVertexContainer * tvxContainer = nullptr;
- result = event.retrieve( tvxContainer, "TruthVertices" );
- if(!result) {
- Error( APP_NAME, "Failed to retrieve TruthVertices on entry %i", static_cast<int>(entry));
- return 1;
- }
-
- int status = matchTool.matchVertices( *vxContainer, *tvxContainer );
- 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<InDetSecVertexTruthMatchUtils::VertexMatchType> matchTypeDecor("VertexMatchType");
-
- //fill the category histogram for all vertices
- for( auto vxit : *vxContainer ) {
- h_cat->Fill( matchTypeDecor( *vxit ) );
- }
-
- }
-
- h_cat->Write();
- fout.Close();
-
- delete h_cat;
-
- return 0;
-}
diff --git a/PhysicsAnalysis/Algorithms/TrackingAnalysisAlgorithms/CMakeLists.txt b/PhysicsAnalysis/Algorithms/TrackingAnalysisAlgorithms/CMakeLists.txt
index 14040b7af0d90b0c3af87827cce12ec1df08eb95..244484d5775a8c5409663356de1013b7ca1dd64b 100644
--- a/PhysicsAnalysis/Algorithms/TrackingAnalysisAlgorithms/CMakeLists.txt
+++ b/PhysicsAnalysis/Algorithms/TrackingAnalysisAlgorithms/CMakeLists.txt
@@ -7,7 +7,7 @@ atlas_subdir( TrackingAnalysisAlgorithms )
atlas_add_library( TrackingAnalysisAlgorithmsLib
TrackingAnalysisAlgorithms/*.h Root/*.cxx
PUBLIC_HEADERS TrackingAnalysisAlgorithms
- LINK_LIBRARIES AnaAlgorithmLib EventBookkeeperToolsLib xAODTracking
+ LINK_LIBRARIES AnaAlgorithmLib EventBookkeeperToolsLib xAODTracking InDetSecVtxTruthMatchToolLib
AsgDataHandlesLib AsgTools AthContainers )
atlas_add_dictionary( TrackingAnalysisAlgorithmsDict
@@ -20,3 +20,6 @@ if( NOT XAOD_STANDALONE )
src/*.h src/*.cxx src/components/*.cxx
LINK_LIBRARIES GaudiKernel TrackingAnalysisAlgorithmsLib )
endif()
+
+atlas_install_python_modules( python/*.py )
+atlas_install_runtime( scripts/*.py POST_BUILD_CMD ${ATLAS_FLAKE8} )
diff --git a/PhysicsAnalysis/Algorithms/TrackingAnalysisAlgorithms/Root/SecVertexTruthMatchAlg.cxx b/PhysicsAnalysis/Algorithms/TrackingAnalysisAlgorithms/Root/SecVertexTruthMatchAlg.cxx
new file mode 100644
index 0000000000000000000000000000000000000000..5dd3c17368d70688cf4bf1de491c3c604018f998
--- /dev/null
+++ b/PhysicsAnalysis/Algorithms/TrackingAnalysisAlgorithms/Root/SecVertexTruthMatchAlg.cxx
@@ -0,0 +1,412 @@
+/*
+ Copyright (C) 2002-2024 CERN for the benefit of the ATLAS collaboration
+*/
+
+#include "TrackingAnalysisAlgorithms/SecVertexTruthMatchAlg.h"
+#include "InDetSecVtxTruthMatchTool/InDetSecVtxTruthMatchTool.h"
+#include "xAODTruth/TruthParticle.h"
+
+#include "TMath.h"
+#include "TH1.h"
+#include "TEfficiency.h"
+
+const float GeV = 1000.;
+
+namespace CP {
+
+ SecVertexTruthMatchAlg::SecVertexTruthMatchAlg( const std::string& name, ISvcLocator* svcLoc )
+ : EL::AnaAlgorithm( name, svcLoc ) {}
+
+ StatusCode SecVertexTruthMatchAlg::initialize() {
+
+ // Initializing Keys
+ ATH_CHECK(m_secVtxContainerKey.initialize());
+ ATH_CHECK(m_truthVtxContainerKey.initialize());
+ ATH_CHECK(m_trackParticleContainerKey.initialize());
+
+ // Retrieving the tool
+ ATH_CHECK(m_matchTool.retrieve());
+
+ if(m_writeHistograms) {
+ std::vector<std::string> recoTypes{"All", "Matched", "Merged", "Fake", "Split", "Other"};
+ std::vector<std::string> truthTypes{"Inclusive", "Reconstructable", "Accepted", "Seeded", "Reconstructed", "ReconstructedSplit"};
+
+ ANA_CHECK (book(TH1F("RecoVertex/matchType", "Vertex Match Type", 8, -0.5, 7.5)));
+
+ for(const auto& recoType : recoTypes) {
+ ANA_CHECK (book(TH1F(("RecoVertex/" + recoType + "_x").c_str(), "Reco vertex x [mm]", 1000, -500, 500)));
+ ANA_CHECK (book(TH1F(("RecoVertex/" + recoType + "_y").c_str(), "Reco vertex y [mm]", 1000, -500, 500)));
+ ANA_CHECK (book(TH1F(("RecoVertex/" + recoType + "_z").c_str(), "Reco vertex z [mm]", 1000, -500, 500)));
+ ANA_CHECK (book(TH1F(("RecoVertex/" + recoType + "_Lxy").c_str(), "Reco vertex L_{xy} [mm]", 500, 0, 500)));
+ ANA_CHECK (book(TH1F(("RecoVertex/" + recoType + "_pT").c_str(), "Reco vertex p_{T} [GeV]", 100, 0, 100)));
+ ANA_CHECK (book(TH1F(("RecoVertex/" + recoType + "_eta").c_str(), "Reco vertex #eta", 100, -5, 5)));
+ ANA_CHECK (book(TH1F(("RecoVertex/" + recoType + "_phi").c_str(), "Reco vertex #phi", 100, -TMath::Pi(), TMath::Pi())));
+ ANA_CHECK (book(TH1F(("RecoVertex/" + recoType + "_mass").c_str(), "Reco vertex mass [GeV]", 500, 0, 100)));
+ ANA_CHECK (book(TH1F(("RecoVertex/" + recoType + "_mu").c_str(), "Reco vertex Red. Mass [GeV]", 500, 0, 100)));
+ ANA_CHECK (book(TH1F(("RecoVertex/" + recoType + "_chi2").c_str(), "Reco vertex recoChi2", 100, 0, 10)));
+ ANA_CHECK (book(TH1F(("RecoVertex/" + recoType + "_dir").c_str(), "Reco vertex recoDirection", 100, -1, 1)));
+ ANA_CHECK (book(TH1F(("RecoVertex/" + recoType + "_charge").c_str(), "Reco vertex recoCharge", 20, -10, 10)));
+ ANA_CHECK (book(TH1F(("RecoVertex/" + recoType + "_H").c_str(), "Reco vertex H [GeV]", 100, 0, 100)));
+ ANA_CHECK (book(TH1F(("RecoVertex/" + recoType + "_HT").c_str(), "Reco vertex Mass [GeV]", 100, 0, 100)));
+ ANA_CHECK (book(TH1F(("RecoVertex/" + recoType + "_minOpAng").c_str(), "Reco vertex minOpAng", 100, -1, 1)));
+ ANA_CHECK (book(TH1F(("RecoVertex/" + recoType + "_maxOpAng").c_str(), "Reco vertex maxOpAng", 100, -1, 1)));
+ ANA_CHECK (book(TH1F(("RecoVertex/" + recoType + "_maxdR").c_str(), "Reco vertex maxDR", 100, 0, 10)));
+ ANA_CHECK (book(TH1F(("RecoVertex/" + recoType + "_mind0").c_str(), "Reco vertex min d0 [mm]", 100, 0, 100)));
+ ANA_CHECK (book(TH1F(("RecoVertex/" + recoType + "_maxd0").c_str(), "Reco vertex max d0 [mm]", 100, 0, 100)));
+ ANA_CHECK (book(TH1F(("RecoVertex/" + recoType + "_ntrk").c_str(), "Reco vertex n tracks", 30, 0, 30)));
+
+ // tracks
+ ANA_CHECK (book(TH1F(("RecoVertex/" + recoType + "_Trk_qOverP").c_str(), "Reco track qOverP ", 100, 0, .01)));
+ ANA_CHECK (book(TH1F(("RecoVertex/" + recoType + "_Trk_theta").c_str(), "Reco track theta ", 64, 0, 3.2)));
+ ANA_CHECK (book(TH1F(("RecoVertex/" + recoType + "_Trk_E").c_str(), "Reco track E ", 100, 0, 100)));
+ ANA_CHECK (book(TH1F(("RecoVertex/" + recoType + "_Trk_M").c_str(), "Reco track M ", 100, 0, 10)));
+ ANA_CHECK (book(TH1F(("RecoVertex/" + recoType + "_Trk_Pt").c_str(), "Reco track Pt ", 100, 0, 100)));
+ ANA_CHECK (book(TH1F(("RecoVertex/" + recoType + "_Trk_Px").c_str(), "Reco track Px ", 100, 0, 100)));
+ ANA_CHECK (book(TH1F(("RecoVertex/" + recoType + "_Trk_Py").c_str(), "Reco track Py ", 100, 0, 100)));
+ ANA_CHECK (book(TH1F(("RecoVertex/" + recoType + "_Trk_Pz").c_str(), "Reco track Pz ", 100, 0, 100)));
+ ANA_CHECK (book(TH1F(("RecoVertex/" + recoType + "_Trk_Eta").c_str(), "Reco track Eta ", 100, -5, 5)));
+ ANA_CHECK (book(TH1F(("RecoVertex/" + recoType + "_Trk_Phi").c_str(), "Reco track Phi ", 63, -3.2, 3.2)));
+ ANA_CHECK (book(TH1F(("RecoVertex/" + recoType + "_Trk_D0").c_str(), "Reco track D0 ", 300, -300, 300)));
+ ANA_CHECK (book(TH1F(("RecoVertex/" + recoType + "_Trk_Z0").c_str(), "Reco track Z0 ", 500, -500, 500)));
+ ANA_CHECK (book(TH1F(("RecoVertex/" + recoType + "_Trk_errD0").c_str(), "Reco track errD0 ", 300, 0, 30)));
+ ANA_CHECK (book(TH1F(("RecoVertex/" + recoType + "_Trk_errZ0").c_str(), "Reco track errZ0 ", 500, 0, 50)));
+ ANA_CHECK (book(TH1F(("RecoVertex/" + recoType + "_Trk_Chi2").c_str(), "Reco track Chi2 ", 100, 0, 10)));
+ ANA_CHECK (book(TH1F(("RecoVertex/" + recoType + "_Trk_nDoF").c_str(), "Reco track nDoF ", 100, 0, 100)));
+ ANA_CHECK (book(TH1F(("RecoVertex/" + recoType + "_Trk_charge").c_str(), "Reco track charge ", 3, -1.5, 1.5)));
+
+ // truth matching -- don't book for non-matched vertices
+ if ( recoType != "Inclusive" and recoType != "Fake" ) {
+ ANA_CHECK (book(TH1F(("RecoVertex/" + recoType + "_positionRes_R").c_str(), "Position resolution for vertices matched to truth decays", 400, -20, 20)));
+ ANA_CHECK (book(TH1F(("RecoVertex/" + recoType + "_positionRes_Z").c_str(), "Position resolution for vertices matched to truth decays", 400, -20, 20)));
+ ANA_CHECK (book(TH1F(("RecoVertex/" + recoType + "_matchScore_weight").c_str(), "Vertex Match Score (weight)", 101, 0, 1.01)));
+ ANA_CHECK (book(TH1F(("RecoVertex/" + recoType + "_matchScore_pt").c_str(), "Vertex Match Score (pT)", 101, 0, 1.01)));
+ ANA_CHECK (book(TH1F(("RecoVertex/" + recoType + "_matchedTruthID").c_str(), "Vertex Truth Match ID", 100, 0, 100)));
+ }
+ }
+ for(const auto& truthType : truthTypes) {
+ ANA_CHECK (book(TH1F(("TruthVertex/" + truthType + "_x").c_str(), "Truth vertex x [mm]", 1000, -500, 500)));
+ ANA_CHECK (book(TH1F(("TruthVertex/" + truthType + "_y").c_str(), "Truth vertex y [mm]", 500, -500, 500)));
+ ANA_CHECK (book(TH1F(("TruthVertex/" + truthType + "_z").c_str(), "Truth vertex z [mm]", 500, -500, 500)));
+ ANA_CHECK (book(TH1F(("TruthVertex/" + truthType + "_R").c_str(), "Truth vertex r [mm]", 6000, 0, 600)));
+ ANA_CHECK (book(TH1F(("TruthVertex/" + truthType + "_distFromPV").c_str(), "Truth vertex distFromPV [mm]", 500, 0, 500)));
+ ANA_CHECK (book(TH1F(("TruthVertex/" + truthType + "_Eta").c_str(), "Truth vertex Eta", 100, -5, 5)));
+ ANA_CHECK (book(TH1F(("TruthVertex/" + truthType + "_Phi").c_str(), "Truth vertex Phi", 64, -3.2, 3.2)));
+ ANA_CHECK (book(TH1F(("TruthVertex/" + truthType + "_Ntrk_out").c_str(), "Truth vertex n outgoing tracks", 100, 0, 100)));
+ ANA_CHECK (book(TH1F(("TruthVertex/" + truthType + "_Parent_E").c_str(), "Reco track E", 100, 0, 100)));
+ ANA_CHECK (book(TH1F(("TruthVertex/" + truthType + "_Parent_M").c_str(), "Reco track M", 500, 0, 500)));
+ ANA_CHECK (book(TH1F(("TruthVertex/" + truthType + "_Parent_Pt").c_str(), "Reco track Pt", 100, 0, 100)));
+ ANA_CHECK (book(TH1F(("TruthVertex/" + truthType + "_Parent_Eta").c_str(), "Reco track Eta", 100, -5, 5)));
+ ANA_CHECK (book(TH1F(("TruthVertex/" + truthType + "_Parent_Phi").c_str(), "Reco track Phi", 63, -3.2, 3.2)));
+ ANA_CHECK (book(TH1F(("TruthVertex/" + truthType + "_Parent_charge").c_str(), "Reco track charge", 3, -1, 1)));
+ ANA_CHECK (book(TH1F(("TruthVertex/" + truthType + "_ParentProdX").c_str(), "truthParentProd vertex x [mm]", 500, -500, 500)));
+ ANA_CHECK (book(TH1F(("TruthVertex/" + truthType + "_ParentProdY").c_str(), "truthParentProd vertex y [mm]", 500, -500, 500)));
+ ANA_CHECK (book(TH1F(("TruthVertex/" + truthType + "_ParentProdZ").c_str(), "truthParentProd vertex z [mm]", 500, -500, 500)));
+ ANA_CHECK (book(TH1F(("TruthVertex/" + truthType + "_ParentProdR").c_str(), "truthParentProd vertex r [mm]", 6000, 0, 600)));
+ ANA_CHECK (book(TH1F(("TruthVertex/" + truthType + "_ParentProddistFromPV").c_str(), "truthParentProd vertex distFromPV [mm]", 500, 0, 500)));
+ ANA_CHECK (book(TH1F(("TruthVertex/" + truthType + "_ParentProdEta").c_str(), "truthParentProd vertex Eta", 100, -5, 5)));
+ ANA_CHECK (book(TH1F(("TruthVertex/" + truthType + "_ParentProdPhi").c_str(), "truthParentProd vertex Phi", 64, -3.2, 3.2)));
+ }
+ // now add the efficiencies
+ Double_t bins[] = {0.0, 0.1, 0.2, 0.3, 0.4, 0.5, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 15, 20, 25, 30, 35, 40, 50, 60, 70, 80, 90, 100, 125, 150, 200, 300, 500};
+ size_t nbins = sizeof(bins)/sizeof(bins[0])-1;
+ ANA_CHECK (book(TEfficiency("Acceptance", "Acceptance", nbins, bins)));
+ ANA_CHECK (book(TEfficiency("eff_seed", "Seed efficiency", nbins, bins)));
+ ANA_CHECK (book(TEfficiency("eff_core", "Core efficiency", nbins, bins)));
+ ANA_CHECK (book(TEfficiency("eff_total", "Total efficiency", nbins, bins)));
+ }
+
+ return StatusCode::SUCCESS;
+ }
+
+ StatusCode SecVertexTruthMatchAlg::execute() {
+
+ //Retrieve the vertices:
+ SG::ReadHandle<xAOD::VertexContainer> recoVertexContainer(m_secVtxContainerKey);
+ SG::ReadHandle<xAOD::TruthVertexContainer> truthVertexContainer(m_truthVtxContainerKey);
+ SG::ReadHandle<xAOD::TrackParticleContainer> trackParticleContainer(m_trackParticleContainerKey);
+
+ std::vector<const xAOD::Vertex*> recoVerticesToMatch;
+ std::vector<const xAOD::TruthVertex*> truthVerticesToMatch;
+
+ for(const auto &recoVertex : *recoVertexContainer) {
+ xAOD::VxType::VertexType vtxType = static_cast<xAOD::VxType::VertexType>( recoVertex->vertexType() );
+
+ if(vtxType != xAOD::VxType::SecVtx ){
+ ATH_MSG_DEBUG("Vertex not labeled as secondary");
+ continue;
+ }
+ recoVerticesToMatch.push_back(recoVertex);
+ }
+
+ for(const auto &truthVertex : *truthVertexContainer) {
+ if(truthVertex->nIncomingParticles() != 1) {
+ continue;
+ }
+ const xAOD::TruthParticle* truthPart = truthVertex->incomingParticle(0);
+ if(not truthPart) {
+ continue;
+ }
+ if(std::find(m_targetPDGIDs.begin(), m_targetPDGIDs.end(), std::abs(truthPart->pdgId())) == m_targetPDGIDs.end()) {
+ continue;
+ }
+ if(truthVertex->nOutgoingParticles() < 2) {
+ continue;
+ }
+ truthVerticesToMatch.push_back(truthVertex);
+ }
+
+ //pass to the tool for decoration:
+ ATH_CHECK( m_matchTool->matchVertices( recoVerticesToMatch, truthVerticesToMatch, trackParticleContainer.cptr() ) );
+
+ if(m_writeHistograms) {
+ static const xAOD::Vertex::Decorator<int> matchTypeDecor("vertexMatchType");
+ for(const auto& secVtx : recoVerticesToMatch) {
+ int matchTypeBitset = matchTypeDecor(*secVtx);
+ hist("RecoVertex/matchType")->Fill(matchTypeBitset);
+
+ if(InDetSecVtxTruthMatchUtils::isMatched(matchTypeBitset)) {
+ fillRecoHistograms(secVtx, "Matched");
+ }
+ if(InDetSecVtxTruthMatchUtils::isMerged(matchTypeBitset)) {
+ fillRecoHistograms(secVtx, "Merged");
+ }
+ if(InDetSecVtxTruthMatchUtils::isFake(matchTypeBitset)) {
+ fillRecoHistograms(secVtx, "Fake");
+ }
+ if(InDetSecVtxTruthMatchUtils::isSplit(matchTypeBitset)) {
+ fillRecoHistograms(secVtx, "Split");
+ }
+ if(InDetSecVtxTruthMatchUtils::isOther(matchTypeBitset)) {
+ fillRecoHistograms(secVtx, "Other");
+ }
+ fillRecoHistograms(secVtx, "All");
+ }
+
+ static const xAOD::Vertex::Decorator<int> truthTypeDecor("truthVertexMatchType");
+ for(const auto& truthVtx : truthVerticesToMatch) {
+ int truthTypeBitset = truthTypeDecor(*truthVtx);
+ if(InDetSecVtxTruthMatchUtils::isReconstructable(truthTypeBitset)) {
+ fillTruthHistograms(truthVtx, "Reconstructable");
+
+ // fill efficiencies
+ efficiency("Acceptance")->Fill(InDetSecVtxTruthMatchUtils::isAccepted(truthTypeBitset), truthVtx->perp());
+ efficiency("eff_total")->Fill(InDetSecVtxTruthMatchUtils::isReconstructed(truthTypeBitset), truthVtx->perp());
+ }
+ if(InDetSecVtxTruthMatchUtils::isAccepted(truthTypeBitset)) {
+ fillTruthHistograms(truthVtx, "Accepted");
+ efficiency("eff_seed")->Fill(InDetSecVtxTruthMatchUtils::isSeeded(truthTypeBitset), truthVtx->perp());
+ }
+ if(InDetSecVtxTruthMatchUtils::isSeeded(truthTypeBitset)) {
+ fillTruthHistograms(truthVtx, "Seeded");
+ efficiency("eff_core")->Fill(InDetSecVtxTruthMatchUtils::isReconstructed(truthTypeBitset), truthVtx->perp());
+ }
+ if(InDetSecVtxTruthMatchUtils::isReconstructed(truthTypeBitset)) {
+ fillTruthHistograms(truthVtx, "Reconstructed");
+ }
+ if(InDetSecVtxTruthMatchUtils::isReconstructedSplit(truthTypeBitset)) {
+ fillTruthHistograms(truthVtx, "ReconstructedSplit");
+ }
+ fillTruthHistograms(truthVtx, "Inclusive");
+
+ }
+ }
+
+ return StatusCode::SUCCESS;
+
+ }
+ void SecVertexTruthMatchAlg::fillRecoHistograms(const xAOD::Vertex* secVtx, std::string matchType) {
+
+ // set of accessors for tracks and weights
+ xAOD::Vertex::ConstAccessor<xAOD::Vertex::TrackParticleLinks_t> trkAcc("trackParticleLinks");
+ xAOD::Vertex::ConstAccessor<std::vector<float> > weightAcc("trackWeights");
+
+ // set of decorators for truth matching info
+ const xAOD::Vertex::Decorator<std::vector<InDetSecVtxTruthMatchUtils::VertexTruthMatchInfo> > matchInfoDecor("truthVertexMatchingInfos");
+
+ TVector3 reco_pos(secVtx->x(), secVtx->y(), secVtx->z());
+ float Lxy = reco_pos.Perp();
+
+ size_t ntracks;
+ const xAOD::Vertex::TrackParticleLinks_t & trkParts = trkAcc( *secVtx );
+ ntracks = trkParts.size();
+
+ TLorentzVector sumP4(0,0,0,0);
+ double H = 0.0;
+ double HT = 0.0;
+ int charge = 0;
+ double minOpAng = -1.0* 1.e10;
+ double maxOpAng = 1.0* 1.e10;
+ double minD0 = 1.0* 1.e10;
+ double maxD0 = 0.0;
+ double maxDR = 0.0;
+
+ xAOD::TrackParticle::ConstAccessor< std::vector< float > > accCovMatrixDiag( "definingParametersCovMatrixDiag" );
+
+ ATH_MSG_DEBUG("Loop over tracks");
+ for(size_t t = 0; t < ntracks; t++){
+ if(!trkParts[t].isValid()){
+ ATH_MSG_DEBUG("Track " << t << " is bad!");
+ continue;
+ }
+ const xAOD::TrackParticle & trk = **trkParts[t];
+
+ double trk_d0 = std::abs(trk.definingParameters()[0]);
+ double trk_z0 = std::abs(trk.definingParameters()[1]);
+
+ if(trk_d0 < minD0){ minD0 = trk_d0; }
+ if(trk_d0 > maxD0){ maxD0 = trk_d0; }
+
+ TLorentzVector vv;
+ // TODO: use values computed w.r.t SV
+ vv.SetPtEtaPhiM(trk.pt(),trk.eta(), trk.phi0(), trk.m());
+ sumP4 += vv;
+ H += vv.Vect().Mag();
+ HT += vv.Pt();
+
+ TLorentzVector v_minus_iv(0,0,0,0);
+ for(size_t j = 0; j < ntracks; j++){
+ if (j == t){ continue; }
+ if(!trkParts[j].isValid()){
+ ATH_MSG_DEBUG("Track " << j << " is bad!");
+ continue;
+ }
+
+ const xAOD::TrackParticle & trk_2 = **trkParts[j];
+
+ TLorentzVector tmp;
+ // TODO: use values computed w.r.t. SV
+ tmp.SetPtEtaPhiM(trk_2.pt(),trk_2.eta(), trk_2.phi0(), trk_2.m());
+ v_minus_iv += tmp;
+
+ if( j > t ) {
+ double tm = vv * tmp / ( vv.Mag() * tmp.Mag() );
+ if( minOpAng < tm ) minOpAng = tm;
+ if( maxOpAng > tm ) maxOpAng = tm;
+ }
+ }
+ double DR = vv.DeltaR(v_minus_iv);
+ if( DR > maxDR ){ maxDR = DR;}
+
+ charge += trk.charge();
+
+ xAOD::TrackParticle::ConstAccessor<float> Trk_Chi2("chiSquared");
+ xAOD::TrackParticle::ConstAccessor<float> Trk_nDoF("numberDoF");
+
+ if ( Trk_Chi2.isAvailable(trk) && Trk_Chi2(trk) && Trk_nDoF.isAvailable(trk) && Trk_nDoF(trk) ) {
+ hist("RecoVertex/" + matchType + "_Trk_Chi2")->Fill(Trk_Chi2(trk) / Trk_nDoF(trk));
+ hist("RecoVertex/" + matchType + "_Trk_nDoF")->Fill(Trk_nDoF(trk));
+ }
+ hist("RecoVertex/" + matchType + "_Trk_D0")->Fill(trk_d0);
+ hist("RecoVertex/" + matchType + "_Trk_Z0")->Fill(trk_z0);
+ hist("RecoVertex/" + matchType + "_Trk_theta")->Fill(trk.definingParameters()[3]);
+ hist("RecoVertex/" + matchType + "_Trk_qOverP")->Fill(trk.definingParameters()[4]);
+ hist("RecoVertex/" + matchType + "_Trk_Eta")->Fill(trk.eta());
+ hist("RecoVertex/" + matchType + "_Trk_Phi")->Fill(trk.phi0());
+ hist("RecoVertex/" + matchType + "_Trk_E")->Fill(trk.e() / GeV);
+ hist("RecoVertex/" + matchType + "_Trk_M")->Fill(trk.m() / GeV);
+ hist("RecoVertex/" + matchType + "_Trk_Pt")->Fill(trk.pt() / GeV);
+ hist("RecoVertex/" + matchType + "_Trk_Px")->Fill(trk.p4().Px() / GeV);
+ hist("RecoVertex/" + matchType + "_Trk_Py")->Fill(trk.p4().Py() / GeV);
+ hist("RecoVertex/" + matchType + "_Trk_Pz")->Fill(trk.p4().Pz() / GeV);
+ hist("RecoVertex/" + matchType + "_Trk_charge")->Fill(trk.charge());
+ hist("RecoVertex/" + matchType + "_Trk_errD0")->Fill(trk.definingParametersCovMatrix()(0,0));
+ hist("RecoVertex/" + matchType + "_Trk_errZ0")->Fill(trk.definingParametersCovMatrix()(1,1));
+ } // end loop over tracks
+
+ const double dir = sumP4.Vect().Dot( reco_pos ) / sumP4.Vect().Mag() / reco_pos.Mag();
+
+ xAOD::Vertex::ConstAccessor<float> Chi2("chiSquared");
+ xAOD::Vertex::ConstAccessor<float> nDoF("numberDoF");
+
+ hist("RecoVertex/" + matchType + "_x")->Fill(secVtx->x());
+ hist("RecoVertex/" + matchType + "_y")->Fill(secVtx->y());
+ hist("RecoVertex/" + matchType + "_z")->Fill(secVtx->z());
+ hist("RecoVertex/" + matchType + "_Lxy")->Fill(Lxy);
+ hist("RecoVertex/" + matchType + "_ntrk")->Fill(ntracks);
+ hist("RecoVertex/" + matchType + "_pT")->Fill(sumP4.Pt() / GeV);
+ hist("RecoVertex/" + matchType + "_eta")->Fill(sumP4.Eta());
+ hist("RecoVertex/" + matchType + "_phi")->Fill(sumP4.Phi());
+ hist("RecoVertex/" + matchType + "_mass")->Fill(sumP4.M() / GeV);
+ hist("RecoVertex/" + matchType + "_mu")->Fill(sumP4.M()/maxDR / GeV);
+ hist("RecoVertex/" + matchType + "_chi2")->Fill(Chi2(*secVtx)/nDoF(*secVtx));
+ hist("RecoVertex/" + matchType + "_dir")->Fill(dir);
+ hist("RecoVertex/" + matchType + "_charge")->Fill(charge);
+ hist("RecoVertex/" + matchType + "_H")->Fill(H / GeV);
+ hist("RecoVertex/" + matchType + "_HT")->Fill(HT / GeV);
+ hist("RecoVertex/" + matchType + "_minOpAng")->Fill(minOpAng);
+ hist("RecoVertex/" + matchType + "_maxOpAng")->Fill(maxOpAng);
+ hist("RecoVertex/" + matchType + "_mind0")->Fill(minD0);
+ hist("RecoVertex/" + matchType + "_maxd0")->Fill(maxD0);
+ hist("RecoVertex/" + matchType + "_maxdR")->Fill(maxDR);
+
+ std::vector<InDetSecVtxTruthMatchUtils::VertexTruthMatchInfo> truthmatchinfo;
+ truthmatchinfo = matchInfoDecor(*secVtx);
+
+ // This includes all matched vertices, including splits
+ if(not truthmatchinfo.empty()){
+ float matchScore_weight = std::get<1>(truthmatchinfo.at(0));
+ float matchScore_pt = std::get<2>(truthmatchinfo.at(0));
+
+ ATH_MSG_DEBUG("Match Score and probability: " << matchScore_weight << " " << matchScore_pt/0.01);
+
+ const ElementLink<xAOD::TruthVertexContainer>& truthVertexLink = std::get<0>(truthmatchinfo.at(0));
+ const xAOD::TruthVertex& truthVtx = **truthVertexLink ;
+
+ hist("RecoVertex/" + matchType + "_positionRes_R")->Fill(Lxy - truthVtx.perp());
+ hist("RecoVertex/" + matchType + "_positionRes_Z")->Fill(secVtx->z() - truthVtx.z());
+ hist("RecoVertex/" + matchType + "_matchScore_weight")->Fill(matchScore_weight);
+ hist("RecoVertex/" + matchType + "_matchScore_pt")->Fill(matchScore_pt);
+ }
+ }
+
+ void SecVertexTruthMatchAlg::fillTruthHistograms(const xAOD::TruthVertex* truthVtx, std::string truthType) {
+
+ hist("TruthVertex/" + truthType + "_x")->Fill(truthVtx->x());
+ hist("TruthVertex/" + truthType + "_y")->Fill(truthVtx->y());
+ hist("TruthVertex/" + truthType + "_z")->Fill(truthVtx->z());
+ hist("TruthVertex/" + truthType + "_R")->Fill(truthVtx->perp());
+ hist("TruthVertex/" + truthType + "_Eta")->Fill(truthVtx->eta());
+ hist("TruthVertex/" + truthType + "_Phi")->Fill(truthVtx->phi());
+ hist("TruthVertex/" + truthType + "_Ntrk_out")->Fill(truthVtx->nOutgoingParticles());
+
+ ATH_MSG_DEBUG("Plotting truth parent");
+ const xAOD::TruthParticle& truthPart = *truthVtx->incomingParticle(0);
+
+ hist("TruthVertex/" + truthType + "_Parent_E")->Fill(truthPart.e() / GeV);
+ hist("TruthVertex/" + truthType + "_Parent_M")->Fill(truthPart.m() / GeV);
+ hist("TruthVertex/" + truthType + "_Parent_Pt")->Fill(truthPart.pt() / GeV);
+ hist("TruthVertex/" + truthType + "_Parent_Phi")->Fill(truthPart.phi());
+ hist("TruthVertex/" + truthType + "_Parent_Eta")->Fill(truthPart.eta());
+ hist("TruthVertex/" + truthType + "_Parent_charge")->Fill(truthPart.charge());
+
+ ATH_MSG_DEBUG("Plotting truth prod vtx");
+ if(truthPart.hasProdVtx()){
+ const xAOD::TruthVertex & vertex = *truthPart.prodVtx();
+
+ hist("TruthVertex/" + truthType + "_ParentProdX")->Fill(vertex.x());
+ hist("TruthVertex/" + truthType + "_ParentProdY")->Fill(vertex.y());
+ hist("TruthVertex/" + truthType + "_ParentProdZ")->Fill(vertex.z());
+ hist("TruthVertex/" + truthType + "_ParentProdR")->Fill(vertex.perp());
+ hist("TruthVertex/" + truthType + "_ParentProdEta")->Fill(vertex.eta());
+ hist("TruthVertex/" + truthType + "_ParentProdPhi")->Fill(vertex.phi());
+ }
+ // m_truthInclusive_r->Fill(truthVtx.perp());
+
+ // if(matchTypeDecor(truthVtx) >= RECONSTRUCTABLE){
+ // m_truthReconstructable_r->Fill(truthVtx.perp());
+ // }
+ // if(matchTypeDecor(truthVtx) >= ACCEPTED){
+ // m_truthAccepted_r->Fill(truthVtx.perp());
+ // }
+ // if(matchTypeDecor(truthVtx) >= SEEDED){
+ // m_truthSeeded_r->Fill(truthVtx.perp());
+ // }
+ // if(matchTypeDecor(truthVtx) >= RECONSTRUCTED){
+ // m_truthReconstructed_r->Fill(truthVtx.perp());
+ // }
+ // if(matchTypeDecor(truthVtx) >= RECONSTRUCTEDSPLIT){
+ // m_truthSplit_r->Fill(truthVtx.perp());
+ // }
+ //
+ }
+} // namespace CP
diff --git a/PhysicsAnalysis/Algorithms/TrackingAnalysisAlgorithms/TrackingAnalysisAlgorithms/SecVertexTruthMatchAlg.h b/PhysicsAnalysis/Algorithms/TrackingAnalysisAlgorithms/TrackingAnalysisAlgorithms/SecVertexTruthMatchAlg.h
new file mode 100644
index 0000000000000000000000000000000000000000..4e70b277bc0f2b301b6c21bc1cc219955ae200c8
--- /dev/null
+++ b/PhysicsAnalysis/Algorithms/TrackingAnalysisAlgorithms/TrackingAnalysisAlgorithms/SecVertexTruthMatchAlg.h
@@ -0,0 +1,61 @@
+/*
+ Copyright (C) 2002-2021 CERN for the benefit of the ATLAS collaboration
+*/
+
+#ifndef TRACKINGANALYSISALGORITHMS_SECVERTEXTRUTHMATCHALG_H
+#define TRACKINGANALYSISALGORITHMS_SECVERTEXTRUTHMATCHALG_H
+
+// Gaudi/Athena include(s):
+#include "AnaAlgorithm/AnaAlgorithm.h"
+#include "AsgTools/ToolHandle.h"
+
+// SG includes
+#include "xAODTracking/VertexContainer.h"
+#include "xAODTracking/TrackParticleContainer.h"
+#include "xAODTruth/TruthVertexContainer.h"
+
+// Local include(s):
+#include "InDetSecVtxTruthMatchTool/IInDetSecVtxTruthMatchTool.h"
+
+// Asg includes
+#include <AsgTools/PropertyWrapper.h>
+#include <AsgDataHandles/ReadHandleKey.h>
+#include <AsgDataHandles/ReadHandle.h>
+
+namespace CP {
+
+ /// Algorithm to perform truth matching on secondary vertices
+ ///
+ /// @author Jackson Burzynski <jackson.carl.burzynski@cern.ch>
+ ///
+ class SecVertexTruthMatchAlg final : public EL::AnaAlgorithm {
+
+ public:
+ /// Regular Algorithm constructor
+ SecVertexTruthMatchAlg( const std::string& name, ISvcLocator* svcLoc );
+ virtual StatusCode initialize() override;
+ virtual StatusCode execute() override;
+
+ private:
+ // Input (reco) Secondary Vertices
+ SG::ReadHandleKey<xAOD::VertexContainer> m_secVtxContainerKey{this, "SecondaryVertexContainer", "VrtSecInclusive_SecondaryVertices",
+ "Secondary vertex container"};
+ // Input Truth Vertices
+ SG::ReadHandleKey<xAOD::TruthVertexContainer> m_truthVtxContainerKey{this, "TruthVertexContainer", "TruthVertices",
+ "Truth vertex container"};
+ // Input Tracks
+ SG::ReadHandleKey<xAOD::TrackParticleContainer> m_trackParticleContainerKey{this, "TrackParticleContainer", "InDetTrackParticles",
+ "Track container"};
+
+ Gaudi::Property<std::vector<int>> m_targetPDGIDs{this, "TargetPDGIDs", {}, "List of PDGIDs to select for matching"};
+
+ Gaudi::Property<bool> m_writeHistograms{this, "WriteHistograms", true, "Write histograms"};
+
+ ToolHandle<IInDetSecVtxTruthMatchTool> m_matchTool{this, "MatchTool", "InDetSecVtxTruthMatchTool"};
+
+ void fillRecoHistograms(const xAOD::Vertex* secVtx, std::string matchType);
+ void fillTruthHistograms(const xAOD::TruthVertex* truthVtx, std::string truthType);
+
+ };
+} // namespace CP
+#endif
diff --git a/PhysicsAnalysis/Algorithms/TrackingAnalysisAlgorithms/TrackingAnalysisAlgorithms/TrackingAnalysisAlgorithmsDict.h b/PhysicsAnalysis/Algorithms/TrackingAnalysisAlgorithms/TrackingAnalysisAlgorithms/TrackingAnalysisAlgorithmsDict.h
index b1168f430529624e4f03d2a3d787269fd712ed3f..40f074f015a9c2b5fa6a18050d2e62d18431baad 100644
--- a/PhysicsAnalysis/Algorithms/TrackingAnalysisAlgorithms/TrackingAnalysisAlgorithms/TrackingAnalysisAlgorithmsDict.h
+++ b/PhysicsAnalysis/Algorithms/TrackingAnalysisAlgorithms/TrackingAnalysisAlgorithms/TrackingAnalysisAlgorithmsDict.h
@@ -1,6 +1,6 @@
// Dear emacs, this is -*- c++ -*-
//
-// Copyright (C) 2002-2023 CERN for the benefit of the ATLAS collaboration
+// Copyright (C) 2002-2024 CERN for the benefit of the ATLAS collaboration
//
#ifndef TRACKINGANALYSISALGORITHMS_TRACKINGANALYSISALGORITHMSDICT_H
#define TRACKINGANALYSISALGORITHMS_TRACKINGANALYSISALGORITHMSDICT_H
@@ -8,5 +8,6 @@
// Local include(s):
#include "TrackingAnalysisAlgorithms/VertexSelectionAlg.h"
#include "TrackingAnalysisAlgorithms/TrackParticleMergerAlg.h"
+#include "TrackingAnalysisAlgorithms/SecVertexTruthMatchAlg.h"
#endif // TRACKINGANALYSISALGORITHMS_TRACKINGANALYSISALGORITHMSDICT_H
diff --git a/PhysicsAnalysis/Algorithms/TrackingAnalysisAlgorithms/TrackingAnalysisAlgorithms/selection.xml b/PhysicsAnalysis/Algorithms/TrackingAnalysisAlgorithms/TrackingAnalysisAlgorithms/selection.xml
index 0713cd73938555d2f016d02620f1c70c5bdfb1e6..a86b756c4f14ce0914f09ad84c22e1e6b70b8be5 100644
--- a/PhysicsAnalysis/Algorithms/TrackingAnalysisAlgorithms/TrackingAnalysisAlgorithms/selection.xml
+++ b/PhysicsAnalysis/Algorithms/TrackingAnalysisAlgorithms/TrackingAnalysisAlgorithms/selection.xml
@@ -1,7 +1,8 @@
-<!-- Copyright (C) 2002-2023 CERN for the benefit of the ATLAS collaboration -->
+<!-- Copyright (C) 2002-2024 CERN for the benefit of the ATLAS collaboration -->
<lcgdict>
<class name="CP::VertexSelectionAlg" />
<class name="CP::TrackParticleMergerAlg" />
+ <class name="CP::SecVertexTruthMatchAlg" />
</lcgdict>
diff --git a/PhysicsAnalysis/Algorithms/TrackingAnalysisAlgorithms/python/TrackingAnalysisAlgorithmsConfig.py b/PhysicsAnalysis/Algorithms/TrackingAnalysisAlgorithms/python/TrackingAnalysisAlgorithmsConfig.py
new file mode 100644
index 0000000000000000000000000000000000000000..d5c56f4f43e069e67641effbdbbf6aa41b95bf0e
--- /dev/null
+++ b/PhysicsAnalysis/Algorithms/TrackingAnalysisAlgorithms/python/TrackingAnalysisAlgorithmsConfig.py
@@ -0,0 +1,31 @@
+# Copyright (C) 2002-2024 CERN for the benefit of the ATLAS collaboration
+
+from AthenaConfiguration.ComponentFactory import CompFactory
+from AthenaConfiguration.ComponentAccumulator import ComponentAccumulator
+
+def InDetSecVtxTruthMatchToolCfg(flags, name="InDetSecVtxTruthMatchTool", **kwargs):
+ acc = ComponentAccumulator()
+ acc.setPrivateTools(CompFactory.InDetSecVtxTruthMatchTool(**kwargs))
+ return acc
+
+def SecVertexTruthMatchAlgCfg(flags, name="SecVertexTruthMatchAlg", useLRTTracks = False, **kwargs):
+
+ acc = ComponentAccumulator()
+
+ if useLRTTracks:
+ from DerivationFrameworkInDet.InDetToolsConfig import InDetLRTMergeCfg
+ acc.merge(InDetLRTMergeCfg(flags))
+ kwargs.setdefault("TrackParticleContainer", "InDetWithLRTTrackParticles")
+
+ kwargs.setdefault("TruthVertexContainer", "TruthVertices")
+ kwargs.setdefault("SecondaryVertexContainer", "VrtSecInclusive_SecondaryVertices")
+ kwargs.setdefault("TargetPDGIDs", [511,521])
+ kwargs.setdefault("MatchTool", acc.popToolsAndMerge(InDetSecVtxTruthMatchToolCfg(flags)))
+
+ truthMatchAlg = CompFactory.CP.SecVertexTruthMatchAlg(name, **kwargs)
+
+ acc.addEventAlgo(truthMatchAlg)
+ acc.addService(CompFactory.THistSvc(Output = [f"ANALYSIS DATAFILE='{flags.Output.HISTFileName}' OPT='RECREATE'"]))
+ acc.setAppProperty("HistogramPersistency","ROOT")
+ return acc
+
diff --git a/PhysicsAnalysis/Algorithms/TrackingAnalysisAlgorithms/scripts/runSecVtxTruthMatching.py b/PhysicsAnalysis/Algorithms/TrackingAnalysisAlgorithms/scripts/runSecVtxTruthMatching.py
new file mode 100755
index 0000000000000000000000000000000000000000..1b022212f93e892b1e99305f4376716af0c54fbc
--- /dev/null
+++ b/PhysicsAnalysis/Algorithms/TrackingAnalysisAlgorithms/scripts/runSecVtxTruthMatching.py
@@ -0,0 +1,57 @@
+#!/usr/bin/env python
+# Copyright (C) 2002-2024 CERN for the benefit of the ATLAS collaboration
+
+from glob import glob
+
+def get_args():
+ from argparse import ArgumentParser
+ parser = ArgumentParser(description='Parser for SecVertexTruthMatching configuration')
+ parser.add_argument("--filesInput", required=True)
+ parser.add_argument("--maxEvents", help="Limit number of events. Default: all input events", default=-1, type=int)
+ parser.add_argument("--skipEvents", help="Skip this number of events. Default: no events are skipped", default=0, type=int)
+ parser.add_argument("--mergeLargeD0Tracks", help='Consider LRT tracks in the matching', action='store_true', default=False)
+ parser.add_argument("--outputFile", help='Name of output file',default="TruthMatchHists.root")
+ parser.add_argument("--pdgIds", help='List of pdgIds to match', nargs='+', type=int, default=[36,51])
+ parser.add_argument("--vertexContainer", help='SG key of secondary vertex container',default='VrtSecInclusive_SecondaryVertices')
+ parser.add_argument("--truthVertexContainer", help='SG key of truth vertex container',default='TruthVertices')
+ return parser.parse_args()
+
+if __name__=='__main__':
+
+ args = get_args()
+
+ from AthenaConfiguration.AllConfigFlags import initConfigFlags
+ flags = initConfigFlags()
+
+ flags.Input.Files = []
+ for path in args.filesInput.split(','):
+ flags.Input.Files += glob(path)
+ flags.Output.HISTFileName = args.outputFile
+
+ flags.Exec.SkipEvents = args.skipEvents
+ flags.Exec.MaxEvents = args.maxEvents
+
+ flags.lock()
+
+ from AthenaConfiguration.MainServicesConfig import MainServicesCfg
+ acc = MainServicesCfg(flags)
+ from AthenaPoolCnvSvc.PoolReadConfig import PoolReadCfg
+ acc.merge(PoolReadCfg(flags))
+
+ from TrackingAnalysisAlgorithms.TrackingAnalysisAlgorithmsConfig import SecVertexTruthMatchAlgCfg
+ acc.merge(SecVertexTruthMatchAlgCfg(flags,
+ useLRTTracks = args.mergeLargeD0Tracks,
+ TargetPDGIDs = args.pdgIds,
+ SecondaryVertexContainer = args.vertexContainer,
+ TruthVertexContainer = args.truthVertexContainer
+ )
+ )
+
+ acc.printConfig(withDetails=True)
+
+ # Execute and finish
+ sc = acc.run()
+
+ # Success should be 0
+ import sys
+ sys.exit(not sc.isSuccess())
diff --git a/PhysicsAnalysis/Algorithms/TrackingAnalysisAlgorithms/src/components/TrackingAnalysisAlgorithms_entries.cxx b/PhysicsAnalysis/Algorithms/TrackingAnalysisAlgorithms/src/components/TrackingAnalysisAlgorithms_entries.cxx
index f254ab815a274580e695e22a7ba2028c973b6f53..2af979ddd3710f0861f907a97bbcc4f8f3909697 100644
--- a/PhysicsAnalysis/Algorithms/TrackingAnalysisAlgorithms/src/components/TrackingAnalysisAlgorithms_entries.cxx
+++ b/PhysicsAnalysis/Algorithms/TrackingAnalysisAlgorithms/src/components/TrackingAnalysisAlgorithms_entries.cxx
@@ -5,9 +5,11 @@
// Local include(s):
#include "TrackingAnalysisAlgorithms/VertexSelectionAlg.h"
#include "TrackingAnalysisAlgorithms/TrackParticleMergerAlg.h"
+#include "TrackingAnalysisAlgorithms/SecVertexTruthMatchAlg.h"
// Declare the component(s) of the package:
DECLARE_COMPONENT( CP::VertexSelectionAlg )
DECLARE_COMPONENT( CP::TrackParticleMergerAlg )
+DECLARE_COMPONENT( CP::SecVertexTruthMatchAlg )
diff --git a/Projects/AnalysisBase/package_filters.txt b/Projects/AnalysisBase/package_filters.txt
index 06a700c85111361ded0a5b92ab6f1c58c349d6da..15a47ef25a8b607c78d0b018c7591280cda44205 100644
--- a/Projects/AnalysisBase/package_filters.txt
+++ b/Projects/AnalysisBase/package_filters.txt
@@ -55,6 +55,7 @@
+ Generators/TruthUtils
+ InnerDetector/InDetRecTools/InDetTrackSelectionTool
+ InnerDetector/InDetRecTools/TrackVertexAssociationTool
++ InnerDetector/InDetRecTools/InDetSecVtxTruthMatchTool
+ MuonSpectrometer/MuonStationIndex
+ PhysicsAnalysis/Algorithms/.*
+ PhysicsAnalysis/AnalysisCommon/AssociationUtils
diff --git a/Projects/AthAnalysis/package_filters.txt b/Projects/AthAnalysis/package_filters.txt
index a4287407eb2104b09d9a966ead6cd92ee7371cf0..206a9134fde7b2f4bcabb0909969cd126cd06610 100644
--- a/Projects/AthAnalysis/package_filters.txt
+++ b/Projects/AthAnalysis/package_filters.txt
@@ -50,6 +50,7 @@
+ Generators/TruthUtils
+ InnerDetector/InDetRecTools/InDetTrackSelectionTool
+ InnerDetector/InDetRecTools/TrackVertexAssociationTool
++ InnerDetector/InDetRecTools/InDetSecVtxTruthMatchTool
+ MuonSpectrometer/MuonStationIndex
- PhysicsAnalysis/Algorithms/StandaloneAnalysisAlgorithms
+ PhysicsAnalysis/Algorithms/.*