Implement Track-Protocluster Match Promoter (PFA0b)

src/algorithms/CMakeLists.txt

@@ -3,6 +3,7 @@ cmake_minimum_required(VERSION 3.16)
 add_subdirectory(interfaces)
 add_subdirectory(calorimetry)
 add_subdirectory(tracking)
+add_subdirectory(particle)
 add_subdirectory(pid)
 add_subdirectory(pid_lut)
 add_subdirectory(digi)

src/algorithms/calorimetry/TrackClusterMergeSplitter.h

@@ -1,101 +1,128 @@
 // SPDX-License-Identifier: LGPL-3.0-or-later
-// Copyright (C) 2024 Derek Anderson
+// Copyright (C) 2025 Derek Anderson, Dmitry Kalinkin
 
 #pragma once
 
 #include <algorithms/algorithm.h>
-#include <algorithms/geo.h>
+#include <edm4eic/CalorimeterHit.h>
+#include <edm4eic/ClusterCollection.h>
+#include <edm4eic/EDM4eicVersion.h>
 #include <edm4eic/ProtoClusterCollection.h>
+#include <edm4eic/TrackClusterMatchCollection.h>
+#if EDM4EIC_VERSION_MAJOR >= 8 && EDM4EIC_VERSION_MINOR >= 4
+#include <edm4eic/TrackProtoClusterMatchCollection.h>
+#endif
+#include <edm4eic/Track.h>
 #include <edm4eic/TrackPoint.h>
 #include <edm4eic/TrackSegmentCollection.h>
-#include <edm4hep/Vector3f.h>
-#include <podio/ObjectID.h>
+#include <algorithm>
 #include <map>
+#include <optional>
 #include <set>
 #include <string>
 #include <string_view>
 #include <vector>
 
 #include "TrackClusterMergeSplitterConfig.h"
