Add HS Vertex Reco/Selection Efficiency and Reso Plots to IDPVM

Attempts to add Figures 8-12 from ATL-PHYS-PUB-2019-015. The scheme for each is:

• Selection efficiency: determine the "best" reco HS vertex using the sumpt2 of reco tracks and using the sum of truth HS track weights (InDetVertexTruthMatchUtils::bestHardScatterMatch?) and see if they match (calculate local PU density using the truth vertex for the sumpt2 best reco HS vertex)
• Reconstruction efficiency: find the truth vertex for every reco vertex and check if it is present in our vector of truth HS vertices for the event - if so, we count this as a "pass"
• HS resolutions: for each vertex in our vector of truth HS vertices (I think we only ever expect 1), calculate the "truth" resolution as the absolute difference between that truth vertex and its matched reco vertex longitudinally and transversely (this is what's in the note) as well as the "reco" resolution using the covariance matrix and error propagation (which I hope I did correctly...)

Additionally, I have also fixed a seg fault which shows up when trying to fill histograms in InDetPerfPlot_VertexTruthMatching.cxx when m_iDetailLevel >= 200 by adding iDetailLevel to the constructor of InDetRttPlots.

Tagging @npetters @goblirsc @vcairo @kostyuk @glee

InnerDetector/InDetValidation/InDetPhysValMonitoring/src/InDetPerfPlot_VertexTruthMatching.cxx

@@ -87,19 +87,24 @@ float InDetPerfPlot_VertexTruthMatching::getLocalPUDensity(const xAOD::TruthVert
     float x0 = vtxOfInterest->x();
     float y0 = vtxOfInterest->y();
     float z0 = vtxOfInterest->z();
+    float radialWindow2 = std::pow(radialWindow, 2);
     int nTracksInWindow = 0;
     float localPUDensity;
-    float radialDiff;
+    float radialDiff2;
     for (const auto& vtx : truthHSVertices) {
-        radialDiff = std::sqrt(std::pow((vtx->x() - x0), 2) + std::pow((vtx->y() - y0), 2) + std::pow((vtx->z() - z0), 2));
-        if ((vtx != vtxOfInterest) && (radialDiff < radialWindow)) {
-            nTracksInWindow += 1;
+        if (vtx != vtxOfInterest) {
+            radialDiff2 = std::pow((vtx->x() - x0), 2) + std::pow((vtx->y() - y0), 2) + std::pow((vtx->z() - z0), 2);
+            if (radialDiff2 < radialWindow2) {
+                nTracksInWindow += 1;
+            }
         }
     }
     for (const auto& vtx : truthPUVertices) {
-        radialDiff = std::sqrt(std::pow((vtx->x() - x0), 2) + std::pow((vtx->y() - y0), 2) + std::pow((vtx->z() - z0), 2));
-        if ((vtx != vtxOfInterest) && (radialDiff < radialWindow)) {
-            nTracksInWindow += 1;
+        if (vtx != vtxOfInterest) {
+            radialDiff2 = std::pow((vtx->x() - x0), 2) + std::pow((vtx->y() - y0), 2) + std::pow((vtx->z() - z0), 2);
+            if (radialDiff2 < radialWindow2) {
+                nTracksInWindow += 1;
+            }
         }
     }
     localPUDensity = (float)(nTracksInWindow) / (2 * radialWindow);
@@ -130,33 +135,25 @@ float InDetPerfPlot_VertexTruthMatching::getTruthTransverseReso(const xAOD::Vert
 
 const xAOD::TruthVertex* InDetPerfPlot_VertexTruthMatching::getTruthVertex(const xAOD::Vertex* recoVtx) const {
     const xAOD::TruthVertex* truthVtx = nullptr;
-    const static xAOD::Vertex::Decorator<std::vector<InDetVertexTruthMatchUtils::VertexTruthMatchInfo>> truthMatchingInfos("TruthEventMatchingInfos");
-    const std::vector<InDetVertexTruthMatchUtils::VertexTruthMatchInfo>& truthInfos = truthMatchingInfos(*recoVtx);
-    if (!truthInfos.empty()) {
-        const InDetVertexTruthMatchUtils::VertexTruthMatchInfo& truthInfo = truthInfos.at(0);
-        const ElementLink<xAOD::TruthEventBaseContainer> truthEventLink = std::get<0>(truthInfo);
-        const xAOD::TruthEvent* truthEvent = nullptr;
-        if (truthEventLink.isValid()) {
-            truthEvent = (const xAOD::TruthEvent*)(*truthEventLink);
-            if (truthEvent) {
-                truthVtx = truthEvent->truthVertex(0);
+    if (recoVtx) {
+        const static xAOD::Vertex::Decorator<std::vector<InDetVertexTruthMatchUtils::VertexTruthMatchInfo>> truthMatchingInfos("TruthEventMatchingInfos");
+        const std::vector<InDetVertexTruthMatchUtils::VertexTruthMatchInfo>& truthInfos = truthMatchingInfos(*recoVtx);
+        if (!truthInfos.empty()) {
+            const InDetVertexTruthMatchUtils::VertexTruthMatchInfo& truthInfo = truthInfos.at(0);
+            const ElementLink<xAOD::TruthEventBaseContainer> truthEventLink = std::get<0>(truthInfo);
+            const xAOD::TruthEvent* truthEvent = nullptr;
+            if (truthEventLink.isValid()) {
+                truthEvent = static_cast<const xAOD::TruthEvent*>(*truthEventLink);
+                if (truthEvent) {
+                    truthVtx = truthEvent->truthVertex(0);
+                }
             }
         }
-    }
-    else {
-        ATH_MSG_WARNING("TruthEventMatchingInfos DECORATOR yields empty vector -- returning nullptr!");
-    }
-    return truthVtx;
-}
-
-template <typename T>
-int InDetPerfPlot_VertexTruthMatching::inIterable(const T* ptrOfInterest, const std::vector<const T*>& iterable) const {
-    for (size_t i = 0; i < iterable.size(); i++) {
-        if (ptrOfInterest == iterable[i]) {
-            return i;
+        else {
+            ATH_MSG_WARNING("TruthEventMatchingInfos DECORATOR yields empty vector -- returning nullptr!");
         }
     }
-    return -1;
+    return truthVtx;
 }
 
 void InDetPerfPlot_VertexTruthMatching::fill(const xAOD::Vertex& vertex) {
@@ -184,8 +181,9 @@ void InDetPerfPlot_VertexTruthMatching::fill(const xAOD::VertexContainer& vertex
 
     // Fill our histograms
     // Inclusive:
-    float nTruthVertices = (float)(truthHSVertices.size() + truthPUVertices.size());
-    fillHisto(m_vx_nReco_vs_nTruth_inclusive, nTruthVertices, (float)vertexContainer.size());
+    int nTruthVertices = (int)(truthHSVertices.size() + truthPUVertices.size());
+    int nRecoVertices = (int)vertexContainer.size();
+    fillHisto(m_vx_nReco_vs_nTruth_inclusive, nTruthVertices, nRecoVertices);
     
     // Let's also plot the vertices by vertex match type:
     const static xAOD::Vertex::Decorator<InDetVertexTruthMatchUtils::VertexMatchType> recoVtxMatchTypeInfo("VertexMatchType");
@@ -209,11 +207,9 @@ void InDetPerfPlot_VertexTruthMatching::fill(const xAOD::VertexContainer& vertex
     localPUDensity = getLocalPUDensity(truthVtx, truthHSVertices, truthPUVertices);
     fillHisto(m_vx_hs_sel_eff, localPUDensity, (bestRecoHSVtx_sumpt2 == bestRecoHSVtx_truth));
 
+    // Iterate over vertices:
     InDetVertexTruthMatchUtils::VertexMatchType matchType;
-    int iTruthVtx;
     std::vector<const xAOD::Vertex*> truthHSVtxRecoMatch = std::vector<const xAOD::Vertex*>(truthHSVertices.size(), nullptr);
-
-    // Iterate over vertices:
     for (const auto& vertex : vertexContainer.stdcont()) {
 
         fill(*vertex);
@@ -224,9 +220,9 @@ void InDetPerfPlot_VertexTruthMatching::fill(const xAOD::VertexContainer& vertex
         // Get the truth vertex for the current reco vertex
         truthVtx = getTruthVertex(vertex);
         if (truthVtx) {
-            iTruthVtx = inIterable(truthVtx, truthHSVertices);
-            if (iTruthVtx != -1) {
-                truthHSVtxRecoMatch[iTruthVtx] = vertex;
+            auto it = std::find(truthHSVertices.begin(), truthHSVertices.end(), truthVtx);
+            if (it != truthHSVertices.end()) {
+                truthHSVtxRecoMatch[std::distance(truthHSVertices.begin(), it)] = vertex;
             }
         }
 
@@ -244,41 +240,40 @@ void InDetPerfPlot_VertexTruthMatching::fill(const xAOD::VertexContainer& vertex
             fillHisto(m_vx_hs_truth_long_reso, localPUDensity, getTruthLongitudinalReso(recoVtx, truthVtx));
             fillHisto(m_vx_hs_reco_trans_reso, localPUDensity, getRecoTransverseReso(recoVtx));
             fillHisto(m_vx_hs_truth_trans_reso, localPUDensity, getTruthTransverseReso(recoVtx, truthVtx));
-
         }
         else {
             fillHisto(m_vx_hs_reco_eff, localPUDensity, 0);
         }
     }
 
-    fillHisto(m_vx_nReco_vs_nTruth_matched, nTruthVertices, (float)breakdown[InDetVertexTruthMatchUtils::VertexMatchType::MATCHED]);
-    fillHisto(m_vx_nReco_vs_nTruth_merged,  nTruthVertices, (float)breakdown[InDetVertexTruthMatchUtils::VertexMatchType::MERGED]);
-    fillHisto(m_vx_nReco_vs_nTruth_split,   nTruthVertices, (float)breakdown[InDetVertexTruthMatchUtils::VertexMatchType::SPLIT]);
-    fillHisto(m_vx_nReco_vs_nTruth_fake,    nTruthVertices, (float)breakdown[InDetVertexTruthMatchUtils::VertexMatchType::FAKE]);
-    fillHisto(m_vx_nReco_vs_nTruth_dummy,   nTruthVertices, (float)breakdown[InDetVertexTruthMatchUtils::VertexMatchType::DUMMY]);
+    fillHisto(m_vx_nReco_vs_nTruth_matched, nTruthVertices, breakdown[InDetVertexTruthMatchUtils::VertexMatchType::MATCHED]);
+    fillHisto(m_vx_nReco_vs_nTruth_merged,  nTruthVertices, breakdown[InDetVertexTruthMatchUtils::VertexMatchType::MERGED]);
+    fillHisto(m_vx_nReco_vs_nTruth_split,   nTruthVertices, breakdown[InDetVertexTruthMatchUtils::VertexMatchType::SPLIT]);
+    fillHisto(m_vx_nReco_vs_nTruth_fake,    nTruthVertices, breakdown[InDetVertexTruthMatchUtils::VertexMatchType::FAKE]);
+    fillHisto(m_vx_nReco_vs_nTruth_dummy,   nTruthVertices, breakdown[InDetVertexTruthMatchUtils::VertexMatchType::DUMMY]);
     
     // And by hardscatter type:
     InDetVertexTruthMatchUtils::HardScatterType hsType = InDetVertexTruthMatchUtils::classifyHardScatter(vertexContainer);
     fillHisto(m_vx_hs_classification, hsType);
     switch (hsType) {
         case InDetVertexTruthMatchUtils::HardScatterType::CLEAN: {
-            fillHisto(m_vx_nReco_vs_nTruth_clean, nTruthVertices, (float)vertexContainer.size());
+            fillHisto(m_vx_nReco_vs_nTruth_clean, nTruthVertices, nRecoVertices);
             break;
         }
         case InDetVertexTruthMatchUtils::HardScatterType::LOWPU: {
-            fillHisto(m_vx_nReco_vs_nTruth_lowpu, nTruthVertices, (float)vertexContainer.size());
+            fillHisto(m_vx_nReco_vs_nTruth_lowpu, nTruthVertices, nRecoVertices);
             break;
         }
         case InDetVertexTruthMatchUtils::HardScatterType::HIGHPU: {
-            fillHisto(m_vx_nReco_vs_nTruth_highpu, nTruthVertices, (float)vertexContainer.size());
+            fillHisto(m_vx_nReco_vs_nTruth_highpu, nTruthVertices, nRecoVertices);
             break;
         }
         case InDetVertexTruthMatchUtils::HardScatterType::HSSPLIT: {
-            fillHisto(m_vx_nReco_vs_nTruth_hssplit, nTruthVertices, (float)vertexContainer.size());
+            fillHisto(m_vx_nReco_vs_nTruth_hssplit, nTruthVertices, nRecoVertices);
             break;
         }
         case InDetVertexTruthMatchUtils::HardScatterType::NONE: {
-            fillHisto(m_vx_nReco_vs_nTruth_none, nTruthVertices, (float)vertexContainer.size());
+            fillHisto(m_vx_nReco_vs_nTruth_none, nTruthVertices, nRecoVertices);
             break;
         }
         default: {