Classes | Public Types | Public Member Functions | Private Member Functions | Private Attributes | List of all members
phot::PDFastSimPAR Class Reference
Inheritance diagram for phot::PDFastSimPAR:
art::EDProducer art::detail::Producer art::detail::LegacyModule art::Modifier art::ModuleBase art::ProductRegistryHelper

Classes

struct  Config
 

Public Types

using Parameters = art::EDProducer::Table< Config >
 
- Public Types inherited from art::EDProducer
using ModuleType = EDProducer
 
using WorkerType = WorkerT< EDProducer >
 
- Public Types inherited from art::detail::Producer
template<typename UserConfig , typename KeysToIgnore = void>
using Table = Modifier::Table< UserConfig, KeysToIgnore >
 
- Public Types inherited from art::Modifier
template<typename UserConfig , typename UserKeysToIgnore = void>
using Table = ProducerTable< UserConfig, detail::ModuleConfig, UserKeysToIgnore >
 

Public Member Functions

 PDFastSimPAR (Parameters const &config)
 
void produce (art::Event &) override
 
- Public Member Functions inherited from art::EDProducer
 EDProducer (fhicl::ParameterSet const &pset)
 
template<typename Config >
 EDProducer (Table< Config > const &config)
 
std::string workerType () const
 
- Public Member Functions inherited from art::detail::Producer
virtual ~Producer () noexcept
 
 Producer (fhicl::ParameterSet const &)
 
 Producer (Producer const &)=delete
 
 Producer (Producer &&)=delete
 
Produceroperator= (Producer const &)=delete
 
Produceroperator= (Producer &&)=delete
 
void doBeginJob (SharedResources const &resources)
 
void doEndJob ()
 
void doRespondToOpenInputFile (FileBlock const &fb)
 
void doRespondToCloseInputFile (FileBlock const &fb)
 
void doRespondToOpenOutputFiles (FileBlock const &fb)
 
void doRespondToCloseOutputFiles (FileBlock const &fb)
 
bool doBeginRun (RunPrincipal &rp, ModuleContext const &mc)
 
bool doEndRun (RunPrincipal &rp, ModuleContext const &mc)
 
bool doBeginSubRun (SubRunPrincipal &srp, ModuleContext const &mc)
 
bool doEndSubRun (SubRunPrincipal &srp, ModuleContext const &mc)
 
bool doEvent (EventPrincipal &ep, ModuleContext const &mc, std::atomic< std::size_t > &counts_run, std::atomic< std::size_t > &counts_passed, std::atomic< std::size_t > &counts_failed)
 
- Public Member Functions inherited from art::Modifier
 ~Modifier () noexcept
 
 Modifier ()
 
 Modifier (Modifier const &)=delete
 
 Modifier (Modifier &&)=delete
 
Modifieroperator= (Modifier const &)=delete
 
Modifieroperator= (Modifier &&)=delete
 
- Public Member Functions inherited from art::ModuleBase
virtual ~ModuleBase () noexcept
 
 ModuleBase ()
 
ModuleDescription const & moduleDescription () const
 
void setModuleDescription (ModuleDescription const &)
 
std::array< std::vector< ProductInfo >, NumBranchTypes > const & getConsumables () const
 
void sortConsumables (std::string const &current_process_name)
 
template<typename T , BranchType BT>
ViewToken< T > consumesView (InputTag const &tag)
 
template<typename T , BranchType BT>
ViewToken< T > mayConsumeView (InputTag const &tag)
 

Private Member Functions

void Initialization ()
 
void detectedNumPhotons (std::map< size_t, int > &DetectedNumPhotons, const std::map< size_t, double > &OpDetVisibilities, const double NumPhotons)
 
void AddOpDetBTR (std::vector< sim::OpDetBacktrackerRecord > &opbtr, std::map< size_t, int > &ChannelMap, sim::OpDetBacktrackerRecord btr)
 
bool isOpDetInSameTPC (geo::Point_t const &ScintPoint, geo::Point_t const &OpDetPoint) const
 

Private Attributes

larg4::ISTPC fISTPC
 
std::unique_ptr< SemiAnalyticalModelfVisibilityModel
 
std::unique_ptr< PropagationTimeModelfPropTimeModel
 
CLHEP::HepRandomEngine & fPhotonEngine
 
std::unique_ptr< CLHEP::RandPoissonQ > fRandPoissPhot
 
