MuonR4::FastReco::GlobalPatternFinder Class Reference

Standalone module to handle global pattern recognition. More...

#include <GlobalPatternFinder.h>

Inheritance diagram for MuonR4::FastReco::GlobalPatternFinder:

Collaboration diagram for MuonR4::FastReco::GlobalPatternFinder:

Classes
struct	Config
	Configuration object. More...
struct	HitPayload
	Hit information stored during pattern building. More...
struct	PatternState
	Pattern state object storing pattern information during construction. More...
struct	LineCompatibilityResult
	: Small struct to encapsulate the checkLineCompatibility result More...

Public Types
using	StIndex = Muon::MuonStationIndex::StIndex
	Type alias for the station index.
using	LayerIndex = Muon::MuonStationIndex::LayerIndex
	Type alias for the station layer index.
using	SpacePointContainerVec = std::vector<const SpacePointContainer*>
	Abrivation for a vector of space-point containers.
using	PatternVec = std::vector<GlobalPattern>
	Abrivation for a vector of global patterns.
using	BucketPerContainer = std::unordered_map<const SpacePointContainer*, std::vector<const SpacePointBucket*>>
	Abrivation for a collection of space-point buckets grouped by their corresponding input container.
using	PatternHitVisualInfo = MuonValR4::IFastRecoVisualizationTool::PatternHitVisualInfo
	Type alias for the visual information of a pattern.
using	PatternHitVisualInfoVec = std::vector<PatternHitVisualInfo>
	Abrivation for a vector of visual information objects.

Public Member Functions
	GlobalPatternFinder (const std::string &name, Config &&config)
	Standard constructor.
PatternVec	findPatterns (const ActsTrk::GeometryContext &gctx, const SpacePointContainerVec &spacepoints, BucketPerContainer &outBuckets) const
	Main methods steering the pattern finding.
bool	msgLvl (const MSG::Level lvl) const
	Test the output level.
MsgStream &	msg () const
	The standard message stream.
MsgStream &	msg (const MSG::Level lvl) const
	The standard message stream.
void	setLevel (MSG::Level lvl)
	Change the current logging level.

Private Types
enum class	SeedCoords : std::uint8_t { eSector , eTheta }
	Abrivation of the seed coordinates. More...
enum class	CompatibilityResult : std::int8_t { eAddHit = 0 , eBranchPattern = 1 , eRejectHit = -1 }
	: Enum for the possible outcomes of the line compatibility test of one pattern against one test hit More...
enum class	LayerOrdering : std::int8_t { eSameLayer , eLowerLayer , eHigherLayer }
	Enum to express the logical measurement layer ordering given two hits. More...
using	SearchTree_t = Acts::KDTree<2, HitPayload, double, std::array, 5>
	Definition of the search tree class.
using	TreeNode = std::pair<SearchTree_t::coordinate_t, HitPayload>
	Type alias for a tree node, formed by a hit payload and its indexing coordinates.
using	PatternStateVec = std::vector<PatternState>

Private Member Functions
SearchTree_t	constructTree (const ActsTrk::GeometryContext &gctx, const SpacePointContainerVec &spacepoints) const
	Method to construct the search tree by filling it up with spacepoints from the given containers.
PatternStateVec	findPatternsInEta (const SearchTree_t &orderedSpacepoints, PatternHitVisualInfoVec *visualInfo=nullptr) const
	Method steering the building of patterns in eta.
void	extendPatterns (PatternStateVec &activePatterns, const HitPayload &test, const HitPayload &seed, const HitPayload &prevCandidate, PatternHitVisualInfoVec *visualInfo=nullptr) const
	Function testing pattern compatibility of a set of active patterns (patterns produced from the same seed hit) against one test hit.
LineCompatibilityResult	checkLineCompatibility (const HitPayload &seed, const HitPayload &test, const PatternState &pat) const
	Method to check the line compatibility of a test hit with a given pattern.
LineCompatibilityResult	computeResidual (const HitPayload &seed, const HitPayload &test, const PatternState &pat, const uint8_t patHitIdx, const Amg::Vector3D &beamSpot) const
	Helper method to compute the residual of a test hit against a reference pattern hit and check whether it is within the acceptance window.
bool	isBetter (const PatternState &a, const PatternState &b) const
	Operator to compare two patterns.
PatternStateVec	resolveOverlaps (PatternStateVec &&toResolve, PatternHitVisualInfoVec *visualInfo=nullptr) const
	Method to remove overlapping patterns.
void	addPhiOnlyHits (const ActsTrk::GeometryContext &gctx, PatternStateVec &patterns) const
	Method to add phi-only measurements to existing PatternStates.
bool	isPhiCompatible (const double testPhi, const PatternState &pattern) const
	Method to check the phi compatibility of a test hit with a given pattern.
GlobalPattern	convertToPattern (const PatternState &candidate) const
	Method to convert a PatternState into a GlobalPattern object.
PatternVec	convertToPattern (const PatternStateVec &candidates) const
	Method to convert a vector of PatternStates into GlobalPattern objects.
void	addVisualInfo (const PatternState &candidate, PatternHitVisualInfo::PatternStatus status, PatternHitVisualInfoVec *visualInfo) const
	Helper function to add visual information of a given pattern (which is usually going to be destroyed) to the final container.
void	initMessaging () const
	Initialize our message level and MessageSvc.

Static Private Member Functions
static LayerOrdering	checkLayerOrdering (const HitPayload &hit1, const HitPayload &hit2)
	Method to check the logical layer ordering of two hits.

Private Attributes
SpacePointPerLayerSorter	m_spSorter {}
	Spacepoint sorter per logical measurement layer.
Config	m_cfg
	Global Pattern Recognition configuration.
std::string	m_nm
	Message source name.
boost::thread_specific_ptr< MsgStream >	m_msg_tls
	MsgStream instance (a std::cout like with print-out levels).
std::atomic< IMessageSvc * >	m_imsg { nullptr }
	MessageSvc pointer.
std::atomic< MSG::Level >	m_lvl { MSG::NIL }
	Current logging level.
std::atomic_flag m_initialized	ATLAS_THREAD_SAFE = ATOMIC_FLAG_INIT
	Messaging initialized (initMessaging).

Detailed Description

Standalone module to handle global pattern recognition.

This tool performs global pattern recognition as the first step of the Phase-2 fast reconstruction stage. It builds global patterns of precision and non-precision hits using space-points created in upstream algorithms. It first builds patterns in eta and then adds compatible phi-only hits to the patterns. The resulting patterns are returned by the main method of the tool.

Definition at line 29 of file GlobalPatternFinder.h.

Member Typedef Documentation

BucketPerContainer

using MuonR4::FastReco::GlobalPatternFinder::BucketPerContainer = std::unordered_map<const SpacePointContainer*, std::vector<const SpacePointBucket*>>

Abrivation for a collection of space-point buckets grouped by their corresponding input container.

Definition at line 40 of file GlobalPatternFinder.h.

LayerIndex

using MuonR4::FastReco::GlobalPatternFinder::LayerIndex = Muon::MuonStationIndex::LayerIndex

Type alias for the station layer index.

Definition at line 34 of file GlobalPatternFinder.h.

PatternHitVisualInfo

using MuonR4::FastReco::GlobalPatternFinder::PatternHitVisualInfo = MuonValR4::IFastRecoVisualizationTool::PatternHitVisualInfo

Type alias for the visual information of a pattern.

Definition at line 42 of file GlobalPatternFinder.h.

PatternHitVisualInfoVec

using MuonR4::FastReco::GlobalPatternFinder::PatternHitVisualInfoVec = std::vector<PatternHitVisualInfo>

Abrivation for a vector of visual information objects.

Definition at line 44 of file GlobalPatternFinder.h.

PatternStateVec

using MuonR4::FastReco::GlobalPatternFinder::PatternStateVec = std::vector<PatternState>

private

Definition at line 274 of file GlobalPatternFinder.h.

PatternVec

using MuonR4::FastReco::GlobalPatternFinder::PatternVec = std::vector<GlobalPattern>

Abrivation for a vector of global patterns.

Definition at line 38 of file GlobalPatternFinder.h.

SearchTree_t

using MuonR4::FastReco::GlobalPatternFinder::SearchTree_t = Acts::KDTree<2, HitPayload, double, std::array, 5>

private

Definition of the search tree class.

Definition at line 164 of file GlobalPatternFinder.h.

SpacePointContainerVec

using MuonR4::FastReco::GlobalPatternFinder::SpacePointContainerVec = std::vector<const SpacePointContainer*>

Abrivation for a vector of space-point containers.

Definition at line 36 of file GlobalPatternFinder.h.

StIndex

