Algorithm to handle segment fits
More...

#include <SegmentFittingAlg.h>

Inheritance diagram for MuonR4::SegmentFittingAlg:

Collaboration diagram for MuonR4::SegmentFittingAlg:

Public Types
using	HitVec = SegmentFitResult::HitVec

Public Member Functions
virtual	~SegmentFittingAlg ()

virtual StatusCode	initialize () override

virtual StatusCode	execute (const EventContext &ctx) const override

	AthReentrantAlgorithm (const std::string &name, ISvcLocator *pSvcLocator)
	Constructor with parameters: More...

virtual StatusCode	sysInitialize () override
	Override sysInitialize. More...

virtual bool	isClonable () const override
	Specify if the algorithm is clonable. More...

virtual unsigned int	cardinality () const override
	Cardinality (Maximum number of clones that can exist) special value 0 means that algorithm is reentrant. More...

virtual StatusCode	sysExecute (const EventContext &ctx) override
	Execute an algorithm. More...

virtual const DataObjIDColl &	extraOutputDeps () const override
	Return the list of extra output dependencies. More...

virtual bool	filterPassed (const EventContext &ctx) const

virtual void	setFilterPassed (bool state, const EventContext &ctx) const

ServiceHandle< StoreGateSvc > &	evtStore ()
	The standard `StoreGateSvc` (event store) Returns (kind of) a pointer to the `StoreGateSvc`. More...

const ServiceHandle< StoreGateSvc > &	evtStore () const
	The standard `StoreGateSvc` (event store) Returns (kind of) a pointer to the `StoreGateSvc`. More...

const ServiceHandle< StoreGateSvc > &	detStore () const
	The standard `StoreGateSvc/DetectorStore` Returns (kind of) a pointer to the `StoreGateSvc`. More...

virtual StatusCode	sysStart () override
	Handle START transition. More...

virtual std::vector< Gaudi::DataHandle * >	inputHandles () const override
	Return this algorithm's input handles. More...

virtual std::vector< Gaudi::DataHandle * >	outputHandles () const override
	Return this algorithm's output handles. More...

Gaudi::Details::PropertyBase &	declareProperty (Gaudi::Property< T > &t)

Gaudi::Details::PropertyBase *	declareProperty (const std::string &name, SG::VarHandleKey &hndl, const std::string &doc, const SG::VarHandleKeyType &)
	Declare a new Gaudi property. More...

Gaudi::Details::PropertyBase *	declareProperty (const std::string &name, SG::VarHandleBase &hndl, const std::string &doc, const SG::VarHandleType &)
	Declare a new Gaudi property. More...

Gaudi::Details::PropertyBase *	declareProperty (const std::string &name, SG::VarHandleKeyArray &hndArr, const std::string &doc, const SG::VarHandleKeyArrayType &)

Gaudi::Details::PropertyBase *	declareProperty (const std::string &name, T &property, const std::string &doc, const SG::NotHandleType &)
	Declare a new Gaudi property. More...

Gaudi::Details::PropertyBase *	declareProperty (const std::string &name, T &property, const std::string &doc="none")
	Declare a new Gaudi property. More...

void	updateVHKA (Gaudi::Details::PropertyBase &)

MsgStream &	msg () const

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

bool	msgLvl (const MSG::Level lvl) const

Protected Member Functions
void	renounceArray (SG::VarHandleKeyArray &handlesArray)
	remove all handles from I/O resolution More...

std::enable_if_t< std::is_void_v< std::result_of_t< decltype(&T::renounce)(T)> > &&!std::is_base_of_v< SG::VarHandleKeyArray, T > &&std::is_base_of_v< Gaudi::DataHandle, T >, void >	renounce (T &h)

void	extraDeps_update_handler (Gaudi::Details::PropertyBase &ExtraDeps)
	Add StoreName to extra input/output deps as needed. More...

Private Types
using	Parameters = SegmentFit::Parameters

typedef ServiceHandle< StoreGateSvc >	StoreGateSvc_t

Private Member Functions
template<class ContainerType >
StatusCode	retrieveContainer (const EventContext &ctx, const SG::ReadHandleKey< ContainerType > &key, const ContainerType *&contToPush) const
	Helper method to fetch data from StoreGate. More...

SegmentFitResult	fitSegmentHits (const EventContext &ctx, const ActsGeometryContext &gctx, const Parameters &startPars, SegmentFitResult::HitVec &&calibHits) const
	Executes the segment fit with start parameters. More...

std::vector< std::unique_ptr< Segment > >	fitSegmentSeed (const EventContext &ctx, const ActsGeometryContext &gctx, const SegmentSeed *seed) const

bool	removeOutliers (const EventContext &ctx, const ActsGeometryContext &gctx, const SegmentSeed &seed, SegmentFitResult &fitResult) const
	Spot hits with large discrepancy from the estimated parameters and remove them from the list. More...

bool	plugHoles (const EventContext &ctx, const ActsGeometryContext &gctx, const SegmentSeed &seed, SegmentFitResult &toRecover) const
	Recovery of missed hits. More...

void	eraseWrongHits (const ActsGeometryContext &gctx, SegmentFitResult &candidate) const
	Removes all hits from the segment which are obvious outliers. More...

void	resolveAmbiguities (const ActsGeometryContext &gctx, std::vector< std::unique_ptr< Segment >> &segmentCandidates) const

Gaudi::Details::PropertyBase &	declareGaudiProperty (Gaudi::Property< T > &hndl, const SG::VarHandleKeyType &)
	specialization for handling Gaudi::Property<SG::VarHandleKey> More...

Gaudi::Details::PropertyBase &	declareGaudiProperty (Gaudi::Property< T > &hndl, const SG::VarHandleKeyArrayType &)
	specialization for handling Gaudi::Property<SG::VarHandleKeyArray> More...

Gaudi::Details::PropertyBase &	declareGaudiProperty (Gaudi::Property< T > &hndl, const SG::VarHandleType &)
	specialization for handling Gaudi::Property<SG::VarHandleBase> More...

Gaudi::Details::PropertyBase &	declareGaudiProperty (Gaudi::Property< T > &t, const SG::NotHandleType &)
	specialization for handling everything that's not a Gaudi::Property<SG::VarHandleKey> or a <SG::VarHandleKeyArray> More...

Static Private Member Functions
static std::unique_ptr< Segment >	convertToSegment (const Amg::Transform3D &locToGlobTrf, const SegmentSeed *parentSeed, SegmentFitResult &&toConvert)
	Converts the fit result into a segment object. More...

Private Attributes
SG::ReadHandleKey< SegmentSeedContainer >	m_seedKey {this, "ReadKey", "MuonHoughStationSegmentSeeds"}
	ReadHandle of the seeds. More...

SG::WriteHandleKey< SegmentContainer >	m_outSegments {this, "MuonSegmentContainer", "R4MuonSegments"}

SG::ReadHandleKey< ActsGeometryContext >	m_geoCtxKey {this, "AlignmentKey", "ActsAlignment", "cond handle key"}

ServiceHandle< Muon::IMuonIdHelperSvc >	m_idHelperSvc {this, "MuonIdHelperSvc", "Muon::MuonIdHelperSvc/MuonIdHelperSvc"}
	IdHelperSvc. More...

ToolHandle< ISpacePointCalibrator >	m_calibTool {this, "Calibrator", "" }
	Handle to the space point calibrator. More...

ToolHandle< MuonValR4::IPatternVisualizationTool >	m_visionTool {this, "VisualizationTool", ""}
	Pattern visualization tool. More...

Gaudi::Property< bool >	m_useMinuit {this, "useMinuit", false}
	Toggle the fitter. More...

Gaudi::Property< bool >	m_doT0Fit {this, "fitSegmentT0", true}

Gaudi::Property< bool >	m_doBeamspotConstraint {this, "doBeamspotConstraint", false}
	Add beamline constraint. More...

Gaudi::Property< double >	m_beamSpotR {this, "BeamSpotRadius", 30.* Gaudi::Units::cm}

Gaudi::Property< double >	m_beamSpotL {this, "BeamSpotLength", 20. * Gaudi::Units::m}

Gaudi::Property< double >	m_seedHitChi2 {this, "ResoSeedHitAssoc", 5. }
	Two mdt seeds are the same if their defining parameters match wihin. More...

Gaudi::Property< bool >	m_recalibSeed {this, "SeedRecalibrate", false}
	Toggle seed recalibration. More...

Gaudi::Property< bool >	m_refineSeed {this, "SeedRefine", false}
	Toggle seed refit. More...

