src/Geometry/Polyhedron.cpp

0001 // This file is part of the Acts project.
0002 //
0003 // Copyright (C) 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/Geometry/Polyhedron.hpp"
0010
0011 #include "Acts/Surfaces/detail/VerticesHelper.hpp"
0012 #include "Acts/Utilities/BinningType.hpp"
0013
0014 #include <algorithm>
0015 #include <cmath>
0016 #include <limits>
0017 #include <utility>
0018
0019 void Acts::Polyhedron::merge(const Acts::Polyhedron& other) {
0020   std::size_t cvert = vertices.size();
0021   vertices.insert(vertices.end(), other.vertices.begin(), other.vertices.end());
0022   /// Add the new faces with offsets
0023   auto join = [&](std::vector<FaceType>& existing,
0024                   const std::vector<FaceType>& additional) -> void {
0025     for (const auto& aface : additional) {
0026       FaceType nface = aface;
0027       std::transform(nface.begin(), nface.end(), nface.begin(),
0028                      [&](std::size_t x) { return (x + cvert); });
0029       existing.push_back(nface);
0030     }
0031   };
0032   // For faces and triangular mesh
0033   join(faces, other.faces);
0034   join(triangularMesh, other.triangularMesh);
0035 }
0036
0037 void Acts::Polyhedron::move(const Transform3& transform) {
0038   for_each(vertices.begin(), vertices.end(),
0039            [&](auto& v) { v = transform * v; });
0040 }
0041
0042 Acts::Extent Acts::Polyhedron::extent(const Transform3& transform) const {
0043   Extent extent;
0044   auto vtxs = vertices;
0045   std::transform(vtxs.begin(), vtxs.end(), vtxs.begin(), [&](auto& v) {
0046     auto vt = (transform * v);
0047     extent.extend(vt);
0048     return (vt);
0049   });
0050
0051   // Special checks of binR for hyper plane surfaces
0052   if (detail::VerticesHelper::onHyperPlane(vtxs)) {
0053     // Check inclusion of origin (i.e. convex around origin)
0054     Vector3 origin = transform * Vector3(0., 0., extent.medium(binZ));
0055     for (const auto& face : faces) {
0056       std::vector<Vector3> tface;
0057       tface.reserve(face.size());
0058       for (auto f : face) {
0059         tface.push_back(vtxs[f]);
0060       }
0061       if (detail::VerticesHelper::isInsidePolygon(origin, tface)) {
0062         extent.range(binR).setMin(0.);
0063         extent.range(binPhi).set(-M_PI, M_PI);
0064         break;
0065       }
0066     }
0067     if (exact) {
0068       // Check for radial extend in 2D
0069       auto radialDistance = [&](const Vector3& pos1,
0070                                 const Vector3& pos2) -> double {
0071         Vector2 O(0, 0);
0072         Vector2 p1p2 = (pos2.block<2, 1>(0, 0) - pos1.block<2, 1>(0, 0));
0073         double L = p1p2.norm();
0074         Vector2 p1O = (O - pos1.block<2, 1>(0, 0));
0075
0076         // Don't try parallel lines
0077         if (L < 1e-7) {
0078           return std::numeric_limits<double>::max();
0079         }
0080         double f = p1p2.dot(p1O) / L;
0081
0082         // Clamp to [0, |p1p2|]
0083         f = std::min(L, std::max(0., f));
0084         Vector2 closest = f * p1p2.normalized() + pos1.block<2, 1>(0, 0);
0085         double dist = (closest - O).norm();
0086         return dist;
0087       };
0088
0089       for (std::size_t iv = 1; iv < vtxs.size() + 1; ++iv) {
0090         std::size_t fpoint = iv < vtxs.size() ? iv : 0;
0091         double testR = radialDistance(vtxs[fpoint], vtxs[iv - 1]);
0092         extent.range(binR).expandMin(testR);
0093       }
0094     }
0095   }
0096   return extent;
0097 }