Plugins/Geant4/Geant4ConvertersTests.cpp

0001 // This file is part of the Acts project.
0002 //
0003 // Copyright (C) 2022 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/unit_test.hpp>
0010
0011 #include "Acts/Definitions/Algebra.hpp"
0012 #include "Acts/Material/ISurfaceMaterial.hpp"
0013 #include "Acts/Material/Material.hpp"
0014 #include "Acts/Material/MaterialSlab.hpp"
0015 #include "Acts/Plugins/Geant4/Geant4Converters.hpp"
0016 #include "Acts/Surfaces/CylinderBounds.hpp"
0017 #include "Acts/Surfaces/LineBounds.hpp"
0018 #include "Acts/Surfaces/RadialBounds.hpp"
0019 #include "Acts/Surfaces/RectangleBounds.hpp"
0020 #include "Acts/Surfaces/Surface.hpp"
0021 #include "Acts/Surfaces/TrapezoidBounds.hpp"
0022 #include "Acts/Tests/CommonHelpers/FloatComparisons.hpp"
0023
0024 #include <array>
0025 #include <cmath>
0026 #include <memory>
0027 #include <stdexcept>
0028 #include <tuple>
0029
0030 #include "G4Box.hh"
0031 #include "G4LogicalVolume.hh"
0032 #include "G4Material.hh"
0033 #include "G4PVPlacement.hh"
0034 #include "G4RotationMatrix.hh"
0035 #include "G4SystemOfUnits.hh"
0036 #include "G4ThreeVector.hh"
0037 #include "G4Trd.hh"
0038 #include "G4Tubs.hh"
0039 #include "G4VPhysicalVolume.hh"
0040
0041 Acts::ActsScalar rho = 1.2345;
0042 G4Material* g4Material = new G4Material("Material", 6., 12., rho);
0043
0044 BOOST_AUTO_TEST_SUITE(Geant4Plugin)
0045
0046 BOOST_AUTO_TEST_CASE(Geant4AlgebraConversion) {
0047   G4ThreeVector g4Translation(10., 20., 30.);
0048
0049   auto translated = Acts::Geant4AlgebraConverter{}.transform(g4Translation);
0050   auto actsTranslation = translated.translation();
0051   BOOST_CHECK_EQUAL(actsTranslation[0], 10.);
0052   BOOST_CHECK_EQUAL(actsTranslation[1], 20.);
0053   BOOST_CHECK_EQUAL(actsTranslation[2], 30.);
0054
0055   auto translatedScaled =
0056       Acts::Geant4AlgebraConverter{10.}.transform(g4Translation);
0057   auto actsTranslationScaled = translatedScaled.translation();
0058   BOOST_CHECK_EQUAL(actsTranslationScaled[0], 100.);
0059   BOOST_CHECK_EQUAL(actsTranslationScaled[1], 200.);
0060   BOOST_CHECK_EQUAL(actsTranslationScaled[2], 300.);
0061 }
0062
0063 BOOST_AUTO_TEST_CASE(Geant4CylinderConversion) {
0064   G4Tubs cylinder("Cylinder", 399., 401., 800., -M_PI * CLHEP::radian,
0065                   2 * M_PI * CLHEP::radian);
0066   auto [bounds, thickness] =
0067       Acts::Geant4ShapeConverter{}.cylinderBounds(cylinder);
0068   CHECK_CLOSE_ABS(bounds->get(Acts::CylinderBounds::BoundValues::eR), 400.,
0069                   10e-10);
0070   CHECK_CLOSE_ABS(bounds->get(Acts::CylinderBounds::BoundValues::eHalfLengthZ),
0071                   800., 10e-10);
0072   CHECK_CLOSE_ABS(
0073       bounds->get(Acts::CylinderBounds::BoundValues::eHalfPhiSector), M_PI,
0074       10e-10);
0075   CHECK_CLOSE_ABS(bounds->get(Acts::CylinderBounds::BoundValues::eAveragePhi),
0076                   0., 10e-10);
0077   CHECK_CLOSE_ABS(thickness, 2., 10e-10);
0078 }
0079
0080 BOOST_AUTO_TEST_CASE(Geant4RadialConversion) {
0081   G4Tubs disc("disc", 40., 400., 2., -M_PI * CLHEP::radian,
0082               2 * M_PI * CLHEP::radian);
0083   auto [bounds, thickness] = Acts::Geant4ShapeConverter{}.radialBounds(disc);
0084   CHECK_CLOSE_ABS(bounds->get(Acts::RadialBounds::BoundValues::eMinR), 40.,
0085                   10e-10);
0086   CHECK_CLOSE_ABS(bounds->get(Acts::RadialBounds::BoundValues::eMaxR), 400.,
0087                   10e-10);
0088   CHECK_CLOSE_ABS(bounds->get(Acts::RadialBounds::BoundValues::eHalfPhiSector),
0089                   M_PI, 10e-10);
0090   CHECK_CLOSE_ABS(bounds->get(Acts::RadialBounds::BoundValues::eAveragePhi), 0.,
0091                   10e-10);
0092   CHECK_CLOSE_ABS(thickness, 4., 10e-10);
0093 }
0094
0095 BOOST_AUTO_TEST_CASE(Geant4LineConversion) {
0096   G4Tubs line("line", 0., 20., 400., 0., 2 * M_PI);
0097   auto bounds = Acts::Geant4ShapeConverter{}.lineBounds(line);
0098   CHECK_CLOSE_ABS(bounds->get(Acts::LineBounds::BoundValues::eR), 20., 10e-10);
0099   CHECK_CLOSE_ABS(bounds->get(Acts::LineBounds::BoundValues::eHalfLengthZ),
0100                   400., 10e-10);
0101 }
0102
0103 BOOST_AUTO_TEST_CASE(Geant4BoxConversion) {
0104   // Test the standard orientations
0105   G4Box sensorXY("SensorXY", 23., 34., 1.);
0106   auto [boundsXY, axesXY, thicknessZ] =
0107       Acts::Geant4ShapeConverter{}.rectangleBounds(sensorXY);
0108   CHECK_CLOSE_ABS(boundsXY->halfLengthX(), 23., 10e-10);
0109   CHECK_CLOSE_ABS(boundsXY->halfLengthY(), 34., 10e-10);
0110   auto refXY = std::array<int, 2u>{0, 1};
0111   BOOST_CHECK(axesXY == refXY);
0112   CHECK_CLOSE_ABS(thicknessZ, 2., 10e-10);
0113
0114   G4Box sensorYZ("SensorYZ", 2., 45., 56.);
0115   auto [boundsYZ, axesYZ, thicknessX] =
0116       Acts::Geant4ShapeConverter{}.rectangleBounds(sensorYZ);
0117   CHECK_CLOSE_ABS(boundsYZ->halfLengthX(), 45., 10e-10);
0118   CHECK_CLOSE_ABS(boundsYZ->halfLengthY(), 56., 10e-10);
0119   auto refYZ = std::array<int, 2u>{1, 2};
0120   BOOST_CHECK(axesYZ == refYZ);
0121   CHECK_CLOSE_ABS(thicknessX, 4., 10e-10);
0122
0123   G4Box sensorZX("SensorZX", 78., 2., 67.);
0124   auto [boundsZX, axesZX, thicknessY] =
0125       Acts::Geant4ShapeConverter{}.rectangleBounds(sensorZX);
0126   CHECK_CLOSE_ABS(boundsZX->halfLengthX(), 67., 10e-10);
0127   CHECK_CLOSE_ABS(boundsZX->halfLengthY(), 78., 10e-10);
0128   auto refZX = std::array<int, 2u>{2, 0};
0129   BOOST_CHECK(axesZX == refZX);
0130   CHECK_CLOSE_ABS(thicknessY, 4., 10e-10);
0131
0132   // Test the flipped axis
0133   G4Box sensorXz("SensorXz", 78., 2., 67.);
0134   auto [boundsXz, axesXz, thicknessY2] =
0135       Acts::Geant4ShapeConverter{1, true}.rectangleBounds(sensorXz);
0136   CHECK_CLOSE_ABS(boundsXz->halfLengthX(), 78., 10e-10);
0137   CHECK_CLOSE_ABS(boundsXz->halfLengthY(), 67., 10e-10);
0138   auto refXz = std::array<int, 2u>{0, -2};
0139   BOOST_CHECK(axesXz == refXz);
0140   CHECK_CLOSE_ABS(thicknessY2, 4., 10e-10);
0141 }
0142
0143 BOOST_AUTO_TEST_CASE(Geant4TrapzoidConversion) {
0144   // Standard TRD: XY are already well defined
0145   G4Trd trdXY("trdXY", 100, 150, 200, 200, 2);
0146   auto [boundsXY, axesXY, thicknessZ] =
0147       Acts::Geant4ShapeConverter{}.trapezoidBounds(trdXY);
0148   CHECK_CLOSE_ABS(
0149       boundsXY->get(Acts::TrapezoidBounds::BoundValues::eHalfLengthXnegY), 100,
0150       10e-10);
0151   CHECK_CLOSE_ABS(
0152       boundsXY->get(Acts::TrapezoidBounds::BoundValues::eHalfLengthXposY), 150,
0153       10e-10);
0154   CHECK_CLOSE_ABS(
0155       boundsXY->get(Acts::TrapezoidBounds::BoundValues::eHalfLengthY), 200,
0156       10e-10);
0157   auto refXY = std::array<int, 2u>{0, 1};
0158   BOOST_CHECK(axesXY == refXY);
0159   CHECK_CLOSE_ABS(thicknessZ, 4., 10e-10);
0160
0161   // Flipped, yX are the coordinates
0162   G4Trd trdyX("trdyX", 200, 200, 100, 150, 2);
0163   auto [boundsyX, axesyX, thicknessZ2] =
0164       Acts::Geant4ShapeConverter{}.trapezoidBounds(trdyX);
0165   CHECK_CLOSE_ABS(
0166       boundsyX->get(Acts::TrapezoidBounds::BoundValues::eHalfLengthXnegY), 100,
0167       10e-10);
0168   CHECK_CLOSE_ABS(
0169       boundsyX->get(Acts::TrapezoidBounds::BoundValues::eHalfLengthXposY), 150,
0170       10e-10);
0171   CHECK_CLOSE_ABS(
0172       boundsyX->get(Acts::TrapezoidBounds::BoundValues::eHalfLengthY), 200,
0173       10e-10);
0174   auto refyX = std::array<int, 2u>{-1, 0};
0175   BOOST_CHECK(axesyX == refyX);
0176   CHECK_CLOSE_ABS(thicknessZ2, 4., 10e-10);
0177
0178   // YZ span the trapezoid
0179   G4Trd trdYZ("trdYZ", 2, 2, 120, 140, 200);
0180   auto [boundsYZ, axesYZ, thicknessX] =
0181       Acts::Geant4ShapeConverter{}.trapezoidBounds(trdYZ);
0182   CHECK_CLOSE_ABS(
0183       boundsYZ->get(Acts::TrapezoidBounds::BoundValues::eHalfLengthXnegY), 120.,
0184       10e-10);
0185   CHECK_CLOSE_ABS(
0186       boundsYZ->get(Acts::TrapezoidBounds::BoundValues::eHalfLengthXposY), 140.,
0187       10e-10);
0188   CHECK_CLOSE_ABS(
0189       boundsYZ->get(Acts::TrapezoidBounds::BoundValues::eHalfLengthY), 200.,
0190       10e-10);
0191   auto refYZ = std::array<int, 2u>{1, 2};
0192   BOOST_CHECK(axesYZ == refYZ);
0193   CHECK_CLOSE_ABS(thicknessX, 4., 10e-10);
0194
0195   // Xz span the trapezoid
0196   G4Trd trdXz("trdXz", 50, 75, 1, 1, 200);
0197   auto [boundsXz, axesXz, thicknessY] =
0198       Acts::Geant4ShapeConverter{}.trapezoidBounds(trdXz);
0199   CHECK_CLOSE_ABS(
0200       boundsXz->get(Acts::TrapezoidBounds::BoundValues::eHalfLengthXnegY), 50.,
0201       10e-10);
0202   CHECK_CLOSE_ABS(
0203       boundsXz->get(Acts::TrapezoidBounds::BoundValues::eHalfLengthXposY), 75.,
0204       10e-10);
0205   CHECK_CLOSE_ABS(
0206       boundsXz->get(Acts::TrapezoidBounds::BoundValues::eHalfLengthY), 200.,
0207       10e-10);
0208   auto refXz = std::array<int, 2u>{0, -2};
0209   BOOST_CHECK(axesXz == refXz);
0210   CHECK_CLOSE_ABS(thicknessY, 2., 10e-10);
0211 }
0212
0213 BOOST_AUTO_TEST_CASE(Geant4PlanarConversion) {
0214   G4Box boxXY("boxXY", 23., 34., 1.);
0215   auto pBoundsBox =
0216       std::get<0u>(Acts::Geant4ShapeConverter{}.planarBounds(boxXY));
0217   auto rBounds = dynamic_cast<const Acts::RectangleBounds*>(pBoundsBox.get());
0218   BOOST_CHECK_NE(rBounds, nullptr);
0219
0220   G4Trd trdXY("trdXY", 100, 150, 200, 200, 2);
0221   auto pBoundsTrd =
0222       std::get<0u>(Acts::Geant4ShapeConverter{}.planarBounds(trdXY));
0223   auto tBounds = dynamic_cast<const Acts::TrapezoidBounds*>(pBoundsTrd.get());
0224   BOOST_CHECK_NE(tBounds, nullptr);
0225 }
0226
0227 BOOST_AUTO_TEST_CASE(Geant4BoxVPhysConversion) {
0228   Acts::ActsScalar thickness = 2.;
0229
0230   G4Box* g4Box = new G4Box("Box", 23., 34., 0.5 * thickness);
0231   G4RotationMatrix* g4Rot = new G4RotationMatrix({0., 0., 1.}, 1.2);
0232   G4LogicalVolume* g4BoxLog = new G4LogicalVolume(g4Box, g4Material, "BoxLog");
0233
0234   G4ThreeVector g4Trans(0., 0., 100.);
0235   G4PVPlacement g4BoxPhys(g4Rot, g4Trans, g4BoxLog, "BoxPhys", nullptr, false,
0236                           1);
0237
0238   auto planeSurface = Acts::Geant4PhysicalVolumeConverter{}.surface(
0239       g4BoxPhys, Acts::Transform3::Identity(), true, thickness);
0240   BOOST_REQUIRE_NE(planeSurface, nullptr);
0241   BOOST_CHECK_EQUAL(planeSurface->type(), Acts::Surface::SurfaceType::Plane);
0242
0243   auto material = planeSurface->surfaceMaterial();
0244   BOOST_REQUIRE_NE(material, nullptr);
0245
0246   auto materialSlab = material->materialSlab(Acts::Vector3{0., 0., 0.});
0247   // Here it should be uncompressed material
0248   CHECK_CLOSE_ABS(materialSlab.material().massDensity(), rho, 0.001);
0249   CHECK_CLOSE_REL(thickness / g4Material->GetRadlen(),
0250                   materialSlab.thicknessInX0(), 0.1);
0251
0252   // Convert with compression
0253   Acts::ActsScalar compression = 4.;
0254   planeSurface = Acts::Geant4PhysicalVolumeConverter{}.surface(
0255       g4BoxPhys, Acts::Transform3::Identity(), true, thickness / compression);
0256   BOOST_REQUIRE_NE(planeSurface, nullptr);
0257   BOOST_CHECK_EQUAL(planeSurface->type(), Acts::Surface::SurfaceType::Plane);
0258
0259   material = planeSurface->surfaceMaterial();
0260   BOOST_REQUIRE_NE(material, nullptr);
0261   materialSlab = material->materialSlab(Acts::Vector3{0., 0., 0.});
0262
0263   // Here it should be uncompressed material
0264   CHECK_CLOSE_ABS(materialSlab.material().massDensity(), compression * rho,
0265                   0.001);
0266   CHECK_CLOSE_REL(thickness / g4Material->GetRadlen(),
0267                   materialSlab.thicknessInX0(), 0.01);
0268
0269   CHECK_CLOSE_ABS(materialSlab.thickness(), thickness / compression, 0.01);
0270   CHECK_CLOSE_REL(materialSlab.material().X0() * compression,
0271                   g4Material->GetRadlen(), 0.01);
0272
0273   delete g4Box;
0274   delete g4Rot;
0275   delete g4BoxLog;
0276 }
0277
0278 BOOST_AUTO_TEST_CASE(Geant4CylVPhysConversion) {
0279   Acts::ActsScalar radius = 45.;
0280   Acts::ActsScalar thickness = 1.;
0281   Acts::ActsScalar halfLengthZ = 200;
0282
0283   G4Tubs* g4Tube = new G4Tubs("Tube", radius, radius + thickness, halfLengthZ,
0284                               -M_PI * CLHEP::radian, 2 * M_PI * CLHEP::radian);
0285
0286   G4RotationMatrix* g4Rot = new G4RotationMatrix({0., 0., 1.}, 0.);
0287   G4LogicalVolume* g4TubeLog =
0288       new G4LogicalVolume(g4Tube, g4Material, "TubeLog");
0289   G4ThreeVector g4Trans(0., 0., 100.);
0290   G4PVPlacement g4CylinderPhys(g4Rot, g4Trans, g4TubeLog, "TubePhys", nullptr,
0291                                false, 1);
0292
0293   auto cylinderSurface = Acts::Geant4PhysicalVolumeConverter{}.surface(
0294       g4CylinderPhys, Acts::Transform3::Identity(), true, thickness);
0295   BOOST_REQUIRE_NE(cylinderSurface, nullptr);
0296   BOOST_CHECK_EQUAL(cylinderSurface->type(),
0297                     Acts::Surface::SurfaceType::Cylinder);
0298
0299   auto material = cylinderSurface->surfaceMaterial();
0300   BOOST_REQUIRE_NE(material, nullptr);
0301
0302   auto materialSlab = material->materialSlab(Acts::Vector3{0., 0., 0.});
0303   CHECK_CLOSE_REL(thickness / g4Material->GetRadlen(),
0304                   materialSlab.thicknessInX0(), 0.1);
0305
0306   // Here it should be uncompressed material
0307   CHECK_CLOSE_ABS(materialSlab.material().massDensity(), rho, 0.001);
0308
0309   /// CHECK exception throwing
0310   BOOST_CHECK_THROW(
0311       Acts::Geant4PhysicalVolumeConverter{Acts::Surface::SurfaceType::Plane}
0312           .surface(g4CylinderPhys, Acts::Transform3::Identity(), true,
0313                    thickness),
0314       std::runtime_error);
0315
0316   delete g4Tube;
0317   delete g4Rot;
0318   delete g4TubeLog;
0319 }
0320
0321 BOOST_AUTO_TEST_CASE(Geant4VDiscVPhysConversion) {
0322   Acts::ActsScalar innerRadius = 45.;
0323   Acts::ActsScalar outerRadius = 75.;
0324   Acts::ActsScalar thickness = 2.;
0325
0326   G4Tubs* g4Tube = new G4Tubs("Disc", innerRadius, outerRadius, 0.5 * thickness,
0327                               -M_PI * CLHEP::radian, 2 * M_PI * CLHEP::radian);
0328
0329   G4RotationMatrix* g4Rot = new G4RotationMatrix({0., 0., 1.}, 0.);
0330   G4LogicalVolume* g4TubeLog =
0331       new G4LogicalVolume(g4Tube, g4Material, "TubeLog");
0332   G4ThreeVector g4Trans(0., 0., 100.);
0333   G4PVPlacement g4discPhys(g4Rot, g4Trans, g4TubeLog, "TubePhys", nullptr,
0334                            false, 1);
0335
0336   auto discSurface = Acts::Geant4PhysicalVolumeConverter{}.surface(
0337       g4discPhys, Acts::Transform3::Identity(), true, thickness);
0338   BOOST_REQUIRE_NE(discSurface, nullptr);
0339   BOOST_CHECK_EQUAL(discSurface->type(), Acts::Surface::SurfaceType::Disc);
0340
0341   auto material = discSurface->surfaceMaterial();
0342   BOOST_REQUIRE_NE(material, nullptr);
0343
0344   auto materialSlab = material->materialSlab(Acts::Vector3{0., 0., 0.});
0345   // Here it should be uncompressed material
0346   CHECK_CLOSE_ABS(materialSlab.material().massDensity(), rho, 0.001);
0347
0348   delete g4Tube;
0349   delete g4Rot;
0350   delete g4TubeLog;
0351 }
0352
0353 BOOST_AUTO_TEST_CASE(Geant4LineVPhysConversion) {
0354   Acts::ActsScalar innerRadius = 0.;
0355   Acts::ActsScalar outerRadius = 20.;
0356   Acts::ActsScalar thickness = 400.;
0357
0358   G4Tubs* g4Tube = new G4Tubs("Line", innerRadius, outerRadius, 0.5 * thickness,
0359                               -M_PI * CLHEP::radian, 2 * M_PI * CLHEP::radian);
0360
0361   G4RotationMatrix* g4Rot = new G4RotationMatrix({0., 0., 1.}, 0.);
0362   G4LogicalVolume* g4TubeLog =
0363       new G4LogicalVolume(g4Tube, g4Material, "LineLog");
0364   G4ThreeVector g4Trans(0., 0., 100.);
0365   G4PVPlacement g4linePhys(g4Rot, g4Trans, g4TubeLog, "LinePhys", nullptr,
0366                            false, 1);
0367
0368   auto lineSurface =
0369       Acts::Geant4PhysicalVolumeConverter{Acts::Surface::SurfaceType::Straw}
0370           .surface(g4linePhys, Acts::Transform3::Identity(), true, thickness);
0371   BOOST_REQUIRE_NE(lineSurface, nullptr);
0372   BOOST_CHECK_EQUAL(lineSurface->type(), Acts::Surface::SurfaceType::Straw);
0373
0374   delete g4Tube;
0375   delete g4Rot;
0376   delete g4TubeLog;
0377 }
0378
0379 BOOST_AUTO_TEST_SUITE_END()