CLHEP::HepRandomEngine & fScintTimeEngine
 
size_t nOpDets
 
std::map< size_t, int > PDChannelToSOCMapDirect
 
std::map< size_t, int > PDChannelToSOCMapReflect
 
std::vector< geo::BoxBoundedGeofActiveVolumes
 
int fNTPC
 
std::vector< geo::Point_tfOpDetCenter
 
art::InputTag simTag
 
bool fDoFastComponent
 
bool fDoSlowComponent
 
bool fDoReflectedLight
 
bool fIncludeAnodeReflections
 
bool fIncludePropTime
 
bool fGeoPropTimeOnly
 
bool fUseLitePhotons
 
bool fOpaqueCathode
 
bool fOnlyActiveVolume
 
bool fOnlyOneCryostat
 
std::unique_ptr< ScintTimefScintTime
 
fhicl::ParameterSet fVUVTimingParams
 
fhicl::ParameterSet fVISTimingParams
 
fhicl::ParameterSet fVUVHitsParams
 
fhicl::ParameterSet fVISHitsParams
 

Additional Inherited Members

- Static Public Member Functions inherited from art::EDProducer
static void commitEvent (EventPrincipal &ep, Event &e)
 
- Protected Member Functions inherited from art::ModuleBase
ConsumesCollectorconsumesCollector ()
 
template<typename T , BranchType = InEvent>
ProductToken< T > consumes (InputTag const &)
 
template<typename Element , BranchType = InEvent>
ViewToken< Element > consumesView (InputTag const &)
 
template<typename T , BranchType = InEvent>
void consumesMany ()
 
template<typename T , BranchType = InEvent>
ProductToken< T > mayConsume (InputTag const &)
 
template<typename Element , BranchType = InEvent>
ViewToken< Element > mayConsumeView (InputTag const &)
 
template<typename T , BranchType = InEvent>
void mayConsumeMany ()
 

Detailed Description

Definition at line 91 of file PDFastSimPAR_module.cc.

Member Typedef Documentation

Definition at line 118 of file PDFastSimPAR_module.cc.

Constructor & Destructor Documentation

phot::PDFastSimPAR::PDFastSimPAR ( Parameters const &  config)
explicit

Definition at line 191 of file PDFastSimPAR_module.cc.

193  , fISTPC{*(lar::providerFrom<geo::Geometry>())}
195  "HepJamesRandom",
196  "photon",
197  config.get_PSet(),
198  "SeedPhoton"))
200  "HepJamesRandom",
201  "scinttime",
202  config.get_PSet(),
203  "SeedScintTime"))
204  , simTag(config().SimulationLabel())
205  , fDoFastComponent(config().DoFastComponent())
206  , fDoSlowComponent(config().DoSlowComponent())
207  , fDoReflectedLight(config().DoReflectedLight())
208  , fIncludeAnodeReflections(config().IncludeAnodeReflections())
209  , fIncludePropTime(config().IncludePropTime())
210  , fGeoPropTimeOnly(config().GeoPropTimeOnly())
211  , fUseLitePhotons(config().UseLitePhotons())
212  , fOpaqueCathode(config().OpaqueCathode())
213  , fOnlyActiveVolume(config().OnlyActiveVolume())
214  , fOnlyOneCryostat(config().OnlyOneCryostat())
215  , fScintTime{art::make_tool<ScintTime>(config().ScintTimeTool.get<fhicl::ParameterSet>())}
216  , fVUVHitsParams(config().VUVHits.get<fhicl::ParameterSet>())
217  {
218 
219  // Validate configuration options
220  if(fIncludePropTime && !config().VUVTiming.get_if_present<fhicl::ParameterSet>(fVUVTimingParams)) {
222  << "Propagation time simulation requested, but VUVTiming not specified." << "\n";
223  }
224 
225  if(fDoReflectedLight && !config().VISHits.get_if_present<fhicl::ParameterSet>(fVISHitsParams)) {
227  << "Reflected light simulation requested, but VisHits not specified." << "\n";
228  }
229 
230  if (fDoReflectedLight && fIncludePropTime && !config().VISTiming.get_if_present<fhicl::ParameterSet>(fVISTimingParams)) {
232  << "Reflected light propagation time simulation requested, but VISTiming not specified." << "\n";
233  }
234 
235  if(fIncludeAnodeReflections && !config().VISHits.get_if_present<fhicl::ParameterSet>(fVISHitsParams)) {
237  << "Anode reflections light simulation requested, but VisHits not specified." << "\n";
238  }
239 
240  Initialization();
241 
242  if (fUseLitePhotons)
243  {
244  mf::LogInfo("PDFastSimPAR") << "Using Lite Photons";
245  produces< std::vector<sim::SimPhotonsLite> >();
246  produces< std::vector<sim::OpDetBacktrackerRecord> >();
247 
249  {
250  mf::LogInfo("PDFastSimPAR") << "Storing Reflected Photons";
251  produces< std::vector<sim::SimPhotonsLite> >("Reflected");
252  produces< std::vector<sim::OpDetBacktrackerRecord> >("Reflected");
253  }
254  }
255  else
256  {
257  mf::LogInfo("PDFastSimPAR") << "Using Sim Photons";
258  produces< std::vector<sim::SimPhotons> >();
260  {
261  mf::LogInfo("PDFastSimPAR") << "Storing Reflected Photons";
262  produces< std::vector<sim::SimPhotons> >("Reflected");
263  }
264  }
265  }
base_engine_t & createEngine(seed_t seed)
fhicl::ParameterSet fVISHitsParams
MaybeLogger_< ELseverityLevel::ELsev_info, false > LogInfo
std::unique_ptr< ScintTime > fScintTime
CLHEP::HepRandomEngine & fScintTimeEngine
fhicl::ParameterSet fVISTimingParams
static Config * config
Definition: config.cpp:1054
fhicl::ParameterSet fVUVTimingParams
cet::coded_exception< errors::ErrorCodes, ExceptionDetail::translate > Exception
Definition: Exception.h:66
CLHEP::HepRandomEngine & fPhotonEngine
fhicl::ParameterSet fVUVHitsParams

