Feature/gp overlay fix

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
                 )

sbncode/CMakeLists.txt

@@ -21,6 +21,7 @@ add_subdirectory(CosmicID)
 add_subdirectory(DetSim)
 add_subdirectory(Cluster3D)
 add_subdirectory(HitFinder)
+add_subdirectory(Utilities)
 
 # Supera
 #

sbncode/DetSim/AdjustSimForTrigger_module.cc

@@ -1,18 +1,17 @@
-////////////////////////////////////////////////////////////////////////
-// 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"
 #include "art/Framework/Principal/Event.h"
 #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 <lardata/DetectorInfoServices/DetectorClocksService.h>
-#include <memory>
-
-class AdjustSimForTrigger;
+#include "lardata/DetectorInfoServices/DetectorClocksService.h"
+#include "sbncode/Utilities/AssnsUtils.h" // sbn::RebindAssociatedProducts()
 
+#include <ios> // std::boolalpha
+#include <limits>
+#include <memory>
+#include <utility> // std::move()
+#include <vector>
+
+
+/**
+ * @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`
+ * &micro;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` &micro;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 &micro;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<double>::max()` and
+ * `std::numeric_limits<double>::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<raw::Trigger>` (`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<raw::OpDetWaveform, icarus::WaveformBaseline>`
+ *   (`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<raw::Trigger>`: 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<sim::BeamGateInfo>`: 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<raw::OpDetWaveform, icarus::WaveformBaseline>`: 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<art::InputTag>("InputTriggerLabel", "undefined")}
+  , fInputTriggerLabel{p.get<art::InputTag>("InputTriggerLabel")}
   , fInitAuxDetSimChannelLabel(p.get<art::InputTag>("InitAuxDetSimChannelLabel", "undefined"))
   , fInitBeamGateInfoLabel{p.get<art::InputTag>("InitBeamGateInfoLabel", "undefined")}
   , fInitSimEnergyDepositLabel{p.get<art::InputTag>("InitSimEnergyDepositLabel", "undefined")}
@@ -78,26 +232,36 @@ AdjustSimForTrigger::AdjustSimForTrigger(fhicl::ParameterSet const& p)
   , fShiftSimEnergyDepositLites{p.get<bool>("ShiftSimEnergyDepositLites", false)}
   , fShiftSimPhotons{p.get<bool>("ShiftSimPhotons", false)}
   , fShiftWaveforms{p.get<bool>("ShiftWaveforms", false)}
+  , fBindWaveformBaselines{p.get<art::InputTag>("BindWaveformBaselines", "")}
   , fAdditionalOffset{p.get<double>("AdditionalOffset", 0.)}
+  , fDropTriggerProduct{p.get<bool>("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<std::vector<raw::Trigger>>();
   if (fShiftAuxDetIDEs) { produces<std::vector<sim::AuxDetSimChannel>>(); }
   if (fShiftBeamGateInfo) { produces<std::vector<sim::BeamGateInfo>>(); }
   if (fShiftSimEnergyDeposits) { produces<std::vector<sim::SimEnergyDeposit>>(); }
   if (fShiftSimEnergyDepositLites) { produces<std::vector<sim::SimEnergyDepositLite>>(); }
   if (fShiftSimPhotons) { produces<std::vector<sim::SimPhotons>>(); }
-  if (fShiftWaveforms) { produces<std::vector<raw::OpDetWaveform>>(); }
+  if (fShiftWaveforms) { 
+    produces<std::vector<raw::OpDetWaveform>>(); 
+    if (!fBindWaveformBaselines.empty())
+      produces<art::Assns<raw::OpDetWaveform, icarus::WaveformBaseline>>(); 
+  }
 }
 
 void AdjustSimForTrigger::produce(art::Event& e)
@@ -123,13 +287,32 @@ void AdjustSimForTrigger::produce(art::Event& e)
     trigger.TriggerTime() <
       (std::numeric_limits<double>::max() - std::numeric_limits<double>::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<std::vector<raw::Trigger>>();
+    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<art::Assns<raw::OpDetWaveform, icarus::WaveformBaseline>>(fBindWaveformBaselines);
+      art::PtrMaker<raw::OpDetWaveform> const makeWaveformPtr{ e };
+
+      auto pWaveformBaselineAssns
+        = std::make_unique<art::Assns<raw::OpDetWaveform, icarus::WaveformBaseline>>
+        (sbn::RebindAssociatedProducts(waveformBaselineAssns, makeWaveformPtr));
+
+      e.put(std::move(pWaveformBaselineAssns));
+    } // if rebinding associations
+
   }
 }
 

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})