-// for algorithm configuration
 #include "algorithms/interfaces/WithPodConfig.h"
 
 namespace eicrecon {
 
-// --------------------------------------------------------------------------
-//! Comparator struct for protoclusters
-// --------------------------------------------------------------------------
-/*! Organizes protoclusters by their ObjectID's in decreasing collection
-   *  ID first, and second by decreasing index second.
-   */
-struct CompareProto {
-
-  bool operator()(const edm4eic::ProtoCluster& lhs, const edm4eic::ProtoCluster& rhs) const {
-    if (lhs.getObjectID().collectionID == rhs.getObjectID().collectionID) {
-      return (lhs.getObjectID().index < rhs.getObjectID().index);
-    } else {
-      return (lhs.getObjectID().collectionID < rhs.getObjectID().collectionID);
-    }
-  }
-
-}; // end CompareObjectID
-
-// --------------------------------------------------------------------------
-//! Convenience types
-// --------------------------------------------------------------------------
-using VecProj       = std::vector<edm4eic::TrackPoint>;
-using VecClust      = std::vector<edm4eic::ProtoCluster>;
-using SetClust      = std::set<edm4eic::ProtoCluster, CompareProto>;
-using MapToVecProj  = std::map<edm4eic::ProtoCluster, VecProj, CompareProto>;
-using MapToVecClust = std::map<edm4eic::ProtoCluster, VecClust, CompareProto>;
-
 // --------------------------------------------------------------------------
 //! Algorithm input/output
 // --------------------------------------------------------------------------
 using TrackClusterMergeSplitterAlgorithm = algorithms::Algorithm<
-    algorithms::Input<edm4eic::ProtoClusterCollection, edm4eic::TrackSegmentCollection>,
+    algorithms::Input<edm4eic::TrackClusterMatchCollection, edm4eic::ClusterCollection,
+                      edm4eic::TrackSegmentCollection>,
+#if EDM4EIC_VERSION_MAJOR >= 8 && EDM4EIC_VERSION_MINOR >= 4
+    algorithms::Output<edm4eic::ProtoClusterCollection,
+                       edm4eic::TrackProtoClusterMatchCollection
+                       >>;
+#else
     algorithms::Output<edm4eic::ProtoClusterCollection>>;
+#endif
 
 // --------------------------------------------------------------------------
 //! Track-Based Cluster Merger/Splitter
 // --------------------------------------------------------------------------
-/*! An algorithm which takes a collection of proto-clusters, matches
-   *  track projections, and then decides to merge or split those proto-
-   *  clusters based on average E/p from simulations.
-   *
-   *  Heavily inspired by Eur. Phys. J. C (2017) 77:466
-   */
+/*! An algorithm which takes a collection of clusters, matches
+ *  track projections, and then decides to merge or split those
+ *  clusters based on average E/p from simulations.
+ *
+ *  Heavily inspired by Eur. Phys. J. C (2017) 77:466
+ */
 class TrackClusterMergeSplitter : public TrackClusterMergeSplitterAlgorithm,
                                   public WithPodConfig<TrackClusterMergeSplitterConfig> {
 
 public:
-  // ctor
+  // --------------------------------------------------------------------------
+  //! Comparator struct for object IDs
+  // --------------------------------------------------------------------------
+  /*! Organizes objects by their ObjectID's in decreasing collection
+   *  ID first, and second by decreasing index second.
+   */
+  template <typename T> struct CompareObjectID {
+    bool operator()(const T& lhs, const T& rhs) const {
+      if (lhs.getObjectID().collectionID == rhs.getObjectID().collectionID) {
+        return (lhs.getObjectID().index < rhs.getObjectID().index);
+      } else {
+        return (lhs.getObjectID().collectionID < rhs.getObjectID().collectionID);
+      }
+    }
+  }; // end CompareObjectID
+
+  ///! Specialization of comparator for clusters
+  using CompareClust = CompareObjectID<edm4eic::Cluster>;
+
+  ///! Specialization of comparator for hits
+  using CompareHit = CompareObjectID<edm4eic::CalorimeterHit>;
+
+  ///! Alias for vectors of track segments
+  using VecSeg = std::vector<edm4eic::TrackSegment>;
+
+  ///! Alias for vectors of mutable protoclusters
+  using VecProto = std::vector<edm4eic::MutableProtoCluster>;
+
+  ///! Alias for vectors of clusters
+  using VecClust = std::vector<edm4eic::Cluster>;
+
+  ///! Alias for sets of clusters
+  using SetClust = std::set<edm4eic::Cluster, CompareClust>;
+
+  ///! Alias for a map of clusters to the segments of matched tracks
+  using MapToVecSeg = std::map<edm4eic::Cluster, VecSeg, CompareClust>;
+
+  ///! Alias for a map of clusters onto clusters to merge
+  using MapToVecClust = std::map<edm4eic::Cluster, VecClust, CompareClust>;
+
+  ///! Alias for a map of hits onto their splitting weights
+  using MapToWeight = std::map<edm4eic::CalorimeterHit, double, CompareHit>;
+
+  ///! Alias for a vector of weight maps
+  using VecWeights = std::vector<MapToWeight>;
+
+  ///! Algorithm constructor
   TrackClusterMergeSplitter(std::string_view name)
       : TrackClusterMergeSplitterAlgorithm{
             name,
-            {"InputProtoClusterCollection", "InputTrackProjections"},
+            {"InputTrackClusterMatches", "InputClusterCollection", "InputTrackProjections"},
+#if EDM4EIC_VERSION_MAJOR >= 8
+            {"OutputProtoClusterCollection", "OutputTrackProtoClusterMatches"},
+#else
             {"OutputProtoClusterCollection"},
-            "Merges or splits clusters based on tracks projected to them."} {}
+#endif
+            "Merges or splits clusters based on tracks matched to them."} {
+  }
 
-  // public methods
-  void init();
+  // public method
   void process(const Input&, const Output&) const final;
 
 private:
-  const algorithms::GeoSvc& m_geo = algorithms::GeoSvc::instance();
-
   // private methods
-  void get_projections(const edm4eic::TrackSegmentCollection* projections,
-                       VecProj& relevant_projects) const;
-  void match_clusters_to_tracks(const edm4eic::ProtoClusterCollection* clusters,
-                                const VecProj& projections, MapToVecProj& matches) const;
-  void merge_and_split_clusters(const VecClust& to_merge, const VecProj& to_split,
-                                edm4eic::ProtoClusterCollection* out_protoclusters) const;
-  float get_cluster_energy(const edm4eic::ProtoCluster& clust) const;
-  edm4hep::Vector3f get_cluster_position(const edm4eic::ProtoCluster& clust) const;
-
-  // calorimeter id
-  unsigned int m_idCalo{0};
+  void merge_and_split_clusters(const VecClust& to_merge, const VecSeg& to_split,
+                                VecProto& new_protos) const;
+  void add_cluster_to_proto(const edm4eic::Cluster& clust, edm4eic::MutableProtoCluster& proto,
+                            std::optional<MapToWeight> split_weights = std::nullopt) const;
 
 }; // end TrackClusterMergeSplitter
 