Gaudi::Property< double >	m_outlierRemovalCut {this, "OutlierRemoval", 5.}
	Cut on the segment chi2 / nDoF to launch the outlier removal. More...

Gaudi::Property< double >	m_recoveryPull {this, "RecoveryPull", 5.}

std::unique_ptr< SegmentAmbiSolver >	m_ambiSolver {}

DataObjIDColl	m_extendedExtraObjects
	Extra output dependency collection, extended by AthAlgorithmDHUpdate to add symlinks. More...

StoreGateSvc_t	m_evtStore
	Pointer to StoreGate (event store by default) More...

StoreGateSvc_t	m_detStore
	Pointer to StoreGate (detector store by default) More...

std::vector< SG::VarHandleKeyArray * >	m_vhka

bool	m_varHandleArraysDeclared

Detailed Description

Algorithm to handle segment fits

This is currently a placeholder to test ideas!

Definition at line 31 of file SegmentFittingAlg.h.

Member Typedef Documentation

◆ HitVec

using MuonR4::SegmentFittingAlg::HitVec = SegmentFitResult::HitVec

Definition at line 38 of file SegmentFittingAlg.h.

◆ Parameters

using MuonR4::SegmentFittingAlg::Parameters = SegmentFit::Parameters

private

Definition at line 41 of file SegmentFittingAlg.h.

◆ StoreGateSvc_t

typedef ServiceHandle<StoreGateSvc> AthCommonDataStore< AthCommonMsg< Gaudi::Algorithm > >::StoreGateSvc_t

privateinherited

Definition at line 388 of file AthCommonDataStore.h.

Constructor & Destructor Documentation

◆ ~SegmentFittingAlg()

MuonR4::SegmentFittingAlg::~SegmentFittingAlg ( )

virtualdefault

Member Function Documentation

◆ AthReentrantAlgorithm()

AthReentrantAlgorithm::AthReentrantAlgorithm

Constructor with parameters:

Definition at line 90 of file AthReentrantAlgorithm.cxx.

                                                                          :
   ::AthCommonDataStore<AthCommonMsg<Gaudi::Algorithm>>   ( name, pSvcLocator )
 {
  
   // Set up to run AthAlgorithmDHUpdate in sysInitialize before
   // merging depedency lists.  This extends the output dependency
   // list with any symlinks implied by inheritance relations.
   m_updateDataHandles =
     std::make_unique<AthenaBaseComps::AthAlgorithmDHUpdate>
     (m_extendedExtraObjects,
      std::move (m_updateDataHandles));
 }

◆ cardinality()

unsigned int AthReentrantAlgorithm::cardinality ( ) const

overridevirtualinherited

Cardinality (Maximum number of clones that can exist) special value 0 means that algorithm is reentrant.

Override this to return 0 for reentrant algorithms.

Definition at line 55 of file AthReentrantAlgorithm.cxx.

 {
   return 0;
 }

◆ convertToSegment()

std::unique_ptr< Segment > MuonR4::SegmentFittingAlg::convertToSegment	(	const Amg::Transform3D &	locToGlobTrf,
		const SegmentSeed *	parentSeed,
		SegmentFitResult &&	toConvert
	)

staticprivate

Converts the fit result into a segment object.

Parameters

locToGlobTrf	Local to global transform to translate the segment parameters into global parameters
parentSeed	Segment seed from which the segment was built
toConvert	Fitted segment that needs conversion

Definition at line 325 of file SegmentFittingAlg.cxx.

                                                                                            {
         const auto [locPos, locDir] = data.makeLine();
         Amg::Vector3D globPos = locToGlob * locPos;
         Amg::Vector3D globDir = locToGlob.linear()* locDir;
  
  
         auto finalSeg = std::make_unique<Segment>(std::move(globPos), std::move(globDir),
                                                   patternSeed, std::move(data.calibMeasurements),
                                                   data.chi2, data.nDoF);
         finalSeg->setCallsToConverge(data.nIter);
         finalSeg->setChi2PerMeasurement(std::move(data.chi2PerMeasurement));
         finalSeg->setParUncertainties(std::move(data.segmentParErrs));
         if (data.timeFit) {
             finalSeg->setSegmentT0(data.segmentPars[toInt(ParamDefs::time)]);
             
         }
  
  
         return finalSeg;
     }

◆ declareGaudiProperty() [1/4]

Gaudi::Details::PropertyBase& AthCommonDataStore< AthCommonMsg< Gaudi::Algorithm > >::declareGaudiProperty	(	Gaudi::Property< T > &	hndl,
		const SG::VarHandleKeyArrayType &
	)

inlineprivateinherited

specialization for handling Gaudi::Property<SG::VarHandleKeyArray>

Definition at line 170 of file AthCommonDataStore.h.

   {
     return *AthCommonDataStore<PBASE>::declareProperty(hndl.name(),
                                                        hndl.value(), 
                                                        hndl.documentation());
  
   }

◆ declareGaudiProperty() [2/4]

Gaudi::Details::PropertyBase& AthCommonDataStore< AthCommonMsg< Gaudi::Algorithm > >::declareGaudiProperty	(	Gaudi::Property< T > &	hndl,
		const SG::VarHandleKeyType &
	)

inlineprivateinherited

specialization for handling Gaudi::Property<SG::VarHandleKey>

Definition at line 156 of file AthCommonDataStore.h.

   {
     return *AthCommonDataStore<PBASE>::declareProperty(hndl.name(),
                                                        hndl.value(), 
                                                        hndl.documentation());
  
   }

◆ declareGaudiProperty() [3/4]

Gaudi::Details::PropertyBase& AthCommonDataStore< AthCommonMsg< Gaudi::Algorithm > >::declareGaudiProperty	(	Gaudi::Property< T > &	hndl,
		const SG::VarHandleType &
	)

inlineprivateinherited

specialization for handling Gaudi::Property<SG::VarHandleBase>

Definition at line 184 of file AthCommonDataStore.h.

   {
     return *AthCommonDataStore<PBASE>::declareProperty(hndl.name(),
                                                        hndl.value(), 
                                                        hndl.documentation());
   }

◆ declareGaudiProperty() [4/4]

Gaudi::Details::PropertyBase& AthCommonDataStore< AthCommonMsg< Gaudi::Algorithm > >::declareGaudiProperty	(	Gaudi::Property< T > &	t,
		const SG::NotHandleType &
	)

inlineprivateinherited

specialization for handling everything that's not a Gaudi::Property<SG::VarHandleKey> or a <SG::VarHandleKeyArray>

Definition at line 199 of file AthCommonDataStore.h.

   {
     return PBASE::declareProperty(t);
   }

◆ declareProperty() [1/6]

Gaudi::Details::PropertyBase* AthCommonDataStore< AthCommonMsg< Gaudi::Algorithm > >::declareProperty	(	const std::string &	name,
		SG::VarHandleBase &	hndl,
		const std::string &	doc,
		const SG::VarHandleType &
	)

inlineinherited

Declare a new Gaudi property.

Parameters

name	Name of the property.
hndl	Object holding the property value.
doc	Documentation string for the property.

This is the version for types that derive from SG::VarHandleBase. The property value object is put on the input and output lists as appropriate; then we forward to the base class.

Definition at line 245 of file AthCommonDataStore.h.

   {
     this->declare(hndl.vhKey());
     hndl.vhKey().setOwner(this);
  
     return PBASE::declareProperty(name,hndl,doc);
   }

◆ declareProperty() [2/6]

Gaudi::Details::PropertyBase* AthCommonDataStore< AthCommonMsg< Gaudi::Algorithm > >::declareProperty	(	const std::string &	name,
		SG::VarHandleKey &	hndl,
		const std::string &	doc,
		const SG::VarHandleKeyType &
	)

inlineinherited

Declare a new Gaudi property.

Parameters

name	Name of the property.
hndl	Object holding the property value.
doc	Documentation string for the property.

This is the version for types that derive from SG::VarHandleKey. The property value object is put on the input and output lists as appropriate; then we forward to the base class.

Definition at line 221 of file AthCommonDataStore.h.

   {
     this->declare(hndl);
     hndl.setOwner(this);
  
     return PBASE::declareProperty(name,hndl,doc);
   }

◆ declareProperty() [3/6]

Gaudi::Details::PropertyBase* AthCommonDataStore< AthCommonMsg< Gaudi::Algorithm > >::declareProperty	(	const std::string &	name,
		SG::VarHandleKeyArray &	hndArr,
		const std::string &	doc,
		const SG::VarHandleKeyArrayType &
	)

inlineinherited

