src/Surfaces/TrapezoidBounds.cpp

0001 // This file is part of the Acts project.
0002 //
0003 // Copyright (C) 2016-2020 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 http://mozilla.org/MPL/2.0/.
0008
0009 #include "Acts/Surfaces/TrapezoidBounds.hpp"
0010
0011 #include "Acts/Definitions/TrackParametrization.hpp"
0012 #include "Acts/Surfaces/ConvexPolygonBounds.hpp"
0013
0014 #include <iomanip>
0015 #include <iostream>
0016
0017 /// Constructor for symmetric Trapezoid
0018 ///
0019 /// @param halfXnegY minimal half length X, definition at negative Y
0020 /// @param halfXposY maximal half length X, definition at positive Y
0021 /// @param halfY half length Y - defined at x=0
0022 /// @param rotAngle: rotation angle of the bounds w.r.t coordinate axes
0023 Acts::TrapezoidBounds::TrapezoidBounds(double halfXnegY, double halfXposY,
0024                                        double halfY,
0025                                        double rotAngle) noexcept(false)
0026     : m_values({halfXnegY, halfXposY, halfY, rotAngle}),
0027       m_boundingBox(std::max(halfXnegY, halfXposY), halfY) {
0028   rotateBoundingBox();
0029   checkConsistency();
0030 }
0031
0032 /// Constructor for symmetric Trapezoid - from fixed size array
0033 ///
0034 /// @param values the values to be stream in
0035 Acts::TrapezoidBounds::TrapezoidBounds(
0036     const std::array<double, eSize>& values) noexcept(false)
0037     : m_values(values),
0038       m_boundingBox(
0039           std::max(values[eHalfLengthXnegY], values[eHalfLengthXposY]),
0040           values[eHalfLengthY]) {
0041   rotateBoundingBox();
0042   checkConsistency();
0043 }
0044
0045 Acts::TrapezoidBounds::~TrapezoidBounds() = default;
0046
0047 Acts::SurfaceBounds::BoundsType Acts::TrapezoidBounds::type() const {
0048   return SurfaceBounds::eTrapezoid;
0049 }
0050
0051 bool Acts::TrapezoidBounds::inside(const Acts::Vector2& lposition,
0052                                    const Acts::BoundaryCheck& bcheck) const {
0053   const double hlXnY = get(TrapezoidBounds::eHalfLengthXnegY);
0054   const double hlXpY = get(TrapezoidBounds::eHalfLengthXposY);
0055   const double hlY = get(TrapezoidBounds::eHalfLengthY);
0056   const double rotAngle = get(TrapezoidBounds::eRotationAngle);
0057
0058   const Acts::Vector2 extPosition = Eigen::Rotation2Dd(rotAngle) * lposition;
0059   const double x = extPosition[0];
0060   const double y = extPosition[1];
0061
0062   if (bcheck.type() == BoundaryCheck::Type::eAbsolute) {
0063     const double tolX = bcheck.tolerance()[eBoundLoc0];
0064     const double tolY = bcheck.tolerance()[eBoundLoc1];
0065
0066     if (std::abs(y) - hlY > tolY) {
0067       // outside y range
0068       return false;
0069     }
0070
0071     if (std::abs(x) - std::max(hlXnY, hlXpY) > tolX) {
0072       // outside x range
0073       return false;
0074     }
0075
0076     if (std::abs(x) - std::min(hlXnY, hlXpY) <= tolX) {
0077       // inside x range
0078       return true;
0079     }
0080   }
0081
0082   // at this stage, the point can only be in the triangles
0083   // run slow-ish polygon check
0084   std::vector<Acts::Vector2> vertices = {
0085       {-hlXnY, -hlY}, {hlXnY, -hlY}, {hlXpY, hlY}, {-hlXpY, hlY}};
0086   return bcheck.isInside(extPosition, vertices);
0087 }
0088
0089 std::vector<Acts::Vector2> Acts::TrapezoidBounds::vertices(
0090     unsigned int /*lseg*/) const {
0091   const double hlXnY = get(TrapezoidBounds::eHalfLengthXnegY);
0092   const double hlXpY = get(TrapezoidBounds::eHalfLengthXposY);
0093   const double hlY = get(TrapezoidBounds::eHalfLengthY);
0094   const double rotAngle = get(TrapezoidBounds::eRotationAngle);
0095
0096   std::vector<Acts::Vector2> vertices = {
0097       {-hlXnY, -hlY}, {hlXnY, -hlY}, {hlXpY, hlY}, {-hlXpY, hlY}};
0098   for (auto& v : vertices) {
0099     v = Eigen::Rotation2Dd(-rotAngle) * v;
0100   }
0101   return vertices;
0102 }
0103
0104 const Acts::RectangleBounds& Acts::TrapezoidBounds::boundingBox() const {
0105   return m_boundingBox;
0106 }
0107
0108 std::ostream& Acts::TrapezoidBounds::toStream(std::ostream& sl) const {
0109   sl << std::setiosflags(std::ios::fixed);
0110   sl << std::setprecision(7);
0111   sl << "Acts::TrapezoidBounds:  (halfXnegY, halfXposY, halfY, rotAngle) = "
0112      << "(" << get(eHalfLengthXnegY) << ", " << get(eHalfLengthXposY) << ", "
0113      << get(eHalfLengthY) << ", " << get(eRotationAngle) << ")";
0114   sl << std::setprecision(-1);
0115   return sl;
0116 }
0117
0118 std::vector<double> Acts::TrapezoidBounds::values() const {
0119   std::vector<double> valvector;
0120   valvector.insert(valvector.begin(), m_values.begin(), m_values.end());
0121   return valvector;
0122 }
0123
0124 void Acts::TrapezoidBounds::rotateBoundingBox() noexcept(false) {
0125   const double rotAngle = get(eRotationAngle);
0126
0127   if (rotAngle != 0.) {
0128     m_boundingBox = ConvexPolygonBounds<4>(vertices()).boundingBox();
0129   }
0130 }
0131
0132 void Acts::TrapezoidBounds::checkConsistency() noexcept(false) {
0133   if (get(eHalfLengthXnegY) <= 0. || get(eHalfLengthXposY) <= 0.) {
0134     throw std::invalid_argument("TrapezoidBounds: invalid local x setup");
0135   }
0136   if (get(eHalfLengthY) <= 0.) {
0137     throw std::invalid_argument("TrapezoidBounds: invalid local y setup");
0138   }
0139 }