using MuonR4::FastReco::GlobalPatternFinder::StIndex = Muon::MuonStationIndex::StIndex

Type alias for the station index.

Definition at line 32 of file GlobalPatternFinder.h.

TreeNode

using MuonR4::FastReco::GlobalPatternFinder::TreeNode = std::pair<SearchTree_t::coordinate_t, HitPayload>

private

Type alias for a tree node, formed by a hit payload and its indexing coordinates.

Definition at line 166 of file GlobalPatternFinder.h.

Member Enumeration Documentation

CompatibilityResult

enum class MuonR4::FastReco::GlobalPatternFinder::CompatibilityResult : std::int8_t

strongprivate

: Enum for the possible outcomes of the line compatibility test of one pattern against one test hit

Enumerator
eAddHit	Test successfull, add the hit to the pattern.
eBranchPattern	Test successfull with multiple pattern hits in the same logical measurement layer, branch the pattern.
eRejectHit	Test failed, discard the hit.

Definition at line 303 of file GlobalPatternFinder.h.

                                           : std::int8_t{
                eAddHit = 0,
                eBranchPattern = 1,
                eRejectHit = -1,
            };

LayerOrdering

enum class MuonR4::FastReco::GlobalPatternFinder::LayerOrdering : std::int8_t

strongprivate

Enum to express the logical measurement layer ordering given two hits.

Enumerator
eSameLayer
eLowerLayer
eHigherLayer

Definition at line 371 of file GlobalPatternFinder.h.

                                     : std::int8_t{
                eSameLayer,
                eLowerLayer,
                eHigherLayer
            };

SeedCoords

enum class MuonR4::FastReco::GlobalPatternFinder::SeedCoords : std::uint8_t

strongprivate

Abrivation of the seed coordinates.

Enumerator
eSector	Expanded sector coordinate of the associated spectrometer sector
eTheta	Global Theta.

Definition at line 168 of file GlobalPatternFinder.h.

                                  : std::uint8_t{
                eSector,
                eTheta
            };

Constructor & Destructor Documentation

GlobalPatternFinder()

MuonR4::FastReco::GlobalPatternFinder::GlobalPatternFinder	(	const std::string &	name,
		Config &&	config )

Standard constructor.

Parameters

name	Name to be printed in the messaging
config	Configuration parameters

Definition at line 52 of file GlobalPatternFinder.cxx.

                                                                               :
    AthMessaging{name},
    m_cfg{config} {
        static_assert(std::is_move_assignable_v<PatternState>);
        static_assert(std::is_move_constructible_v<PatternState>);
        static_assert(std::is_copy_assignable_v<PatternState>);
        static_assert(std::is_copy_constructible_v<PatternState>);
        static_assert(std::is_nothrow_move_constructible_v<PatternState>);
        static_assert(std::is_nothrow_move_assignable_v<PatternState>);
};

Member Function Documentation

addPhiOnlyHits()

void MuonR4::FastReco::GlobalPatternFinder::addPhiOnlyHits ( const ActsTrk::GeometryContext & gctx, PatternStateVec & patterns ) const

private

Method to add phi-only measurements to existing PatternStates.

Parameters

gctx	Geometry context
patterns	Vector of pattern states to which to add phi-only hits

Returns

: void

Struct to model the projection of the pattern line onto the phi strip in a certain station

Parameters

refLay	Reference layer coordinate (e.g. R or Z) increasing with the layer number
refStrip	Reference strip coordinate (normal to the layer coordinate, along the strip direction)
invSlope	Inverse slope, defined as deltaStrip / deltaLay

We look for phi-only hits in the buckets associated with the pattern

Definition at line 569 of file GlobalPatternFinder.cxx.

                                                                          {
    constexpr auto covIdxEta {Acts::toUnderlying(SpacePoint::CovIdx::etaCov)};
    struct PhiStripProjectionModel {
        StIndex station{};
        double refLay{0.};
        double refStrip{0.};
        double invSlope{0.};
        bool isBarrel{false};
        bool isValid{false};
 
        double project(const Amg::Vector3D& pos) const {
            const double layCoord   = isBarrel ? pos.perp() : pos.z();
            return refStrip + (layCoord - refLay) * invSlope;
        }
        double residual(const Amg::Vector3D& pos) const {
            const double stripCoord = isBarrel ? pos.z() : pos.perp();
            return std::abs(project(pos) - stripCoord);
        }
    };
    auto makeProjectionModel = [this](const PatternState& pat, const StIndex station) {
        PhiStripProjectionModel result{};
        result.station = station;
        const std::vector<HitPayload>& hits {pat.hitsPerStation.at(station)};
        if (hits.empty()) return result;
        const bool isBarrel {hits.front().hit->msSector()->barrel()};
        // Define the coordinate across the layers, i.e. orthogonal to the strip
        const auto layCoord = [isBarrel](const HitPayload& hit) {
            return (isBarrel) ? hit.R : hit.Z;
        };
        // Define the coordinate along the strip, i.e. orthogonal to the measureent layers
        const auto stripCoord = [isBarrel](const HitPayload& hit) {
            return (isBarrel) ? hit.Z : hit.R;
        };
 
        const HitPayload* sp1 {nullptr};
        const HitPayload* sp2 {nullptr};
        if (hits.size() > 1) {
            // if we have two eta hits in the station, we use the furthestmost to define the pattern line
            for (const HitPayload& hit : hits) {
                if (!hit->measuresEta()) continue;
                if (!sp1 || layCoord(hit) < layCoord(*sp1)) {
                    sp1 = &hit;
                }
                if (!sp2 || layCoord(hit) > layCoord(*sp2)) {
                    sp2 = &hit;
                }
            }
        } else {
            // if we have only one eta hit, we use it and look for the closest eta hit in the closest station to define the pattern line
            sp1 = &hits.front();
            auto layDistanceFromSp1 = [&layCoord,&sp1](const HitPayload& hit) {
                return std::abs(layCoord(hit) - layCoord(*sp1));
            };
            const std::vector<HitPayload>& closestSt {std::ranges::min_element(pat.hitsPerStation, 
                [&layDistanceFromSp1, &station](const auto& st1, const auto& st2){ 
                    if (st1.first == station) return false;
                    if (st2.first == station) return true;
                    return layDistanceFromSp1(st1.second.front()) < layDistanceFromSp1(st2.second.front()); 
            })->second};
            for (const HitPayload& hit : closestSt) {
                if (!hit->measuresEta()) continue;
                if (!sp2 || layDistanceFromSp1(hit) < layDistanceFromSp1(*sp2)) {
                    sp2 = &hit;
                }
            }
        }
        if (!sp1 || !sp2 ) return result;
        const double deltaLay  {layCoord(*sp2) - layCoord(*sp1)};
        if (std::abs(deltaLay) < m_cfg.minLayerSeparation) return result;
 
        result.refLay = layCoord(*sp1);
        result.refStrip = stripCoord(*sp1);
        result.invSlope = (stripCoord(*sp2) - stripCoord(*sp1)) / deltaLay;
        result.isBarrel = isBarrel;
        result.isValid = true;
        return result;
    };
 
    PatternStateVec survivingPatterns{};
    survivingPatterns.reserve(patterns.size());
    for (PatternState& pat : patterns) {
        ATH_MSG_VERBOSE(__func__<<"() Search for phi-only hits for pattern: " << pat);
        // Projection model of pattern line onto a given phi strip
        std::optional<PhiStripProjectionModel> patProjOnStrip{};
        for (const auto& [spc, bucketVec] : pat.bucketsPerContainer) {
            for (const SpacePointBucket* bucket : bucketVec) {
                const Amg::Transform3D& localToGlobal {bucket->msSector()->localToGlobalTransform(gctx)};
                const StIndex station {m_cfg.idHelperSvc->stationIndex(bucket->front()->identify())};
                // If the projection model is not valid, we will use the pattern theta. We cache the local Y in glob frame
                const Amg::Vector3D locY {localToGlobal.linear() * Amg::Vector3D::UnitY()};
                
                for (const auto& hit : *bucket) {
                    // We are looking for phi-only hits
                    if (hit->measuresEta()){
                        continue;
                    }
                    ATH_MSG_VERBOSE(__func__<<"() *** Test phi-only hit "<<*hit);
                    // Check phi compatibility
                    const Amg::Vector3D locPosTest {hit->localPosition()};
                    const Amg::Vector3D globPosTest {localToGlobal * locPosTest};
                    const double globPhi {globPosTest.phi()};
                    if (!isPhiCompatible(globPhi, pat)) {
                        ATH_MSG_VERBOSE(__func__<<"() Phi-only hit not compatible");
                        continue;
                    }
                    // Check there are not other phi hits in the same layer
                    const uint8_t layNum = m_spSorter.sectorLayerNum(*hit);
                    if (auto it {pat.hitsPerStation.find(station)}; it != pat.hitsPerStation.end()) {
                        if (std::ranges::any_of(it->second, [&](const HitPayload& h){
                                return h->measuresPhi() && hit->msSector() == h->msSector() && layNum == h.layerNum; })) {
                            ATH_MSG_VERBOSE(__func__<<"() The pattern already has a phi hit in the same layer - skip this test hit.");
                            continue;
                        }
                    }
                    // Check eta compatibility. Try first to use the pattern line projection model
                    bool isEtaCompatible {false};
                    const double sigmaEta {std::sqrt(hit->covariance()[covIdxEta])};
                    if (!patProjOnStrip.has_value() || patProjOnStrip->station != station) {
                        patProjOnStrip = makeProjectionModel(pat, station);
                    }
                    if(patProjOnStrip->isValid) {
                        isEtaCompatible = patProjOnStrip->residual(globPosTest) <= 1.1*sigmaEta;
                        ATH_MSG_VERBOSE(__func__<<"() Distance pattern line from strip center: "<<patProjOnStrip->residual(globPosTest)<<", strip half-length: "<<sigmaEta<<", isCompatible: "<<isEtaCompatible);
                    } else {
                        // Check pattern theta against global theta at the boundaries of the phi strip
                        double thetaMin {(globPosTest - sigmaEta * locY).theta()};
                        double thetaMax {(globPosTest + sigmaEta * locY).theta()};
                        if (thetaMax < thetaMin) {
                            std::swap(thetaMin, thetaMax);
                        }
                        isEtaCompatible = std::max(pat.theta - thetaMax, thetaMin - pat.theta) < 0.5 * m_cfg.thetaSearchWindow;
                        ATH_MSG_VERBOSE(__func__<<"() Pattern theta "<<inDegrees(pat.theta)<<", strip theta window: ["<<inDegrees(thetaMin)<<", "<<inDegrees(thetaMax)<<"]");
                        
                    }
                    if (!isEtaCompatible) {
                        ATH_MSG_VERBOSE(__func__<<"() The pattern falls outside the test hit strip in eta - skip this test hit.");
                        continue;
                    }
                    // Create the hit payload and add the hit to the pattern
                    ATH_MSG_VERBOSE(__func__<<"() Phi-only hit compatible - add it to the pattern.");
                    pat.addHit(HitPayload{hit.get(), station, layNum, globPhi}, 0., 0.);
                }
            }
        }
        if (pat.nPhiHits < m_cfg.minPhiHits) {
            ATH_MSG_VERBOSE(__func__<<"() Pattern "<< pat<<" has only "<<pat.nPhiHits<<" phi hits, below the minimum required - reject this pattern.");
            continue;
        }
        pat.finalizePatternPhi();
        survivingPatterns.push_back(std::move(pat));
    }
    std::swap(patterns, survivingPatterns);
}

