tests/algorithms_test/digi_PulseGeneration.cc

0001 // SPDX-License-Identifier: LGPL-3.0-or-later
0002 // Copyright (C) 2025, Dmitry Kalinkin, Simon Gardner
0003
0004 #include <catch2/catch_test_macros.hpp>
0005 #include <catch2/generators/catch_generators.hpp>
0006 #include <edm4eic/EDM4eicVersion.h>
0007 #include <edm4eic/unit_system.h>
0008 #include <edm4hep/SimTrackerHitCollection.h>
0009 #if EDM4EIC_VERSION_MAJOR > 8 || (EDM4EIC_VERSION_MAJOR == 8 && EDM4EIC_VERSION_MINOR >= 1)
0010 #include <edm4eic/SimPulseCollection.h>
0011 #else
0012 #include <edm4hep/TimeSeriesCollection.h>
0013 #endif
0014 #include <podio/RelationRange.h>
0015 #include <cmath>
0016 #include <cstddef>
0017 #include <memory>
0018 #include <string>
0019 #include <tuple>
0020 #include <vector>
0021
0022 #include "algorithms/digi/PulseGeneration.h"
0023 #include "algorithms/digi/PulseGenerationConfig.h"
0024
0025 #if EDM4EIC_VERSION_MAJOR > 8 || (EDM4EIC_VERSION_MAJOR == 8 && EDM4EIC_VERSION_MINOR >= 1)
0026 using PulseType = edm4eic::SimPulse;
0027 #else
0028 using PulseType = edm4hep::TimeSeries;
0029 #endif
0030
0031 TEST_CASE("PulseGeneration generates correct number of pulses", "[PulseGeneration]") {
0032
0033   eicrecon::PulseGeneration<edm4hep::SimTrackerHit> algo("PulseGeneration");
0034   eicrecon::PulseGenerationConfig cfg;
0035   cfg.pulse_shape_function = "LandauPulse"; // Example pulse shape
0036   cfg.pulse_shape_params   = {1.0, 1.0};    // Example parameters for the pulse shape
0037   cfg.ignore_thres         = 1;
0038
0039   algo.applyConfig(cfg);
0040   algo.init();
0041
0042   SECTION("from different hits") {
0043     std::size_t nHits = GENERATE(1, 2, 3);
0044
0045     edm4hep::SimTrackerHitCollection hits_coll;
0046
0047     for (std::size_t i = 0; i < nHits; i++) {
0048       hits_coll.create(12345 + i, 10.0, 5.0); // cellID, charge, time
0049     }
0050
0051     auto pulses = std::make_unique<PulseType::collection_type>();
0052
0053     auto input  = std::make_tuple(&hits_coll);
0054     auto output = std::make_tuple(pulses.get());
0055
0056     algo.process(input, output);
0057
0058     REQUIRE(pulses->size() == nHits);
0059     REQUIRE((*pulses)[0].getCellID() == 12345);
0060     if (nHits > 1) {
0061       REQUIRE((*pulses)[1].getCellID() == 12346);
0062     }
0063     if (nHits > 2) {
0064       REQUIRE((*pulses)[2].getCellID() == 12347);
0065     }
0066   }
0067 }
0068
0069 TEST_CASE("Test the EvaluatorSvc pulse generation with a square pulse", "[PulseGeneration]") {
0070
0071   eicrecon::PulseGeneration<edm4hep::SimTrackerHit> algo("PulseGeneration");
0072   eicrecon::PulseGenerationConfig cfg;
0073
0074   // Square wave expression
0075   std::string expression = "(time >= param0 && time < param1) ? charge : 0";
0076
0077   double startTime      = 0.0 * edm4eic::unit::ns;
0078   double endTime        = 1.0 * edm4eic::unit::ns;
0079   std::size_t nTimeBins = 10;
0080   double timeStep       = (endTime - startTime) / nTimeBins;
0081
0082   cfg.pulse_shape_function = expression;
0083   cfg.pulse_shape_params   = {startTime, endTime}; // Example parameters for the square pulse
0084   cfg.ignore_thres         = 1;
0085   cfg.timestep             = timeStep;
0086   cfg.min_sampling_time    = startTime + timeStep;
0087
0088   algo.applyConfig(cfg);
0089   algo.init();
0090
0091   double charge      = 10.0 * cfg.ignore_thres;
0092   double time        = GENERATE_COPY(0.0, 0.5 * timeStep, timeStep);
0093   float rounded_time = std::floor(time / timeStep) * timeStep;
0094
0095   edm4hep::SimTrackerHitCollection hits_coll;
0096   hits_coll.create(12345, charge, time); // cellID, charge, time
0097
0098   auto pulses = std::make_unique<PulseType::collection_type>();
0099
0100   auto input  = std::make_tuple(&hits_coll);
0101   auto output = std::make_tuple(pulses.get());
0102
0103   algo.process(input, output);
0104
0105   REQUIRE(pulses->size() == 1);
0106   REQUIRE((*pulses)[0].getCellID() == 12345);
0107   REQUIRE((*pulses)[0].getTime() == rounded_time);
0108   auto amplitudes = (*pulses)[0].getAmplitude();
0109   REQUIRE(amplitudes.size() == nTimeBins); // Two time bins for the square pulse
0110   for (auto amplitude : amplitudes) {
0111     REQUIRE(amplitude == charge); // All time bins should be zero
0112   }
0113 }