Core/Material/MaterialGridHelperTests.cpp

0001 // This file is part of the Acts project.
0002 //
0003 // Copyright (C) 2019 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/AccumulatedVolumeMaterial.hpp"
0013 #include "Acts/Material/Material.hpp"
0014 #include "Acts/Material/MaterialGridHelper.hpp"
0015 #include "Acts/Material/MaterialSlab.hpp"
0016 #include "Acts/Tests/CommonHelpers/FloatComparisons.hpp"
0017 #include "Acts/Utilities/BinUtility.hpp"
0018 #include "Acts/Utilities/BinningType.hpp"
0019 #include "Acts/Utilities/Grid.hpp"
0020 #include "Acts/Utilities/detail/Axis.hpp"
0021 #include "Acts/Utilities/detail/AxisFwd.hpp"
0022
0023 #include <cmath>
0024 #include <functional>
0025 #include <memory>
0026 #include <utility>
0027 #include <vector>
0028
0029 namespace Acts::Test {
0030
0031 using EAxis = Acts::detail::EquidistantAxis;
0032 using Grid2D = Acts::Grid<Acts::AccumulatedVolumeMaterial, EAxis, EAxis>;
0033 using Grid3D = Acts::Grid<Acts::AccumulatedVolumeMaterial, EAxis, EAxis, EAxis>;
0034 using MaterialGrid2D =
0035     Acts::Grid<Acts::Material::ParametersVector, EAxis, EAxis>;
0036 using MaterialGrid3D =
0037     Acts::Grid<Acts::Material::ParametersVector, EAxis, EAxis, EAxis>;
0038
0039 /// @brief Various test for the Material in the case of a Cuboid volume and 2D
0040 /// Grid
0041 BOOST_AUTO_TEST_CASE(Square_Grid_test) {
0042   BinUtility bu(7, -3., 3., open, binX);
0043   bu += BinUtility(3, -2., 2., open, binY);
0044   auto bd = bu.binningData();
0045   std::function<Acts::Vector2(Acts::Vector3)> transfoGlobalToLocal;
0046
0047   Grid2D Grid = createGrid2D(bu, transfoGlobalToLocal);
0048
0049   // Test Global To Local transform
0050   Acts::Vector3 pos(1., 2., 3.);
0051   Acts::Vector2 pos_2d(1., 2.);
0052   BOOST_CHECK_EQUAL(pos_2d, transfoGlobalToLocal(pos));
0053
0054   // Test Grid
0055   BOOST_CHECK_EQUAL(Grid.dimensions(), 2);
0056
0057   BOOST_CHECK_EQUAL(Grid.numLocalBins()[0], bd[0].bins());
0058   BOOST_CHECK_EQUAL(Grid.numLocalBins()[1], bd[1].bins());
0059
0060   BOOST_CHECK_EQUAL(Grid.minPosition()[0], bd[0].min);
0061   BOOST_CHECK_EQUAL(Grid.minPosition()[1], bd[1].min);
0062
0063   float max1 =
0064       bd[0].max + std::fabs(bd[0].max - bd[0].min) / (bd[0].bins() - 1);
0065   float max2 =
0066       bd[1].max + std::fabs(bd[1].max - bd[1].min) / (bd[1].bins() - 1);
0067
0068   BOOST_CHECK_EQUAL(Grid.maxPosition()[0], max1);
0069   BOOST_CHECK_EQUAL(Grid.maxPosition()[1], max2);
0070
0071   // Test pos to index
0072   Grid2D::index_t index1 = {1, 1};
0073   Grid2D::index_t index2 = {7, 2};
0074   Grid2D::index_t index3 = {1, 3};
0075
0076   Acts::Vector3 pos1 = {-2.6, -1.5, -0.7};
0077   Acts::Vector3 pos2 = {2.8, 0, 0.2};
0078   Acts::Vector3 pos3 = {-2.7, 1.8, 0.8};
0079
0080   for (int i = 0; i < 2; i++) {
0081     BOOST_CHECK_EQUAL(
0082         Grid.localBinsFromLowerLeftEdge(transfoGlobalToLocal(pos1))[i],
0083         index1[i]);
0084     BOOST_CHECK_EQUAL(
0085         Grid.localBinsFromLowerLeftEdge(transfoGlobalToLocal(pos2))[i],
0086         index2[i]);
0087     BOOST_CHECK_EQUAL(
0088         Grid.localBinsFromLowerLeftEdge(transfoGlobalToLocal(pos3))[i],
0089         index3[i]);
0090   }
0091   // Test material mapping
0092
0093   std::vector<Acts::Vector3> vectPos1;
0094   vectPos1.push_back(pos1);
0095   std::vector<Acts::Vector3> vectPos2;
0096   vectPos2.push_back(pos2);
0097   std::vector<Acts::Vector3> vectPos3;
0098   vectPos3.push_back(pos3);
0099
0100   std::vector<std::pair<MaterialSlab, std::vector<Vector3>>> matRecord;
0101   Material mat1 = Material::fromMolarDensity(1., 2., 3., 4., 5.);
0102   Material mat2 = Material::fromMolarDensity(6., 7., 8., 9., 10.);
0103   Material vacuum;
0104
0105   MaterialSlab matprop1(mat1, 1);
0106   MaterialSlab matprop2(mat2, 1);
0107
0108   matRecord.clear();
0109   matRecord.push_back(std::make_pair(matprop1, vectPos1));
0110   matRecord.push_back(std::make_pair(matprop2, vectPos2));
0111
0112   // Walk over each property
0113   for (const auto& rm : matRecord) {
0114     // Walk over each point associated with the properties
0115     for (const auto& point : rm.second) {
0116       // Search for fitting grid point and accumulate
0117       Acts::Grid2D::index_t index =
0118           Grid.localBinsFromLowerLeftEdge(transfoGlobalToLocal(point));
0119       Grid.atLocalBins(index).accumulate(rm.first);
0120     }
0121   }
0122
0123   MaterialGrid2D matMap = mapMaterialPoints(Grid);
0124
0125   CHECK_CLOSE_REL(matMap.atLocalBins(index1), mat1.parameters(), 1e-4);
0126   CHECK_CLOSE_REL(matMap.atLocalBins(index2), mat2.parameters(), 1e-4);
0127   BOOST_CHECK_EQUAL(matMap.atLocalBins(index3), vacuum.parameters());
0128 }
0129
0130 /// @brief Various test for the Material in the case of a Cylindrical volume
0131 /// with a 2D grid
0132 BOOST_AUTO_TEST_CASE(PhiZ_Grid_test) {
0133   BinUtility bu(2, -2., 2., open, binZ);
0134   bu += BinUtility(3, -M_PI, M_PI, closed, binPhi);
0135   auto bd = bu.binningData();
0136   std::function<Acts::Vector2(Acts::Vector3)> transfoGlobalToLocal;
0137
0138   Grid2D Grid = createGrid2D(bu, transfoGlobalToLocal);
0139
0140   // Test Global To Local transform
0141   Acts::Vector3 pos(1., 2., 3.);
0142
0143   CHECK_CLOSE_REL(transfoGlobalToLocal(pos)[1], atan2(2, 1), 1e-4);
0144   CHECK_CLOSE_REL(transfoGlobalToLocal(pos)[0], 3, 1e-4);
0145
0146   // Test Grid
0147   BOOST_CHECK_EQUAL(Grid.dimensions(), 2);
0148
0149   BOOST_CHECK_EQUAL(Grid.numLocalBins()[0], bd[0].bins());
0150   BOOST_CHECK_EQUAL(Grid.numLocalBins()[1], bd[1].bins());
0151
0152   BOOST_CHECK_EQUAL(Grid.minPosition()[0], bd[0].min);
0153   BOOST_CHECK_EQUAL(Grid.minPosition()[1], bd[1].min);
0154
0155   float max1 =
0156       bd[0].max + std::fabs(bd[0].max - bd[0].min) / (bd[0].bins() - 1);
0157   float max2 =
0158       bd[1].max + std::fabs(bd[1].max - bd[1].min) / (bd[1].bins() - 1);
0159
0160   BOOST_CHECK_EQUAL(Grid.maxPosition()[0], max1);
0161   BOOST_CHECK_EQUAL(Grid.maxPosition()[1], max2);
0162
0163   // Test pos to index
0164   Grid2D::index_t index1 = {1, 1};
0165   Grid2D::index_t index2 = {1, 2};
0166   Grid2D::index_t index3 = {2, 3};
0167
0168   Acts::Vector3 pos1 = {-0.2, -1, -1};
0169   Acts::Vector3 pos2 = {3.6, 0., -1.5};
0170   Acts::Vector3 pos3 = {-1, 0.3, 0.8};
0171
0172   for (int i = 0; i < 2; i++) {
0173     BOOST_CHECK_EQUAL(
0174         Grid.localBinsFromLowerLeftEdge(transfoGlobalToLocal(pos1))[i],
0175         index1[i]);
0176     BOOST_CHECK_EQUAL(
0177         Grid.localBinsFromLowerLeftEdge(transfoGlobalToLocal(pos2))[i],
0178         index2[i]);
0179     BOOST_CHECK_EQUAL(
0180         Grid.localBinsFromLowerLeftEdge(transfoGlobalToLocal(pos3))[i],
0181         index3[i]);
0182   }
0183
0184   // Test material mapping
0185   std::vector<Acts::Vector3> vectPos1;
0186   vectPos1.push_back(pos1);
0187   std::vector<Acts::Vector3> vectPos2;
0188   vectPos2.push_back(pos2);
0189   std::vector<Acts::Vector3> vectPos3;
0190   vectPos3.push_back(pos3);
0191
0192   std::vector<std::pair<MaterialSlab, std::vector<Vector3>>> matRecord;
0193   Material mat1 = Material::fromMolarDensity(1., 2., 3., 4., 5.);
0194   Material mat2 = Material::fromMolarDensity(6., 7., 8., 9., 10.);
0195   Material vacuum;
0196
0197   MaterialSlab matprop1(mat1, 1);
0198   MaterialSlab matprop2(mat2, 1);
0199
0200   matRecord.clear();
0201   matRecord.push_back(std::make_pair(matprop1, vectPos1));
0202   matRecord.push_back(std::make_pair(matprop2, vectPos2));
0203
0204   // Walk over each property
0205   for (const auto& rm : matRecord) {
0206     // Walk over each point associated with the properties
0207     for (const auto& point : rm.second) {
0208       // Search for fitting grid point and accumulate
0209       Acts::Grid2D::index_t index =
0210           Grid.localBinsFromLowerLeftEdge(transfoGlobalToLocal(point));
0211       Grid.atLocalBins(index).accumulate(rm.first);
0212     }
0213   }
0214
0215   MaterialGrid2D matMap = mapMaterialPoints(Grid);
0216
0217   CHECK_CLOSE_REL(matMap.atLocalBins(index1), mat1.parameters(), 1e-4);
0218   CHECK_CLOSE_REL(matMap.atLocalBins(index2), mat2.parameters(), 1e-4);
0219   BOOST_CHECK_EQUAL(matMap.atLocalBins(index3), vacuum.parameters());
0220 }
0221
0222 /// @brief Various test for the Material in the case of a Cuboid volume
0223 BOOST_AUTO_TEST_CASE(Cubic_Grid_test) {
0224   BinUtility bu(7, -3., 3., open, binX);
0225   bu += BinUtility(3, -2., 2., open, binY);
0226   bu += BinUtility(2, -1., 1., open, binZ);
0227   auto bd = bu.binningData();
0228   std::function<Acts::Vector3(Acts::Vector3)> transfoGlobalToLocal;
0229
0230   Grid3D Grid = createGrid3D(bu, transfoGlobalToLocal);
0231
0232   // Test Global To Local transform
0233   Acts::Vector3 pos(1., 2., 3.);
0234   BOOST_CHECK_EQUAL(pos, transfoGlobalToLocal(pos));
0235
0236   // Test Grid
0237   BOOST_CHECK_EQUAL(Grid.dimensions(), 3);
0238
0239   BOOST_CHECK_EQUAL(Grid.numLocalBins()[0], bd[0].bins());
0240   BOOST_CHECK_EQUAL(Grid.numLocalBins()[1], bd[1].bins());
0241   BOOST_CHECK_EQUAL(Grid.numLocalBins()[2], bd[2].bins());
0242
0243   BOOST_CHECK_EQUAL(Grid.minPosition()[0], bd[0].min);
0244   BOOST_CHECK_EQUAL(Grid.minPosition()[1], bd[1].min);
0245   BOOST_CHECK_EQUAL(Grid.minPosition()[2], bd[2].min);
0246
0247   float max1 =
0248       bd[0].max + std::fabs(bd[0].max - bd[0].min) / (bd[0].bins() - 1);
0249   float max2 =
0250       bd[1].max + std::fabs(bd[1].max - bd[1].min) / (bd[1].bins() - 1);
0251   float max3 =
0252       bd[2].max + std::fabs(bd[2].max - bd[2].min) / (bd[2].bins() - 1);
0253
0254   BOOST_CHECK_EQUAL(Grid.maxPosition()[0], max1);
0255   BOOST_CHECK_EQUAL(Grid.maxPosition()[1], max2);
0256   BOOST_CHECK_EQUAL(Grid.maxPosition()[2], max3);
0257
0258   // Test pos to index
0259   Grid3D::index_t index1 = {1, 1, 1};
0260   Grid3D::index_t index2 = {7, 2, 2};
0261   Grid3D::index_t index3 = {1, 3, 2};
0262
0263   Acts::Vector3 pos1 = {-2.6, -1.5, -0.7};
0264   Acts::Vector3 pos2 = {2.8, 0, 0.2};
0265   Acts::Vector3 pos3 = {-2.7, 1.8, 0.8};
0266
0267   for (int i = 0; i < 3; i++) {
0268     BOOST_CHECK_EQUAL(
0269         Grid.localBinsFromLowerLeftEdge(transfoGlobalToLocal(pos1))[i],
0270         index1[i]);
0271     BOOST_CHECK_EQUAL(
0272         Grid.localBinsFromLowerLeftEdge(transfoGlobalToLocal(pos2))[i],
0273         index2[i]);
0274     BOOST_CHECK_EQUAL(
0275         Grid.localBinsFromLowerLeftEdge(transfoGlobalToLocal(pos3))[i],
0276         index3[i]);
0277   }
0278   // Test material mapping
0279   std::vector<Acts::Vector3> vectPos1;
0280   vectPos1.push_back(pos1);
0281   std::vector<Acts::Vector3> vectPos2;
0282   vectPos2.push_back(pos2);
0283   std::vector<Acts::Vector3> vectPos3;
0284   vectPos3.push_back(pos3);
0285
0286   std::vector<std::pair<MaterialSlab, std::vector<Vector3>>> matRecord;
0287   Material mat1 = Material::fromMolarDensity(1., 2., 3., 4., 5.);
0288   Material mat2 = Material::fromMolarDensity(6., 7., 8., 9., 10.);
0289   Material vacuum;
0290
0291   MaterialSlab matprop1(mat1, 1);
0292   MaterialSlab matprop2(mat2, 1);
0293
0294   matRecord.clear();
0295   matRecord.push_back(std::make_pair(matprop1, vectPos1));
0296   matRecord.push_back(std::make_pair(matprop2, vectPos2));
0297
0298   // Walk over each property
0299   for (const auto& rm : matRecord) {
0300     // Walk over each point associated with the properties
0301     for (const auto& point : rm.second) {
0302       // Search for fitting grid point and accumulate
0303       Acts::Grid3D::index_t index =
0304           Grid.localBinsFromLowerLeftEdge(transfoGlobalToLocal(point));
0305       Grid.atLocalBins(index).accumulate(rm.first);
0306     }
0307   }
0308
0309   MaterialGrid3D matMap = mapMaterialPoints(Grid);
0310
0311   CHECK_CLOSE_REL(matMap.atLocalBins(index1), mat1.parameters(), 1e-4);
0312   CHECK_CLOSE_REL(matMap.atLocalBins(index2), mat2.parameters(), 1e-4);
0313   BOOST_CHECK_EQUAL(matMap.atLocalBins(index3), vacuum.parameters());
0314 }
0315
0316 /// @brief Various test for the Material in the case of a Cylindrical volume
0317 BOOST_AUTO_TEST_CASE(Cylindrical_Grid_test) {
0318   BinUtility bu(4, 1., 4., open, binR);
0319   bu += BinUtility(3, -M_PI, M_PI, closed, binPhi);
0320   bu += BinUtility(2, -2., 2., open, binZ);
0321   auto bd = bu.binningData();
0322   std::function<Acts::Vector3(Acts::Vector3)> transfoGlobalToLocal;
0323
0324   Grid3D Grid = createGrid3D(bu, transfoGlobalToLocal);
0325
0326   // Test Global To Local transform
0327   Acts::Vector3 pos(1., 2., 3.);
0328
0329   CHECK_CLOSE_REL(transfoGlobalToLocal(pos)[0], sqrt(5), 1e-4);
0330   CHECK_CLOSE_REL(transfoGlobalToLocal(pos)[1], atan2(2, 1), 1e-4);
0331   CHECK_CLOSE_REL(transfoGlobalToLocal(pos)[2], 3, 1e-4);
0332
0333   // Test Grid
0334   BOOST_CHECK_EQUAL(Grid.dimensions(), 3);
0335
0336   BOOST_CHECK_EQUAL(Grid.numLocalBins()[0], bd[0].bins());
0337   BOOST_CHECK_EQUAL(Grid.numLocalBins()[1], bd[1].bins());
0338   BOOST_CHECK_EQUAL(Grid.numLocalBins()[2], bd[2].bins());
0339
0340   BOOST_CHECK_EQUAL(Grid.minPosition()[0], bd[0].min);
0341   BOOST_CHECK_EQUAL(Grid.minPosition()[1], bd[1].min);
0342   BOOST_CHECK_EQUAL(Grid.minPosition()[2], bd[2].min);
0343
0344   float max1 =
0345       bd[0].max + std::fabs(bd[0].max - bd[0].min) / (bd[0].bins() - 1);
0346   float max2 =
0347       bd[1].max + std::fabs(bd[1].max - bd[1].min) / (bd[1].bins() - 1);
0348   float max3 =
0349       bd[2].max + std::fabs(bd[2].max - bd[2].min) / (bd[2].bins() - 1);
0350
0351   BOOST_CHECK_EQUAL(Grid.maxPosition()[0], max1);
0352   BOOST_CHECK_EQUAL(Grid.maxPosition()[1], max2);
0353   BOOST_CHECK_EQUAL(Grid.maxPosition()[2], max3);
0354
0355   // Test pos to index
0356   Grid3D::index_t index1 = {1, 1, 1};
0357   Grid3D::index_t index2 = {4, 2, 1};
0358   Grid3D::index_t index3 = {1, 3, 2};
0359
0360   Acts::Vector3 pos1 = {-0.2, -1, -1};
0361   Acts::Vector3 pos2 = {3.6, 0., -1.5};
0362   Acts::Vector3 pos3 = {-1, 0.3, 0.8};
0363
0364   for (int i = 0; i < 3; i++) {
0365     BOOST_CHECK_EQUAL(
0366         Grid.localBinsFromLowerLeftEdge(transfoGlobalToLocal(pos1))[i],
0367         index1[i]);
0368     BOOST_CHECK_EQUAL(
0369         Grid.localBinsFromLowerLeftEdge(transfoGlobalToLocal(pos2))[i],
0370         index2[i]);
0371     BOOST_CHECK_EQUAL(
0372         Grid.localBinsFromLowerLeftEdge(transfoGlobalToLocal(pos3))[i],
0373         index3[i]);
0374   }
0375
0376   // Test material mapping
0377   std::vector<Acts::Vector3> vectPos1;
0378   vectPos1.push_back(pos1);
0379   std::vector<Acts::Vector3> vectPos2;
0380   vectPos2.push_back(pos2);
0381   std::vector<Acts::Vector3> vectPos3;
0382   vectPos3.push_back(pos3);
0383
0384   std::vector<std::pair<MaterialSlab, std::vector<Vector3>>> matRecord;
0385   Material mat1 = Material::fromMolarDensity(1., 2., 3., 4., 5.);
0386   Material mat2 = Material::fromMolarDensity(6., 7., 8., 9., 10.);
0387   Material vacuum;
0388
0389   MaterialSlab matprop1(mat1, 1);
0390   MaterialSlab matprop2(mat2, 1);
0391
0392   matRecord.clear();
0393   matRecord.push_back(std::make_pair(matprop1, vectPos1));
0394   matRecord.push_back(std::make_pair(matprop2, vectPos2));
0395
0396   // Walk over each property
0397   for (const auto& rm : matRecord) {
0398     // Walk over each point associated with the properties
0399     for (const auto& point : rm.second) {
0400       // Search for fitting grid point and accumulate
0401       Acts::Grid3D::index_t index =
0402           Grid.localBinsFromLowerLeftEdge(transfoGlobalToLocal(point));
0403       Grid.atLocalBins(index).accumulate(rm.first);
0404     }
0405   }
0406
0407   MaterialGrid3D matMap = mapMaterialPoints(Grid);
0408
0409   CHECK_CLOSE_REL(matMap.atLocalBins(index1), mat1.parameters(), 1e-4);
0410   CHECK_CLOSE_REL(matMap.atLocalBins(index2), mat2.parameters(), 1e-4);
0411   BOOST_CHECK_EQUAL(matMap.atLocalBins(index3), vacuum.parameters());
0412 }
0413
0414 }  // namespace Acts::Test