ATLAS Offline Software
DCMathSegmentMaker.cxx
Go to the documentation of this file.
1 /*
2  Copyright (C) 2002-2024 CERN for the benefit of the ATLAS collaboration
3 */
4 
5 #include "DCMathSegmentMaker.h"
6 
7 #include <cassert>
8 #include <functional>
9 #include <iostream>
10 
11 #include "AthenaKernel/Timeout.h"
12 #include "EventPrimitives/EventPrimitivesHelpers.h"
13 #include "EventPrimitives/EventPrimitivesToStringConverter.h"
14 #include "FourMomUtils/xAODP4Helpers.h"
15 #include "GeoPrimitives/GeoPrimitivesToStringConverter.h"
16 #include "MuonCompetingRIOsOnTrack/CompetingMuonClustersOnTrack.h"
17 #include "MuonPrepRawData/MdtPrepData.h"
18 #include "MuonPrepRawData/RpcPrepData.h"
19 #include "MuonPrepRawData/TgcPrepData.h"
20 #include "MuonRIO_OnTrack/MdtDriftCircleOnTrack.h"
21 #include "MuonRIO_OnTrack/RpcClusterOnTrack.h"
22 #include "MuonRIO_OnTrack/TgcClusterOnTrack.h"
23 #include "MuonReadoutGeometry/MdtReadoutElement.h"
24 #include "MuonReadoutGeometry/MuonDetectorManager.h"
25 #include "MuonReadoutGeometry/RpcReadoutElement.h"
26 #include "MuonSegment/MuonSegmentQuality.h"
27 #include "TrkDriftCircleMath/ClusterId.h"
28 #include "TrkDriftCircleMath/DriftCircle.h"
29 #include "TrkDriftCircleMath/MdtMultiChamberGeometry.h"
30 #include "TrkDriftCircleMath/ResidualWithSegment.h"
31 #include "TrkDriftCircleMath/Road.h"
32 #include "TrkDriftCircleMath/Segment.h"
33 #include "TrkDriftCircleMath/SegmentFinder.h"
34 #include "TrkEventPrimitives/DefinedParameter.h"
35 #include "TrkEventPrimitives/FitQuality.h"
36 #include "TrkEventPrimitives/JacobianCotThetaPtToThetaP.h"
37 #include "TrkEventPrimitives/LocalDirection.h"
38 #include "TrkEventPrimitives/LocalParameters.h"
39 #include "TrkEventPrimitives/ParamDefs.h"
40 #include "TrkParameters/TrackParameters.h"
41 #include "TrkRIO_OnTrack/RIO_OnTrack.h"
42 #include "TrkRoad/TrackRoad.h"
43 #include "TrkSurfaces/PlaneSurface.h"
44 #include "TrkSurfaces/StraightLineSurface.h"
45 #include "TrkTrack/Track.h"
46 
47 namespace {
48 
49  double cot(double x) {
51  constexpr double dX = std::numeric_limits<float>::epsilon();
52  if (std::abs(x) < dX || std::abs(x - M_PI) < dX) return std::numeric_limits<float>::max();
54  return std::tan(M_PI_2 - x);
55  }
56  template <typename T> constexpr T absmax(const T& a, const T& b) {
57  return std::abs(a) > std::abs(b) ? a : b;}
58 
59 } // namespace
60 namespace Muon {
61 
62  DCMathSegmentMaker::DCMathSegmentMaker(const std::string& t, const std::string& n, const IInterface* p) :
63  AthAlgTool(t, n, p) {
64  declareInterface<IMuonSegmentMaker>(this);
65  }
66 
67  StatusCode DCMathSegmentMaker::initialize() {
68  // retrieve MuonDetectorManager
69  ATH_CHECK(m_DetectorManagerKey.initialize());
70  ATH_CHECK(m_mdtCreator.retrieve());
71  ATH_CHECK(m_mdtCreatorT0.retrieve(DisableTool{m_mdtCreatorT0.empty()}));
72  ATH_CHECK(m_clusterCreator.retrieve());
73  ATH_CHECK(m_compClusterCreator.retrieve());
74  ATH_CHECK(m_idHelperSvc.retrieve());
75  ATH_CHECK(m_printer.retrieve());
76  ATH_CHECK(m_edmHelperSvc.retrieve());
77  ATH_CHECK(m_segmentFinder.retrieve());
78  ATH_CHECK(m_segmentSelectionTool.retrieve());
79  ATH_CHECK(m_segmentFitter.retrieve(DisableTool{!m_refitParameters}));
80  ATH_CHECK(m_dcslFitProvider.retrieve(DisableTool{m_dcslFitProvider.empty()}));
81 
82  // initialise for data handles
83  ATH_CHECK(m_rpcKey.initialize());
84  ATH_CHECK(m_tgcKey.initialize());
85  ATH_CHECK(m_mdtKey.initialize());
86  ATH_CHECK(m_chamberGeoKey.initialize());
87 
88  return StatusCode::SUCCESS;
89  }
90 
91  void DCMathSegmentMaker::find(const Amg::Vector3D& roadpos, const Amg::Vector3D& roaddir,
92  const std::vector<const MdtDriftCircleOnTrack*>& mdts,
93  const std::vector<const MuonClusterOnTrack*>& clusters, bool hasPhiMeasurements,
94  Trk::SegmentCollection* segColl, double momentum, double sinAngleCut) const {
95  const EventContext& ctx = Gaudi::Hive::currentContext();
96  if (m_doTimeOutChecks && Athena::Timeout::instance().reached()) {
97  ATH_MSG_DEBUG("Timeout reached. Aborting sequence.");
98  return;
99  }
100 
101  ATH_MSG_DEBUG("In find, passed " << mdts.size() << " MDTs & "<<clusters.size()<<" clusters");
102 
103  if (mdts.size() < 3) return;
104 
105  const MdtDriftCircleOnTrack* firstRot = findFirstRotInChamberWithMostHits(mdts);
106 
107  if (!firstRot) { return; }
108 
109  const MuonGM::MdtReadoutElement* detEl = firstRot->detectorElement();
110 
111  if (!detEl) {
112  ATH_MSG_WARNING(" no MdtReadoutElement found, returning 0 ");
113  return;
114  }
115 
116  // identifier
117  Identifier chid = firstRot->identify();
118 
119  // endcap or barrel
120  bool isEndcap = m_idHelperSvc->mdtIdHelper().isEndcap(chid);
121 
122  // global to local transformation for chamber
123  Amg::Transform3D gToStation = detEl->GlobalToAmdbLRSTransform();
124 
125  // define axis of chamber in global coordinates
126  Amg::Transform3D amdbToGlobal = detEl->AmdbLRSToGlobalTransform();
127 
128  // transform nominal pointing chamber position into surface frame
129  Amg::Vector3D globalDirCh{detEl->center()};
130  Amg::Vector3D dirCh(gToStation.linear() * globalDirCh);
131  double chamber_angleYZ = std::atan2(dirCh.z(), dirCh.y());
132 
133  Amg::Vector3D roaddir2 = roaddir;
134  double dotprod = globalDirCh.perp() * std::sin(roaddir2.theta()) + globalDirCh.z() * std::cos(roaddir2.theta());
135  if (dotprod < 0) roaddir2 = -roaddir2;
136 
137  // transform the global direction into the surface frame
138  Amg::Vector3D d(gToStation.linear() * roaddir2);
139  // calculate the local road angles in the surface frame
140  double road_angleXZ = std::atan2(d.z(), d.x());
141  double road_angleYZ = std::atan2(d.z(), d.y());
142 
143  if (!hasPhiMeasurements) road_angleXZ = M_PI; // if no phi, take phi perpendicular to plane
144  ATH_MSG_VERBOSE("global road pos "<<Amg::toString(roadpos)<<", global road dir " << Amg::toString(roaddir2) << " XZ " << road_angleXZ << " YZ " << road_angleYZ << " isEndcap "
145  << isEndcap << " central phi " << detEl->center().phi() << " r " << detEl->center().perp()
146  << " z " << detEl->center().z());
147 
148  // rescale errors for low momentum
149  double errorScale = errorScaleFactor(chid, momentum, hasPhiMeasurements);
150 
151  /* ***** create cluster hits ******** */
152  ATH_MSG_DEBUG(" adding clusters " << clusters.size());
153  ClusterVecPair spVecs;
154  if (m_doSpacePoints)
155  spVecs = create2DClusters(clusters);
156  else
157  spVecs = create1DClusters(clusters);
158  TrkDriftCircleMath::CLVec cls = createClusterVec(chid, spVecs.first, gToStation);
159 
160  /* ***** create MDT hits ************ */
161  if (msgLvl(MSG::VERBOSE)) {
162  std::stringstream sstr{};
163  for (const MdtDriftCircleOnTrack* mdt : mdts)
164  sstr<<m_printer->print(*mdt)<<std::endl;
165  ATH_MSG_VERBOSE(" adding mdts " << mdts.size()<<std::endl<<sstr.str());
166  }
167 
168  // set to get Identifiers of chambers with hits
169  std::set<Identifier> chamberSet;
170  double phimin{-9999}, phimax{9999};
171  TrkDriftCircleMath::DCStatistics dcStatistics; // statistics on chamber occupancy
172  TrkDriftCircleMath::DCVec dcs = createDCVec(mdts, errorScale, chamberSet, phimin, phimax, dcStatistics, gToStation, amdbToGlobal);
173 
174  // create geometry
175  std::shared_ptr<const TrkDriftCircleMath::ChamberGeometry> multiGeo;
176  if (m_doGeometry) {
177  ATH_MSG_VERBOSE(" using chamber geometry with #chambers " << chamberSet.size());
178  // vector to store chamber geometries
179  std::vector<TrkDriftCircleMath::MdtChamberGeometry> geos{};
180 
181  // loop over chambers
182  geos.reserve(chamberSet.size());
183  for (const Identifier& id : chamberSet) {
184  ATH_MSG_VERBOSE("Chamber: "<<m_idHelperSvc->toStringChamber(id));
185  geos.push_back(createChamberGeometry(id, gToStation));
186  }
187  // create new geometry
188  multiGeo = std::make_unique<TrkDriftCircleMath::MdtMultiChamberGeometry>(geos);
189 
190 
191  }
192 
193  double angle = m_sinAngleCut;
194  if (sinAngleCut > 0) angle = sinAngleCut;
195  TrkDriftCircleMath::Road road(TrkDriftCircleMath::LocVec2D{0.,0.}, road_angleYZ, chamber_angleYZ, angle);
196 
197  // call segment finder
198  TrkDriftCircleMath::SegVec segs = m_segmentFinder->findSegments(dcs, cls, road, dcStatistics, multiGeo.get());
199 
200  if (msgLvl(MSG::VERBOSE)) {
201  std::stringstream sstr{};
202  unsigned int seg_n{0};
203  for (const TrkDriftCircleMath::Segment& seg: segs) {
204  constexpr double toDeg = 1./Gaudi::Units::degree;
205  sstr<<"Segment number "<<seg_n<<" is at ("<<seg.line().x0()<<","<<seg.line().y0()<<") pointing to "<<seg.line().phi()*toDeg<<" chi2: "<<
206  (seg.chi2()/seg.ndof())<<"("<<seg.ndof()<<")"<<std::endl;
207  sstr<<"Mdt measurements: "<<seg.dcs().size()<<std::endl;
208  for (const TrkDriftCircleMath::DCOnTrack & mdt_meas : seg.dcs()){
209  sstr<<" **** "<<m_printer->print(*mdt_meas.rot());
210  sstr<<" ("<<mdt_meas.state()<<")"<<std::endl;
211  }
212  sstr<<"Cluster measurements "<<seg.clusters().size()<<std::endl;
213  for (const TrkDriftCircleMath::Cluster& clus: seg.clusters()) {
214  sstr<<" ---- "<<m_printer->print(*clus.rot())<<std::endl;
215  }
216  sstr<<std::endl;
217  ++seg_n;
218  }
219  ATH_MSG_VERBOSE("Found " << segs.size() << " segments "<<std::endl<<sstr.str());
220  }
221 
222  // return
223  if (segs.empty()) { return; }
224 
225  // loop over segments
226  segmentCreationInfo sInfo(spVecs, multiGeo.get(), gToStation, amdbToGlobal, phimin, phimax);
227  for (TrkDriftCircleMath::Segment& seg : segs) {
228  std::unique_ptr<MuonSegment> segment = createSegment(ctx, seg, chid, roadpos, roaddir2, mdts, hasPhiMeasurements, sInfo);
229  if (segment) segColl->push_back(segment.release());
230  }
231  ATH_MSG_DEBUG(" Done ");
232  }
233 
234  std::unique_ptr<MuonSegment> DCMathSegmentMaker::createSegment(const EventContext& ctx, TrkDriftCircleMath::Segment& segment, const Identifier& chid,
235  const Amg::Vector3D& roadpos, const Amg::Vector3D& roaddir2,
236  const std::vector<const MdtDriftCircleOnTrack*>& mdts, bool hasPhiMeasurements,
237  segmentCreationInfo& sInfo) const {
238  bool isEndcap = m_idHelperSvc->isEndcap(chid);
239  // find all curved segments
240  MuonStationIndex::ChIndex chIndex = m_idHelperSvc->chamberIndex(chid);
241 
242  static constexpr std::array<Muon::MuonStationIndex::ChIndex ,4> statWithField{MuonStationIndex::BIL, MuonStationIndex::BML,
243  MuonStationIndex::BMS, MuonStationIndex::BOL};
244  const bool isCurvedSegment = segment.hasCurvatureParameters() &&
245  std::find(statWithField.begin(), statWithField.end(), chIndex) != statWithField.end();
246 
247  // remove segments with too few hits
248  if (segment.hitsOnTrack() < 3) return nullptr;
249 
250  // convert segment parameters + x position from road
251  const TrkDriftCircleMath::Line& line = segment.line();
252 
253  ATH_MSG_DEBUG("New segment: chi2 " << segment.chi2() << " ndof " << segment.ndof() << " line " << line.position().x() << ","
254  << line.position().y() << " phi " << line.phi() << " associated clusters "
255  << segment.clusters().size());
256 
257  // local position along x from road
258  Amg::Vector3D lroadpos = sInfo.globalTrans * roadpos;
259  Amg::Vector3D lroaddir = sInfo.globalTrans.linear() * roaddir2;
260 
261  // local x position of first tube used if no phi measurement is present
262  double lxroad = 0.;
263 
264  if (hasPhiMeasurements) {
265  // calculate local position of segment along tube using the road
266  // calculate intersect pattern measurement plane
267  double sphi = 0.;
268  double cphi = lroaddir.x();
269  // swap local y and z in the endcaps
270  if (isEndcap) {
271  sphi = lroaddir.y();
272  lxroad = lroadpos.x() + (-lroadpos.y() + line.position().x()) * cphi / absmax(sphi, std::numeric_limits<double>::min());
273  } else {
274  sphi = lroaddir.z();
275  lxroad = lroadpos.x() + (-lroadpos.z() + line.position().y()) * cphi / absmax(sphi, std::numeric_limits<double>::min());
276  }
277 
278  double shortestTubeLen = 1e9;
279  // loop over hits and get the shortest tube on the segment
280  for (const TrkDriftCircleMath::DCOnTrack& driftCircle : segment.dcs()) {
281  if (driftCircle.state() != TrkDriftCircleMath::DCOnTrack::OnTrack) continue;
282 
283  const MdtDriftCircleOnTrack* riodc{driftCircle.rot()};
284  int lay = m_idHelperSvc->mdtIdHelper().tubeLayer(riodc->identify());
285  int tube = m_idHelperSvc->mdtIdHelper().tube(riodc->identify());
286  double tubelen = 0.5 * riodc->prepRawData()->detectorElement()->getActiveTubeLength(lay, tube);
287  if (tubelen < shortestTubeLen) shortestTubeLen = tubelen;
288  }
289  // if the predicted position lies outside the chamber move it back inside
290  if (std::abs(lxroad) > shortestTubeLen) {
291  ATH_MSG_DEBUG("coordinates far outside chamber! using global position of first hit ");
292  if (lxroad < 0.) shortestTubeLen *= -1.;
293  lxroad = shortestTubeLen;
294  }
295  } else {
296  lxroad = (sInfo.globalTrans * mdts[0]->prepRawData()->detectorElement()->surface(mdts[0]->identify()).center()).x();
297  }
298 
299  // calculate local direction vector
300  Amg::Vector3D lpos(lxroad, line.position().x(), line.position().y());
301 
302  // global position segment
303  Amg::Vector3D gpos = sInfo.amdbTrans * lpos;
304 
305  // create new surface
306  Amg::Transform3D surfaceTransform(sInfo.amdbTrans.rotation());
307  surfaceTransform.pretranslate(gpos);
308  double surfDim = 500.;
309  std::unique_ptr<Trk::PlaneSurface> surf = std::make_unique<Trk::PlaneSurface>(surfaceTransform, surfDim, surfDim);
310 
311  // measurements
312  Amg::Vector2D segLocPos{Amg::Vector2D::Zero()};
313  double linephi = line.phi();
314 
315  // now update the global direction using the local precision angle of the segment and the global phi angle of the
316  // road.
317  Amg::Vector3D gdir = updateDirection(linephi, *surf, roaddir2, isCurvedSegment);
318 
319  // extract RIO_OnTracks
320  std::vector<std::pair<double, std::unique_ptr<const Trk::MeasurementBase>> > rioDistVec; // vector to store the distance of a ROT to the segment
321 
322  // associate MDT hits to segment
323  std::set<Identifier> deltaVec;
324  std::set<Identifier> outoftimeVec;
325 
326  associateMDTsToSegment(gdir, segment, sInfo.geom, sInfo.globalTrans, sInfo.amdbTrans, deltaVec, outoftimeVec, rioDistVec);
327  std::vector<std::pair<double, std::unique_ptr<const Trk::MeasurementBase>>> garbage_collector;
328 
329  TrkDriftCircleMath::DCSLHitSelector hitSelector;
330 
331  if (m_redo2DFit && !isCurvedSegment) {
332  // refit segment after recalibration
333  TrkDriftCircleMath::DCSLFitter defaultFitter;
334  TrkDriftCircleMath::Segment result(TrkDriftCircleMath::Line{0.,0.,0.}, TrkDriftCircleMath::DCOnTrackVec());
335  bool goodFit = defaultFitter.fit(result, line, segment.dcs(), hitSelector.selectHitsOnTrack(segment.dcs()));
336  if (goodFit) {
337  if (std::abs(xAOD::P4Helpers::deltaPhi(segment.line().phi(), result.line().phi())) > 0.01 ||
338  std::abs(segment.line().x0() - result.line().x0()) > 0.01 ||
339  std::abs(segment.line().y0() - result.line().y0()) > 0.01) {
340  // update local position and global
341  linephi = result.line().phi();
342  lpos[1] = result.line().position().x();
343  lpos[2] = result.line().position().y();
344  gpos = sInfo.amdbTrans * lpos;
345 
346  // recreate surface
347  surfaceTransform = Amg::Transform3D(sInfo.amdbTrans.rotation());
348  surfaceTransform.pretranslate(gpos);
349  surf = std::make_unique<Trk::PlaneSurface>(surfaceTransform, surfDim, surfDim);
350 
351  // finally update global direction
352  gdir = updateDirection(linephi, *surf, roaddir2, isCurvedSegment);
353  }
354  }
355  }
356 
357  // create local segment direction
358  Trk::LocalDirection segLocDir;
359  surf->globalToLocalDirection(gdir, segLocDir);
360  if (segLocDir.angleYZ() == 0 && segLocDir.angleXZ() == 0) {
361  ATH_MSG_DEBUG("invalid local direction");
362  return nullptr;
363  }
364 
365  // sanity checks
366  const double diff_phi = xAOD::P4Helpers::deltaPhi(roaddir2.phi(), gdir.phi());
367  const double diff_prec = xAOD::P4Helpers::deltaPhi(linephi, segLocDir.angleYZ());
369  if (std::min(std::abs(diff_phi), std::abs( std::abs(diff_phi) - M_PI)) > 1.e-3 ||
370  std::min(std::abs(diff_prec), std::abs(std::abs(diff_prec) - M_PI)) > 1.e-3) {
371  ATH_MSG_WARNING(" ALARM updated angles wrong: diff phi " << diff_phi << " prec " << diff_prec << " phi rdir " << roaddir2.phi()
372  << " gdir " << gdir.phi() << " lphi " << linephi << " seg "
373  << segLocDir.angleYZ());
374  }
375 
376  // associate Clusters to segment, uses spVecs to get clusters
377  std::pair<std::pair<int, int>, bool> netaPhiHits =
378  associateClustersToSegment(segment, chid, sInfo.globalTrans, sInfo.clusters, sInfo.phimin, sInfo.phimax, rioDistVec);
379 
380  if (rioDistVec.empty()){
381  ATH_MSG_VERBOSE("No measurements were collected.");
382  return nullptr;
383 
384  }
386  auto meas_for_fit = [&rioDistVec] () {
387  std::vector<const Trk::MeasurementBase*> out{};
388  out.reserve(rioDistVec.size());
389  std::sort(rioDistVec.begin(), rioDistVec.end(), SortByDistanceToSegment());
390  for (const std::pair<double, std::unique_ptr<const Trk::MeasurementBase>>& ele : rioDistVec) out.push_back(ele.second.get());
391  return out;
392  };
393 
394 
395  double dlocx{1000.}, dangleXZ{1000.}, qoverp{-99999.}, dqoverp{-99999.};
396  bool hasMeasuredCoordinate = false;
397  if (m_refitParameters && netaPhiHits.second) {
398  ATH_MSG_DEBUG(" distance between first and last phi hit sufficient to perform 4D fit: phi " << gdir.phi() << " theta "
399  << gdir.theta());
400 
401  std::unique_ptr<Trk::Track> track{m_segmentFitter->fit(gpos, gdir, *surf, meas_for_fit())};
402 
403  if (track) {
404  if (isCurvedSegment && track->perigeeParameters() && track->perigeeParameters()->covariance()) {
405  qoverp = track->perigeeParameters()->parameters()[Trk::qOverP];
406  dqoverp = Amg::error(*track->perigeeParameters()->covariance(), Trk::qOverP);
407  }
408  hasMeasuredCoordinate = true;
409  // hack to update the second coordinate errors
410  Amg::MatrixX updatedCov(5, 5);
411  updatedCov.setZero();
412  m_segmentFitter->updateSegmentParameters(*track, *surf, segLocPos, segLocDir, updatedCov);
413  if (Amg::error(updatedCov, Trk::locX) > 0 && Amg::error(updatedCov, Trk::phi) > 0.) {
414  dlocx = Amg::error(updatedCov, Trk::locX);
415  dangleXZ = Amg::error(updatedCov, Trk::phi); // hack (2): phi not always angleXZ
416  } else {
417  ATH_MSG_WARNING(" Corrupt error matrix returned from fitter " << Amg::toString(updatedCov));
418  }
419 
421  surf->localToGlobal(segLocPos, gdir, gpos);
422  surf->localToGlobalDirection(segLocDir, gdir);
423 
424  if (track->measurementsOnTrack() && rioDistVec.size() != track->measurementsOnTrack()->size()) {
425  if (track->measurementsOnTrack()->empty()) {
426  ATH_MSG_DEBUG("No measurements survived");
427  return nullptr;
428  }
429  ATH_MSG_DEBUG(" ROT vector size changed after fit, updating ");
430  garbage_collector = std::move(rioDistVec);
431  rioDistVec.reserve(track->measurementsOnTrack()->size());
432  const Trk::TrackParameters* firstPars = nullptr;
433  for (const Trk::TrackStateOnSurface* tsit : *track->trackStateOnSurfaces()) {
434  const Trk::TrackParameters* pars = tsit->trackParameters();
435  if (!pars) continue;
436  if (!firstPars) firstPars = pars;
437 
438  // check whether state is a measurement, skip outliers if they are not MDT
439  const Trk::MeasurementBase* meas = tsit->measurementOnTrack();
440  if (!meas) continue;
441  if (tsit->type(Trk::TrackStateOnSurface::Outlier) && !dynamic_cast<const MdtDriftCircleOnTrack*>(meas)) continue;
442  double dist = (pars->position() - firstPars->position()).dot(firstPars->momentum().unit());
443  rioDistVec.emplace_back(dist, meas->uniqueClone());
444  }
445  }
446  } else {
447  ATH_MSG_DEBUG(" refit of segment failed!! ");
448  netaPhiHits.second = false;
449  }
450  }
451 
452  // optional update of phi position and direction, only performed if the segment was not refitted and there are phi
453  // hits
454  if (m_updatePhiUsingPhiHits && !netaPhiHits.second) {
455  if (updateSegmentPhi(gpos, gdir, segLocPos, segLocDir, *surf, meas_for_fit(), sInfo.phimin, sInfo.phimax)) {
456  surf->localToGlobal(segLocPos, gpos, gpos);
457  surf->localToGlobalDirection(segLocDir, gdir);
458  hasMeasuredCoordinate = true;
459  dlocx = 100.;
460  dangleXZ = 0.1;
461  }
462  }
463 
464  if (msgLvl(MSG::DEBUG)) {
465  std::vector<const Trk::MeasurementBase*> debug_meas = meas_for_fit();
466  ATH_MSG_DEBUG(" number of hits " << debug_meas.size() << " of which trigger " << netaPhiHits.first.first << " eta and "
467  << netaPhiHits.first.second << " phi ");
468  for (const Trk::MeasurementBase* mit : debug_meas) {
469  const Trk::RIO_OnTrack* rot = dynamic_cast<const Trk::RIO_OnTrack*>(mit);
470  if (rot) {
471  ATH_MSG_DEBUG(m_idHelperSvc->toString(rot->identify()));
472  const MdtDriftCircleOnTrack* mdt = dynamic_cast<const MdtDriftCircleOnTrack*>(rot);
473  if (mdt)
474  ATH_MSG_DEBUG(std::setprecision(4)
475  << " radius " << std::setw(6) << mdt->driftRadius() << " time " << std::setw(6) << mdt->driftTime());
476  continue;
477  }
478  const CompetingMuonClustersOnTrack* crot = dynamic_cast<const CompetingMuonClustersOnTrack*>(mit);
479  if (crot) {
480  ATH_MSG_DEBUG(m_idHelperSvc->toString(crot->rioOnTrack(0).identify())
481  << " comp rot with hits " << crot->containedROTs().size());
482  continue;
483  }
484  ATH_MSG_WARNING("failed to dynamic_cast to ROT ");
485  }
486  }
487  // recalculate holes
488  std::vector<Identifier> holeVec = calculateHoles(ctx, chid, gpos, gdir, hasMeasuredCoordinate, deltaVec, outoftimeVec, rioDistVec);
489 
490  // currently not taking into account masked channels
491  if (!outoftimeVec.empty()) holeVec.insert(holeVec.end(), std::make_move_iterator(outoftimeVec.begin()),
492  std::make_move_iterator(outoftimeVec.end()));
493  MuonSegmentQuality* quality = new MuonSegmentQuality(segment.chi2(), segment.ndof(), holeVec);
494 
495  const TrkDriftCircleMath::DCSLFitter* dcslFitter = m_dcslFitProvider->getFitter();
496  TrkDriftCircleMath::Segment result(TrkDriftCircleMath::Line{0., 0., 0.}, TrkDriftCircleMath::DCOnTrackVec());
497  if (dcslFitter && !segment.hasT0Shift() && m_outputFittedT0) {
498  if (!dcslFitter->fit(result, segment.line(), segment.dcs(), hitSelector.selectHitsOnTrack(segment.dcs()))) {
499  ATH_MSG_DEBUG(" T0 refit failed ");
500  } else {
501  ATH_MSG_DEBUG(" Fitted T0 " << result.t0Shift()<<" is valid "<<result.hasT0Shift());
502  }
503  }
504  bool hasFittedT0 = false;
505  double fittedT0{0}, errorFittedT0{1.};
506  if (m_outputFittedT0 && (segment.hasT0Shift() || (dcslFitter && result.hasT0Shift()))) {
507  hasFittedT0 = true;
508  if (segment.hasT0Shift()) {
509  fittedT0 = segment.t0Shift();
510  errorFittedT0 = segment.t0Error();
511  } else if (dcslFitter && result.hasT0Shift()) {
512  fittedT0 = result.t0Shift();
513  errorFittedT0 = result.t0Error();
514  } else {
515  ATH_MSG_WARNING(" Failed to access fitted t0 ");
516  hasFittedT0 = false;
517  }
518  }
519  // create new segment
520  std::unique_ptr<MuonSegment> msegment;
521  if (isCurvedSegment) { // curved segments
522  if (qoverp == -99999.) {
523  double charge = gpos.z() * std::tan(gdir.theta());
524  charge = charge / std::abs(charge);
525  // if the curved segment was not refit, then use a momentum estimate
526  constexpr double BILALPHA(28.4366), BMLALPHA(62.8267), BMSALPHA(53.1259), BOLALPHA(29.7554);
527  if (chIndex == MuonStationIndex::BIL) {
528  qoverp = (charge * segment.deltaAlpha()) / BILALPHA;
529  dqoverp = M_SQRT2 * segment.dtheta() / BILALPHA;
530  } else if (chIndex == MuonStationIndex::BML) {
531  qoverp = (charge * segment.deltaAlpha()) / BMLALPHA;
532  dqoverp = M_SQRT2 * segment.dtheta() / BMLALPHA;
533  } else if (chIndex == MuonStationIndex::BMS) {
534  qoverp = (charge * segment.deltaAlpha()) / BMSALPHA;
535  dqoverp = M_SQRT2 * segment.dtheta() / BMSALPHA;
536  } else if (chIndex == MuonStationIndex::BOL) {
537  qoverp = (charge * segment.deltaAlpha()) / BOLALPHA;
538  dqoverp = M_SQRT2 * segment.dtheta() / BOLALPHA;
539  }
540  }
541  Amg::MatrixX covMatrix(5, 5);
542  covMatrix.setIdentity();
543  covMatrix(0, 0) = dlocx * dlocx;
544  covMatrix(1, 1) = segment.dy0() * segment.dy0();
545  covMatrix(2, 2) = dangleXZ * dangleXZ;
546  covMatrix(3, 3) = segment.dtheta() * segment.dtheta();
547  covMatrix(4, 4) = dqoverp * dqoverp;
548 
549  std::vector<Trk::DefinedParameter> defPars;
550  defPars.emplace_back(segLocPos[Trk::loc1], Trk::loc1);
551  defPars.emplace_back(segLocPos[Trk::loc2], Trk::loc2);
552  defPars.emplace_back(gdir.phi(), Trk::phi);
553  defPars.emplace_back(gdir.theta(), Trk::theta);
554  defPars.emplace_back(qoverp, Trk::qOverP);
555  Trk::LocalParameters segLocPar(defPars);
556  msegment = std::make_unique<MuonSegment>(
557  std::move(segLocPar),
558  std::move(covMatrix),
559  surf.release(),
560  createROTVec(rioDistVec),
561  quality,
562  Trk::Segment::DCMathSegmentMakerCurved);
563  } else { // straight segments
564  // errors (for now no correlations)
565  Amg::MatrixX covMatrix(4, 4);
566  covMatrix.setIdentity();
567  covMatrix(0, 0) = dlocx * dlocx;
568  covMatrix(1, 1) = segment.dy0() * segment.dy0();
569  covMatrix(2, 2) = dangleXZ * dangleXZ;
570  covMatrix(3, 3) = segment.dtheta() * segment.dtheta();
571  msegment =
572  std::make_unique<MuonSegment>(segLocPos,
573  segLocDir,
574  std::move(covMatrix),
575  surf.release(),
576  createROTVec(rioDistVec),
577  quality,
578  Trk::Segment::DCMathSegmentMaker);
579  }
580 
581  if (hasFittedT0) msegment->setT0Error(fittedT0, errorFittedT0);
582 
583  // check whether segment satisfies minimum quality criteria
584  int segmentQuality = m_segmentSelectionTool->quality(*msegment);
585 
586  if (msgLvl(MSG::DEBUG)) {
587  ATH_MSG_DEBUG(m_printer->print(*msegment) << " quality " << segmentQuality);
588  if (segmentQuality < 0) ATH_MSG_DEBUG(" BAD segment ");
589  if (hasFittedT0) ATH_MSG_DEBUG(" T0 " << fittedT0);
590  if (isCurvedSegment) ATH_MSG_DEBUG(" Curved " << fittedT0);
591  }
592  if (segmentQuality < 0) { return nullptr; }
593  return msegment;
594  }
595 
596  void DCMathSegmentMaker::find(const std::vector<const Trk::RIO_OnTrack*>& rios, Trk::SegmentCollection* segColl) const {
597  std::vector<const MdtDriftCircleOnTrack*> mdts;
598  std::vector<const MuonClusterOnTrack*> clusters;
599 
600  for (const Trk::RIO_OnTrack* it : rios) {
601  Identifier id = it->identify();
602  if (m_idHelperSvc->isMdt(id)) {
603  const MdtDriftCircleOnTrack* mdt = dynamic_cast<const MdtDriftCircleOnTrack*>(it);
604  if (!mdt) { ATH_MSG_WARNING("failed dynamic_cast, not a MDT but hit has MDT id!!!"); }
605  mdts.push_back(mdt);
606  } else if (m_idHelperSvc->isTrigger(id)) {
607  const MuonClusterOnTrack* clus = dynamic_cast<const MuonClusterOnTrack*>(it);
608  if (!clus) { ATH_MSG_WARNING("failed dynamic_cast, not a cluster but hit has RPC/TGC id!!!"); }
609  clusters.push_back(clus);
610  }
611  }
612  find(mdts, clusters, segColl);
613  }
614 
615  void DCMathSegmentMaker::find(const std::vector<const MdtDriftCircleOnTrack*>& mdts,
616  const std::vector<const MuonClusterOnTrack*>& clusters, Trk::SegmentCollection* segColl) const {
617  if (mdts.empty()) return;
618 
619  const MdtDriftCircleOnTrack* mdt = mdts.front();
620  if (!mdt) return;
621 
622  bool hasPhiMeasurements = false;
623  Amg::Vector3D gpos = mdt->globalPosition();
624  Amg::Vector3D gdir = gpos.unit();
625  find(gpos, gdir, mdts, clusters, hasPhiMeasurements, segColl);
626  }
627 
628  void DCMathSegmentMaker::find(const std::vector<const Trk::RIO_OnTrack*>& rios1,
629  const std::vector<const Trk::RIO_OnTrack*>& rios2) const {
630  std::vector<const Trk::RIO_OnTrack*> rios = rios1;
631  rios.insert(rios.end(), rios2.begin(), rios2.end());
632  find(rios);
633  }
634 
635  void DCMathSegmentMaker::find(const Trk::TrackRoad& road, const std::vector<std::vector<const MdtDriftCircleOnTrack*> >& mdts,
636  const std::vector<std::vector<const MuonClusterOnTrack*> >& clusters, Trk::SegmentCollection* segColl,
637  bool hasPhiMeasurements, double momentum) const {
638  // copy all mdt hits into one vector
639  std::vector<const MdtDriftCircleOnTrack*> all_mdts;
640  for (const std::vector<const MdtDriftCircleOnTrack*>& circle_vec : mdts) { std::copy(circle_vec.begin(), circle_vec.end(), std::back_inserter(all_mdts)); }
641 
642  // copy all clusters into one vector
643  std::vector<const MuonClusterOnTrack*> all_clus;
644  for ( const std::vector<const MuonClusterOnTrack*>& clus_vec : clusters) { std::copy(clus_vec.begin(), clus_vec.end(), std::back_inserter(all_clus)); }
645 
646  const Amg::Vector3D& gpos = road.globalPosition();
647  const Amg::Vector3D& gdir = road.globalDirection();
648  find(gpos, gdir, all_mdts, all_clus, hasPhiMeasurements, segColl, momentum, road.deltaEta());
649  }
650 
651  double DCMathSegmentMaker::errorScaleFactor(const Identifier& id, double curvature, bool hasPhiMeasurements) const {
652  double scale = 1.;
653  if (!m_curvedErrorScaling) return scale;
654 
655  if (!errorScalingRegion(id)) return scale;
656 
657  double scaleMax = 5.;
658  if (m_curvedErrorScaling && curvature > 2) {
659  scale = std::min(scaleMax, 1. + curvature / 10000);
660  ATH_MSG_DEBUG(" rescaled errors " << scale << " curvature " << curvature);
661  }
662  scale *= 2;
663 
664  // rescale errors is no phi measurement was found
665  if (!hasPhiMeasurements) {
666  double phiScale = 1.;
667  // rescale errors
668  int stRegion = m_idHelperSvc->mdtIdHelper().stationRegion(id);
669  if (stRegion == 0)
670  phiScale = 2.; // inner
671  else if (stRegion == 1)
672  phiScale = 2.5; // extended
673  else if (stRegion == 2)
674  phiScale = 2.5; // middle
675  else
676  phiScale = 3.; // outer
677  scale = std::sqrt(scale*scale + phiScale*phiScale);
678  ATH_MSG_DEBUG(" rescaled error for missing phi road " << scale);
679  }
680 
681  return scale;
682  }
683 
684  bool DCMathSegmentMaker::errorScalingRegion(const Identifier& id) const
685 
686  {
687  // simple division of MuonSpectrometer in regions using barrel/endcap seperation plus
688  // inner/middle/outer seperation
689 
690  bool isEndcap = m_idHelperSvc->isEndcap(id);
691 
692  if (isEndcap) {
693  std::string stName = m_idHelperSvc->mdtIdHelper().stationNameString(m_idHelperSvc->mdtIdHelper().stationName(id));
694  if (stName[1] == 'I') return true;
695 
696  } else {
697  return true;
698  }
699  return false;
700  }
701 
702  DCMathSegmentMaker::ClusterVecPair DCMathSegmentMaker::create1DClusters(const std::vector<const MuonClusterOnTrack*>& clusters) const {
703  // if empty return
704  if (clusters.empty()) return {};
705  // some useful typedefs...
706 
707  // create a vector to hold the clusters
708  ClusterVec clVec;
709  ClusterVec phiVec;
710  clVec.reserve(clusters.size());
711 
712  for (const MuonClusterOnTrack* clust : clusters) {
713  const Identifier id = clust->identify();
714  const Identifier gasGapId = m_idHelperSvc->gasGapId(id);
715 
716  if (m_idHelperSvc->measuresPhi(id)) {
717  phiVec.push_back(createSpacePoint(gasGapId, nullptr, clust));
718  if (phiVec.back().corrupt()) phiVec.pop_back();
719  } else {
720  clVec.push_back(createSpacePoint(gasGapId, clust, nullptr));
721  if (clVec.back().corrupt()) clVec.pop_back();
722  }
723  }
724 
725  return ClusterVecPair(clVec, phiVec);
726  }
727 
728  DCMathSegmentMaker::ClusterVecPair DCMathSegmentMaker::create2DClusters(const std::vector<const MuonClusterOnTrack*>& clusters) const {
729  // if empty return
730  if (clusters.empty()) return {};
731 
732  ChIdHitMap gasGapHitMap;
733  for (const MuonClusterOnTrack* clus: clusters) {
734  const Identifier id = clus->identify();
735  ATH_MSG_VERBOSE(" new trigger hit " << m_idHelperSvc->toString(id));
736 
737  const Identifier chId = m_idHelperSvc->chamberId(id);
738  const Identifier gasGapId = m_idHelperSvc->gasGapId(id);
739  // eta hits first than phi hits
740  if (!m_idHelperSvc->measuresPhi(id))
741  gasGapHitMap[chId][gasGapId].first.push_back(clus);
742  else
743  gasGapHitMap[chId][gasGapId].second.push_back(clus);
744  }
745 
746  return createSpacePoints(gasGapHitMap);
747  }
748 
749  DCMathSegmentMaker::ClusterVecPair DCMathSegmentMaker::createSpacePoints(const DCMathSegmentMaker::ChIdHitMap& chIdHitMap) const {
750  // vector to store output
751  ClusterVec spacePoints{}, phiVec{};
752  spacePoints.reserve(20);
753  phiVec.reserve(20);
754  // loop over chambers
755  for ( const auto& [id, gasGapHits] : chIdHitMap) {
756  // create clusters per chamber and copy them in to result vector
757  ClusterVecPair cls = createSpacePoints(gasGapHits);
758  std::copy(std::make_move_iterator(cls.first.begin()),
759  std::make_move_iterator(cls.first.end()), std::back_inserter(spacePoints));
760  std::copy(std::make_move_iterator(cls.second.begin()),
761  std::make_move_iterator(cls.second.end()), std::back_inserter(phiVec));
762  }
763 
764  return std::make_pair(std::move(spacePoints), std::move(phiVec));
765  }
766 
767  DCMathSegmentMaker::ClusterVecPair DCMathSegmentMaker::createSpacePoints(const DCMathSegmentMaker::IdHitMap& gasGapHitMap) const {
768  ClusterVec spacePoints{}, phiVec{};
769  bool isEndcap = m_idHelperSvc->isEndcap((*(gasGapHitMap.begin())).first);
770 
771  ATH_MSG_VERBOSE(" creating Space points for " << gasGapHitMap.size() << " gas gaps ");
772 
773  for (const auto& [gasGapId, etaPhiHits] : gasGapHitMap) {
774  // flag whether phi hits are matched with a eta hit
775  std::vector<bool> flagPhihit(etaPhiHits.second.size(), 0);
776 
777  // store Identifier of previous hit to remove duplicates
778  Identifier prevEtaId;
779 
780  ATH_MSG_VERBOSE(" New gasgap " << m_idHelperSvc->toString(gasGapId) << " neta " << etaPhiHits.first.size() << " nphi "
781  << etaPhiHits.second.size());
782 
783  for (const MuonClusterOnTrack* etaHit : etaPhiHits.first) {
784  // check whether we are not dealing with a duplicate hit
785  if (etaHit->identify() == prevEtaId) continue;
786  prevEtaId = etaHit->identify();
787 
788  ATH_MSG_VERBOSE(" Eta hit " << m_idHelperSvc->toString(etaHit->identify()));
789 
790  if (isEndcap) {
791  // check whether match with phi hits was found
792  bool foundSP = false;
793  Identifier prevPhiId;
794  int phi_idx{-1};
795  for (const MuonClusterOnTrack* phiHit : etaPhiHits.second) {
796  // check for duplicate phi hits
797  ++phi_idx;
798  if (phiHit->identify() == prevPhiId) continue;
799  prevPhiId = phiHit->identify();
800 
801  ATH_MSG_VERBOSE(" Phi hit " << m_idHelperSvc->toString(phiHit->identify()));
802 
803  Cluster2D sp = createTgcSpacePoint(gasGapId, etaHit, phiHit);
804  if (sp.corrupt()) continue;
805  spacePoints.push_back(std::move(sp));
806  // mark as used
807  foundSP = true;
808  flagPhihit[phi_idx] = true;
809  }
810 
811  // add single eta hit if not matching phi hit was found
812  if (!foundSP) {
813  Cluster2D sp = createSpacePoint(gasGapId, etaHit, nullptr);
814  if (sp.corrupt()) continue;
815  spacePoints.push_back(std::move(sp));
816  }
817  } else {
818  Cluster2D sp = createRpcSpacePoint(gasGapId, etaHit, etaPhiHits.second);
819  if (sp.corrupt()) continue;
820  // flag all phi hits, not very elegant, but works
821  flagPhihit = std::vector<bool>(etaPhiHits.second.size(), 1);
822  spacePoints.push_back(std::move(sp));
823  }
824  }
825  if (isEndcap) {
826  // loop over flag vector and add unmatched phi hits to phiVec;
827  Identifier prevPhiId;
828  for (unsigned int i = 0; i < flagPhihit.size(); ++i) {
829  if (flagPhihit[i]) continue;
830 
831  // check for duplicate phi hits
832  if (etaPhiHits.second[i]->identify() == prevPhiId) continue;
833  prevPhiId = etaPhiHits.second[i]->identify();
834 
835  Cluster2D sp = createTgcSpacePoint(gasGapId, nullptr, etaPhiHits.second[i]);
836  if (sp.corrupt()) continue;
837  phiVec.push_back(std::move(sp));
838  }
839  } else if (etaPhiHits.first.empty() && !etaPhiHits.second.empty()) {
840  // if there were no eta hits create one phi spacePoint of all phi hits in gasgap
841  Cluster2D sp = createRpcSpacePoint(gasGapId, nullptr, etaPhiHits.second);
842  if (sp.corrupt()) continue;
843  phiVec.push_back(std::move(sp));
844  }
845  }
846 
847  ATH_MSG_VERBOSE(" Creating space points, number of gas-gaps " << gasGapHitMap.size() << " space points " << spacePoints.size());
848 
849  return std::make_pair(std::move(spacePoints), std::move(phiVec));
850  }
851 
852  DCMathSegmentMaker::Cluster2D DCMathSegmentMaker::createSpacePoint(const Identifier& gasGapId, const MuonClusterOnTrack* etaHit,
853  const MuonClusterOnTrack* phiHit) const {
854  bool isEndcap = m_idHelperSvc->isEndcap(gasGapId);
855  double error{1.}, lpx{0.}, lpy{0.};
856  // case one hit missing. Take position and error of the available hit
857  if (!etaHit) {
858  if (!phiHit) {
859  ATH_MSG_WARNING("Both eta and phi hits missing");
860  error = 0;
861  } else {
862  lpx = phiHit->localParameters()[Trk::locX];
863  error = Amg::error(phiHit->localCovariance(), Trk::locX);
864  }
865  } else if (!phiHit) {
866  lpx = etaHit->localParameters()[Trk::locX];
867  error = Amg::error(etaHit->localCovariance(), Trk::locX);
868  } else if (etaHit && phiHit) {
869  if (isEndcap) {
870  return createTgcSpacePoint(gasGapId, etaHit, phiHit);
871  } else {
872  std::vector<const MuonClusterOnTrack*> phiVec{phiHit};
873  return createRpcSpacePoint(gasGapId, etaHit, phiVec);
874  }
875  }
876  Identifier detElId = m_idHelperSvc->detElId(gasGapId);
877  if (std::abs(error) < 0.001) {
878  ATH_MSG_WARNING(" Unphysical error assigned for gasgap " << m_idHelperSvc->toString(gasGapId));
879  error = 0.;
880  }
881  return Cluster2D(detElId, gasGapId, Amg::Vector2D(lpx, lpy), error, etaHit, phiHit);
882  }
883 
884  DCMathSegmentMaker::Cluster2D DCMathSegmentMaker::createTgcSpacePoint(const Identifier& gasGapId,
885  const MuonClusterOnTrack* etaHit,
886  const MuonClusterOnTrack* phiHit) const {
887  double error{1.}, lpx{0.}, lpy{0.};
888  Identifier detElId = m_idHelperSvc->detElId(gasGapId);
889  // case one hit missing. Take position and error of the available hit
890  if (!etaHit) {
891  lpx = phiHit->localParameters()[Trk::locX];
892  error = Amg::error(phiHit->localCovariance(), Trk::locX);
893  } else if (!phiHit) {
894  lpx = etaHit->localParameters()[Trk::locX];
895  error = Amg::error(etaHit->localCovariance(), Trk::locX);
896  } else if (etaHit && phiHit) {
897  // get orientation angle of strip to rotate back from local frame to strip
898  // copy code from ROT creator
899  const MuonGM::TgcReadoutElement* detEl = dynamic_cast<const MuonGM::TgcReadoutElement*>(etaHit->detectorElement());
900  const Amg::Vector3D lSpacePoint = Amg::getRotateZ3D(-90 * Gaudi::Units::deg) * detEl->localSpacePoint(phiHit->identify(),
901  etaHit->globalPosition(),
902  phiHit->globalPosition());
903 
904 
905  lpx = lSpacePoint.x();
906  lpy = lSpacePoint.y();
907  error = Amg::error(etaHit->localCovariance(), Trk::locX);
908  if (error <= std::numeric_limits<double>::epsilon()) {
909  ATH_MSG_WARNING(" Unphysical error assigned for " << m_idHelperSvc->toString(etaHit->identify()));
910  if (etaHit->prepRawData())
911  ATH_MSG_WARNING(" PRD error " << Amg::error(etaHit->prepRawData()->localCovariance(), Trk::locX));
912  }
913  ATH_MSG_DEBUG(" TGC space point: error " << error << " stripWith " << error * M_SQRT2 << std::endl
914  << " " << m_idHelperSvc->toString(etaHit->identify()) << std::endl
915  << " " << m_idHelperSvc->toString(phiHit->identify()));
916  }
917  if (std::abs(error) < 0.001) {
918  ATH_MSG_WARNING(" Unphysical error assigned for gasgap " << m_idHelperSvc->toString(gasGapId));
919  error = 1.;
920  }
921  ATH_MSG_VERBOSE("New space point for "<<m_idHelperSvc->toStringGasGap(gasGapId)<<" at ("<<lpx<<","<<lpy<<")");
922  return Cluster2D(detElId, gasGapId, Amg::Vector2D(lpx, lpy), error, etaHit, phiHit);
923  }
924 
925  DCMathSegmentMaker::Cluster2D DCMathSegmentMaker::createRpcSpacePoint(const Identifier& gasGapId, const MuonClusterOnTrack* etaHit,
926  const std::vector<const MuonClusterOnTrack*>& phiHits) const {
927  // create vector to store phi hits after removal of duplicate hits
928  std::vector<const MuonClusterOnTrack*> cleanPhihits;
929  cleanPhihits.reserve(phiHits.size());
930 
931  double error{1.}, lpx{0.}, lpy{0.};
932  // case one hit missing. Take position and error of the available hit
933  if (!etaHit) {
934  lpx = phiHits.front()->localParameters()[Trk::locX];
935  error = Amg::error(phiHits.front()->localCovariance(), Trk::locX);
936  // loop over phi hits, remove duplicate phi hits
937  Identifier prevId;
938  for (const MuonClusterOnTrack* clus : phiHits) {
939  // remove duplicate phi hits
940  if (clus->identify() == prevId) continue;
941  prevId = clus->identify();
942  cleanPhihits.push_back(clus);
943  }
944  } else if (phiHits.empty()) {
945  lpx = etaHit->localParameters()[Trk::locX];
946  error = Amg::error(etaHit->localCovariance(), Trk::locX);
947  } else if (etaHit && !phiHits.empty()) {
948  lpx = etaHit->localParameters()[Trk::locX];
949  error = Amg::error(etaHit->localCovariance(), Trk::locX);
950 
951  ATH_MSG_DEBUG(" RPC space point: error " << error << " stripWith " << error * M_SQRT2 << std::endl
952  << " " << m_idHelperSvc->toString(etaHit->identify()));
953 
954  double minPos{1e9}, maxPos{-1e9};
955  Identifier prevId;
956 
957  // loop over phi hits, calculate average position + cluster width, remove duplicate phi hits
958  for (const MuonClusterOnTrack* phiHit : phiHits) {
959  // remove duplicate phi hits
960  if (phiHit->identify() == prevId) continue;
961  prevId = phiHit->identify();
962 
963  // calculate phi hit position in local eta hit reference frame
964  Amg::Vector2D phiLocPos{Amg::Vector2D::Zero()};
965  if (etaHit->associatedSurface().globalToLocal(phiHit->globalPosition(), phiHit->globalPosition(), phiLocPos)) {
966  lpy = phiLocPos[Trk::locY];
967  minPos = std::min(minPos, lpy);
968  maxPos = std::max(maxPos, lpy);
969  ATH_MSG_DEBUG(" " << m_idHelperSvc->toString(phiHit->identify()));
970  cleanPhihits.push_back(phiHit);
971  }
972  }
973  if (cleanPhihits.size() > 1)
974  ATH_MSG_DEBUG(" multiple phi hits: nhits " << cleanPhihits.size() << " cl width " << maxPos - minPos);
975  } else {
976  ATH_MSG_DEBUG(" ARRRGGG got two empty pointers!!! ");
977  }
978  Identifier detElId = m_idHelperSvc->detElId(gasGapId);
979  if (std::abs(error) < 0.001) {
980  ATH_MSG_WARNING(" Unphysical error assigned for gasgap " << m_idHelperSvc->toString(gasGapId));
981  error = 1.;
982  }
983  return Cluster2D(detElId, gasGapId, Amg::Vector2D(lpx, lpy), error, etaHit, !cleanPhihits.empty() ? cleanPhihits : phiHits);
984  }
985 
986  TrkDriftCircleMath::CLVec DCMathSegmentMaker::createClusterVec(const Identifier& chid, ClusterVec& spVec,
987  const Amg::Transform3D& gToStation) const {
988  TrkDriftCircleMath::CLVec cls;
989 
990  const int chPhi = m_idHelperSvc->mdtIdHelper().stationPhi(chid);
991 
992  // loop over clusters
993  int index{-1};
994  cls.reserve(spVec.size());
995  for (const Cluster2D& clust : spVec) {
996  ++index;
997  const MuonClusterOnTrack* meas = clust.etaHit ? clust.etaHit : clust.phiHit;
998  // construct cluster id
999  const Identifier id = meas->identify();
1000  const int measuresPhi = m_idHelperSvc->measuresPhi(id);
1001  const int eta = m_idHelperSvc->stationEta(id);
1002  const int phi = m_idHelperSvc->stationPhi(id);
1003  const int isTgc = m_idHelperSvc->isTgc(id);
1004  const int name = isTgc ? m_idHelperSvc->tgcIdHelper().stationName(id) : m_idHelperSvc->rpcIdHelper().stationName(id);
1005  if (!isTgc) {
1006  if (chPhi != phi) {
1007  ATH_MSG_VERBOSE(" Discarding cluster, wrong station phi " << m_idHelperSvc->toString(id));
1008  continue;
1009  }
1010  }
1011  TrkDriftCircleMath::ClusterId clid{name, eta, phi, isTgc, measuresPhi};
1012 
1013  // calculate local cluster position
1014  Amg::Vector3D locPos = gToStation * clust.globalPos;
1015  TrkDriftCircleMath::LocVec2D lp(locPos.y(), locPos.z());
1016 
1017  if (std::abs(lp.y()) > m_maxAssociateClusterDistance) {
1018  ATH_MSG_VERBOSE(" Discarding cluster with large distance from chamber " << m_idHelperSvc->toString(id));
1019  continue;
1020  }
1021  ATH_MSG_VERBOSE(" " << m_idHelperSvc->toString(id) << " clid: " << clid.id() << " central phi "
1022  << meas->detectorElement()->center().phi() << " index " << index);
1023  cls.emplace_back(lp, clust.error, clid, meas, index);
1024  }
1025  return cls;
1026  }
1027 
1028  TrkDriftCircleMath::DCVec DCMathSegmentMaker::createDCVec(const std::vector<const MdtDriftCircleOnTrack*>& mdts, double errorScale,
1029  std::set<Identifier>& chamberSet, double& phimin, double& phimax,
1030  TrkDriftCircleMath::DCStatistics& dcStatistics, const Amg::Transform3D& gToStation,
1031  const Amg::Transform3D& amdbToGlobal) const {
1032  TrkDriftCircleMath::DCVec dcs;
1033  dcs.reserve(mdts.size());
1034  /* ******** Mdt hits ******** */
1035  bool firstMdt = true;
1036 
1037  for (const MdtDriftCircleOnTrack* rot : mdts) {
1038 
1039  Identifier id = rot->identify();
1040  Identifier elId = m_idHelperSvc->mdtIdHelper().elementID(id);
1041 
1042  // calculate local AMDB position
1043  Amg::Vector3D locPos = gToStation * rot->prepRawData()->globalPosition();
1044  TrkDriftCircleMath::LocVec2D lpos(locPos.y(), locPos.z());
1045 
1046  double r = rot->localParameters()[Trk::locR];
1047  double dr = Amg::error(rot->localCovariance(), Trk::locR) * errorScale;
1048 
1049  // create identifier
1050  TrkDriftCircleMath::MdtId mdtid(m_idHelperSvc->mdtIdHelper().isBarrel(id), m_idHelperSvc->mdtIdHelper().multilayer(id) - 1,
1051  m_idHelperSvc->mdtIdHelper().tubeLayer(id) - 1, m_idHelperSvc->mdtIdHelper().tube(id) - 1);
1052 
1053  //
1054  double preciseError = dr;
1055  if (m_usePreciseError) { preciseError = m_preciseErrorScale * (0.23 * std::exp(-std::abs(r) / 6.06) + 0.0362); }
1056  // create new DriftCircle
1057  TrkDriftCircleMath::DriftCircle dc(lpos, r, dr, preciseError, TrkDriftCircleMath::DriftCircle::InTime, mdtid, rot);
1058 
1059  TubeEnds tubeEnds = localTubeEnds(*rot, gToStation, amdbToGlobal);
1060  if (firstMdt) {
1061  phimin = tubeEnds.phimin;
1062  phimax = tubeEnds.phimax;
1063  firstMdt = false;
1064  } else {
1065  updatePhiRanges(tubeEnds.phimin, tubeEnds.phimax, phimin, phimax);
1066  }
1067 
1068  ATH_MSG_VERBOSE(" new MDT hit " << m_idHelperSvc->toString(id) << " x " << lpos.x() << " y " << lpos.y() << " time "
1069  << rot->driftTime() << " r " << r << " dr " << dr << " phi range " << tubeEnds.phimin << " "
1070  << tubeEnds.phimax<<" precise error "<<preciseError);
1071  dcs.push_back(std::move(dc));
1072 
1073  chamberSet.insert(elId);
1074 
1075  ++dcStatistics[rot->prepRawData()->detectorElement()];
1076  }
1077 
1078  return dcs;
1079  }
1080 
1081  TrkDriftCircleMath::MdtChamberGeometry DCMathSegmentMaker::createChamberGeometry(const Identifier& chid,
1082  const Amg::Transform3D& gToStation) const {
1083  /* calculate chamber geometry
1084  it takes as input:
1085  distance between the first and second tube in the chamber within a layer along the tube layer (tube distance)
1086  distance between the first tube in the first layer and the first tube in the second layer along the tube layer
1087  (tube stagering) distance between the first and second layer perpendicular to the tube layers (layer distance)
1088  position of the first hit in ml 0 and ml 1 (2D in plane)
1089  total number of multilayers
1090  total number of layers
1091  total number of tubes per layer for each multilayer
1092  an identifier uniquely identifying the chamber
1093  */
1094 
1095  Amg::Vector3D firstTubeMl0{Amg::Vector3D::Zero()};
1096  Amg::Vector3D firstTubeMl1{Amg::Vector3D::Zero()};
1097 
1098  // get id
1099  int eta = m_idHelperSvc->mdtIdHelper().stationEta(chid);
1100  int phi = m_idHelperSvc->mdtIdHelper().stationPhi(chid);
1101  int name = m_idHelperSvc->mdtIdHelper().stationName(chid);
1102 
1103  SG::ReadCondHandle<MuonGM::MuonDetectorManager> DetectorManagerHandle{m_DetectorManagerKey};
1104  const MuonGM::MuonDetectorManager* MuonDetMgr{*DetectorManagerHandle};
1105  if (!MuonDetMgr) {
1106  ATH_MSG_ERROR("Null pointer to the read MuonDetectorManager conditions object");
1107  return {};
1108  }
1109 
1110  // get detEL for first ml (always there)
1111  const MuonGM::MdtReadoutElement* detEl1 =
1112  MuonDetMgr->getMdtReadoutElement(m_idHelperSvc->mdtIdHelper().channelID(name, eta, phi, 1, 1, 1));
1113  const MuonGM::MdtReadoutElement* detEl2 = nullptr;
1114  int ntube2 = 0;
1115  // number of multilayers in chamber
1116  int nml = detEl1->nMDTinStation();
1117 
1118  // treament of chambers with two ml
1119  if (nml == 2) {
1120  Identifier firstIdml1 = m_idHelperSvc->mdtIdHelper().channelID(name, eta, phi, 2, 1, 1);
1121  detEl2 = MuonDetMgr->getMdtReadoutElement(firstIdml1);
1122  firstTubeMl1 = gToStation * (detEl2->surface(firstIdml1).center());
1123  ntube2 = detEl2->getNtubesperlayer();
1124  }
1125 
1126  // number of layers and tubes
1127  int nlay = detEl1->getNLayers();
1128  int ntube1 = detEl1->getNtubesperlayer();
1129 
1130  // position first tube in ml 0
1131  Identifier firstIdml0 = m_idHelperSvc->mdtIdHelper().channelID(name, eta, phi, 1, 1, 1);
1132  firstTubeMl0 = gToStation * (detEl1->surface(firstIdml0).center());
1133 
1134  // position second tube in ml 0
1135  Identifier secondIdml0 = m_idHelperSvc->mdtIdHelper().channelID(name, eta, phi, 1, 1, 2);
1136  Amg::Vector3D secondTubeMl0 = gToStation * (detEl1->surface(secondIdml0).center());
1137 
1138  TrkDriftCircleMath::LocVec2D firstTube0(firstTubeMl0.y(), firstTubeMl0.z());
1139  TrkDriftCircleMath::LocVec2D firstTube1(firstTubeMl1.y(), firstTubeMl1.z());
1140 
1141  // position first tube ml 0 and 1
1142  Identifier firstIdml0lay1 = m_idHelperSvc->mdtIdHelper().channelID(name, eta, phi, 1, 2, 1);
1143  Amg::Vector3D firstTubeMl0lay1 = gToStation * (detEl1->surface(firstIdml0lay1).center());
1144 
1145  double tubeDist = (secondTubeMl0 - firstTubeMl0).y(); // distance between tube in a given layer
1146  double tubeStage = (firstTubeMl0lay1 - firstTubeMl0).y(); // tube stagering distance
1147  double layDist = (firstTubeMl0lay1 - firstTubeMl0).z(); // distance between layers
1148 
1149  TrkDriftCircleMath::MdtChamberGeometry mdtgeo(chid, m_idHelperSvc.get(), nml, nlay, ntube1, ntube2, firstTube0, firstTube1, tubeDist, tubeStage,
1150  layDist, detEl1->surface().center().theta());
1151 
1152  if (msgLvl(MSG::VERBOSE)) mdtgeo.print(msgStream());
1153 
1154  return mdtgeo;
1155  }
1156 
1157  void DCMathSegmentMaker::associateMDTsToSegment(
1158  const Amg::Vector3D& gdir, TrkDriftCircleMath::Segment& segment,
1159  const TrkDriftCircleMath::ChamberGeometry* multiGeo, const Amg::Transform3D& gToStation, const Amg::Transform3D& amdbToGlobal,
1160  std::set<Identifier>& deltaVec, std::set<Identifier>& outoftimeVec,
1161  std::vector<std::pair<double, std::unique_ptr<const Trk::MeasurementBase>> >& rioDistVec) const {
1162  // clear result vectors
1163 
1164  // convert segment parameters + x position from road
1165  const TrkDriftCircleMath::Line& line = segment.line();
1166  TrkDriftCircleMath::TransformToLine toLineml1(line);
1167  TrkDriftCircleMath::TransformToLine toLineml2(line);
1168  if (segment.hasCurvatureParameters()) {
1169  // ml2 segment direction
1170  double ml2phi = line.phi() - segment.deltaAlpha();
1171  TrkDriftCircleMath::LocVec2D ml2dir(std::cos(ml2phi), std::sin(ml2phi));
1172  // ml2 segment position
1173  const TrkDriftCircleMath::LocVec2D ml1LocPos = multiGeo->tubePosition(0, multiGeo->nlay(), 0);
1174  const TrkDriftCircleMath::LocVec2D ml2LocPos = multiGeo->tubePosition(1, 1, 0);
1175  double chamberMidPtY = (ml1LocPos.y() + ml2LocPos.y()) / 2.0;
1176  TrkDriftCircleMath::LocVec2D ml2pos(segment.deltab(), chamberMidPtY);
1177  // construct the new ml2 segment line & transform
1178  const TrkDriftCircleMath::Line ml2line(ml2pos, ml2dir);
1179  TrkDriftCircleMath::TransformToLine tmptoLine(ml2line);
1180  // set the ml2 line
1181  toLineml2 = tmptoLine;
1182  }
1183 
1184  for (TrkDriftCircleMath::DCOnTrack& dcit : segment.dcs()) {
1185  if (dcit.state() == TrkDriftCircleMath::DCOnTrack::Delta) { deltaVec.insert(dcit.rot()->identify()); }
1186 
1187  if (dcit.state() == TrkDriftCircleMath::DCOnTrack::OutOfTime) { outoftimeVec.insert(dcit.rot()->identify()); }
1188 
1189  if (dcit.state() != TrkDriftCircleMath::DCOnTrack::OnTrack) continue;
1190 
1191  const MdtDriftCircleOnTrack* riodc{dcit.rot()};
1192 
1193  // choose which line to use (ml1 or ml2)
1194  TrkDriftCircleMath::TransformToLine toLine = toLineml1;
1195  if (m_idHelperSvc->mdtIdHelper().multilayer(riodc->identify()) == 2) toLine = toLineml2;
1196  // calculate position of hit in line frame
1197  TrkDriftCircleMath::LocVec2D pointOnHit = toLine.toLine(dcit.position());
1198 
1199  // calculate position of hit on line in line frame
1200  TrkDriftCircleMath::LocVec2D pointOnLine(pointOnHit.x(), 0.);
1201 
1202  // transform back to local AMDB coordinates
1203  TrkDriftCircleMath::LocVec2D pointOnLineAMDB = toLine.toLocal(pointOnLine);
1204 
1205  // get position along wire from ROT
1206  Amg::Vector3D posAlong = gToStation * riodc->globalPosition();
1207 
1208  // set yz components
1209  posAlong[1] = pointOnLineAMDB.x();
1210  posAlong[2] = pointOnLineAMDB.y();
1211 
1212  // back to global
1213  Amg::Vector3D mdtGP = amdbToGlobal * posAlong;
1214 
1215  const Trk::StraightLineSurface* surf = dynamic_cast<const Trk::StraightLineSurface*>(&riodc->associatedSurface());
1216  if (!surf) {
1217  ATH_MSG_WARNING(" dynamic cast to StraightLineSurface failed for mdt!!! ");
1218  continue;
1219  }
1220 
1221  // calculate Amg::Vector2D using surf to obtain sign
1222  Amg::Vector2D locPos{Amg::Vector2D::Zero()};
1223  if (!surf->globalToLocal(mdtGP, gdir, locPos)) ATH_MSG_WARNING(" globalToLocal failed ");
1224 
1225  // calculate side
1226  Trk::DriftCircleSide side = locPos[Trk::driftRadius] < 0 ? Trk::LEFT : Trk::RIGHT;
1227 
1228  std::unique_ptr<MdtDriftCircleOnTrack> nonconstDC;
1229  bool hasT0 = segment.hasT0Shift();
1230  if (!hasT0 || m_mdtCreatorT0.empty()) {
1231  // ATH_MSG_VERBOSE(" recalibrate MDT hit");
1232  nonconstDC.reset(m_mdtCreator->createRIO_OnTrack(*riodc->prepRawData(), mdtGP, &gdir));
1233  if (hasT0) ATH_MSG_WARNING("Attempted to change t0 without a properly configured MDT creator tool. ");
1234  } else {
1235  ATH_MSG_VERBOSE(" recalibrate MDT hit with shift " << segment.t0Shift()<<" "<<m_printer->print(*riodc));
1236  nonconstDC.reset(m_mdtCreatorT0->createRIO_OnTrack(*riodc->prepRawData(), mdtGP, &gdir, segment.t0Shift()));
1237  }
1238 
1239  if (!nonconstDC) {
1240  dcit.state(TrkDriftCircleMath::DCOnTrack::OutOfTime);
1241  continue;
1242  }
1243  ATH_MSG_VERBOSE("Post calibrated hit "<<m_printer->print(*nonconstDC));
1244 
1245  // update the drift radius after recalibration, keep error
1246  TrkDriftCircleMath::DriftCircle new_dc(dcit.position(), std::abs(nonconstDC->driftRadius()), dcit.dr(), dcit.drPrecise(),
1247  dcit.driftState(), dcit.id(),
1248  nonconstDC.get());
1249  TrkDriftCircleMath::DCOnTrack new_dc_on_track(new_dc, dcit.residual(), dcit.errorTrack());
1250  dcit = std::move(new_dc_on_track);
1251 
1252  if (hasT0) {
1253  if (msgLvl(MSG::VERBOSE)) {
1254  double shift = riodc->driftTime() - nonconstDC->driftTime();
1255  ATH_MSG_VERBOSE(" t0 shift " << segment.t0Shift() << " from hit " << shift << " recal " << nonconstDC->driftRadius()
1256  << " t " << nonconstDC->driftTime() << " from fit " << dcit.r() << " old "
1257  << riodc->driftRadius() << " t " << riodc->driftTime());
1258  if (std::abs(std::abs(nonconstDC->driftRadius()) - std::abs(dcit.r())) > 0.1 && nonconstDC->driftRadius() < 19. &&
1259  nonconstDC->driftRadius() > 1.) {
1260  ATH_MSG_WARNING("Detected invalid recalibration after T0 shift");
1261  }
1262  }
1263  }
1264  m_mdtCreator->updateSign(*nonconstDC, side);
1265  double dist = pointOnHit.x();
1266  rioDistVec.emplace_back(dist, std::move(nonconstDC));
1267  }
1268  }
1269 
1270  template <class T> struct IdDataVec {
1271  using Entry = T;
1272  using EntryVec = std::vector<Entry>;
1273 
1274  IdDataVec() = default;
1275  explicit IdDataVec(const Identifier& i) : id(i) {}
1276 
1277  Identifier id{0};
1278  EntryVec data{};
1279  };
1280 
1281  template <class T> struct SortIdDataVec {
1282  bool operator()(const IdDataVec<T>& d1, const IdDataVec<T>& d2) { return d1.id < d2.id; }
1283  };
1284 
1285  struct SortClByPull {
1286  bool operator()(const std::pair<double, DCMathSegmentMaker::Cluster2D>& d1,
1287  const std::pair<double, DCMathSegmentMaker::Cluster2D>& d2) {
1288  return std::abs(d1.first) < std::abs(d2.first);
1289  }
1290  };
1291 
1292  std::pair<std::pair<int, int>, bool> DCMathSegmentMaker::associateClustersToSegment(
1293  const TrkDriftCircleMath::Segment& segment, const Identifier& chid, const Amg::Transform3D& gToStation, ClusterVecPair& spVecs,
1294  double phimin, double phimax, std::vector<std::pair<double, std::unique_ptr<const Trk::MeasurementBase>> >& rioDistVec) const {
1295  typedef IdDataVec<std::pair<double, Cluster2D> > GasGapData;
1296  typedef IdDataVec<GasGapData> ChamberData;
1297  typedef std::vector<ChamberData> ChamberDataVec;
1298  ChamberDataVec chamberDataVec;
1299  bool isEndcap = m_idHelperSvc->isEndcap(chid);
1300 
1301  // keep track of the number of eta/phi hits on the segment
1302  bool refit = false;
1303  std::pair<std::pair<int, int>, bool> netaPhiHits(std::make_pair(0, 0), false);
1304  if (segment.clusters().empty()) return netaPhiHits;
1305 
1306  std::vector<const Trk::MeasurementBase*> phiHits;
1307 
1308  // only refit if there are sufficient phi hits and no multiple phi hits per surface
1309  refit = true;
1310 
1311  // keep track of detector elements which space points added to the track
1312  std::set<Identifier> detElOnSegments;
1313  std::set<MuonStationIndex::PhiIndex> phiIndices;
1314 
1315  ATH_MSG_DEBUG(" Associating clusters: " << segment.clusters().size() << " number of space points " << spVecs.first.size());
1316 
1317  // associate space points and sort them per detector element and gas gap
1318  for (const TrkDriftCircleMath::Cluster& clust : segment.clusters()) {
1319  ATH_MSG_VERBOSE(" accessing cluster: " << clust.index());
1320  const Cluster2D& spacePoint = spVecs.first[clust.index()];
1321 
1322  // skip corrupt space points
1323  if (spacePoint.corrupt()) {
1324  ATH_MSG_DEBUG(" Found corrupt space point: index " << clust.index());
1325  continue;
1326  }
1327  // reject TGC clusters that are not 2D
1328  if (m_reject1DTgcSpacePoints && !spacePoint.is2D() && m_idHelperSvc->isTgc(spacePoint.identify())) {
1329  ATH_MSG_DEBUG(" Rejecting 1D tgc space point " << m_idHelperSvc->toString(spacePoint.identify()));
1330  continue;
1331  }
1332  if (m_assumePointingPhi && spacePoint.is2D() && !checkPhiConsistency(spacePoint.globalPos.phi(), phimin, phimax)) {
1333  ATH_MSG_DEBUG(" Inconsistent phi angle, dropping space point: phi " << spacePoint.globalPos.phi() << " range " << phimin
1334  << " " << phimax);
1335  continue;
1336  }
1337 
1338  std::pair<double, double> resPull = residualAndPullWithSegment(segment, spacePoint, gToStation);
1339 
1340  // if empty or new chamber, add chamber
1341  if (chamberDataVec.empty() || chamberDataVec.back().id != spacePoint.detElId) {
1342  detElOnSegments.insert(spacePoint.detElId);
1343  chamberDataVec.emplace_back(spacePoint.detElId);
1344  MuonStationIndex::PhiIndex phiIndex = m_idHelperSvc->phiIndex(spacePoint.detElId);
1345  phiIndices.insert(phiIndex);
1346  }
1347 
1348  // reference to current chamber data
1349  ChamberData& chamber = chamberDataVec.back();
1350 
1351  // if same detector element
1352  if (spacePoint.detElId == chamber.id) {
1353  // if chamber empty or new gas gap, add gasp gap
1354  if (chamber.data.empty() || chamber.data.back().id != spacePoint.gasGapId) {
1355  chamber.data.emplace_back(spacePoint.gasGapId);
1356  }
1357  }
1358 
1359  // reference to current gas gap data
1360  GasGapData& gasGap = chamber.data.back();
1361  gasGap.data.emplace_back(resPull.second, spacePoint);
1362  }
1363 
1364  // calculate the distance between the first and last station, use r in barrel and z in endcaps
1365  double posFirstPhiStation{FLT_MAX}, posLastPhiStation{0.};
1366 
1367  // loop over chambers and create competing ROTs per chamber
1368  for (ChamberData& chamb : chamberDataVec) {
1369  // select best clusters per gas gap in chamber
1370  std::list<const Trk::PrepRawData*> etaClusterVec{}, phiClusterVec{};
1371  std::set<Identifier> etaIds;
1372  // loop over gas gaps
1373  for (GasGapData& gasGap : chamb.data) {
1374  // sort space points by their pull with the segment
1375  std::sort(gasGap.data.begin(), gasGap.data.end(), SortClByPull());
1376 
1377  // select all space points with a pull that is within 1 of the best pull
1378  double bestPull = std::abs(gasGap.data.front().first);
1379 
1380  // count number of associated clusters in gas gap
1381  unsigned int nassociatedSp = 0;
1382  GasGapData::EntryVec::const_iterator cl_it = gasGap.data.begin();
1383  while (cl_it != gasGap.data.end() && std::abs(cl_it->first) - bestPull < 1.) {
1384  const Cluster2D& sp = cl_it->second;
1386  double dist = distanceToSegment(segment, sp.globalPos, gToStation);
1387  ATH_MSG_VERBOSE(" selected space point: " << m_idHelperSvc->toString(sp.identify()) << " pull "
1388  << std::abs(cl_it->first) << " distance to segment " << dist << " phi "
1389  << sp.globalPos.phi());
1390 
1391  // here keep open the option not to create CompetingMuonClustersOnTrack
1392  if (sp.etaHit) {
1393  if (!etaIds.count(sp.etaHit->identify())) {
1394  etaIds.insert(sp.etaHit->identify());
1395 
1396  if (m_createCompetingROTsEta)
1397  etaClusterVec.push_back(sp.etaHit->prepRawData());
1398  else {
1399  rioDistVec.emplace_back(dist, sp.etaHit->uniqueClone());
1400  ++netaPhiHits.first.first;
1401  }
1402  }
1403  }
1404  if (!sp.phiHits.empty()) {
1405  if (m_createCompetingROTsPhi) {
1406  // can have multiple phi hits per cluster, loop over phi hits and add them
1407  std::transform(sp.phiHits.begin(), sp.phiHits.end(), std::back_inserter(phiClusterVec),
1408  [](const Muon::MuonClusterOnTrack* clus){
1409  return clus->prepRawData();
1410  });
1411  } else {
1412  // can have multiple phi hits per cluster, loop over phi hits and add them
1413  for (const MuonClusterOnTrack* phi_hit : sp.phiHits) {
1414  rioDistVec.emplace_back(dist, phi_hit->uniqueClone());
1415  ++netaPhiHits.first.second;
1416  phiHits.push_back(phi_hit);
1417 
1418  // calculate position
1419  double phiPos = isEndcap ? std::abs(phi_hit->globalPosition().z()) : phi_hit->globalPosition().perp();
1420  posFirstPhiStation = std::min(phiPos, posFirstPhiStation);
1421  posLastPhiStation = std::max(phiPos, posLastPhiStation);
1422  }
1423  if (sp.phiHits.size() > 1) refit = false;
1424  }
1425  }
1426  ++nassociatedSp;
1427  ++cl_it;
1428  }
1429  // multiple clusters in same gas gap, don't refit
1430  if (!m_createCompetingROTsPhi && nassociatedSp > 1) refit = false;
1431  }
1432 
1433  if (m_createCompetingROTsEta) {
1434  // create competing ROT for eta hits
1435  if (!etaClusterVec.empty()) {
1436  std::unique_ptr<const CompetingMuonClustersOnTrack> etaCompCluster = m_compClusterCreator->createBroadCluster(etaClusterVec, 0.);
1437  if (!etaCompCluster) {
1438  ATH_MSG_DEBUG(" failed to create competing ETA ROT " << etaClusterVec.size());
1439  } else {
1440  double dist = distanceToSegment(segment, etaCompCluster->globalPosition(), gToStation);
1441  ++netaPhiHits.first.first;
1442  if (msgLvl(MSG::VERBOSE)) {
1443  ATH_MSG_VERBOSE(" selected cluster: " << m_idHelperSvc->toString(etaClusterVec.front()->identify()));
1444  for (unsigned int i = 0; i < etaCompCluster->containedROTs().size(); ++i) {
1445  ATH_MSG_VERBOSE(
1446  " content: " << m_idHelperSvc->toString(etaCompCluster->containedROTs()[i]->identify()));
1447  }
1448  }
1449  rioDistVec.emplace_back(dist, std::move(etaCompCluster));
1450  }
1451  }
1452  }
1453 
1454  if (m_createCompetingROTsPhi) {
1455  // create competing ROT for phi hits
1456  if (!phiClusterVec.empty()) {
1457  std::unique_ptr<const CompetingMuonClustersOnTrack> phiCompCluster = m_compClusterCreator->createBroadCluster(phiClusterVec, 0.);
1458  if (!phiCompCluster) {
1459  ATH_MSG_DEBUG(" failed to create competing PHI ROT " << phiClusterVec.size());
1460  } else {
1461  double dist = distanceToSegment(segment, phiCompCluster->globalPosition(), gToStation);
1462  phiHits.push_back(phiCompCluster.get());
1463 
1464  ++netaPhiHits.first.second;
1465 
1466  if (msgLvl(MSG::VERBOSE)) {
1467  ATH_MSG_VERBOSE(" selected cluster: " << m_idHelperSvc->toString(phiClusterVec.front()->identify()));
1468  for (unsigned int i = 0; i < phiCompCluster->containedROTs().size(); ++i) {
1469  ATH_MSG_VERBOSE(
1470  " content: " << m_idHelperSvc->toString(phiCompCluster->containedROTs()[i]->identify()));
1471  }
1472  }
1473 
1474 
1475  // calculate position
1476  double phiPos = isEndcap ? std::abs(phiCompCluster->globalPosition().z()) :
1477  phiCompCluster->globalPosition().perp();
1478  posFirstPhiStation = std::min(phiPos,posFirstPhiStation);
1479  posLastPhiStation = std::max(phiPos,posLastPhiStation);
1480  rioDistVec.emplace_back(dist, std::move(phiCompCluster));
1481 
1482  }
1483  }
1484  }
1485  }
1486 
1487  // add phi hits that were not associated with an eta hit (only in barrel)
1488  if ((!spVecs.second.empty() || m_recoverBadRpcCabling) && m_addUnassociatedPhiHits && !isEndcap) {
1489  TrkDriftCircleMath::ResidualWithSegment resWithSegment(segment);
1490 
1491  std::map<Identifier, std::list<const Trk::PrepRawData*> > phiClusterMap;
1492 
1493  std::set<const MuonClusterOnTrack*> selectedClusters;
1494  std::vector<const Cluster2D*> phiClusters;
1495  phiClusters.reserve(spVecs.second.size());
1496 
1497  // create lists of PrepRawData per detector element
1498  for (const Cluster2D& phi_clus :spVecs.second) {
1499  if (!phi_clus.phiHit || phi_clus.corrupt()) {
1500  ATH_MSG_WARNING(" phi clusters without phi hit!!");
1501  continue;
1502  }
1503  phiClusters.push_back(&phi_clus);
1504  selectedClusters.insert(phi_clus.phiHit);
1505  }
1506 
1507  unsigned int recoveredUnassociatedPhiHits(0);
1508  if (m_recoverBadRpcCabling) {
1509  // now loop over 2D space points and add the phi hits to the list if the detEl is not yet added to the
1510  // segment
1511  for (const Cluster2D& rpc_clust : spVecs.first) {
1512  // skip clusters without phi hit
1513  if (!rpc_clust.phiHit || rpc_clust.corrupt()) continue;
1514 
1515  // skip clusters in detector elements that are already associated (ok as this is only done for RPC)
1516  if (detElOnSegments.count(rpc_clust.detElId)) continue;
1517 
1518  MuonStationIndex::PhiIndex phiIndex = m_idHelperSvc->phiIndex(rpc_clust.detElId);
1519  // skip clusters in detector layer
1520  if (phiIndices.count(phiIndex)) continue;
1521 
1522  bool wasFound = false;
1523  for (const MuonClusterOnTrack* phi_hit : rpc_clust.phiHits) {
1524  // now to avoid duplicate also skip if the given ROT is already in the list
1525  if (!selectedClusters.insert(phi_hit).second) {
1526  // flag as found
1527  wasFound = true;
1528 
1529  // just because I'm paranoid, remove the hits from this cluster that were already inserted up to
1530  // here
1531  for (const MuonClusterOnTrack* erase_me : rpc_clust.phiHits) {
1532  if (erase_me == phi_hit) break;
1533  selectedClusters.erase(erase_me);
1534  }
1535  break;
1536  }
1537  }
1538  if (wasFound) continue;
1539 
1540  // if we get here we should add the hit
1541  phiClusters.push_back(&rpc_clust);
1542  ++recoveredUnassociatedPhiHits;
1543  }
1544  }
1545 
1546  unsigned int addedPhiHits(0);
1547  for (const Cluster2D* phi_clus : phiClusters) {
1548  const Identifier& detElId = phi_clus->detElId;
1549 
1550  // check that detector element is not already added to segment
1551  if (detElOnSegments.count(detElId)) continue;
1552 
1553  MuonStationIndex::PhiIndex phiIndex = m_idHelperSvc->phiIndex(detElId);
1554  // skip clusters in detector layer
1555  if (phiIndices.count(phiIndex)) continue;
1556 
1557  // calculate local cluster position
1558  Amg::Vector3D locPos = gToStation * phi_clus->globalPos;
1559 
1560  // calculate intersect of segment with cluster
1561  TrkDriftCircleMath::Cluster cl(TrkDriftCircleMath::LocVec2D(locPos.y(), locPos.z()), 1.);
1562  double residual = resWithSegment.residual(cl);
1563  double segError = std::sqrt(resWithSegment.trackError2(cl));
1564  const MuonGM::RpcReadoutElement* detEl = dynamic_cast<const MuonGM::RpcReadoutElement*>(phi_clus->phiHit->detectorElement());
1565  if (!detEl) {
1566  ATH_MSG_WARNING("found RPC prd without RpcReadoutElement");
1567  continue;
1568  }
1569 
1570  // perform bound check
1571  double stripLength = detEl->StripLength(1);
1572  bool inBounds = std::abs(residual) < 0.5 * stripLength + 2. + segError;
1573  if (msgLvl(MSG::DEBUG)) {
1574  ATH_MSG_DEBUG(" Unassociated " << m_idHelperSvc->toString(phi_clus->phiHit->identify()) << " pos x " << cl.position().x()
1575  << " pos y " << cl.position().y() << " : residual " << residual << " strip half length "
1576  << 0.5 * stripLength << " segment error " << segError);
1577  if (inBounds)
1578  ATH_MSG_DEBUG(" inBounds");
1579  else
1580  ATH_MSG_DEBUG(" outBounds");
1581  }
1582  if (inBounds) {
1583  // can have multiple phi hits per cluster, loop over phi hits and add them
1584  std::list<const Trk::PrepRawData*>& cham_hits{phiClusterMap[detElId]};
1585  std::transform(phi_clus->phiHits.begin(), phi_clus->phiHits.end(), std::back_inserter(cham_hits),
1586  [](const MuonClusterOnTrack* clus){
1587  return clus->prepRawData();
1588  });
1589  }
1590  }
1591 
1592  // loop over detector elements and created competing ROTs
1593  for (const auto& [phi_id, prds] : phiClusterMap) {
1594  if (prds.empty()) {
1595  ATH_MSG_WARNING(" chamber without phi hits " << m_idHelperSvc->toString(phi_id));
1596  continue;
1597  }
1598 
1599  std::unique_ptr<const CompetingMuonClustersOnTrack> phiCompCluster = m_compClusterCreator->createBroadCluster(prds, 0.);
1600  if (!phiCompCluster) {
1601  ATH_MSG_DEBUG(" failed to create competing PHI ROT " << prds.size());
1602  } else {
1603  double dist = distanceToSegment(segment, phiCompCluster->globalPosition(), gToStation);
1604 
1605  if (std::abs(dist) > m_maxAssociateClusterDistance) {
1606 
1607  ATH_MSG_VERBOSE(" rejected unassociated cluster: " << m_idHelperSvc->toString(prds.front()->identify())
1608  << " distance to segment " << dist);
1609  continue;
1610  }
1611  phiHits.push_back(phiCompCluster.get());
1612  ++netaPhiHits.first.second;
1613  ++addedPhiHits;
1614  if (msgLvl(MSG::VERBOSE)) {
1615  ATH_MSG_VERBOSE(" selected unassociated cluster: " << m_idHelperSvc->toString(prds.front()->identify())
1616  << " distance to segment " << dist);
1617  for (unsigned int i = 0; i < phiCompCluster->containedROTs().size(); ++i) {
1618  ATH_MSG_VERBOSE(
1619  " content: " << m_idHelperSvc->toString(phiCompCluster->containedROTs()[i]->identify()));
1620  }
1621  }
1622  rioDistVec.emplace_back(dist, std::move(phiCompCluster));
1623  }
1624  }
1625  ATH_MSG_VERBOSE("Added " << addedPhiHits << " unass phi hits out of " << spVecs.second.size()
1626  << " phi hits without eta hit and " << recoveredUnassociatedPhiHits << " with unassociated eta hit ");
1627  }
1628 
1629  // calculate distance between first and last phi trigger hit, refit if there is a good constraint on phi
1630  double phiDistanceMax = posLastPhiStation - posFirstPhiStation;
1631  if (isEndcap && phiDistanceMax < 1000.)
1632  refit = false;
1633  else if (phiDistanceMax < 400.)
1634  refit = false;
1635 
1636  netaPhiHits.second = refit;
1637  return netaPhiHits;
1638  }
1639 
1640  double DCMathSegmentMaker::distanceToSegment(const TrkDriftCircleMath::Segment& segment, const Amg::Vector3D& hitPos,
1641  const Amg::Transform3D& gToStation) {
1642  const TrkDriftCircleMath::Line& line = segment.line();
1643  TrkDriftCircleMath::TransformToLine toLine(line);
1644  double cos_sinLine = cot(line.phi());
1645 
1646  // calculate local AMDB position
1647  Amg::Vector3D locPos = gToStation * hitPos;
1648 
1649  TrkDriftCircleMath::LocVec2D lpos(locPos.y(), locPos.z());
1650 
1651  // calculate distance of segment to measurement surface
1652  double delta_y = lpos.y() - line.position().y();
1653 
1654  // calculate position of hit in line frame
1655  TrkDriftCircleMath::LocVec2D lineSurfaceIntersect(delta_y * cos_sinLine + line.position().x(), lpos.y());
1656 
1657  // calculate position of hit in line frame
1658  TrkDriftCircleMath::LocVec2D pointOnHit = toLine.toLine(lineSurfaceIntersect);
1659 
1660  return pointOnHit.x();
1661  }
1662 
1663  DataVector<const Trk::MeasurementBase> DCMathSegmentMaker::createROTVec(
1664  std::vector<std::pair<double, std::unique_ptr<const Trk::MeasurementBase>> >& rioDistVec) {
1665  // sort hits according to they distance to the segment position
1666  std::sort(rioDistVec.begin(), rioDistVec.end(), SortByDistanceToSegment());
1667 
1668  DataVector<const Trk::MeasurementBase> rioVec{SG::OWN_ELEMENTS};
1669  rioVec.reserve(rioDistVec.size());
1670  for (std::pair<double, std::unique_ptr<const Trk::MeasurementBase>>& rdit : rioDistVec) { rioVec.push_back(std::move(rdit.second)); }
1671  rioDistVec.clear();
1672  return rioVec;
1673  }
1674 
1675  std::pair<double, double> DCMathSegmentMaker::residualAndPullWithSegment(const TrkDriftCircleMath::Segment& segment,
1676  const Cluster2D& spacePoint,
1677  const Amg::Transform3D& gToStation) {
1678  const TrkDriftCircleMath::Line& line = segment.line();
1679  double cos_sinLine = cot(line.phi());
1680 
1681  // calculate sp postion in AMDB reference frame
1682  Amg::Vector3D locPos = gToStation * spacePoint.globalPos;
1683  TrkDriftCircleMath::LocVec2D lpos(locPos.y(), locPos.z());
1684 
1685  // calculate distance of segment to measurement surface
1686  double delta_y = lpos.y() - line.position().y();
1687 
1688  // calculate position of hit in line frame
1689  TrkDriftCircleMath::LocVec2D lineSurfaceIntersect(delta_y * cos_sinLine + line.position().x(), lpos.y());
1690 
1691  // calculate position of hit in line frame
1692  double residual = lpos.x() - lineSurfaceIntersect.x();
1693  double pull = residual / spacePoint.error;
1694  return std::make_pair(residual, pull);
1695  }
1696 
1697  std::vector<Identifier> DCMathSegmentMaker::calculateHoles(const EventContext& ctx,
1698  Identifier chid, const Amg::Vector3D& gpos, const Amg::Vector3D& gdir, bool hasMeasuredCoordinate, std::set<Identifier>& deltaVec,
1699  std::set<Identifier>& outoftimeVec, const std::vector<std::pair<double, std::unique_ptr<const Trk::MeasurementBase>> >& rioDistVec) const {
1700  // calculate crossed tubes
1701  SG::ReadCondHandle<Muon::MuonIntersectGeoData> InterSectSvc{m_chamberGeoKey, ctx};
1702  if (!InterSectSvc.isValid()) {
1703  ATH_MSG_ERROR("Null pointer to the read MuonDetectorManager conditions object");
1704  return {};
1705  }
1706  const MuonStationIntersect intersect = InterSectSvc->tubesCrossedByTrack(chid, gpos, gdir);
1707  const MuonGM::MuonDetectorManager* MuonDetMgr = InterSectSvc->detMgr();
1708 
1709  // set to identify the hit on the segment
1710  std::set<Identifier> hitsOnSegment, chambersOnSegment;
1711  int firstLayer{-1}, lastLayer{-1};
1712  for (const std::pair<double, std::unique_ptr<const Trk::MeasurementBase>>& rdit : rioDistVec) {
1713  const MdtDriftCircleOnTrack* mdt = dynamic_cast<const MdtDriftCircleOnTrack*>(rdit.second.get());
1714  if (!mdt) continue;
1715  const Identifier id = mdt->identify();
1716  int layer = (m_idHelperSvc->mdtIdHelper().tubeLayer(id) - 1) + 4 * (m_idHelperSvc->mdtIdHelper().multilayer(id) - 1);
1717  if (firstLayer == -1)
1718  firstLayer = layer;
1719  else
1720  lastLayer = layer;
1721 
1722  hitsOnSegment.insert(id);
1723  chambersOnSegment.insert(m_idHelperSvc->chamberId(id));
1724  }
1725 
1726  // cross check for cosmic case
1727  if (firstLayer > lastLayer) { std::swap(firstLayer, lastLayer); }
1728  ATH_MSG_VERBOSE(" Tube layer ranges: " << firstLayer << " -- " << lastLayer << " crossed tubes "
1729  << intersect.tubeIntersects().size());
1730  // clear hole vector
1731  std::vector<Identifier> holeVec;
1732  for (const MuonTubeIntersect& tint : intersect.tubeIntersects()) {
1733  if (!chambersOnSegment.count(m_idHelperSvc->chamberId(tint.tubeId))) {
1734  ATH_MSG_VERBOSE(" chamber not on segment, not counting tube " << tint.rIntersect << " l " << tint.xIntersect << " "
1735  << m_idHelperSvc->toString(tint.tubeId));
1736  continue;
1737  }
1738 
1739  const Identifier& id = tint.tubeId;
1740  int layer = (m_idHelperSvc->mdtIdHelper().tubeLayer(id) - 1) + 4 * (m_idHelperSvc->mdtIdHelper().multilayer(id) - 1);
1741 
1742  bool notBetweenHits = layer < firstLayer || layer > lastLayer;
1743  double distanceCut = hasMeasuredCoordinate ? -20 : -200.;
1744  double innerRadius = MuonDetMgr->getMdtReadoutElement(id)->innerTubeRadius();
1745  if (notBetweenHits && (std::abs(tint.rIntersect) > innerRadius || (!m_allMdtHoles && tint.xIntersect > distanceCut))) {
1746  ATH_MSG_VERBOSE(" not counting tube: distance to wire " << tint.rIntersect << " dist to tube end " << tint.xIntersect
1747  << " " << m_idHelperSvc->toString(tint.tubeId));
1748  continue;
1749  }
1750  // check whether there is a hit in this tube
1751  if (hitsOnSegment.count(tint.tubeId)) {
1752  ATH_MSG_VERBOSE(" tube on segment: distance to wire " << tint.rIntersect << " dist to tube end " << tint.xIntersect
1753  << " " << m_idHelperSvc->toString(tint.tubeId));
1754  continue;
1755  }
1756  // check whether there is a delta electron in this tube
1757  if (m_removeDeltas) {
1758  if (deltaVec.count(tint.tubeId)) {
1759  ATH_MSG_VERBOSE(" removing delta, distance to wire " << tint.rIntersect << " dist to tube end " << tint.xIntersect
1760  << " " << m_idHelperSvc->toString(tint.tubeId));
1761  continue;
1762  }
1763 
1764  const MdtPrepData* prd = findMdt(ctx, id);
1765  if (prd && std::abs(prd->localPosition()[Trk::locR]) < std::abs(tint.rIntersect)) {
1766  ATH_MSG_VERBOSE(" found and removed delta, distance to wire " << tint.rIntersect << " dist to tube end "
1767  << tint.xIntersect << " "
1768  << m_idHelperSvc->toString(tint.tubeId));
1769  continue;
1770  }
1771  }
1772  ATH_MSG_VERBOSE((outoftimeVec.count(tint.tubeId) ? "Out-of-time" : "hole") << " distance to wire "
1773  << tint.rIntersect << " dist to tube end " << tint.xIntersect << " "
1774  << m_idHelperSvc->toString(tint.tubeId)<<(notBetweenHits ? "outside hits" : "between hits"));
1775 
1776  holeVec.push_back(tint.tubeId);
1777  }
1778  return holeVec;
1779  }
1780 
1781  const MdtPrepData* DCMathSegmentMaker::findMdt(const EventContext& ctx, const Identifier& id) const {
1782  IdentifierHash colHash;
1783  if (m_idHelperSvc->mdtIdHelper().get_module_hash(m_idHelperSvc->chamberId(id), colHash)){
1784  ATH_MSG_VERBOSE("Invalid Mdt identifier "<<m_idHelperSvc->toString(id));
1785  return nullptr;
1786  }
1787  SG::ReadHandle<Muon::MdtPrepDataContainer> MdtCont{m_mdtKey, ctx};
1788  if (!MdtCont.isValid()) {
1789  ATH_MSG_WARNING("Cannot retrieve MdtPrepDataContainer ");
1790  return nullptr;
1791  }
1792  const MdtPrepDataCollection* collptr = MdtCont->indexFindPtr(colHash);
1793  if (!collptr) return nullptr;
1794  for (const MdtPrepData* prd : *collptr) {
1795  if (prd->identify() == id) return prd;
1796  }
1797  return nullptr;
1798  }
1799 
1800  const MdtDriftCircleOnTrack* DCMathSegmentMaker::findFirstRotInChamberWithMostHits(
1801  const std::vector<const MdtDriftCircleOnTrack*>& mdts) const {
1802  int hitsInChamberWithMostHits = 0;
1803  std::map<Identifier, int> hitsPerChamber;
1804  int currentSector = -1;
1805  const MdtDriftCircleOnTrack* rotInChamberWithMostHits = nullptr;
1806 
1807  // loop over all MDTs and count number of MDTs per chamber
1808  for (const MdtDriftCircleOnTrack* rot : mdts) {
1809  if (!rot) {
1810  ATH_MSG_WARNING(" rot not a MdtDriftCircleOnTrack ");
1811  continue;
1812  }
1813  Identifier chId = m_idHelperSvc->chamberId(rot->identify());
1814  int sector = m_idHelperSvc->sector(chId);
1815  if (currentSector == -1) {
1816  currentSector = sector;
1817  } else if (sector != currentSector) {
1818  return nullptr;
1819  }
1820  int& hitsInCh = hitsPerChamber[chId];
1821  ++hitsInCh;
1822  if (hitsInCh > hitsInChamberWithMostHits) {
1823  hitsInChamberWithMostHits = hitsInCh;
1824  rotInChamberWithMostHits = rot;
1825  }
1826  }
1827  return rotInChamberWithMostHits;
1828  }
1829 
1830  bool DCMathSegmentMaker::checkBoundsInXZ(double xline, double zline, double dXdZ,
1831  const std::vector<DCMathSegmentMaker::HitInXZ>& hits) const {
1832  bool ok = true;
1833 
1834  // look over hits and check whether all are in bounds
1835  for (const HitInXZ& hit : hits) {
1836  bool outBounds = false;
1837  double locExX = xline + dXdZ * (hit.z - zline);
1838  if (hit.isMdt && (locExX < hit.xmin - 1. || locExX > hit.xmax + 1.)) {
1839  ok = false;
1840  outBounds = true;
1841  if (!msgLvl(MSG::DEBUG)) break;
1842  }
1843 
1844  if (outBounds && msgLvl(MSG::DEBUG)) {
1845  ATH_MSG_DEBUG(" " << std::setw(65) << m_idHelperSvc->toString(hit.id) << " pos (" << std::setw(6) << (int)hit.x << ","
1846  << std::setw(6) << (int)hit.z << ") ex pos " << std::setw(6) << (int)locExX << " min " << std::setw(6)
1847  << (int)hit.xmin << " max " << std::setw(6) << (int)hit.xmax << " phimin " << std::setw(6)
1848  << hit.phimin << " phimax " << std::setw(6) << hit.phimax << " outBounds, cross-check");
1849  }
1850  }
1851  return ok;
1852  }
1853 
1854  bool DCMathSegmentMaker::updateSegmentPhi(const Amg::Vector3D& gpos, const Amg::Vector3D& gdir, Amg::Vector2D& segLocPos,
1855  Trk::LocalDirection& segLocDir, Trk::PlaneSurface& surf,
1856  const std::vector<const Trk::MeasurementBase*>& rots, double seg_phimin,
1857  double seg_phimax) const {
1858  bool hasUpdated = false;
1859 
1860  const Amg::Transform3D& segmentToGlobal = surf.transform();
1861  Amg::Transform3D gToSegment = surf.transform().inverse();
1862  Amg::Vector3D ldir = gToSegment * gdir;
1863 
1864  // ensure that we can extrapolate
1865  if (ldir.z() < 0.0001) return false;
1866 
1867  double dXdZ = ldir.x() / ldir.z();
1868  Amg::Vector3D lsegPos = gToSegment * gpos;
1869  double xline = lsegPos.x();
1870  double zline = lsegPos.z();
1871  ATH_MSG_VERBOSE(" Associated hits " << rots.size() << " angleXZ " << 90. * segLocDir.angleXZ() / (M_PI_2) << " dXdZ " << dXdZ
1872  << " seg Pos (" << xline << " " << zline << ") " << segLocPos);
1873 
1874  std::vector<HitInXZ> hits;
1875  hits.reserve(rots.size());
1876 
1877  unsigned int nphiHits(0);
1878  const HitInXZ* firstPhiHit{nullptr}, *lastPhiHit{nullptr};
1879 
1880  for (const Trk::MeasurementBase* meas : rots) {
1881  Identifier id = m_edmHelperSvc->getIdentifier(*meas);
1882  if (!id.is_valid()) continue;
1883  Amg::Vector3D lpos{Amg::Vector3D::Zero()};
1884  double lxmin{0}, lxmax{0}, phimin{0}, phimax{0};
1885  bool isMdt = m_idHelperSvc->isMdt(id);
1886  bool measuresPhi = m_idHelperSvc->measuresPhi(id);
1887  if (isMdt) {
1888  lpos.setZero();
1889  const MdtDriftCircleOnTrack* mdt = static_cast<const MdtDriftCircleOnTrack*>(meas);
1890  TubeEnds tubeEnds = localTubeEnds(*mdt, gToSegment, segmentToGlobal);
1891 
1892  lxmin = tubeEnds.lxmin;
1893  lxmax = tubeEnds.lxmax;
1894  phimin = tubeEnds.phimin;
1895  phimax = tubeEnds.phimax;
1896  } else {
1897  lpos = gToSegment * meas->globalPosition();
1898  lxmin = lpos.x() - 5 * Amg::error(meas->localCovariance(), Trk::locX);
1899  lxmax = lpos.x() + 5 * Amg::error(meas->localCovariance(), Trk::locX);
1900 
1901  const CompetingMuonClustersOnTrack* crot = dynamic_cast<const CompetingMuonClustersOnTrack*>(meas);
1902  if (!measuresPhi) {
1903  if (crot) {
1904  const MuonGM::RpcReadoutElement* detEl =
1905  dynamic_cast<const MuonGM::RpcReadoutElement*>(crot->containedROTs().front()->prepRawData()->detectorElement());
1906  if (detEl) {
1907  // perform bound check
1908  double stripLength = detEl->StripLength(0);
1909  lxmin = lpos.x() - 0.5 * stripLength;
1910  lxmax = lpos.x() + 0.5 * stripLength;
1911  }
1912  }
1913  Amg::Vector3D locPosition = lpos;
1914  locPosition[0] = lxmin;
1915  Amg::Vector3D globalPos = segmentToGlobal * locPosition;
1916  double phi1 = globalPos.phi();
1917 
1918  locPosition[0] = lxmax;
1919  globalPos = segmentToGlobal * locPosition;
1920  double phi2 = globalPos.phi();
1921  phimin = std::min(phi1, phi2);
1922  phimax = std::max(phi1, phi2);
1923 
1924  } else {
1925  if (m_idHelperSvc->isTgc(id)) {
1926  // need some special tricks for TGC phi hits as their reference plane can be rotated
1927  // with respect to the MDT frame
1928 
1929  // get orientation angle of strip to rotate back from local frame to strip
1930  // copy code from ROT creator
1931  int stripNo = m_idHelperSvc->tgcIdHelper().channel(id);
1932  int gasGap = m_idHelperSvc->tgcIdHelper().gasGap(id);
1933  if (!crot) {
1934  ATH_MSG_WARNING("dynamic cast failed for CompetingMuonClustersOnTrack");
1935  continue;
1936  }
1937  auto detEl = dynamic_cast<const MuonGM::TgcReadoutElement*>(crot->containedROTs().front()->prepRawData()->detectorElement());
1938 
1939  // transform the two points inth
1940  const Amg::Vector3D segFrame_StripDir = gToSegment.linear()* detEl->stripDir(gasGap, stripNo);
1941  const Amg::Vector3D segFrame_stripPos = gToSegment * detEl->channelPos(id);
1942 
1943  lpos = segFrame_stripPos +
1944  Amg::intersect<3>(lsegPos, ldir, segFrame_stripPos, segFrame_StripDir).value_or(0) * segFrame_StripDir;
1945 
1946  ATH_MSG_VERBOSE(" In seg frame: phi pos " << Amg::toString(lsegPos)
1947  << " shifted pos " << Amg::toString(segFrame_StripDir)
1948  << " intersect with segment " << Amg::toString(lpos));
1949  }
1950  Amg::Vector3D globalPos = segmentToGlobal * lpos;
1951  phimin = globalPos.phi();
1952  phimax = phimin;
1953 
1954  // check whether phi is consistent with segment phi range
1955  bool phiOk = checkPhiConsistency(phimin, seg_phimin, seg_phimax);
1956  if (!phiOk) {
1957  ATH_MSG_DEBUG(" Inconsistent phi " << phimin << " range " << seg_phimin << " " << seg_phimax);
1958  }
1959  }
1960  }
1961 
1962  hits.emplace_back(id, isMdt, measuresPhi, lpos.x(), lpos.z(), lxmin, lxmax, phimin, phimax);
1963  if (measuresPhi) {
1964  ++nphiHits;
1965  if (!firstPhiHit)
1966  firstPhiHit = &hits.back();
1967  else {
1968  double distPhiHits = std::abs(firstPhiHit->z - hits.back().z);
1969  if (distPhiHits > 500.) {
1970  lastPhiHit = &hits.back();
1971  } else {
1972  // not count this phi hit
1973  --nphiHits;
1974  ATH_MSG_DEBUG(" close phi hits, distance " << distPhiHits);
1975  }
1976  }
1977  }
1978  if (msgLvl(MSG::VERBOSE)) {
1979  double locExX = xline + dXdZ * (lpos.z() - zline);
1980  ATH_MSG_VERBOSE(" " << std::setw(65) << m_idHelperSvc->toString(id) << " pos (" << std::setw(6) << (int)lpos.x() << ","
1981  << std::setw(6) << (int)lpos.z() << ") ex pos " << std::setw(6) << (int)locExX << " min "
1982  << std::setw(6) << (int)lxmin << " max " << std::setw(6) << (int)lxmax << " phimin " << std::setw(6)
1983  << phimin << " phimax " << std::setw(6) << phimax);
1984  if (lpos.x() < lxmin || lpos.x() > lxmax) ATH_MSG_VERBOSE(" outBounds");
1985  }
1986  }
1987 
1988  if (nphiHits == 1) {
1989  if (!firstPhiHit) {
1990  ATH_MSG_WARNING(" Pointer to first phi hit not set, this should not happen! ");
1991  } else {
1992  if (xline != firstPhiHit->x) {
1993  hasUpdated = true;
1994 
1995  // use phi position of the phi hit
1996  xline = firstPhiHit->x;
1997  zline = firstPhiHit->z;
1998 
1999  if (m_assumePointingPhi) {
2000  Amg::Vector3D ipLocPos = gToSegment.translation();
2001  ATH_MSG_VERBOSE(" IP position in local frame " << ipLocPos);
2002 
2003  double dz = ipLocPos.z() - zline;
2004  if (std::abs(dz) > 0.001) {
2005  ATH_MSG_VERBOSE(" hit (" << xline << "," << zline << ") IP (" << ipLocPos.x() << "," << ipLocPos.z()
2006  << ") dXdZ " << (ipLocPos.x() - xline) / dz << " old " << dXdZ);
2007  dXdZ = (ipLocPos.x() - xline) / dz;
2008  }
2009  }
2010  }
2011  }
2012  } else if (nphiHits == 2) {
2013  if (!firstPhiHit || !lastPhiHit) {
2014  ATH_MSG_WARNING(" Pointer to one of the two phi hit not set, this should not happen! ");
2015  } else {
2016  double dz = lastPhiHit->z - firstPhiHit->z;
2017  // use phi position of the first hit
2018  xline = firstPhiHit->x;
2019  zline = firstPhiHit->z;
2020  if (std::abs(dz) > 300.) {
2021  double dx = lastPhiHit->x - firstPhiHit->x;
2022  hasUpdated = true;
2023 
2024  // if the two hits are far enough apart, also use the direction of the line connecting the two hits.
2025  dXdZ = dx / dz;
2026  }
2027  }
2028  } else {
2029  // in all other cases just rotate until the MDTs are ok
2030  }
2031 
2032  if (hasUpdated) {
2033  // move segment to position of phi hit
2034  double segX = xline - dXdZ * zline;
2035 
2036  // finally check whether now everything is in bounds
2037  bool ok = checkBoundsInXZ(segX, 0., dXdZ, hits);
2038  if (!ok) {
2039  // give WARNING and give up for now
2040  ATH_MSG_DEBUG("still several out of bounds hits after rotation: posx(" << segX << ") dXdZ " << dXdZ
2041  << " keeping old result ");
2042  }
2043 
2044  // update segment parameters
2045  double alphaYZ = segLocDir.angleYZ();
2046  double alphaXZ = std::atan2(1, dXdZ);
2047 
2048  segLocPos[Trk::locX] = segX;
2049  segLocDir = Trk::LocalDirection(alphaXZ, alphaYZ);
2050  }
2051  return hasUpdated;
2052  }
2053 
2054  DCMathSegmentMaker::TubeEnds DCMathSegmentMaker::localTubeEnds(const MdtDriftCircleOnTrack& mdt, const Amg::Transform3D& gToSegment,
2055  const Amg::Transform3D& segmentToG) const {
2056  TubeEnds tubeEnds;
2057  const Identifier& id = mdt.identify();
2058  Amg::Vector3D lpos = gToSegment * mdt.prepRawData()->globalPosition();
2059 
2060  // use readout and hv side as the surface frame is not that of the chamber
2061  Amg::Vector3D lropos = gToSegment * mdt.prepRawData()->detectorElement()->ROPos(id);
2062  Amg::Vector3D lhvpos = lpos + (lpos - lropos);
2063 
2064  // rescale to correctly take into account active tube length
2065  double tubeLen = (lropos - lhvpos).mag();
2066  double activeTubeLen =
2067  mdt.detectorElement()->getActiveTubeLength(m_idHelperSvc->mdtIdHelper().tubeLayer(id), m_idHelperSvc->mdtIdHelper().tube(id));
2068  double scaleFactor = activeTubeLen / tubeLen;
2069  lropos[0] = scaleFactor * lropos.x();
2070  lhvpos[0] = scaleFactor * lhvpos.x();
2071 
2072  tubeEnds.lxmin = std::min(lropos.x(), lhvpos.x());
2073  tubeEnds.lxmax = std::max(lropos.x(), lhvpos.x());
2074 
2075  Amg::Vector3D ropos = segmentToG * lropos;
2076  Amg::Vector3D hvpos = segmentToG * lhvpos;
2077  const double phiRO = ropos.phi();
2078  const double phiHV = hvpos.phi();
2079  tubeEnds.phimin = std::min(phiRO, phiHV);
2080  tubeEnds.phimax = std::max(phiRO, phiHV);
2081  return tubeEnds;
2082  }
2083 
2084  void DCMathSegmentMaker::updatePhiRanges(double phiminNew, double phimaxNew, double& phiminRef, double& phimaxRef) {
2085  // check whether we are at the boundary where phi changes sign
2086  if (phiminRef * phimaxRef < 0.) {
2087  if (phiminRef < -1.1) {
2088  if (phiminRef > phiminNew) phiminRef = phiminNew;
2089  if (phimaxRef < phimaxNew) phimaxRef = phimaxNew;
2090  } else {
2091  if (phiminRef < phiminNew) phiminRef = phiminNew;
2092  if (phimaxRef > phimaxNew) phimaxRef = phimaxNew;
2093  }
2094  } else {
2095  // if not life is 'easy'
2096  if (phiminRef < 0.) {
2097  if (phiminRef < phiminNew) phiminRef = phiminNew;
2098  if (phimaxRef > phimaxNew) phimaxRef = phimaxNew;
2099  } else {
2100  if (phiminRef > phiminNew) phiminRef = phiminNew;
2101  if (phimaxRef < phimaxNew) phimaxRef = phimaxNew;
2102  }
2103  }
2104  }
2105 
2106  bool DCMathSegmentMaker::checkPhiConsistency(double phi, double phimin, double phimax) const {
2107  // only if assuming pointing phi
2108  if (!m_assumePointingPhi) return true;
2109 
2110  bool phiOk = true;
2111  double offset = 0.05;
2112  if (phimin * phimax < 0.) {
2113  if (phi < 0.) {
2114  if (phi > -1.1) {
2115  if (phi < phimin - offset) phiOk = false;
2116  } else {
2117  if (phi > phimin + offset) phiOk = false;
2118  }
2119  } else {
2120  if (phi > 1.1) {
2121  if (phi < phimax - offset) phiOk = false;
2122  } else {
2123  if (phi > phimax + offset) phiOk = false;
2124  }
2125  }
2126  } else {
2127  if (phi < phimin - offset || phi > phimax + offset) phiOk = false;
2128  }
2129  return phiOk;
2130  }
2131 
2132  Amg::Vector3D DCMathSegmentMaker::updateDirection(double linephi, const Trk::PlaneSurface& surf, const Amg::Vector3D& roaddir,
2133  bool isCurvedSegment) const {
2134  // Local direction along precision measurement (0, dy, dz)
2135  Trk::LocalDirection segLocDirs(M_PI_2, linephi);
2136  Amg::Vector3D gdirs;
2137  surf.localToGlobalDirection(segLocDirs, gdirs);
2138  // Local direction in plane (1,0,0)
2139  Trk::LocalDirection segLocDiro(0., M_PI_2);
2140  Amg::Vector3D gdiro;
2141  surf.localToGlobalDirection(segLocDiro, gdiro);
2142 
2143  // recalculate the value of the local XZ angle for the give YZ angle of the segment such that the global phi
2144  // direction remains unchanged
2145  double dx = std::sin(gdirs.theta()) * std::cos(gdirs.phi());
2146  double dy = std::sin(gdirs.theta()) * std::sin(gdirs.phi());
2147  double dz = std::cos(gdirs.theta());
2148 
2149  // vector gdiro
2150 
2151  double dxo = std::sin(gdiro.theta()) * std::cos(gdiro.phi());
2152  double dyo = std::sin(gdiro.theta()) * std::sin(gdiro.phi());
2153  double dzo = std::cos(gdiro.theta());
2154 
2155  // solve system real direction = A * gdir + B * gdiro
2156  // phi global constraint: (1)*sin(phi road) - (2)*cos(phi road) = 0
2157  // ( A * dx + B * dxo ) sin (phi ) - (A * dy + B *dyo ) cos (phi) = 0
2158  // A ( dx sin - dy cos ) + B (dx0 sin -dy0 cos) = A * a0 + B * b0 = 0
2159  // psi = atan (-b0 , a0)
2160 
2161  double a0 = dx * std::sin(roaddir.phi()) - dy * std::cos(roaddir.phi());
2162  double b0 = dxo * std::sin(roaddir.phi()) - dyo * std::cos(roaddir.phi());
2163  if (b0 < 1e-8 && b0 > 0) b0 = 1e-8;
2164  if (b0 > -1e-8 && b0 < 0) b0 = -1e-8;
2165  double dxn = dx - a0 * dxo / b0;
2166  double dyn = dy - a0 * dyo / b0;
2167  double dzn = dz - a0 * dzo / b0;
2168  double norm = std::sqrt(dxn * dxn + dyn * dyn + dzn * dzn);
2169 
2170  // flip the sign if the direction NOT parallel to road
2171  if (m_assumePointingPhi) {
2172  if (dxn * roaddir.x() + dyn * roaddir.y() + dzn * roaddir.z() < 0.) { norm = -norm; }
2173  } else {
2174  if (dxn * roaddir.x() + dyn * roaddir.y() < 0.) { norm = -norm; }
2175  }
2176 
2177  if (isCurvedSegment) norm = norm / 2.;
2178 
2179  //
2180  // Follow segment fit direction
2181  //
2182  dxn = dxn / norm;
2183  dyn = dyn / norm;
2184  dzn = dzn / norm;
2185 
2186  return Amg::Vector3D(dxn, dyn, dzn);
2187  }
2188 } // namespace Muon
MdtMultiChamberGeometry.h
Muon::MuonSegmentQuality
Definition: MuonSegmentQuality.h:34
Trk::TrackRoad
Encapsulates the information required by the find() method of the muon segment makers.
Definition: TrackRoad.h:21
Trk::RIO_OnTrack::uniqueClone
std::unique_ptr< RIO_OnTrack > uniqueClone() const
NVI clone returning unique_ptr.
Definition: RIO_OnTrack.h:97
Muon::DCMathSegmentMaker::m_reject1DTgcSpacePoints
Gaudi::Property< bool > m_reject1DTgcSpacePoints
Definition: DCMathSegmentMaker.h:423
make_hlt_rep.pars
pars
Definition: make_hlt_rep.py:90
Muon::DCMathSegmentMaker::residualAndPullWithSegment
static std::pair< double, double > residualAndPullWithSegment(const TrkDriftCircleMath::Segment &segment, const Cluster2D &spacePoint, const Amg::Transform3D &gToStation)
Definition: DCMathSegmentMaker.cxx:1675
beamspotman.r
def r
Definition: beamspotman.py:676
Trk::LocalParameters
Definition: LocalParameters.h:98
Muon::DCMathSegmentMaker::errorScaleFactor
double errorScaleFactor(const Identifier &id, double curvature, bool hasPhiMeasurements) const
calculate error scaling factor
Definition: DCMathSegmentMaker.cxx:651
Muon::DCMathSegmentMaker::m_addUnassociatedPhiHits
Gaudi::Property< bool > m_addUnassociatedPhiHits
Definition: DCMathSegmentMaker.h:416
Trk::PlaneSurface::localToGlobalDirection
void localToGlobalDirection(const Trk::LocalDirection &locdir, Amg::Vector3D &globdir) const
This method transforms a local direction wrt the plane to a global direction.
Definition: PlaneSurface.cxx:242
MdtReadoutElement.h
TrkDriftCircleMath::MdtChamberGeometry::print
void print(MsgStream &msg) const override
Definition: MdtChamberGeometry.cxx:178
Muon::DCMathSegmentMaker::m_createCompetingROTsEta
Gaudi::Property< bool > m_createCompetingROTsEta
Definition: DCMathSegmentMaker.h:413
Muon::DCMathSegmentMaker::ClusterVec
std::vector< Cluster2D > ClusterVec
Definition: DCMathSegmentMaker.h:164
Trk::Segment::DCMathSegmentMaker
@ DCMathSegmentMaker
Definition: TrkEvent/TrkSegment/TrkSegment/Segment.h:68
dumpTgcDigiDeadChambers.gasGap
list gasGap
Definition: dumpTgcDigiDeadChambers.py:33
PlotCalibFromCool.norm
norm
Definition: PlotCalibFromCool.py:100
Muon::DCMathSegmentMaker::associateClustersToSegment
std::pair< std::pair< int, int >, bool > associateClustersToSegment(const TrkDriftCircleMath::Segment &segment, const Identifier &chid, const Amg::Transform3D &gToStation, ClusterVecPair &spVecs, double phimin, double phimax, std::vector< std::pair< double, std::unique_ptr< const Trk::MeasurementBase >> > &rioDistVec) const
Definition: DCMathSegmentMaker.cxx:1292
MuonGM::MdtReadoutElement::getNLayers
int getNLayers() const
Returns the number of tube layers inside the multilayer.
MuonGM::MuonReadoutElement::AmdbLRSToGlobalTransform
virtual Amg::Transform3D AmdbLRSToGlobalTransform() const
Definition: MuonDetDescr/MuonReadoutGeometry/src/MuonReadoutElement.cxx:145
Muon::DCMathSegmentMaker::DCMathSegmentMaker
DCMathSegmentMaker(const std::string &, const std::string &, const IInterface *)
Definition: DCMathSegmentMaker.cxx:62
Muon::DCMathSegmentMaker::ClusterVecPair
std::pair< ClusterVec, ClusterVec > ClusterVecPair
Definition: DCMathSegmentMaker.h:165
checkFileSG.line
line
Definition: checkFileSG.py:75
Muon::DCMathSegmentMaker::Cluster2D::phiHit
const MuonClusterOnTrack * phiHit
Definition: DCMathSegmentMaker.h:146
Amg::MatrixX
Eigen::Matrix< double, Eigen::Dynamic, Eigen::Dynamic > MatrixX
Dynamic Matrix - dynamic allocation.
Definition: EventPrimitives.h:29
get_generator_info.result
result
Definition: get_generator_info.py:21
TrkDriftCircleMath::DCOnTrackVec
std::vector< DCOnTrack > DCOnTrackVec
Definition: DCOnTrack.h:59
python.PerfMonSerializer.p
def p
Definition: PerfMonSerializer.py:743
max
#define max(a, b)
Definition: cfImp.cxx:41
StraightLineSurface.h
TrackParameters.h
Muon::DCMathSegmentMaker::segmentCreationInfo::phimax
double phimax
Definition: DCMathSegmentMaker.h:189
Muon::IdDataVec::id
Identifier id
Definition: DCMathSegmentMaker.cxx:1277
MuonGM::RpcReadoutElement::StripLength
double StripLength(bool measphi) const
returns the strip length for the phi or eta plane
Muon::MuonClusterOnTrack::globalPosition
virtual const Amg::Vector3D & globalPosition() const override
Returns global position.
Definition: MuonClusterOnTrack.cxx:93
phi
Scalar phi() const
phi method
Definition: AmgMatrixBasePlugin.h:64
Muon::MuonTubeIntersect
Definition: MuonTubeIntersect.h:12
SG::ReadCondHandle
Definition: ReadCondHandle.h:44
Muon::MdtDriftCircleOnTrack::prepRawData
virtual const MdtPrepData * prepRawData() const override final
Returns the PrepRawData used to create this corrected measurement.
Definition: MdtDriftCircleOnTrack.h:257
Trk::locX
@ locX
Definition: ParamDefs.h:43
find
std::string find(const std::string &s)
return a remapped string
Definition: hcg.cxx:135
TrkDriftCircleMath::DCOnTrack::Delta
@ Delta
used in fit
Definition: DCOnTrack.h:21
TrkDriftCircleMath::MdtId
Definition: MdtId.h:14
TrkDriftCircleMath::ClusterId
Definition: ClusterId.h:12
Trk::locY
@ locY
local cartesian
Definition: ParamDefs.h:44
MuonGM::MdtReadoutElement::center
virtual const Amg::Vector3D & center(const Identifier &) const override final
Return the center of the surface associated with this identifier In the case of silicon it returns th...
Muon::DCMathSegmentMaker::Cluster2D::gasGapId
Identifier gasGapId
Definition: DCMathSegmentMaker.h:142
Muon::DCMathSegmentMaker::m_segmentSelectionTool
ToolHandle< IMuonSegmentSelectionTool > m_segmentSelectionTool
Definition: DCMathSegmentMaker.h:398
MuonGM::MdtReadoutElement::innerTubeRadius
double innerTubeRadius() const
Returns the inner tube radius excluding the aluminium walls.
MuonGM::MuonClusterReadoutElement::center
virtual const Amg::Vector3D & center() const override
Return the center of the element.
Definition: MuonClusterReadoutElement.h:125
SegmentFinder.h
Muon::DCMathSegmentMaker::create1DClusters
ClusterVecPair create1DClusters(const std::vector< const MuonClusterOnTrack * > &clusters) const
Definition: DCMathSegmentMaker.cxx:702
ClusterSeg::residual
@ residual
Definition: ClusterNtuple.h:20
TrkDriftCircleMath::DCSLFitter::fit
virtual bool fit(Segment &result, const Line &line, const DCOnTrackVec &dcs, double t0Seed=-99999.) const
Definition: Tracking/TrkUtilityPackages/TrkDriftCircleMath/TrkDriftCircleMath/DCSLFitter.h:38
calibdata.chamber
chamber
Definition: calibdata.py:32
eta
Scalar eta() const
pseudorapidity method
Definition: AmgMatrixBasePlugin.h:79
Trk::ParametersBase::position
const Amg::Vector3D & position() const
Access method for the position.
Amg::Vector2D
Eigen::Matrix< double, 2, 1 > Vector2D
Definition: GeoPrimitives.h:48
SG::ReadHandle
Definition: StoreGate/StoreGate/ReadHandle.h:70
index
Definition: index.py:1
TrkDriftCircleMath::ChamberGeometry
Definition: ChamberGeometry.h:17
xAODP4Helpers.h
Muon::DCMathSegmentMaker::find
void find(const std::vector< const Trk::RIO_OnTrack * > &rios, Trk::SegmentCollection *segColl=nullptr) const
find segments starting from a list of RIO_OnTrack objects, implementation of IMuonSegmentMaker interf...
Definition: DCMathSegmentMaker.cxx:596
hist_file_dump.d
d
Definition: hist_file_dump.py:137
Muon::MuonClusterOnTrack::prepRawData
virtual const MuonCluster * prepRawData() const override=0
Returns the Trk::PrepRawData - is a MuonCluster in this scope.
Muon::SortByDistanceToSegment
Function object to sort pairs containing a double and a pointer to a MuonClusterOnTrack.
Definition: DCMathSegmentMaker.h:88
TrkDriftCircleMath::CLVec
std::vector< Cluster > CLVec
Definition: Tracking/TrkUtilityPackages/TrkDriftCircleMath/TrkDriftCircleMath/Cluster.h:70
Muon::MuonStationIndex::BML
@ BML
Definition: MuonStationIndex.h:17
Muon::IdDataVec::data
EntryVec data
Definition: DCMathSegmentMaker.cxx:1278
EventPrimitivesHelpers.h
Muon::DCMathSegmentMaker::m_removeDeltas
Gaudi::Property< bool > m_removeDeltas
Definition: DCMathSegmentMaker.h:422
MuonGM::MuonReadoutElement::GlobalToAmdbLRSTransform
virtual Amg::Transform3D GlobalToAmdbLRSTransform() const
Definition: MuonDetDescr/MuonReadoutGeometry/src/MuonReadoutElement.cxx:153
Trk::PrepRawData::localCovariance
const Amg::MatrixX & localCovariance() const
return const ref to the error matrix
Muon::DCMathSegmentMaker::createDCVec
TrkDriftCircleMath::DCVec createDCVec(const std::vector< const MdtDriftCircleOnTrack * > &mdts, double errorScale, std::set< Identifier > &chamberSet, double &phimin, double &phimax, TrkDriftCircleMath::DCStatistics &dcStatistics, const Amg::Transform3D &gToStation, const Amg::Transform3D &amdbToGlobal) const
Definition: DCMathSegmentMaker.cxx:1028
Muon::DCMathSegmentMaker::m_printer
PublicToolHandle< MuonEDMPrinterTool > m_printer
Definition: DCMathSegmentMaker.h:383
skel.it
it
Definition: skel.GENtoEVGEN.py:423
Muon::SortClByPull
Definition: DCMathSegmentMaker.cxx:1285
python.AthDsoLogger.out
out
Definition: AthDsoLogger.py:71
TrkDriftCircleMath::DCSLFitter
Definition: Tracking/TrkUtilityPackages/TrkDriftCircleMath/TrkDriftCircleMath/DCSLFitter.h:17
CompetingMuonClustersOnTrack.h
M_PI
#define M_PI
Definition: ActiveFraction.h:11
TrkDriftCircleMath::ResidualWithLine::residual
double residual(const LocVec2D &pos) const
Definition: ResidualWithLine.h:21
Muon::DCMathSegmentMaker::m_segmentFinder
ToolHandle< IMdtSegmentFinder > m_segmentFinder
Definition: DCMathSegmentMaker.h:388
Muon::DCMathSegmentMaker::Cluster2D::error
double error
Definition: DCMathSegmentMaker.h:144
CaloClusterListBadChannel.cls
cls
Definition: CaloClusterListBadChannel.py:8
deg
#define deg
Definition: SbPolyhedron.cxx:17
Muon::DCMathSegmentMaker::ChIdHitMap
std::map< Identifier, IdHitMap > ChIdHitMap
Definition: DCMathSegmentMaker.h:108
AthCommonMsg< AlgTool >::msgLvl
bool msgLvl(const MSG::Level lvl) const
Definition: AthCommonMsg.h:30
dq_defect_virtual_defect_validation.d1
d1
Definition: dq_defect_virtual_defect_validation.py:79
xAOD::P4Helpers::deltaPhi
double deltaPhi(double phiA, double phiB)
delta Phi in range [-pi,pi[
Definition: xAODP4Helpers.h:69
MuonGM::RpcReadoutElement
An RpcReadoutElement corresponds to a single RPC module; therefore typicaly a barrel muon station con...
Definition: MuonDetDescr/MuonReadoutGeometry/MuonReadoutGeometry/RpcReadoutElement.h:54
Trk::loc2
@ loc2
generic first and second local coordinate
Definition: ParamDefs.h:41
GasGapData
Definition: RPCDQUtils.h:98
DriftCircle.h
Trk::RIO_OnTrack
Definition: RIO_OnTrack.h:70
Muon::DCMathSegmentMaker::Cluster2D::identify
Identifier identify() const
Definition: DCMathSegmentMaker.h:154
TrkDriftCircleMath::LocVec2D::x
double x() const
Returns the x coordinate of the vector.
Definition: LocVec2D.h:27
MdtPrepData.h
python.TurnDataReader.dr
dr
Definition: TurnDataReader.py:112
TrkDriftCircleMath::DCOnTrack::OutOfTime
@ OutOfTime
delta electron
Definition: DCOnTrack.h:22
Muon::DCMathSegmentMaker::updateSegmentPhi
bool updateSegmentPhi(const Amg::Vector3D &gpos, const Amg::Vector3D &gdir, Amg::Vector2D &segLocPos, Trk::LocalDirection &segLocDir, Trk::PlaneSurface &surf, const std::vector< const Trk::MeasurementBase * > &rots, double phimin, double phimax) const
Definition: DCMathSegmentMaker.cxx:1854
MuonSegmentQuality.h
Trk::locR
@ locR
Definition: ParamDefs.h:50
TrkDriftCircleMath::DriftCircle
This class represents a drift time measurement.
Definition: DriftCircle.h:22
read_hist_ntuple.t
t
Definition: read_hist_ntuple.py:5
drawFromPickle.cos
cos
Definition: drawFromPickle.py:36
Muon::MuonStationIndex::BMS
@ BMS
Definition: MuonStationIndex.h:17
Muon::DCMathSegmentMaker::Cluster2D::phiHits
std::vector< const MuonClusterOnTrack * > phiHits
Definition: DCMathSegmentMaker.h:147
Muon::MdtDriftCircleOnTrack::detectorElement
virtual const MuonGM::MdtReadoutElement * detectorElement() const override final
Returns the detector element, assoicated with the PRD of this class.
Definition: MdtDriftCircleOnTrack.h:268
ATH_MSG_VERBOSE
#define ATH_MSG_VERBOSE(x)
Definition: AthMsgStreamMacros.h:28
Muon::CompetingMuonClustersOnTrack
Definition: CompetingMuonClustersOnTrack.h:54
MdtDriftCircleOnTrack.h
Muon::DCMathSegmentMaker::TubeEnds::lxmax
double lxmax
Definition: DCMathSegmentMaker.h:170
Trk::Surface::center
const Amg::Vector3D & center() const
Returns the center position of the Surface.
Muon
This class provides conversion from CSC RDO data to CSC Digits.
Definition: TrackSystemController.h:49
drawFromPickle.exp
exp
Definition: drawFromPickle.py:36
Muon::IdDataVec::Entry
T Entry
Definition: DCMathSegmentMaker.cxx:1271
MuonGM::TgcReadoutElement::stripDir
Amg::Vector3D stripDir(int gasGap, int strip) const
Returns the direction of a strip.
yodamerge_tmp.scale
scale
Definition: yodamerge_tmp.py:138
Muon::DCMathSegmentMaker::updatePhiRanges
static void updatePhiRanges(double phiminNew, double phimaxNew, double &phiminRef, double &phimaxRef)
update phi ranges
Definition: DCMathSegmentMaker.cxx:2084
x
#define x
ParamDefs.h
Muon::MuonStationIndex::PhiIndex
PhiIndex
enum to classify the different phi layers in the muon spectrometer
Definition: MuonStationIndex.h:31
TgcPrepData.h
Muon::DCMathSegmentMaker::m_DetectorManagerKey
SG::ReadCondHandleKey< MuonGM::MuonDetectorManager > m_DetectorManagerKey
pointers to IdHelpers
Definition: DCMathSegmentMaker.h:344
TrkDriftCircleMath::TransformToLine
Definition: TransformToLine.h:15
Muon::DCMathSegmentMaker::createSpacePoints
ClusterVecPair createSpacePoints(const ChIdHitMap &chIdHitMap) const
Definition: DCMathSegmentMaker.cxx:749
Muon::MdtDriftCircleOnTrack::globalPosition
virtual const Amg::Vector3D & globalPosition() const override final
Returns the global Position.
Definition: MdtDriftCircleOnTrack.cxx:108
Muon::MdtDriftCircleOnTrack::driftTime
double driftTime() const
Returns the value of the drift time used to obtain the drift radius.
Definition: MdtDriftCircleOnTrack.h:280
Muon::DCMathSegmentMaker::checkBoundsInXZ
bool checkBoundsInXZ(double xline, double zline, double dXdZ, const std::vector< HitInXZ > &hits) const
check whether all hits are in bounds in the XZ plane
Definition: DCMathSegmentMaker.cxx:1830
TrkDriftCircleMath::Segment
Definition: TrkUtilityPackages/TrkDriftCircleMath/TrkDriftCircleMath/Segment.h:18
TrkDriftCircleMath::LocVec2D
Implementation of 2 dimensional vector class.
Definition: LocVec2D.h:16
TrkDriftCircleMath::DCVec
std::vector< DriftCircle > DCVec
Definition: DriftCircle.h:117
Muon::DCMathSegmentMaker::m_curvedErrorScaling
Gaudi::Property< bool > m_curvedErrorScaling
Definition: DCMathSegmentMaker.h:411
Amg::getRotateZ3D
Amg::Transform3D getRotateZ3D(double angle)
get a rotation transformation around Z-axis
Definition: GeoPrimitivesHelpers.h:270
Trk::LocalDirection::angleYZ
double angleYZ() const
access method for angle of local YZ projection
Definition: LocalDirection.h:106
Muon::IdDataVec
Definition: DCMathSegmentMaker.cxx:1270
Muon::DCMathSegmentMaker::m_mdtCreator
ToolHandle< IMdtDriftCircleOnTrackCreator > m_mdtCreator
Definition: DCMathSegmentMaker.h:363
Trk::TrackStateOnSurface::Outlier
@ Outlier
This TSoS contains an outlier, that is, it contains a MeasurementBase/RIO_OnTrack which was not used ...
Definition: TrackStateOnSurface.h:122
Muon::DCMathSegmentMaker::segmentCreationInfo::globalTrans
Amg::Transform3D globalTrans
Definition: DCMathSegmentMaker.h:186
TRT::Hit::side
@ side
Definition: HitInfo.h:83
Track.h
TgcClusterOnTrack.h
Trk::StraightLineSurface::globalToLocal
virtual bool globalToLocal(const Amg::Vector3D &glob, const Amg::Vector3D &mom, Amg::Vector2D &loc) const override final
Specified for StraightLineSurface: GlobalToLocal method without dynamic memory allocation This method...
Definition: StraightLineSurface.cxx:164
TrkDriftCircleMath::Line
Definition: Line.h:17
Muon::MuonClusterOnTrack::detectorElement
virtual const MuonGM::MuonClusterReadoutElement * detectorElement() const override=0
Returns the detector element, associated with the PRD of this class.
Muon::DCMathSegmentMaker::Cluster2D::detElId
Identifier detElId
Definition: DCMathSegmentMaker.h:141
python.setupRTTAlg.size
int size
Definition: setupRTTAlg.py:39
TrkDriftCircleMath::Road
TrkDriftCircleMath::Road - encodes the road given to the segment finder in station coordinates.
Definition: Road.h:15
Trk::TrackRoad::globalDirection
const Amg::Vector3D & globalDirection() const
Get the global direction of the road.
Definition: TrackRoad.h:143
Muon::DCMathSegmentMaker::m_idHelperSvc
ServiceHandle< Muon::IMuonIdHelperSvc > m_idHelperSvc
Definition: DCMathSegmentMaker.h:351
DCMathSegmentMaker.h
Muon::DCMathSegmentMaker::m_compClusterCreator
ToolHandle< IMuonCompetingClustersOnTrackCreator > m_compClusterCreator
Definition: DCMathSegmentMaker.h:378
RpcClusterOnTrack.h
Muon::DCMathSegmentMaker::createROTVec
static DataVector< const Trk::MeasurementBase > createROTVec(std::vector< std::pair< double, std::unique_ptr< const Trk::MeasurementBase >> > &rioDistVec)
Definition: DCMathSegmentMaker.cxx:1663
Muon::IdDataVec::IdDataVec
IdDataVec()=default
Amg::toString
std::string toString(const Translation3D &translation, int precision=4)
GeoPrimitvesToStringConverter.
Definition: GeoPrimitivesToStringConverter.h:40
Muon::DCMathSegmentMaker::Cluster2D::etaHit
const MuonClusterOnTrack * etaHit
Definition: DCMathSegmentMaker.h:145
ATH_MSG_ERROR
#define ATH_MSG_ERROR(x)
Definition: AthMsgStreamMacros.h:33
Muon::DCMathSegmentMaker::createClusterVec
TrkDriftCircleMath::CLVec createClusterVec(const Identifier &chid, ClusterVec &spVec, const Amg::Transform3D &gToStation) const
Definition: DCMathSegmentMaker.cxx:986
DefinedParameter.h
Muon::DCMathSegmentMaker::m_usePreciseError
Gaudi::Property< bool > m_usePreciseError
Definition: DCMathSegmentMaker.h:424
MuonGM::MuonDetectorManager::getMdtReadoutElement
const MdtReadoutElement * getMdtReadoutElement(const Identifier &id) const
access via extended identifier (requires unpacking)
Definition: MuonDetDescr/MuonReadoutGeometry/src/MuonDetectorManager.cxx:204
EventPrimitivesToStringConverter.h
ParticleGun_EoverP_Config.momentum
momentum
Definition: ParticleGun_EoverP_Config.py:63
MuonGM::MdtReadoutElement
Definition: MuonDetDescr/MuonReadoutGeometry/MuonReadoutGeometry/MdtReadoutElement.h:50
Muon::DCMathSegmentMaker::TubeEnds::phimin
double phimin
Definition: DCMathSegmentMaker.h:171
Muon::DCMathSegmentMaker::associateMDTsToSegment
void associateMDTsToSegment(const Amg::Vector3D &gdir, TrkDriftCircleMath::Segment &segment, const TrkDriftCircleMath::ChamberGeometry *multiGeo, const Amg::Transform3D &gToStation, const Amg::Transform3D &amdbToGlobal, std::set< Identifier > &deltaVec, std::set< Identifier > &outoftimeVec, std::vector< std::pair< double, std::unique_ptr< const Trk::MeasurementBase >> > &rioDistVec) const
Definition: DCMathSegmentMaker.cxx:1157
lumiFormat.i
int i
Definition: lumiFormat.py:92
SG::OWN_ELEMENTS
@ OWN_ELEMENTS
this data object owns its elements
Definition: OwnershipPolicy.h:17
z
#define z
Muon::DCMathSegmentMaker::errorScalingRegion
bool errorScalingRegion(const Identifier &id) const
apply error scaling for low mometum tracks
Definition: DCMathSegmentMaker.cxx:684
Trk::TrackRoad::globalPosition
const Amg::Vector3D & globalPosition() const
Get the global position of the road.
Definition: TrackRoad.h:137
TrkDriftCircleMath::DCSLHitSelector
Definition: DCSLHitSelector.h:13
Muon::DCMathSegmentMaker::m_segmentFitter
ToolHandle< IMuonSegmentFittingTool > m_segmentFitter
Definition: DCMathSegmentMaker.h:393
Identifier
Definition: DetectorDescription/Identifier/Identifier/Identifier.h:32
beamspotman.n
n
Definition: beamspotman.py:731
Trk::theta
@ theta
Definition: ParamDefs.h:72
Trk::MeasurementBase::uniqueClone
std::unique_ptr< MeasurementBase > uniqueClone() const
NVI Clone giving up unique pointer.
Definition: MeasurementBase.h:77
EL::StatusCode
::StatusCode StatusCode
StatusCode definition for legacy code.
Definition: PhysicsAnalysis/D3PDTools/EventLoop/EventLoop/StatusCode.h:22
ATH_MSG_DEBUG
#define ATH_MSG_DEBUG(x)
Definition: AthMsgStreamMacros.h:29
TrkDriftCircleMath::ChamberGeometry::nlay
virtual unsigned int nlay() const =0
angle
double angle(const GeoTrf::Vector2D &a, const GeoTrf::Vector2D &b)
Definition: TRTDetectorFactory_Full.cxx:73
xAOD::covMatrix
covMatrix
Definition: TrackMeasurement_v1.cxx:19
LArG4ShowerLibProcessing.hits
hits
Definition: LArG4ShowerLibProcessing.py:136
TRT::Hit::layer
@ layer
Definition: HitInfo.h:79
Muon::DCMathSegmentMaker::Cluster2D
Definition: DCMathSegmentMaker.h:125
Trk::driftRadius
@ driftRadius
trt, straws
Definition: ParamDefs.h:59
MuonGM::MdtReadoutElement::getActiveTubeLength
double getActiveTubeLength(const int tubeLayer, const int tube) const
Muon::DCMathSegmentMaker::segmentCreationInfo::phimin
double phimin
Definition: DCMathSegmentMaker.h:188
TrkDriftCircleMath::MdtChamberGeometry
Definition: MdtChamberGeometry.h:18
Muon::DCMathSegmentMaker::m_edmHelperSvc
ServiceHandle< IMuonEDMHelperSvc > m_edmHelperSvc
Definition: DCMathSegmentMaker.h:356
Amg::Transform3D
Eigen::Affine3d Transform3D
Definition: GeoPrimitives.h:46
Trk::RIGHT
@ RIGHT
the drift radius is positive (see Trk::AtaStraightLine)
Definition: DriftCircleSide.h:22
Muon::DCMathSegmentMaker::m_rpcKey
SG::ReadHandleKey< Muon::RpcPrepDataContainer > m_rpcKey
Definition: DCMathSegmentMaker.h:437
Amg::transform
Amg::Vector3D transform(Amg::Vector3D &v, Amg::Transform3D &tr)
Transform a point from a Trasformation3D.
Definition: GeoPrimitivesHelpers.h:156
Muon::DCMathSegmentMaker::m_chamberGeoKey
SG::ReadCondHandleKey< Muon::MuonIntersectGeoData > m_chamberGeoKey
Definition: DCMathSegmentMaker.h:441
MuonGM::TgcReadoutElement
A TgcReadoutElement corresponds to a single TGC chamber; therefore typically a TGC station contains s...
Definition: MuonDetDescr/MuonReadoutGeometry/MuonReadoutGeometry/TgcReadoutElement.h:42
Muon::DCMathSegmentMaker::Cluster2D::corrupt
bool corrupt() const
Definition: DCMathSegmentMaker.h:162
Muon::SortIdDataVec::operator()
bool operator()(const IdDataVec< T > &d1, const IdDataVec< T > &d2)
Definition: DCMathSegmentMaker.cxx:1282
TrkDriftCircleMath::LocVec2D::y
double y() const
Returns the y coordinate of the vector.
Definition: LocVec2D.h:29
Muon::DCMathSegmentMaker::TubeEnds::phimax
double phimax
Definition: DCMathSegmentMaker.h:172
Muon::DCMathSegmentMaker::TubeEnds::lxmin
double lxmin
Definition: DCMathSegmentMaker.h:169
Muon::MdtPrepData::globalPosition
virtual const Amg::Vector3D & globalPosition() const
Returns the global position of the CENTER of the drift tube (i.e.
Definition: MdtPrepData.h:149
python.StandardJetMods.pull
pull
Definition: StandardJetMods.py:264
ClusterId.h
Muon::DCMathSegmentMaker::createSpacePoint
Cluster2D createSpacePoint(const Identifier &gasGapId, const MuonClusterOnTrack *etaHit, const MuonClusterOnTrack *phiHit) const
Definition: DCMathSegmentMaker.cxx:852
ATH_CHECK
#define ATH_CHECK
Definition: AthCheckMacros.h:40
MuonGM::TgcReadoutElement::localSpacePoint
Amg::Vector3D localSpacePoint(const Identifier &stripId, const Amg::Vector3D &etaHitPos, const Amg::Vector3D &phiHitPos) const
Definition: MuonDetDescr/MuonReadoutGeometry/src/TgcReadoutElement.cxx:315
Muon::IdDataVec::IdDataVec
IdDataVec(const Identifier &i)
Definition: DCMathSegmentMaker.cxx:1275
MuonGM::MdtReadoutElement::nMDTinStation
unsigned int nMDTinStation() const
How many MDT chambers are in the station.
Definition: MuonDetDescr/MuonReadoutGeometry/MuonReadoutGeometry/MdtReadoutElement.h:61
Muon::DCMathSegmentMaker::Cluster2D::is2D
bool is2D() const
Definition: DCMathSegmentMaker.h:161
Muon::DCMathSegmentMaker::segmentCreationInfo::clusters
ClusterVecPair clusters
Definition: DCMathSegmentMaker.h:184
Muon::DCMathSegmentMaker::m_refitParameters
Gaudi::Property< bool > m_refitParameters
Definition: DCMathSegmentMaker.h:415
Trk::PlaneSurface::globalToLocalDirection
void globalToLocalDirection(const Amg::Vector3D &glodir, Trk::LocalDirection &locdir) const
This method transforms the global direction to a local direction wrt the plane.
Definition: PlaneSurface.cxx:260
Muon::DCMathSegmentMaker::m_allMdtHoles
Gaudi::Property< bool > m_allMdtHoles
Definition: DCMathSegmentMaker.h:421
Muon::DCMathSegmentMaker::createSegment
std::unique_ptr< MuonSegment > createSegment(const EventContext &ctx, TrkDriftCircleMath::Segment &segment, const Identifier &chid, const Amg::Vector3D &roadpos, const Amg::Vector3D &roaddir2, const std::vector< const MdtDriftCircleOnTrack * > &mdts, bool hasPhiMeasurements, segmentCreationInfo &sInfo) const
Definition: DCMathSegmentMaker.cxx:234
drawFromPickle.tan
tan
Definition: drawFromPickle.py:36
Trk::ParametersBase
Definition: ParametersBase.h:55
Muon::DCMathSegmentMaker::initialize
virtual StatusCode initialize()
Definition: DCMathSegmentMaker.cxx:67
Muon::DCMathSegmentMaker::m_createCompetingROTsPhi
Gaudi::Property< bool > m_createCompetingROTsPhi
Definition: DCMathSegmentMaker.h:414
Muon::DCMathSegmentMaker::m_recoverBadRpcCabling
Gaudi::Property< bool > m_recoverBadRpcCabling
Definition: DCMathSegmentMaker.h:429
Muon::MdtDriftCircleOnTrack::driftRadius
double driftRadius() const
Returns the value of the drift radius.
Definition: MdtDriftCircleOnTrack.h:277
Athena::Timeout::instance
static Timeout & instance()
Get reference to Timeout singleton.
Definition: Timeout.h:64
a0
double a0
Definition: globals.cxx:27
TrkDriftCircleMath::TransformToLine::toLocal
LocVec2D toLocal(const LocVec2D &pos) const
Definition: TransformToLine.h:26
Muon::MuonPrepDataCollection
Template to hold collections of MuonPrepRawData objects.
Definition: MuonPrepDataCollection.h:46
DataVector< Trk::Segment >
WriteCalibToCool.swap
swap
Definition: WriteCalibToCool.py:94
Trk::LocalDirection
represents the three-dimensional global direction with respect to a planar surface frame.
Definition: LocalDirection.h:81
Muon::DCMathSegmentMaker::IdHitMap
std::map< Identifier, EtaPhiHitsPair > IdHitMap
Definition: DCMathSegmentMaker.h:107
TrkDriftCircleMath::DCOnTrack::OnTrack
@ OnTrack
Definition: DCOnTrack.h:20
dot.dot
def dot(G, fn, nodesToHighlight=[])
Definition: dot.py:5
Trk::MeasurementBase::localCovariance
const Amg::MatrixX & localCovariance() const
Interface method to get the localError.
Definition: MeasurementBase.h:138
Trk::Segment::DCMathSegmentMakerCurved
@ DCMathSegmentMakerCurved
Definition: TrkEvent/TrkSegment/TrkSegment/Segment.h:75
min
#define min(a, b)
Definition: cfImp.cxx:40
Road.h
Muon::DCMathSegmentMaker::localTubeEnds
TubeEnds localTubeEnds(const MdtDriftCircleOnTrack &mdt, const Amg::Transform3D &gToSegment, const Amg::Transform3D &segmentToG) const
calculate positions of tube ends
Definition: DCMathSegmentMaker.cxx:2054
Trk::MeasurementBase
Definition: MeasurementBase.h:58
TrkDriftCircleMath::TransformToLine::toLine
LocVec2D toLine(const LocVec2D &pos) const
Definition: TransformToLine.h:20
Trk::TrackStateOnSurface
represents the track state (measurement, material, fit parameters and quality) at a surface.
Definition: TrackStateOnSurface.h:71
name
std::string name
Definition: Control/AthContainers/Root/debug.cxx:192
Muon::DCMathSegmentMaker::calculateHoles
std::vector< Identifier > calculateHoles(const EventContext &ctx, Identifier chid, const Amg::Vector3D &gpos, const Amg::Vector3D &gdir, bool hasMeasuredCoordinate, std::set< Identifier > &deltaVec, std::set< Identifier > &outoftimeVec, const std::vector< std::pair< double, std::unique_ptr< const Trk::MeasurementBase >> > &rioDistVec) const
Definition: DCMathSegmentMaker.cxx:1697
RIO_OnTrack.h
Muon::DCMathSegmentMaker::m_mdtKey
SG::ReadHandleKey< Muon::MdtPrepDataContainer > m_mdtKey
Definition: DCMathSegmentMaker.h:439
TrkDriftCircleMath::Cluster
class representing a cluster meaurement
Definition: Tracking/TrkUtilityPackages/TrkDriftCircleMath/TrkDriftCircleMath/Cluster.h:22
Muon::DCMathSegmentMaker::create2DClusters
ClusterVecPair create2DClusters(const std::vector< const MuonClusterOnTrack * > &clusters) const
Definition: DCMathSegmentMaker.cxx:728
plotBeamSpotMon.b
b
Definition: plotBeamSpotMon.py:77
Amg::error
double error(const Amg::MatrixX &mat, int index)
return diagonal error of the matrix caller should ensure the matrix is symmetric and the index is in ...
Definition: EventPrimitivesHelpers.h:40
Trk::PrepRawData::localPosition
const Amg::Vector2D & localPosition() const
return the local position reference
Muon::MdtDriftCircleOnTrack
This class represents the corrected MDT measurements, where the corrections include the effects of wi...
Definition: MdtDriftCircleOnTrack.h:37
charge
double charge(const T &p)
Definition: AtlasPID.h:494
Muon::DCMathSegmentMaker::findFirstRotInChamberWithMostHits
const MdtDriftCircleOnTrack * findFirstRotInChamberWithMostHits(const std::vector< const MdtDriftCircleOnTrack * > &mdts) const
Definition: DCMathSegmentMaker.cxx:1800
DataVector::push_back
value_type push_back(value_type pElem)
Add an element to the end of the collection.
SG::CondHandleKey::initialize
StatusCode initialize(bool used=true)
TrkDriftCircleMath::ResidualWithSegment
class to calculate residual of a hit with a segment and calculate the local track errors
Definition: ResidualWithSegment.h:20
Muon::MdtPrepData
Class to represent measurements from the Monitored Drift Tubes.
Definition: MdtPrepData.h:37
Muon::CompetingMuonClustersOnTrack::containedROTs
const std::vector< const MuonClusterOnTrack * > & containedROTs() const
returns the vector of SCT_ClusterOnTrack objects .
Definition: CompetingMuonClustersOnTrack.h:184
Amg::Vector3D
Eigen::Matrix< double, 3, 1 > Vector3D
Definition: GeoPrimitives.h:47
LocalParameters.h
MuonDetectorManager.h
TrkDriftCircleMath::ChamberGeometry::tubePosition
virtual LocVec2D tubePosition(unsigned int ml, unsigned int lay, unsigned int tube) const =0
Muon::DCMathSegmentMaker::m_dcslFitProvider
ToolHandle< IDCSLFitProvider > m_dcslFitProvider
Definition: DCMathSegmentMaker.h:403
MuonGM::MdtReadoutElement::ROPos
Amg::Vector3D ROPos(const int tubelayer, const int tube) const
Definition: MuonDetDescr/MuonReadoutGeometry/src/MdtReadoutElement.cxx:288
Trk::MeasurementBase::localParameters
const LocalParameters & localParameters() const
Interface method to get the LocalParameters.
Definition: MeasurementBase.h:132
Muon::DCMathSegmentMaker::m_doSpacePoints
Gaudi::Property< bool > m_doSpacePoints
Definition: DCMathSegmentMaker.h:412
Trk::DriftCircleSide
DriftCircleSide
Enumerates the 'side' of the wire on which the tracks passed (i.e.
Definition: DriftCircleSide.h:16
Muon::SortIdDataVec
Definition: DCMathSegmentMaker.cxx:1281
eflowRec::phiIndex
unsigned int phiIndex(float phi, float binsize)
calculate phi index for a given phi
Definition: EtaPhiLUT.cxx:23
Timeout.h
Timeout singleton.
makeTRTBarrelCans.dy
tuple dy
Definition: makeTRTBarrelCans.py:21
Muon::MuonStationIntersect
Definition: MuonStationIntersect.h:12
Muon::DCMathSegmentMaker::createTgcSpacePoint
Cluster2D createTgcSpacePoint(const Identifier &gasGapId, const MuonClusterOnTrack *etaHit, const MuonClusterOnTrack *phiHit) const
Definition: DCMathSegmentMaker.cxx:884
Trk::ParametersBase::momentum
const Amg::Vector3D & momentum() const
Access method for the momentum.
DeMoScan.index
string index
Definition: DeMoScan.py:362
Muon::CompetingMuonClustersOnTrack::rioOnTrack
const MuonClusterOnTrack & rioOnTrack(unsigned int) const
returns the RIO_OnTrack (also known as ROT) objects depending on the integer
Definition: CompetingMuonClustersOnTrack.h:190
Muon::DCMathSegmentMaker::m_doTimeOutChecks
Gaudi::Property< bool > m_doTimeOutChecks
Definition: DCMathSegmentMaker.h:427
Muon::DCMathSegmentMaker::m_maxAssociateClusterDistance
Gaudi::Property< double > m_maxAssociateClusterDistance
Definition: DCMathSegmentMaker.h:419
DiTauMassTools::MaxHistStrategyV2::e
e
Definition: PhysicsAnalysis/TauID/DiTauMassTools/DiTauMassTools/HelperFunctions.h:26
a
TList * a
Definition: liststreamerinfos.cxx:10
Segment.h
y
#define y
Trk::Surface::globalToLocal
virtual bool globalToLocal(const Amg::Vector3D &glob, const Amg::Vector3D &mom, Amg::Vector2D &loc) const =0
Specified by each surface type: GlobalToLocal method without dynamic memory allocation - boolean chec...
Muon::DCMathSegmentMaker::m_sinAngleCut
Gaudi::Property< double > m_sinAngleCut
Definition: DCMathSegmentMaker.h:409
Muon::DCMathSegmentMaker::m_redo2DFit
Gaudi::Property< bool > m_redo2DFit
Definition: DCMathSegmentMaker.h:432
Muon::DCMathSegmentMaker::m_doGeometry
Gaudi::Property< bool > m_doGeometry
Definition: DCMathSegmentMaker.h:410
Amg::intersect
std::optional< double > intersect(const AmgVector(N)&posA, const AmgVector(N)&dirA, const AmgVector(N)&posB, const AmgVector(N)&dirB)
Calculates the closest approach of two lines.
Definition: GeoPrimitivesHelpers.h:302
MuonGM::MuonDetectorManager
The MuonDetectorManager stores the transient representation of the Muon Spectrometer geometry and pro...
Definition: MuonDetDescr/MuonReadoutGeometry/MuonReadoutGeometry/MuonDetectorManager.h:49
ATH_MSG_WARNING
#define ATH_MSG_WARNING(x)
Definition: AthMsgStreamMacros.h:32
Muon::MuonStationIndex::BIL
@ BIL
Definition: MuonStationIndex.h:17
Trk::PlaneSurface
Definition: PlaneSurface.h:64
TrkDriftCircleMath::DCSLHitSelector::selectHitsOnTrack
const HitSelection selectHitsOnTrack(const DCOnTrackVec &dcs) const
Definition: DCSLHitSelector.h:20
Muon::DCMathSegmentMaker::segmentCreationInfo::geom
const TrkDriftCircleMath::ChamberGeometry * geom
Definition: DCMathSegmentMaker.h:185
PlaneSurface.h
Muon::SortClByPull::operator()
bool operator()(const std::pair< double, DCMathSegmentMaker::Cluster2D > &d1, const std::pair< double, DCMathSegmentMaker::Cluster2D > &d2)
Definition: DCMathSegmentMaker.cxx:1286
Muon::DCMathSegmentMaker::findMdt
const MdtPrepData * findMdt(const EventContext &ctx, const Identifier &id) const
Definition: DCMathSegmentMaker.cxx:1781
ResidualWithSegment.h
JacobianCotThetaPtToThetaP.h
Muon::DCMathSegmentMaker::m_assumePointingPhi
Gaudi::Property< bool > m_assumePointingPhi
Definition: DCMathSegmentMaker.h:431
dq_defect_virtual_defect_validation.d2
d2
Definition: dq_defect_virtual_defect_validation.py:81
DEBUG
#define DEBUG
Definition: page_access.h:11
makeTRTBarrelCans.dx
tuple dx
Definition: makeTRTBarrelCans.py:20
Trk::qOverP
@ qOverP
perigee
Definition: ParamDefs.h:73
RunTileMonitoring.clusters
clusters
Definition: RunTileMonitoring.py:133
MuonGM::MdtReadoutElement::surface
virtual const Trk::Surface & surface() const override final
Return surface associated with this detector element.
Definition: MuonDetDescr/MuonReadoutGeometry/src/MdtReadoutElement.cxx:875
Trk::RIO_OnTrack::identify
virtual Identifier identify() const final
return the identifier -extends MeasurementBase
Definition: RIO_OnTrack.h:155
GeoPrimitivesToStringConverter.h
convertTimingResiduals.offset
offset
Definition: convertTimingResiduals.py:71
Muon::DCMathSegmentMaker::m_updatePhiUsingPhiHits
Gaudi::Property< bool > m_updatePhiUsingPhiHits
Definition: DCMathSegmentMaker.h:430
TrkDriftCircleMath::DCOnTrack
class representing a drift circle meaurement on segment
Definition: DCOnTrack.h:16
TrkDriftCircleMath::SegVec
std::vector< Segment > SegVec
Definition: TrkUtilityPackages/TrkDriftCircleMath/TrkDriftCircleMath/Segment.h:122
Muon::DCMathSegmentMaker::m_clusterCreator
ToolHandle< IMuonClusterOnTrackCreator > m_clusterCreator
Definition: DCMathSegmentMaker.h:373
Muon::MuonStationIndex::ChIndex
ChIndex
enum to classify the different chamber layers in the muon spectrometer
Definition: MuonStationIndex.h:15
Muon::DCMathSegmentMaker::distanceToSegment
static double distanceToSegment(const TrkDriftCircleMath::Segment &segment, const Amg::Vector3D &hitPos, const Amg::Transform3D &gToStation)
Definition: DCMathSegmentMaker.cxx:1640
Muon::DCMathSegmentMaker::createRpcSpacePoint
Cluster2D createRpcSpacePoint(const Identifier &gasGapId, const MuonClusterOnTrack *etaHit, const std::vector< const MuonClusterOnTrack * > &phiHits) const
Definition: DCMathSegmentMaker.cxx:925
Muon::DCMathSegmentMaker::updateDirection
Amg::Vector3D updateDirection(double linephi, const Trk::PlaneSurface &surf, const Amg::Vector3D &roaddir, bool isCurvedSegment) const
update the global direction, keeping the phi of the input road direction but using the local angle YZ
Definition: DCMathSegmentMaker.cxx:2132
Muon::DCMathSegmentMaker::segmentCreationInfo::amdbTrans
Amg::Transform3D amdbTrans
Definition: DCMathSegmentMaker.h:187
Trk::LEFT
@ LEFT
the drift radius is negative (see Trk::AtaStraightLine)
Definition: DriftCircleSide.h:20
Trk::phi
@ phi
Definition: ParamDefs.h:81
Trk::RIO_OnTrack::associatedSurface
virtual const Surface & associatedSurface() const override=0
returns the surface for the local to global transformation
LocalDirection.h
xAOD::track
@ track
Definition: TrackingPrimitives.h:512
calibdata.copy
bool copy
Definition: calibdata.py:27
Muon::DCMathSegmentMaker::Cluster2D::globalPos
Amg::Vector3D globalPos
Definition: DCMathSegmentMaker.h:160
Muon::DCMathSegmentMaker::m_outputFittedT0
Gaudi::Property< bool > m_outputFittedT0
Definition: DCMathSegmentMaker.h:425
Muon::DCMathSegmentMaker::m_tgcKey
SG::ReadHandleKey< Muon::TgcPrepDataContainer > m_tgcKey
Definition: DCMathSegmentMaker.h:438
python.Constants.VERBOSE
int VERBOSE
Definition: Control/AthenaCommon/python/Constants.py:14
TrkDriftCircleMath::ResidualWithSegment::trackError2
double trackError2(const DriftCircle &dc) const
calculate the track error at the position of a drift circle
Definition: ResidualWithSegment.h:29
drawFromPickle.sin
sin
Definition: drawFromPickle.py:36
AthAlgTool
Definition: AthAlgTool.h:26
Trk::LocalDirection::angleXZ
double angleXZ() const
access method for angle of local XZ projection
Definition: LocalDirection.h:103
RpcPrepData.h
IdentifierHash
Definition: IdentifierHash.h:38
Trk::PlaneSurface::localToGlobal
virtual void localToGlobal(const Amg::Vector2D &locp, const Amg::Vector3D &mom, Amg::Vector3D &glob) const override final
Specified for PlaneSurface: LocalToGlobal method without dynamic memory allocation.
Definition: PlaneSurface.cxx:204
FitQuality.h
MuonGM::MdtReadoutElement::getNtubesperlayer
int getNtubesperlayer() const
Returns the number of tubes in each tube layer.
Muon::DCMathSegmentMaker::HitInXZ
Definition: DCMathSegmentMaker.h:111
Trk::loc1
@ loc1
Definition: ParamDefs.h:40
Muon::DCMathSegmentMaker::m_mdtCreatorT0
ToolHandle< IMdtDriftCircleOnTrackCreator > m_mdtCreatorT0
Definition: DCMathSegmentMaker.h:368
error
Definition: IImpactPoint3dEstimator.h:70
Trk::Surface::transform
const Amg::Transform3D & transform() const
Returns HepGeom::Transform3D by reference.
Muon::DCMathSegmentMaker::TubeEnds
Definition: DCMathSegmentMaker.h:167
Muon::MuonStationIndex::BOL
@ BOL
Definition: MuonStationIndex.h:17
dq_make_web_display.cl
cl
print [x.__class__ for x in toList(dqregion.getSubRegions()) ]
Definition: dq_make_web_display.py:26
Muon::MdtPrepData::detectorElement
virtual const MuonGM::MdtReadoutElement * detectorElement() const override
Returns the detector element corresponding to this PRD.
Definition: MdtPrepData.h:156
mag
Scalar mag() const
mag method
Definition: AmgMatrixBasePlugin.h:25
TrackRoad.h
calibdata.tube
tube
Definition: calibdata.py:31
python.SystemOfUnits.degree
tuple degree
Definition: SystemOfUnits.py:106
Trk::driftCircle
@ driftCircle
Definition: MeasurementType.h:25
Muon::DCMathSegmentMaker::createChamberGeometry
TrkDriftCircleMath::MdtChamberGeometry createChamberGeometry(const Identifier &chid, const Amg::Transform3D &gToStation) const
Definition: DCMathSegmentMaker.cxx:1081
Muon::DCMathSegmentMaker::checkPhiConsistency
bool checkPhiConsistency(double phi, double phimin, double phimax) const
check whether phi is consistent with segment phi
Definition: DCMathSegmentMaker.cxx:2106
TSU::T
unsigned long long T
Definition: L1TopoDataTypes.h:35
Muon::DCMathSegmentMaker::m_preciseErrorScale
Gaudi::Property< double > m_preciseErrorScale
Definition: DCMathSegmentMaker.h:426
Trk::StraightLineSurface
Definition: StraightLineSurface.h:51
Muon::MuonClusterOnTrack
Base class for Muon cluster RIO_OnTracks.
Definition: MuonClusterOnTrack.h:34
Muon::DCMathSegmentMaker::segmentCreationInfo
Definition: DCMathSegmentMaker.h:175
generate::Zero
void Zero(TH1D *hin)
Definition: generate.cxx:32
NSWL1::PadTriggerAdapter::segment
Muon::NSW_PadTriggerSegment segment(const NSWL1::PadTrigger &data)
Definition: PadTriggerAdapter.cxx:5
Muon::IdDataVec::EntryVec
std::vector< Entry > EntryVec
Definition: DCMathSegmentMaker.cxx:1272
TrkDriftCircleMath::DriftCircle::InTime
@ InTime
drift time too small to be compatible with drift spectrum
Definition: DriftCircle.h:27
Trk::TrackRoad::deltaEta
double deltaEta() const
Get the width of the road in the eta direction.
Definition: TrackRoad.h:149
RpcReadoutElement.h
TrkDriftCircleMath::DCStatistics
This class offers no functionality, but to define a standard device for the maker to transfer to the ...
Definition: DCStatistics.h:19
Generated on Sun Jun 30 2024 21:13:59 for ATLAS Offline Software by   doxygen 1.8.18