Member Function Documentation

void phot::PDFastSimPAR::AddOpDetBTR ( std::vector< sim::OpDetBacktrackerRecord > &  opbtr,
std::map< size_t, int > &  ChannelMap,
sim::OpDetBacktrackerRecord  btr 
)
private

Definition at line 500 of file PDFastSimPAR_module.cc.

504  {
505  size_t iChan = btr.OpDetNum();
506  auto channelPosition = ChannelMap.find(iChan);
507 
508  if (channelPosition == ChannelMap.end()) {
509  ChannelMap[iChan] = opbtr.size();
510  opbtr.emplace_back(std::move(btr));
511  }
512  else {
513  unsigned int idtest = channelPosition->second;
514  auto const& timePDclockSDPsMap = btr.timePDclockSDPsMap();
515 
516  for (auto const& timePDclockSDP : timePDclockSDPsMap) {
517  for (auto const& sdp : timePDclockSDP.second) {
518  double xyz[3] = {sdp.x, sdp.y, sdp.z};
519  opbtr.at(idtest).AddScintillationPhotons(
520  sdp.trackID, timePDclockSDP.first, sdp.numPhotons, xyz, sdp.energy);
521  }
522  }
523  }
524  }
int OpDetNum() const
Returns the readout Optical Detector this object describes.
def move(depos, offset)
Definition: depos.py:107
timePDclockSDPs_t const & timePDclockSDPsMap() const
Returns all the deposited energy information as stored.
void phot::PDFastSimPAR::detectedNumPhotons ( std::map< size_t, int > &  DetectedNumPhotons,
const std::map< size_t, double > &  OpDetVisibilities,
const double  NumPhotons 
)
private

Definition at line 584 of file PDFastSimPAR_module.cc.

585  {
586  for (auto const& x : OpDetVisibilities)
587  {
588  DetectedNumPhotons[x.first] = fRandPoissPhot->fire(x.second * NumPhotons);
589  }
590  }
std::unique_ptr< CLHEP::RandPoissonQ > fRandPoissPhot
list x
Definition: train.py:276
void phot::PDFastSimPAR::Initialization ( )
private

Definition at line 528 of file PDFastSimPAR_module.cc.

