Core/Surfaces/SurfaceArrayTests.cpp

0001 // This file is part of the Acts project.
0002 //
0003 // Copyright (C) 2017-2018 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 <boost/test/data/test_case.hpp>
0010 #include <boost/test/tools/output_test_stream.hpp>
0011 #include <boost/test/unit_test.hpp>
0012
0013 #include "Acts/Definitions/Algebra.hpp"
0014 #include "Acts/Geometry/GeometryContext.hpp"
0015 #include "Acts/Surfaces/PlanarBounds.hpp"
0016 #include "Acts/Surfaces/PlaneSurface.hpp"
0017 #include "Acts/Surfaces/RectangleBounds.hpp"
0018 #include "Acts/Surfaces/Surface.hpp"
0019 #include "Acts/Surfaces/SurfaceArray.hpp"
0020 #include "Acts/Surfaces/SurfaceBounds.hpp"
0021 #include "Acts/Utilities/BinningType.hpp"
0022 #include "Acts/Utilities/Grid.hpp"
0023 #include "Acts/Utilities/Helpers.hpp"
0024 #include "Acts/Utilities/detail/Axis.hpp"
0025 #include "Acts/Utilities/detail/AxisFwd.hpp"
0026 #include "Acts/Utilities/detail/grid_helper.hpp"
0027
0028 #include <cmath>
0029 #include <cstddef>
0030 #include <fstream>
0031 #include <iomanip>
0032 #include <iostream>
0033 #include <memory>
0034 #include <string>
0035 #include <tuple>
0036 #include <utility>
0037 #include <vector>
0038
0039 #include <boost/format.hpp>
0040
0041 using Acts::VectorHelpers::phi;
0042
0043 namespace Acts::Test {
0044
0045 // Create a test context
0046 GeometryContext tgContext = GeometryContext();
0047
0048 using SrfVec = std::vector<std::shared_ptr<const Surface>>;
0049 struct SurfaceArrayFixture {
0050   std::vector<std::shared_ptr<const Surface>> m_surfaces;
0051
0052   SurfaceArrayFixture() { BOOST_TEST_MESSAGE("setup fixture"); }
0053   ~SurfaceArrayFixture() { BOOST_TEST_MESSAGE("teardown fixture"); }
0054
0055   SrfVec fullPhiTestSurfacesEC(std::size_t n = 10, double shift = 0,
0056                                double zbase = 0, double r = 10) {
0057     SrfVec res;
0058
0059     double phiStep = 2 * M_PI / n;
0060     for (std::size_t i = 0; i < n; ++i) {
0061       double z = zbase + ((i % 2 == 0) ? 1 : -1) * 0.2;
0062
0063       Transform3 trans;
0064       trans.setIdentity();
0065       trans.rotate(Eigen::AngleAxisd(i * phiStep + shift, Vector3(0, 0, 1)));
0066       trans.translate(Vector3(r, 0, z));
0067
0068       auto bounds = std::make_shared<const RectangleBounds>(2, 1);
0069       std::shared_ptr<const Surface> srf =
0070           Surface::makeShared<PlaneSurface>(trans, bounds);
0071
0072       res.push_back(srf);
0073       m_surfaces.push_back(
0074           std::move(srf));  // keep shared, will get destroyed at the end
0075     }
0076
0077     return res;
0078   }
0079
0080   SrfVec fullPhiTestSurfacesBRL(int n = 10, double shift = 0, double zbase = 0,
0081                                 double incl = M_PI / 9., double w = 2,
0082                                 double h = 1.5) {
0083     SrfVec res;
0084
0085     double phiStep = 2 * M_PI / n;
0086     for (int i = 0; i < n; ++i) {
0087       double z = zbase;
0088
0089       Transform3 trans;
0090       trans.setIdentity();
0091       trans.rotate(Eigen::AngleAxisd(i * phiStep + shift, Vector3(0, 0, 1)));
0092       trans.translate(Vector3(10, 0, z));
0093       trans.rotate(Eigen::AngleAxisd(incl, Vector3(0, 0, 1)));
0094       trans.rotate(Eigen::AngleAxisd(M_PI / 2., Vector3(0, 1, 0)));
0095
0096       auto bounds = std::make_shared<const RectangleBounds>(w, h);
0097       std::shared_ptr<const Surface> srf =
0098           Surface::makeShared<PlaneSurface>(trans, bounds);
0099
0100       res.push_back(srf);
0101       m_surfaces.push_back(
0102           std::move(srf));  // keep shared, will get destroyed at the end
0103     }
0104
0105     return res;
0106   }
0107
0108   SrfVec straightLineSurfaces(
0109       std::size_t n = 10., double step = 3, const Vector3& origin = {0, 0, 1.5},
0110       const Transform3& pretrans = Transform3::Identity(),
0111       const Vector3& dir = {0, 0, 1}) {
0112     SrfVec res;
0113     for (std::size_t i = 0; i < n; ++i) {
0114       Transform3 trans;
0115       trans.setIdentity();
0116       trans.translate(origin + dir * step * i);
0117       // trans.rotate(AngleAxis3(M_PI/9., Vector3(0, 0, 1)));
0118       trans.rotate(AngleAxis3(M_PI / 2., Vector3(1, 0, 0)));
0119       trans = trans * pretrans;
0120
0121       auto bounds = std::make_shared<const RectangleBounds>(2, 1.5);
0122
0123       std::shared_ptr<const Surface> srf =
0124           Surface::makeShared<PlaneSurface>(trans, bounds);
0125
0126       res.push_back(srf);
0127       m_surfaces.push_back(
0128           std::move(srf));  // keep shared, will get destroyed at the end
0129     }
0130
0131     return res;
0132   }
0133
0134   SrfVec makeBarrel(int nPhi, int nZ, double w, double h) {
0135     double z0 = -(nZ - 1) * w;
0136     SrfVec res;
0137
0138     for (int i = 0; i < nZ; i++) {
0139       double z = i * w * 2 + z0;
0140       // std::cout << "z=" << z << std::endl;
0141       SrfVec ring = fullPhiTestSurfacesBRL(nPhi, 0, z, M_PI / 9., w, h);
0142       res.insert(res.end(), ring.begin(), ring.end());
0143     }
0144
0145     return res;
0146   }
0147
0148   void draw_surfaces(const SrfVec& surfaces, const std::string& fname) {
0149     std::ofstream os;
0150     os.open(fname);
0151
0152     os << std::fixed << std::setprecision(4);
0153
0154     std::size_t nVtx = 0;
0155     for (const auto& srfx : surfaces) {
0156       std::shared_ptr<const PlaneSurface> srf =
0157           std::dynamic_pointer_cast<const PlaneSurface>(srfx);
0158       const PlanarBounds* bounds =
0159           dynamic_cast<const PlanarBounds*>(&srf->bounds());
0160
0161       for (const auto& vtxloc : bounds->vertices()) {
0162         Vector3 vtx =
0163             srf->transform(tgContext) * Vector3(vtxloc.x(), vtxloc.y(), 0);
0164         os << "v " << vtx.x() << " " << vtx.y() << " " << vtx.z() << "\n";
0165       }
0166
0167       // connect them
0168       os << "f";
0169       for (std::size_t i = 1; i <= bounds->vertices().size(); ++i) {
0170         os << " " << nVtx + i;
0171       }
0172       os << "\n";
0173
0174       nVtx += bounds->vertices().size();
0175     }
0176
0177     os.close();
0178   }
0179 };
0180
0181 BOOST_AUTO_TEST_SUITE(Surfaces)
0182
0183 BOOST_FIXTURE_TEST_CASE(SurfaceArray_create, SurfaceArrayFixture) {
0184   GeometryContext tgContext = GeometryContext();
0185
0186   SrfVec brl = makeBarrel(30, 7, 2, 1);
0187   std::vector<const Surface*> brlRaw = unpack_shared_vector(brl);
0188   draw_surfaces(brl, "SurfaceArray_create_BRL_1.obj");
0189
0190   detail::Axis<detail::AxisType::Equidistant, detail::AxisBoundaryType::Closed>
0191       phiAxis(-M_PI, M_PI, 30u);
0192   detail::Axis<detail::AxisType::Equidistant, detail::AxisBoundaryType::Bound>
0193       zAxis(-14, 14, 7u);
0194
0195   double angleShift = 2 * M_PI / 30. / 2.;
0196   auto transform = [angleShift](const Vector3& pos) {
0197     return Vector2(phi(pos) + angleShift, pos.z());
0198   };
0199   double R = 10;
0200   auto itransform = [angleShift, R](const Vector2& loc) {
0201     return Vector3(R * std::cos(loc[0] - angleShift),
0202                    R * std::sin(loc[0] - angleShift), loc[1]);
0203   };
0204
0205   auto sl = std::make_unique<
0206       SurfaceArray::SurfaceGridLookup<decltype(phiAxis), decltype(zAxis)>>(
0207       transform, itransform,
0208       std::make_tuple(std::move(phiAxis), std::move(zAxis)));
0209   sl->fill(tgContext, brlRaw);
0210   SurfaceArray sa(std::move(sl), brl);
0211
0212   // let's see if we can access all surfaces
0213   sa.toStream(tgContext, std::cout);
0214
0215   for (const auto& srf : brl) {
0216     Vector3 ctr = srf->binningPosition(tgContext, binR);
0217     std::vector<const Surface*> binContent = sa.at(ctr);
0218
0219     BOOST_CHECK_EQUAL(binContent.size(), 1u);
0220     BOOST_CHECK_EQUAL(srf.get(), binContent.at(0));
0221   }
0222
0223   std::vector<const Surface*> neighbors =
0224       sa.neighbors(itransform(Vector2(0, 0)));
0225   BOOST_CHECK_EQUAL(neighbors.size(), 9u);
0226
0227   auto sl2 = std::make_unique<
0228       SurfaceArray::SurfaceGridLookup<decltype(phiAxis), decltype(zAxis)>>(
0229       transform, itransform,
0230       std::make_tuple(std::move(phiAxis), std::move(zAxis)));
0231   // do NOT fill, only completebinning
0232   sl2->completeBinning(tgContext, brlRaw);
0233   SurfaceArray sa2(std::move(sl2), brl);
0234   sa.toStream(tgContext, std::cout);
0235   for (const auto& srf : brl) {
0236     Vector3 ctr = srf->binningPosition(tgContext, binR);
0237     std::vector<const Surface*> binContent = sa2.at(ctr);
0238
0239     BOOST_CHECK_EQUAL(binContent.size(), 1u);
0240     BOOST_CHECK_EQUAL(srf.get(), binContent.at(0));
0241   }
0242 }
0243
0244 BOOST_AUTO_TEST_CASE(SurfaceArray_singleElement) {
0245   double w = 3, h = 4;
0246   auto bounds = std::make_shared<const RectangleBounds>(w, h);
0247   auto srf = Surface::makeShared<PlaneSurface>(Transform3::Identity(), bounds);
0248
0249   SurfaceArray sa(srf);
0250
0251   auto binContent = sa.at(Vector3(42, 42, 42));
0252   BOOST_CHECK_EQUAL(binContent.size(), 1u);
0253   BOOST_CHECK_EQUAL(binContent.at(0), srf.get());
0254   BOOST_CHECK_EQUAL(sa.surfaces().size(), 1u);
0255   BOOST_CHECK_EQUAL(sa.surfaces().at(0), srf.get());
0256 }
0257
0258 BOOST_AUTO_TEST_CASE(SurfaceArray_manyElementsSingleLookup) {
0259   double w = 3, h = 4;
0260   auto bounds = std::make_shared<const RectangleBounds>(w, h);
0261   auto srf0 = Surface::makeShared<PlaneSurface>(Transform3::Identity(), bounds);
0262   auto srf1 = Surface::makeShared<PlaneSurface>(Transform3::Identity(), bounds);
0263
0264   std::vector<const Surface*> sfPointers = {srf0.get(), srf1.get()};
0265   std::vector<std::shared_ptr<const Surface>> surfaces = {srf0, srf1};
0266
0267   auto singleLookUp =
0268       std::make_unique<Acts::SurfaceArray::SingleElementLookup>(sfPointers);
0269
0270   SurfaceArray sa(std::move(singleLookUp), surfaces);
0271
0272   auto binContent = sa.at(Vector3(42, 42, 42));
0273   BOOST_CHECK_EQUAL(binContent.size(), 2u);
0274   BOOST_CHECK_EQUAL(sa.surfaces().size(), 2u);
0275 }
0276
0277 BOOST_AUTO_TEST_SUITE_END()
0278 }  // namespace Acts::Test