diff --git a/Core/include/ACTS/Surfaces/TrapezoidBounds.hpp b/Core/include/ACTS/Surfaces/TrapezoidBounds.hpp
index 3453d8a395d5a1a876e3ff10e92dd6478e08868d..fa3e6e92040f4e4a822366a88ef042dbc6f005c9 100644
--- a/Core/include/ACTS/Surfaces/TrapezoidBounds.hpp
+++ b/Core/include/ACTS/Surfaces/TrapezoidBounds.hpp
@@ -50,14 +50,6 @@ public:
/// @param haley half length Y - defined at x=0
TrapezoidBounds(double minhalex, double maxhalex, double haley);
- /// Constructor for arbitrary Trapezoid
- ///
- /// @param minhalex minimal half lenght X, definition at negative halflength Y
- /// @param haley half length Y - defined at x=0
- /// @param alpha opening angle at @todo check
- /// @param beta opentin angle at @todo check
- TrapezoidBounds(double minhalex, double haley, double alpha, double beta);
-
/// Copy constructor
///
/// @param trabo are the source bounds for assignment
@@ -99,16 +91,6 @@ public:
double
halflengthY() const;
- /// This method returns the opening angle alpha in point A
- // (negative local phi)
- double
- alpha() const;
-
- /// This method returns the opening angle beta in point B
- /// (positive local phi)
- double
- beta() const;
-
/// The orientation of the Trapezoid is according to the figure above,
/// in words: the shorter of the two parallel sides of the trapezoid
/// intersects
@@ -205,54 +187,6 @@ public:
dump(std::ostream& sl) const final override;
private:
- /// private helper method for inside check
- ///
- ///
- /// @param lpos Local position (assumed to be in right surface frame)
- /// @param tol0 absulote tolerance parameter on the first coordinate
- /// @param tol1 absulote tolerance parameter on the second coordinate
- ///
- /// @return boolean indicator for the success of this operation
- bool
- inside(const Vector2D& lpos, double tol0, double tol1) const;
-
- /// private helper method inside() method for a full symmetric trapezoid
- ///
- /// @param lpos Local position (assumed to be in right surface frame)
- /// @param tol0 absulote tolerance parameter on the first coordinate
- /// @param tol1 absulote tolerance parameter on the second coordinate
- ///
- /// @return boolean indicator for the success of this operation
- bool
- insideFull(const Vector2D& lpos, double tol0 = 0., double tol1 = 0.) const;
-
- /// private inside() method for the triangular exclude
- /// area for an arbitrary trapezoid
- ///
- /// @param lpos Local position (assumed to be in right surface frame)
- /// @param tol0 absulote tolerance parameter on the first coordinate
- /// @param tol1 absulote tolerance parameter on the second coordinate
- ///
- /// @return boolean indicator for the success of this operation
- bool
- insideExclude(const Vector2D& lpos, double tol0 = 0., double tol1 = 0.) const;
-
- /// private isAbove() method for checking whether a point
- /// lies above or under a straight line
- ///
- /// @param lpos Local position (assumed to be in right surface frame)
- /// @param tol0 absulote tolerance parameter on the first coordinate
- /// @param tol1 absulote tolerance parameter on the second coordinate
- /// @param k is the first parameter of the parametric line equation
- /// @param d is the second parameter of the parameteric line equation
- ///
- /// @return boolean indicator for the success of this operation
- bool
- isAbove(const Vector2D& lpos, double tol0, double tol1, double k, double d)
- const;
-
- TDD_real_t m_alpha; ///< private cache of angle alpha
- TDD_real_t m_beta; ///< private cache of angle beta
RectangleBounds m_boundingBox; ///< internal bounding box cache
};
@@ -265,136 +199,35 @@ TrapezoidBounds::clone() const
inline double
TrapezoidBounds::minHalflengthX() const
{
- return m_valueStore.at(TrapezoidBounds::bv_minHalfX);
+ return m_valueStore[TrapezoidBounds::bv_minHalfX];
}
inline double
TrapezoidBounds::maxHalflengthX() const
{
- return m_valueStore.at(TrapezoidBounds::bv_maxHalfX);
+ return m_valueStore[TrapezoidBounds::bv_maxHalfX];
}
inline double
TrapezoidBounds::halflengthY() const
{
- return m_valueStore.at(TrapezoidBounds::bv_halfY);
-}
-
-inline double
-TrapezoidBounds::alpha() const
-{
- return m_alpha;
-}
-
-inline double
-TrapezoidBounds::beta() const
-{
- return m_beta;
+ return m_valueStore[TrapezoidBounds::bv_halfY];
}
inline bool
TrapezoidBounds::inside(const Vector2D& lpos, const BoundaryCheck& bcheck) const
{
- if (bcheck.bcType == 0)
- return inside(lpos, bcheck.toleranceLoc0, bcheck.toleranceLoc1);
-
- // a fast FALSE
- double fabsY = std::abs(lpos[Acts::eLOC_Y]);
- double max_ell = bcheck.lCovariance(0, 0) > bcheck.lCovariance(1, 1)
- ? bcheck.lCovariance(0, 0)
- : bcheck.lCovariance(1, 1);
- double limit = bcheck.nSigmas * sqrt(max_ell);
- if (fabsY > (m_valueStore.at(TrapezoidBounds::bv_halfY) + limit))
- return false;
- // a fast FALSE
- double fabsX = std::abs(lpos[Acts::eLOC_X]);
- if (fabsX > (m_valueStore.at(TrapezoidBounds::bv_maxHalfX) + limit))
- return false;
- // a fast TRUE
- double min_ell = bcheck.lCovariance(0, 0) < bcheck.lCovariance(1, 1)
- ? bcheck.lCovariance(0, 0)
- : bcheck.lCovariance(1, 1);
- limit = bcheck.nSigmas * sqrt(min_ell);
- if (fabsX < (m_valueStore.at(TrapezoidBounds::bv_minHalfX) + limit)
- && fabsY < (m_valueStore.at(TrapezoidBounds::bv_halfY) + limit))
- return true;
-
- // compute KDOP and axes for surface polygon
- std::vector<KDOP> elementKDOP(3);
- std::vector<Vector2D> elementP(4);
- float theta
- = (bcheck.lCovariance(1, 0) != 0
- && (bcheck.lCovariance(1, 1) - bcheck.lCovariance(0, 0)) != 0)
- ? .5 * std::atan(
- 2 * bcheck.lCovariance(1, 0)
- / (bcheck.lCovariance(1, 1) - bcheck.lCovariance(0, 0)))
- : 0.;
- auto rotMatrix = Eigen::Rotation2D<double>(theta).toRotationMatrix();
- ActsMatrixD<2, 2> normal;
- normal << 0, -1, 1, 0;
- // ellipse is always at (0,0), surface is moved to ellipse position and then
- // rotated
- Vector2D p;
- p << m_valueStore.at(TrapezoidBounds::bv_minHalfX),
- -m_valueStore.at(TrapezoidBounds::bv_halfY);
- elementP.at(0) = (rotMatrix * (p - lpos));
- p << -m_valueStore.at(TrapezoidBounds::bv_minHalfX),
- -m_valueStore.at(TrapezoidBounds::bv_halfY);
- elementP.at(1) = (rotMatrix * (p - lpos));
- p << m_valueStore.at(TrapezoidBounds::bv_minHalfX)
- + (2. * m_valueStore.at(TrapezoidBounds::bv_halfY))
- * std::tan(m_beta),
- m_valueStore.at(TrapezoidBounds::bv_halfY);
- elementP.at(2) = (rotMatrix * (p - lpos));
- p << -(m_valueStore.at(TrapezoidBounds::bv_minHalfX)
- + (2. * m_valueStore[TrapezoidBounds::bv_halfY]) * std::tan(m_alpha)),
- m_valueStore.at(TrapezoidBounds::bv_halfY);
- elementP.at(3) = (rotMatrix * (p - lpos));
- std::vector<Vector2D> axis = {normal * (elementP.at(1) - elementP.at(0)),
- normal * (elementP.at(3) - elementP.at(1)),
- normal * (elementP.at(2) - elementP.at(0))};
- bcheck.ComputeKDOP(elementP, axis, elementKDOP);
- // compute KDOP for error ellipse
- std::vector<KDOP> errelipseKDOP(3);
- bcheck.ComputeKDOP(bcheck.EllipseToPoly(3), axis, errelipseKDOP);
- // check if KDOPs overlap and return result
- return bcheck.TestKDOPKDOP(elementKDOP, errelipseKDOP);
-}
-
-inline bool
-TrapezoidBounds::insideLoc0(const Vector2D& lpos, double tol0) const
-{
- return (std::abs(lpos[Acts::eLOC_X])
- < m_valueStore.at(TrapezoidBounds::bv_maxHalfX) + tol0);
-}
-
-inline bool
-TrapezoidBounds::insideLoc1(const Vector2D& lpos, double tol1) const
-{
- return (std::abs(lpos[Acts::eLOC_Y])
- < m_valueStore.at(TrapezoidBounds::bv_halfY) + tol1);
+ return bcheck.isInsidePolygon(lpos, vertices());
}
inline const std::vector<Vector2D>
TrapezoidBounds::vertices() const
{
- // create the return vector
- std::vector<Vector2D> vertices;
- // fill the vertices
- vertices.reserve(4);
- vertices.push_back(
- Vector2D(m_valueStore.at(TrapezoidBounds::bv_minHalfX),
- -m_valueStore.at(TrapezoidBounds::bv_halfY))); // [0]
- vertices.push_back(
- Vector2D(m_valueStore.at(TrapezoidBounds::bv_maxHalfX),
- m_valueStore.at(TrapezoidBounds::bv_halfY))); // [1]
- vertices.push_back(
- Vector2D(-m_valueStore.at(TrapezoidBounds::bv_maxHalfX),
- m_valueStore.at(TrapezoidBounds::bv_halfY))); // [1]
- vertices.push_back(
- Vector2D(-m_valueStore.at(TrapezoidBounds::bv_minHalfX),
- -m_valueStore.at(TrapezoidBounds::bv_halfY))); // [3]
- return vertices;
+ // counter-clockwise from bottom-right corner
+ return {{minHalflengthX(), -halflengthY()},
+ {maxHalflengthX(), halflengthY()},
+ {-maxHalflengthX(), halflengthY()},
+ {-minHalflengthX(), -halflengthY()}};
}
inline const RectangleBounds&
diff --git a/Core/src/Surfaces/TrapezoidBounds.cpp b/Core/src/Surfaces/TrapezoidBounds.cpp
index e9bb3a425bea456ffbcb7ab5378ee57bc4cc2bca..100d2548c7d83defd76cd437083051dc90325ef4 100644
--- a/Core/src/Surfaces/TrapezoidBounds.cpp
+++ b/Core/src/Surfaces/TrapezoidBounds.cpp
@@ -18,45 +18,12 @@
Acts::TrapezoidBounds::TrapezoidBounds(double minhalex,
double maxhalex,
double haley)
- : PlanarBounds(TrapezoidBounds::bv_length)
- , m_alpha(0.)
- , m_beta(0.)
- , m_boundingBox(0., 0.)
+ : PlanarBounds(TrapezoidBounds::bv_length), m_boundingBox(0.0, 0.0)
{
- m_valueStore.at(TrapezoidBounds::bv_minHalfX) = std::abs(minhalex);
- m_valueStore.at(TrapezoidBounds::bv_maxHalfX) = std::abs(maxhalex);
- m_valueStore.at(TrapezoidBounds::bv_halfY) = std::abs(haley);
- // find the maximum at for the bounding box
- double mx = m_valueStore.at(TrapezoidBounds::bv_minHalfX)
- > m_valueStore.at(TrapezoidBounds::bv_maxHalfX)
- ? m_valueStore.at(TrapezoidBounds::bv_minHalfX)
- : m_valueStore.at(TrapezoidBounds::bv_maxHalfX);
- m_boundingBox
- = RectangleBounds(mx, m_valueStore.at(TrapezoidBounds::bv_halfY));
-}
-
-Acts::TrapezoidBounds::TrapezoidBounds(double minhalex,
- double haley,
- double alpha,
- double beta)
- : PlanarBounds(TrapezoidBounds::bv_length)
- , m_alpha(alpha)
- , m_beta(beta)
- , m_boundingBox(0., 0.)
-{
- double gamma = (alpha > beta) ? (alpha - 0.5 * M_PI) : (beta - 0.5 * M_PI);
- // now fill them
- m_valueStore.at(TrapezoidBounds::bv_minHalfX) = std::abs(minhalex);
- m_valueStore.at(TrapezoidBounds::bv_maxHalfX) = minhalex
- + (2. * m_valueStore.at(TrapezoidBounds::bv_halfY)) * tan(gamma);
- m_valueStore.at(TrapezoidBounds::bv_halfY) = std::abs(haley);
- // find the maximum for the bounding box
- double mx = m_valueStore.at(TrapezoidBounds::bv_minHalfX)
- > m_valueStore.at(TrapezoidBounds::bv_maxHalfX)
- ? m_valueStore.at(TrapezoidBounds::bv_minHalfX)
- : m_valueStore.at(TrapezoidBounds::bv_maxHalfX);
- m_boundingBox
- = RectangleBounds(mx, m_valueStore.at(TrapezoidBounds::bv_halfY));
+ m_valueStore[TrapezoidBounds::bv_minHalfX] = std::abs(minhalex);
+ m_valueStore[TrapezoidBounds::bv_maxHalfX] = std::abs(maxhalex);
+ m_valueStore[TrapezoidBounds::bv_halfY] = std::abs(haley);
+ m_boundingBox = RectangleBounds(std::max(minhalex, maxhalex), haley);
}
Acts::TrapezoidBounds::~TrapezoidBounds()
@@ -68,140 +35,15 @@ Acts::TrapezoidBounds::operator=(const TrapezoidBounds& trabo)
{
if (this != &trabo) {
PlanarBounds::operator=(trabo);
- m_alpha = trabo.m_alpha;
- m_beta = trabo.m_beta;
m_boundingBox = trabo.m_boundingBox;
}
return *this;
}
-bool
-Acts::TrapezoidBounds::inside(const Acts::Vector2D& lpos,
- double tol0,
- double tol1) const
-{
- if (m_alpha == 0.) return insideFull(lpos, tol0, tol1);
- return (insideFull(lpos, tol0, tol1) && !insideExclude(lpos, tol0, tol1));
-}
-
-bool
-Acts::TrapezoidBounds::insideFull(const Acts::Vector2D& lpos,
- double tol0,
- double tol1) const
-{
- // the cases:
- double absX = std::abs(lpos[Acts::eLOC_X]);
- double absY = std::abs(lpos[Acts::eLOC_Y]);
- // (1) a fast FALSE
- if (absY > (m_valueStore.at(TrapezoidBounds::bv_halfY) + tol1)) return false;
- // (2) a fast FALSE
- if (absX > (m_valueStore.at(TrapezoidBounds::bv_maxHalfX) + tol0))
- return false;
- // (3) a fast TRUE
- if (absX < (m_valueStore.at(TrapezoidBounds::bv_minHalfX) - tol0))
- return true;
- // (4) particular case - a rectangle
- if (m_valueStore.at(TrapezoidBounds::bv_maxHalfX)
- == m_valueStore.at(TrapezoidBounds::bv_minHalfX))
- return true;
- // (5) /** use basic calculation of a straight line */
- double k = 2.0 * m_valueStore.at(TrapezoidBounds::bv_halfY)
- / (m_valueStore.at(TrapezoidBounds::bv_maxHalfX)
- - m_valueStore.at(TrapezoidBounds::bv_minHalfX))
- * ((lpos[Acts::eLOC_X] > 0.) ? 1.0 : -1.0);
- double d
- = -std::abs(k) * 0.5 * (m_valueStore.at(TrapezoidBounds::bv_maxHalfX)
- + m_valueStore.at(TrapezoidBounds::bv_minHalfX));
- return (isAbove(lpos, tol0, tol1, k, d));
-}
-
-bool
-Acts::TrapezoidBounds::insideExclude(const Acts::Vector2D& lpos,
- double tol0,
- double tol1) const
-{
- // line a
- bool alphaBiggerBeta(m_alpha > m_beta);
- double ka = alphaBiggerBeta ? tan(M_PI - m_alpha) : tan(m_alpha);
- double kb = alphaBiggerBeta ? tan(M_PI - m_beta) : tan(m_beta);
- double sign = alphaBiggerBeta ? -1. : 1.;
- double da = -m_valueStore.at(TrapezoidBounds::bv_halfY)
- + sign * ka * m_valueStore.at(TrapezoidBounds::bv_minHalfX);
- double db = -m_valueStore.at(TrapezoidBounds::bv_halfY)
- + sign * kb * m_valueStore.at(TrapezoidBounds::bv_minHalfX);
-
- return (isAbove(lpos, tol0, tol1, ka, da)
- && isAbove(lpos, tol0, tol1, kb, db));
-}
-
-bool
-Acts::TrapezoidBounds::isAbove(const Acts::Vector2D& lpos,
- double tol0,
- double tol1,
- double k,
- double d) const
-{
- // the most tolerant approach for tol0 and tol1
- double sign = k > 0. ? -1. : +1.;
- return (lpos[Acts::eLOC_Y] + tol1
- > (k * (lpos[Acts::eLOC_X] + sign * tol0) + d));
-}
-
double
Acts::TrapezoidBounds::distanceToBoundary(const Acts::Vector2D& lpos) const
{
- const int Np = 4;
-
- double xl = -m_valueStore.at(TrapezoidBounds::bv_maxHalfX);
- double xr = m_valueStore.at(TrapezoidBounds::bv_maxHalfX);
- if (m_alpha != 0.) {
- xl = -m_valueStore.at(TrapezoidBounds::bv_minHalfX)
- - 2. * tan(m_alpha) * m_valueStore.at(TrapezoidBounds::bv_halfY);
- } else if (m_beta != 0.) {
- xr = m_valueStore.at(TrapezoidBounds::bv_minHalfX)
- + 2. * tan(m_beta) * m_valueStore.at(TrapezoidBounds::bv_halfY);
- }
- double X[4] = {-m_valueStore.at(TrapezoidBounds::bv_minHalfX),
- xl,
- xr,
- m_valueStore.at(TrapezoidBounds::bv_minHalfX)};
- double Y[4] = {-m_valueStore.at(TrapezoidBounds::bv_halfY),
- m_valueStore.at(TrapezoidBounds::bv_halfY),
- m_valueStore.at(TrapezoidBounds::bv_halfY),
- -m_valueStore.at(TrapezoidBounds::bv_halfY)};
-
- double dm = 1.e+20;
- double Ao = 0.;
- bool in = true;
-
- for (int i = 0; i != Np; ++i) {
- int j = i + 1;
- if (j == Np) j = 0;
-
- double x = X[i] - lpos[0];
- double y = Y[i] - lpos[1];
- double dx = X[j] - X[i];
- double dy = Y[j] - Y[i];
- double A = x * dy - y * dx;
- double S = -(x * dx + y * dy);
-
- if (S <= 0.) {
- double d = x * x + y * y;
- if (d < dm) dm = d;
- } else {
- double a = dx * dx + dy * dy;
- if (S <= a) {
- double d = (A * A) / a;
- if (d < dm) dm = d;
- }
- }
- if (i && in && Ao * A < 0.) in = false;
- Ao = A;
- }
- if (in)
- return -sqrt(dm);
- else
- return sqrt(dm);
+ return BoundaryCheck(true).distanceToBoundary(lpos, vertices());
}
std::ostream&
diff --git a/Tests/Surfaces/TrapezoidBoundsTests.cpp b/Tests/Surfaces/TrapezoidBoundsTests.cpp
index 56a36c99577155d12e35985ea0cae1fb13e302f4..2f69d7a2f15bbcb4210c84b43baa5d2a4029477d 100644
--- a/Tests/Surfaces/TrapezoidBoundsTests.cpp
+++ b/Tests/Surfaces/TrapezoidBoundsTests.cpp
@@ -36,7 +36,6 @@ namespace Test {
BOOST_AUTO_TEST_CASE(TrapezoidBoundsConstruction)
{
double minHalfX(1.), maxHalfX(6.), halfY(2.);
- double alpha(std::atan(5. / 4.)), beta(std::atan(1.));
//
// default construction deleted
// TrapezoidBounds defaultConstructedTrapezoidBounds;
@@ -44,10 +43,6 @@ namespace Test {
/// Test construction with defining half lengths
BOOST_TEST(TrapezoidBounds(minHalfX, maxHalfX, halfY).type()
== SurfaceBounds::Trapezoid);
- //
- /// Test construction with lengths and angles
- BOOST_TEST(TrapezoidBounds(minHalfX, halfY, alpha, beta).type()
- == SurfaceBounds::Trapezoid);
/// Copy constructor
TrapezoidBounds original(minHalfX, maxHalfX, halfY);
TrapezoidBounds copied(original);
@@ -58,11 +53,8 @@ namespace Test {
BOOST_AUTO_TEST_CASE(TrapezoidBoundsProperties, *utf::expected_failures(3))
{
double minHalfX(1.), maxHalfX(6.), halfY(2.);
- double alpha(std::atan(5. / 4.)), beta(std::atan(1.));
//
TrapezoidBounds trapezoidBoundsObject(minHalfX, maxHalfX, halfY);
- TrapezoidBounds asymmetricTrapezoidBoundsObject(
- minHalfX, halfY, alpha, beta);
/// Test clone
auto pClonedTrapezoidBounds = trapezoidBoundsObject.clone();
BOOST_TEST(bool(pClonedTrapezoidBounds));
@@ -77,24 +69,9 @@ namespace Test {
/// Test maxHalfLengthX
BOOST_TEST(trapezoidBoundsObject.maxHalflengthX() == maxHalfX);
//
- /// Test maxHalfLengthX for asymmetric trapezoid, fails
- BOOST_TEST(asymmetricTrapezoidBoundsObject.maxHalflengthX() == maxHalfX);
- //
/// Test halflengthY
BOOST_TEST(trapezoidBoundsObject.halflengthY() == halfY);
//
- /// Test alpha for symmetric trapezoid, will fail
- BOOST_TEST(trapezoidBoundsObject.alpha() == std::atan(5. / 4.));
- //
- /// Test alpha for asymmetric trapezoid
- BOOST_TEST(asymmetricTrapezoidBoundsObject.alpha() == alpha);
- //
- /// Test beta for symmetric trapezoid, will fail
- BOOST_TEST(trapezoidBoundsObject.beta() == std::atan(5. / 4.));
- //
- /// Test beta for asymmetric trapezoid
- BOOST_TEST(asymmetricTrapezoidBoundsObject.beta() == beta);
- //
/// Test distanceToBoundary
Vector2D origin(0., 0.);
Vector2D outside(30., 0.);
@@ -147,7 +124,7 @@ namespace Test {
//
/// Test assignment
TrapezoidBounds assignedTrapezoidBoundsObject(
- NaN, NaN, NaN, NaN); // invalid object, in some sense
+ NaN, NaN, NaN); // invalid object, in some sense
assignedTrapezoidBoundsObject = trapezoidBoundsObject;
BOOST_TEST(assignedTrapezoidBoundsObject == trapezoidBoundsObject);
}