LArSoft  v09_90_00
Liquid Argon Software toolkit - https://larsoft.org/
OpHitFinder_module.cc
Go to the documentation of this file.
1 // -*- mode: c++; c-basic-offset: 2; -*-
2 // Gleb Sinev, Duke, 2016
3 //
4 // This module finds periods of time-localized activity
5 // on each channel of the optical system and creates OpHits.
6 //
7 // Split from OpFlashFinder_module.cc
8 // by Ben Jones, MIT, 2013
9 //
10 
11 // LArSoft includes
18 #include "larcore/CoreUtils/ServiceUtil.h" // lar::providerFrom()
27 
28 // Framework includes
35 #include "fhiclcpp/ParameterSet.h"
37 
38 // ROOT includes
39 
40 // C++ Includes
41 #include <map>
42 #include <memory>
43 #include <string>
44 
45 namespace {
46 
80  fhicl::ParameterSet makeAlgoToolConfig(fhicl::ParameterSet const& baseConfig,
81  std::string const& configKey,
82  std::string const& algoClassPrefix = "",
83  std::string const& algoNameKey = "Name")
84  {
85 
86  fhicl::ParameterSet toolConfig;
87 
88  // add algo configuration
89  fhicl::ParameterSet algoConfig = baseConfig.get<fhicl::ParameterSet>(configKey);
90 
91  if (auto toolType = algoConfig.get_if_present<std::string>("tool_type")) {
92  toolConfig.put("tool_type", *toolType);
93  algoConfig.erase("tool_type");
94  }
95  else if (auto algoName = algoConfig.get_if_present<std::string>(algoNameKey)) {
96  toolConfig.put("tool_type", algoClassPrefix + *algoName + "Maker");
97  // we leave the algorithm Name in its configuration for compatibility
98  }
99  toolConfig.put(configKey, std::move(algoConfig));
100 
101  return toolConfig;
102  } // makeAlgoToolConfig()
103 
114  fhicl::ParameterSet makePedAlgoToolConfig(fhicl::ParameterSet const& baseConfig)
115  {
116  return makeAlgoToolConfig(baseConfig, "PedAlgoPset", "PedAlgo");
117  }
118 
140  fhicl::ParameterSet makeHitAlgoToolConfig(fhicl::ParameterSet const& baseConfig)
141  {
142  fhicl::ParameterSet toolConfig = makeAlgoToolConfig(baseConfig, "HitAlgoPset", "Algo");
143 
144  // add rise time calculator configuration ("no configuration" is ok)
145  if (auto riseCalcCfg = baseConfig.get_if_present<fhicl::ParameterSet>("RiseTimeCalculator")) {
146  toolConfig.put("RiseTimeCalculator", std::move(*riseCalcCfg));
147  }
148 
149  return toolConfig;
150  } // makeHitAlgoToolConfig()
151 
152 }
153 
154 namespace opdet {
155 
156  class OpHitFinder : public art::EDProducer {
157  public:
158  // Standard constructor and destructor for an ART module.
159  explicit OpHitFinder(const fhicl::ParameterSet&);
160 
161  // The producer routine, called once per event.
162  void produce(art::Event&);
163 
164  private:
165  std::map<int, int> GetChannelMap();
166  std::vector<double> GetSPEScales();
167  std::vector<double> GetSPEShifts();
168 
169  // The parameters we'll read from the .fcl file.
170  std::string fInputModule; // Input tag for OpDetWaveform collection
171  std::string fGenModule;
172  std::vector<std::string> fInputLabels;
173  std::set<unsigned int> fChannelMasks;
174 
176  std::unique_ptr<pmtana::PMTPulseRecoBase> const fThreshAlg;
177  std::unique_ptr<pmtana::PMTPedestalBase> const fPedAlg;
178 
179  Float_t fHitThreshold;
180  unsigned int fMaxOpChannel;
182 
184  };
185 
186 }
187 
188 namespace opdet {
190 }
191 
192 namespace opdet {
193 
194  //----------------------------------------------------------------------------
195  // Constructor
197  : EDProducer{pset}
198  , fPulseRecoMgr()
199  , fThreshAlg{art::make_tool<opdet::IHitAlgoMakerTool>(makeHitAlgoToolConfig(pset))->makeAlgo()}
200  , fPedAlg{art::make_tool<opdet::IPedAlgoMakerTool>(makePedAlgoToolConfig(pset))->makeAlgo()}
201  {
202  // Indicate that the Input Module comes from .fcl
203  fInputModule = pset.get<std::string>("InputModule");
204  fGenModule = pset.get<std::string>("GenModule");
205  fInputLabels = pset.get<std::vector<std::string>>("InputLabels");
206  fUseStartTime = pset.get<bool>("UseStartTime", false);
207 
208  for (auto const& ch :
209  pset.get<std::vector<unsigned int>>("ChannelMasks", std::vector<unsigned int>()))
210  fChannelMasks.insert(ch);
211 
212  fHitThreshold = pset.get<float>("HitThreshold");
213  bool useCalibrator = pset.get<bool>("UseCalibrator", false);
214 
215  auto const& geometry(*lar::providerFrom<geo::Geometry>());
216  fMaxOpChannel = geometry.MaxOpChannel();
217 
218  if (useCalibrator) {
219  // If useCalibrator, get it from ART
220  fCalib = lar::providerFrom<calib::IPhotonCalibratorService>();
221  }
222  else {
223  // If not useCalibrator, make an internal one based
224  // on fhicl settings to hit finder.
225  bool areaToPE = pset.get<bool>("AreaToPE");
226  float SPEArea = pset.get<float>("SPEArea");
227  float SPEShift = pset.get<float>("SPEShift", 0.);
228 
229  // Reproduce behavior from GetSPEScales()
230  if (!areaToPE) SPEArea = 20;
231 
232  // Delete and replace if we are reconfiguring
233  if (fCalib) { delete fCalib; }
234 
235  fCalib = new calib::PhotonCalibratorStandard(SPEArea, SPEShift, areaToPE);
236  }
237 
238  produces<std::vector<recob::OpHit>>();
239 
242 
243  // show the algorithm selection on screen
244  mf::LogInfo{"OpHitFinder"} << "Pulse finder algorithm: '" << fThreshAlg->Name() << "'"
245  << "\nPedestal algorithm: '" << fPedAlg->Name() << "'";
246  }
247 
248  //----------------------------------------------------------------------------
250  {
251 
252  // These is the storage pointer we will put in the event
253  std::unique_ptr<std::vector<recob::OpHit>> HitPtr(new std::vector<recob::OpHit>);
254 
255  std::vector<const sim::BeamGateInfo*> beamGateArray;
256  try {
257  evt.getView(fGenModule, beamGateArray);
258  }
259  catch (art::Exception const& err) {
260  if (err.categoryCode() != art::errors::ProductNotFound) throw;
261  }
262 
263  auto const& geometry(*lar::providerFrom<geo::Geometry>());
264  auto const clock_data =
266  auto const& calibrator(*fCalib);
267  //
268  // Get the pulses from the event
269  //
270 
271  // Load pulses into WaveformVector
272  if (fChannelMasks.empty() && fInputLabels.size() < 2) {
274  if (fInputLabels.empty())
275  evt.getByLabel(fInputModule, wfHandle);
276  else
277  evt.getByLabel(fInputModule, fInputLabels.front(), wfHandle);
278  assert(wfHandle.isValid());
279  RunHitFinder(*wfHandle,
280  *HitPtr,
282  *fThreshAlg,
283  geometry,
285  clock_data,
286  calibrator,
287  fUseStartTime);
288  }
289  else {
290 
291  // Reserve a large enough array
292  int totalsize = 0;
293  for (auto label : fInputLabels) {
295  evt.getByLabel(fInputModule, label, wfHandle);
296  if (!wfHandle.isValid()) continue; // Skip non-existent collections
297  totalsize += wfHandle->size();
298  }
299 
300  std::vector<raw::OpDetWaveform> WaveformVector;
301  WaveformVector.reserve(totalsize);
302 
303  for (auto label : fInputLabels) {
305  evt.getByLabel(fInputModule, label, wfHandle);
306  if (!wfHandle.isValid()) continue; // Skip non-existent collections
307 
308  //WaveformVector.insert(WaveformVector.end(),
309  // wfHandle->begin(), wfHandle->end());
310  for (auto const& wf : *wfHandle) {
311  if (fChannelMasks.find(wf.ChannelNumber()) != fChannelMasks.end()) continue;
312  WaveformVector.push_back(wf);
313  }
314  }
315 
316  RunHitFinder(WaveformVector,
317  *HitPtr,
319  *fThreshAlg,
320  geometry,
322  clock_data,
323  calibrator,
324  fUseStartTime);
325  }
326  // Store results into the event
327  evt.put(std::move(HitPtr));
328  }
329 
330 } // namespace opdet
Utilities related to art service access.
Tool interface for creating a pedestal estimator algorithm.
std::unique_ptr< pmtana::PMTPedestalBase > const fPedAlg
std::map< int, int > GetChannelMap()
EDProducer(fhicl::ParameterSet const &pset)
Definition: EDProducer.cc:6
void AddRecoAlgo(pmtana::PMTPulseRecoBase *algo, PMTPedestalBase *ped_algo=nullptr)
A method to set pulse reconstruction algorithm.
OpHitFinder(const fhicl::ParameterSet &)
Class definition file of PMTPedestalBase.
std::set< unsigned int > fChannelMasks
bool isValid() const noexcept
Definition: Handle.h:203
PutHandle< PROD > put(std::unique_ptr< PROD > &&edp, std::string const &instance={})
Definition: Event.h:77
#define DEFINE_ART_MODULE(klass)
Definition: ModuleMacros.h:65
T get(std::string const &key) const
Definition: ParameterSet.h:314
std::vector< double > GetSPEShifts()
std::size_t getView(std::string const &moduleLabel, std::string const &productInstanceName, std::string const &processName, std::vector< ELEMENT const * > &result) const
void produce(art::Event &)
cet::coded_exception< errors::ErrorCodes, ExceptionDetail::translate > Exception
Definition: Exception.h:66
calib::IPhotonCalibrator const * fCalib
bool getByLabel(std::string const &label, std::string const &instance, Handle< PROD > &result) const
Tool interface for creating a hit finder algorithm.
void RunHitFinder(std::vector< raw::OpDetWaveform > const &opDetWaveformVector, std::vector< recob::OpHit > &hitVector, pmtana::PulseRecoManager const &pulseRecoMgr, pmtana::PMTPulseRecoBase const &threshAlg, geo::GeometryCore const &geometry, float hitThreshold, detinfo::DetectorClocksData const &clocksData, calib::IPhotonCalibrator const &calibrator, bool use_start_time)
Definition: OpHitAlg.cxx:29
pmtana::PulseRecoManager fPulseRecoMgr
void SetDefaultPedAlgo(pmtana::PMTPedestalBase *algo)
A method to set a choice of pedestal estimation method.
std::optional< T > get_if_present(std::string const &key) const
Definition: ParameterSet.h:267
std::vector< double > GetSPEScales()
bool erase(std::string const &key)
Class definition file of PMTPulseRecoBase.
TCEvent evt
Definition: DataStructs.cxx:8
std::vector< std::string > fInputLabels
void put(std::string const &key)
art framework interface to geometry description
std::unique_ptr< pmtana::PMTPulseRecoBase > const fThreshAlg
Class definition file of PulseRecoManager.