20 #include "Acts/Utilities/RangeXD.hpp"
26 double yCross = (
y0 + loc.
location().z() * tanTheta);
27 return (loc.
minY() < yCross && yCross < loc. maxY());
44 const std::array<double, 2>& chambEdges) {
46 if (chambEdges[0] <= seedEdges[0] && chambEdges[1] >= seedEdges[1]) {
50 else if (chambEdges[0]<= seedEdges[0]) {
51 return (chambEdges[1] - seedEdges[0]) / (seedEdges[1] - seedEdges[0]);
54 else if (chambEdges[1] >= seedEdges[1]) {
55 return (seedEdges[1] - chambEdges[0]) / (seedEdges[1] - seedEdges[0]);
58 else if (seedEdges[0] <= chambEdges[0] && seedEdges[1] >= chambEdges[1]) {
59 return (chambEdges[1] - chambEdges[0]) / (seedEdges[1] - seedEdges[0]);
70 return StatusCode::SUCCESS;
80 ATH_CHECK(writeMaxima.
record(std::make_unique<EtaHoughMaxContainer>()));
93 for (
auto& [station, stationHoughBuckets] :
data.houghSetups) {
95 for (
auto& bucket : stationHoughBuckets) {
99 writeMaxima->push_back(std::make_unique<HoughMaximum>(std::move(
max)));
103 std::stable_sort(writeMaxima->begin(), writeMaxima->end(),
105 return (*a->parentBucket()) < (*b->parentBucket());
107 return StatusCode::SUCCESS;
119 std::vector<HoughSetupForBucket>& buckets =
data.houghSetups[bucket->front()->msSector()];
122 Amg::Vector3D leftSide = globToLoc.translation() - (
hs.bucket->coveredMin() * Amg::Vector3D::UnitY());
123 Amg::Vector3D rightSide = globToLoc.translation() - (
hs.bucket->coveredMax() * Amg::Vector3D::UnitY());
127 for (
const std::shared_ptr<MuonR4::SpacePoint> & sp : *bucket) {
136 hs.searchWindowTanAngle = {tanThetaLeft, tanThetaRight};
141 for (
const std::shared_ptr<MuonR4::SpacePoint> & hit : *bucket){
143 double y0l = hit->positionInChamber().y() - hit->positionInChamber().z() * tanThetaLeft;
144 double y0r = hit->positionInChamber().y() - hit->positionInChamber().z() * tanThetaRight;
153 switch (hit->
type()){
163 return meas->
channelType() == sTgcIdHelper::sTgcChannelTypes::Strip;
171 HoughPlaneConfig
cfg;
178 data.houghPlane = std::make_unique<HoughPlane>(
cfg);
179 data.peakFinder = std::make_unique<ActsPeakFinderForMuon>(peakFinderCfg);
185 int expectedPrecisionChambers{0}, seenPrecisionChambers{0};
186 bool hasTrig =
false;
188 std::unordered_set<const MuonGMR4::MuonReadoutElement*> seenChambers{};
189 std::set<std::pair<int,int>> seenLayers;
190 using enum Acts::TrapezoidVolumeBounds::BoundValues;
193 std::array<double, 2> tubeExtend{halfX, -halfX};
195 auto addSeenHit = [&tubeExtend, &seenLayers, &seenChambers](
const SpacePoint& sp,
197 const double sensorL,
198 const int mL,
const int layer){
200 seenChambers.insert(
re);
201 tubeExtend[0] =
std::min(tubeExtend[0], sp.positionInChamber().x() - sensorL);
202 tubeExtend[1] =
std::max(tubeExtend[1], sp.positionInChamber().x() + sensorL);
204 for (
const SpacePoint* SP : maximum.hitIdentifiers){
205 ATH_MSG_VERBOSE(__func__<<
"() - "<<__LINE__<<
" Maximum has associated hit in "
212 addSeenHit(*SP,
re, 0.5*
re->activeTubeLength(dc->measurementHash()),
213 re->multilayer(), dc->tubeLayer());
217 addSeenHit(*SP,
re, 0.5*
re->stripLayer(clust->measurementHash()).design().stripLength(clust->channelNumber()),
218 re->multilayer(), clust->gasGap());
222 addSeenHit(*SP,
re, 0.5*
re->stripLayer(clust->measurementHash()).design().stripLength(clust->channelNumber()),
223 re->multilayer(), clust->gasGap());
233 ATH_MSG_VERBOSE(__func__<<
"() - "<<__LINE__<<
" maximum does not cross "
234 <<
m_idHelperSvc->toStringDetEl(muonChamber.readoutEle()->identify()));
238 <<
m_idHelperSvc->toStringDetEl(muonChamber.readoutEle()->identify())<<
", "
239 <<(*muonChamber.bounds())<<
", @: "<<
Amg::toString(muonChamber.location()));
244 const bool hasHit = seenChambers.count(muonChamber.readoutEle());
250 ++seenPrecisionChambers;
255 }
else if (precTech) {
258 const double lowL = muonChamber.width(maximum.y + muonChamber.location().z() * maximum.x);
259 const std::array<double, 2> chambEdges{muonChamber.location().x() - lowL,
260 muonChamber.location().x() + lowL};
263 if (coverage < 0.95){
264 ATH_MSG_VERBOSE(__func__<<
"() - "<<__LINE__<<
" Reject chamber due to partial coverage: "<<coverage
265 <<
", chamb extend: "<<chambEdges[0]<<
"-"<<chambEdges[1]
266 <<
", tube extend: "<<tubeExtend[0]<<
"-"<<tubeExtend[1]
267 <<
", "<<(*muonChamber.readoutEle()->msSector()->bounds()));
273 ++expectedPrecisionChambers;
277 int minLayers = seenLayers.size() / 2 + 1;
279 int minSeenPrecisionChambers = (expectedPrecisionChambers > 1) + 1;
283 minSeenPrecisionChambers = 1;
286 ", seen prec: "<<seenPrecisionChambers<<
", required: "<<minSeenPrecisionChambers
287 <<
" -- layers: "<<seenLayers.size()<<
", required: "<<
minLayers);
288 return seenPrecisionChambers >= minSeenPrecisionChambers && (
int)seenLayers.size() >=
minLayers;
319 data.currAxisRanges = Acts::HoughTransformUtils::HoughAxisRanges{
320 searchStartTanTheta, searchEndTanTheta, searchStart, searchEnd};
322 data.houghPlane->reset();
323 for (
const SpacePointBucket::value_type& hit : *(bucket.
bucket)) {
326 auto maxima =
data.peakFinder->findPeaks(*(
data.houghPlane),
data.currAxisRanges);
329 "#eta Hough accumulator");
331 if (maxima.empty()) {
333 <<
":\n Mean tanTheta was "<<tanThetaMean
334 <<
" and my intercept "<<chamberCenter
337 <<
". The bucket found a search range of ("
342 <<
") , and my final search range is ["
343 <<searchStartTanTheta<<
" - "<<searchEndTanTheta
344 <<
"] and ["<<searchStart<<
" - "<<searchEnd
352 std::set<HoughHitType> seenHits;
355 for (
const auto&
max : maxima) {
358 unsigned int nPrec{0};
359 auto toBins = [&
data](
double x,
double y){
360 return std::make_pair(
365 auto accumulatorBins = toBins(
max.x,
max.y);
366 for (
const HoughHitType& hit :
data.houghPlane->hitIds(accumulatorBins.first, accumulatorBins.second)) {
367 auto res = seenHits.emplace(hit);
378 std::vector<HoughHitType> hitList{
max.hitIdentifiers.begin(),
max.hitIdentifiers.end()};
394 std::abs(eta), eta > 0?
'A' :
'C',
m_idHelperSvc->stationPhi(detId));
395 switch (
chamber.readoutEle()->detectorType()) {
404 primitives.push_back(
MuonValR4::drawLabel(
std::format(
"Missed seed - score {}, layer score {}, comprising {} measurements ",
data.houghPlane->nHits(accumulatorBins.first, accumulatorBins.second),
data.houghPlane->nLayers(accumulatorBins.first, accumulatorBins.second),hitList.size()),0.05,0.03,12));
406 primitivesForAcc.push_back(
MuonValR4::drawLabel(
std::format(
"Missed seed - score {}, layer score {}, comprising {} measurements ",
data.houghPlane->nHits(accumulatorBins.first, accumulatorBins.second),
data.houghPlane->nLayers(accumulatorBins.first, accumulatorBins.second),hitList.size()),0.05,0.03,12));
408 "MissedAccumulator", std::move(primitivesForAcc));
409 m_visionTool->visualizeSeed(ctx, seed,
"Missed seed",std::move(primitives));
415 size_t nHits = hitList.size();
420 std::ranges::stable_sort(hitList,
sorter);
432 primitives.push_back(
MuonValR4::drawLabel(
std::format(
"score {}, layer score {}, comprising {} measurements. wx = {:.2f}, wy = {:.1f} ",
data.houghPlane->nHits(accumulatorBins.first, accumulatorBins.second),
data.houghPlane->nLayers(accumulatorBins.first, accumulatorBins.second),hitList.size(),
max.wx,
max.wy),0.05,0.03,12));
433 primitivesForAcc.push_back(
MuonValR4::drawLabel(
std::format(
"score {}, layer score {}, comprising {} measurements. wx = {:.2f}, wy = {:.1f} ",
data.houghPlane->nHits(accumulatorBins.first, accumulatorBins.second),
data.houghPlane->nLayers(accumulatorBins.first, accumulatorBins.second),hitList.size(),
max.wx /
m_targetResoTanTheta,
max.wy /
m_targetResoIntercept),0.05,0.03,12));
435 m_visionTool->visualizeAccumulator(ctx, *
data.houghPlane,
data.currAxisRanges, {max},
"#eta Hough accumulator", std::move(primitivesForAcc));
436 m_visionTool->visualizeSeed(ctx, seed,
"#eta-HoughSeed", std::move(primitives));
442 using namespace std::placeholders;
457 const unsigned precisionLayerIndex = (dc->readoutElement()->multilayer() * 10 + dc->tubeLayer());
472 for (
const SpacePointBucket::value_type& hit : *bucket.
bucket) {
473 if (!hit->measuresEta()) {
474 hitList.push_back(hit.get());