Definition at line 259 of file AthCommonDataStore.h.

   {
  
     // std::ostringstream ost;
     // ost << Algorithm::name() << " VHKA declareProp: " << name 
     //     << " size: " << hndArr.keys().size() 
     //     << " mode: " << hndArr.mode() 
     //     << "  vhka size: " << m_vhka.size()
     //     << "\n";
     // debug() << ost.str() << endmsg;
  
     hndArr.setOwner(this);
     m_vhka.push_back(&hndArr);
  
     Gaudi::Details::PropertyBase* p =  PBASE::declareProperty(name, hndArr, doc);
     if (p != 0) {
       p->declareUpdateHandler(&AthCommonDataStore<PBASE>::updateVHKA, this);
     } else {
       ATH_MSG_ERROR("unable to call declareProperty on VarHandleKeyArray " 
                     << name);
     }
  
     return p;
  
   }

◆ declareProperty() [4/6]

Gaudi::Details::PropertyBase* AthCommonDataStore< AthCommonMsg< Gaudi::Algorithm > >::declareProperty	(	const std::string &	name,
		T &	property,
		const std::string &	doc,
		const SG::NotHandleType &
	)

inlineinherited

Declare a new Gaudi property.

Parameters

name	Name of the property.
property	Object holding the property value.
doc	Documentation string for the property.

This is the generic version, for types that do not derive from SG::VarHandleKey. It just forwards to the base class version of declareProperty.

Definition at line 333 of file AthCommonDataStore.h.

   {
     return PBASE::declareProperty(name, property, doc);
   }

◆ declareProperty() [5/6]

Gaudi::Details::PropertyBase* AthCommonDataStore< AthCommonMsg< Gaudi::Algorithm > >::declareProperty	(	const std::string &	name,
		T &	property,
		const std::string &	doc = `"none"`
	)

inlineinherited

Declare a new Gaudi property.

Parameters

name	Name of the property.
property	Object holding the property value.
doc	Documentation string for the property.

This dispatches to either the generic declareProperty or the one for VarHandle/Key/KeyArray.

Definition at line 352 of file AthCommonDataStore.h.

   {
     typedef typename SG::HandleClassifier<T>::type htype;
     return declareProperty (name, property, doc, htype());
   }

◆ declareProperty() [6/6]

Gaudi::Details::PropertyBase& AthCommonDataStore< AthCommonMsg< Gaudi::Algorithm > >::declareProperty ( Gaudi::Property< T > & t )

inlineinherited

Definition at line 145 of file AthCommonDataStore.h.

                                                                  {
     typedef typename SG::HandleClassifier<T>::type htype;
     return AthCommonDataStore<PBASE>::declareGaudiProperty(t, htype());
   }

◆ detStore()

const ServiceHandle<StoreGateSvc>& AthCommonDataStore< AthCommonMsg< Gaudi::Algorithm > >::detStore ( ) const

inlineinherited

The standard StoreGateSvc/DetectorStore Returns (kind of) a pointer to the StoreGateSvc.

Definition at line 95 of file AthCommonDataStore.h.

95 { return m_detStore; }

◆ eraseWrongHits()

void MuonR4::SegmentFittingAlg::eraseWrongHits	(	const ActsGeometryContext &	gctx,
		SegmentFitResult &	candidate
	)		const

private

Removes all hits from the segment which are obvious outliers.

E.g. tubes which cannot be crossed by the segment.

Parameters

gctx	Geometry context needed to calculate the nominal time of arrival, if the time fit is activated
candidate	Reference of the segment candidate to prune.

The segment has never crossed the tube

Definition at line 529 of file SegmentFittingAlg.cxx.

                                                                                                              {
         auto [segPos, segDir] = makeLine(candidate.segmentPars);
         candidate.calibMeasurements.erase(std::remove_if(candidate.calibMeasurements.begin(), candidate.calibMeasurements.end(),
                                                 [&segPos, &segDir](const HitVec::value_type& hit){
                                                 if (hit->fitState() != CalibratedSpacePoint::State::Outlier) {
                                                     return false;
                                                 }
                                                 if (hit->type() == xAOD::UncalibMeasType::MdtDriftCircleType) {
                                                     const double dist = Amg::lineDistance(segPos, segDir, hit->positionInChamber(), hit->directionInChamber());
                                                     const auto* dc = static_cast<const xAOD::MdtDriftCircle*>(hit->spacePoint()->primaryMeasurement());
                                                     return dist >= dc->readoutElement()->innerTubeRadius();
                                                 }
                                                 return true;
                                                 }), candidate.calibMeasurements.end());
         const MuonGMR4::SpectrometerSector* chamber{nullptr};
         for (const auto& hit : candidate.calibMeasurements) {
             if (hit->type() != xAOD::UncalibMeasType::Other){
                 chamber = hit->spacePoint()->msSector();
                 break;
             }
         }
         std::optional<double> timeOfFlight = candidate.timeFit ? std::make_optional<double>((chamber->localToGlobalTrans(gctx)*segPos).mag() * inv_c) 
                                                                 : std::nullopt;
         auto [updatedMeasChi2, updateChi2] = SegmentFitHelpers::postFitChi2PerMas(candidate.segmentPars, timeOfFlight, 
                                                                                   candidate.calibMeasurements, msgStream());
         ATH_MSG_VERBOSE("The measurements before "<<candidate.chi2PerMeasurement<<", after: "<<updatedMeasChi2<<", chi2: "<<updateChi2);
         candidate.chi2PerMeasurement = std::move(updatedMeasChi2);
     }

◆ evtStore() [1/2]

ServiceHandle<StoreGateSvc>& AthCommonDataStore< AthCommonMsg< Gaudi::Algorithm > >::evtStore ( )

inlineinherited

The standard StoreGateSvc (event store) Returns (kind of) a pointer to the StoreGateSvc.

Definition at line 85 of file AthCommonDataStore.h.

85 { return m_evtStore; }

◆ evtStore() [2/2]

const ServiceHandle<StoreGateSvc>& AthCommonDataStore< AthCommonMsg< Gaudi::Algorithm > >::evtStore ( ) const

inlineinherited

The standard StoreGateSvc (event store) Returns (kind of) a pointer to the StoreGateSvc.

Definition at line 90 of file AthCommonDataStore.h.

90 { return m_evtStore; }

◆ execute()

StatusCode MuonR4::SegmentFittingAlg::execute ( const EventContext & ctx ) const

overridevirtual

Definition at line 74 of file SegmentFittingAlg.cxx.

                                                                        {
     
         const ActsGeometryContext* gctx{nullptr};
         ATH_CHECK(retrieveContainer(ctx, m_geoCtxKey, gctx));
         const SegmentSeedContainer* segmentSeeds=nullptr; 
         ATH_CHECK(retrieveContainer(ctx, m_seedKey, segmentSeeds));
     
         SG::WriteHandle writeSegments{m_outSegments, ctx};
         ATH_CHECK(writeSegments.record(std::make_unique<SegmentContainer>()));
         std::vector<std::unique_ptr<Segment>> allSegments{};
         for (const SegmentSeed* seed : *segmentSeeds) {
             std::vector<std::unique_ptr<Segment>> segments = fitSegmentSeed(ctx, *gctx, seed);
             if (m_visionTool.isEnabled() && segments.size() > 1) {
                 auto drawFinalReco = [this, &segments, &gctx, &ctx,&seed](const std::string& nameTag) {
                     PrimitiveVec segmentLines{};
                     double yLegend{0.85};
                     segmentLines.push_back(drawLabel(std::format("# segments: {:d}",  segments.size()), 0.12, yLegend, 14));
                     yLegend-=0.04;
                     for (const std::unique_ptr<Segment>& seg : segments) {
                         const Parameters pars = localSegmentPars(*gctx, *seg);
                         const auto [locPos, locDir] = makeLine(pars);
  
                         segmentLines.emplace_back(drawLine(pars, -Gaudi::Units::m, Gaudi::Units::m, kRed));
                         std::stringstream signStream{};
                         signStream<<"#chi^{2}/nDoF: "<<std::format("{:.2f}", seg->chi2() / seg->nDoF())<<"("<<seg->nDoF()<<"), ";
                         signStream<<"y_{0}="<<std::format("{:.2f}",pars[toInt(ParamDefs::y0)])<<", ";
                         signStream<<std::format("#theta={:.2f}", pars[toInt(ParamDefs::theta)]/ Gaudi::Units::deg )<<", ";
                         for (const Segment::MeasType& m : seg->measurements()) {
                             if (m->type() == xAOD::UncalibMeasType::MdtDriftCircleType && m->fitState() == CalibratedSpacePoint::State::Valid) {
                                 signStream<<(SegmentFitHelpers::driftSign(locPos, locDir, *m, msgStream()) == -1 ? "L" : "R");
                             }
                         }
                         segmentLines.push_back(drawLabel(signStream.str(), 0.12, yLegend, 12));
                         yLegend-=0.03;
                     }
  
                     m_visionTool->visualizeBucket(ctx, *seed->parentBucket(), nameTag, std::move(segmentLines));
                 };
                 drawFinalReco("all segments");
                 const unsigned int nBeforeAmbi = segments.size();
                 segments = m_ambiSolver->resolveAmbiguity(*gctx, std::move(segments));
                 if (nBeforeAmbi != segments.size()) {
                     drawFinalReco("post ambiguity");
                 }
             }
             allSegments.insert(allSegments.end(), std::make_move_iterator(segments.begin()),
                                                   std::make_move_iterator(segments.end()));
         }
         resolveAmbiguities(*gctx, allSegments);
         writeSegments->insert(writeSegments->end(),
                                   std::make_move_iterator(allSegments.begin()),
                                   std::make_move_iterator(allSegments.end()));
         ATH_MSG_VERBOSE("Found in total "<<writeSegments->size()<<" segments. ");
         return StatusCode::SUCCESS; 
     }