addVisualInfo()

void MuonR4::FastReco::GlobalPatternFinder::addVisualInfo ( const PatternState & candidate, PatternHitVisualInfo::PatternStatus status, PatternHitVisualInfoVec * visualInfo ) const

private

Helper function to add visual information of a given pattern (which is usually going to be destroyed) to the final container.

Parameters

candidate	PatternState whome visual information is to be added
status	Status of the pattern (e.g. successfull, failed or overlap)
visualInfo	Final vector of visual information to store the visual information

Definition at line 1044 of file GlobalPatternFinder.cxx.

                                                                                   {
    if (!visualInfo) {
        return;
    }
    GlobalPattern pattern {convertToPattern(cache)};
    // Check whether the visual info about this pattern is already in the container
    if (auto it =std::ranges::find_if(*visualInfo, [&pattern](const auto& v){ 
        return v.patternCopy && *v.patternCopy == pattern; }); it != visualInfo->end()) {
        it->status = status; // Update the status if the pattern is already in the container
        return;
    }
    visualInfo->push_back(*cache.visualInfo);
    auto& bucketVec {visualInfo->back().parentBuckets};
    for (const auto& [_, buckets] : cache.bucketsPerContainer) {
        bucketVec.insert(bucketVec.end(), buckets.begin(), buckets.end());
    }
    visualInfo->back().patternCopy = std::make_unique<GlobalPattern>(std::move(pattern));
    visualInfo->back().status = status;
}

checkLayerOrdering()

GlobalPatternFinder::LayerOrdering MuonR4::FastReco::GlobalPatternFinder::checkLayerOrdering ( const HitPayload & hit1, const HitPayload & hit2 )

staticprivate

Method to check the logical layer ordering of two hits.

Parameters

hit1	first hit
hit2	second hit

Returns

: the logical measurement layer ordering of the two hits

Hits in the same spectrometer sector

Hits in the same station and different sectors. We can have this case for hits in the overlap region of two adjacent sectors.

Hit in different stations but same station layer. Expected to happen only for Inner and Middle

If both hits are in the middle layer, the one in the barrel comes first

If both hits are in the inner layer, we use the global R, since in large sector BI comes first, while in small sector EI comes first.

If we have one hit in Inner layer for sure it comes first

If we have one hit in Outer layer for sure it comes last

If we have one hit in BarrelExtended and the other in the Middle layer, the former comes first

If we have one hit in Extended (EE) layer and the other in the Middle layer, it depends if the latter is endcap or barrel

Definition at line 845 of file GlobalPatternFinder.cxx.

                                                                {
    auto getLayerOrdering = [](const bool isLayer1Lower) {
        return isLayer1Lower ? eLowerLayer : eHigherLayer;
    };
    if (hit1 == hit2) {
        return eSameLayer;
    }
    if (hit1->msSector() == hit2->msSector()) {
        if (hit1.layerNum == hit2.layerNum) {
            return eSameLayer;
        } else {
            return getLayerOrdering(hit1.layerNum < hit2.layerNum);
        }
    }
    StIndex st1 {hit1.station};
    StIndex st2 {hit2.station};
    if (st1 == st2) {
        return getLayerOrdering(hit1->msSector()->barrel() ? hit1.R < hit2.R : std::abs(hit1.Z) < std::abs(hit2.Z));
    }
    LayerIndex layer1 {toLayerIndex(st1)};
    LayerIndex layer2 {toLayerIndex(st2)};
    if (layer1 == layer2) {
        if (layer1 == LayerIndex::Middle) {
            return getLayerOrdering(st1 == StIndex::BM);
        } else if (layer1 == LayerIndex::Inner) {
            return getLayerOrdering(hit1.R < hit2.R);
        } else {
            throw std::runtime_error("Unexpected to have two pattern-compatible hits one in BO and the other in EO.");
        }
    }
    if (layer1 == LayerIndex::Inner || layer2 == LayerIndex::Inner) {
        return getLayerOrdering(layer1 == LayerIndex::Inner);
    }
    if (layer1 == LayerIndex::Outer || layer2 == LayerIndex::Outer) {
        return getLayerOrdering(layer2 == LayerIndex::Outer);
    }
    if (layer1 == LayerIndex::BarrelExtended || layer2 == LayerIndex::BarrelExtended) {
        return getLayerOrdering(layer1 == LayerIndex::BarrelExtended);
    }
    if (layer1 == LayerIndex::Extended) {
        return getLayerOrdering(st2 == StIndex::EM);
    }
    return getLayerOrdering(st1 == StIndex::BM);
}

checkLineCompatibility()

GlobalPatternFinder::LineCompatibilityResult MuonR4::FastReco::GlobalPatternFinder::checkLineCompatibility ( const HitPayload & seed, const HitPayload & test, const PatternState & pat ) const

private

Method to check the line compatibility of a test hit with a given pattern.

Parameters

seed	seed hit information
test	test hit information
pattern	pattern to be extended

Returns

: result of the test, including the computed line residual and acceptance window

Safety check that the test hit is not on the same layer as the seed

Here we can retrieve the beamspot if desired

Helper function to make the result

Parameters

patHitIdx	Pattern hit index to use to compute the pattern line
isNewLayer	Flag indicating if the test hit is in a new layer

Returns

: Line compatibility result

Fetch the last inserted hit

Test hit coincides with the last inserted hit

We add the test hit to the pattern if they are consecutive MDT hits

sTGCs: if the test hit is a trigger hit (Pad) and the last inserted is precision (strip), we keep the precision hit

If they are not consecutive MDT hits && we have no inserted hits yet (last hit = seed), we give priority to the seed and reject the test hit

Find first previously inserted hit on a different layer from test, starting counting from the second-to-last hit