529  {
530  std::cout << "PDFastSimPAR Initialization" << std::endl;
531  std::cout << "Initializing the geometry of the detector." << std::endl;
532  std::cout << "Simulate using semi-analytic model for number of hits." << std::endl;
533 
534  fRandPoissPhot = std::make_unique<CLHEP::RandPoissonQ>(fPhotonEngine);
535  geo::GeometryCore const& geom = *(lar::providerFrom<geo::Geometry>());
536 
537  // photo-detector visibility model (semi-analytical model)
539 
540  // propagation time model
542 
543  // Store info from the Geometry service
544  nOpDets = geom.NOpDets();
546  fNTPC = geom.NTPC();
547 
548  {
549  auto log = mf::LogTrace("PDFastSimPAR") << "PDFastSimPAR: active volume boundaries from "
550  << fActiveVolumes.size() << " volumes:";
551  for (auto const& [iCryo, box] : util::enumerate(fActiveVolumes)) {
552  log << "\n - C:" << iCryo << ": " << box.Min() << " -- " << box.Max() << " cm";
553  }
554  } // local scope
555 
556  if (geom.Ncryostats() > 1U) {
557  if (fOnlyOneCryostat) {
558  mf::LogWarning("PDFastSimPAR")
559  << std::string(80, '=') << "\nA detector with " << geom.Ncryostats()
560  << " cryostats is configured"
561  << " , and semi-analytic model is requested for scintillation photon propagation."
562  << " THIS CONFIGURATION IS NOT SUPPORTED and it is open to bugs"
563  << " (e.g. scintillation may be detected only in cryostat #0)."
564  << "\nThis would be normally a fatal error, but it has been forcibly overridden."
565  << "\n"
566  << std::string(80, '=');
567  }
568  else {
570  << "Photon propagation via semi-analytic model is not supported yet"
571  << " on detectors with more than one cryostat.";
572  }
573  }
574 
575  for (size_t const i : util::counter(nOpDets)) {
576  geo::OpDetGeo const& opDet = geom.OpDetGeoFromOpDet(i);
577  fOpDetCenter.push_back(opDet.GetCenter());
578  }
579  }
std::vector< geo::Point_t > fOpDetCenter
fhicl::ParameterSet fVISHitsParams
static std::vector< geo::BoxBoundedGeo > extractActiveLArVolume(geo::GeometryCore const &geom)
Definition: ISTPC.cxx:49
std::string string
Definition: nybbler.cc:12
std::unique_ptr< PropagationTimeModel > fPropTimeModel
std::vector< geo::BoxBoundedGeo > fActiveVolumes
void GetCenter(double *xyz, double localz=0.0) const
Definition: OpDetGeo.cxx:40
auto enumerate(Iterables &&...iterables)
Range-for loop helper tracking the number of iteration.
Definition: enumerate.h:69
std::unique_ptr< CLHEP::RandPoissonQ > fRandPoissPhot
CLHEP::HepRandomEngine & fScintTimeEngine
fhicl::ParameterSet fVISTimingParams
auto counter(T begin, T end)
Returns an object to iterate values from begin to end in a range-for loop.
Definition: counter.h:285
std::unique_ptr< SemiAnalyticalModel > fVisibilityModel
fhicl::ParameterSet fVUVTimingParams
Description of geometry of one entire detector.
cet::coded_exception< errors::ErrorCodes, ExceptionDetail::translate > Exception
Definition: Exception.h:66
CLHEP::HepRandomEngine & fPhotonEngine
MaybeLogger_< ELseverityLevel::ELsev_warning, false > LogWarning
MaybeLogger_< ELseverityLevel::ELsev_success, true > LogTrace
QTextStream & endl(QTextStream &s)
fhicl::ParameterSet fVUVHitsParams
bool phot::PDFastSimPAR::isOpDetInSameTPC ( geo::Point_t const &  ScintPoint,
geo::Point_t const &  OpDetPoint 
) const
private

Definition at line 595 of file PDFastSimPAR_module.cc.

597  {
598  // check optical channel is in same TPC as scintillation light, if not doesn't see light
599  // temporary method working for SBND, uBooNE, DUNE 1x2x6; to be replaced to work in full DUNE geometry
600  // check x coordinate has same sign or is close to zero
601  if (((ScintPoint.X() < 0.) != (OpDetPoint.X() < 0.)) &&
602  std::abs(OpDetPoint.X()) > 10. && fNTPC == 2) { // TODO: replace with geometry service method
603  return false;
604  }
605  return true;
606  }
T abs(T value)
void phot::PDFastSimPAR::produce ( art::Event event)
overridevirtual

Implements art::EDProducer.