◆ extraDeps_update_handler()

void AthCommonDataStore< AthCommonMsg< Gaudi::Algorithm > >::extraDeps_update_handler ( Gaudi::Details::PropertyBase & ExtraDeps )

protectedinherited

Add StoreName to extra input/output deps as needed.

use the logic of the VarHandleKey to parse the DataObjID keys supplied via the ExtraInputs and ExtraOuputs Properties to add the StoreName if it's not explicitly given

◆ extraOutputDeps()

const DataObjIDColl & AthReentrantAlgorithm::extraOutputDeps ( ) const

overridevirtualinherited

Return the list of extra output dependencies.

This list is extended to include symlinks implied by inheritance relations.

Definition at line 79 of file AthReentrantAlgorithm.cxx.

 {
   // If we didn't find any symlinks to add, just return the collection
   // from the base class.  Otherwise, return the extended collection.
   if (!m_extendedExtraObjects.empty()) {
     return m_extendedExtraObjects;
   }
   return Algorithm::extraOutputDeps();
 }

◆ filterPassed()

virtual bool AthReentrantAlgorithm::filterPassed ( const EventContext & ctx ) const

inlinevirtualinherited

Definition at line 135 of file AthReentrantAlgorithm.h.

                                                            {
     return execState( ctx ).filterPassed();
   }

◆ fitSegmentHits()

SegmentFitResult MuonR4::SegmentFittingAlg::fitSegmentHits	(	const EventContext &	ctx,
		const ActsGeometryContext &	gctx,
		const Parameters &	startPars,
		SegmentFitResult::HitVec &&	calibHits
	)		const

private

Executes the segment fit with start parameters.

The returned fit result indicates whether the fit was a success and all relevant output parameters

ctx: Event context needed to access the calibration constants of the hits

gctx: Geometry context needed to place the segment globally within ATLAS to calculate the nominal time of arrival

startPars: Segment parameters preestimated from the SegmentSeed (either hough pattern or two drift circle seed)

calibHits: Vector of strip & mdt hits to consider for the fit

Add beamspot if neccessary

placeholder for a very generous beam spot: 300mm in X,Y (tracking volume), 20000 along Z

Configure the minimizer

Assumption that the particle travels at the speed of light

Execute fit

Definition at line 145 of file SegmentFittingAlg.cxx.

                                                                                                  {
         SegmentFitResult data{};
         if (calibHits.empty()) {
             ATH_MSG_WARNING("No hits, no segment");
             return data;
         }
  
         if (m_doBeamspotConstraint) {
             unsigned int numPhi = std::accumulate(calibHits.begin(), calibHits.end(),0,
                                                 [](unsigned int n, const HitVec::value_type& hit){
                                                     return n + (hit->fitState() == CalibratedSpacePoint::State::Valid &&
                                                                 hit->measuresPhi());
                                                 });
             if (numPhi) {
                 const Amg::Transform3D globToLoc{calibHits[0]->spacePoint()->msSector()->globalToLocalTrans(gctx)};
                 Amg::Vector3D beamSpot{globToLoc.translation()};
                 AmgSymMatrix(3) covariance{AmgSymMatrix(3)::Identity()}; 
                 covariance(0,0) = std::pow(m_beamSpotR, 2);
                 covariance(1,1) = std::pow(m_beamSpotR, 2);
                 covariance(2,2) = std::pow(m_beamSpotL, 2);
                 AmgSymMatrix(3) jacobian =  globToLoc.linear();
                 covariance = jacobian * covariance * jacobian.transpose(); 
                 AmgSymMatrix(2) beamSpotCov{covariance.block<2,2>(0,0)};
                 auto beamSpotSP = std::make_unique<CalibratedSpacePoint>(nullptr, std::move(beamSpot), Amg::Vector3D::Zero());
                 beamSpotSP->setCovariance<2>(std::move(beamSpotCov));
                 ATH_MSG_VERBOSE("Beam spot constraint "<<Amg::toString(beamSpotSP->positionInChamber())<<", "<<toString(beamSpotSP->covariance()));
                 calibHits.push_back(std::move(beamSpotSP));
             }
         }
  
         const Amg::Transform3D& locToGlob{calibHits[0]->spacePoint()->msSector()->localToGlobalTrans(gctx)};
  
         if (!m_useMinuit) {
             MdtSegmentFitter::Config fitCfg{};
             fitCfg.calibrator = m_calibTool.get();
             fitCfg.doTimeFit = m_doT0Fit;
  
             MdtSegmentFitter fitter{name(), std::move(fitCfg)};
             return fitter.fitSegment(ctx, std::move(calibHits), startPars, locToGlob);
         }
  
         data.segmentPars = startPars;
  
         CalibSegmentChi2Minimizer c2f{name(), ctx, locToGlob, copy(calibHits), m_calibTool.get(), m_doT0Fit};
         data.hasPhi = c2f.hasPhiMeas();
         data.timeFit = c2f.doTimeFit();
         data.nDoF = c2f.nDoF();
         if (data.nDoF <= 0) {
             ATH_MSG_DEBUG("Reject fit due to 0 degrees of freedom");
             return data;
         }
         ROOT::Minuit2::Minuit2Minimizer minimizer((name() + std::to_string(ctx.eventID().event_number())).c_str());
         minimizer.SetMaxFunctionCalls(100000);
         minimizer.SetTolerance(0.0001);
         minimizer.SetPrintLevel(-1);
         minimizer.SetStrategy(1);
  
         minimizer.SetVariable(toInt(ParamDefs::y0), "y0", startPars[toInt(ParamDefs::y0)], 1.e-5);
         minimizer.SetVariable(toInt(ParamDefs::theta), "theta", startPars[toInt(ParamDefs::theta)], 1.e-5);
         minimizer.SetVariableLimits(toInt(ParamDefs::y0), 
                                     startPars[toInt(ParamDefs::y0)] - 60. *Gaudi::Units::cm, 
                                     startPars[toInt(ParamDefs::y0)] + 60. *Gaudi::Units::cm);
         minimizer.SetVariableLimits(toInt(ParamDefs::theta),
                                     startPars[toInt(ParamDefs::theta)] - 0.6, 
                                     startPars[toInt(ParamDefs::theta)] + 0.6);
         
         if (data.hasPhi) {
             minimizer.SetVariable(toInt(ParamDefs::x0), "x0", startPars[toInt(ParamDefs::x0)], 1.e-5);
             minimizer.SetVariable(toInt(ParamDefs::phi), "phi", startPars[toInt(ParamDefs::phi)], 1.e-5);
             minimizer.SetVariableLimits(toInt(ParamDefs::x0), 
                                         startPars[toInt(ParamDefs::x0)] - 600, 
                                         startPars[toInt(ParamDefs::x0)] + 600);
             minimizer.SetVariableLimits(toInt(ParamDefs::phi), 
                                         startPars[toInt(ParamDefs::phi)] - 0.6, 
                                         startPars[toInt(ParamDefs::phi)] + 0.6);
         } else {
             minimizer.SetFixedVariable(toInt(ParamDefs::x0), "x0", 0.);
             minimizer.SetFixedVariable(toInt(ParamDefs::phi), "phi", 90.*Gaudi::Units::deg);
         }
         if (data.timeFit) {
             minimizer.SetVariable(toInt(ParamDefs::time), "t0", startPars[toInt(ParamDefs::time)] , 1.);
             minimizer.SetVariableLimits(toInt(ParamDefs::time), 
                                         startPars[toInt(ParamDefs::time)] - 50,
                                         startPars[toInt(ParamDefs::time)] + 50);
         } else{
             minimizer.SetFixedVariable(toInt(ParamDefs::time), "t0", 0.);
         }
         minimizer.SetFunction(c2f);
         if (!minimizer.Minimize() || !minimizer.Hesse()) {
             data.calibMeasurements = std::move(calibHits);
         } else {
             const double* xs = minimizer.X();
             const double* errs = minimizer.Errors();
  
             for (unsigned int p = 0; p < toInt(ParamDefs::nPars); ++p) {
                 data.segmentPars[p] = xs[p];
                 data.segmentParErrs(p,p) = errs[p];
             }
             std::optional<double> timeOfArrival{std::nullopt};
             const auto [locPos, locDir] = data.makeLine();
             if (data.timeFit) {
                 timeOfArrival = std::make_optional<double>((locToGlob*locPos).mag() * inv_c);
             }
             data.nIter = minimizer.NCalls();
             data.calibMeasurements = c2f.release(xs);
  
             auto [chiPerMeas, finalChi2] = SegmentFitHelpers::postFitChi2PerMas(data.segmentPars, timeOfArrival, 
                                                                                 data.calibMeasurements, msgStream());
             data.chi2PerMeasurement = std::move(chiPerMeas);
             data.chi2 = finalChi2;
             data.converged = true;
         }
         return data;
     }

