Back to home page

EIC code displayed by LXR

 
 

    


File indexing completed on 2026-07-26 08:22:26

0001 /** TRACCC library, part of the ACTS project (R&D line)
0002  *
0003  * (c) 2023-2026 CERN for the benefit of the ACTS project
0004  *
0005  * Mozilla Public License Version 2.0
0006  */
0007 
0008 // Project include(s).
0009 #include "traccc/bfield/construct_const_bfield.hpp"
0010 #include "traccc/bfield/magnetic_field_types.hpp"
0011 #include "traccc/cuda/finding/combinatorial_kalman_filter_algorithm.hpp"
0012 #include "traccc/finding/combinatorial_kalman_filter_algorithm.hpp"
0013 #include "traccc/io/read_detector.hpp"
0014 #include "traccc/io/read_measurements.hpp"
0015 #include "traccc/io/utils.hpp"
0016 #include "traccc/performance/container_comparator.hpp"
0017 #include "traccc/simulation/event_generators.hpp"
0018 #include "traccc/simulation/simulator.hpp"
0019 #include "traccc/utils/event_data.hpp"
0020 #include "traccc/utils/ranges.hpp"
0021 #include "traccc/utils/seed_generator.hpp"
0022 
0023 // Test include(s).
0024 #include "tests/ckf_toy_detector_test.hpp"
0025 
0026 // VecMem include(s).
0027 #include <vecmem/memory/cuda/device_memory_resource.hpp>
0028 #include <vecmem/memory/cuda/managed_memory_resource.hpp>
0029 #include <vecmem/memory/host_memory_resource.hpp>
0030 #include <vecmem/utils/cuda/async_copy.hpp>
0031 #include <vecmem/utils/cuda/stream_wrapper.hpp>
0032 
0033 // GTest include(s).
0034 #include <gtest/gtest.h>
0035 
0036 // System include(s).
0037 #include <filesystem>
0038 #include <string>
0039 
0040 using namespace traccc;
0041 
0042 TEST_P(CkfToyDetectorTests, Run) {
0043   // Get the parameters
0044   const std::string name = std::get<0>(GetParam());
0045   const traccc::pdg_particle<scalar> ptc = std::get<6>(GetParam());
0046   const unsigned int n_truth_tracks = std::get<7>(GetParam());
0047   const unsigned int n_events = std::get<8>(GetParam());
0048   const bool random_charge = std::get<9>(GetParam());
0049 
0050   constexpr bool use_material_maps = false;
0051   WriteDetector(use_material_maps, name);
0052 
0053   /*****************************
0054    * Build a toy detector
0055    *****************************/
0056 
0057   // Memory resources used by the application.
0058   vecmem::host_memory_resource host_mr;
0059   vecmem::cuda::device_memory_resource device_mr;
0060   traccc::memory_resource mr{device_mr, &host_mr};
0061   vecmem::cuda::managed_memory_resource mng_mr;
0062   vecmem::copy host_copy;
0063 
0064   // Read back detector file
0065   const std::string path = name + "/";
0066   detray::io::detector_reader_config reader_cfg{};
0067   reader_cfg.add_file(path + "toy_detector_geometry.json")
0068       .add_file(path + "toy_detector_surface_grids.json")
0069       .add_file(path + "toy_detector_homogeneous_material.json")
0070       .do_check(true);
0071 
0072   auto [io_det, names] =
0073       detray::io::read_detector<traccc::default_detector::host>(host_mr,
0074                                                                 reader_cfg);
0075   traccc::host_detector detector{};
0076   detector.template set<
0077       traccc::detector_traits<traccc::default_detector::host::metadata>>(
0078       std::move(io_det));
0079 
0080   traccc::detector_buffer detector_buffer =
0081       traccc::buffer_from_host_detector(detector, mng_mr, host_copy);
0082 
0083   const auto field = traccc::construct_const_bfield(B);
0084 
0085   /***************************
0086    * Generate simulation data
0087    ***************************/
0088 
0089   // Track generator
0090   using generator_type =
0091       detray::random_track_generator<traccc::free_track_parameters<>,
0092                                      uniform_gen_t>;
0093   generator_type::configuration gen_cfg{};
0094   gen_cfg.n_tracks(n_truth_tracks);
0095   gen_cfg.origin(std::get<1>(GetParam()));
0096   gen_cfg.origin_stddev(std::get<2>(GetParam()));
0097   gen_cfg.phi_range(std::get<5>(GetParam()));
0098   gen_cfg.eta_range(std::get<4>(GetParam()));
0099   gen_cfg.mom_range(std::get<3>(GetParam()));
0100   gen_cfg.randomize_charge(random_charge);
0101   gen_cfg.seed(42);
0102   generator_type generator(gen_cfg);
0103 
0104   // Smearing value for measurements
0105   traccc::measurement_smearer<traccc::default_algebra> meas_smearer(
0106       smearing[0], smearing[1]);
0107 
0108   using writer_type = traccc::smearing_writer<
0109       traccc::measurement_smearer<traccc::default_algebra>>;
0110 
0111   typename writer_type::config smearer_writer_cfg{meas_smearer};
0112   traccc::seed_generator<host_detector_type>::config seed_cfg{};
0113   seed_cfg.initial_sigmas = stddevs;
0114 
0115   // Run simulator
0116   const std::string full_path = io::data_directory() + path;
0117   std::filesystem::create_directories(full_path);
0118   auto sim = traccc::simulator<host_detector_type, b_field_t, generator_type,
0119                                writer_type>(
0120       ptc, n_events, detector.as<detector_traits>(),
0121       field.as_field<traccc::const_bfield_backend_t<traccc::scalar>>(),
0122       std::move(generator), std::move(smearer_writer_cfg), full_path);
0123   sim.get_config().propagation.navigation.search_window = search_window;
0124   sim.run();
0125 
0126   /*****************************
0127    * Do the reconstruction
0128    *****************************/
0129 
0130   // Stream object
0131   vecmem::cuda::stream_wrapper vecmem_stream;
0132   traccc::cuda::stream_wrapper stream{vecmem_stream.stream()};
0133 
0134   // Copy objects
0135   vecmem::cuda::async_copy copy{stream.cudaStream()};
0136 
0137   // Seed generator
0138   seed_generator<host_detector_type> sg(detector.as<detector_traits>(),
0139                                         seed_cfg);
0140 
0141   // Finding algorithm configuration
0142   typename traccc::cuda::combinatorial_kalman_filter_algorithm::config_type cfg;
0143   cfg.ptc_hypothesis = ptc;
0144   cfg.max_num_branches_per_seed = 500;
0145   cfg.max_num_branches_per_surface = 2;
0146   cfg.chi2_max = 10.f;
0147   cfg.propagation.navigation.search_window = search_window;
0148   cfg.run_smoother = smoother_type::e_none;
0149 
0150   // Finding algorithm object
0151   traccc::host::combinatorial_kalman_filter_algorithm host_finding(cfg,
0152                                                                    host_mr);
0153 
0154   // Finding algorithm object
0155   traccc::cuda::combinatorial_kalman_filter_algorithm device_finding(
0156       cfg, mr, copy, stream);
0157 
0158   // Iterate over events
0159   for (std::size_t i_evt = 0; i_evt < n_events; i_evt++) {
0160     // Truth Track Candidates
0161     traccc::event_data evt_data(path, i_evt, host_mr);
0162 
0163     traccc::edm::measurement_collection::host truth_measurements{host_mr};
0164     traccc::edm::track_container<traccc::default_algebra>::host
0165         truth_track_candidates{host_mr};
0166     evt_data.generate_truth_candidates(truth_track_candidates,
0167                                        truth_measurements, sg, host_mr);
0168     truth_track_candidates.measurements = vecmem::get_data(truth_measurements);
0169 
0170     ASSERT_EQ(truth_track_candidates.tracks.size(), n_truth_tracks);
0171 
0172     // Prepare truth seeds
0173     traccc::bound_track_parameters_collection_types::host seeds(&host_mr);
0174     for (unsigned int i_trk = 0; i_trk < n_truth_tracks; i_trk++) {
0175       seeds.push_back(truth_track_candidates.tracks.at(i_trk).params());
0176     }
0177     ASSERT_EQ(seeds.size(), n_truth_tracks);
0178 
0179     traccc::bound_track_parameters_collection_types::buffer seeds_buffer{
0180         static_cast<unsigned int>(seeds.size()), mr.main};
0181     copy.setup(seeds_buffer)->wait();
0182     copy(vecmem::get_data(seeds), seeds_buffer,
0183          vecmem::copy::type::host_to_device)
0184         ->wait();
0185 
0186     // Prepare the measurements
0187     traccc::edm::measurement_collection::host measurements_per_event{host_mr};
0188     traccc::io::read_measurements(measurements_per_event, i_evt, path);
0189 
0190     traccc::edm::measurement_collection::buffer measurements_buffer(
0191         static_cast<unsigned int>(measurements_per_event.size()), mr.main);
0192     copy.setup(measurements_buffer)->wait();
0193     copy(vecmem::get_data(measurements_per_event), measurements_buffer)->wait();
0194 
0195     // Run host finding
0196     auto track_candidates =
0197         host_finding(detector, field, vecmem::get_data(measurements_per_event),
0198                      vecmem::get_data(seeds));
0199 
0200     // Run device finding
0201     traccc::edm::track_container<traccc::default_algebra>::buffer
0202         track_candidates_cuda_buffer = device_finding(
0203             detector_buffer, field, measurements_buffer, seeds_buffer);
0204 
0205     traccc::edm::track_collection<traccc::default_algebra>::host
0206         track_candidates_cuda{host_mr};
0207     copy(track_candidates_cuda_buffer.tracks, track_candidates_cuda)->wait();
0208 
0209     // Simple check
0210     ASSERT_GE(track_candidates.tracks.size(), n_truth_tracks)
0211         << "No. tracks (host): " << track_candidates.tracks.size() << "/"
0212         << n_truth_tracks;
0213     ASSERT_LE(static_cast<double>(
0214                   std::llabs(static_cast<long>(track_candidates.tracks.size()) -
0215                              static_cast<long>(track_candidates_cuda.size()))) /
0216                   static_cast<double>(track_candidates.tracks.size()),
0217               0.001f)
0218         << "No. tracks (host): " << track_candidates.tracks.size() << "/"
0219         << n_truth_tracks
0220         << "\nNo. tracks (device): " << track_candidates_cuda.size() << "/"
0221         << n_truth_tracks;
0222 
0223     // Make sure that the outputs from cpu and cuda CKF are equivalent
0224     unsigned int n_matches = 0u;
0225     for (unsigned int i = 0u; i < track_candidates.tracks.size(); i++) {
0226       traccc::details::is_same_object<traccc::edm::track_collection<
0227           traccc::default_algebra>::host::const_proxy_type>
0228           iso{track_candidates.measurements, track_candidates.measurements,
0229               vecmem::get_data(track_candidates.states),
0230               vecmem::get_data(track_candidates.states),
0231               track_candidates.tracks.at(i)};
0232 
0233       for (unsigned int j = 0u; j < track_candidates_cuda.size(); j++) {
0234         if (iso(track_candidates_cuda.at(j))) {
0235           n_matches++;
0236           break;
0237         }
0238       }
0239     }
0240 
0241     float matching_rate = float(n_matches) / static_cast<float>(std::max(
0242                                                  track_candidates.tracks.size(),
0243                                                  track_candidates_cuda.size()));
0244     EXPECT_GE(matching_rate, 0.998f);
0245   }
0246 }
0247 
0248 INSTANTIATE_TEST_SUITE_P(
0249     CUDACkfToyDetectorValidation, CkfToyDetectorTests,
0250     ::testing::Values(
0251         std::make_tuple("toy_n_particles_1",
0252                         std::array<scalar, 3u>{0.f, 0.f, 0.f},
0253                         std::array<scalar, 3u>{0.f, 0.f, 0.f},
0254                         std::array<scalar, 2u>{1.f, 100.f},
0255                         std::array<scalar, 2u>{-4.f, 4.f},
0256                         std::array<scalar, 2u>{-traccc::constant<scalar>::pi,
0257                                                traccc::constant<scalar>::pi},
0258                         traccc::muon<scalar>(), 1, 1, false),
0259         std::make_tuple("toy_n_particles_10000",
0260                         std::array<scalar, 3u>{0.f, 0.f, 0.f},
0261                         std::array<scalar, 3u>{0.f, 0.f, 0.f},
0262                         std::array<scalar, 2u>{1.f, 100.f},
0263                         std::array<scalar, 2u>{-4.f, 4.f},
0264                         std::array<scalar, 2u>{-traccc::constant<scalar>::pi,
0265                                                traccc::constant<scalar>::pi},
0266                         traccc::muon<scalar>(), 10000, 1, false),
0267         std::make_tuple("toy_n_particles_10000_random_charge",
0268                         std::array<scalar, 3u>{0.f, 0.f, 0.f},
0269                         std::array<scalar, 3u>{0.f, 0.f, 0.f},
0270                         std::array<scalar, 2u>{1.f, 100.f},
0271                         std::array<scalar, 2u>{-4.f, 4.f},
0272                         std::array<scalar, 2u>{-traccc::constant<scalar>::pi,
0273                                                traccc::constant<scalar>::pi},
0274                         traccc::muon<scalar>(), 10000, 1, true)));