Definition at line 269 of file PDFastSimPAR_module.cc.

270  {
271  mf::LogTrace("PDFastSimPAR") << "PDFastSimPAR Module Producer"
272  << "EventID: " << event.event();
273 
274  auto phot = std::make_unique<std::vector<sim::SimPhotons>>();
275  auto phlit = std::make_unique<std::vector<sim::SimPhotonsLite>>();
276  auto opbtr = std::make_unique<std::vector<sim::OpDetBacktrackerRecord>>();
277 
278  auto phot_ref = std::make_unique<std::vector<sim::SimPhotons>>();
279  auto phlit_ref = std::make_unique<std::vector<sim::SimPhotonsLite>>();
280  auto opbtr_ref = std::make_unique<std::vector<sim::OpDetBacktrackerRecord>>();
281 
282  auto& dir_photcol(*phot);
283  auto& ref_photcol(*phot_ref);
284  auto& dir_phlitcol(*phlit);
285  auto& ref_phlitcol(*phlit_ref);
286  dir_photcol.resize(nOpDets);
287  ref_photcol.resize(nOpDets);
288  dir_phlitcol.resize(nOpDets);
289  ref_phlitcol.resize(nOpDets);
290  for (unsigned int i = 0; i < nOpDets; i ++)
291  {
292  dir_photcol[i].fOpChannel = i;
293  ref_photcol[i].fOpChannel = i;
294  dir_phlitcol[i].OpChannel = i;
295  ref_phlitcol[i].OpChannel = i;
296  }
297 
299  if (!event.getByLabel(simTag, edepHandle)) {
300  mf::LogError("PDFastSimPAR") << "PDFastSimPAR Module Cannot getByLabel: " << simTag;
301  return;
302  }
303 
304  auto const& edeps = edepHandle;
305 
306  int num_points = 0;
307  int num_fastph = 0;
308  int num_slowph = 0;
309  int num_fastdp = 0;
310  int num_slowdp = 0;
311 
312  for (auto const& edepi : *edeps) {
313  num_points++;
314 
315  int trackID = edepi.TrackID();
316  double nphot = edepi.NumPhotons();
317  double edeposit = edepi.Energy() / nphot;
318  double pos[3] = {edepi.MidPointX(), edepi.MidPointY(), edepi.MidPointZ()};
319  geo::Point_t const ScintPoint = {pos[0], pos[1], pos[2]};
320 
321  if (fOnlyActiveVolume && !fISTPC.isScintInActiveVolume(ScintPoint)) continue;
322 
323  double nphot_fast = edepi.NumFPhotons();
324  double nphot_slow = edepi.NumSPhotons();
325 
326  num_fastph += nphot_fast;
327  num_slowph += nphot_slow;
328 
329  // direct light
330  std::map<size_t, int> DetectedNumFast;
331  std::map<size_t, int> DetectedNumSlow;
332 
333  bool needHits = (nphot_fast > 0 && fDoFastComponent) || (nphot_slow > 0 && fDoSlowComponent);
334  if ( needHits ) {
335  std::map<size_t, double> OpDetVisibilities;
336  fVisibilityModel->detectedDirectVisibilities(OpDetVisibilities, ScintPoint);
337  detectedNumPhotons(DetectedNumFast, OpDetVisibilities, nphot_fast);
338  detectedNumPhotons(DetectedNumSlow, OpDetVisibilities, nphot_slow);
339 
340  if ( fIncludeAnodeReflections ) {
341  std::map<size_t, int> AnodeDetectedNumFast;
342  std::map<size_t, int> AnodeDetectedNumSlow;
343 
344  std::map<size_t, double> OpDetVisibilitiesAnode;
345  fVisibilityModel->detectedReflectedVisibilities(OpDetVisibilitiesAnode, ScintPoint, true);
346  detectedNumPhotons(AnodeDetectedNumFast, OpDetVisibilitiesAnode, nphot_fast);
347  detectedNumPhotons(AnodeDetectedNumSlow, OpDetVisibilitiesAnode, nphot_slow);
348 
349  // add to existing count
350  for (auto const& x : AnodeDetectedNumFast) DetectedNumFast[x.first] += x.second;
351  for (auto const& x : AnodeDetectedNumSlow) DetectedNumSlow[x.first] += x.second;
352  }
353  }
354 
355  // reflected light, if enabled
356  std::map<size_t, int> ReflDetectedNumFast;
357  std::map<size_t, int> ReflDetectedNumSlow;
358 
359  if (fDoReflectedLight && needHits) {
360  std::map<size_t, double> OpDetVisibilitiesRefl;
361  fVisibilityModel->detectedReflectedVisibilities(OpDetVisibilitiesRefl, ScintPoint, false);
362  detectedNumPhotons(ReflDetectedNumFast, OpDetVisibilitiesRefl, nphot_fast);
363  detectedNumPhotons(ReflDetectedNumSlow, OpDetVisibilitiesRefl, nphot_slow);
364  }
365 
366  // propagation time
367  std::vector<double> transport_time;
368 
369  // loop through direct photons then reflected photons cases
370  for (size_t Reflected = 0; Reflected <= 1; ++Reflected) {
371 
372  // only do the reflected loop if including reflected light
373  if (Reflected && !fDoReflectedLight) continue;
374 
375  for (size_t channel = 0; channel < nOpDets; channel++) {
376 
377  if (fOpaqueCathode && !isOpDetInSameTPC(ScintPoint, fOpDetCenter[channel])) continue;
378 
379  int ndetected_fast = DetectedNumFast[channel];
380  int ndetected_slow = DetectedNumSlow[channel];
381  if (Reflected) {
382  ndetected_fast = ReflDetectedNumFast[channel];
383  ndetected_slow = ReflDetectedNumSlow[channel];
384  }
385 
386  // calculate propagation time, does not matter whether fast or slow photon
387  if (fIncludePropTime && needHits) {
388  transport_time.resize(ndetected_fast + ndetected_slow);
389  fPropTimeModel->propagationTime(transport_time, ScintPoint, channel, Reflected);
390  }
391 
392  // SimPhotonsLite case
393  if (fUseLitePhotons) {
394 
395  sim::OpDetBacktrackerRecord tmpbtr(channel);
396 
397  if (ndetected_fast > 0 && fDoFastComponent) {
398  int n = ndetected_fast;
399  num_fastdp += n;
400  for (long i = 0; i < n; ++i) {
401  // calculates the time at which the photon was produced
402  fScintTime->GenScintTime(true, fScintTimeEngine);
403  int time;
404  if (fIncludePropTime) time = static_cast<int>(edepi.StartT() + fScintTime->GetScintTime() + transport_time[i]);
405  else time = static_cast<int>(edepi.StartT() + fScintTime->GetScintTime());
406  if (Reflected) ++ref_phlitcol[channel].DetectedPhotons[time];
407  else ++dir_phlitcol[channel].DetectedPhotons[time];
408  tmpbtr.AddScintillationPhotons(trackID, time, 1, pos, edeposit);
409  }
410  }
411 
412  if (ndetected_slow > 0 && fDoSlowComponent) {
413  int n = ndetected_slow;
414  num_slowdp += n;
415  for (long i = 0; i < n; ++i) {
416  fScintTime->GenScintTime(false, fScintTimeEngine);
417  int time;
418  if (fIncludePropTime) time = static_cast<int>(edepi.StartT() + fScintTime->GetScintTime() + transport_time[ndetected_fast + i]);
419  else time = static_cast<int>(edepi.StartT() + fScintTime->GetScintTime());
420  if (Reflected) ++ref_phlitcol[channel].DetectedPhotons[time];
421  else ++dir_phlitcol[channel].DetectedPhotons[time];
422  tmpbtr.AddScintillationPhotons(trackID, time, 1, pos, edeposit);
423  }
424  }
425 
426  if (Reflected) AddOpDetBTR(*opbtr_ref, PDChannelToSOCMapReflect, tmpbtr);
427  else AddOpDetBTR(*opbtr, PDChannelToSOCMapDirect, tmpbtr);
428  }
429  // SimPhotons case
430  else {
431 
432  sim::OnePhoton photon;
433  photon.SetInSD = false;
434  photon.InitialPosition = edepi.End();
435  if (Reflected) photon.Energy = 2.9 * CLHEP::eV; // 430 nm
436  else photon.Energy = 9.7 * CLHEP::eV; // 128 nm
437 
438  if (ndetected_fast > 0 && fDoFastComponent) {
439  int n = ndetected_fast;
440  num_fastdp += n;
441  for (long i = 0; i < n; ++i) {
442  // calculates the time at which the photon was produced
443  fScintTime->GenScintTime(true, fScintTimeEngine);
444  int time;
445  if (fIncludePropTime) time = static_cast<int>(edepi.StartT() + fScintTime->GetScintTime() + transport_time[i]);
446  else time = static_cast<int>(edepi.StartT() + fScintTime->GetScintTime());
447  photon.Time = time;
448  if(Reflected) ref_photcol[channel].insert(ref_photcol[channel].end(), 1, photon);
449  else dir_photcol[channel].insert(dir_photcol[channel].end(), 1, photon);
450  }
451  }
452 
453  if (ndetected_slow > 0 && fDoSlowComponent) {
454  int n = ndetected_slow;
455  num_slowdp += n;
456  for (long i = 0; i < n; ++i) {
457  fScintTime->GenScintTime(false, fScintTimeEngine);
458  int time;
459  if (fIncludePropTime) time = static_cast<int>(edepi.StartT() + fScintTime->GetScintTime() + transport_time[ndetected_fast + i]);
460  else time = static_cast<int>(edepi.StartT() + fScintTime->GetScintTime());
461  photon.Time = time;
462  if(Reflected) ref_photcol[channel].insert(ref_photcol[channel].end(), 1, photon);
463  else dir_photcol[channel].insert(dir_photcol[channel].end(), 1, photon);
464  }
465  }
466  }
467  }
468  }
469  }
470 
471  mf::LogTrace("PDFastSimPAR") << "Total points: " << num_points
472  << ", total fast photons: " << num_fastph
473  << ", total slow photons: " << num_slowph
474  << "\ndetected fast photons: " << num_fastdp
475  << ", detected slow photons: " << num_slowdp;
476 
477  PDChannelToSOCMapDirect.clear();
478  PDChannelToSOCMapReflect.clear();
479 
480  if (fUseLitePhotons) {
481  event.put(move(phlit));
482  event.put(move(opbtr));
483  if (fDoReflectedLight) {
484  event.put(move(phlit_ref), "Reflected");
485  event.put(move(opbtr_ref), "Reflected");
486  }
487  }
488  else {
489  event.put(move(phot));
490  if (fDoReflectedLight) {
491  event.put(move(phot_ref), "Reflected");
492  }
493  }
494 
495  return;
496  }
std::map< size_t, int > PDChannelToSOCMapDirect
end
while True: pbar.update(maxval-len(onlies[E][S])) #print iS, "/", len(onlies[E][S]) found = False for...
std::vector< geo::Point_t > fOpDetCenter
bool isOpDetInSameTPC(geo::Point_t const &ScintPoint, geo::Point_t const &OpDetPoint) const
void detectedNumPhotons(std::map< size_t, int > &DetectedNumPhotons, const std::map< size_t, double > &OpDetVisibilities, const double NumPhotons)
std::unique_ptr< PropagationTimeModel > fPropTimeModel
All information of a photon entering the sensitive optical detector volume.
Definition: SimPhotons.h:64
std::unique_ptr< ScintTime > fScintTime
uint8_t channel
Definition: CRTFragment.hh:201
MaybeLogger_< ELseverityLevel::ELsev_error, false > LogError
Energy deposited on a readout Optical Detector by simulated tracks.
geo::Point_t InitialPosition
Scintillation position in world coordinates [cm].
Definition: SimPhotons.h:67
CLHEP::HepRandomEngine & fScintTimeEngine
bool getByLabel(std::string const &label, std::string const &instance, Handle< PROD > &result) const
Definition: DataViewImpl.h:633
std::map< size_t, int > PDChannelToSOCMapReflect
static constexpr double eV
Definition: Units.h:127
std::unique_ptr< SemiAnalyticalModel > fVisibilityModel
std::void_t< T > n
def move(depos, offset)
Definition: depos.py:107
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:184
General LArSoft Utilities.
bool SetInSD
Whether the photon reaches the sensitive detector.
Definition: SimPhotons.h:88
list x
Definition: train.py:276
MaybeLogger_< ELseverityLevel::ELsev_success, true > LogTrace
float Energy
Scintillation photon energy [GeV].
Definition: SimPhotons.h:82
bool isScintInActiveVolume(geo::Point_t const &ScintPoint)
Definition: ISTPC.cxx:41
void AddOpDetBTR(std::vector< sim::OpDetBacktrackerRecord > &opbtr, std::map< size_t, int > &ChannelMap, sim::OpDetBacktrackerRecord btr)

