Add MinP and MinPT properties to PrMatchNN

Pr/PrAlgorithms/src/PrMatchNN.cpp

@@ -128,6 +128,8 @@ namespace LHCb::Pr::MatchNN {
     Gaudi::Property<float> m_maxDistY{this, "MaxDistY", 250 * Gaudi::Units::mm};
     Gaudi::Property<float> m_maxDSlope{this, "MaxDSlope", 1.5};
     Gaudi::Property<float> m_maxDSlopeY{this, "MaxDSlopeY", 0.15};
+    Gaudi::Property<float> m_minP{this, "MinP", 0.0 * Gaudi::Units::MeV};
+    Gaudi::Property<float> m_minPt{this, "MinPt", 0.0 * Gaudi::Units::MeV};
 
     // -- Counters
     mutable Gaudi::Accumulators::MsgCounter<MSG::WARNING>     m_momentum_failed{this, "momentum determination failed!"};
@@ -360,30 +362,17 @@ namespace LHCb::Pr::MatchNN {
     auto const seediter = seeds.simd();
     auto const veloiter = velos.simd();
 
+    const bool          apply_min_p_cut = m_minP.value() > 0;
+    const simd::float_v maxQOvP =
+        ( apply_min_p_cut ? 1.0 / m_minP.value() : std::numeric_limits<simd::float_v>::max() );
+    const bool apply_min_pt_cut = m_minPt.value() > 0;
+
     for ( auto const& match : matches.simd() ) {
-      auto       loopMask = match.loop_mask();
-      auto const oTrack   = result.compress_back<SIMDWrapper::InstructionSet::Best>( loopMask );
 
-      oTrack.field<TracksTag::trackVP>().set( match.get<Tag::veloIndex>() );
-      oTrack.field<TracksTag::trackUT>().set( -1 );
-      oTrack.field<TracksTag::trackSeed>().set( match.get<Tag::seedIndex>() );
+      auto loopMask = match.loop_mask();
 
       auto const seed_track = seediter.gather( match.get<Tag::seedIndex>(), loopMask );
       auto const velo_track = veloiter.gather( match.get<Tag::veloIndex>(), loopMask );
-      auto const n_fthits   = seed_track.nHits();
-      auto const n_vphits   = velo_track.nHits();
-      oTrack.field<TracksTag::VPHits>().resize( n_vphits );
-      oTrack.field<TracksTag::UTHits>().resize( 0 );
-      oTrack.field<TracksTag::FTHits>().resize( n_fthits );
-
-      for ( auto idx{0}; idx < n_vphits.hmax( loopMask ); ++idx ) {
-        oTrack.field<TracksTag::VPHits>()[idx].template field<TracksTag::Index>().set( velo_track.vp_index( idx ) );
-        oTrack.field<TracksTag::VPHits>()[idx].template field<TracksTag::LHCbID>().set( velo_track.vp_lhcbID( idx ) );
-      }
-      for ( auto idx{0}; idx < n_fthits.hmax( loopMask ); ++idx ) {
-        oTrack.field<TracksTag::FTHits>()[idx].template field<TracksTag::Index>().set( seed_track.ft_index( idx ) );
-        oTrack.field<TracksTag::FTHits>()[idx].template field<TracksTag::LHCbID>().set( seed_track.ft_lhcbID( idx ) );
-      }
 
       //== get Velo and T states at the usual pattern reco positions
       auto state_endvelo = velo_track.get<VeloTag::States>( 1 );
@@ -401,6 +390,42 @@ namespace LHCb::Pr::MatchNN {
           } )
           .ignore();
 
+      // If enabled filter on min track P, Pt
+      if ( apply_min_p_cut ) {
+        auto pMask = abs( qOverP ) < maxQOvP;
+        loopMask   = loopMask & pMask;
+      }
+      if ( apply_min_pt_cut ) {
+        const auto momentum = abs( 1.0 / qOverP );
+        const auto tx2      = state_beam.tx() * state_beam.tx();
+        const auto ty2      = state_beam.ty() * state_beam.ty();
+        const auto sinTrack = sqrt( 1.0 - ( 1.0 / ( 1.0 + tx2 + ty2 ) ) );
+        const auto pt       = sinTrack * momentum;
+        auto       ptMask   = m_minPt.value() < pt;
+        loopMask            = loopMask & ptMask;
+      }
+
+      auto const oTrack = result.compress_back<SIMDWrapper::InstructionSet::Best>( loopMask );
+
+      oTrack.field<TracksTag::trackVP>().set( match.get<Tag::veloIndex>() );
+      oTrack.field<TracksTag::trackUT>().set( -1 );
+      oTrack.field<TracksTag::trackSeed>().set( match.get<Tag::seedIndex>() );
+
+      auto const n_fthits = seed_track.nHits();
+      auto const n_vphits = velo_track.nHits();
+      oTrack.field<TracksTag::VPHits>().resize( n_vphits );
+      oTrack.field<TracksTag::UTHits>().resize( 0 );
+      oTrack.field<TracksTag::FTHits>().resize( n_fthits );
+
+      for ( auto idx{0}; idx < n_vphits.hmax( loopMask ); ++idx ) {
+        oTrack.field<TracksTag::VPHits>()[idx].template field<TracksTag::Index>().set( velo_track.vp_index( idx ) );
+        oTrack.field<TracksTag::VPHits>()[idx].template field<TracksTag::LHCbID>().set( velo_track.vp_lhcbID( idx ) );
+      }
+      for ( auto idx{0}; idx < n_fthits.hmax( loopMask ); ++idx ) {
+        oTrack.field<TracksTag::FTHits>()[idx].template field<TracksTag::Index>().set( seed_track.ft_index( idx ) );
+        oTrack.field<TracksTag::FTHits>()[idx].template field<TracksTag::LHCbID>().set( seed_track.ft_lhcbID( idx ) );
+      }
+
       // store end of VELO state
       oTrack.field<TracksTag::States>( 0 ).setPosition( state_endvelo.x(), state_endvelo.y(), state_endvelo.z() );
       oTrack.field<TracksTag::States>( 0 ).setDirection( state_endvelo.tx(), state_endvelo.ty() );