ATLAS Offline Software
SegmentFittingAlg.cxx
Go to the documentation of this file.
1 /*
2  Copyright (C) 2002-2025 CERN for the benefit of the ATLAS collaboration
3 */
4 
5 #include "SegmentFittingAlg.h"
6 
7 #include <TrkEventPrimitives/ParticleHypothesis.h>
8 #include <GeoPrimitives/GeoPrimitivesHelpers.h>
9 
10 #include <MuonPatternHelpers/SegmentFitHelperFunctions.h>
11 #include <MuonPatternHelpers/MdtSegmentSeedGenerator.h>
12 #include <MuonPatternHelpers/MdtSegmentFitter.h>
13 #include "xAODMuonPrepData/MdtDriftCircleContainer.h"
14 #include "xAODMuonPrepData/RpcStripContainer.h"
15 #include <MuonSpacePoint/SpacePointPerLayerSplitter.h>
16 #include <MuonSpacePoint/UtilFunctions.h>
17 
18 #include <xAODMuonPrepData/RpcMeasurement.h>
19 #include <xAODMuonPrepData/TgcStrip.h>
20 #include "xAODMuonSimHit/MuonSimHitContainer.h"
21 
22 #include <GaudiKernel/PhysicalConstants.h>
23 #include <Minuit2/Minuit2Minimizer.h>
24 #include <Math/Minimizer.h>
25 
26 #include <MuonVisualizationHelpersR4/VisualizationHelpers.h>
27 
28 #include <format>
29 
30 namespace MuonR4 {
31  using namespace SegmentFit;
32  using namespace MuonValR4;
33 
34  MuonR4::SegmentFitResult::HitVec copy(const MuonR4::SegmentFitResult::HitVec& hits) {
35  MuonR4::SegmentFitResult::HitVec copied{};
36  copied.reserve(hits.size());
37  std::ranges::transform(hits, std::back_inserter(copied),
38  [](const auto& hit) { return std::make_unique<MuonR4::CalibratedSpacePoint>(*hit);});
39  return copied;
40  }
41  MuonR4::SegmentFitResult copy(const MuonR4::SegmentFitResult& toCopy) {
42  MuonR4::SegmentFitResult toRet{};
43  toRet.segmentPars = toCopy.segmentPars;
44  toRet.calibMeasurements = copy(toCopy.calibMeasurements);
45  toRet.chi2 = toCopy.chi2;
46  toRet.nDoF = toCopy.nDoF;
47  toRet.timeFit = toCopy.timeFit;
48  return toRet;
49  }
50  bool removeBeamSpot(MuonR4::SegmentFitResult::HitVec& hits){
51  std::size_t before = hits.size();
52  hits.erase(std::remove_if(hits.begin(), hits.end(),
53  [](const auto& a){
54  return a->type() == xAOD::UncalibMeasType::Other;
55  }), hits.end());
56  return before != hits.size();
57  }
58 
59  using PrimitiveVec = MuonValR4::IPatternVisualizationTool::PrimitiveVec;
60 
61  SegmentFittingAlg::~SegmentFittingAlg() = default;
62  StatusCode SegmentFittingAlg::initialize() {
63  ATH_CHECK(m_geoCtxKey.initialize());
64  ATH_CHECK(m_seedKey.initialize());
65  ATH_CHECK(m_outSegments.initialize());
66  ATH_CHECK(m_calibTool.retrieve());
67  ATH_CHECK(m_idHelperSvc.retrieve());
68  ATH_CHECK(m_visionTool.retrieve(EnableTool{!m_visionTool.empty()}));
69  SegmentAmbiSolver::Config cfg{};
70  m_ambiSolver = std::make_unique<SegmentAmbiSolver>(name(), std::move(cfg));
71  return StatusCode::SUCCESS;
72  }
73  StatusCode SegmentFittingAlg::execute(const EventContext& ctx) const {
74  const ActsGeometryContext* gctx{nullptr};
75  ATH_CHECK(SG::get(gctx, m_geoCtxKey, ctx));
76  const SegmentSeedContainer* segmentSeeds=nullptr;
77  ATH_CHECK(SG::get(segmentSeeds, m_seedKey, ctx));
78 
79  SG::WriteHandle writeSegments{m_outSegments, ctx};
80  ATH_CHECK(writeSegments.record(std::make_unique<SegmentContainer>()));
81 
82  std::vector<std::unique_ptr<Segment>> allSegments{};
83  for (const SegmentSeed* seed : *segmentSeeds) {
84  std::vector<std::unique_ptr<Segment>> segments = fitSegmentSeed(ctx, *gctx, seed);
85  if (m_visionTool.isEnabled() && segments.size() > 1) {
86  auto drawFinalReco = [this, &segments, &gctx, &ctx,&seed](const std::string& nameTag) {
87  PrimitiveVec segmentLines{};
88  double yLegend{0.85};
89  segmentLines.push_back(drawLabel(std::format("# segments: {:d}", segments.size()), 0.2, yLegend, 14));
90  yLegend-=0.04;
91  for (const std::unique_ptr<Segment>& seg : segments) {
92  const Parameters pars = localSegmentPars(*gctx, *seg);
93  const auto [locPos, locDir] = makeLine(pars);
94 
95  segmentLines.emplace_back(drawLine(pars, -Gaudi::Units::m, Gaudi::Units::m, kRed));
96  std::stringstream signStream{};
97  signStream<<std::format("#chi^{{2}}/nDoF: {:.2f} ({:}), ", seg->chi2() / seg->nDoF(), seg->nDoF());
98  signStream<<std::format("y_{{0}}={:.2f}",pars[toInt(ParamDefs::y0)])<<", ";
99  signStream<<std::format("#theta={:.2f}^{{#circ}}", pars[toInt(ParamDefs::theta)]/ Gaudi::Units::deg )<<", ";
100  for (const Segment::MeasType& m : seg->measurements()) {
101  if (m->type() == xAOD::UncalibMeasType::MdtDriftCircleType && m->fitState() == CalibratedSpacePoint::State::Valid) {
102  signStream<<(SegmentFitHelpers::driftSign(locPos, locDir, *m, msgStream()) == -1 ? "L" : "R");
103  }
104  }
105  segmentLines.push_back(drawLabel(signStream.str(), 0.2, yLegend, 13));
106  yLegend-=0.03;
107  }
108 
109  m_visionTool->visualizeBucket(ctx, *seed->parentBucket(), nameTag, std::move(segmentLines));
110  };
111  drawFinalReco("all segments");
112  const unsigned int nBeforeAmbi = segments.size();
113  segments = m_ambiSolver->resolveAmbiguity(*gctx, std::move(segments));
114  if (nBeforeAmbi != segments.size()) {
115  drawFinalReco("post ambiguity");
116  }
117  } else if (m_visionTool.isEnabled() && segments.empty() &&
118  std::ranges::count_if(seed->getHitsInMax(),[this](const SpacePoint* hit){
119  return m_visionTool->isLabeled(*hit);
120  })) {
121  m_visionTool->visualizeSeed(ctx, *seed, "Failed fit");
122  }
123  allSegments.insert(allSegments.end(), std::make_move_iterator(segments.begin()),
124  std::make_move_iterator(segments.end()));
125  }
126  resolveAmbiguities(*gctx, allSegments);
127  writeSegments->insert(writeSegments->end(),
128  std::make_move_iterator(allSegments.begin()),
129  std::make_move_iterator(allSegments.end()));
130  ATH_MSG_VERBOSE("Found in total "<<writeSegments->size()<<" segments. ");
131  return StatusCode::SUCCESS;
132  }
133 
134  SegmentFitResult SegmentFittingAlg::fitSegmentHits(const EventContext& ctx,
135  const ActsGeometryContext& gctx,
136  const Parameters& startPars,
137  SegmentFitResult::HitVec&& calibHits) const {
138  SegmentFitResult data{};
139  if (calibHits.empty()) {
140  ATH_MSG_WARNING("No hits, no segment");
141  return data;
142  }
143 
145  if (m_doBeamspotConstraint) {
146  unsigned int numPhi = std::accumulate(calibHits.begin(), calibHits.end(),0,
147  [](unsigned int n, const HitVec::value_type& hit){
148  return n + (hit->fitState() == CalibratedSpacePoint::State::Valid &&
149  hit->measuresPhi());
150  });
151  if (numPhi) {
152  const Amg::Transform3D globToLoc{calibHits[0]->spacePoint()->msSector()->globalToLocalTrans(gctx)};
153  Amg::Vector3D beamSpot{globToLoc.translation()};
154  AmgSymMatrix(3) covariance{AmgSymMatrix(3)::Identity()};
156  covariance(0,0) = std::pow(m_beamSpotR, 2);
157  covariance(1,1) = std::pow(m_beamSpotR, 2);
158  covariance(2,2) = std::pow(m_beamSpotL, 2);
159  AmgSymMatrix(3) jacobian = globToLoc.linear();
160  covariance = jacobian * covariance * jacobian.transpose();
161  AmgSymMatrix(2) beamSpotCov{covariance.block<2,2>(0,0)};
162  auto beamSpotSP = std::make_unique<CalibratedSpacePoint>(nullptr, std::move(beamSpot), Amg::Vector3D::Zero());
163  beamSpotSP->setCovariance<2>(std::move(beamSpotCov));
164  ATH_MSG_VERBOSE("Beam spot constraint "<<Amg::toString(beamSpotSP->positionInChamber())<<", "<<toString(beamSpotSP->covariance()));
165  calibHits.push_back(std::move(beamSpotSP));
166  }
167  }
168 
169  const Amg::Transform3D& locToGlob{calibHits[0]->spacePoint()->msSector()->localToGlobalTrans(gctx)};
170 
171  MdtSegmentFitter::Config fitCfg{};
172  fitCfg.calibrator = m_calibTool.get();
173  fitCfg.doTimeFit = m_doT0Fit;
174  fitCfg.reCalibrate = m_recalibInFit;
175 
176  MdtSegmentFitter fitter{name(), std::move(fitCfg)};
177  return fitter.fitSegment(ctx, std::move(calibHits), startPars, locToGlob);
178  }
179  std::vector<std::unique_ptr<Segment>>
180  SegmentFittingAlg::fitSegmentSeed(const EventContext& ctx,
181  const ActsGeometryContext& gctx,
182  const SegmentSeed* patternSeed) const {
183 
184  const Amg::Transform3D& locToGlob{patternSeed->msSector()->localToGlobalTrans(gctx)};
185  std::vector<std::unique_ptr<Segment>> segments{};
186 
187  MdtSegmentSeedGenerator::Config genCfg{};
188  genCfg.hitPullCut = m_seedHitChi2;
189  genCfg.recalibSeedCircles = m_recalibSeed;
190  genCfg.fastSeedFit = m_refineSeed;
191  genCfg.fastSegFitWithT0 = m_doT0Fit;
192  genCfg.calibrator = m_calibTool.get();
193 
194 
196  genCfg.busyLayerLimit = 2 + 2*(patternSeed->parameters()[toInt(ParamDefs::theta)] > 50 * Gaudi::Units::deg);
198  if (m_visionTool.isEnabled()) {
199  PrimitiveVec seedLines{};
200  MdtSegmentSeedGenerator drawMe{name(), patternSeed, genCfg};
201  while(auto s = drawMe.nextSeed(ctx)) {
202  seedLines.push_back(drawLine(s->parameters, -Gaudi::Units::m, Gaudi::Units::m, kViolet));
203  }
204  seedLines.push_back(drawLabel(std::format("possible seeds: {:d}", drawMe.numGenerated()), 0.2, 0.85, 14));
205  m_visionTool->visualizeSeed(ctx, *patternSeed, "pattern", std::move(seedLines));
206  }
207 
208  MdtSegmentSeedGenerator seedGen{name(), patternSeed, std::move(genCfg)};
209  ATH_MSG_VERBOSE("fitSegmentHits() - Start segment seed search");
210  while (auto seed = seedGen.nextSeed(ctx)) {
211  SegmentFitResult data{};
212  data.segmentPars = seed->parameters;
213  data.calibMeasurements = std::move(seed->measurements);
215  if (m_visionTool.isEnabled()) {
216  auto seedCopy = convertToSegment(locToGlob, patternSeed, copy(data));
217  m_visionTool->visualizeSegment(ctx, *seedCopy, std::format("Pre fit {:d}", seedGen.numGenerated()));
218  }
219  data = fitSegmentHits(ctx, gctx, seed->parameters, std::move(data.calibMeasurements));
220  data.nIter += seed->nIter;
221  if (m_visionTool.isEnabled() && data.converged) {
222  auto seedCopy = convertToSegment(locToGlob, patternSeed, copy(data));
223  m_visionTool->visualizeSegment(ctx, *seedCopy, std::format("Intermediate fit {:d}", seedGen.numGenerated()));
224  }
225  if (!removeOutliers(ctx, gctx, *patternSeed, data)) {
226  continue;
227  }
228  if (!plugHoles(ctx, gctx, *patternSeed, data)) {
229  continue;
230  }
231  if (m_visionTool.isEnabled()) {
232  auto seedCopy = convertToSegment(locToGlob, patternSeed, copy(data));
233  m_visionTool->visualizeSegment(ctx, *seedCopy, std::format("Final fit {:d}", seedGen.numGenerated()));
234  }
235  segments.push_back(convertToSegment(locToGlob, patternSeed, std::move(data)));
236  }
237  ATH_MSG_VERBOSE("fitSegmentHits() - In total "<<segments.size()<<" segment were constructed ");
238 
239  return segments;
240  }
241  std::unique_ptr<Segment> SegmentFittingAlg::convertToSegment(const Amg::Transform3D& locToGlob,
242  const SegmentSeed* patternSeed,
243  SegmentFitResult&& data) {
244  const auto [locPos, locDir] = data.makeLine();
245  Amg::Vector3D globPos = locToGlob * locPos;
246  Amg::Vector3D globDir = locToGlob.linear()* locDir;
247 
248 
249  auto finalSeg = std::make_unique<Segment>(std::move(globPos), std::move(globDir),
250  patternSeed, std::move(data.calibMeasurements),
251  data.chi2, data.nDoF);
252  finalSeg->setCallsToConverge(data.nIter);
253  finalSeg->setParUncertainties(std::move(data.segmentParErrs));
254  if (data.timeFit) {
255  finalSeg->setSegmentT0(data.segmentPars[toInt(ParamDefs::time)]);
256  }
257  return finalSeg;
258  }
259 
260  bool SegmentFittingAlg::removeOutliers(const EventContext& ctx,
261  const ActsGeometryContext& gctx,
262  const SegmentSeed& seed,
263  SegmentFitResult& data) const {
264 
266  if (data.nDoF<=0 || data.calibMeasurements.empty() || data.nPrecMeas < m_precHitCut) {
267  ATH_MSG_VERBOSE("No degree of freedom available. What shall be removed?!. nDoF: "
268  <<data.nDoF<<", n-meas: "<<data.calibMeasurements);
269  return false;
270  }
271 
272  const auto [segPos, segDir] = data.makeLine();
273 
274  if (data.converged && data.chi2 / data.nDoF < m_outlierRemovalCut) {
275  ATH_MSG_VERBOSE("The segment "<<Amg::toString(segPos)<<" + "<<Amg::toString(segDir)
276  <<" is already of good quality "<<data.chi2 / std::max(data.nDoF, 1)
277  <<". Don't remove outliers");
278  return true;
279  }
280  ATH_MSG_VERBOSE("Segment "<<toString(data.segmentPars)<<" is of badish quality.");
283  if (m_doBeamspotConstraint && removeBeamSpot(data.calibMeasurements)) {
284  data.nDoF-=2;
285  data.nPhiMeas-=1;
286  }
287 
289  std::ranges::sort(data.calibMeasurements,
290  [&, this](const HitVec::value_type& a, const HitVec::value_type& b){
291  return SegmentFitHelpers::chiSqTerm(segPos, segDir, data.segmentPars[toInt(ParamDefs::time)], std::nullopt, *a, msgStream()) <
292  SegmentFitHelpers::chiSqTerm(segPos, segDir, data.segmentPars[toInt(ParamDefs::time)], std::nullopt, *b, msgStream());
293  });
294 
296  data.calibMeasurements.back()->setFitState(CalibratedSpacePoint::State::Outlier);
297  data.nDoF -= data.calibMeasurements.back()->measuresEta();
298  data.nDoF -= data.calibMeasurements.back()->measuresPhi();
299  if (m_doT0Fit && data.calibMeasurements.back()->measuresTime() &&
300  data.calibMeasurements.back()->type() != xAOD::UncalibMeasType::MdtDriftCircleType) {
301  --data.nDoF;
302  }
304  std::vector<HoughHitType> uncalib{};
305  for (const HitVec::value_type& calib : data.calibMeasurements) {
306  uncalib.push_back(calib->spacePoint());
307  }
308  SegmentFitResult newAttempt = fitSegmentHits(ctx, gctx, data.segmentPars, std::move(data.calibMeasurements));
309  if (newAttempt.converged) {
310  newAttempt.nIter+=data.nIter;
311  data = std::move(newAttempt);
312  if (m_visionTool.isEnabled()) {
313  auto seedCopy = convertToSegment(seed.msSector()->localToGlobalTrans(gctx), &seed, copy(data));
314  m_visionTool->visualizeSegment(ctx, *seedCopy, "Bad fit recovery");
315  }
316  } else {
317  data.calibMeasurements = std::move(newAttempt.calibMeasurements);
318  }
319  return removeOutliers(ctx, gctx, seed, data);
320  }
321 
322  bool SegmentFittingAlg::plugHoles(const EventContext& ctx,
323  const ActsGeometryContext& gctx,
324  const SegmentSeed& seed,
325  SegmentFitResult& beforeRecov) const {
327  ATH_MSG_VERBOSE("plugHoles() -- segment "<<toString(beforeRecov.segmentPars)
328  <<", chi2: "<<beforeRecov.chi2<<", nDoF: "<<beforeRecov.nDoF<<", "
329  <<beforeRecov.chi2 /std::max(beforeRecov.nDoF, 1) );
331  std::unordered_set<const SpacePoint*> usedSpacePoint{};
332  for (const HitVec::value_type& hit : beforeRecov.calibMeasurements) {
333  usedSpacePoint.insert(hit->spacePoint());
334  }
335 
336  HitVec candidateHits{};
337  SpacePointPerLayerSplitter hitLayers{*seed.parentBucket()};
338  bool hasCandidate{false};
339  const auto [locPos, locDir] = beforeRecov.makeLine();
340  for (const std::vector<HoughHitType>& mdtLayer : hitLayers.mdtHits()) {
341  for (const SpacePoint* mdtHit: mdtLayer) {
343  if (usedSpacePoint.count(mdtHit)) {
344  continue;
345  }
346  const double dist = Amg::lineDistance(locPos, locDir, mdtHit->positionInChamber(), mdtHit->directionInChamber());
347  const auto* dc = static_cast<const xAOD::MdtDriftCircle*>(mdtHit->primaryMeasurement());
348  if (dist >= dc->readoutElement()->innerTubeRadius()) {
349  continue;
350  }
351  HitVec::value_type calibHit{m_calibTool->calibrate(ctx, mdtHit, locPos, locDir, beforeRecov.segmentPars[toInt(ParamDefs::time)])};
352  const double pull = std::sqrt(SegmentFitHelpers::chiSqTermMdt(locPos, locDir, *calibHit, msgStream()));
353  ATH_MSG_VERBOSE(__func__<<"() :"<<__LINE__<<" Candidate hit for recovery "<<m_idHelperSvc->toString(mdtHit->identify())<<", chi2: "<<pull);
354  if (pull <= m_recoveryPull) {
355  hasCandidate |= calibHit->fitState() == CalibratedSpacePoint::State::Valid;
356  candidateHits.push_back(std::move(calibHit));
357  } else {
358  calibHit->setFitState(CalibratedSpacePoint::State::Outlier);
359  candidateHits.push_back(std::move(calibHit));
360  }
361  }
362  }
363  if (!hasCandidate) {
364  ATH_MSG_VERBOSE("No space point candidates for recovery were found");
365  beforeRecov.calibMeasurements.insert(beforeRecov.calibMeasurements.end(),
366  std::make_move_iterator(candidateHits.begin()),
367  std::make_move_iterator(candidateHits.end()));
368  eraseWrongHits(beforeRecov);
369  return beforeRecov.nDoF > 0;
370  }
371 
372  HitVec copied = copy(beforeRecov.calibMeasurements), copiedCandidates = copy(candidateHits);
374  if (m_doBeamspotConstraint) {
375  removeBeamSpot(copied);
376  }
377 
378  candidateHits.insert(candidateHits.end(), std::make_move_iterator(copied.begin()), std::make_move_iterator(copied.end()));
379 
380  SegmentFitResult recovered = fitSegmentHits(ctx, gctx, beforeRecov.segmentPars, std::move(candidateHits));
381  if (!recovered.converged) {
382  return false;
383  }
385  if (recovered.nDoF + recovered.timeFit <= beforeRecov.nDoF + beforeRecov.timeFit) {
386  for (HitVec::value_type& hit : copiedCandidates) {
387  hit->setFitState(CalibratedSpacePoint::State::Outlier);
388  beforeRecov.calibMeasurements.push_back(std::move(hit));
389  }
390  eraseWrongHits(beforeRecov);
391  return true;
392  }
393  ATH_MSG_VERBOSE("Chi2, nDOF before:"<<beforeRecov.chi2<<", "<<beforeRecov.nDoF<<" after recovery: "<<recovered.chi2<<", "<<recovered.nDoF);
394  double redChi2 = recovered.chi2 / std::max(recovered.nDoF,1);
396  if (redChi2 < m_outlierRemovalCut || (beforeRecov.nDoF == 0) || redChi2 < beforeRecov.chi2 / beforeRecov.nDoF) {
397  ATH_MSG_VERBOSE("Accept segment with recovered "<<(recovered.nDoF + recovered.timeFit) - (beforeRecov.nDoF + beforeRecov.timeFit)<<" hits.");
398  recovered.nIter += beforeRecov.nIter;
399  beforeRecov = std::move(recovered);
401  while (true) {
402  bool runAnotherTrial = false;
403  copied = copy(beforeRecov.calibMeasurements);
404  if (m_doBeamspotConstraint) {
405  removeBeamSpot(copied);
406  }
407  const auto [beforePos, beforeDir] = beforeRecov.makeLine();
408  for (HitVec::value_type& copyHit : copied) {
409  if (copyHit->fitState() != CalibratedSpacePoint::State::Outlier) {
410  continue;
411  }
412  copyHit->setFitState(CalibratedSpacePoint::State::Valid);
413  bool append = std::sqrt(SegmentFitHelpers::chiSqTerm(beforePos, beforeDir,
414  beforeRecov.segmentPars[toInt(ParamDefs::time)],
415  std::nullopt, *copyHit, msgStream())) < m_recoveryPull;
416  if (append) {
417  runAnotherTrial = true;
418  } else {
419  copyHit->setFitState(CalibratedSpacePoint::State::Outlier);
420  }
421  }
422  if (!runAnotherTrial) {
423  break;
424  }
425  recovered = fitSegmentHits(ctx, gctx, beforeRecov.segmentPars, std::move(copied));
426  if (!recovered.converged) {
427  break;
428  }
429  if (recovered.nDoF + recovered.timeFit <= beforeRecov.nDoF + beforeRecov.timeFit) {
430  break;
431  }
432  redChi2 = recovered.chi2 / std::max(recovered.nDoF, 1);
433  if (redChi2 < m_outlierRemovalCut || redChi2 < beforeRecov.chi2 / beforeRecov.nDoF) {
434  recovered.nIter += beforeRecov.nIter;
435  beforeRecov = std::move(recovered);
436  } else {
437  break;
438  }
439  }
441  eraseWrongHits(beforeRecov);
442  } else{
443  for (HitVec::value_type& hit : copiedCandidates) {
444  hit->setFitState(CalibratedSpacePoint::State::Outlier);
445  beforeRecov.calibMeasurements.push_back(std::move(hit));
446  }
447  }
448  return true;
449  }
450  void SegmentFittingAlg::eraseWrongHits(SegmentFitResult& candidate) const {
451  auto [segPos, segDir] = makeLine(candidate.segmentPars);
452  candidate.calibMeasurements.erase(std::remove_if(candidate.calibMeasurements.begin(), candidate.calibMeasurements.end(),
453  [&segPos, &segDir](const HitVec::value_type& hit){
454  if (hit->fitState() != CalibratedSpacePoint::State::Outlier) {
455  return false;
456  }
458  if (hit->type() == xAOD::UncalibMeasType::MdtDriftCircleType) {
459  const double dist = Amg::lineDistance(segPos, segDir, hit->positionInChamber(), hit->directionInChamber());
460  const auto* dc = static_cast<const xAOD::MdtDriftCircle*>(hit->spacePoint()->primaryMeasurement());
461  return dist >= dc->readoutElement()->innerTubeRadius();
462  }
463  return true;
464  }), candidate.calibMeasurements.end());
465 
466  std::ranges::sort(candidate.calibMeasurements, [](const Segment::MeasType& a, const Segment::MeasType& b){
467  return a->positionInChamber().z() < b->positionInChamber().z();
468  });
469  }
470  void SegmentFittingAlg::resolveAmbiguities(const ActsGeometryContext& gctx,
471  std::vector<std::unique_ptr<Segment>>& segmentCandidates) const {
472  using SegmentVec = std::vector<std::unique_ptr<Segment>>;
473  ATH_MSG_VERBOSE("Resolve ambiguities amongst "<<segmentCandidates.size()<<" segment candidates. ");
474  std::unordered_map<const MuonGMR4::SpectrometerSector*, SegmentVec> candidatesPerChamber{};
475 
476  for (std::unique_ptr<Segment>& sortMe : segmentCandidates) {
477  const MuonGMR4::SpectrometerSector* chamb = sortMe->msSector();
478  candidatesPerChamber[chamb].push_back(std::move(sortMe));
479  }
480  segmentCandidates.clear();
481  for (auto& [chamber, resolveMe] : candidatesPerChamber) {
482  SegmentVec resolvedSegments = m_ambiSolver->resolveAmbiguity(gctx, std::move(resolveMe));
483  segmentCandidates.insert(segmentCandidates.end(),
484  std::make_move_iterator(resolvedSegments.begin()),
485  std::make_move_iterator(resolvedSegments.end()));
486  }
487  }
488 }
MuonR4::SegmentFitResult::nIter
unsigned nIter
Number of iterations called to reach the minimum.
Definition: SegmentFitterEventData.h:78
MuonR4::PrimitiveVec
MuonValR4::IPatternVisualizationTool::PrimitiveVec PrimitiveVec
Definition: SegmentFittingAlg.cxx:59
make_hlt_rep.pars
pars
Definition: make_hlt_rep.py:90
LArSamples::FitterData::fitter
const ShapeFitter * fitter
Definition: ShapeFitter.cxx:23
MuonR4::MdtSegmentFitter
Definition: MdtSegmentFitter.h:19
data
char data[hepevt_bytes_allocation_ATLAS]
Definition: HepEvt.cxx:11
MuonR4::SegmentFittingAlg::convertToSegment
static std::unique_ptr< Segment > convertToSegment(const Amg::Transform3D &locToGlobTrf, const SegmentSeed *parentSeed, SegmentFitResult &&toConvert)
Converts the fit result into a segment object.
Definition: SegmentFittingAlg.cxx:241
TRTCalib_Extractor.hits
hits
Definition: TRTCalib_Extractor.py:35
MuonR4::SpacePointPerLayerSplitter
The SpacePointPerLayerSplitter takes a set of spacepoints already sorted by layer Identifier (see Muo...
Definition: SpacePointPerLayerSplitter.h:16
MuonGMR4::SpectrometerSector
A spectrometer sector forms the envelope of all chambers that are placed in the same MS sector & laye...
Definition: SpectrometerSector.h:40
vtune_athena.format
format
Definition: vtune_athena.py:14
calibdata.chamber
chamber
Definition: calibdata.py:31
VisualizationHelpers.h
max
constexpr double max()
Definition: ap_fixedTest.cxx:33
MuonR4::SegmentFittingAlg::fitSegmentHits
SegmentFitResult fitSegmentHits(const EventContext &ctx, const ActsGeometryContext &gctx, const Parameters &startPars, SegmentFitResult::HitVec &&calibHits) const
Executes the segment fit with start parameters.
Definition: SegmentFittingAlg.cxx:134
MuonR4::SegmentFitResult::calibMeasurements
HitVec calibMeasurements
Calibrated measurements used in the fit.
Definition: SegmentFitterEventData.h:64
SegmentFittingAlg.h
deg
#define deg
Definition: SbPolyhedron.cxx:17
MuonR4::SegmentFittingAlg::initialize
virtual StatusCode initialize() override
Definition: SegmentFittingAlg.cxx:62
MuonR4::SegmentFittingAlg::HitVec
SegmentFitResult::HitVec HitVec
Definition: SegmentFittingAlg.h:38
MuonR4::CalibratedSpacePoint::State::Outlier
@ Outlier
ATH_MSG_VERBOSE
#define ATH_MSG_VERBOSE(x)
Definition: AthMsgStreamMacros.h:28
MuonR4::SegmentAmbiSolver::Config
Definition: SegmentAmbiSolver.h:18
MuonR4::SegmentFittingAlg::Parameters
SegmentFit::Parameters Parameters
Definition: SegmentFittingAlg.h:41
dumpHVPathFromNtuple.append
bool append
Definition: dumpHVPathFromNtuple.py:91
MuonR4::SegmentFittingAlg::resolveAmbiguities
void resolveAmbiguities(const ActsGeometryContext &gctx, std::vector< std::unique_ptr< Segment >> &segmentCandidates) const
Definition: SegmentFittingAlg.cxx:470
MuonR4::SegmentFitHelpers::chiSqTerm
double chiSqTerm(const Amg::Vector3D &posInChamber, const Amg::Vector3D &dirInChamber, const SpacePoint &measurement, MsgStream &msg)
Calculates the chi2 contribuation to a linear segment line from an uncalibrated measurement.
Definition: SegmentFitHelperFunctions.cxx:28
MuonR4::SegmentFit::makeLine
std::pair< Amg::Vector3D, Amg::Vector3D > makeLine(const Parameters &pars)
Returns the parsed parameters into an Eigen line parametrization.
Definition: SegmentFitterEventData.cxx:26
MuonValR4::drawLabel
std::unique_ptr< TLatex > drawLabel(const std::string &text, const double xPos, const double yPos, const unsigned int fontSize=18)
Create a TLatex label,.
Definition: VisualizationHelpers.cxx:32
MuonR4::SegmentFitResult::HitVec
std::vector< HitType > HitVec
Definition: SegmentFitterEventData.h:53
MuonR4::MdtSegmentSeedGenerator
Helper class to generate valid seeds for the segment fit.
Definition: MdtSegmentSeedGenerator.h:24
RpcStripContainer.h
MuonR4::MdtSegmentFitter::Config::calibrator
const ISpacePointCalibrator * calibrator
Pointer to the calibrator tool.
Definition: MdtSegmentFitter.h:40
MuonR4::SegmentFitResult::nDoF
int nDoF
degrees of freedom
Definition: SegmentFitterEventData.h:68
MuonR4::SegmentFitHelpers::driftSign
int driftSign(const Amg::Vector3D &posInChamber, const Amg::Vector3D &dirInChamber, const SpacePoint &uncalibHit, MsgStream &msg)
Calculates whether a segement line travereses the tube measurement on the left (-1) or right (1) side...
Definition: SegmentFitHelperFunctions.cxx:239
MuonR4::Segment::MeasType
std::unique_ptr< CalibratedSpacePoint > MeasType
Calibrated space point type.
Definition: MuonSpectrometer/MuonPhaseII/Event/MuonPatternEvent/MuonPatternEvent/Segment.h:23
MuonR4::copy
MuonR4::SegmentFitResult::HitVec copy(const MuonR4::SegmentFitResult::HitVec &hits)
Definition: SegmentFittingAlg.cxx:34
MuonR4::MdtSegmentSeedGenerator::Config
Configuration switches of the module
Definition: MdtSegmentSeedGenerator.h:29
MdtSegmentSeedGenerator.h
SpacePointPerLayerSplitter.h
MuonR4::SegmentFittingAlg::fitSegmentSeed
std::vector< std::unique_ptr< Segment > > fitSegmentSeed(const EventContext &ctx, const ActsGeometryContext &gctx, const SegmentSeed *seed) const
Definition: SegmentFittingAlg.cxx:180
Amg::toString
std::string toString(const Translation3D &translation, int precision=4)
GeoPrimitvesToStringConverter.
Definition: GeoPrimitivesToStringConverter.h:40
MuonR4::SegmentFittingAlg::plugHoles
bool plugHoles(const EventContext &ctx, const ActsGeometryContext &gctx, const SegmentSeed &seed, SegmentFitResult &toRecover) const
Recovery of missed hits.
Definition: SegmentFittingAlg.cxx:322
MuonR4::SegmentFitHelpers::chiSqTermMdt
double chiSqTermMdt(const Amg::Vector3D &posInChamber, const Amg::Vector3D &dirInChamber, const SpacePoint &measurement, MsgStream &msg)
Calculates the chi2 contribuation to a linear segment line from an uncalibrated Mdt measurement.
Definition: SegmentFitHelperFunctions.cxx:52
MuonR4::SegmentFittingAlg::removeOutliers
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.
Definition: SegmentFittingAlg.cxx:260
SG::get
const T * get(const ReadCondHandleKey< T > &key, const EventContext &ctx)
Convenience function to retrieve an object given a ReadCondHandleKey.
Definition: ReadCondHandle.h:287
MuonR4::SegmentFitResult::converged
bool converged
Is the fit converged.
Definition: SegmentFitterEventData.h:76
beamspotman.n
n
Definition: beamspotman.py:729
EL::StatusCode
::StatusCode StatusCode
StatusCode definition for legacy code.
Definition: PhysicsAnalysis/D3PDTools/EventLoop/EventLoop/StatusCode.h:22
UtilFunctions.h
SegmentFitHelperFunctions.h
Amg::Transform3D
Eigen::Affine3d Transform3D
Definition: GeoPrimitives.h:46
Amg::transform
Amg::Vector3D transform(Amg::Vector3D &v, Amg::Transform3D &tr)
Transform a point from a Trasformation3D.
Definition: GeoPrimitivesHelpers.h:156
MuonR4::MdtSegmentFitter::Config
Definition: MdtSegmentFitter.h:26
python.StandardJetMods.pull
pull
Definition: StandardJetMods.py:304
MuonSimHitContainer.h
ATH_CHECK
#define ATH_CHECK
Definition: AthCheckMacros.h:40
MuonR4::SegmentFitResult::segmentPars
Parameters segmentPars
Final segment parameters.
Definition: SegmentFitterEventData.h:60
MuonR4::SegmentFittingAlg::eraseWrongHits
void eraseWrongHits(SegmentFitResult &candidate) const
Removes all hits from the segment which are obvious outliers.
Definition: SegmentFittingAlg.cxx:450
MuonR4::SegmentFit::toString
std::string toString(const Parameters &pars)
Definition: SegmentFitterEventData.cxx:59
MuonR4::removeBeamSpot
bool removeBeamSpot(MuonR4::SegmentFitResult::HitVec &hits)
Definition: SegmentFittingAlg.cxx:50
MuonR4::MdtSegmentSeedGenerator::Config::hitPullCut
double hitPullCut
Upper cut on the hit chi2 w.r.t.
Definition: MdtSegmentSeedGenerator.h:35
DataVector
Derived DataVector<T>.
Definition: DataVector.h:794
ParticleHypothesis.h
ActsGeometryContext
Include the GeoPrimitives which need to be put first.
Definition: ActsGeometryContext.h:27
MuonR4::SegmentFit::ParamDefs::time
@ time
MuonR4::CalibratedSpacePoint::State::Valid
@ Valid
MuonR4::SpacePoint
The muon space point is the combination of two uncalibrated measurements one of them measures the eta...
Definition: MuonSpectrometer/MuonPhaseII/Event/MuonSpacePoint/MuonSpacePoint/SpacePoint.h:19
MuonValR4
Lightweight algorithm to read xAOD MDT sim hits and (fast-digitised) drift circles from SG and fill a...
Definition: IPatternVisualizationTool.h:23
MdtDriftCircleContainer.h
TgcStrip.h
MuonR4::SegmentFit::ParamDefs::y0
@ y0
MuonValR4::IPatternVisualizationTool::PrimitiveVec
std::vector< PrimitivePtr > PrimitiveVec
Definition: IPatternVisualizationTool.h:31
name
std::string name
Definition: Control/AthContainers/Root/debug.cxx:240
plotBeamSpotMon.b
b
Definition: plotBeamSpotMon.py:76
MuonR4::SegmentFit::toInt
constexpr int toInt(const ParamDefs p)
Definition: MuonHoughDefs.h:42
RpcMeasurement.h
WriteCaloSwCorrections.cfg
cfg
Definition: WriteCaloSwCorrections.py:23
Amg::lineDistance
double lineDistance(const AmgVector(N)&posA, const AmgVector(N)&dirA, const AmgVector(N)&posB, const AmgVector(N)&dirB)
: Calculates the shortest distance between two lines
Definition: GeoPrimitivesHelpers.h:308
Amg::Vector3D
Eigen::Matrix< double, 3, 1 > Vector3D
Definition: GeoPrimitives.h:47
python.copyTCTOutput.locDir
locDir
Definition: copyTCTOutput.py:112
MuonR4
This header ties the generic definitions in this package.
Definition: HoughEventData.h:16
PlotSFuncertainty.calib
calib
Definition: PlotSFuncertainty.py:110
SG::WriteHandle
Definition: StoreGate/StoreGate/WriteHandle.h:73
MuonR4::SegmentVec
SegmentAmbiSolver::SegmentVec SegmentVec
Definition: SegmentAmbiSolver.cxx:10
runIDAlign.accumulate
accumulate
Update flags based on parser line args.
Definition: runIDAlign.py:60
a
TList * a
Definition: liststreamerinfos.cxx:10
MdtSegmentFitter.h
ATH_MSG_WARNING
#define ATH_MSG_WARNING(x)
Definition: AthMsgStreamMacros.h:32
MuonGMR4::SpectrometerSector::localToGlobalTrans
const Amg::Transform3D & localToGlobalTrans(const ActsGeometryContext &gctx) const
Returns the local -> global tarnsformation from the sector.
Definition: SpectrometerSector.cxx:75
GeoPrimitivesHelpers.h
MuonR4::SegmentFitResult::timeFit
bool timeFit
Was the time fitted.
Definition: SegmentFitterEventData.h:56
MuonR4::SegmentSeed
Representation of a segment seed (a fully processed hough maximum) produced by the hough transform.
Definition: SegmentSeed.h:14
if
if(febId1==febId2)
Definition: LArRodBlockPhysicsV0.cxx:567
MuonR4::SegmentSeed::parameters
const Parameters & parameters() const
Returns the parameter array.
Definition: SegmentSeed.cxx:35
python.SystemOfUnits.s
float s
Definition: SystemOfUnits.py:147
MuonR4::SegmentFit::localSegmentPars
Parameters localSegmentPars(const xAOD::MuonSegment &seg)
Returns the localSegPars decoration from a xAODMuon::Segment.
Definition: SegmentFitterEventData.cxx:32
python.BuildSignatureFlags.beamSpot
AthConfigFlags beamSpot(AthConfigFlags flags, str instanceName, str recoMode)
Definition: BuildSignatureFlags.py:475
MuonR4::SegmentFitResult
Definition: SegmentFitterEventData.h:46
MuonR4::SegmentFitResult::chi2
double chi2
chi2 of the fit
Definition: SegmentFitterEventData.h:66
drawLine
void drawLine(std::ostream &os)
Definition: PrintMC.cxx:16
xAOD::MdtDriftCircle_v1
https://gitlab.cern.ch/atlas/athena/-/blob/master/MuonSpectrometer/MuonReconstruction/MuonRecEvent/Mu...
Definition: MdtDriftCircle_v1.h:21
MuonR4::AmgSymMatrix
const AmgSymMatrix(2) &SpacePoint
Definition: MuonSpectrometer/MuonPhaseII/Event/MuonSpacePoint/src/SpacePoint.cxx:90
pow
constexpr int pow(int base, int exp) noexcept
Definition: ap_fixedTest.cxx:15
MuonR4::SegmentFittingAlg::execute
virtual StatusCode execute(const EventContext &ctx) const override
Definition: SegmentFittingAlg.cxx:73
xAOD::UncalibMeasType::MdtDriftCircleType
@ MdtDriftCircleType
MuonR4::SegmentFittingAlg::~SegmentFittingAlg
virtual ~SegmentFittingAlg()
MuonR4::SegmentFit::ParamDefs::theta
@ theta
MuonR4::SegmentSeed::msSector
const MuonGMR4::SpectrometerSector * msSector() const
Returns the associated chamber.
Definition: SegmentSeed.cxx:39
python.SystemOfUnits.m
float m
Definition: SystemOfUnits.py:106
generate::Zero
void Zero(TH1D *hin)
Definition: generate.cxx:32
MuonR4::SegmentFitResult::makeLine
std::pair< Amg::Vector3D, Amg::Vector3D > makeLine() const
Returns the defining parameters as a pair of Amg::Vector3D The first part is the position expressed a...
Definition: SegmentFitterEventData.h:83
Generated on Sun Aug 31 2025 21:18:41 for ATLAS Offline Software by   doxygen 1.8.18