Member Data Documentation

std::vector<geo::BoxBoundedGeo> phot::PDFastSimPAR::fActiveVolumes
private

Definition at line 158 of file PDFastSimPAR_module.cc.

bool phot::PDFastSimPAR::fDoFastComponent
private

Definition at line 170 of file PDFastSimPAR_module.cc.

bool phot::PDFastSimPAR::fDoReflectedLight
private

Definition at line 172 of file PDFastSimPAR_module.cc.

bool phot::PDFastSimPAR::fDoSlowComponent
private

Definition at line 171 of file PDFastSimPAR_module.cc.

bool phot::PDFastSimPAR::fGeoPropTimeOnly
private

Definition at line 175 of file PDFastSimPAR_module.cc.

bool phot::PDFastSimPAR::fIncludeAnodeReflections
private

Definition at line 173 of file PDFastSimPAR_module.cc.

bool phot::PDFastSimPAR::fIncludePropTime
private

Definition at line 174 of file PDFastSimPAR_module.cc.

larg4::ISTPC phot::PDFastSimPAR::fISTPC
private

Definition at line 139 of file PDFastSimPAR_module.cc.

int phot::PDFastSimPAR::fNTPC
private

Definition at line 159 of file PDFastSimPAR_module.cc.

bool phot::PDFastSimPAR::fOnlyActiveVolume
private

