diff --git a/sbncode/CAFMaker/RecoUtils/CMakeLists.txt b/sbncode/CAFMaker/RecoUtils/CMakeLists.txt index a992d7c63..295847117 100644 --- a/sbncode/CAFMaker/RecoUtils/CMakeLists.txt +++ b/sbncode/CAFMaker/RecoUtils/CMakeLists.txt @@ -11,4 +11,5 @@ art_make_library( LIBRARY_NAME caf_RecoUtils larsim::MCCheater_BackTrackerService_service larsim::MCCheater_ParticleInventoryService_service larcorealg::Geometry + sbnanaobj::StandardRecord ) diff --git a/sbncode/CMakeLists.txt b/sbncode/CMakeLists.txt index daed2590a..ab6eb0aab 100644 --- a/sbncode/CMakeLists.txt +++ b/sbncode/CMakeLists.txt @@ -21,6 +21,7 @@ add_subdirectory(CosmicID) add_subdirectory(DetSim) add_subdirectory(Cluster3D) add_subdirectory(HitFinder) +add_subdirectory(Utilities) # Supera # diff --git a/sbncode/DetSim/AdjustSimForTrigger_module.cc b/sbncode/DetSim/AdjustSimForTrigger_module.cc index 5f6c5eee9..624ba1c3e 100644 --- a/sbncode/DetSim/AdjustSimForTrigger_module.cc +++ b/sbncode/DetSim/AdjustSimForTrigger_module.cc @@ -1,11 +1,8 @@ -//////////////////////////////////////////////////////////////////////// -// Class: AdjustSimForTrigger -// Plugin Type: producer (Unknown Unknown) -// File: AdjustSimForTrigger_module.cc -// -// Generated at December 2023 by Bruce Howard (howard@fnal.gov) using -// cetskelgen. -//////////////////////////////////////////////////////////////////////// +/** + * @file sbncode/DetSim/AdjustSimForTrigger_module.cc + * @author Bruce Howard (howard@fnal.gov) + * @date December 2023 + */ #include "art/Framework/Core/EDProducer.h" #include "art/Framework/Core/ModuleMacros.h" @@ -13,6 +10,8 @@ #include "art/Framework/Principal/Handle.h" #include "art/Framework/Principal/Run.h" #include "art/Framework/Principal/SubRun.h" +#include "art/Persistency/Common/PtrMaker.h" +#include "canvas/Persistency/Common/Assns.h" #include "canvas/Utilities/InputTag.h" #include "fhiclcpp/ParameterSet.h" #include "messagefacility/MessageLogger/MessageLogger.h" @@ -24,12 +23,165 @@ #include "lardataobj/Simulation/SimEnergyDeposit.h" #include "lardataobj/Simulation/SimEnergyDepositLite.h" #include "lardataobj/Simulation/SimPhotons.h" +#include "sbnobj/ICARUS/PMT/Data/WaveformBaseline.h" -#include -#include - -class AdjustSimForTrigger; +#include "lardata/DetectorInfoServices/DetectorClocksService.h" +#include "sbncode/Utilities/AssnsUtils.h" // sbn::RebindAssociatedProducts() +#include // std::boolalpha +#include +#include +#include // std::move() +#include + + +/** + * @brief Applies a time shift to selected data products for simulation. + * + * This module produces a selection of simulation data products shifting their + * time reference to a new "event" (in the sense of something happening, not in + * the DAQ sense). + * + * The time point of the new reference event is specified in + * @ref DetectorClocksElectronicsTime "electronics time scale". + * Technically, this module applies a time shift so that the time of that + * reference event becomes the previous reference value in electronics time scale. + * + * More practically, we start from an electronics time scale where the reference + * time is set at `DetectorClocksData::TriggerTime()` (typically `1500.0` + * µs). + * That time reference may represent a hardware trigger, as the function name + * suggest, or anything else. At this point we introduce a new trigger time + * (from the simulation; e.g. `1502.0` µs) and we want that time to become + * the new reference (and "hardware trigger") of the electronics time scale. + * In this example, this module will add a shift of -2 µs to all the + * selected data products, so that the time of the new trigger is `1500.0` + * and all follows. + * + * @note The new trigger time is going to be whatever value is in + * `DetectorClocksData::TriggerTime()`, rather than an hard-coded value + * like `1500.0`. Therefore some care needs to be taken in the workflow + * to make sure that at the time of the execution of this module + * `DetectorClocksService` is yielding the desired new-reference value. + * Also, after this module is called `DetectorClocksService` itself may + * report obsolete values. Specifically, `DetectorClocksServiceStandard` + * will report the old `BeamGateTime()`, which was from ether a + * configuration parameter or from a trigger data product. + * + * This module reads the new reference time from a `raw::Trigger` object + * (assumed to hold a time on the current electronics scale), and, unless + * configured not to (`DropTriggerProduct`) it produces a new trigger data + * product with the shifted trigger information. This data product is suitable + * for configuring services like `DetectorClocksServiceStandard` for the + * following stages of the workflow. + * + * If the reference time is not valid, no shift is performed at all. + * A reference time is valid if it is neither the maximum nor the minimum value + * of a double (`std::numeric_limits::max()` and + * `std::numeric_limits::min()`). + * + * + * ### Optional shift + * + * It is possible to specify a fixed additional shift to the time reference. + * This may be useful for example if the simulation of the trigger yields + * exactly the time at which the triggering activity is detected, but the + * trigger hardware would take still some time in order to tag that time. + * Note, however, that the additional shift is also applied to the beam gate + * time. + * This additional time is specified via the `AdditionalOffset` configuration + * parameter. + * If the new reference time is invalid, no shift is applied and this offset is + * also ignored. + * + * + * Input + * ------ + * + * * `std::vector` (`InputTriggerLabel`): the first of the + * triggers in the collection will be used as a new reference. + * If the trigger time value is not valid, no shift at all will be performed. + * However, the beam gate time is still expected to be valid. + * An empty collection is not allowed, even when there is no valid trigger. + * + * * `art::Assns` + * (`BindWaveformBaselines`): the association of the original optical detector + * waveforms to their baselines. + * + * + * Output + * ------- + * + * For each enabled data product to be shifted, the corresponding shifted data + * product is produced, with elements in the same order as in the original + * collections. Normally, no associations are ported on. + * + * * `std::vector`: a collection of shifted triggers is produced; + * the first one is the shifted version of the reference trigger, which can + * then be used as new trigger data product e.g. for + * `DetectorClocksServiceStandard`. If the reference trigger time is not + * valid, the trigger time will be overwritten: the trigger time will be set + * to the value from `DetectorClocksService`, and the beam gate time will + * be the same as the input trigger object. This collection is produced by + * default, but it can be disabled via `DropTriggerProduct` configuration + * parameter. + * * `std::vector`: the beam gate used by the event + * generator, shifted. Generator times and particles from the detector + * simulation (GEANT4) are not shifted, but pretty much everything else is, + * including scintillation photons and energy depositions. Depending on which + * aspect of the simulation is being investigated, either the unshifted + * (input) or shifted (output of this module) gate needs to be used. + * * `art::Assns`: enabled only + * if the optical waveforms are being shifted _and_ an association data + * product name is specified (`BindWaveformBaselines`), it rebinds the + * original waveform baselines to the shifted waveforms. Note that neither the + * baseline nor the order nor the content of the waveforms change: only their + * start time (and implicitly the end time) does. + * + * + * Configuration parameters + * ------------------------- + * + * * `InputTriggerLabel` (input tag, mandatory): the tag of the trigger data + * product with the new reference time. It must be available and not empty. + * * `ShiftAuxDetIDE` (bool, default: `false`): enable the shifting of auxiliary + * detector simulation data product at `InitAuxDetSimChannelLabel`. + * * `InitAuxDetSimChannelLabel` (input tag): tag of the auxiliary detector + * simulation `sim::AuxDetSimChannel` product to be shifted. + * * `ShiftBeamGateInfo` (bool, default: `false`): enable the shifting of + * beam gate data product at `InitBeamGateInfoLabel`. + * * `InitBeamGateInfoLabel` (input tag): tag of the simulation beam gate data + * product to be shifted. This can be produced by a LArSoft generation module + * or by a trigger module. + * * `ShiftSimEnergyDeposits` (bool, default: `false`): enable the shifting of + * full energy deposit data product at `InitSimEnergyDepositLabel`. + * * `InitSimEnergyDepositLabel` (input tag): tag of the simulated energy + * deposition data product to be shifted. + * * `ShiftSimEnergyDepositLites` (bool, default: `false`): enable the shifting + * of lightweight energy deposit data product at + * `InitSimEnergyDepositLiteLabel`. + * * `InitSimEnergyDepositLiteLabel` (input tag): tag of the lightweight + * simulated energy deposition data product to be shifted. This reduced + * version is typically kept around for tracking back to truth information. + * * `ShiftSimPhotons` (bool, default: `false`): enable the shifting of + * scintillation photon data product at `InitSimPhotonsLabel`. + * * `InitSimPhotonsLabel` (input tag): tag of the simulated scintillation photon + * data product to be shifted. + * * `ShiftWaveforms` (bool, default: `false`): enable the shifting of optical + * detector waveform data product at `InitWaveformLabel`. + * * `InitWaveformLabel` (input tag): tag of the simulated optical detector + * waveform data product to be shifted. These waveforms may already be + * available if the worflow extracted the trigger time (new reference) out of + * them. + * * `BindWaveformBaselines` (input tag, default: `""`): tag of the association + * between the optical detector waveforms being shifted and their baselines. + * If empty (default), baseline associations will not be produced. + * * `AdditionalOffset` (real value, default: `0`): additional offset in + * microseconds to be added to the new reference trigger. + * * `DropTriggerProduct` (flag, default: `false`): if set, no shifted trigger + * data product will be produced. + * + */ class AdjustSimForTrigger : public art::EDProducer { public: explicit AdjustSimForTrigger(fhicl::ParameterSet const& p); @@ -59,13 +211,15 @@ class AdjustSimForTrigger : public art::EDProducer { bool fShiftSimEnergyDepositLites; bool fShiftSimPhotons; bool fShiftWaveforms; + art::InputTag fBindWaveformBaselines; ///< Tag of OpDetWaveform-baseline associations to be rebound. double fAdditionalOffset; + bool fDropTriggerProduct; ///< Do not put the shifted trigger data product into the event. static constexpr auto& kModuleName = "AdjustSimForTrigger"; }; AdjustSimForTrigger::AdjustSimForTrigger(fhicl::ParameterSet const& p) : EDProducer{p} - , fInputTriggerLabel{p.get("InputTriggerLabel", "undefined")} + , fInputTriggerLabel{p.get("InputTriggerLabel")} , fInitAuxDetSimChannelLabel(p.get("InitAuxDetSimChannelLabel", "undefined")) , fInitBeamGateInfoLabel{p.get("InitBeamGateInfoLabel", "undefined")} , fInitSimEnergyDepositLabel{p.get("InitSimEnergyDepositLabel", "undefined")} @@ -78,26 +232,36 @@ AdjustSimForTrigger::AdjustSimForTrigger(fhicl::ParameterSet const& p) , fShiftSimEnergyDepositLites{p.get("ShiftSimEnergyDepositLites", false)} , fShiftSimPhotons{p.get("ShiftSimPhotons", false)} , fShiftWaveforms{p.get("ShiftWaveforms", false)} + , fBindWaveformBaselines{p.get("BindWaveformBaselines", "")} , fAdditionalOffset{p.get("AdditionalOffset", 0.)} + , fDropTriggerProduct{p.get("DropTriggerProduct", false)} { if (!(fShiftSimEnergyDeposits || fShiftSimPhotons || fShiftWaveforms || fShiftAuxDetIDEs || fShiftBeamGateInfo || fShiftSimEnergyDepositLites)) { throw art::Exception(art::errors::EventProcessorFailure) << kModuleName << ": NO SHIFTS ENABLED!\n"; } + bool const doWaveformBaselines = fShiftWaveforms && !fBindWaveformBaselines.empty(); mf::LogInfo(kModuleName) << std::boolalpha << "SHIFTING AUXDETIDES? " << fShiftAuxDetIDEs << '\n' << "SHIFTING BEAMGATEINFO? " << fShiftBeamGateInfo << '\n' << "SHIFTING SIMENERGYDEPOSITS? " << fShiftSimEnergyDeposits << '\n' << "SHIFTING SIMENERGYDEPOSITLITES? " << fShiftSimEnergyDepositLites << '\n' << "SHIFTING SIMPHOTONS? " << fShiftSimPhotons << '\n' - << "SHIFTING OPDETWAVEFORMS? " << fShiftWaveforms; + << "SHIFTING OPDETWAVEFORMS? " << fShiftWaveforms << '\n' + << " ASSNS OPDETWAVEFORM-BASELINES? " << doWaveformBaselines + << (doWaveformBaselines? (" ('" + fBindWaveformBaselines.encode() + ")"): ""); + if (!fDropTriggerProduct) produces>(); if (fShiftAuxDetIDEs) { produces>(); } if (fShiftBeamGateInfo) { produces>(); } if (fShiftSimEnergyDeposits) { produces>(); } if (fShiftSimEnergyDepositLites) { produces>(); } if (fShiftSimPhotons) { produces>(); } - if (fShiftWaveforms) { produces>(); } + if (fShiftWaveforms) { + produces>(); + if (!fBindWaveformBaselines.empty()) + produces>(); + } } void AdjustSimForTrigger::produce(art::Event& e) @@ -123,13 +287,32 @@ void AdjustSimForTrigger::produce(art::Event& e) trigger.TriggerTime() < (std::numeric_limits::max() - std::numeric_limits::epsilon()); - const double timeShiftForTrigger_us = - hasValidTriggerTime ? clock_data.TriggerTime() - trigger.TriggerTime() + fAdditionalOffset : 0.; + const double newReferenceTime = + hasValidTriggerTime ? trigger.TriggerTime() - fAdditionalOffset : clock_data.TriggerTime(); + const double timeShiftForTrigger_us = clock_data.TriggerTime() - newReferenceTime; const double timeShiftForTrigger_ns = 1000. * timeShiftForTrigger_us; mf::LogInfo(kModuleName) << "FOR THIS EVENT THE TIME SHIFT BEING ASSUMED IS " << timeShiftForTrigger_ns << " ns ..."; + // Shifted trigger (beam gate info, optional, is later) + // sbn::ExtraTriggerInfo and raw::ExternalTrigger are not shifted (so far); + // it's debatable if they should be, since they only hold absolute timestamps + if (!fDropTriggerProduct) { + auto pShiftedTriggers = std::make_unique>(); + for (raw::Trigger const& unshiftedTrigger: triggers) { // ok, we required there is just one + pShiftedTriggers->emplace_back( + unshiftedTrigger.TriggerNumber(), + hasValidTriggerTime // trigger_time + ? (unshiftedTrigger.TriggerTime() + timeShiftForTrigger_us) + : clock_data.TriggerTime(), + unshiftedTrigger.BeamGateTime() + timeShiftForTrigger_us,// beamgate_time + unshiftedTrigger.TriggerBits() + ); + } // for all triggers + e.put(std::move(pShiftedTriggers)); + } // if produce shifted trigger + // Loop over the sim::AuxDetIDE and shift time BACK by the TRIGGER if (fShiftAuxDetIDEs) { auto const& simChannels = @@ -253,6 +436,21 @@ void AdjustSimForTrigger::produce(art::Event& e) waveform.SetTimeStamp(waveform.TimeStamp() + timeShiftForTrigger_us); } e.put(std::move(pWaveforms)); + + if (!fBindWaveformBaselines.empty()) { + // given that the shifting is one-to-one, rebinding is just replacing + // each existing waveform pointer with one to the new waveform in the same position + auto const& waveformBaselineAssns + = e.getProduct>(fBindWaveformBaselines); + art::PtrMaker const makeWaveformPtr{ e }; + + auto pWaveformBaselineAssns + = std::make_unique> + (sbn::RebindAssociatedProducts(waveformBaselineAssns, makeWaveformPtr)); + + e.put(std::move(pWaveformBaselineAssns)); + } // if rebinding associations + } } diff --git a/sbncode/DetSim/CMakeLists.txt b/sbncode/DetSim/CMakeLists.txt index 113aa375a..813562b60 100644 --- a/sbncode/DetSim/CMakeLists.txt +++ b/sbncode/DetSim/CMakeLists.txt @@ -1,34 +1,11 @@ set( MODULE_LIBRARIES + lardata::DetectorClocksService + larcore::headers + sbncode::Utilities larcorealg::Geometry - larcore::Geometry_Geometry_service - larsim::Simulation - nug4::ParticleNavigation - lardataobj::Simulation - lardata::Utilities - lardataalg::DetectorInfo - larevt::Filters lardataobj::RawData - larevt::CalibrationDBI_Providers - nurandom::RandomUtils_NuRandomService_service - art::Framework_Core - art::Framework_Principal - art::Framework_Services_Registry - art_root_io::tfile_support ROOT::Core - art::Framework_Services_Optional_RandomNumberGenerator_service - art_root_io::TFileService_service - art::Persistency_Common - art::Persistency_Provenance - art::Utilities - canvas::canvas - messagefacility::MF_MessageLogger - messagefacility::headers - fhiclcpp::fhiclcpp - cetlib::cetlib - CLHEP::Random - ROOT::Geom - ROOT::XMLIO - ROOT::Gdml - FFTW3::FFTW3 + lardataobj::Simulation + sbnobj::ICARUS_PMT_Data ) cet_build_plugin(AdjustSimForTrigger art::module LIBRARIES ${MODULE_LIBRARIES}) diff --git a/sbncode/Utilities/AssnsUtils.h b/sbncode/Utilities/AssnsUtils.h new file mode 100644 index 000000000..46be1e4c5 --- /dev/null +++ b/sbncode/Utilities/AssnsUtils.h @@ -0,0 +1,144 @@ +/** + * @file sbncode/Utilities/AssnsUtils.h + * @brief Framework-level utilities dealing with associations. + * @date April 2, 2026 + * @author Gianluca Petrillo (petrillo@slac.stanford.edu) + * + * This library is currently header-only. + */ + +#ifndef SBNCODE_UTILITIES_ASSNSUTILS_H +#define SBNCODE_UTILITIES_ASSNSUTILS_H + +// framework libraries +#include "art/Persistency/Common/PtrMaker.h" +#include "canvas/Persistency/Common/Ptr.h" +#include "canvas/Persistency/Common/Assns.h" + +// C/C++ libraries +#include // std::size_t +#include // std::get() +#include // std::is_same_v, ... + + +namespace sbn { + /** + * @brief Creates an association of a new product mirroring an existing one. + * @tparam DiscardData (default: `false`) whether to omit copying the data + * @tparam L type of left-side data in the association + * @tparam R type of right-side data in the association + * @tparam Data type of metadata in the association + * @tparam LorR type whose pointer is being rebound + * @param assns the existing association + * @param reboundPtrMaker tool to create a _art_ pointer to the new product + * @return the newly created association + * + * The existing association `assns` is replicated replacing all the pointers + * in the rebinding side (left or right) with their corresponding pointers + * from another collection of the same type. + * + * The application is for one-to-one processing of data product elements, + * where the element resulting from the processing is still associated to the + * same object as the processed element. + * For example, let's have a reprocessing that modifies the signal on a + * optical detector waveform (`raw::OpDetWaveform`), leaving the baseline + * unchanged. + * ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~{.cpp} + * using Assns_t = art::Assns; + * + * Assns_t const& oldAssns = event.getProduct(assnsTag); + * art::PtrMaker const makeWaveformPtr{ event }; + * + * auto newAssns + * = std::make_unique(sbn::RebindAssociatedProducts(oldAssns, makeWaveformPtr)); + * + * event.put(std::move(pWaveformBaselineAssns)); + * ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ + * The alias `Assns_t` was used only to make the example more readable. + * + * As shown in the example, which side to rebind is determined by the type + * of the `art::PtrMaker` object. If its operand type is neither of the + * association types, a compilation error will occur. + * + * Note that the invalid pointers are also rebound, but their rebound pointer + * should be as invalid as the old one. + * + * If metadata is present in the associations, it is _copied_ element by + * element into the new one. It is possible to omit the metadata by calling + * with `DiscardData` set to `true` (`RebindAssociatedProducts(...)`), + * in which case the returned type will always be an association without data. + * + */ + template < + bool DiscardData = false, + typename L, typename R, typename Data, + typename LorR + > + art::Assns> + RebindAssociatedProducts( + art::Assns const& assns, + art::PtrMaker const& reboundPtrMaker + ); + +} // namespace sbn + + +// ----------------------------------------------------------------------------- +// --- Template implementation +// ----------------------------------------------------------------------------- +template < + bool DiscardData /* = false */, + typename L, typename R, typename Data, + typename LorR + > +art::Assns> +sbn::RebindAssociatedProducts( + art::Assns const& assns, + art::PtrMaker const& reboundPtrMaker +) { + static_assert(std::is_same_v || std::is_same_v, + "The PtrMaker tool must operate on either of the association data types."); + + constexpr bool RebindLeft = std::is_same_v; + constexpr bool hasMetadata = !std::is_void_v; + constexpr bool wantsMetadata = hasMetadata && !DiscardData; + + using DestData_t = std::conditional_t; + + /* + * How to deal with invalid pointers? (there should be none though!) + * a) preserve them (drawback: different product IDs are mixed in) + * b) replace them with default-constructed (which is also invalid) + * c) rebind it (drawback: rely on PtrMaker doing the right thing with them) + * Here I chose (c). + */ + auto rebind = [&reboundPtrMaker](auto& ptr) + { /* if (ptr) */ ptr = reboundPtrMaker(ptr.key()); }; + + art::Assns newAssns; + for (auto elements: assns) { + + art::Ptr leftPtr = std::move(std::get<0>(elements)); + art::Ptr rightPtr = std::move(std::get<1>(elements)); + + if constexpr(RebindLeft) { rebind(leftPtr); } + else { rebind(rightPtr); } + + if constexpr(wantsMetadata) { + newAssns.addSingle + (std::move(leftPtr), std::move(rightPtr), std::move(std::get<2>(elements))); + } + else { + newAssns.addSingle(std::move(leftPtr), std::move(rightPtr)); + } + + } // for assn pairs + + return newAssns; + +} // RebindAssociatedProducts() + + +// ----------------------------------------------------------------------------- + +#endif // SBNCODE_UTILITIES_ASSNSUTILS_H diff --git a/sbncode/Utilities/CMakeLists.txt b/sbncode/Utilities/CMakeLists.txt new file mode 100644 index 000000000..a94a8b369 --- /dev/null +++ b/sbncode/Utilities/CMakeLists.txt @@ -0,0 +1,11 @@ +cet_make_library( + INTERFACE + HEADERS_TARGET_ONLY + LIBRARIES + art::Persistency_Common + canvas::canvas + ) + +install_headers() +install_source() +