refactor+fix: use collections and multi-factory in `RichTrack_factory…

…`; fix potential leaks in `TrackPropagation` (#656)

### Briefly, what does this PR introduce?
- return `Collection` of `TrackSegments` in
`TrackPropagation::propagateToSurfaceList`, rather than a (potentially
leaky) new `TrackSegment` pointer
- the loop over input trajectories is moved from `RichTrack_factory` to
`TrackPropagation::propagateToSurfaceList`, which makes the PR diff
appear to be more complex than it actually is
- fix additional possible memory leaks by replacing any `return new
raw_pointer` calls with `return std::make_unique<...>(...)`, and update
any code which depends on this change
- improve the configuration with `RichTrackConfig`, which configures the
_factory_, not the algorithm, since the factory is used to connect to
the `richgeo` service (and not the more general `TrackPropagation`
algorithm)


### What kind of change does this PR introduce?
- [x] Bug fix (memory leaks)
- [x] New feature (issue #__)
- [ ] Documentation update
- [ ] Other: __

### Please check if this PR fulfills the following:
- [x] Tests for the changes have been added: tested in `irt-algo` branch
- [ ] Documentation has been added / updated
- [ ] Changes have been communicated to collaborators

### Does this PR introduce breaking changes? What changes might users
need to make to their code?
no, unless the `TrackPropagation` algorithm is called by anything
external of EICrecon

### Does this PR change default behavior?
no

src/algorithms/tracking/TrackPropagation.cc

@@ -1,5 +1,5 @@
 // SPDX-License-Identifier: LGPL-3.0-or-later
-// Copyright (C) 2022 Wenqing Fan, Barak Schmookler, Whitney Armstrong, Sylvester Joosten, Dmitry Romanov
+// Copyright (C) 2022, 2023 Wenqing Fan, Barak Schmookler, Whitney Armstrong, Sylvester Joosten, Dmitry Romanov, Christopher Dilks
 
 #include <cmath>
 #include <algorithm>
@@ -51,11 +51,11 @@ namespace eicrecon {
 
 
 
-    std::vector<edm4eic::TrackPoint *>
+    std::vector<std::unique_ptr<edm4eic::TrackPoint>>
     TrackPropagation::propagateMany(std::vector<const eicrecon::TrackingResultTrajectory *> trajectories,
                                     const std::shared_ptr<const Acts::Surface> &targetSurf) {
         // output collection
-        std::vector<edm4eic::TrackPoint *> track_points;
+        std::vector<std::unique_ptr<edm4eic::TrackPoint>> track_points;
         m_log->trace("Track propagation evnet process. Num of input trajectories: {}", std::size(trajectories));
 
         // Loop over the trajectories
@@ -67,66 +67,82 @@ namespace eicrecon {
             if(!result) continue;
 
             // Add to output collection
-            track_points.push_back(result);
+            track_points.push_back(std::move(result));
         }
 
         return track_points;
     }
 
 
 
-    edm4eic::TrackSegment* TrackPropagation::propagateToSurfaceList(const eicrecon::TrackingResultTrajectory *traj,
-                                                                    std::vector<std::shared_ptr<Acts::Surface>> targetSurfaces) {
-      // start a mutable TrackSegment
-      edm4eic::MutableTrackSegment track_segment;
-      decltype(edm4eic::TrackSegmentData::length)      length       = 0;
-      decltype(edm4eic::TrackSegmentData::lengthError) length_error = 0;
-
-      // loop over projection-target surfaces
-      for(const auto& targetSurf : targetSurfaces) {
-
-        // project the trajectory `traj` to this surface
-        edm4eic::TrackPoint *point;
-        try {
-          point = propagate(traj, targetSurf);
-        } catch(std::exception &e) {
-          m_log->warn("<> Exception in TrackPropagation::propagateToSurfaceList: {}; skip this TrackPoint and surface", e.what());
-        }
-        if(!point) {
-          m_log->trace("<> Failed to propagate trajectory to this plane");
-          continue;
-        }
-
-        // logging
-        m_log->trace("<> trajectory: x=( {:>10.2f} {:>10.2f} {:>10.2f} )",
-            point->position.x, point->position.y, point->position.z);
-        m_log->trace("               p=( {:>10.2f} {:>10.2f} {:>10.2f} )",
-            point->momentum.x, point->momentum.y, point->momentum.z);
-
-        // update the `TrackSegment` length
-        // FIXME: `length` and `length_error` are currently not used by any callers, and may not be correctly calculated here
-        if(track_segment.points_size()>0) {
-          auto pos0 = point->position;
-          auto pos1 = std::prev(track_segment.points_end())->position;
-          auto dist = edm4eic::magnitude(pos0-pos1);
-          length += dist;
-          m_log->trace("               dist to previous point: {}", dist);
-        }
-
-        // add the `TrackPoint` to the `TrackSegment`
-        track_segment.addToPoints(*point);
-
-      } // end `targetSurfaces` loop
-
-      // output
-      track_segment.setLength(length);
-      track_segment.setLengthError(length_error);
-      return new edm4eic::TrackSegment(track_segment);
+    std::unique_ptr<edm4eic::TrackSegmentCollection> TrackPropagation::propagateToSurfaceList(
+        std::vector<const eicrecon::TrackingResultTrajectory*> trajectories,
+        std::vector<std::shared_ptr<Acts::Surface>> targetSurfaces
+        )
+    {
+      // logging
+      m_log->trace("Propagate trajectories: --------------------");
+      m_log->trace("number of trajectories: {}",trajectories.size());
+
+      // start output collection
+      auto track_segments = std::make_unique<edm4eic::TrackSegmentCollection>();
+
+      // loop over input trajectories
+      for(const auto& traj : trajectories) {
+
+        // start a mutable TrackSegment
+        auto track_segment = track_segments->create();
+        decltype(edm4eic::TrackSegmentData::length)      length       = 0;
+        decltype(edm4eic::TrackSegmentData::lengthError) length_error = 0;
+
+        // loop over projection-target surfaces
+        for(const auto& targetSurf : targetSurfaces) {
+
+          // project the trajectory `traj` to this surface
+          std::unique_ptr<edm4eic::TrackPoint> point;
+          try {
+            point = propagate(traj, targetSurf);
+          } catch(std::exception &e) {
+            m_log->warn("<> Exception in TrackPropagation::propagateToSurfaceList: {}; skip this TrackPoint and surface", e.what());
+          }
+          if(!point) {
+            m_log->trace("<> Failed to propagate trajectory to this plane");
+            continue;
+          }
+
+          // logging
+          m_log->trace("<> trajectory: x=( {:>10.2f} {:>10.2f} {:>10.2f} )",
+              point->position.x, point->position.y, point->position.z);
+          m_log->trace("               p=( {:>10.2f} {:>10.2f} {:>10.2f} )",
+              point->momentum.x, point->momentum.y, point->momentum.z);
+
+          // update the `TrackSegment` length
+          // FIXME: `length` and `length_error` are currently not used by any callers, and may not be correctly calculated here
+          if(track_segment.points_size()>0) {
+            auto pos0 = point->position;
+            auto pos1 = std::prev(track_segment.points_end())->position;
+            auto dist = edm4eic::magnitude(pos0-pos1);
+            length += dist;
+            m_log->trace("               dist to previous point: {}", dist);
+          }
+
+          // add the `TrackPoint` to the `TrackSegment`
+          track_segment.addToPoints(*point);
+
+        } // end `targetSurfaces` loop
+
+        // set final length and length error
+        track_segment.setLength(length);
+        track_segment.setLengthError(length_error);
+
+      } // end loop over input trajectories
+
+      return track_segments;
     }
 
 
 
-    edm4eic::TrackPoint *TrackPropagation::propagate(const eicrecon::TrackingResultTrajectory *traj,
+    std::unique_ptr<edm4eic::TrackPoint> TrackPropagation::propagate(const eicrecon::TrackingResultTrajectory *traj,
                                                      const std::shared_ptr<const Acts::Surface> &targetSurf) {
         // Get the entry index for the single trajectory
         // The trajectory entry indices and the multiTrajectory
@@ -263,7 +279,7 @@ namespace eicrecon {
           float pathlength{}; ///< Pathlength from the origin to this point
           float pathlengthError{}; ///< Error on the pathlenght
          */
-        return new edm4eic::TrackPoint({
+        return std::make_unique<edm4eic::TrackPoint>(edm4eic::TrackPoint{
                                                position,
                                                positionError,
                                                momentum,

src/algorithms/tracking/TrackPropagation.h

@@ -39,17 +39,19 @@ namespace eicrecon {
         void init(std::shared_ptr<const ActsGeometryProvider> geo_svc, std::shared_ptr<spdlog::logger> logger);
 
         /** Propagates a single trajectory to a given surface */
-        edm4eic::TrackPoint * propagate(const eicrecon::TrackingResultTrajectory *, const std::shared_ptr<const Acts::Surface>& targetSurf);
+        std::unique_ptr<edm4eic::TrackPoint> propagate(const eicrecon::TrackingResultTrajectory *, const std::shared_ptr<const Acts::Surface>& targetSurf);
 
         /** Propagates a collection of trajectories to a given surface
          * @remark: being a simple wrapper of propagate(...) this method is more sutable for factories */
-        std::vector<edm4eic::TrackPoint *> propagateMany(std::vector<const eicrecon::TrackingResultTrajectory *> trajectories,
+        std::vector<std::unique_ptr<edm4eic::TrackPoint>> propagateMany(std::vector<const eicrecon::TrackingResultTrajectory *> trajectories,
                                                          const std::shared_ptr<const Acts::Surface> &targetSurf);
 
-        /** Propagates a trajectory to a list of surfaces, and returns the full `TrackSegment`
+        /** Propagates a collection of trajectories to a list of surfaces, and returns the full `TrackSegment`
          * @remark: being a simple wrapper of propagate(...) this method is more sutable for factories */
-        edm4eic::TrackSegment* propagateToSurfaceList(const eicrecon::TrackingResultTrajectory *traj,
-                                                      std::vector<std::shared_ptr<Acts::Surface>> targetSurfaces);
+        std::unique_ptr<edm4eic::TrackSegmentCollection> propagateToSurfaceList(
+            std::vector<const eicrecon::TrackingResultTrajectory*> trajectories,
+            std::vector<std::shared_ptr<Acts::Surface>> targetSurfaces
+            );
 
     private:
 

src/detectors/DRICH/DRICH.cc

@@ -14,6 +14,7 @@
 
 // algorithm configurations
 #include <algorithms/digi/PhotoMultiplierHitDigiConfig.h>
+#include <global/pid/RichTrackConfig.h>
 
 
 extern "C" {
@@ -54,6 +55,11 @@ extern "C" {
     digi_cfg.quantumEfficiency.push_back({850, 0.06});
     digi_cfg.quantumEfficiency.push_back({900, 0.04});
 
+    // track propagation to each radiator
+    RichTrackConfig track_cfg;
+    track_cfg.numPlanes.insert({ "Aerogel", 5  });
+    track_cfg.numPlanes.insert({ "Gas",     10 });
+
 
     // wiring between factories and data ///////////////////////////////////////
     // clang-format off
@@ -68,13 +74,12 @@ extern "C" {
           ));
 
     // charged particle tracks
-    app->Add(new JChainFactoryGeneratorT<RichTrack_factory>(
-          {"CentralCKFTrajectories"},
-          "DRICHAerogelTracks"
-          ));
-    app->Add(new JChainFactoryGeneratorT<RichTrack_factory>(
+    app->Add(new JChainMultifactoryGeneratorT<RichTrack_factory>(
+          "DRICHTracks",
           {"CentralCKFTrajectories"},
-          "DRICHGasTracks"
+          {"DRICHAerogelTracks", "DRICHGasTracks"},
+          track_cfg,
+          app
           ));
 
     // clang-format on

src/global/pid/RichTrackConfig.h

@@ -0,0 +1,36 @@
+// Copyright 2023, Christopher Dilks
+// Subject to the terms in the LICENSE file found in the top-level directory.
+
+#pragma once
+
+#include <spdlog/spdlog.h>
+
+namespace eicrecon {
+
+  class RichTrackConfig {
+    public:
+
+      /////////////////////////////////////////////////////
+      // CONFIGURATION PARAMETERS
+      //   NOTE: some defaults are hard-coded here; override externally
+
+      std::map <std::string,unsigned> numPlanes; // number of xy-planes for track projections (for each radiator)
+
+      //
+      /////////////////////////////////////////////////////
+
+
+      // print all parameters
+      void Print(
+          std::shared_ptr<spdlog::logger> m_log,
+          spdlog::level::level_enum lvl=spdlog::level::debug
+          )
+      {
+        m_log->log(lvl, "{:=^60}"," RichTrackConfig Settings ");
+        for(const auto& [rad,val] : numPlanes)
+          m_log->log(lvl, "  {:>20} = {:<}", fmt::format("{} numPlanes", rad), val);
+        m_log->log(lvl, "{:=^60}","");
+      }
+
+  };
+}

src/global/pid/RichTrack_factory.cc

@@ -1,53 +1,53 @@
-// Copyright 2022, Christopher Dilks
+// Copyright (C) 2022, 2023, Christopher Dilks
 // Subject to the terms in the LICENSE file found in the top-level directory.
 
 #include "RichTrack_factory.h"
 
-#include <edm4eic/vector_utils.h>
-
 //-----------------------------------------------------------------------------
 void eicrecon::RichTrack_factory::Init() {
   auto app = GetApplication();
 
   // input tags
   auto detector_name = eicrecon::str::ReplaceAll(GetPluginName(), ".so", "");
-  auto param_prefix = detector_name + ":" + GetTag();
-  InitDataTags(param_prefix);
-  m_radiatorID = richgeo::ParseRadiatorName(GetTag());
+  auto param_prefix  = detector_name + GetPrefix();
 
   // services
   m_richGeoSvc = app->GetService<RichGeo_service>();
   m_actsSvc    = app->GetService<ACTSGeo_service>();
   InitLogger(param_prefix, "info");
   m_propagation_algo.init(m_actsSvc->actsGeoProvider(), m_log);
-  m_log->debug("detector_name='{}'  param_prefix='{}'  m_radiatorID={}", detector_name, param_prefix, m_radiatorID);
+  m_log->debug("detector_name='{}'  param_prefix='{}'", detector_name, param_prefix);
 
-  // default configuration parameters
-  /* FIXME: eventually we will move this factory to an independent algorithm
-   * and be able to use a Config class, overridable in {D,PF}RICH.cc; for now,
-   * a quick way to set the defaults without reco_flags.py, to tide us over
-   * until we have a proper config file front end
-   */
-  m_numPlanes = 5;
-  if(param_prefix.rfind("DRICH")!=std::string::npos || param_prefix.rfind("PFRICH")!=std::string::npos) {
-    if(param_prefix.rfind("Aerogel")!=std::string::npos) m_numPlanes = 5;
-    else if(param_prefix.rfind("Gas")!=std::string::npos) m_numPlanes = 10;
+  // get list of radiators
+  std::map<int,std::string> radiator_list;
+  for(auto& output_tag : GetOutputTags()) {
+    auto radiator_id = richgeo::ParseRadiatorName(output_tag, m_log);
+    radiator_list.insert({
+        radiator_id,
+        richgeo::RadiatorName(radiator_id, m_log)
+        });
   }
 
   // configuration parameters
-  auto set_param = [&param_prefix, &app, this] (std::string name, auto &val, std::string description) {
-    app->SetDefaultParameter(param_prefix+":"+name, val, description);
+  auto cfg = GetDefaultConfig();
+  auto set_param = [&param_prefix, &app] (std::string name, auto &val, std::string description) {
+    name = param_prefix + ":" + name;
+    app->SetDefaultParameter(name, val, description);
   };
-  set_param("numPlanes", m_numPlanes, "number of track-projection planes");
-  m_log->debug("numPlanes = {}",m_numPlanes);
+  for(auto& [radiator_id, radiator_name] : radiator_list)
+    set_param(radiator_name+":numPlanes", cfg.numPlanes[radiator_name], "");
+  cfg.Print(m_log, spdlog::level::debug);
 
-  // get RICH geometry for track projection
+  // get RICH geometry for track projection, for each radiator
   m_actsGeo = m_richGeoSvc->GetActsGeo(detector_name);
-  m_trackingPlanes = m_actsGeo->TrackingPlanes(m_radiatorID, m_numPlanes);
+  for(auto& [radiator_id, radiator_name] : radiator_list)
+    m_tracking_planes.push_back(
+        m_actsGeo->TrackingPlanes(radiator_id, cfg.numPlanes.at(radiator_name))
+        );
 }
 
 //-----------------------------------------------------------------------------
-void eicrecon::RichTrack_factory::ChangeRun(const std::shared_ptr<const JEvent> &event) {
+void eicrecon::RichTrack_factory::BeginRun(const std::shared_ptr<const JEvent> &event) {
 }
 
 //-----------------------------------------------------------------------------
@@ -65,15 +65,16 @@ void eicrecon::RichTrack_factory::Process(const std::shared_ptr<const JEvent> &e
     }
   }
 
-  // result will be `TrackSegments`s for each trajectory
-  std::vector<edm4eic::TrackSegment*> result;
-
-  // loop over trajectories
-  m_log->debug("Propagate trajectories: --------------------");
-  m_log->debug("number of trajectories: {}",trajectories.size());
-  for(const auto& traj : trajectories)
-    result.push_back(m_propagation_algo.propagateToSurfaceList( traj, m_trackingPlanes ));
-
-  // output
-  Set(std::move(result));
+  // run algorithm, for each radiator
+  for(unsigned i=0; i<m_tracking_planes.size(); i++) {
+    auto& radiator_tracking_planes = m_tracking_planes.at(i);
+    auto& out_collection_name      = GetOutputTags().at(i);
+    try {
+      auto result = m_propagation_algo.propagateToSurfaceList(trajectories, radiator_tracking_planes);
+      SetCollection<edm4eic::TrackSegment>(out_collection_name, std::move(result));
+    }
+    catch(std::exception &e) {
+      m_log->warn("Exception in underlying algorithm: {}. Event data will be skipped", e.what());
+    }
+  }
 }