Add new TrackParticleAnalysisAlg to monitor track particles.

The particles are assumed to be created from Acts tracks and
additional monitoring is performed using the original Acts
tracks.

Tracking/Acts/ActsConfig/python/ActsAnalysisConfig.py

@@ -38,7 +38,67 @@ def ActsTrackAnalysisAlgCfg(flags,
 
     acc.merge(helper.result())
     return acc
-    
+
+def ActsTrackParticleAnalysisAlgCfg(flags,
+                                    name: str = 'ActsTrackParticleAnalysisAlg',
+                                    **kwargs) -> ComponentAccumulator:
+    acc = ComponentAccumulator()
+    kwargs.setdefault('MonitorTrackStateCounts',True)
+    kwargs.setdefault('TrackParticleLocation', 'ActsTrackParticles')
+
+    kwargs.setdefault('MonGroupName', kwargs['TrackParticleLocation'])
+
+    from AthenaMonitoring import AthMonitorCfgHelper
+    helper = AthMonitorCfgHelper(flags, 'TrackParticleAnalysisAlgCfg')
+    monitoringAlgorithm = helper.addAlgorithm(CompFactory.ActsTrk.TrackParticleAnalysisAlg, name, **kwargs)
+    monitoringGroup = helper.addGroup(monitoringAlgorithm, kwargs['MonGroupName'], '/'+name+'/')
+
+    monitoringGroup.defineHistogram('pt', title='TrackParticle pt;pt (MeV);Entries', type='TH1F', path=kwargs['MonGroupName'],
+                                    xbins=100, xmin=0, xmax=10e3)
+    monitoringGroup.defineHistogram('eta', title='TrackParticle eta;eta;Entries', type='TH1F', path=kwargs['MonGroupName'],
+                                    xbins=50, xmin=-4, xmax=4)
+
+    # hit counts
+    monitoringGroup.defineHistogram('pixelHits', title='Number of pixel hits;N;Entries', type='TH1I', path=kwargs['MonGroupName'],
+                                    xbins=50, xmin=-1, xmax=49)
+    monitoringGroup.defineHistogram('innermostHits', title='Number of innermost pixel hits;N;Entries', type='TH1I', path=kwargs['MonGroupName'],
+                                    xbins=6, xmin=-1, xmax=5)
+    monitoringGroup.defineHistogram('nextToInnermostHits', title='Number of next-to-innermost pixel hits;N;Entries', type='TH1I',
+                                    path=kwargs['MonGroupName'], xbins=6, xmin=-1, xmax=5)
+    monitoringGroup.defineHistogram('expectInnermostHit', title='Innermost pixel hit expected;N;Entries', type='TH1I', path=kwargs['MonGroupName'],
+                                    xbins=3, xmin=-1, xmax=2)
+    monitoringGroup.defineHistogram('expectNextToInnermostHit', title='Next-to-innermost pixel hit expected;N;Entries', type='TH1I',
+                                    path=kwargs['MonGroupName'],xbins=3, xmin=-1, xmax=2)
+
+    if kwargs['MonitorTrackStateCounts'] :
+        # have to add artifical dependency, because the ActsTracks are created by
+        # a special reader algorithm from the various component branches, but the
+        # element links pointing to the ActsTracks would not add this dependency
+        # by themselves.
+        if 'ExtraInputs' not in kwargs :
+            tracks_name = kwargs['TrackParticleLocation']
+            if tracks_name[-len('TrackParticles'):] == 'TrackParticles' :
+                kwargs.setdefault('ExtraInputs',{('ActsTrk::TrackContainer',tracks_name[:-len('TrackParticles')]+'Tracks')})
+
+        monitoringGroup.defineHistogram('States', title='Number of states on track;N;Entries', type='TH1I', path=kwargs['MonGroupName'],
+                                        xbins=50, xmin=0, xmax=50)
+        monitoringGroup.defineHistogram('Measurements', title='Number of measurements on track;N;Entries', type='TH1I', path=kwargs['MonGroupName'],
+                                        xbins=50, xmin=0, xmax=50)
+        monitoringGroup.defineHistogram('Parameters', title='Number of parameters on track;N;Entries', type='TH1I', path=kwargs['MonGroupName'],
+                                        xbins=50, xmin=0, xmax=50)
+        monitoringGroup.defineHistogram('Outliers', title='Number of outliers on track;N;Entries', type='TH1I', path=kwargs['MonGroupName'],
+                                        xbins=50, xmin=0, xmax=50)
+        monitoringGroup.defineHistogram('Holes', title='Number of holes on track;N;Entries', type='TH1I', path=kwargs['MonGroupName'],
+                                        xbins=50, xmin=0, xmax=50)
+        monitoringGroup.defineHistogram('SharedHits', title='Number of shared hits on track;N;Entries', type='TH1I', path=kwargs['MonGroupName'],
+                                        xbins=50, xmin=0, xmax=50)
+        monitoringGroup.defineHistogram('MaterialStates', title='Number of material states on track;N;Entries', type='TH1I', path=kwargs['MonGroupName'],
+                                        xbins=50, xmin=0, xmax=50)
+    acc.merge(helper.result())
+
+    acc.addEventAlgo( CompFactory.ActsTrk.TrackParticleAnalysisAlg(name, **kwargs) )
+    return acc
+
 def ActsHgtdClusterAnalysisAlgCfg(flags,
                                   name: str = "ActsHgtdClusterAnalysisAlg",
                                   **kwargs) -> ComponentAccumulator:

Tracking/Acts/ActsConfig/test/ActsPersistifyEDM.sh

@@ -29,7 +29,9 @@ ActsReadEDM.py \
    readClusters=True \
    readSpacePoints=True \
    readTracks=True \
-   tracks="ActsTracks"
+   tracks="ActsTracks" \
+   readTrackParticles=True \
+   trackParticles="ActsCombinedTracksParticlesAlt"
 
 rc=$?
 if [ $rc != 0 ]; then

Tracking/Acts/ActsMonitoring/src/TrackParticleAnalysisAlg.cxx

+/*
+  Copyright (C) 2002-2024 CERN for the benefit of the ATLAS collaboration
+*/
+
+#include "TrackParticleAnalysisAlg.h"
+#include "ActsEvent/TrackContainer.h"
+#include "SGTools/CurrentEventStore.h"
+#include <iostream>
+#include <sstream>
+#include <iomanip>
+
+namespace {
+   int getSummaryValue(const xAOD::TrackParticle &track_particle, xAOD::SummaryType summary_type) {
+      uint8_t tmp;
+      return track_particle.summaryValue(tmp, summary_type) ? static_cast<int>(tmp) : -1;
+   }
+}
+
+
+namespace ActsTrk {
+
+  TrackParticleAnalysisAlg::TrackParticleAnalysisAlg(const std::string& name, ISvcLocator* pSvcLocator)
+    : AthMonitorAlgorithm(name, pSvcLocator)
+  {}
+
+  StatusCode TrackParticleAnalysisAlg::initialize() {
+    ATH_MSG_DEBUG("Initializing " << name() << " ...");
+
+    ATH_MSG_VERBOSE("Properties:");
+    ATH_MSG_VERBOSE(m_tracksKey);
+
+    ATH_CHECK(m_tracksKey.initialize());
+
+    ATH_MSG_VERBOSE("Monitoring settings ...");
+    ATH_MSG_VERBOSE(m_monGroupName);
+
+    return AthMonitorAlgorithm::initialize();
+  }
+
+  StatusCode TrackParticleAnalysisAlg::fillHistograms(const EventContext& ctx) const {
+    ATH_MSG_DEBUG( "Filling Histograms for " << name() << " ... " );
+
+    // Retrieve the tracks
+    SG::ReadHandle<xAOD::TrackParticleContainer> trackParticleskHandle = SG::makeHandle(m_tracksKey, ctx);
+    ATH_CHECK(trackParticleskHandle.isValid());
+    const xAOD::TrackParticleContainer *track_particles = trackParticleskHandle.cptr();
+
+    if (m_monitorTrackStateCounts) {
+       for (const xAOD::TrackParticle *track_particle : *track_particles) {
+          monitorTrackStateCounts(*track_particle);
+       }
+    }
+    auto monitor_pt = Monitored::Collection("pt", *track_particles,
+                                            [](const xAOD::TrackParticle *tp){return static_cast<double>(tp->pt()); } );
+    auto monitor_eta = Monitored::Collection("eta", *track_particles,
+                                             [](const xAOD::TrackParticle *tp){return static_cast<double>(tp->eta()); } );
+
+    auto monitor_pixelHits = Monitored::Collection("pixelHits", *track_particles,
+                                                   [](const xAOD::TrackParticle *tp){
+                                                      return getSummaryValue(*tp,xAOD::numberOfPixelHits); });
+
+    auto monitor_innermostHits = Monitored::Collection("innermostHits", *track_particles,
+                                                       [](const xAOD::TrackParticle *tp){
+                                                          return getSummaryValue(*tp,xAOD::numberOfInnermostPixelLayerHits); });
+
+    auto monitor_nextToInnermostHits = Monitored::Collection("nextToInnermostHits", *track_particles,
+                                                       [](const xAOD::TrackParticle *tp){
+                                                          return getSummaryValue(*tp,xAOD::numberOfNextToInnermostPixelLayerHits);});
+
+    auto monitor_expectInnermost = Monitored::Collection("expectInnermostHit", *track_particles,
+                                             [](const xAOD::TrackParticle *tp){
+                                                return getSummaryValue(*tp,xAOD::expectInnermostPixelLayerHit); });
+
+    auto monitor_expectNextToInnermost = Monitored::Collection("expectNextToInnermostHit", *track_particles,
+                                             [](const xAOD::TrackParticle *tp){
+                                                return getSummaryValue(*tp,xAOD::expectNextToInnermostPixelLayerHit); });
+
+    fill(m_monGroupName.value(),monitor_pt,monitor_eta, monitor_pixelHits,  monitor_innermostHits, monitor_nextToInnermostHits,
+         monitor_expectInnermost, monitor_expectNextToInnermost);
+    return StatusCode::SUCCESS;
+  }
+
+  void TrackParticleAnalysisAlg::monitorTrackStateCounts(const xAOD::TrackParticle &track_particle) const {
+     static const SG::AuxElement::ConstAccessor<ElementLink<ActsTrk::TrackContainer> > actsTrackLink("actsTrack");
+
+     ElementLink<ActsTrk::TrackContainer> link_to_track = actsTrackLink(track_particle);
+     // to ensure that the code does not suggest something stupid (i.e. creating an unnecessary copy)
+     static_assert( std::is_same<ElementLink<ActsTrk::TrackContainer>::ElementConstReference,
+                    std::optional<ActsTrk::TrackContainer::ConstTrackProxy> >::value);
+     std::optional<ActsTrk::TrackContainer::ConstTrackProxy> optional_track = *link_to_track;
+     if (optional_track.has_value()) {
+        ActsTrk::TrackContainer::ConstTrackProxy track = optional_track.value();
+        std::array<uint8_t, Acts::NumTrackStateFlags+1> counts{};
+
+        const ActsTrk::TrackContainer::ConstTrackProxy::IndexType
+           lastMeasurementIndex = track.tipIndex();
+
+        track.container().trackStateContainer().visitBackwards(
+                 lastMeasurementIndex,
+                 [&counts
+                  ](const typename ActsTrk::TrackStateBackend::ConstTrackStateProxy &state) -> void
+                 {
+                    Acts::ConstTrackStateType flag = state.typeFlags();
+                    ++counts[Acts::NumTrackStateFlags];
+                    for (unsigned int flag_i=0; flag_i<Acts::NumTrackStateFlags; ++flag_i) {
+                       if (flag.test(flag_i)) {
+                          if (flag_i == Acts::TrackStateFlag::HoleFlag) {
+                             if (!state.hasReferenceSurface() || !state.referenceSurface().associatedDetectorElement()) continue;
+                          }
+                          ++counts[flag_i];
+                       }
+                    }
+                 });
+        auto monitor_states = Monitored::Scalar<int>("States",counts[Acts::TrackStateFlag::NumTrackStateFlags]);
+        auto monitor_measurement = Monitored::Scalar<int>("Measurements",counts[Acts::TrackStateFlag::MeasurementFlag]);
+        auto monitor_parameter = Monitored::Scalar<int>("Parameters",counts[Acts::TrackStateFlag::ParameterFlag]);
+        auto monitor_outlier = Monitored::Scalar<int>("Outliers",counts[Acts::TrackStateFlag::OutlierFlag]);
+        auto monitor_hole = Monitored::Scalar<int>("Holes",counts[Acts::TrackStateFlag::HoleFlag]);
+        auto monitor_material = Monitored::Scalar<int>("MaterialStates",counts[Acts::TrackStateFlag::MaterialFlag]);
+        auto monitor_sharedHit = Monitored::Scalar<int>("SharedHits",counts[Acts::TrackStateFlag::SharedHitFlag]);
+        fill(m_monGroupName.value(), monitor_states, monitor_measurement, monitor_parameter, monitor_outlier, monitor_hole,
+             monitor_material, monitor_sharedHit);
+     }
+     else {
+        ATH_MSG_WARNING("Invalid track link for particle  " << track_particle.index());
+     }
+  }
+
+
+}

Tracking/Acts/ActsMonitoring/src/TrackParticleAnalysisAlg.h

+/*
+  Copyright (C) 2002-2024 CERN for the benefit of the ATLAS collaboration
+*/
+
+#ifndef ACTSTRKANALYSIS_TRACKPARTICLEANALYSISALG_H
+#define ACTSTRKANALYSIS_TRACKPARTICLEANALYSISALG_H
+
+#include "AthenaMonitoring/AthMonitorAlgorithm.h"
+#include "xAODTracking/TrackParticleContainer.h"
+#include "StoreGate/ReadHandleKey.h"
+
+namespace ActsTrk {
+
+  class TrackParticleAnalysisAlg final :
+    public AthMonitorAlgorithm {
+  public:
+    TrackParticleAnalysisAlg(const std::string& name, ISvcLocator* pSvcLocator);
+    virtual ~TrackParticleAnalysisAlg() override = default;
+
+    virtual StatusCode initialize() override;
+    virtual StatusCode fillHistograms(const EventContext& ctx) const override;
+
+  private:
+    void monitorTrackStateCounts(const xAOD::TrackParticle &track_particle) const;
+
+    SG::ReadHandleKey<xAOD::TrackParticleContainer> m_tracksKey {this, "TrackParticleLocation", "",
+        "Input track particle collection"};
+    Gaudi::Property< std::string > m_monGroupName
+      {this, "MonGroupName", "TrackParticleAnalysisAlg"};
+
+    Gaudi::Property< bool > m_monitorTrackStateCounts
+      {this, "MonitorTrackStateCounts", true};
+  };
+
+}
+
+#endif

Tracking/Acts/ActsMonitoring/src/components/ActsMonitoring_entries.cxx

 /*
-  Copyright (C) 2002-2023 CERN for the benefit of the ATLAS collaboration
+  Copyright (C) 2002-2024 CERN for the benefit of the ATLAS collaboration
 */
 
 // Algs
@@ -11,6 +11,7 @@
 #include "src/EstimatedTrackParamsAnalysisAlg.h"
 #include "src/SeedingAlgorithmAnalysisAlg.h"
 #include "src/TrackAnalysisAlg.h"
+#include "src/TrackParticleAnalysisAlg.h"
 // Tools
 #include "src/PhysValTool.h"
 
@@ -23,5 +24,6 @@ DECLARE_COMPONENT( ActsTrk::SeedAnalysisAlg )
 DECLARE_COMPONENT( ActsTrk::SeedingAlgorithmAnalysisAlg )
 DECLARE_COMPONENT( ActsTrk::EstimatedTrackParamsAnalysisAlg )
 DECLARE_COMPONENT( ActsTrk::TrackAnalysisAlg )
+DECLARE_COMPONENT( ActsTrk::TrackParticleAnalysisAlg )
 // Tools
 DECLARE_COMPONENT( ActsTrk::PhysValTool )

Tracking/Acts/ActsMonitoring/test/ActsReadEDM.py

@@ -3,6 +3,7 @@
 # Copyright (C) 2002-2023 CERN for the benefit of the ATLAS collaboration
 
 if __name__ == "__main__":
+    import re
     from AthenaConfiguration.AllConfigFlags import initConfigFlags
     flags = initConfigFlags()
 
@@ -20,6 +21,8 @@ if __name__ == "__main__":
     flags.addFlag("readSpacePoints", False)
     flags.addFlag("readTracks", False)
     flags.addFlag("tracks", "")
+    flags.addFlag("readTrackParticles", False)
+    flags.addFlag("trackParticles", "ActsCombinedTracksParticlesAlt")
     flags.fillFromArgs()
     
     flags.lock()
@@ -58,6 +61,16 @@ if __name__ == "__main__":
                                               name=f"{track}AnalysisAlg",
                                               OutputLevel=2,
                                               TracksLocation=track))
+    if flags.readTrackParticles:
+        from ActsConfig.ActsAnalysisConfig import ActsTrackParticleAnalysisAlgCfg
+        for tp in flags.trackParticles.split(','):
+            src_track_name = re.search(r'^(.*)ParticlesAlt.*',tp).group(1)
+            acc.merge(ActsTrackParticleAnalysisAlgCfg(flags,
+                                                      name=f"{tp}AnalysisAlg",
+                                                      OutputLevel=2,
+                                                      TrackParticleLocation=tp,
+                                                      ExtraInputs={('ActsTrk::TrackContainer',src_track_name)}, # ensure scheduled after reader
+                                                      MonGroupName=f"{tp}Analysis"))
             
     acc.printConfig()
     status = acc.run()