◆ fitSegmentSeed()

std::vector< std::unique_ptr< Segment > > MuonR4::SegmentFittingAlg::fitSegmentSeed	(	const EventContext &	ctx,
		const ActsGeometryContext &	gctx,
		const SegmentSeed *	seed
	)		const

private

At very high inclanation angles, the muon may traverse 3 hits in the same layer (E.g. BEE)

Draw the pattern with all possible seeds

Draw the prefit

Definition at line 268 of file SegmentFittingAlg.cxx.

                                                                                  {
  
         const Amg::Transform3D& locToGlob{patternSeed->msSector()->localToGlobalTrans(gctx)};
         std::vector<std::unique_ptr<Segment>> segments{};
  
         MdtSegmentSeedGenerator::Config genCfg{};
         genCfg.hitPullCut = m_seedHitChi2;
         genCfg.recalibSeedCircles = m_recalibSeed;
         genCfg.fastSeedFit = m_refineSeed;
         genCfg.calibrator = m_calibTool.get();
         
  
         genCfg.busyLayerLimit = 2 + 2*(patternSeed->parameters()[toInt(ParamDefs::theta)] > 50 * Gaudi::Units::deg);
         if (m_visionTool.isEnabled()) { 
             PrimitiveVec seedLines{};
             MdtSegmentSeedGenerator drawMe{name(), patternSeed, genCfg};
             while(auto s = drawMe.nextSeed(ctx)) {
                 seedLines.push_back(drawLine(s->parameters, -Gaudi::Units::m, Gaudi::Units::m, kViolet));
             }
             seedLines.push_back(drawLabel(std::format("possible seeds: {:d}",  drawMe.numGenerated()), 0.15, 0.85, 14));
             m_visionTool->visualizeSeed(ctx, *patternSeed, "pattern", std::move(seedLines));
         }
  
         MdtSegmentSeedGenerator seedGen{name(), patternSeed, std::move(genCfg)};
         while (auto seed = seedGen.nextSeed(ctx)) {
             SegmentFitResult data{};
             data.segmentPars = seed->parameters;
             data.calibMeasurements = std::move(seed->measurements);            
             if (m_visionTool.isEnabled()) {
                 auto seedCopy = convertToSegment(locToGlob, patternSeed, copy(data));
                 m_visionTool->visualizeSegment(ctx, *seedCopy, std::format("Pre fit {:d}", seedGen.numGenerated()));
             } 
             data = fitSegmentHits(ctx, gctx, seed->parameters, std::move(data.calibMeasurements));
             data.nIter +=  seed->nIter;
             if (m_visionTool.isEnabled() && data.converged) {
                 auto seedCopy = convertToSegment(locToGlob, patternSeed, copy(data));
                 m_visionTool->visualizeSegment(ctx, *seedCopy, std::format("Intermediate fit {:d}", seedGen.numGenerated()));
             }
             if (!removeOutliers(ctx, gctx, *patternSeed, data)) {
                 continue;
             }          
             if (!plugHoles(ctx, gctx, *patternSeed, data)) {
                 continue;
             }
             if (m_visionTool.isEnabled()) {
                 auto seedCopy = convertToSegment(locToGlob, patternSeed, copy(data));
                 m_visionTool->visualizeSegment(ctx, *seedCopy, std::format("Final fit {:d}", seedGen.numGenerated()));
             }
             segments.push_back(convertToSegment(locToGlob, patternSeed, std::move(data)));
         }
         return segments;
     }

◆ initialize()

StatusCode MuonR4::SegmentFittingAlg::initialize ( )

overridevirtual

Definition at line 63 of file SegmentFittingAlg.cxx.

                                              {
         ATH_CHECK(m_geoCtxKey.initialize());
         ATH_CHECK(m_seedKey.initialize());
         ATH_CHECK(m_outSegments.initialize());  
         ATH_CHECK(m_calibTool.retrieve());
         ATH_CHECK(m_idHelperSvc.retrieve());
         ATH_CHECK(m_visionTool.retrieve(EnableTool{!m_visionTool.empty()}));
         SegmentAmbiSolver::Config cfg{};
         m_ambiSolver = std::make_unique<SegmentAmbiSolver>(name(), std::move(cfg));
         return StatusCode::SUCCESS;
     }

◆ inputHandles()

virtual std::vector<Gaudi::DataHandle*> AthCommonDataStore< AthCommonMsg< Gaudi::Algorithm > >::inputHandles ( ) const

overridevirtualinherited

Return this algorithm's input handles.

We override this to include handle instances from key arrays if they have not yet been declared. See comments on updateVHKA.

◆ isClonable()

bool AthReentrantAlgorithm::isClonable ( ) const

overridevirtualinherited

Specify if the algorithm is clonable.

Reentrant algorithms are clonable.

Reimplemented in Simulation::BeamEffectsAlg, InDet::SiTrackerSpacePointFinder, InDet::SCT_Clusterization, InDet::SiSPSeededTrackFinder, SCTRawDataProvider, InDet::GNNSeedingTrackMaker, SCT_PrepDataToxAOD, RoIBResultToxAOD, SCT_CablingCondAlgFromCoraCool, SCT_ReadCalibDataTestAlg, SCT_CablingCondAlgFromText, InDet::SiSPGNNTrackMaker, SCT_ReadCalibChipDataTestAlg, SCT_TestCablingAlg, SCT_ConfigurationConditionsTestAlg, ITkPixelCablingAlg, ITkStripCablingAlg, SCTEventFlagWriter, SCT_ConditionsSummaryTestAlg, SCT_ModuleVetoTestAlg, SCT_MonitorConditionsTestAlg, SCT_LinkMaskingTestAlg, SCT_MajorityConditionsTestAlg, SCT_RODVetoTestAlg, SCT_SensorsTestAlg, SCT_TdaqEnabledTestAlg, SCT_SiliconConditionsTestAlg, SCTSiLorentzAngleTestAlg, SCT_ByteStreamErrorsTestAlg, SCT_ConditionsParameterTestAlg, SCT_FlaggedConditionTestAlg, SCT_StripVetoTestAlg, SCT_RawDataToxAOD, and SCTSiPropertiesTestAlg.

Definition at line 44 of file AthReentrantAlgorithm.cxx.

 {
   // Reentrant algorithms are clonable.
   return true;
 }

◆ msg() [1/2]

MsgStream& AthCommonMsg< Gaudi::Algorithm >::msg ( ) const

inlineinherited

Definition at line 24 of file AthCommonMsg.h.

                                 {
     return this->msgStream();
   }

◆ msg() [2/2]

MsgStream& AthCommonMsg< Gaudi::Algorithm >::msg ( const MSG::Level lvl ) const

inlineinherited

Definition at line 27 of file AthCommonMsg.h.

                                                   {
     return this->msgStream(lvl);
   }

◆ msgLvl()

bool AthCommonMsg< Gaudi::Algorithm >::msgLvl ( const MSG::Level lvl ) const

inlineinherited

Definition at line 30 of file AthCommonMsg.h.

                                                {
     return this->msgLevel(lvl);
   }

◆ outputHandles()

virtual std::vector<Gaudi::DataHandle*> AthCommonDataStore< AthCommonMsg< Gaudi::Algorithm > >::outputHandles ( ) const