Definition at line 337 of file GlobalPatternFinder.cxx.

                                                                           {
    // We test hits in the same **expanded** sector, so we need just to compare hit's phi with the pattern's phi, if both available
    if (test->measuresPhi() && pat.nPhiHits && 
        std::abs(CxxUtils::deltaPhi(pat.phi, test.phi)) > m_cfg.phiTolerance) {
        ATH_MSG_VERBOSE(__func__<<"() The pattern with phi = "<<inDegrees(pat.phi)
            <<" is not compatible with the test hit with phi "<<inDegrees(test.phi));
        return LineCompatibilityResult{};
    }
    if (checkLayerOrdering(test, seed) == eSameLayer) {
        ATH_MSG_VERBOSE(__func__<<"() Test hit is on the same layer as the seed - reject.");
        return LineCompatibilityResult{};
    }
    
    const Amg::Vector3D beamSpot{Amg::Vector3D::Zero()};

    auto makeResult = [&](const uint8_t patHitIdx,
                          const bool isNewLayer) -> LineCompatibilityResult {
        LineCompatibilityResult res {computeResidual(seed, test, pat, patHitIdx, beamSpot)};
        if (res.residual < res.accWindow) {
            res.result = isNewLayer ? CompatibilityResult::eAddHit : CompatibilityResult::eBranchPattern;
        }
        return res;
    };
    
    const HitPayload& lastPatHit {pat.getNthLastHit(1u)};
    const int nMeas {pat.nBendingTriggerHits + pat.nPrecisionHits};
 
    /*************** Test hit is on a new layer — draw line from seed to lastHit */
    if(checkLayerOrdering(test, lastPatHit) != eSameLayer) {
        return makeResult(1u, /*isNewLayer*/ true);
    }
    /*************** Test hit is on the same layer as the last inserted hit ***************/
    if (test == lastPatHit) {
        ATH_MSG_VERBOSE(__func__<<"() Test hit is the same as last inserted hit - rejecting.");
        return LineCompatibilityResult{};
    }
    if (areConsecutiveMdt(*test, *lastPatHit)) {
        ATH_MSG_VERBOSE(__func__<<"() Consecutive MDT hits on the same layer - accepting.");
        return LineCompatibilityResult{CompatibilityResult::eAddHit, pat.lastResidual, pat.lastAcceptWindow};
    }
    if (!test.isPrecision && lastPatHit.isPrecision) {
        ATH_MSG_VERBOSE(__func__<<"() Test hit is a trigger hit and last inserted hit is precision on the same layer - keep the precision hit.");
        return LineCompatibilityResult{};
    }
    if (nMeas == 1) {
        ATH_MSG_VERBOSE(__func__<<"() Test hit on same layer as seed with no prior hits, and they are not consecutive MDT hits - rejecting.");
        return LineCompatibilityResult{};
    }
    uint8_t n {2u};
    while (n < nMeas &&
           checkLayerOrdering(test, pat.getNthLastHit(n)) == eSameLayer) {
        ++n;
    }
    return makeResult(n, /*isNewLayer*/ false);
}

computeResidual()

GlobalPatternFinder::LineCompatibilityResult MuonR4::FastReco::GlobalPatternFinder::computeResidual ( const HitPayload & seed, const HitPayload & test, const PatternState & pat, const uint8_t patHitIdx, const Amg::Vector3D & beamSpot ) const

private

Helper method to compute the residual of a test hit against a reference pattern hit and check whether it is within the acceptance window.

The residual is defined against the line passing through seed and the reference pattern hit

Parameters

seed	seed hit information
test	test hit information
pattern	pattern to be extended
patHitIdx	Index of the pattern hit to use for residual computation
beamSpot	Needed to draw the pattern line when the provided pattern hit is too close to the seed

Returns

: result of the test, including the computed line residual and acceptance window

The dynamic acceptance window is defined using the error propagation law for the residual. We have two contributions: the line extrapolation distance and the uncertainty on the line slope

Alpha parameter determining the propagation/ line extrapolation magnitude

Loosen the window when using the beamspot as reference, as the line slope is less well defined

Apply LR correction if the pat hit is a straw

Definition at line 407 of file GlobalPatternFinder.cxx.

                                                                        {
    const HitPayload& patHit {pat.getNthLastHit(patHitIdx)};
    const int nMeas {pat.nBendingTriggerHits + pat.nPrecisionHits};
 
    ATH_MSG_VERBOSE(__func__<<"() Start residual check given pat hit: " << *patHit << ", st: " << stName(patHit.station));
    // Check whether we have to use the beamspot instead of the last pattern hit to draw the line with the seed.
    const double patDeltaZ {std::abs(patHit.Z - seed.Z)};
    const bool useBeamspot {patHit == seed || patDeltaZ < Acts::s_epsilon || 
        (patHit.station == seed.station && (patHit->msSector()->barrel() ? std::abs(patHit.R - seed.R) : patDeltaZ) <= m_cfg.minLayerSeparation)};
    
    if (useBeamspot && patHitIdx + 1 < nMeas) {
        ATH_MSG_VERBOSE(__func__<<"() Distance seed to pat hit is very small - try with next pat hit as reference.");
        return computeResidual(seed, test, pat, patHitIdx + 1u, beamSpot);
    }
 
    /* If the pat hit is a straw, find the consecutive (MDT) hits on the same layer and return use middle one
       for next computations — gives a more central reference for the line direction */
    const HitPayload& centralPatHit { patHit->isStraw() ? [&]() -> const HitPayload& {
        uint8_t nSameLayer{1u};
        while (patHitIdx + nSameLayer < nMeas &&
               checkLayerOrdering(patHit, pat.getNthLastHit(patHitIdx + nSameLayer)) == eSameLayer) {
            ++nSameLayer;
        }
        return nSameLayer > 2u ? pat.getNthLastHit(patHitIdx + (nSameLayer - 1u) / 2u) : patHit;
    }() : patHit};
 
    ATH_MSG_VERBOSE(__func__<<"() Distance seed to pat hit - deltaR_pat= " << (centralPatHit.R - seed.R) << ", deltaZ_pat= " << (centralPatHit.Z - seed.Z) << ". Use beamspot: " << useBeamspot);
    LineCompatibilityResult res {};
 
    // Determine the coordinates to use as reference for the line computation
    const double refR {useBeamspot ? beamSpot.perp() : centralPatHit.R};
    const double refZ {useBeamspot ? beamSpot.z()    : centralPatHit.Z};
 
    // Compute the line slope. 
    const double dZ_slope {refZ - seed.Z};
    const double dR_slope {refR - seed.R};
    assert(dZ_slope != 0);
    const double lineSlope {dR_slope / dZ_slope};
    ATH_MSG_VERBOSE(__func__<<"() Slope computation - deltaR_slope= " << dR_slope << ", deltaZ_slope= " << dZ_slope << ", slope= " << lineSlope);
 
    // Compute the residual
    const double dZ_res {test.Z - seed.Z};
    res.residual = (test.R - seed.R) - lineSlope * dZ_res;
    ATH_MSG_VERBOSE(__func__<<"() Residual computation - deltaR_res: " << (test.R - seed.R) << ", deltaZ_res: " << dZ_res << ", signed residual: " << res.residual);

    const double geometricalFactor {1.+ Acts::square(lineSlope)};
    const double alpha {dZ_res / dZ_slope};  
    const double propagationFactor {1. - alpha + Acts::square(alpha)};
    res.accWindow = m_cfg.baseRWindow * std::sqrt(2*geometricalFactor * propagationFactor);
    if (useBeamspot) res.accWindow *= 1.5; 
    ATH_MSG_VERBOSE(__func__<<"() Window computation - Geometrical Factor: " << std::sqrt(geometricalFactor) << " alpha: " << alpha
                    << ", Propagation Factor: " << std::sqrt(propagationFactor) << ", Computed Window: " << res.accWindow);
        
    if (!useBeamspot && centralPatHit->isStraw()) {
        const double LRcorrection {dZ_res * geometricalFactor * centralPatHit->driftRadius()/ Acts::fastHypot(dZ_slope, dR_slope)};
        res.residual = std::min(std::abs(res.residual-LRcorrection), std::abs(res.residual+LRcorrection));
        ATH_MSG_VERBOSE(__func__<<"() Apply LR correction: " << LRcorrection << ", corrected residual: " << res.residual);
    } else {
        res.residual = std::abs(res.residual);
    }
    if (pat.visualInfo) {
        pat.visualInfo->hitLineInfo[test.hit] = std::make_pair(lineSlope, res.accWindow);
    }
    return res;
}

constructTree()

