src/Geometry/CutoutCylinderVolumeBounds.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/CutoutCylinderVolumeBounds.hpp"
0010
0011 #include "Acts/Definitions/Algebra.hpp"
0012 #include "Acts/Definitions/Direction.hpp"
0013 #include "Acts/Definitions/Tolerance.hpp"
0014 #include "Acts/Geometry/BoundarySurfaceFace.hpp"
0015 #include "Acts/Geometry/Volume.hpp"
0016 #include "Acts/Geometry/VolumeBounds.hpp"
0017 #include "Acts/Surfaces/CylinderBounds.hpp"
0018 #include "Acts/Surfaces/CylinderSurface.hpp"
0019 #include "Acts/Surfaces/DiscSurface.hpp"
0020 #include "Acts/Surfaces/RadialBounds.hpp"
0021 #include "Acts/Surfaces/Surface.hpp"
0022 #include "Acts/Utilities/BoundingBox.hpp"
0023
0024 #include <memory>
0025 #include <ostream>
0026 #include <stdexcept>
0027 #include <utility>
0028
0029 namespace Acts {
0030
0031 std::vector<double> CutoutCylinderVolumeBounds::values() const {
0032   return {m_values.begin(), m_values.end()};
0033 }
0034
0035 void CutoutCylinderVolumeBounds::checkConsistency() noexcept(false) {
0036   if (get(eMinR) < 0. || get(eMedR) <= 0. || get(eMaxR) <= 0. ||
0037       get(eMinR) >= get(eMedR) || get(eMinR) >= get(eMaxR) ||
0038       get(eMedR) >= get(eMaxR)) {
0039     throw std::invalid_argument(
0040         "CutoutCylinderVolumeBounds: invalid radial input.");
0041   }
0042   if (get(eHalfLengthZ) <= 0 || get(eHalfLengthZcutout) <= 0. ||
0043       get(eHalfLengthZcutout) > get(eHalfLengthZ)) {
0044     throw std::invalid_argument(
0045         "CutoutCylinderVolumeBounds: invalid longitudinal input.");
0046   }
0047 }
0048
0049 bool CutoutCylinderVolumeBounds::inside(const Vector3& gpos, double tol) const {
0050   // first check whether we are in the outer envelope at all (ignore r_med)
0051   using VectorHelpers::perp;
0052   using VectorHelpers::phi;
0053   double ros = perp(gpos);
0054
0055   bool insideR = (ros >= get(eMinR) - tol) && (ros <= get(eMaxR) + tol);
0056   bool insideZ = std::abs(gpos.z()) <= get(eHalfLengthZ) + tol;
0057
0058   if (!insideR || !insideZ) {
0059     return false;
0060   }
0061
0062   // we're inside the outer volume, but we might be in inside the
0063   // cutout section in the middle
0064   bool insideRInner = ros <= get(eMedR) - tol;
0065   bool insideZInner = std::abs(gpos.z()) < get(eHalfLengthZcutout) - tol;
0066
0067   return !insideRInner || !insideZInner;  // we are not, inside bounds
0068 }
0069
0070 std::vector<OrientedSurface> CutoutCylinderVolumeBounds::orientedSurfaces(
0071     const Transform3& transform) const {
0072   std::vector<OrientedSurface> oSurfaces;
0073
0074   if (get(eMinR) == 0.) {
0075     oSurfaces.resize(6);  // exactly six surfaces (no choke inner cover)
0076   } else {
0077     oSurfaces.resize(8);  // exactly eight surfaces
0078   }
0079
0080   // Outer cylinder envelope
0081   auto outer =
0082       Surface::makeShared<CylinderSurface>(transform, m_outerCylinderBounds);
0083   oSurfaces.at(tubeOuterCover) = {std::move(outer),
0084                                   Direction::OppositeNormal()};
0085
0086   // Inner (cutout) cylinder envelope
0087   auto cutoutInner =
0088       Surface::makeShared<CylinderSurface>(transform, m_cutoutCylinderBounds);
0089   oSurfaces.at(tubeInnerCover) = {std::move(cutoutInner),
0090                                   Direction::AlongNormal()};
0091
0092   // z position of the pos and neg choke points
0093   double hlChoke = (get(eHalfLengthZ) - get(eHalfLengthZcutout)) * 0.5;
0094   double zChoke = get(eHalfLengthZcutout) + hlChoke;
0095
0096   if (m_innerCylinderBounds != nullptr) {
0097     auto posChokeTrf = transform * Translation3(Vector3(0, 0, zChoke));
0098     auto posInner = Surface::makeShared<CylinderSurface>(posChokeTrf,
0099                                                          m_innerCylinderBounds);
0100     oSurfaces.at(index7) = {std::move(posInner), Direction::AlongNormal()};
0101
0102     auto negChokeTrf = transform * Translation3(Vector3(0, 0, -zChoke));
0103     auto negInner = Surface::makeShared<CylinderSurface>(negChokeTrf,
0104                                                          m_innerCylinderBounds);
0105     oSurfaces.at(index6) = {std::move(negInner), Direction::AlongNormal()};
0106   }
0107
0108   // Two Outer disks
0109   auto posOutDiscTrf =
0110       transform * Translation3(Vector3(0, 0, get(eHalfLengthZ)));
0111   auto posOutDisc =
0112       Surface::makeShared<DiscSurface>(posOutDiscTrf, m_outerDiscBounds);
0113   oSurfaces.at(positiveFaceXY) = {std::move(posOutDisc),
0114                                   Direction::OppositeNormal()};
0115
0116   auto negOutDiscTrf =
0117       transform * Translation3(Vector3(0, 0, -get(eHalfLengthZ)));
0118   auto negOutDisc =
0119       Surface::makeShared<DiscSurface>(negOutDiscTrf, m_outerDiscBounds);
0120   oSurfaces.at(negativeFaceXY) = {std::move(negOutDisc),
0121                                   Direction::AlongNormal()};
0122
0123   // Two Inner disks
0124   auto posInDiscTrf =
0125       transform * Translation3(Vector3(0, 0, get(eHalfLengthZcutout)));
0126   auto posInDisc =
0127       Surface::makeShared<DiscSurface>(posInDiscTrf, m_innerDiscBounds);
0128   oSurfaces.at(index5) = {std::move(posInDisc), Direction::AlongNormal()};
0129
0130   auto negInDiscTrf =
0131       transform * Translation3(Vector3(0, 0, -get(eHalfLengthZcutout)));
0132   auto negInDisc =
0133       Surface::makeShared<DiscSurface>(negInDiscTrf, m_innerDiscBounds);
0134   oSurfaces.at(index4) = {std::move(negInDisc), Direction::OppositeNormal()};
0135
0136   return oSurfaces;
0137 }
0138
0139 Volume::BoundingBox CutoutCylinderVolumeBounds::boundingBox(
0140     const Transform3* trf, const Vector3& envelope,
0141     const Volume* entity) const {
0142   Vector3 vmin, vmax;
0143
0144   // no phi sector is possible, so this is just the outer size of
0145   // the cylinder
0146
0147   vmax = {get(eMaxR), get(eMaxR), get(eHalfLengthZ)};
0148   vmin = {-get(eMaxR), -get(eMaxR), -get(eHalfLengthZ)};
0149
0150   Volume::BoundingBox box(entity, vmin - envelope, vmax + envelope);
0151   // transform at the very end, if required
0152   return trf == nullptr ? box : box.transformed(*trf);
0153 }
0154
0155 std::ostream& CutoutCylinderVolumeBounds::toStream(std::ostream& sl) const {
0156   sl << "CutoutCylinderVolumeBounds(\n";
0157   sl << "rmin = " << get(eMinR) << " rmed = " << get(eMedR)
0158      << " rmax = " << get(eMaxR) << "\n";
0159   sl << "dz1 = " << get(eHalfLengthZ) << " dz2 = " << get(eHalfLengthZcutout);
0160   return sl;
0161 }
0162
0163 void CutoutCylinderVolumeBounds::buildSurfaceBounds() {
0164   if (get(eMinR) > s_epsilon) {
0165     double hlChoke = (get(eHalfLengthZ) - get(eHalfLengthZcutout)) * 0.5;
0166     m_innerCylinderBounds =
0167         std::make_shared<CylinderBounds>(get(eMinR), hlChoke);
0168   }
0169
0170   m_cutoutCylinderBounds =
0171       std::make_shared<CylinderBounds>(get(eMedR), get(eHalfLengthZcutout));
0172
0173   m_outerCylinderBounds =
0174       std::make_shared<CylinderBounds>(get(eMaxR), get(eHalfLengthZ));
0175
0176   m_innerDiscBounds = std::make_shared<RadialBounds>(get(eMinR), get(eMedR));
0177
0178   m_outerDiscBounds = std::make_shared<RadialBounds>(get(eMinR), get(eMaxR));
0179 }
0180
0181 }  // namespace Acts