Definition at line 178 of file PDFastSimPAR_module.cc.

bool phot::PDFastSimPAR::fOnlyOneCryostat
private

Definition at line 179 of file PDFastSimPAR_module.cc.

bool phot::PDFastSimPAR::fOpaqueCathode
private

Definition at line 177 of file PDFastSimPAR_module.cc.

std::vector<geo::Point_t> phot::PDFastSimPAR::fOpDetCenter
private

Definition at line 162 of file PDFastSimPAR_module.cc.

CLHEP::HepRandomEngine& phot::PDFastSimPAR::fPhotonEngine
private

Definition at line 148 of file PDFastSimPAR_module.cc.

std::unique_ptr<PropagationTimeModel> phot::PDFastSimPAR::fPropTimeModel
private

Definition at line 145 of file PDFastSimPAR_module.cc.

std::unique_ptr<CLHEP::RandPoissonQ> phot::PDFastSimPAR::fRandPoissPhot
private

Definition at line 149 of file PDFastSimPAR_module.cc.

std::unique_ptr<ScintTime> phot::PDFastSimPAR::fScintTime
private

Definition at line 180 of file PDFastSimPAR_module.cc.

CLHEP::HepRandomEngine& phot::PDFastSimPAR::fScintTimeEngine
private

Definition at line 150 of file PDFastSimPAR_module.cc.