overridevirtualinherited

Return this algorithm's output handles.

We override this to include handle instances from key arrays if they have not yet been declared. See comments on updateVHKA.

◆ plugHoles()

bool MuonR4::SegmentFittingAlg::plugHoles	(	const EventContext &	ctx,
		const ActsGeometryContext &	gctx,
		const SegmentSeed &	seed,
		SegmentFitResult &	toRecover
	)		const

private

Recovery of missed hits.

Hits in the space point bucket that are maximally <RecoveryPull> away from the fitted segment are put onto the segment candidate and the candidate is refitted. If the refitted candidate has a chi2/nDoF < <OutlierRemoval> the canidate is automatically choosen otherwise, its chi needs to be better.

Parameters

ctx	EventContext needed for calibration of the hits
gctx	Geometry context needed if the beamspot constaint will be packed onto the segment
seed	Segment seed from the pattern recognition to access the underlying bucket
toRecover	Fit result with parameters & hits to recover.

We've the first estimator of the segment fit

Setup a map to replace space points if they better suite

Hit is already used in the segment fit

Remove the beamspot constraint measurement

Nothing has been recovered. Just bail out

If the chi2 is less than 5, no outlier rejection is launched. So also accept any recovered segment below that threshold

Next check whether the recovery made measurements marked as outlier feasable for the hole recovery

Finally remove all hits from the calib measurements which are obvious outliers

Definition at line 409 of file SegmentFittingAlg.cxx.

                                                                            {
         ATH_MSG_VERBOSE("plugHoles() -- segment "<<toString(beforeRecov.segmentPars)
                         <<", chi2: "<<beforeRecov.chi2<<", nDoF: "<<beforeRecov.nDoF<<", "
                         <<beforeRecov.chi2 /std::max(beforeRecov.nDoF, 1) );
         std::unordered_set<const SpacePoint*> usedSpacePoint{};
         for (const HitVec::value_type& hit : beforeRecov.calibMeasurements) {
             usedSpacePoint.insert(hit->spacePoint());
         }
  
         HitVec candidateHits{};
         SpacePointPerLayerSorter hitLayers{*seed.parentBucket()};
         bool hasCandidate{false};
         const auto [locPos, locDir] = beforeRecov.makeLine();
         for (const std::vector<HoughHitType>& mdtLayer : hitLayers.mdtHits()) {
             for (const SpacePoint* mdtHit: mdtLayer) {
                 if (usedSpacePoint.count(mdtHit)) {
                     continue;
                 }
                 const double dist = Amg::lineDistance(locPos, locDir, mdtHit->positionInChamber(), mdtHit->directionInChamber());
                 const auto* dc = static_cast<const xAOD::MdtDriftCircle*>(mdtHit->primaryMeasurement());
                 if (dist >= dc->readoutElement()->innerTubeRadius()) {
                     continue;
                 }
                 HitVec::value_type calibHit{m_calibTool->calibrate(ctx, mdtHit, locPos, locDir, beforeRecov.segmentPars[toInt(ParamDefs::time)])};
                 const double pull = std::sqrt(SegmentFitHelpers::chiSqTermMdt(locPos, locDir, *calibHit, msgStream()));
                 ATH_MSG_VERBOSE(__func__<<"() :"<<__LINE__<<" Candidate hit for recovery "<<m_idHelperSvc->toString(mdtHit->identify())<<", chi2: "<<pull);
                 if (pull <= m_recoveryPull) {
                     hasCandidate |= calibHit->fitState() == CalibratedSpacePoint::State::Valid;
                     candidateHits.push_back(std::move(calibHit));
                 } else {
                     calibHit->setFitState(CalibratedSpacePoint::State::Outlier);
                     candidateHits.push_back(std::move(calibHit));
                 }
             }
         }
         if (!hasCandidate) {
             ATH_MSG_VERBOSE("No space point candidates for recovery were found");
             beforeRecov.calibMeasurements.insert(beforeRecov.calibMeasurements.end(), 
                                                  std::make_move_iterator(candidateHits.begin()),
                                                  std::make_move_iterator(candidateHits.end()));
             eraseWrongHits(gctx, beforeRecov);
             return beforeRecov.nDoF > 0;
         }
    
         HitVec copied = copy(beforeRecov.calibMeasurements), copiedCandidates = copy(candidateHits);
         if (m_doBeamspotConstraint) {
             removeBeamSpot(copied);
         }
  
         candidateHits.insert(candidateHits.end(), std::make_move_iterator(copied.begin()), std::make_move_iterator(copied.end()));
  
         SegmentFitResult recovered = fitSegmentHits(ctx, gctx, beforeRecov.segmentPars, std::move(candidateHits));
         if (!recovered.converged) {
             return false;
         }
         if (recovered.nDoF + recovered.timeFit <= beforeRecov.nDoF + beforeRecov.timeFit) {
             for (HitVec::value_type& hit : copiedCandidates) {
                 hit->setFitState(CalibratedSpacePoint::State::Outlier);
                 beforeRecov.calibMeasurements.push_back(std::move(hit));
             }
             eraseWrongHits(gctx, beforeRecov);
             return true;
         }
         ATH_MSG_VERBOSE("Chi2, nDOF before:"<<beforeRecov.chi2<<", "<<beforeRecov.nDoF<<" after recovery: "<<recovered.chi2<<", "<<recovered.nDoF);
         double redChi2 = recovered.chi2 / std::max(recovered.nDoF,1);
         if (redChi2 < m_outlierRemovalCut || (beforeRecov.nDoF == 0) || redChi2 < beforeRecov.chi2 / beforeRecov.nDoF) {
             ATH_MSG_VERBOSE("Accept segment with recovered "<<(recovered.nDoF + recovered.timeFit) - (beforeRecov.nDoF + beforeRecov.timeFit)<<" hits.");
             recovered.nIter += beforeRecov.nIter;
             beforeRecov = std::move(recovered);
             while (true) {
                 bool runAnotherTrial = false;
                 copied = copy(beforeRecov.calibMeasurements);
                 if (m_doBeamspotConstraint) {
                     removeBeamSpot(copied);
                 }
                 for (unsigned int m = 0; m < beforeRecov.calibMeasurements.size(); ++m) {
                     if (beforeRecov.calibMeasurements[m]->fitState() == CalibratedSpacePoint::State::Outlier && 
                         std::sqrt(beforeRecov.chi2PerMeasurement[m]) < m_recoveryPull) {
                         copied[m]->setFitState(CalibratedSpacePoint::State::Valid);
                         runAnotherTrial = true;
                     }
                 }
                 if (!runAnotherTrial) {
                     break;
                 }
                 recovered = fitSegmentHits(ctx, gctx, beforeRecov.segmentPars, std::move(copied));
                 if (!recovered.converged) {
                     break;
                 }
                 if (recovered.nDoF + recovered.timeFit <= beforeRecov.nDoF + beforeRecov.timeFit) {
                     break;
                 }
                 redChi2 = recovered.chi2 / std::max(recovered.nDoF, 1);
                 if (redChi2 < m_outlierRemovalCut || redChi2 < beforeRecov.chi2 / beforeRecov.nDoF) {
                     recovered.nIter += beforeRecov.nIter;
                     beforeRecov = std::move(recovered);
                 } else {
                     break;
                 }
             }
             eraseWrongHits(gctx, beforeRecov);
         } else{
             for (HitVec::value_type& hit : copiedCandidates) {
                 hit->setFitState(CalibratedSpacePoint::State::Outlier);
                 beforeRecov.calibMeasurements.push_back(std::move(hit));
             }
         }
         return true;
     }

◆ removeOutliers()

bool MuonR4::SegmentFittingAlg::removeOutliers	(	const EventContext &	ctx,
		const ActsGeometryContext &	gctx,
		const SegmentSeed &	seed,
		SegmentFitResult &	fitResult
	)		const

private

Spot hits with large discrepancy from the estimated parameters and remove them from the list.

Parameters

ctx	EventContext needed for hit calibration
gctx	Geometry context needed if the beamspot constaint will be packed onto the segment
fitResult	Data structure carrying the fit parameters & calibrated hits

If no degree of freedom is in the segment fit then try to plug the holes

Remove a priori the beamspot constaint as it never should pose any problem and another one will be added anyway in the next iteration

Next sort the measurements by chi2

Declare the hit with the largest chi2 as outlier.

Remove the last measurement as it has the largest discrepancy

