Plugins/Root/TGeoTrd2ConversionTests.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 <boost/test/unit_test.hpp>
0010
0011 #include "Acts/Definitions/Algebra.hpp"
0012 #include "Acts/Definitions/Tolerance.hpp"
0013 #include "Acts/Geometry/GeometryContext.hpp"
0014 #include "Acts/Surfaces/Surface.hpp"
0015 #include "Acts/Surfaces/SurfaceBounds.hpp"
0016 #include "Acts/Surfaces/TrapezoidBounds.hpp"
0017 #include "Acts/Visualization/GeometryView3D.hpp"
0018 #include "Acts/Visualization/ObjVisualization3D.hpp"
0019 #include "Acts/Visualization/ViewConfig.hpp"
0020 #include "ActsPlugins/Root/TGeoSurfaceConverter.hpp"
0021 #include "ActsTests/CommonHelpers/FloatComparisons.hpp"
0022
0023 #include <algorithm>
0024 #include <cstddef>
0025 #include <memory>
0026 #include <stdexcept>
0027 #include <string>
0028 #include <utility>
0029 #include <vector>
0030
0031 #include "TGeoManager.h"
0032 #include "TGeoMaterial.h"
0033 #include "TGeoMatrix.h"
0034 #include "TGeoMedium.h"
0035 #include "TGeoTrd1.h"
0036 #include "TGeoVolume.h"
0037 #include "TView.h"
0038
0039 using namespace Acts;
0040 using namespace ActsPlugins;
0041
0042 namespace ActsTests {
0043
0044 GeometryContext tgContext = GeometryContext();
0045
0046 ViewConfig red{.color = {200, 0, 0}};
0047 ViewConfig green{.color = {0, 200, 0}};
0048 ViewConfig blue{.color = {0, 0, 200}};
0049
0050 BOOST_AUTO_TEST_SUITE(RootSuite)
0051
0052 /// @brief Unit test to convert a TGeoTrd2 into a Plane
0053 ///
0054 /// * The TGeoTrd2 has x/z orientation
0055 BOOST_AUTO_TEST_CASE(TGeoTrd2_xz_to_PlaneSurface) {
0056   ObjVisualization3D objVis;
0057
0058   double hxmin = 10.;
0059   double hxmax = 30.;
0060   double ht = 1.;  // this is the half thickness
0061   double hy = 40.;
0062
0063   new TGeoManager("trd1", "poza9");
0064   TGeoMaterial *mat = new TGeoMaterial("Al", 26.98, 13, 2.7);
0065   TGeoMedium *med = new TGeoMedium("MED", 1, mat);
0066   TGeoVolume *top = gGeoManager->MakeBox("TOP", med, 100, 100, 100);
0067   gGeoManager->SetTopVolume(top);
0068   TGeoVolume *vol =
0069       gGeoManager->MakeTrd2("Trd2", med, hxmin, hxmax, ht, ht, hy);
0070   gGeoManager->CloseGeometry();
0071
0072   // Check the 4 possible ways
0073   std::vector<std::string> axesTypes = {"XZ*", "xZ*", "xz*", "Xz*"};
0074
0075   std::size_t itrd = 0;
0076   for (const auto &axes : axesTypes) {
0077     auto [plane, thickness] = TGeoSurfaceConverter::toSurface(
0078         *vol->GetShape(), *gGeoIdentity, axes, 1);
0079     BOOST_REQUIRE_NE(plane, nullptr);
0080     BOOST_CHECK_EQUAL(plane->type(), Surface::Plane);
0081     CHECK_CLOSE_ABS(thickness, 2 * ht, s_epsilon);
0082
0083     auto bounds = dynamic_cast<const TrapezoidBounds *>(&(plane->bounds()));
0084     BOOST_REQUIRE_NE(bounds, nullptr);
0085     double hXminY = bounds->get(TrapezoidBounds::eHalfLengthXnegY);
0086     double hXmaxY = bounds->get(TrapezoidBounds::eHalfLengthXposY);
0087     double hY = bounds->get(TrapezoidBounds::eHalfLengthY);
0088
0089     CHECK_CLOSE_ABS(hxmin, std::min(hXminY, hXmaxY), s_epsilon);
0090     CHECK_CLOSE_ABS(hxmax, std::max(hXminY, hXmaxY), s_epsilon);
0091     CHECK_CLOSE_ABS(hy, hY, s_epsilon);
0092
0093     // Check if the surface is the (negative) identity
0094     auto transform = plane->transform(tgContext);
0095     auto rotation = transform.rotation();
0096     const Vector3 offset{(-5.5 + (itrd++) * 2.5) * hxmax, 0., 0.};
0097     GeometryView3D::drawSurface(objVis, *plane, tgContext,
0098                                 Translation3{offset} * Transform3::Identity());
0099     const Vector3 center = plane->center(tgContext) + offset;
0100     GeometryView3D::drawArrowForward(
0101         objVis, center, center + 1.2 * (hXminY + hXmaxY) * rotation.col(0), 4.,
0102         2.5, red);
0103     GeometryView3D::drawArrowForward(
0104         objVis, center, center + 1.2 * hY * rotation.col(1), 4., 2.5, green);
0105     GeometryView3D::drawArrowForward(
0106         objVis, center, center + 2 * rotation.col(2), 4., 2.5, blue);
0107   }
0108   objVis.write("TGeoConversion_TGeoTrd2_xz_PlaneSurface");
0109
0110   // Check exceptions for not allowed axis definition
0111   std::vector<std::string> notAllowed = {"XY*", "xy*", "Xy*", "xY*"};
0112   for (const auto &naxis : notAllowed) {
0113     BOOST_CHECK_THROW(TGeoSurfaceConverter::toSurface(*vol->GetShape(),
0114                                                       *gGeoIdentity, naxis, 1),
0115                       std::invalid_argument);
0116   }
0117 }
0118
0119 /// @brief Unit test to convert a TGeoTrd2 into a Plane
0120 ///
0121 /// * The TGeoTrd2 has y/z orientation
0122 BOOST_AUTO_TEST_CASE(TGeoTrd2_yz_to_PlaneSurface) {
0123   ObjVisualization3D objVis;
0124
0125   double hxmin = 10.;
0126   double hxmax = 30.;
0127   double ht = 1.;  // this is the half thickness
0128   double hy = 40.;
0129
0130   new TGeoManager("trd1", "poza9");
0131   TGeoMaterial *mat = new TGeoMaterial("Al", 26.98, 13, 2.7);
0132   TGeoMedium *med = new TGeoMedium("MED", 1, mat);
0133   TGeoVolume *top = gGeoManager->MakeBox("TOP", med, 100, 100, 100);
0134   gGeoManager->SetTopVolume(top);
0135   TGeoVolume *vol =
0136       gGeoManager->MakeTrd2("Trd2", med, ht, ht, hxmin, hxmax, hy);
0137   gGeoManager->CloseGeometry();
0138
0139   // Check the 4 possible ways
0140   std::vector<std::string> axesTypes = {"YZ*", "yZ*", "yz*", "Yz*"};
0141
0142   std::size_t itrd = 0;
0143   for (const auto &axes : axesTypes) {
0144     auto [plane, thickness] = TGeoSurfaceConverter::toSurface(
0145         *vol->GetShape(), *gGeoIdentity, axes, 1);
0146     BOOST_REQUIRE_NE(plane, nullptr);
0147     BOOST_CHECK_EQUAL(plane->type(), Surface::Plane);
0148     CHECK_CLOSE_ABS(thickness, 2 * ht, s_epsilon);
0149
0150     auto bounds = dynamic_cast<const TrapezoidBounds *>(&(plane->bounds()));
0151     BOOST_REQUIRE_NE(bounds, nullptr);
0152     double hXminY = bounds->get(TrapezoidBounds::eHalfLengthXnegY);
0153     double hXmaxY = bounds->get(TrapezoidBounds::eHalfLengthXposY);
0154     double hY = bounds->get(TrapezoidBounds::eHalfLengthY);
0155
0156     CHECK_CLOSE_ABS(hxmin, std::min(hXminY, hXmaxY), s_epsilon);
0157     CHECK_CLOSE_ABS(hxmax, std::max(hXminY, hXmaxY), s_epsilon);
0158     CHECK_CLOSE_ABS(hy, hY, s_epsilon);
0159
0160     // Check if the surface is the (negative) identity
0161     auto transform = plane->transform(tgContext);
0162     auto rotation = transform.rotation();
0163     const Vector3 offset{(-5.5 + (itrd++) * 2.5) * hxmax, 0., 0.};
0164     GeometryView3D::drawSurface(objVis, *plane, tgContext,
0165                                 Translation3{offset} * Transform3::Identity());
0166     const Vector3 center = plane->center(tgContext) + offset;
0167     GeometryView3D::drawArrowForward(
0168         objVis, center, center + 1.2 * (hXminY + hXmaxY) * rotation.col(0), 4.,
0169         2.5, red);
0170     GeometryView3D::drawArrowForward(
0171         objVis, center, center + 1.2 * hY * rotation.col(1), 4., 2.5, green);
0172     GeometryView3D::drawArrowForward(
0173         objVis, center, center + 2 * rotation.col(2), 4., 2.5, blue);
0174   }
0175   objVis.write("TGeoConversion_TGeoTrd2_yz_PlaneSurface");
0176
0177   // Check exceptions for not allowed axis definition
0178   std::vector<std::string> notAllowed = {"YX*", "yx*", "yX*", "Yx*"};
0179   for (const auto &naxis : notAllowed) {
0180     BOOST_CHECK_THROW(TGeoSurfaceConverter::toSurface(*vol->GetShape(),
0181                                                       *gGeoIdentity, naxis, 1),
0182                       std::invalid_argument);
0183   }
0184 }
0185
0186 BOOST_AUTO_TEST_SUITE_END()
0187
0188 }  // namespace ActsTests