src/algorithms/calorimetry/TrackClusterMergeSplitterConfig.h

@@ -1,23 +1,28 @@
 // SPDX-License-Identifier: LGPL-3.0-or-later
-// Copyright (C) 2024 Derek Anderson
+// Copyright (C) 2025 Derek Anderson, Dmitry Kalinkin
 
 #pragma once
 
-#include <string>
-
 namespace eicrecon {
 
 struct TrackClusterMergeSplitterConfig {
 
-  std::string idCalo = "HcalBarrel_ID"; // id of calorimeter to match projections to
+  ///! any clusters below this will be merged
+  double minSigCut = -1.0;
+
+  ///! average of E/p distribution
+  double avgEP = 1.0;
+
+  ///! RMS of E/p distribution
+  double sigEP = 1.0;
+
+  ///! window to merge clusters
+  double drAdd = 0.4;
 
-  double minSigCut = -1.; // min significance
-  double avgEP     = 1.0; // mean E/p
-  double sigEP     = 1.0; // rms of E/p
-  double drAdd     = 0.4; // window to add clusters
-  double sampFrac  = 1.0; // allows for sampling fraction correction
+  ///! index of surface to use for track projections
+  uint64_t surfaceToUse = 1;
 
-  // scale for hit-track distance
+  ///! scale for hit-track distance
   double transverseEnergyProfileScale = 1.0;
 
 }; // end TrackClusterMergeSplitterConfig

src/algorithms/particle/CMakeLists.txt

@@ -0,0 +1,18 @@
+cmake_minimum_required(VERSION 3.16)
+
+set(PLUGIN_NAME "algorithms_particle")
+
+# Function creates ${PLUGIN_NAME}_plugin and ${PLUGIN_NAME}_library targets
+# Setting default includes, libraries and installation paths
+plugin_add(${PLUGIN_NAME} WITH_SHARED_LIBRARY WITHOUT_PLUGIN)
+
+# The macro grabs sources as *.cc *.cpp *.c and headers as *.h *.hh *.hpp Then
+# correctly sets sources for ${_name}_plugin and ${_name}_library targets Adds
+# headers to the correct installation directory
+plugin_glob_all(${PLUGIN_NAME})
+
+# Find dependencies
+plugin_add_algorithms(${PLUGIN_NAME})
+plugin_add_dd4hep(${PLUGIN_NAME})
+plugin_add_event_model(${PLUGIN_NAME})
+plugin_add_eigen3(${PLUGIN_NAME})

src/algorithms/particle/TrackProtoClusterMatchPromoter.cc

@@ -0,0 +1,67 @@
+// SPDX-License-Identifier: LGPL-3.0-or-later
+// Copyright (C) 2025 Derek Anderson
+
+#include <optional>
+
+#include "TrackProtoClusterMatchPromoter.h"
+
+namespace eicrecon {
+
+// ----------------------------------------------------------------------------
+//! Process inputs
+// ----------------------------------------------------------------------------
+/*! For each track-protocluster match, algorithm
+ *  will create a corresponding track-cluster match.
+ *  Input protocluster and cluster collections are
+ *  assumed to be 1-to-1, i.e. assuming that the Nth 
+ *  cluster was reconstructed from the Nth
+ *  protocluster.
+ */
+void TrackProtoClusterMatchPromoter::process(const TrackProtoClusterMatchPromoter::Input& input,
+                                             const TrackProtoClusterMatchPromoter::Output& output) const {
+
+  // grab inputs/outputs
+  const auto [in_matches, in_protos, in_clusts] = input;
+  auto [out_matches] = output;
+
+  // exit if no matches in input collection
+  if (in_matches->size() == 0) {
+    debug("No track-protocluster matches in collection.");
+    return;
+  }
+
+  // exit if protocluster/cluster collection
+  // sizes are different
+  //   --> 1-to-1 ordering can't be assumed!
+  if (in_protos->size() != in_clusts->size()) {
+    error("Number of input protoclusters ({}) not the same as number of input clusters ({})",
+          in_protos->size(), in_clusts->size());
+    return;
+  }
+
+  // loop through protoclusters
+  //   --> if a protocluster has a match, create a
+  //       match for the corresponding cluster
+  for (std::size_t icl = 0; const auto& proto : *in_protos) {
+
+    // NB track-protocluster matches are
+    //   - FROM track
+    //   - TO   protocluster
+    std::optional<edm4eic::TrackProtoClusterMatch> match;
+    for (const auto& pr_match : *in_matches) {
+      if (pr_match.getTo() == proto) {
+        match = pr_match;
+        break;
+      }
+    }
+
+    if (match.has_value()) {
+      edm4eic::MutableTrackClusterMatch cl_match = out_matches->create();
+      cl_match.setCluster((*in_clusts)[icl]);
+      cl_match.setTrack(match.value().getFrom());
+    }
+    ++icl;
+  }
+}  // end 'process(Input&, Output&)'
+
+} // namespace eicrecon

src/algorithms/particle/TrackProtoClusterMatchPromoter.h

@@ -0,0 +1,57 @@
+// SPDX-License-Identifier: LGPL-3.0-or-later
+// Copyright (C) 2025 Derek Anderson
+
+#pragma once
+
+#include <algorithms/algorithm.h>
+#include <edm4eic/ClusterCollection.h>
+#include <edm4eic/ProtoClusterCollection.h>
+#include <edm4eic/TrackClusterMatchCollection.h>
+#include <edm4eic/TrackProtoClusterMatchCollection.h>
+#include <spdlog/spdlog.h>
+#include <string>
+#include <string_view>
+
+#include "services/log/Log_service.h"
+#include "algorithms/interfaces/WithPodConfig.h"
+
+namespace eicrecon {
+
+// ----------------------------------------------------------------------------
+//! Algorithm input/output
+// ----------------------------------------------------------------------------
+using TrackProtoClusterMatchPromoterAlgorithm =
+    algorithms::Algorithm<typename algorithms::Input<
+                              edm4eic::TrackProtoClusterMatchCollection,
+                              edm4eic::ProtoClusterCollection,
+                              edm4eic::ClusterCollection>,
+                          typename algorithms::Output<
+                              edm4eic::TrackClusterMatchCollection>>;
+
+// ============================================================================
+//! Track-Protocluster Match Promoter
+// ============================================================================
+/*! An algorithm to promote track-protocluster matches to
+ *  track-cluster matches. Assumes input cluster,
+ *  protocluster collections are 1-to-1.
+ */
+class TrackProtoClusterMatchPromoter : public TrackProtoClusterMatchPromoterAlgorithm
+                                     , public WithPodConfig<NoConfig> {
+
+public:
+  ///! Algorithm constructor
+  TrackProtoClusterMatchPromoter(std::string_view name)
+    : TrackProtoClusterMatchPromoterAlgorithm{
+          name,
+          {"inputTrackProtoclusterMatches",
+           "inputProtoclusters",
+           "inputClusters"},
+          {"outputTrackClusterMatches"},
+          "Copies track-protocluster matches onto track-cluster matches"} {}
+
+  // public method
+  void process(const Input&, const Output&) const final;
+
+}; // end TrackProtoClusterMatchPromoter
+
+} // namespace eicrecon