GlobalPatternFinder::SearchTree_t MuonR4::FastReco::GlobalPatternFinder::constructTree ( const ActsTrk::GeometryContext & gctx, const SpacePointContainerVec & spacepoints ) const

private

Method to construct the search tree by filling it up with spacepoints from the given containers.

The tree does not contain only-phi hits.

Parameters

gctx	Geometry context
spacepoints	Vector of space point containers

Returns

: Constructed search tree

Try to duplicate the hit in the neighboring sectors if it is close to the sector border. This ensures that we can find patterns crossing the sector borders.

Check whether the hit belongs to the left or right sector as well

Blow-up the number of sectors by a factor of 2. The even numbers represent the segments expressed @ the sector centre. The odd numbers represent the overlap region between two adjacent sectors. For sector 16, the right overlap region is mapped to 1

Definition at line 782 of file GlobalPatternFinder.cxx.

                                                                                    {
    SearchTree_t::vector_t rawData{};
    using SectorProjector = MsTrackSeeder::SectorProjector;
    // Before the loops: estimate the total number of hits 
    size_t totalHits = 0;
    for (const SpacePointContainer* spc : spacepoints) {
        for (const SpacePointBucket* bucket : *spc) {
            totalHits += bucket->size();
        }
    }
    // We can have up to 3 entries per hit (when the hit does not measure phi).
    rawData.reserve(3 * totalHits);
 
    for (const SpacePointContainer* spc : spacepoints) {
        ATH_MSG_VERBOSE("Adding to the search tree "<<spc->size()<<" space point buckets");
        for (const SpacePointBucket* bucket : *spc) {
            ATH_MSG_VERBOSE("Processing " << bucket->size() << " spacepoint...");
            const Amg::Transform3D& localToGlobal {bucket->msSector()->localToGlobalTransform(gctx)};
            const StIndex bucketStation {m_cfg.idHelperSvc->stationIndex(bucket->front()->identify())};
            const int sector = bucket->msSector()->sector();
 
            for (const auto& hit : *bucket) {
                // Ignore only-phi hits and MDT hits if desired
                if (!hit->measuresEta() || (!m_cfg.useMdtHits && hit->isStraw())) {
                    continue;
                }
                const Amg::Vector3D globalPos {localToGlobal * hit->localPosition()};
                const double globalTheta {globalPos.theta()};
                const HitPayload newHit {hit.get(), bucket, spc, bucketStation,
                    static_cast<uint8_t>(m_spSorter.sectorLayerNum(*hit)),globalPos.perp(), globalPos.z(), globalPos.phi()};
                 
                for (const auto proj : {SectorProjector::leftOverlap, SectorProjector::center, SectorProjector::rightOverlap}) {
                    const int projSector = MsTrackSeeder::ringSector(sector + Acts::toUnderlying(proj));
                    if (hit->measuresPhi() && proj != SectorProjector::center && 
                        !sectorMap.insideSector(projSector, newHit.phi)) {
                        ATH_MSG_VERBOSE(__func__<<"() Hit @"<< *hit<<"\nwith globPhi "<<inDegrees(newHit.phi)
                            <<" is not in sector "<<projSector<<" ["<<inDegrees(sectorMap.sectorPhi(projSector)-sectorMap.sectorWidth(projSector))
                            <<", "<<inDegrees(sectorMap.sectorPhi(projSector)+sectorMap.sectorWidth(projSector))
                            <<"] which is "<<MsTrackSeeder::to_string(proj) <<" to "<< sector);
                        continue;
                    }
                    std::array<double, 2> coords{};
                    coords[Acts::toUnderlying(SeedCoords::eTheta)] = globalTheta;

                    coords[Acts::toUnderlying(SeedCoords::eSector)] = MsTrackSeeder::ringOverlap(2*sector + Acts::toUnderlying(proj));
                    ATH_MSG_VERBOSE(__func__<<"() Add hit @"<< *hit 
                                            <<"\nwith global position "<<Amg::toString(globalPos) <<" and coordinates "<<coords<<" to the search tree");
                    rawData.emplace_back(std::move(coords), newHit);
                }  
            }
        }
    }
    ATH_MSG_VERBOSE("Create a new tree with "<<rawData.size()<<" entries. ");
    return SearchTree_t{std::move(rawData)};
}

convertToPattern() [1/2]

GlobalPattern MuonR4::FastReco::GlobalPatternFinder::convertToPattern ( const PatternState & candidate ) const

private

Method to convert a PatternState into a GlobalPattern object.

Parameters

candidate

PatternState to be converted

Returns

: Converted GlobalPattern

Definition at line 748 of file GlobalPatternFinder.cxx.

                                                                                   {
    GlobalPattern::HitCollection hitPerStation{};
    for (const auto& [station, hits] : cache.hitsPerStation) {
        auto& out = hitPerStation[station];
        out.reserve(hits.size());
        std::ranges::transform(hits, std::back_inserter(out), [](const HitPayload& h){ return h.hit;});
    }
    GlobalPattern pattern{std::move(hitPerStation)};
    pattern.setTheta(cache.theta);
    pattern.setPhi(cache.phi);
    // Set the pattern sector(s) and theta.
    pattern.setSector(cache.sector1);
    pattern.setSecondarySector(cache.sector2);
    // Set pattern quality information.
    pattern.setNPrecisionHits(cache.nPrecisionHits);
    pattern.setNEtaNonPrecisionHits(cache.nBendingTriggerHits);
    pattern.setNPhiHits(cache.nPhiHits);
    pattern.setMeanNormResidual2(cache.meanNormResidual2);
    return pattern;
}

convertToPattern() [2/2]

GlobalPatternFinder::PatternVec MuonR4::FastReco::GlobalPatternFinder::convertToPattern ( const PatternStateVec & candidates ) const

private

Method to convert a vector of PatternStates into GlobalPattern objects.

Parameters

candidates

PatternStates to be converted

Returns

: Vector of converted GlobalPatterns

Definition at line 770 of file GlobalPatternFinder.cxx.

                                                                        {
    PatternVec patterns{};
    patterns.reserve(cache.size());
    std::transform(cache.begin(), cache.end(), std::back_inserter(patterns), [this](const PatternState& cacheEntry) {
        GlobalPattern pattern {convertToPattern(cacheEntry)};
        ATH_MSG_VERBOSE("convertToPattern() Converted new pattern "<<pattern);
        return pattern;
    });
    return patterns;
}

extendPatterns()

void MuonR4::FastReco::GlobalPatternFinder::extendPatterns ( PatternStateVec & activePatterns, const HitPayload & test, const HitPayload & seed, const HitPayload & prevCandidate, PatternHitVisualInfoVec * visualInfo = nullptr ) const

private

Function testing pattern compatibility of a set of active patterns (patterns produced from the same seed hit) against one test hit.

At the end, activePatterns contains the surviving patterns.

Parameters

activePatterns	Vector of active patterns to be extended
test	Hit to be tested
seed	Seed hit information
prevCandidate	Previous candidate hit information
visualInfo	Pointer to visual information for pattern visualization (nullptr if the VisualizationTool is disabled)

Check angular compatibility of the test hit and the pattern

Add the test hit to the pattern

Branch the pattern: we clone it and overwrite the existing hit with the test hit

Update visual information of the original pattern

Add the new pattern to the list of next patterns

