LArSoft  v10_04_05
Liquid Argon Software toolkit - https://larsoft.org/
PDFastSimPAR_module.cc
Go to the documentation of this file.
1 // Class: PDFastSimPAR
3 // Plugin Type: producer
4 // File: PDFastSimPAR_module.cc
5 // Description:
6 // - acts on sim::SimEnergyDeposit from LArG4Main,
7 // - simulate (fast, photon visibility service) the OpDet response to optical
8 // photons Input: 'sim::SimEnergyDeposit' Output: 'sim::OpDetBacktrackerRecord'
9 // Fast simulation of propagating the photons created from SimEnergyDeposits.
10 
11 // This module does a fast simulation of propagating the photons created from
12 // SimEnergyDeposits, This simulation is done using the Semi-Analytical model, which
13 // stores the visibilities of each optical channel with respect to each optical
14 // voxel in the TPC volume, to avoid propagating single photons using Geant4. At
15 // the end of this module a collection of the propagated photons either as
16 //'sim::OpDetBacktrackerRecord' are placed into the art event.
17 
18 // The steps this module takes are:
19 // - to take number of photon and the vertex information from
20 // 'sim::SimEnergyDeposits',
21 // - use the Semi-Analytical model (visibilities) to determine the amount of visible
22 // photons at each optical channel,
23 // - visible photons: the number of photons times the visibility at the middle
24 // of the Geant4 step for a given optical channel.
25 // - other photon information is got from 'sim::SimEnergyDeposits'
26 // - add 'sim::OpDetBacktrackerRecord' to event
27 // Aug. 19 by Mu Wei
28 // Restructured Nov. 21 by P. Green
30 
31 // LArSoft libraries
44 
46 
47 // Art libraries
56 #include "cetlib_except/exception.h"
57 #include "fhiclcpp/ParameterSet.h"
58 #include "fhiclcpp/types/Atom.h"
59 #include "fhiclcpp/types/Comment.h"
61 #include "fhiclcpp/types/Name.h"
64 
65 // Random numbers
66 #include "CLHEP/Random/RandPoissonQ.h"
67 
68 #include "range/v3/view/enumerate.hpp"
69 
70 #include <cmath>
71 #include <ctime>
72 #include <map>
73 #include <memory>
74 #include <vector>
75 
76 namespace phot {
77  class PDFastSimPAR : public art::EDProducer {
78  public:
79  // Define the fhicl configuration
80  struct Config {
81  using Name = fhicl::Name;
85 
87  Comment("SimEnergyDeposit label.")};
89  Comment("Simulate slow scintillation light, default true"),
90  true};
92  Comment("Simulate slow scintillation light")};
94  Comment("Simulate reflected visible light")};
96  Name("IncludeAnodeReflections"),
97  Comment("Simulate anode reflections, default false"),
98  false};
100  Comment("Simulate light propagation time")};
102  Name("GeoPropTimeOnly"),
103  Comment("Simulate light propagation time geometric approximation, default false"),
104  false};
106  Name("UseLitePhotons"),
107  Comment("Store SimPhotonsLite/OpDetBTRs instead of SimPhotons")};
109  Comment("Photons cannot cross the cathode")};
111  Name("OnlyActiveVolume"),
112  Comment("PAR fast sim usually only for active volume, default true"),
113  true};
115  Comment("Set to true if light is only supported in C:1")};
116  DP ScintTimeTool{Name("ScintTimeTool"),
117  Comment("Tool describing scintillation time structure")};
119  Name("UseXeAbsorption"),
120  Comment("Use Xe absorption length instead of Ar, default false"),
121  false};
122  ODP VUVTiming{Name("VUVTiming"), Comment("Configuration for UV timing parameterization")};
123  ODP VISTiming{Name("VISTiming"),
124  Comment("Configuration for visible timing parameterization")};
125  DP VUVHits{Name("VUVHits"), Comment("Configuration for UV visibility parameterization")};
126  ODP VISHits{Name("VISHits"),
127  Comment("Configuration for visibile visibility parameterization")};
128  };
130 
131  explicit PDFastSimPAR(Parameters const& config);
132  void produce(art::Event&) override;
133 
134  private:
135  void detectedNumPhotons(std::vector<int>& DetectedNumPhotons,
136  const std::vector<double>& OpDetVisibilities,
137  const int NumPhotons) const;
138 
139  void AddOpDetBTR(std::vector<sim::OpDetBacktrackerRecord>& opbtr,
140  std::vector<int>& ChannelMap,
141  const sim::OpDetBacktrackerRecord& btr) const;
142 
143  bool isOpDetInSameTPC(geo::Point_t const& ScintPoint, geo::Point_t const& OpDetPoint) const;
144  std::vector<geo::Point_t> opDetCenters() const;
145 
146  // semi-analytical model
147  std::unique_ptr<SemiAnalyticalModel> fVisibilityModel;
148 
149  // propagation time model
150  std::unique_ptr<PropagationTimeModel> fPropTimeModel;
151 
152  // random numbers
153  CLHEP::HepRandomEngine& fPhotonEngine;
154  std::unique_ptr<CLHEP::RandPoissonQ> fRandPoissPhot;
155  CLHEP::HepRandomEngine& fScintTimeEngine;
156  std::unique_ptr<ScintTime> fScintTime; // Tool to retrieve timinig of scintillation
157 
158  // geometry properties
161  const size_t fNOpChannels;
162  const std::vector<geo::BoxBoundedGeo> fActiveVolumes;
163  const int fNTPC;
165  const std::vector<geo::Point_t> fOpDetCenter;
166 
167  // Module behavior
169  const bool fDoFastComponent;
170  const bool fDoSlowComponent;
171  const bool fDoReflectedLight;
173  const bool fIncludePropTime;
174  const bool fGeoPropTimeOnly;
175  const bool fUseLitePhotons;
176  const bool fOpaqueCathode;
177  const bool fOnlyActiveVolume;
178  const bool fOnlyOneCryostat;
179  const bool fUseXeAbsorption;
180  };
181 
182  //......................................................................
184  : art::EDProducer{config}
185  , fPhotonEngine(art::ServiceHandle<rndm::NuRandomService>()->registerAndSeedEngine(
186  createEngine(0, "HepJamesRandom", "photon"),
187  "HepJamesRandom",
188  "photon",
189  config.get_PSet(),
190  "SeedPhoton"))
191  , fRandPoissPhot(std::make_unique<CLHEP::RandPoissonQ>(fPhotonEngine))
193  createEngine(0, "HepJamesRandom", "scinttime"),
194  "HepJamesRandom",
195  "scinttime",
196  config.get_PSet(),
197  "SeedScintTime"))
198  , fScintTime{art::make_tool<phot::ScintTime>(config().ScintTimeTool.get<fhicl::ParameterSet>())}
199  , fGeom(*(lar::providerFrom<geo::Geometry>()))
200  , fISTPC{fGeom}
203  , fNTPC(fGeom.NTPC())
205  , fSimTag(config().SimulationLabel())
206  , fDoFastComponent(config().DoFastComponent())
207  , fDoSlowComponent(config().DoSlowComponent())
208  , fDoReflectedLight(config().DoReflectedLight())
209  , fIncludeAnodeReflections(config().IncludeAnodeReflections())
210  , fIncludePropTime(config().IncludePropTime())
211  , fGeoPropTimeOnly(config().GeoPropTimeOnly())
212  , fUseLitePhotons(config().UseLitePhotons())
213  , fOpaqueCathode(config().OpaqueCathode())
214  , fOnlyActiveVolume(config().OnlyActiveVolume())
215  , fOnlyOneCryostat(config().OnlyOneCryostat())
216  , fUseXeAbsorption(config().UseXeAbsorption())
217  {
218  mf::LogInfo("PDFastSimPAR") << "Initializing PDFastSimPAR." << std::endl;
219 
220  // Parameterized Simulation
221  fhicl::ParameterSet VUVHitsParams = config().VUVHits.get<fhicl::ParameterSet>();
222  fhicl::ParameterSet VUVTimingParams;
223  fhicl::ParameterSet VISHitsParams;
224  fhicl::ParameterSet VISTimingParams;
225 
226  // Validate configuration options
227  if (fIncludePropTime &&
228  !config().VUVTiming.get_if_present<fhicl::ParameterSet>(VUVTimingParams)) {
230  << "Propagation time simulation requested, but VUVTiming not specified."
231  << "\n";
232  }
234  !config().VISHits.get_if_present<fhicl::ParameterSet>(VISHitsParams)) {
236  << "Reflected light or anode reflections simulation requested, but VisHits not specified."
237  << "\n";
238  }
240  !config().VISTiming.get_if_present<fhicl::ParameterSet>(VISTimingParams)) {
242  << "Reflected light propagation time simulation requested, but VISTiming not specified."
243  << "\n";
244  }
245  if (fGeom.Ncryostats() > 1U) {
246  if (fOnlyOneCryostat) {
247  mf::LogWarning("PDFastSimPAR")
248  << std::string(80, '=') << "\nA detector with " << fGeom.Ncryostats()
249  << " cryostats is configured"
250  << " , and semi-analytic model is requested for scintillation photon propagation."
251  << " THIS CONFIGURATION IS NOT SUPPORTED and it is open to bugs"
252  << " (e.g. scintillation may be detected only in cryostat #0)."
253  << "\nThis would be normally a fatal error, but it has been forcibly overridden."
254  << "\n"
255  << std::string(80, '=');
256  }
257  else {
259  << "Photon propagation via semi-analytic model is not supported yet"
260  << " on detectors with more than one cryostat.";
261  }
262  }
263 
264  // Initialise the Scintillation Time
265  fScintTime->initRand(fScintTimeEngine);
266 
267  // photo-detector visibility model (semi-analytical model)
268  fVisibilityModel = std::make_unique<SemiAnalyticalModel>(
269  VUVHitsParams, VISHitsParams, fDoReflectedLight, fIncludeAnodeReflections, fUseXeAbsorption);
270 
271  // propagation time model
272  if (fIncludePropTime)
273  fPropTimeModel = std::make_unique<PropagationTimeModel>(
274  VUVTimingParams, VISTimingParams, fScintTimeEngine, fDoReflectedLight, fGeoPropTimeOnly);
275 
276  {
277  auto log = mf::LogTrace("PDFastSimPAR") << "PDFastSimPAR: active volume boundaries from "
278  << fActiveVolumes.size() << " volumes:";
279  for (auto const& [iCryo, box] : ::ranges::views::enumerate(fActiveVolumes)) {
280  log << "\n - C:" << iCryo << ": " << box.Min() << " -- " << box.Max() << " cm";
281  }
282  }
283 
284  if (fUseLitePhotons) {
285  mf::LogInfo("PDFastSimPAR") << "Using Lite Photons";
286  produces<std::vector<sim::SimPhotonsLite>>();
287  produces<std::vector<sim::OpDetBacktrackerRecord>>();
288  if (fDoReflectedLight) {
289  mf::LogInfo("PDFastSimPAR") << "Storing Reflected Photons";
290  produces<std::vector<sim::SimPhotonsLite>>("Reflected");
291  produces<std::vector<sim::OpDetBacktrackerRecord>>("Reflected");
292  }
293  }
294  else {
295  mf::LogInfo("PDFastSimPAR") << "Using Sim Photons";
296  produces<std::vector<sim::SimPhotons>>();
297  if (fDoReflectedLight) {
298  mf::LogInfo("PDFastSimPAR") << "Storing Reflected Photons";
299  produces<std::vector<sim::SimPhotons>>("Reflected");
300  }
301  }
302 
303  // determine drift distance
305  // for multiple TPCs, use second TPC to skip small volume at edges of detector (DUNE)
306  if (fNTPC > 1 && fDriftDistance < 50)
307  fDriftDistance = fGeom.TPC(geo::TPCID{0, 1}).DriftDistance();
308 
309  mf::LogInfo("PDFastSimPAR") << "PDFastSimPAR Initialization finish.\n"
310  << "Simulate using semi-analytic model for number of hits."
311  << std::endl;
312  }
313 
314  //......................................................................
316  {
317  mf::LogTrace("PDFastSimPAR") << "PDFastSimPAR Module Producer"
318  << "EventID: " << event.event();
319 
320  std::vector<int> PDChannelToSOCMapDirect(fNOpChannels, -1); // Where each OpChan is.
321  std::vector<int> PDChannelToSOCMapReflect(fNOpChannels, -1); // Where each OpChan is.
322 
323  // SimPhotonsLite
324  auto phlit = std::make_unique<std::vector<sim::SimPhotonsLite>>();
325  auto opbtr = std::make_unique<std::vector<sim::OpDetBacktrackerRecord>>();
326  auto phlit_ref = std::make_unique<std::vector<sim::SimPhotonsLite>>();
327  auto opbtr_ref = std::make_unique<std::vector<sim::OpDetBacktrackerRecord>>();
328  auto& dir_phlitcol(*phlit);
329  auto& ref_phlitcol(*phlit_ref);
330  // SimPhotons
331  auto phot = std::make_unique<std::vector<sim::SimPhotons>>();
332  auto phot_ref = std::make_unique<std::vector<sim::SimPhotons>>();
333  auto& dir_photcol(*phot);
334  auto& ref_photcol(*phot_ref);
335  if (fUseLitePhotons) {
336  dir_phlitcol.resize(fNOpChannels);
337  ref_phlitcol.resize(fNOpChannels);
338  for (unsigned int i = 0; i < fNOpChannels; i++) {
339  dir_phlitcol[i].OpChannel = i;
340  ref_phlitcol[i].OpChannel = i;
341  }
342  }
343  else { // SimPhotons
344  dir_photcol.resize(fNOpChannels);
345  ref_photcol.resize(fNOpChannels);
346  for (unsigned int i = 0; i < fNOpChannels; i++) {
347  dir_photcol[i].fOpChannel = i;
348  ref_photcol[i].fOpChannel = i;
349  }
350  }
351 
353  if (!event.getByLabel(fSimTag, edepHandle)) {
354  mf::LogError("PDFastSimPAR") << "PDFastSimPAR Module Cannot getByLabel: " << fSimTag;
355  return;
356  }
357 
358  auto const& edeps = edepHandle;
359 
360  int num_points = 0;
361  int num_fastph = 0;
362  int num_slowph = 0;
363  int num_fastdp = 0;
364  int num_slowdp = 0;
365 
366  for (auto const& edepi : *edeps) {
367  num_points++;
368 
369  int nphot_fast = edepi.NumFPhotons();
370  int nphot_slow = edepi.NumSPhotons();
371 
372  num_fastph += nphot_fast;
373  num_slowph += nphot_slow;
374 
375  if (!((nphot_fast > 0 && fDoFastComponent) || (nphot_slow > 0 && fDoSlowComponent))) continue;
376 
377  int trackID = edepi.TrackID();
378  int nphot = edepi.NumPhotons();
379  double edeposit = edepi.Energy() / nphot;
380  double pos[3] = {edepi.MidPointX(), edepi.MidPointY(), edepi.MidPointZ()};
381  geo::Point_t const ScintPoint = {pos[0], pos[1], pos[2]};
382 
383  if (fOnlyActiveVolume && !fISTPC.isScintInActiveVolume(ScintPoint)) continue;
384 
385  // direct light
386  std::vector<int> DetectedNumFast(fNOpChannels);
387  std::vector<int> DetectedNumSlow(fNOpChannels);
388 
389  std::vector<double> OpDetVisibilities;
390  fVisibilityModel->detectedDirectVisibilities(OpDetVisibilities, ScintPoint);
391  detectedNumPhotons(DetectedNumFast, OpDetVisibilities, nphot_fast);
392  detectedNumPhotons(DetectedNumSlow, OpDetVisibilities, nphot_slow);
393 
395  std::vector<int> AnodeDetectedNumFast(fNOpChannels);
396  std::vector<int> AnodeDetectedNumSlow(fNOpChannels);
397 
398  std::vector<double> OpDetVisibilitiesAnode;
399  fVisibilityModel->detectedReflectedVisibilities(OpDetVisibilitiesAnode, ScintPoint, true);
400  detectedNumPhotons(AnodeDetectedNumFast, OpDetVisibilitiesAnode, nphot_fast);
401  detectedNumPhotons(AnodeDetectedNumSlow, OpDetVisibilitiesAnode, nphot_slow);
402 
403  // add to existing count
404  for (size_t i = 0; i < AnodeDetectedNumFast.size(); ++i) {
405  DetectedNumFast[i] += AnodeDetectedNumFast[i];
406  }
407  for (size_t i = 0; i < AnodeDetectedNumSlow.size(); ++i) {
408  DetectedNumSlow[i] += AnodeDetectedNumSlow[i];
409  }
410  }
411 
412  // reflected light, if enabled
413  std::vector<int> ReflDetectedNumFast(fNOpChannels);
414  std::vector<int> ReflDetectedNumSlow(fNOpChannels);
415  if (fDoReflectedLight) {
416  std::vector<double> OpDetVisibilitiesRefl;
417  fVisibilityModel->detectedReflectedVisibilities(OpDetVisibilitiesRefl, ScintPoint, false);
418  detectedNumPhotons(ReflDetectedNumFast, OpDetVisibilitiesRefl, nphot_fast);
419  detectedNumPhotons(ReflDetectedNumSlow, OpDetVisibilitiesRefl, nphot_slow);
420  }
421 
422  // loop through direct photons then reflected photons cases
423  size_t DoReflected = (fDoReflectedLight) ? 1 : 0;
424  for (size_t Reflected = 0; Reflected <= DoReflected; ++Reflected) {
425  for (size_t channel = 0; channel < fNOpChannels; channel++) {
426 
427  if (fOpaqueCathode && !isOpDetInSameTPC(ScintPoint, fOpDetCenter[channel])) continue;
428 
429  int ndetected_fast = DetectedNumFast[channel];
430  int ndetected_slow = DetectedNumSlow[channel];
431  if (Reflected) {
432  ndetected_fast = ReflDetectedNumFast[channel];
433  ndetected_slow = ReflDetectedNumSlow[channel];
434  }
435  if (!((ndetected_fast > 0 && fDoFastComponent) ||
436  (ndetected_slow > 0 && fDoSlowComponent)))
437  continue;
438 
439  // calculate propagation time, does not matter whether fast or slow photon
440  std::vector<double> transport_time;
441  if (fIncludePropTime) {
442  transport_time.resize(ndetected_fast + ndetected_slow);
443  fPropTimeModel->propagationTime(transport_time, ScintPoint, channel, Reflected);
444  }
445 
446  // SimPhotonsLite case
447  if (fUseLitePhotons) {
448  sim::OpDetBacktrackerRecord tmpbtr(channel);
449  if (ndetected_fast > 0 && fDoFastComponent) {
450  int n = ndetected_fast;
451  num_fastdp += n;
452  for (int i = 0; i < n; ++i) {
453  // calculates the time at which the photon was produced
454  double dtime = edepi.StartT() + fScintTime->fastScintTime();
455  if (fIncludePropTime) dtime += transport_time[i];
456  int time = static_cast<int>(std::round(dtime));
457  if (Reflected)
458  ++ref_phlitcol[channel].DetectedPhotons[time];
459  else
460  ++dir_phlitcol[channel].DetectedPhotons[time];
461  tmpbtr.AddScintillationPhotons(trackID, time, 1, pos, edeposit);
462  }
463  }
464  if (ndetected_slow > 0 && fDoSlowComponent) {
465  int n = ndetected_slow;
466  num_slowdp += n;
467  for (int i = 0; i < n; ++i) {
468  double dtime = edepi.StartT() + fScintTime->slowScintTime();
469  if (fIncludePropTime) dtime += transport_time[ndetected_fast + i];
470  int time = static_cast<int>(std::round(dtime));
471  if (Reflected)
472  ++ref_phlitcol[channel].DetectedPhotons[time];
473  else
474  ++dir_phlitcol[channel].DetectedPhotons[time];
475  tmpbtr.AddScintillationPhotons(trackID, time, 1, pos, edeposit);
476  }
477  }
478  if (Reflected)
479  AddOpDetBTR(*opbtr_ref, PDChannelToSOCMapReflect, tmpbtr);
480  else
481  AddOpDetBTR(*opbtr, PDChannelToSOCMapDirect, tmpbtr);
482  }
483  // SimPhotons case
484  else {
485  sim::OnePhoton photon;
486  photon.SetInSD = false;
487  photon.InitialPosition = edepi.End();
488  photon.MotherTrackID = edepi.TrackID();
489  if (Reflected)
490  photon.Energy = 2.9 * CLHEP::eV; // 430 nm
491  else
492  photon.Energy = 9.7 * CLHEP::eV; // 128 nm
493  // TODO: un-hardcode and add another energy for Xe scintillation
494  if (ndetected_fast > 0 && fDoFastComponent) {
495  int n = ndetected_fast;
496  num_fastdp += n;
497  for (int i = 0; i < n; ++i) {
498  // calculates the time at which the photon was produced
499  double dtime = edepi.StartT() + fScintTime->fastScintTime();
500  if (fIncludePropTime) dtime += transport_time[i];
501  int time = static_cast<int>(std::round(dtime));
502  photon.Time = time;
503  if (Reflected)
504  ref_photcol[channel].insert(ref_photcol[channel].end(), 1, photon);
505  else
506  dir_photcol[channel].insert(dir_photcol[channel].end(), 1, photon);
507  }
508  }
509  if (ndetected_slow > 0 && fDoSlowComponent) {
510  int n = ndetected_slow;
511  num_slowdp += n;
512  for (int i = 0; i < n; ++i) {
513  double dtime = edepi.StartT() + fScintTime->slowScintTime();
514  if (fIncludePropTime) dtime += transport_time[ndetected_fast + i];
515  int time = static_cast<int>(std::round(dtime));
516  photon.Time = time;
517  if (Reflected)
518  ref_photcol[channel].insert(ref_photcol[channel].end(), 1, photon);
519  else
520  dir_photcol[channel].insert(dir_photcol[channel].end(), 1, photon);
521  }
522  }
523  }
524  }
525  }
526  }
527 
528  mf::LogTrace("PDFastSimPAR") << "Total points: " << num_points
529  << ", total fast photons: " << num_fastph
530  << ", total slow photons: " << num_slowph
531  << "\ndetected fast photons: " << num_fastdp
532  << ", detected slow photons: " << num_slowdp;
533 
534  if (fUseLitePhotons) {
535  event.put(move(phlit));
536  event.put(move(opbtr));
537  if (fDoReflectedLight) {
538  event.put(move(phlit_ref), "Reflected");
539  event.put(move(opbtr_ref), "Reflected");
540  }
541  }
542  else {
543  event.put(move(phot));
544  if (fDoReflectedLight) { event.put(move(phot_ref), "Reflected"); }
545  }
546 
547  return;
548  }
549 
550  //......................................................................
551  void PDFastSimPAR::AddOpDetBTR(std::vector<sim::OpDetBacktrackerRecord>& opbtr,
552  std::vector<int>& ChannelMap,
553  const sim::OpDetBacktrackerRecord& btr) const
554  {
555  int iChan = btr.OpDetNum();
556  if (ChannelMap[iChan] < 0) {
557  ChannelMap[iChan] = opbtr.size();
558  opbtr.emplace_back(std::move(btr));
559  }
560  else {
561  size_t idtest = ChannelMap[iChan];
562  auto const& timePDclockSDPsMap = btr.timePDclockSDPsMap();
563  for (auto const& timePDclockSDP : timePDclockSDPsMap) {
564  for (auto const& sdp : timePDclockSDP.second) {
565  double xyz[3] = {sdp.x, sdp.y, sdp.z};
566  opbtr.at(idtest).AddScintillationPhotons(
567  sdp.trackID, timePDclockSDP.first, sdp.numPhotons, xyz, sdp.energy);
568  }
569  }
570  }
571  }
572 
573  //......................................................................
574  // calculates number of photons detected given visibility and emitted number of photons
575  void PDFastSimPAR::detectedNumPhotons(std::vector<int>& DetectedNumPhotons,
576  const std::vector<double>& OpDetVisibilities,
577  const int NumPhotons) const
578  {
579  for (size_t i = 0; i < OpDetVisibilities.size(); ++i) {
580  DetectedNumPhotons[i] = fRandPoissPhot->fire(OpDetVisibilities[i] * NumPhotons);
581  }
582  }
583 
584  //......................................................................
585  // checks whether photo-detector is able to see the emitted light scintillation
587  geo::Point_t const& OpDetPoint) const
588  {
589  // check optical channel is in same TPC as scintillation light, if not doesn't see light
590  // temporary method, needs to be replaced with geometry service
591  // working for SBND, uBooNE, DUNE HD 1x2x6, DUNE HD 10kt and DUNE VD subset
592 
593  // special case for SBND = 2 TPCs
594  // check x coordinate has same sign or is close to zero
595  if (fNTPC == 2 && ((ScintPoint.X() < 0.) != (OpDetPoint.X() < 0.)) &&
596  std::abs(OpDetPoint.X()) > 10.) {
597  return false;
598  }
599 
600  // special case for DUNE-HD 10kt = 300 TPCs
601  // check whether distance in drift direction > 1 drift distance
602  if (fNTPC == 300 && std::abs(ScintPoint.X() - OpDetPoint.X()) > fDriftDistance) {
603  return false;
604  }
605  // not needed for DUNE HD 1x2x6, DUNE VD subset, uBooNE
606 
607  return true;
608  }
609 
610  std::vector<geo::Point_t> PDFastSimPAR::opDetCenters() const
611  {
612  std::vector<geo::Point_t> opDetCenter;
613  for (size_t const i : ::ranges::views::ints(size_t(0), fNOpChannels)) {
614  geo::OpDetGeo const& opDet = fGeom.OpDetGeoFromOpDet(i);
615  opDetCenter.push_back(opDet.GetCenter());
616  }
617  return opDetCenter;
618  }
619  // ---------------------------------------------------------------------------
620 
621 } // namespace phot
622 
const std::vector< geo::BoxBoundedGeo > fActiveVolumes
base_engine_t & createEngine(seed_t seed)
bool isOpDetInSameTPC(geo::Point_t const &ScintPoint, geo::Point_t const &OpDetPoint) const
void AddOpDetBTR(std::vector< sim::OpDetBacktrackerRecord > &opbtr, std::vector< int > &ChannelMap, const sim::OpDetBacktrackerRecord &btr) const
static std::vector< geo::BoxBoundedGeo > extractActiveLArVolume(geo::GeometryCore const &geom)
Definition: ISTPC.cxx:46
Utilities related to art service access.
MaybeLogger_< ELseverityLevel::ELsev_info, false > LogInfo
std::unique_ptr< PropagationTimeModel > fPropTimeModel
unsigned int NTPC(CryostatID const &cryoid=details::cryostat_zero) const
Returns the total number of TPCs in the specified cryostat.
Definition: GeometryCore.h:416
fhicl::Atom< bool > DoSlowComponent
Point_t const & GetCenter() const
Definition: OpDetGeo.h:71
EDProducer(fhicl::ParameterSet const &pset)
Definition: EDProducer.cc:6
PDFastSimPAR(Parameters const &config)
All information of a photon entering the sensitive optical detector volume.
Definition: SimPhotons.h:60
fhicl::Atom< bool > DoReflectedLight
constexpr auto abs(T v)
Returns the absolute value of the argument.
std::unique_ptr< ScintTime > fScintTime
MaybeLogger_< ELseverityLevel::ELsev_error, false > LogError
unsigned int Ncryostats() const
Returns the number of cryostats in the detector.
Definition: GeometryCore.h:303
Energy deposited on a readout Optical Detector by simulated tracks.
geo::Point_t InitialPosition
Scintillation position in world coordinates [cm].
Definition: SimPhotons.h:63
fhicl::Atom< bool > UseLitePhotons
decltype(auto) constexpr end(T &&obj)
ADL-aware version of std::end.
Definition: StdUtils.h:77
auto enumerate(Iterables &&...iterables)
Range-for loop helper tracking the number of iteration.
Definition: enumerate.h:49
std::unique_ptr< CLHEP::RandPoissonQ > fRandPoissPhot
CLHEP::HepRandomEngine & fScintTimeEngine
Simulation objects for optical detectors.
int OpDetNum() const
Returns the readout Optical Detector this object describes.
#define DEFINE_ART_MODULE(klass)
Definition: ModuleMacros.h:65
const std::vector< geo::Point_t > fOpDetCenter
Definitions of geometry vector data types.
std::unique_ptr< SemiAnalyticalModel > fVisibilityModel
fhicl::Atom< bool > OpaqueCathode
T get(std::string const &key) const
Definition: ParameterSet.h:314
fhicl::Atom< bool > IncludeAnodeReflections
void produce(art::Event &) override
The data type to uniquely identify a TPC.
Definition: geo_types.h:306
Description of the physical geometry of one entire detector.
Definition: GeometryCore.h:91
const art::InputTag fSimTag
An art service to assist in the distribution of guaranteed unique seeds to all engines within an art ...
unsigned int NOpDets() const
Number of OpDets in the whole detector.
Provides a base class aware of world box coordinates.
void AddScintillationPhotons(TrackID_t trackID, timePDclock_t timePDclock, double numberPhotons, double const *xyz, double energy)
Add scintillation photons and energy to this OpticalDetector.
Encapsulate the geometry of an optical detector.
cet::coded_exception< errors::ErrorCodes, ExceptionDetail::translate > Exception
Definition: Exception.h:66
ROOT::Math::PositionVector3D< ROOT::Math::Cartesian3D< double >, ROOT::Math::GlobalCoordinateSystemTag > Point_t
Type for representation of position in physical 3D space.
Definition: geo_vectors.h:180
fhicl::Atom< bool > UseXeAbsorption
General LArSoft Utilities.
fhicl::Atom< bool > GeoPropTimeOnly
CLHEP::HepRandomEngine & fPhotonEngine
double DriftDistance() const
Drift distance is defined as the distance between the anode and the cathode, in centimeters.
Definition: TPCGeo.cxx:165
bool getByLabel(std::string const &label, std::string const &instance, Handle< PROD > &result) const
std::vector< geo::Point_t > opDetCenters() const
bool SetInSD
Whether the photon reaches the sensitive detector.
Definition: SimPhotons.h:84
contains information for a single step in the detector simulation
void detectedNumPhotons(std::vector< int > &DetectedNumPhotons, const std::vector< double > &OpDetVisibilities, const int NumPhotons) const
MaybeLogger_< ELseverityLevel::ELsev_warning, false > LogWarning
fhicl::Atom< art::InputTag > SimulationLabel
Definition: MVAAlg.h:12
Char_t n[5]
MaybeLogger_< ELseverityLevel::ELsev_success, true > LogTrace
fhicl::Atom< bool > IncludePropTime
float Energy
Scintillation photon energy [GeV].
Definition: SimPhotons.h:78
fhicl::Atom< bool > OnlyActiveVolume
TPCGeo const & TPC(TPCID const &tpcid=details::tpc_zero) const
Returns the specified TPC.
Definition: GeometryCore.h:448
bool isScintInActiveVolume(geo::Point_t const &ScintPoint)
Definition: ISTPC.cxx:37
const bool fIncludeAnodeReflections
timePDclockSDPs_t const & timePDclockSDPsMap() const
Returns all the deposited energy information as stored.
fhicl::Atom< bool > OnlyOneCryostat
art framework interface to geometry description
geo::GeometryCore const & fGeom
fhicl::Atom< bool > DoFastComponent
OpDetGeo const & OpDetGeoFromOpDet(unsigned int OpDet) const
Returns the geo::OpDetGeo object for the given detector number.
Event finding and building.
int MotherTrackID
ID of the GEANT4 track causing the scintillation.
Definition: SimPhotons.h:81