Definition at line 348 of file SegmentFittingAlg.cxx.

                                                                          {
         
         if (data.nDoF<=0 || data.calibMeasurements.empty()) {
             ATH_MSG_VERBOSE("No degree of freedom available. What shall be removed?!");
             return false;
         }
  
         const auto [segPos, segDir] = data.makeLine();
  
         if (data.converged && data.chi2 / data.nDoF < m_outlierRemovalCut) {
             ATH_MSG_VERBOSE("The segment "<<Amg::toString(segPos)<<" + "<<Amg::toString(segDir)
                             <<" is already of good quality "<<data.chi2 / std::max(data.nDoF, 1)
                             <<". Don't remove outliers");
             return true;
         }
         ATH_MSG_VERBOSE("Segment "<<toString(data.segmentPars)<<" is of badish quality.");
         if (m_doBeamspotConstraint && removeBeamSpot(data.calibMeasurements)) {
             data.nDoF-=2;
             data.nPhiMeas-=1;
         }
  
         std::sort(data.calibMeasurements.begin(), data.calibMeasurements.end(),
                   [&, this](const HitVec::value_type& a, const HitVec::value_type& b){
                     return SegmentFitHelpers::chiSqTerm(segPos, segDir, data.segmentPars[toInt(ParamDefs::time)], std::nullopt, *a, msgStream()) <
                            SegmentFitHelpers::chiSqTerm(segPos, segDir, data.segmentPars[toInt(ParamDefs::time)], std::nullopt, *b, msgStream());
                   });
         
         data.calibMeasurements.back()->setFitState(CalibratedSpacePoint::State::Outlier);
         data.nDoF -= data.calibMeasurements.back()->measuresEta();
         data.nDoF -= data.calibMeasurements.back()->measuresPhi();
         if (m_doT0Fit && data.calibMeasurements.back()->measuresTime() && 
                          data.calibMeasurements.back()->type() != xAOD::UncalibMeasType::MdtDriftCircleType) {
             --data.nDoF;
         }
         std::vector<HoughHitType> uncalib{};
         for (const HitVec::value_type& calib : data.calibMeasurements) {
            uncalib.push_back(calib->spacePoint());
         }
         SegmentFitResult newAttempt = fitSegmentHits(ctx, gctx, data.segmentPars, std::move(data.calibMeasurements));
         if (newAttempt.converged) {
             newAttempt.nIter+=data.nIter;
             data = std::move(newAttempt);
             if (m_visionTool.isEnabled()) {
                 auto seedCopy = convertToSegment(seed.msSector()->localToGlobalTrans(gctx), &seed, copy(data));
                 m_visionTool->visualizeSegment(ctx, *seedCopy, "Bad fit recovery");
             }
         } else {
             data.calibMeasurements = std::move(newAttempt.calibMeasurements);
         }
         return removeOutliers(ctx, gctx, seed, data);
     }

◆ renounce()

std::enable_if_t<std::is_void_v<std::result_of_t<decltype(&T::renounce)(T)> > && !std::is_base_of_v<SG::VarHandleKeyArray, T> && std::is_base_of_v<Gaudi::DataHandle, T>, void> AthCommonDataStore< AthCommonMsg< Gaudi::Algorithm > >::renounce ( T & h )

inlineprotectedinherited

Definition at line 380 of file AthCommonDataStore.h.

   {
     h.renounce();
     PBASE::renounce (h);
   }

◆ renounceArray()

void AthCommonDataStore< AthCommonMsg< Gaudi::Algorithm > >::renounceArray ( SG::VarHandleKeyArray & handlesArray )

inlineprotectedinherited

remove all handles from I/O resolution

Definition at line 364 of file AthCommonDataStore.h.

                                                           {
     handlesArray.renounce();
   }

◆ resolveAmbiguities()

void MuonR4::SegmentFittingAlg::resolveAmbiguities	(	const ActsGeometryContext &	gctx,
		std::vector< std::unique_ptr< Segment >> &	segmentCandidates
	)		const

private

Definition at line 558 of file SegmentFittingAlg.cxx.

                                                                                                            {
         using SegmentVec = std::vector<std::unique_ptr<Segment>>;
         ATH_MSG_VERBOSE("Resolve ambiguities amongst "<<segmentCandidates.size()<<" segment candidates. ");
         std::unordered_map<const MuonGMR4::SpectrometerSector*, SegmentVec> candidatesPerChamber{};
         
         for (std::unique_ptr<Segment>& sortMe : segmentCandidates) {
             const MuonGMR4::SpectrometerSector* chamb = sortMe->msSector();
             candidatesPerChamber[chamb].push_back(std::move(sortMe));
         }
         segmentCandidates.clear();
         for (auto& [chamber, resolveMe] : candidatesPerChamber) {
             SegmentVec resolvedSegments = m_ambiSolver->resolveAmbiguity(gctx, std::move(resolveMe));
             segmentCandidates.insert(segmentCandidates.end(), 
                                      std::make_move_iterator(resolvedSegments.begin()),
                                      std::make_move_iterator(resolvedSegments.end()));
         }
     }

◆ retrieveContainer()

template<class ContainerType >

StatusCode MuonR4::SegmentFittingAlg::retrieveContainer	(	const EventContext &	ctx,
		const SG::ReadHandleKey< ContainerType > &	key,
		const ContainerType *&	contToPush
	)		const

private

Helper method to fetch data from StoreGate.

If the key is empty, a nullptr is assigned to the container ptr Failure is returned in cases, of non-empty keys and failed retrieval

Definition at line 131 of file SegmentFittingAlg.cxx.

                                                                                                 {
             contToPush = nullptr;
             if (key.empty()) {
                 ATH_MSG_VERBOSE("No key has been parsed for object "<< typeid(ContainerType).name());
                 return StatusCode::SUCCESS;
             }
             SG::ReadHandle readHandle{key, ctx};
             ATH_CHECK(readHandle.isPresent());
             contToPush = readHandle.cptr();
             return StatusCode::SUCCESS;
         }

◆ setFilterPassed()

virtual void AthReentrantAlgorithm::setFilterPassed	(	bool	state,
		const EventContext &	ctx
	)		const

inlinevirtualinherited

Definition at line 139 of file AthReentrantAlgorithm.h.

                                                                             {
     execState( ctx ).setFilterPassed( state );
   }

◆ sysExecute()

StatusCode AthReentrantAlgorithm::sysExecute ( const EventContext & ctx )

overridevirtualinherited

Execute an algorithm.

We override this in order to work around an issue with the Algorithm base class storing the event context in a member variable that can cause crashes in MT jobs.

Definition at line 67 of file AthReentrantAlgorithm.cxx.

 {
   return Gaudi::Algorithm::sysExecute (ctx);
 }

◆ sysInitialize()

StatusCode AthReentrantAlgorithm::sysInitialize ( )

overridevirtualinherited

Override sysInitialize.

Override sysInitialize from the base class.

Loop through all output handles, and if they're WriteCondHandles, automatically register them and this Algorithm with the CondSvc

Scan through all outputHandles, and if they're WriteCondHandles, register them with the CondSvc

Reimplemented from AthCommonDataStore< AthCommonMsg< Gaudi::Algorithm > >.

Reimplemented in InputMakerBase, and HypoBase.

Definition at line 96 of file AthReentrantAlgorithm.cxx.

                                                 {
   StatusCode sc=AthCommonDataStore<AthCommonMsg<Gaudi::Algorithm>>::sysInitialize();
  
   if (sc.isFailure()) {
     return sc;
   }
   
   ServiceHandle<ICondSvc> cs("CondSvc",name());
   for (auto h : outputHandles()) {
     if (h->isCondition() && h->mode() == Gaudi::DataHandle::Writer) {
       // do this inside the loop so we don't create the CondSvc until needed
       if ( cs.retrieve().isFailure() ) {
         ATH_MSG_WARNING("no CondSvc found: won't autoreg WriteCondHandles");
         return StatusCode::SUCCESS;
       }      
       if (cs->regHandle(this,*h).isFailure()) {
         sc = StatusCode::FAILURE;
         ATH_MSG_ERROR("unable to register WriteCondHandle " << h->fullKey()
                       << " with CondSvc");
       }
     }
   }
   return sc;  
 }

◆ sysStart()

virtual StatusCode AthCommonDataStore< AthCommonMsg< Gaudi::Algorithm > >::sysStart ( )

overridevirtualinherited

Handle START transition.

We override this in order to make sure that conditions handle keys can cache a pointer to the conditions container.

◆ updateVHKA()

void AthCommonDataStore< AthCommonMsg< Gaudi::Algorithm > >::updateVHKA ( Gaudi::Details::PropertyBase & )

inlineinherited

