src/Geometry/Extent.cpp

0001 // This file is part of the ACTS project.
0002 //
0003 // Copyright (C) 2016 CERN for the benefit of the ACTS project
0004 //
0005 // This Source Code Form is subject to the terms of the Mozilla Public
0006 // License, v. 2.0. If a copy of the MPL was not distributed with this
0007 // file, You can obtain one at https://mozilla.org/MPL/2.0/.
0008
0009 #include "Acts/Geometry/Extent.hpp"
0010
0011 #include "Acts/Utilities/VectorHelpers.hpp"
0012
0013 #include <algorithm>
0014 #include <cstddef>
0015 #include <iomanip>
0016 #include <limits>
0017 #include <numbers>
0018
0019 namespace Acts {
0020
0021 Extent::Extent(const ExtentEnvelope& envelope)
0022     : m_constrains(0), m_envelope(envelope) {
0023   m_range[toUnderlying(AxisDirection::AxisR)] =
0024       Range1D<double>(0., std::numeric_limits<double>::max());
0025   m_range[toUnderlying(AxisDirection::AxisPhi)] =
0026       Range1D<double>(-std::numbers::pi, std::numbers::pi);
0027   m_range[toUnderlying(AxisDirection::AxisRPhi)] =
0028       Range1D<double>(0., std::numeric_limits<double>::max());
0029   m_range[toUnderlying(AxisDirection::AxisMag)] =
0030       Range1D<double>(0., std::numeric_limits<double>::max());
0031 }
0032
0033 void Extent::extend(const Vector3& vtx, const std::vector<AxisDirection>& aDirs,
0034                     bool applyEnv, bool fillHistograms) {
0035   for (auto aDir : aDirs) {
0036     // Get the casted value given the binning value description
0037     double cValue = VectorHelpers::cast(vtx, aDir);
0038     if (fillHistograms) {
0039       m_valueHistograms[toUnderlying(aDir)].push_back(cValue);
0040     }
0041     // Apply envelope as suggested
0042     double lEnv = applyEnv ? m_envelope[aDir][0] : 0.;
0043     double hEnv = applyEnv ? m_envelope[aDir][1] : 0.;
0044     double mValue = cValue - lEnv;
0045     // Special protection for radial value
0046     if (aDir == AxisDirection::AxisR && mValue < 0.) {
0047       mValue = std::max(mValue, 0.);
0048     }
0049     if (constrains(aDir)) {
0050       m_range[toUnderlying(aDir)].expand(mValue, cValue + hEnv);
0051     } else {
0052       m_range[toUnderlying(aDir)].shrink(mValue, cValue + hEnv);
0053     }
0054     m_constrains.set(toUnderlying(aDir));
0055   }
0056 }
0057
0058 void Extent::extend(const Extent& rhs, const std::vector<AxisDirection>& aDirs,
0059                     bool applyEnv) {
0060   for (auto aDir : aDirs) {
0061     // The value is constraint, envelope can be optional
0062     if (rhs.constrains(aDir)) {
0063       double lEnv = applyEnv ? m_envelope[aDir][0] : 0.;
0064       double hEnv = applyEnv ? m_envelope[aDir][1] : 0.;
0065       if (constrains(aDir)) {
0066         m_range[toUnderlying(aDir)].expand(
0067             rhs.range()[toUnderlying(aDir)].min() - lEnv,
0068             rhs.range()[toUnderlying(aDir)].max() + hEnv);
0069       } else {
0070         m_range[toUnderlying(aDir)].shrink(
0071             rhs.range()[toUnderlying(aDir)].min() - lEnv,
0072             rhs.range()[toUnderlying(aDir)].max() + hEnv);
0073       }
0074       m_constrains.set(toUnderlying(aDir));
0075     } else if (rhs.envelope()[aDir] != zeroEnvelope) {
0076       // Only an envelope given, but value is not constraint -> apply envelope
0077       m_range[toUnderlying(aDir)].expand(
0078           m_range[toUnderlying(aDir)].min() - rhs.envelope()[aDir][0],
0079           m_range[toUnderlying(aDir)].max() + rhs.envelope()[aDir][1]);
0080       m_constrains.set(toUnderlying(aDir));
0081     }
0082   }
0083 }
0084
0085 void Extent::addConstrain(const Extent& rhs, const ExtentEnvelope& envelope) {
0086   for (const auto& aDir : allAxisDirections()) {
0087     if (rhs.constrains(aDir) && !constrains(aDir)) {
0088       const auto& cRange = rhs.range(aDir);
0089       m_range[toUnderlying(aDir)].setMin(cRange.min() - envelope[aDir][0u]);
0090       m_range[toUnderlying(aDir)].setMax(cRange.max() + envelope[aDir][1u]);
0091       m_constrains.set(toUnderlying(aDir));
0092     }
0093   }
0094 }
0095
0096 void Extent::set(AxisDirection aDir, double min, double max) {
0097   double minval = min;
0098   if (aDir == AxisDirection::AxisR && minval < 0.) {
0099     minval = 0.;
0100   }
0101   m_range[toUnderlying(aDir)] = Range1D<double>{minval, max};
0102   m_constrains.set(toUnderlying(aDir));
0103 }
0104
0105 void Extent::setMin(AxisDirection aDir, double min) {
0106   double minval = min;
0107   if (aDir == AxisDirection::AxisR && minval < 0.) {
0108     minval = 0.;
0109   }
0110   m_range[toUnderlying(aDir)].setMin(0u, minval);
0111   m_constrains.set(toUnderlying(aDir));
0112 }
0113
0114 void Extent::setMax(AxisDirection aDir, double max) {
0115   m_range[toUnderlying(aDir)].setMax(0u, max);
0116   m_constrains.set(toUnderlying(aDir));
0117 }
0118
0119 void Extent::setEnvelope(const ExtentEnvelope& envelope) {
0120   m_envelope = envelope;
0121 }
0122
0123 bool Extent::contains(const Vector3& vtx) const {
0124   Extent checkExtent;
0125   for (const auto& bv : allAxisDirections()) {
0126     if (constrains(bv)) {
0127       double vtxVal = VectorHelpers::cast(vtx, bv);
0128       checkExtent.set(bv, vtxVal, vtxVal);
0129     }
0130   }
0131   return contains(checkExtent);
0132 }
0133
0134 bool Extent::contains(const Extent& rhs,
0135                       std::optional<AxisDirection> aDir) const {
0136   // Helper to check including a constraint bit set check
0137   auto checkContainment = [&](AxisDirection bvc) -> bool {
0138     if (!constrains(bvc)) {
0139       return true;
0140     }
0141     return (rhs.range()[toUnderlying(bvc)] <= m_range[toUnderlying(bvc)]);
0142   };
0143
0144   // Check all
0145   if (!aDir.has_value()) {
0146     return std::ranges::all_of(allAxisDirections(), checkContainment);
0147   }
0148   // Check specific
0149   return checkContainment(aDir.value());
0150 }
0151
0152 bool Extent::intersects(const Extent& rhs,
0153                         std::optional<AxisDirection> aDir) const {
0154   // Helper to check including a constraint bit set check
0155   auto checkIntersect = [&](AxisDirection bvc) -> bool {
0156     if (!constrains(bvc) || !rhs.constrains(bvc)) {
0157       return false;
0158     }
0159     return (m_range[toUnderlying(bvc)] && rhs.range()[toUnderlying(bvc)]);
0160   };
0161
0162   // Check all
0163   if (!aDir.has_value()) {
0164     return std::ranges::any_of(allAxisDirections(), checkIntersect);
0165   }
0166   // Check specific
0167   return checkIntersect(aDir.value());
0168 }
0169
0170 bool Extent::constrains(AxisDirection aDir) const {
0171   return m_constrains.test(static_cast<std::size_t>(aDir));
0172 }
0173
0174 bool Extent::constrains() const {
0175   return m_constrains.count() > 0;
0176 }
0177
0178 bool Extent::operator==(const Extent& e) const {
0179   if (m_constrains != e.m_constrains) {
0180     return false;
0181   }
0182   if (m_envelope != e.m_envelope) {
0183     return false;
0184   }
0185   if (!(m_range == e.m_range)) {
0186     return false;
0187   }
0188   if (m_valueHistograms != e.m_valueHistograms) {
0189     return false;
0190   }
0191   return true;
0192 }
0193
0194 std::string Extent::toString(const std::string& indent) const {
0195   std::stringstream sl;
0196   sl << indent << "Extent in space :" << std::endl;
0197   for (const auto& bv : allAxisDirections()) {
0198     if (constrains(bv)) {
0199       sl << indent << "  - value :" << std::setw(10) << axisDirectionName(bv)
0200          << " | range = [" << m_range[toUnderlying(bv)].min() << ", "
0201          << m_range[toUnderlying(bv)].max() << "]" << std::endl;
0202     }
0203   }
0204   return sl.str();
0205 }
0206
0207 // Overload of << operator for std::ostream for debug output
0208 std::ostream& operator<<(std::ostream& sl, const Extent& rhs) {
0209   sl << rhs.toString();
0210   return sl;
0211 }
0212
0213 }  // namespace Acts