Definition at line 240 of file GlobalPatternFinder.cxx.

                                                                                    {
    // Filter patterns: deduplicate groups sharing the same last hit, but only when the next 
    // test hit is on a different layer — if it's on the same layer we might branch again and need both
    if (activePatterns.size() > 1) {
        // Sort pointers to avoid moving pattern states during deduplication
        std::vector<PatternState*> patPtrs(activePatterns.size());
        std::ranges::transform(activePatterns, patPtrs.begin(), [](PatternState& p){ return &p; });
        // Sort patteres by last inserted hit to find groups of patterns sharing the same last hit
        std::ranges::sort(patPtrs, {}, [](const PatternState* p){ return p->lastInsertedHit; });
 
        PatternStateVec deduplicated{};
        deduplicated.reserve(activePatterns.size());
        for (auto it = patPtrs.begin(); it != patPtrs.end(); ) {
            // Find the group of patterns sharing the same last hit
            const SpacePoint* groupLastHit {(*it)->lastInsertedHit};
            auto groupEnd = std::ranges::find_if(it, patPtrs.end(), [&](const PatternState* p) {
                return p->lastInsertedHit != groupLastHit;
            });
            if (checkLayerOrdering(test, (*it)->getNthLastHit(1u)) != eSameLayer) {
                // If last hit of the group is on a different layer than the test hit, we deduplicate and pick the best pattern in the group
                deduplicated.push_back(std::move(**std::ranges::max_element(it, groupEnd, 
                    [this](const PatternState* a, const PatternState* b){ return isBetter(*b, *a); })));
            } else {
                // If the last hit of the group is on the same layer as the test hit, we keep all patterns in the group, as we might branch again with the next hits.
                std::ranges::transform(it, groupEnd, std::back_inserter(deduplicated), 
                    [](PatternState* p){ return std::move(*p); });
            }
            it = groupEnd;
        }
        std::swap(activePatterns, deduplicated);
    }
    PatternStateVec nextPatterns{};
    nextPatterns.reserve(activePatterns.size()* 2);
 
    const unsigned refMissedLayerHits {std::ranges::min_element(activePatterns,
            [](const auto& a, const auto& b) { return a.nMissedLayerHits < b.nMissedLayerHits; }
        )->nMissedLayerHits};
    
    for (PatternState& pattern : activePatterns) {
 
        if (pattern.nMissedLayerHits >= m_cfg.maxMissedLayerHits && pattern.nMissedLayerHits > refMissedLayerHits) {
            ATH_MSG_VERBOSE(__func__<<"() Pattern " << pattern << " has missed " << pattern.nMissedLayerHits << " hits in different layers, above the maximum allowed - abort this pattern.");
            addVisualInfo(pattern, PatternHitVisualInfo::PatternStatus::eFailed, visualInfo);
            continue;
        }

        const auto [result, residual, accWindow] {checkLineCompatibility(seed, test, pattern)};
        switch (result) {
            case CompatibilityResult::eAddHit: {
                pattern.addHit(test, residual, accWindow);
                ATH_MSG_VERBOSE(__func__<<"() Hit is compatible, add it to the pattern. Updated pattern: " << pattern);
                break;
            }
            case CompatibilityResult::eBranchPattern: {
                /* Check first if the branched pattern already exists*/
                if (std::ranges::any_of(nextPatterns, [&test, &pattern](const PatternState& p) {
                        return p.lastInsertedHit == test.hit && p.prevLayerHit == pattern.prevLayerHit; })) {
                    ATH_MSG_VERBOSE(__func__<<"() Hit is compatible & on same layer of last added hit, but branched pattern already exists");
                    break;
                }
                PatternState newPattern {pattern};
                const HitPayload& lastPatHit {pattern.getNthLastHit(1u)};
                newPattern.overWriteHit(lastPatHit, test, residual, accWindow);
                ATH_MSG_VERBOSE(__func__<<"() Hit is compatible & on same layer of last added hit - new branched pattern: " << newPattern);

                if (newPattern.visualInfo && pattern.visualInfo) {
                    pattern.visualInfo->replacedHits.push_back(test.hit);
                    newPattern.visualInfo->replacedHits.push_back(lastPatHit.hit);
                }
                nextPatterns.push_back(std::move(newPattern));
                break;
            }
            case CompatibilityResult::eRejectHit: {
                ATH_MSG_VERBOSE(__func__<<"() Hit is not compatible with the pattern.");
                if (pattern.visualInfo) {
                    pattern.visualInfo->discardedHits.push_back(test.hit);
                }
                if (checkLayerOrdering(test, prevCandidate) != eSameLayer){
                    pattern.nMissedLayerHits++;
                }
                break;
            }
        }
        nextPatterns.push_back(std::move(pattern));
    }
    std::swap(activePatterns, nextPatterns);
};

findPatterns()

GlobalPatternFinder::PatternVec MuonR4::FastReco::GlobalPatternFinder::findPatterns	(	const ActsTrk::GeometryContext &	gctx,
		const SpacePointContainerVec &	spacepoints,
		BucketPerContainer &	outBuckets ) const

Main methods steering the pattern finding.

Given the space-point containers, it creates the search tree,
builds patterns in eta, attach compatible only-phi measurements, and convert PatternStates into GlobalPatterns

Parameters

gctx	Geometry context
spacepoints	Vector of space point containers
outBuckets	output collection of space-point buckets grouped by their corresponding input container to be written into StoreGate.

Returns

: Vector of found patterns

Create the search tree by ordering hits in theta and expanded spectrometer sector and find patterns in eta

Add phi-only hits to the patterns

Fill the output buckets & visual info

Add the successfull pattern to visual info, as we won't touch it again

Plot patterns

Definition at line 65 of file GlobalPatternFinder.cxx.

                                                                        {
    auto visualInfo {m_cfg.visionTool ? std::make_unique<PatternHitVisualInfoVec>() : nullptr};
    PatternStateVec patterns{findPatternsInEta(constructTree(gctx, spacepoints), visualInfo.get())};

    addPhiOnlyHits(gctx, patterns);

    for (const PatternState& pat : patterns) {
        addVisualInfo(pat, PatternHitVisualInfo::PatternStatus::eSuccessful, visualInfo.get());
 
        for (const auto& [spc, bucketVec] : pat.bucketsPerContainer) {
            if (outBuckets.find(spc) == outBuckets.end()) {
                throw std::runtime_error("The space point container associated to the pattern is not present in the output bucket map.");
            }
            for (const SpacePointBucket* bucket : bucketVec) {
                auto& outBucketVec = outBuckets[spc];
                if (std::ranges::find(outBucketVec, bucket) == outBucketVec.end()) {
                    outBucketVec.push_back(bucket);
                }
            }
        }
    }
    if (visualInfo) {
        m_cfg.visionTool->plotPatternBuckets(Gaudi::Hive::currentContext(), "GlobPatFind_", std::move(*visualInfo));
    }
    return convertToPattern(patterns);
}

findPatternsInEta()

GlobalPatternFinder::PatternStateVec MuonR4::FastReco::GlobalPatternFinder::findPatternsInEta ( const SearchTree_t & orderedSpacepoints, PatternHitVisualInfoVec * visualInfo = nullptr ) const

private

Method steering the building of patterns in eta.

Parameters

orderedSpacepoints	Search tree with spacepoints ordered by their corresponding coordinates
visualInfo	Pointer to visual information for pattern visualization (nullptr if the VisualizationTool is disabled). Needed to add visual info about pattern candidates discarded during the building stage.

Returns

: resulting vector of PatternStates successfully built

We try to build a pattern in eta starting from every hit in the three

Check the seed is in the current seeding layer, and if seeding from MDT hits is enabled

check how many existing (and overlapping) patterns contain this hit

Define the search range.

Search hits with same expanded sector

Define theta window size. While a middle-layer seed already constrains the track direction more tightly (the line must connect to hits on both sides), inner- and outer-layer seeds are more loosely constraining, and we need a larger theta window with double size to achieve the same angular acceptance as inner/outer seeds

Search for compatible spacepoints with the seed and check if there are enough to build a pattern

Check that the candidate hits extend at least in two layers

Sort the compatible spacepoints by global logical layer

If one hit is precision and one trigger, put first the trigger one. Relevant for STGCs

If the two hits are in the same layer, sort them by local y coordinate.

Initialize a pattern from the seed hit

Helper function to process a new candidate hit. During pattern building, we can have pattern branching when the initial pattern is compatible with multiple hits in the same layer. These patterns will be stored in activePatterns. For each new hit, extendPatterns
will try to extend every active pattern and remove the ones not meeting continuation criteria.

Parameters

testPair

The new candidate hit to process.

Helper function ensuring we have at the end one pattern for seed hit. If multiple survived at the end, we keep the best.

First search for compatible hits from the seed layer onwards

Then try to proceed toward the innermost layer

Check that the pattern meets the minimum requirements for trigger and precision hits

Definition at line 100 of file GlobalPatternFinder.cxx.

                                                                                  {
    constexpr auto thetaIdx {Acts::toUnderlying(SeedCoords::eTheta)};
    constexpr auto sectorIdx {Acts::toUnderlying(SeedCoords::eSector)};
    PatternStateVec patterns{};
    using enum SeedCoords;
    for (const auto seedingLayer : m_cfg.layerSeedings) {
        ATH_MSG_VERBOSE(__func__<<"() Start pattern search in eta with seeding layer "<< layerName(seedingLayer));
        for (const auto& [seedCoords, seed] : orderedSpacepoints) {
            const LayerIndex seedLayer {toLayerIndex(seed.station)};
            if (seedLayer != seedingLayer || (seed->isStraw() && !m_cfg.seedFromMdt)) {
                continue;
            }
            ATH_MSG_VERBOSE(__func__<<"() New seed hit "<<*seed<<", coordinates "<<seedCoords);
            auto nExistingPatterns {std::ranges::count_if(patterns, [&seed, &seedCoords, this](const PatternState& pattern){
                if (std::abs(pattern.theta - seedCoords[thetaIdx]) > m_cfg.thetaSearchWindow ||
                    pattern.sectorCoord != static_cast<int>(seedCoords[sectorIdx])) {
                    return false;
                }
                return pattern.isInPattern(seed);
            })};
            if (nExistingPatterns >= m_cfg.maxSeedAttempts) {
                ATH_MSG_VERBOSE(__func__<<"() Seed has already been used in "<<nExistingPatterns<<" patterns, which is above the limit - skip this seed.");
                continue;
            }    
            SearchTree_t::range_t selectRange{};
            selectRange[sectorIdx].shrink(seedCoords[sectorIdx] - 0.1, seedCoords[sectorIdx] + 0.1); 
            const double thetaWindow {(seedLayer == LayerIndex::Inner || seedLayer == LayerIndex::Outer) 
                                    ? m_cfg.thetaSearchWindow : 0.5*m_cfg.thetaSearchWindow};
            selectRange[thetaIdx].shrink(seedCoords[thetaIdx] - thetaWindow, seedCoords[thetaIdx] + thetaWindow);
            std::vector<TreeNode> candidateHits = orderedSpacepoints.rangeSearchWithKey(selectRange);
            if (candidateHits.size() < m_cfg.minBendingTriggerHits + m_cfg.minBendingPrecisionHits) {
                ATH_MSG_VERBOSE(__func__<<"() Found "<<candidateHits.size()<<" candidate hits, below the minimum required - skip this seed.");
                continue;
            }
            if (std::ranges::none_of(candidateHits, [this, seedLayer](const TreeNode& c){
                    return m_cfg.idHelperSvc->layerIndex(c.second->identify()) != seedLayer; }) ) {
                ATH_MSG_VERBOSE(__func__<<"() All candidates are in the same station layer, and we need at least two - skip this seed.");
                continue;
            }
            std::ranges::sort(candidateHits, [](const TreeNode& c1, const TreeNode& c2){
                const HitPayload& hit1 {c1.second};
                const HitPayload& hit2 {c2.second};
                LayerOrdering ordering {checkLayerOrdering(hit1, hit2)};
                if (ordering == eSameLayer) {
                    if (hit1.isPrecision != hit2.isPrecision) {
                        return hit2.isPrecision;
                    }
                    return hit1->localPosition().y() < hit2->localPosition().y();
                }
                return ordering == eLowerLayer;
            });
            if (msgLvl(MSG::VERBOSE)) {
                ATH_MSG_VERBOSE(__func__<<"() Found "<< candidateHits.size()<<" candidate hits: ");
                for (const auto& [coords, hit] : candidateHits) {
                    ATH_MSG_VERBOSE(__func__<<"() \t**"<<*hit<<", coords: "<<coords << ", glob Z/R/phi: "<<hit.Z<<" / "<<hit.R<<" / "<<inDegrees(hit.phi) 
                        << ", st/layer: " << stName(hit.station) << " / "<< hit.layerNum);
                }
            }
            const auto seedItr {std::ranges::find_if(candidateHits,
                [&seed](const TreeNode& c){ return c.second == seed; })};
            assert(seedItr != candidateHits.end());
            std::vector<PatternState> activePatterns{};
            activePatterns.emplace_back(seed, static_cast<int>(seedCoords[sectorIdx]), seedCoords[thetaIdx]);
            if (visualInfo) {
                activePatterns.back().visualInfo = std::make_unique<PatternHitVisualInfo>(
                    seed.hit, seedCoords[thetaIdx] - thetaWindow, seedCoords[thetaIdx] + thetaWindow);
            }
 
            const HitPayload* prevCandidate {&seed};
            auto processNewHit = [&](const TreeNode& testPair){
                const HitPayload& test {testPair.second};
                ATH_MSG_VERBOSE("processNewHit() *** Test "<<*test<<" against " << activePatterns.size() << " active patterns.");
                extendPatterns(activePatterns, test, seed, *prevCandidate, visualInfo);
                prevCandidate = &test;
            };
            auto ensureOnePattern = [&activePatterns, &visualInfo, this]() {
                if (activePatterns.size() > 1) {
                    ATH_MSG_VERBOSE("ensureOnePattern() Found "<<activePatterns.size()<<" patterns in current search stage - keep the best one.");
                    activePatterns = resolveOverlaps(std::move(activePatterns), visualInfo);
                    resetVisualToOverlap(visualInfo);
                    assert(activePatterns.size() == 1);
                }
            };

            ATH_MSG_VERBOSE(__func__<<"() Search between " << candidateHits.size() <<" candidates for compatible hits...");
            for (auto it = std::next(seedItr); it != candidateHits.end(); ++it) {
                processNewHit(*it);
            }
            ensureOnePattern();
            for (auto it = std::reverse_iterator(seedItr); it != candidateHits.rend(); ++it) {
                processNewHit(*it);
            }
            ensureOnePattern();
            PatternState& newPattern {activePatterns.back()};
            if (newPattern.nBendingTriggerHits < m_cfg.minBendingTriggerHits || 
                newPattern.nPrecisionHits < m_cfg.minBendingPrecisionHits ||
                std::ranges::count_if(newPattern.hitsPerStation, [](const auto& st) { return st.second.size()>= 3; }) < 2) {
                ATH_MSG_VERBOSE(__func__<<"() Pattern " << newPattern << "\ndoes not meet minimum hit requirements - reject.");
                addVisualInfo(newPattern, PatternHitVisualInfo::PatternStatus::eFailed, visualInfo);
                continue;
            }
            newPattern.meanNormResidual2 /= (newPattern.nPrecisionHits + newPattern.nBendingTriggerHits);
            if (newPattern.meanNormResidual2 > m_cfg.meanNormRes2Cut) {
                ATH_MSG_VERBOSE(__func__<<"() Pattern " << newPattern << "\ndoes not meet the mean norm residual2 cut - reject.");
                addVisualInfo(newPattern, PatternHitVisualInfo::PatternStatus::eFailed, visualInfo);
                continue;
            }
            newPattern.isFinalized = true;
            ATH_MSG_VERBOSE(__func__<<"() Add new pattern "<<newPattern);
            patterns.push_back(std::move(newPattern));
        }
    }
    ATH_MSG_VERBOSE("Found in total "<<patterns.size()<<" patterns in eta before overlap removal");
    return resolveOverlaps(std::move(patterns), visualInfo);
}

initMessaging()

void AthMessaging::initMessaging ( ) const

privateinherited

Initialize our message level and MessageSvc.

This method should only be called once.

Definition at line 39 of file AthMessaging.cxx.

{
  m_imsg = Athena::getMessageSvc();
  // If user did not set an explicit level, set a default
  if (m_lvl == MSG::NIL) {
    m_lvl = m_imsg ?
      static_cast<MSG::Level>( m_imsg.load()->outputLevel(m_nm) ) :
      MSG::INFO;
  }
}

isBetter()

bool MuonR4::FastReco::GlobalPatternFinder::isBetter ( const PatternState & a, const PatternState & b ) const

private

Operator to compare two patterns.

Parameters

a	First pattern to compare
b	Second pattern to compare

Returns

: true if pattern a is better than pattern b, false otherwise

Function computing the score for a pattern

Definition at line 479 of file GlobalPatternFinder.cxx.

                                                                {
    auto score = [this](const PatternState& p) {
        const long nStations {std::ranges::count_if(p.hitsPerStation, [](const auto& pair){ return pair.second.size() >= 4u; })};
        const double etaHitCount  {1.*p.nBendingTriggerHits + m_cfg.precisionWeight * p.nPrecisionHits};
        const double etaHitScore  {std::tanh(etaHitCount / (std::max(nStations, 1L) * m_cfg.hitScoreSaturation))};
        const double stationScore {std::tanh(1.*nStations / 2.)};
        double residualRatio      {p.meanNormResidual2 / (1.5*m_cfg.meanNormRes2Cut)};
        if (!p.isFinalized) residualRatio /= (p.nBendingTriggerHits + p.nPrecisionHits);
        const double resPenalty {m_cfg.residualPenalty * residualRatio  / (1. + residualRatio)};
        const double etaScore {stationScore * etaHitScore * (1. - resPenalty)};
        const double phiScore {std::tanh(1.*p.nPhiHits / m_cfg.phiBonusSaturation)};
        ATH_MSG_VERBOSE("score() stationScore: "<<stationScore<<", etaHitScore: "<<etaHitScore<<", resPenalty: "<<resPenalty<<", etaScore: "<<etaScore<<", phiScore: "<<phiScore);
        return 0.9 * etaScore + 0.1 * phiScore;
    };
    const double scoreA {score(a)};
    const double scoreB {score(b)};
    ATH_MSG_VERBOSE(__func__<<"() Pattern " << a << " with score " << scoreA << " is " <<
        (scoreA > scoreB ? "BETTER" : "WORSE") << " than " << b << " with score " << scoreB);
    return scoreA > scoreB;                              
}