Definition at line 308 of file AthCommonDataStore.h.

                                                {
     // debug() << "updateVHKA for property " << p.name() << " " << p.toString() 
     //         << "  size: " << m_vhka.size() << endmsg;
     for (auto &a : m_vhka) {
       std::vector<SG::VarHandleKey*> keys = a->keys();
       for (auto k : keys) {
         k->setOwner(this);
       }
     }
   }

Member Data Documentation

◆ m_ambiSolver

std::unique_ptr<SegmentAmbiSolver> MuonR4::SegmentFittingAlg::m_ambiSolver {}

private

Definition at line 139 of file SegmentFittingAlg.h.

◆ m_beamSpotL

Gaudi::Property<double> MuonR4::SegmentFittingAlg::m_beamSpotL {this, "BeamSpotLength", 20. * Gaudi::Units::m}

private

Definition at line 124 of file SegmentFittingAlg.h.

◆ m_beamSpotR

Gaudi::Property<double> MuonR4::SegmentFittingAlg::m_beamSpotR {this, "BeamSpotRadius", 30.* Gaudi::Units::cm}

private

Definition at line 123 of file SegmentFittingAlg.h.

◆ m_calibTool

ToolHandle<ISpacePointCalibrator> MuonR4::SegmentFittingAlg::m_calibTool {this, "Calibrator", "" }

private

Handle to the space point calibrator.

Definition at line 113 of file SegmentFittingAlg.h.

◆ m_detStore

StoreGateSvc_t AthCommonDataStore< AthCommonMsg< Gaudi::Algorithm > >::m_detStore

privateinherited

Pointer to StoreGate (detector store by default)

Definition at line 393 of file AthCommonDataStore.h.

◆ m_doBeamspotConstraint

Gaudi::Property<bool> MuonR4::SegmentFittingAlg::m_doBeamspotConstraint {this, "doBeamspotConstraint", false}

private

Add beamline constraint.

Definition at line 122 of file SegmentFittingAlg.h.

◆ m_doT0Fit

Gaudi::Property<bool> MuonR4::SegmentFittingAlg::m_doT0Fit {this, "fitSegmentT0", true}

private

Definition at line 120 of file SegmentFittingAlg.h.

◆ m_evtStore

StoreGateSvc_t AthCommonDataStore< AthCommonMsg< Gaudi::Algorithm > >::m_evtStore

privateinherited

Pointer to StoreGate (event store by default)

Definition at line 390 of file AthCommonDataStore.h.

◆ m_extendedExtraObjects

DataObjIDColl AthReentrantAlgorithm::m_extendedExtraObjects

privateinherited

Extra output dependency collection, extended by AthAlgorithmDHUpdate to add symlinks.

Empty if no symlinks were found.

Definition at line 153 of file AthReentrantAlgorithm.h.

◆ m_geoCtxKey

SG::ReadHandleKey<ActsGeometryContext> MuonR4::SegmentFittingAlg::m_geoCtxKey {this, "AlignmentKey", "ActsAlignment", "cond handle key"}

private

Definition at line 109 of file SegmentFittingAlg.h.

◆ m_idHelperSvc

ServiceHandle<Muon::IMuonIdHelperSvc> MuonR4::SegmentFittingAlg::m_idHelperSvc {this, "MuonIdHelperSvc", "Muon::MuonIdHelperSvc/MuonIdHelperSvc"}

private

IdHelperSvc.

Definition at line 111 of file SegmentFittingAlg.h.

◆ m_outlierRemovalCut

Gaudi::Property<double> MuonR4::SegmentFittingAlg::m_outlierRemovalCut {this, "OutlierRemoval", 5.}

private

Cut on the segment chi2 / nDoF to launch the outlier removal.

Definition at line 136 of file SegmentFittingAlg.h.

◆ m_outSegments

SG::WriteHandleKey<SegmentContainer> MuonR4::SegmentFittingAlg::m_outSegments {this, "MuonSegmentContainer", "R4MuonSegments"}

private

Definition at line 107 of file SegmentFittingAlg.h.

◆ m_recalibSeed

Gaudi::Property<bool> MuonR4::SegmentFittingAlg::m_recalibSeed {this, "SeedRecalibrate", false}

private

Toggle seed recalibration.

The two seed circles are recalibrated using the initial seed

Definition at line 131 of file SegmentFittingAlg.h.

◆ m_recoveryPull

Gaudi::Property<double> MuonR4::SegmentFittingAlg::m_recoveryPull {this, "RecoveryPull", 5.}

private

Definition at line 137 of file SegmentFittingAlg.h.

◆ m_refineSeed

Gaudi::Property<bool> MuonR4::SegmentFittingAlg::m_refineSeed {this, "SeedRefine", false}

private

Toggle seed refit.

The segment seed is fastly refitted using the collected seed drift circles

Definition at line 134 of file SegmentFittingAlg.h.

◆ m_seedHitChi2

Gaudi::Property<double> MuonR4::SegmentFittingAlg::m_seedHitChi2 {this, "ResoSeedHitAssoc", 5. }

private

Two mdt seeds are the same if their defining parameters match wihin.

Definition at line 128 of file SegmentFittingAlg.h.

◆ m_seedKey

SG::ReadHandleKey<SegmentSeedContainer> MuonR4::SegmentFittingAlg::m_seedKey {this, "ReadKey", "MuonHoughStationSegmentSeeds"}

private

ReadHandle of the seeds.

Definition at line 105 of file SegmentFittingAlg.h.

◆ m_useMinuit

Gaudi::Property<bool> MuonR4::SegmentFittingAlg::m_useMinuit {this, "useMinuit", false}

private

Toggle the fitter.

Definition at line 118 of file SegmentFittingAlg.h.

◆ m_varHandleArraysDeclared

bool AthCommonDataStore< AthCommonMsg< Gaudi::Algorithm > >::m_varHandleArraysDeclared

privateinherited

Definition at line 399 of file AthCommonDataStore.h.

◆ m_vhka

std::vector<SG::VarHandleKeyArray*> AthCommonDataStore< AthCommonMsg< Gaudi::Algorithm > >::m_vhka

privateinherited

Definition at line 398 of file AthCommonDataStore.h.

◆ m_visionTool

ToolHandle<MuonValR4::IPatternVisualizationTool> MuonR4::SegmentFittingAlg::m_visionTool {this, "VisualizationTool", ""}

private

Pattern visualization tool.

Definition at line 115 of file SegmentFittingAlg.h.

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

Public Types

Public Member Functions

Protected Member Functions

Private Types

Private Member Functions

Static Private Member Functions

Private Attributes

Detailed Description

Member Typedef Documentation

◆ HitVec

◆ Parameters

◆ StoreGateSvc_t

Constructor & Destructor Documentation

◆ ~SegmentFittingAlg()

Member Function Documentation

◆ AthReentrantAlgorithm()

◆ cardinality()

◆ convertToSegment()

◆ declareGaudiProperty() [1/4]

◆ declareGaudiProperty() [2/4]

◆ declareGaudiProperty() [3/4]

◆ declareGaudiProperty() [4/4]

◆ declareProperty() [1/6]

◆ declareProperty() [2/6]

◆ declareProperty() [3/6]

◆ declareProperty() [4/6]

◆ declareProperty() [5/6]

◆ declareProperty() [6/6]

◆ detStore()

◆ eraseWrongHits()

◆ evtStore() [1/2]

◆ evtStore() [2/2]

◆ execute()

◆ extraDeps_update_handler()

◆ extraOutputDeps()

◆ filterPassed()

◆ fitSegmentHits()

◆ fitSegmentSeed()

◆ initialize()

◆ inputHandles()

◆ isClonable()

◆ msg() [1/2]

◆ msg() [2/2]

◆ msgLvl()

◆ outputHandles()

◆ plugHoles()

◆ removeOutliers()

◆ renounce()

◆ renounceArray()

◆ resolveAmbiguities()

◆ retrieveContainer()

◆ setFilterPassed()

◆ sysExecute()

◆ sysInitialize()

◆ sysStart()

◆ updateVHKA()

Member Data Documentation

◆ m_ambiSolver

◆ m_beamSpotL

◆ m_beamSpotR

◆ m_calibTool

◆ m_detStore

◆ m_doBeamspotConstraint

◆ m_doT0Fit

◆ m_evtStore

◆ m_extendedExtraObjects

◆ m_geoCtxKey

◆ m_idHelperSvc

◆ m_outlierRemovalCut

◆ m_outSegments

◆ m_recalibSeed

◆ m_recoveryPull

◆ m_refineSeed

◆ m_seedHitChi2

◆ m_seedKey

◆ m_useMinuit

◆ m_varHandleArraysDeclared

◆ m_vhka

◆ m_visionTool