bool phot::PDFastSimPAR::fUseLitePhotons
private

Definition at line 176 of file PDFastSimPAR_module.cc.

fhicl::ParameterSet phot::PDFastSimPAR::fVISHitsParams
private

Definition at line 186 of file PDFastSimPAR_module.cc.

std::unique_ptr<SemiAnalyticalModel> phot::PDFastSimPAR::fVisibilityModel
private

Definition at line 142 of file PDFastSimPAR_module.cc.

fhicl::ParameterSet phot::PDFastSimPAR::fVISTimingParams
private

Definition at line 184 of file PDFastSimPAR_module.cc.

fhicl::ParameterSet phot::PDFastSimPAR::fVUVHitsParams
private

Definition at line 185 of file PDFastSimPAR_module.cc.

fhicl::ParameterSet phot::PDFastSimPAR::fVUVTimingParams
private

Definition at line 183 of file PDFastSimPAR_module.cc.

size_t phot::PDFastSimPAR::nOpDets
private

Definition at line 152 of file PDFastSimPAR_module.cc.

std::map<size_t, int> phot::PDFastSimPAR::PDChannelToSOCMapDirect
private

Definition at line 154 of file PDFastSimPAR_module.cc.

std::map<size_t, int> phot::PDFastSimPAR::PDChannelToSOCMapReflect
private

Definition at line 155 of file PDFastSimPAR_module.cc.

art::InputTag phot::PDFastSimPAR::simTag
private

Definition at line 169 of file PDFastSimPAR_module.cc.


The documentation for this class was generated from the following file: