Tests/Physics/ScatteringTests.cpp

0001 // This file is part of the Acts project.
0002 //
0003 // Copyright (C) 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 ///  Boost include(s)
0010 #define BOOST_TEST_MODULE AbortList Tests
0011
0012 #include <boost/test/included/unit_test.hpp>
0013 // leave blank line
0014
0015 #include <boost/test/data/test_case.hpp>
0016 // leave blank line
0017
0018 #include <boost/test/output_test_stream.hpp>
0019 // leave blank line
0020
0021 #include "Acts/Material/Material.hpp"
0022 #include "Acts/Material/MaterialProperties.hpp"
0023
0024 #include "Fatras/Kernel/PhysicsList.hpp"
0025 #include "Fatras/Kernel/Process.hpp"
0026 #include "Fatras/Physics/Scattering/GaussianMixture.hpp"
0027 #include "Fatras/Physics/Scattering/GeneralMixture.hpp"
0028 #include "Fatras/Physics/Scattering/Highland.hpp"
0029 #include "Fatras/Physics/Scattering/Scattering.hpp"
0030 #include "Particle.hpp"
0031 #include <fstream>
0032 #include <random>
0033
0034 namespace bdata = boost::unit_test::data;
0035 namespace tt = boost::test_tools;
0036
0037 namespace Fatras {
0038
0039 namespace Test {
0040
0041 // the generator
0042 typedef std::mt19937 Generator;
0043
0044 // standard generator
0045 Generator generator;
0046
0047 /// The selector
0048 struct Selector {
0049
0050   /// call operator
0051   template <typename detector_t, typename particle_t>
0052   bool operator()(const detector_t &, const particle_t &) const {
0053     return true;
0054   }
0055 };
0056
0057 // some material
0058 Acts::Material berilium = Acts::Material(352.8, 407., 9.012, 4., 1.848e-3);
0059
0060 bool write_csv = true;
0061
0062 std::ofstream os("ScatteringAngles.csv",
0063                  std::ofstream::out | std::ofstream::trunc);
0064
0065 /// Test the scattering implementation
0066 BOOST_DATA_TEST_CASE(
0067     HighlandScattering_test_,
0068     bdata::random(
0069         (bdata::seed = 20,
0070          bdata::distribution = std::uniform_real_distribution<>(0., 1.))) ^
0071         bdata::random(
0072             (bdata::seed = 21,
0073              bdata::distribution = std::uniform_real_distribution<>(0., 1.))) ^
0074         bdata::random(
0075             (bdata::seed = 22,
0076              bdata::distribution = std::uniform_real_distribution<>(0., 1.))) ^
0077         bdata::random((bdata::seed = 23,
0078                        bdata::distribution =
0079                            std::uniform_real_distribution<>(0.5, 10.5))) ^
0080         bdata::xrange(10000),
0081     x, y, z, p, index) {
0082
0083   // a detector with 1 mm Be
0084   Acts::MaterialProperties detector(berilium, 1. * Acts::units::_mm);
0085
0086   // create the particle and set the momentum
0087   /// position at 0.,0.,0
0088   Acts::Vector3D position{0., 0., 0.};
0089   // p of 1 GeV
0090   Acts::Vector3D momentum =
0091       p * Acts::units::_GeV * Acts::Vector3D(x, y, z).normalized();
0092   // positively charged
0093   double q = -1.;
0094   double m = 105.658367 * Acts::units::_MeV; // muon mass
0095
0096   // create the particle
0097   Particle particle(position, momentum, m, q, 13, 1);
0098
0099   // make the highland scatterer
0100   Highland hscat;
0101   GaussianMixture gamscat;
0102   GeneralMixture genscat;
0103
0104   double hsr = hscat(generator, detector, particle);
0105   double gamr = gamscat(generator, detector, particle);
0106   double genr = genscat(generator, detector, particle);
0107
0108   BOOST_CHECK(hsr != 0.);
0109
0110   // Run the Scatterer as a plug-in process
0111   Scattering<Highland> hsScattering;
0112   auto out = hsScattering(generator, detector, particle);
0113   BOOST_CHECK(!out.size());
0114
0115   // Run the Scatterer as a physics list
0116   typedef Selector All;
0117   std::vector<Particle> outgoing;
0118   typedef Process<Scattering<Highland>, All, All, All> HighlandProcess;
0119   PhysicsList<HighlandProcess> hsPhysicsList;
0120   hsPhysicsList(generator, detector, particle, outgoing);
0121   BOOST_CHECK(!outgoing.size());
0122
0123   // write out a csv file
0124   if (write_csv) {
0125     if (!index)
0126       os << "p,highland,gaussian_mixture,general_mixture" << '\n';
0127     os << particle.p() << "," << hsr << "," << gamr << "," << genr << '\n';
0128   }
0129 }
0130
0131 } // namespace Test
0132
0133 } // namespace Fatras