isPhiCompatible()

bool MuonR4::FastReco::GlobalPatternFinder::isPhiCompatible ( const double testPhi, const PatternState & pattern ) const

private

Method to check the phi compatibility of a test hit with a given pattern.

Parameters

testPhi	test global phi
pattern	pattern to be extended

Returns

: true if the test hit is phi compatible with the pattern, false otherwise

We check that the test hit is compatible with the pattern phi, if available, which is given by the first phi measurement in the pattern. If the pattern doesn't have a phi yet, we check that the test hit is in the same pattern sector(s)

Definition at line 728 of file GlobalPatternFinder.cxx.

                                                                             {
    if (pattern.nPhiHits) {
        if (std::abs(CxxUtils::deltaPhi(pattern.phi, testPhi)) > m_cfg.phiTolerance) {
            ATH_MSG_VERBOSE(__func__<<"() The pattern with phi = "<<inDegrees(pattern.phi)<<" is not compatible with the test hit with phi "<<inDegrees(testPhi));
            return false;
        }
    } else {
        const bool isCompatible {pattern.sector1 == pattern.sector2 ? sectorMap.insideSector(pattern.sector1, testPhi) : 
                sectorMap.insideSector(pattern.sector1, testPhi) && sectorMap.insideSector(pattern.sector2, testPhi)};
        if (!isCompatible) {
            ATH_MSG_VERBOSE(__func__<<"() The test hit with phi = "<<inDegrees(testPhi)<<" is not inside the pattern sectors: "<<pattern.sector1<<" and "<<pattern.sector2);
            return false;
        }            
    }
    return true;
}

msg() [1/2]

MsgStream & AthMessaging::msg ( ) const

inlineinherited

The standard message stream.

Returns a reference to the default message stream May not be invoked before sysInitialize() has been invoked.

Definition at line 167 of file AthMessaging.h.

{
  MsgStream* ms = m_msg_tls.get();
  if (!ms) {
    if (!m_initialized.test_and_set()) initMessaging();
    ms = new MsgStream(m_imsg,m_nm);
    m_msg_tls.reset( ms );
  }
 
  ms->setLevel (m_lvl);
  return *ms;
}

msg() [2/2]

MsgStream & AthMessaging::msg ( const MSG::Level lvl ) const

inlineinherited

The standard message stream.

Returns a reference to the default message stream May not be invoked before sysInitialize() has been invoked.

Definition at line 182 of file AthMessaging.h.

{ return msg() << lvl; }

msgLvl()

bool AthMessaging::msgLvl ( const MSG::Level lvl ) const

inlineinherited

Test the output level.

Parameters

lvl	The message level to test against

Returns

boolean Indicating if messages at given level will be printed

Return values

true	Messages at level "lvl" will be printed

Definition at line 151 of file AthMessaging.h.

{
  // If user did not set explicit message level we have to initialize
  // the messaging and retrieve the default via the MessageSvc.
  if (m_lvl==MSG::NIL && !m_initialized.test_and_set()) initMessaging();
 
  if (m_lvl <= lvl) {
    msg() << lvl;
    return true;
  } else {
    return false;
  }
}

resolveOverlaps()

GlobalPatternFinder::PatternStateVec MuonR4::FastReco::GlobalPatternFinder::resolveOverlaps ( PatternStateVec && toResolve, PatternHitVisualInfoVec * visualInfo = nullptr ) const

private

Method to remove overlapping patterns.

Parameters

toResolve	pattern to be resolved
visualInfo	Pointer to visual information for pattern visualization (nullptr if the VisualizationTool is disabled)

Returns

: resolved patterns

Check if two patterns overlap in space

Check first the geometrical overlap

If we reach here, the patterns can overlap geometrically, so check the hit content

Overlap if more than half the hits of the smaller pattern are shared

Definition at line 502 of file GlobalPatternFinder.cxx.

                                                                                {
    PatternStateVec outputPatterns{};
    outputPatterns.reserve(toResolve.size());
    auto areOverlapping = [this](const PatternState& a, const PatternState& b) {
        if(std::abs(a.sectorCoord - b.sectorCoord) > 1) {
            return false;
        }
        if (a.nPhiHits > 0 && b.nPhiHits > 0) {
            if (std::abs(CxxUtils::deltaPhi(a.phi, b.phi)) > 2.*m_cfg.phiTolerance) return false;
        } else if (a.nPhiHits > 0) {
            if (!sectorMap.insideSector(b.sector1, a.phi) || !sectorMap.insideSector(b.sector2, a.phi)) {
                return false;
            }
        } else if (b.nPhiHits > 0) {
            if (!sectorMap.insideSector(a.sector1, b.phi) || !sectorMap.insideSector(a.sector2, b.phi)) {
                return false;
            }
        }
        if (std::abs(a.theta - b.theta) > 2.*m_cfg.thetaSearchWindow) {
            return false;
        }
        int nSharedHits{0};
        for (const auto& [stA, hitsA] : a.hitsPerStation) {
            auto itB = b.hitsPerStation.find(stA);
            if (itB == b.hitsPerStation.end()) {
                continue;
            }
            nSharedHits += std::ranges::count_if(hitsA, [&](const HitPayload& hitA){
                return std::ranges::find(itB->second, hitA) != itB->second.end();
            });
        }
        const int minHits {std::min(a.nBendingTriggerHits + a.nPrecisionHits, b.nBendingTriggerHits + b.nPrecisionHits)};
        return nSharedHits > minHits / 2;
    };
    for (auto it = toResolve.begin(); it != toResolve.end(); ++it) {
        if (it->isOverlap) {
            // If already marked as overlap, add to visual info, and discard the pattern
            addVisualInfo(*it, PatternHitVisualInfo::PatternStatus::eOverlap, visualInfo);
            continue;
        }
        for (auto jt = std::next(it); jt != toResolve.end(); ++jt) {
            if (jt->isOverlap || !areOverlapping(*it, *jt)) {
                continue;
            }
            if (isBetter(*it, *jt)) {
                jt->isOverlap = true;
            } else {
                it->isOverlap = true;
                break;
            }
        }
        if (!it->isOverlap) {
            it->finalizePatternEta();
            outputPatterns.push_back( std::move(*it));
        } else {
            // If overlap, add to visual info, as the pattern will be discarded
            addVisualInfo(*it, PatternHitVisualInfo::PatternStatus::eOverlap, visualInfo);
        }
    }
    ATH_MSG_VERBOSE(__func__<<"() Patterns surviving overlap removal: "<< outputPatterns.size());
    return outputPatterns;
}

setLevel()

void AthMessaging::setLevel ( MSG::Level lvl )

inherited

Change the current logging level.

Use this rather than msg().setLevel() for proper operation with MT.

Definition at line 28 of file AthMessaging.cxx.

{
  m_lvl = lvl;
}

Member Data Documentation

ATLAS_THREAD_SAFE

std::atomic_flag m_initialized AthMessaging::ATLAS_THREAD_SAFE = ATOMIC_FLAG_INIT

mutableprivateinherited

Messaging initialized (initMessaging).

Definition at line 141 of file AthMessaging.h.

m_cfg

Config MuonR4::FastReco::GlobalPatternFinder::m_cfg

private

Global Pattern Recognition configuration.

Definition at line 394 of file GlobalPatternFinder.h.

m_imsg

std::atomic<IMessageSvc*> AthMessaging::m_imsg { nullptr }

mutableprivateinherited

MessageSvc pointer.

Definition at line 135 of file AthMessaging.h.

{ nullptr };

m_lvl

std::atomic<MSG::Level> AthMessaging::m_lvl { MSG::NIL }

mutableprivateinherited

Current logging level.

Definition at line 138 of file AthMessaging.h.

{ MSG::NIL };

m_msg_tls

boost::thread_specific_ptr<MsgStream> AthMessaging::m_msg_tls

mutableprivateinherited

MsgStream instance (a std::cout like with print-out levels).

Definition at line 132 of file AthMessaging.h.

m_nm

std::string AthMessaging::m_nm

privateinherited

Message source name.

Definition at line 129 of file AthMessaging.h.

m_spSorter

SpacePointPerLayerSorter MuonR4::FastReco::GlobalPatternFinder::m_spSorter {}

private

Spacepoint sorter per logical measurement layer.

Definition at line 392 of file GlobalPatternFinder.h.

{};

---

The documentation for this class was generated from the following files:

• GlobalPatternFinder.h
• GlobalPatternFinder.cxx