LArSoft  v10_04_05
Liquid Argon Software toolkit - https://larsoft.org/
lar::example::AnalysisExample Class Reference

Example analyzer module. More...

Inheritance diagram for lar::example::AnalysisExample:
art::EDAnalyzer art::detail::Analyzer art::detail::LegacyModule art::Observer art::ModuleBase

Classes

struct  Config
 

Public Types

using Parameters = art::EDAnalyzer::Table< Config >
 
using ModuleType = EDAnalyzer
 

Public Member Functions

 AnalysisExample (Parameters const &config)
 Constructor: configures the module (see the Config structure above) More...
 
void beginJob () override
 
void beginRun (const art::Run &run) override
 
void analyze (const art::Event &event) override
 
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)
 
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)
 
std::unique_ptr< Worker > makeWorker (WorkerParams const &wp)
 
template<typename T , BranchType BT>
ViewToken< T > consumesView (InputTag const &tag)
 
template<typename T , BranchType BT>
ViewToken< T > mayConsumeView (InputTag const &tag)
 

Protected Member Functions

std::string const & processName () const
 
bool wantAllEvents () const noexcept
 
bool wantEvent (ScheduleID id, Event const &e) const
 
Handle< TriggerResults > getTriggerResults (Event const &e) const
 
ConsumesCollector & consumesCollector ()
 
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 ()
 

Private Attributes

art::InputTag fSimulationProducerLabel
 
art::InputTag fHitProducerLabel
 The name of the producer that created hits. More...
 
art::InputTag fClusterProducerLabel
 
int fSelectedPDG
 PDG code of particle we'll focus on. More...
 
double fBinSize
 For dE/dx work: the value of dx. More...
 
TH1D * fPDGCodeHist
 PDG code of all particles. More...
 
TH1D * fMomentumHist
 momentum [GeV] of all selected particles More...
 
TH1D * fTrackLengthHist
 true length [cm] of all selected particles More...
 
TTree * fSimulationNtuple
 tuple with simulated data More...
 
TTree * fReconstructionNtuple
 tuple with reconstructed data More...
 
geo::GeometryCore const * fGeometryService
 pointer to Geometry provider More...
 
double fElectronsToGeV
 conversion factor More...
 
int fTriggerOffset
 (units of ticks) time of expected neutrino event More...
 
The variables that will go into both n-tuples.
int fEvent
 number of the event being processed More...
 
int fRun
 number of the run being processed More...
 
int fSubRun
 
The variables that will go into the simulation n-tuple.
int fSimPDG
 PDG ID of the particle being processed. More...
 
int fSimTrackID
 GEANT ID of the particle being processed. More...
 
double fStartXYZT [4]
 (x,y,z,t) of the true start of the particle More...
 
double fEndXYZT [4]
 (x,y,z,t) of the true end of the particle More...
 
double fStartPE [4]
 (Px,Py,Pz,E) at the true start of the particle More...
 
double fEndPE [4]
 (Px,Py,Pz,E) at the true end of the particle More...
 
int fSimNdEdxBins
 Number of dE/dx bins in a given track. More...
 
std::vector< double > fSimdEdxBins
 The vector that will be used to accumulate dE/dx values as a function of range. More...
 
Variables used in the reconstruction n-tuple
int fRecoPDG
 PDG ID of the particle being processed. More...
 
int fRecoTrackID
 GEANT ID of the particle being processed. More...
 
int fRecoNdEdxBins
 Number of dE/dx bins in a given track. More...
 
std::vector< double > fRecodEdxBins
 The vector that will be used to accumulate dE/dx values as a function of range. More...
 

Detailed Description

Example analyzer module.

See also
analysis module example overview

This class extracts information from the generated and reconstructed particles.

It produces histograms for the simulated particles in the input file:

  • PDG ID (flavor) of all particles
  • momentum of the primary particles selected to have a specific PDG ID
  • length of the selected particle trajectory

It also produces two ROOT trees.

The first ROOT tree contains information on the simulated particles, including "dEdx", a binned histogram of collected charge as function of track range.

The second ROOT tree contains information on the reconstructed hits.

Configuration parameters

  • SimulationLabel (string, default: "largeant"): tag of the input data product with the detector simulation information (typically an instance of the LArG4 module)
  • HitLabel (string, mandatory): tag of the input data product with reconstructed hits
  • ClusterLabel (string, mandatory): tag of the input data product with reconstructed clusters
  • PDGcode (integer, mandatory): particle type (PDG ID) to be selected; only primary particles of this type will be considered
  • BinSize (real, mandatory): $ dx $ [cm] used for the $ dE/dx $ calculation

Definition at line 161 of file AnalysisExample_module.cc.

Member Typedef Documentation

Definition at line 22 of file EDAnalyzer.h.

Constructor & Destructor Documentation

lar::example::AnalysisExample::AnalysisExample ( Parameters const &  config)
explicit

Constructor: configures the module (see the Config structure above)

Definition at line 324 of file AnalysisExample_module.cc.

References fClusterProducerLabel, fGeometryService, fHitProducerLabel, fSimulationProducerLabel, fTriggerOffset, and detinfo::trigger_offset().

325  : EDAnalyzer(config)
326  , fSimulationProducerLabel(config().SimulationLabel())
327  , fHitProducerLabel(config().HitLabel())
328  , fClusterProducerLabel(config().ClusterLabel())
329  , fSelectedPDG(config().PDGcode())
330  , fBinSize(config().BinSize())
331  {
332  // Get a pointer to the geometry service provider.
333  fGeometryService = lar::providerFrom<geo::Geometry>();
334 
335  auto const clock_data =
337  fTriggerOffset = trigger_offset(clock_data);
338 
339  // Since art 2.8, you can and should tell beforehand, here in the constructor, all the
340  // data the module is going to read ("consumes") or might read
341  // ("may_consume"). Diligence here will in the future help the framework execute
342  // modules in parallel, making sure the order is correct.
343  consumes<std::vector<simb::MCParticle>>(fSimulationProducerLabel);
344  consumes<std::vector<sim::SimChannel>>(fSimulationProducerLabel);
345  consumes<art::Assns<simb::MCTruth, simb::MCParticle>>(fSimulationProducerLabel);
346  consumes<std::vector<recob::Hit>>(fHitProducerLabel);
347  consumes<std::vector<recob::Cluster>>(fClusterProducerLabel);
348  consumes<art::Assns<recob::Cluster, recob::Hit>>(fHitProducerLabel);
349  }
art::InputTag fHitProducerLabel
The name of the producer that created hits.
geo::GeometryCore const * fGeometryService
pointer to Geometry provider
EDAnalyzer(fhicl::ParameterSet const &pset)
Definition: EDAnalyzer.cc:6
int trigger_offset(DetectorClocksData const &data)
int fTriggerOffset
(units of ticks) time of expected neutrino event
double fBinSize
For dE/dx work: the value of dx.
int fSelectedPDG
PDG code of particle we&#39;ll focus on.

Member Function Documentation

void lar::example::AnalysisExample::analyze ( const art::Event event)
override

Definition at line 421 of file AnalysisExample_module.cc.

References bin, tca::dEdx(), DEFINE_ART_MODULE, energy, fBinSize, fClusterProducerLabel, fElectronsToGeV, fEndPE, fEndXYZT, fEvent, fHitProducerLabel, fMomentumHist, fPDGCodeHist, fRecodEdxBins, fRecoNdEdxBins, fReconstructionNtuple, fRecoPDG, fRecoTrackID, fRun, fSelectedPDG, fSimdEdxBins, fSimNdEdxBins, fSimPDG, fSimTrackID, fSimulationNtuple, fSimulationProducerLabel, fStartPE, fStartXYZT, fSubRun, fTrackLengthHist, Get, art::ProductRetriever::getByLabel(), hits(), geo::kCollection, geo::TPCGeo::Length(), MF_LOG_DEBUG, simb::MCParticle::PdgCode(), simb::MCParticle::Position(), simb::MCParticle::Process(), and geo::GeometryCore::TPC().

422  {
423  // Start by fetching some basic event information for our n-tuple.
424  fEvent = event.id().event();
425  fRun = event.run();
426  fSubRun = event.subRun();
427 
428  // This is the standard method of reading multiple objects associated with the same
429  // event; see
430  // <https://cdcvs.fnal.gov/redmine/projects/larsoft/wiki/Using_art_in_LArSoft> for
431  // more information.
432  //
433  // Define a "handle" to point to a vector of the objects.
435 
436  // Then tell the event to fill the vector with all the objects of that type produced
437  // by a particular producer.
438  //
439  // Note that if I don't test whether this is successful, and there aren't any
440  // simb::MCParticle objects, then the first time we access particleHandle, art will
441  // display a "ProductNotFound" exception message and, depending on art settings, it
442  // may skip all processing for the rest of this event (including any subsequent
443  // analysis steps) or stop the execution.
444  if (!event.getByLabel(fSimulationProducerLabel, particleHandle)) {
445  // If we have no MCParticles at all in an event, then we're in big trouble. Throw an
446  // exception to force this module to stop. Try to include enough information for the
447  // user to figure out what's going on. Note that when we throw a cet::exception, the
448  // run and event number will automatically be displayed.
449  //
450  // __LINE__ and __FILE__ are values computed by the compiler.
451 
452  throw cet::exception("AnalysisExample")
453  << " No simb::MCParticle objects in this event - "
454  << " Line " << __LINE__ << " in file " << __FILE__ << std::endl;
455  }
456 
457  // Get all the simulated channels for the event. These channels include the energy
458  // deposited for each simulated track.
459  //
460  // Here we use a different method to access objects: art::ValidHandle. Using this
461  // method, if there aren't any objects of the given type (sim::SimChannel in this
462  // case) in the input file for this art::Event, art will throw a ProductNotFound
463  // exception.
464 
465  auto simChannelHandle =
466  event.getValidHandle<std::vector<sim::SimChannel>>(fSimulationProducerLabel);
467 
468  //
469  // Let's compute the variables for the simulation n-tuple first.
470  //
471 
472  // The MCParticle objects are not necessarily in any particular order. Since we may
473  // have to search the list of particles, let's put them into a map, a sorted container
474  // that will make searching fast and easy. To save both space and time, the map will
475  // not contain a copy of the MCParticle, but a pointer to it.
476  std::map<int, const simb::MCParticle*> particleMap;
477 
478  auto const& wireReadoutGeom = art::ServiceHandle<geo::WireReadout const>()->Get();
479  for (auto const& particle : *particleHandle) {
480  // For the methods you can call to get particle information, see
481  // ${NUSIMDATA_INC}/nusimdata/SimulationBase/MCParticle.h.
482  fSimTrackID = particle.TrackId();
483 
484  // Add the address of the MCParticle to the map, with the track ID as the key.
485  particleMap[fSimTrackID] = &particle;
486 
487  // Histogram the PDG code of every particle in the event.
488  fSimPDG = particle.PdgCode();
489  fPDGCodeHist->Fill(fSimPDG);
490 
491  // For this example, we want to fill the n-tuples and histograms only with
492  // information from the primary particles in the event, whose PDG codes match a
493  // value supplied in the .fcl file.
494  if (particle.Process() != "primary" || fSimPDG != fSelectedPDG) continue;
495 
496  // A particle has a trajectory, consisting of a set of 4-positions and 4-mommenta.
497  const size_t numberTrajectoryPoints = particle.NumberTrajectoryPoints();
498 
499  // For trajectories, as for vectors and arrays, the first point is #0, not #1.
500  const int last = numberTrajectoryPoints - 1;
501  const TLorentzVector& positionStart = particle.Position(0);
502  const TLorentzVector& positionEnd = particle.Position(last);
503  const TLorentzVector& momentumStart = particle.Momentum(0);
504  const TLorentzVector& momentumEnd = particle.Momentum(last);
505 
506  // Make a histogram of the starting momentum.
507  fMomentumHist->Fill(momentumStart.P());
508 
509  // Fill arrays with the 4-values. (Don't be fooled by the name of the method; it
510  // just puts the numbers from the 4-vector into the array.)
511  positionStart.GetXYZT(fStartXYZT);
512  positionEnd.GetXYZT(fEndXYZT);
513  momentumStart.GetXYZT(fStartPE);
514  momentumEnd.GetXYZT(fEndPE);
515 
516  // Use a polar-coordinate view of the 4-vectors to get the track length.
517  const double trackLength = (positionEnd - positionStart).Rho();
518 
519  // Let's print some debug information in the job output to see that everything is
520  // fine. MF_LOG_DEBUG() is a messagefacility macro that prints stuff when the
521  // message level is set to standard_debug in the .fcl file.
522  MF_LOG_DEBUG("AnalysisExample") << "Track length: " << trackLength << " cm";
523 
524  // Fill a histogram of the track length.
525  fTrackLengthHist->Fill(trackLength);
526 
527  MF_LOG_DEBUG("AnalysisExample")
528  << "track ID=" << fSimTrackID << " (PDG ID: " << fSimPDG << ") " << trackLength
529  << " cm long, momentum " << momentumStart.P() << " GeV/c, has " << numberTrajectoryPoints
530  << " trajectory points";
531 
532  // Determine the number of dE/dx bins for the n-tuple.
533  fSimNdEdxBins = int(trackLength / fBinSize) + 1;
534  // Initialize the vector of dE/dx bins to be empty.
535  fSimdEdxBins.clear();
536 
537  // To look at the energy deposited by this particle's track, we loop over the
538  // SimChannel objects in the event.
539  for (auto const& channel : *simChannelHandle) {
540  // Get the numeric ID associated with this channel. (The channel number is a
541  // 32-bit unsigned int, which normally requires a special data type. Let's use
542  // "auto" so we don't have to remember "raw::ChannelID_t".
543  auto const channelNumber = channel.Channel();
544 
545  // A little care: There is more than one plane that reacts to charge in the
546  // TPC. We only want to include the energy from the collection plane. Note:
547  // geo::kCollection is defined in
548  // ${LARCOREOBJ_INC}/larcoreobj/SimpleTypesAndConstants/geo_types.h
549  if (wireReadoutGeom.SignalType(channelNumber) != geo::kCollection) continue;
550 
551  // Each channel has a map inside it that connects a time slice to energy deposits
552  // in the detector. We'll use "auto", but it's worth noting that the full type of
553  // this map is std::map<unsigned short, std::vector<sim::IDE>>
554  auto const& timeSlices = channel.TDCIDEMap();
555 
556  // For every time slice in this channel:
557  for (auto const& timeSlice : timeSlices) {
558  // Each entry in a map is a pair<first,second>. For the timeSlices map, the
559  // 'first' is a time slice number. The 'second' is a vector of IDE objects.
560  auto const& energyDeposits = timeSlice.second;
561 
562  // Loop over the energy deposits. An "energy deposit" object is something that
563  // knows how much charge/energy was deposited in a small volume, by which
564  // particle, and where. The type of 'energyDeposit' will be sim::IDE, which is
565  // defined in ${LARDATAOBJ_INC}/lardataobj/Simulation/SimChannel.h.
566  for (auto const& energyDeposit : energyDeposits) {
567  // Check if the track that deposited the energy matches the track of the
568  // particle.
569  if (energyDeposit.trackID != fSimTrackID) continue;
570  // Get the (x,y,z) of the energy deposit.
571  TVector3 location(energyDeposit.x, energyDeposit.y, energyDeposit.z);
572 
573  // Distance from the start of the track.
574  const double distance = (location - positionStart.Vect()).Mag();
575 
576  // Into which bin of the dE/dx array do we add the energy?
577  const unsigned int bin = (unsigned int)(distance / fBinSize);
578 
579  // Is the dE/dx array big enough to include this bin?
580  if (fSimdEdxBins.size() < bin + 1) {
581  // Increase the array size to accomodate the new bin, padding it with
582  // zeros.
583  fSimdEdxBins.resize(bin + 1, 0.);
584  }
585 
586  // Add the energy deposit to that bin. (If you look at the definition of
587  // sim::IDE, you'll see that there's another way to get the energy. Are the
588  // two methods equivalent? Compare the results and see!)
589  fSimdEdxBins[bin] += energyDeposit.numElectrons * fElectronsToGeV;
590 
591  } // For each energy deposit
592  } // For each time slice
593  } // For each SimChannel
594 
595  // At this point we've filled in the values of all the variables and arrays we want
596  // to write to the n-tuple for this particle. The following command actually writes
597  // those values.
598  fSimulationNtuple->Fill();
599 
600  } // loop over all particles in the event.
601 
602  //
603  // Reconstruction n-tuple
604  //
605 
606  // All of the above is based on objects entirely produced by the simulation. Let's try
607  // to do something based on reconstructed objects. A Hit (see
608  // ${LARDATAOBJ_INC}/lardataobj/RecoBase/Hit.h) is a 2D object in a plane.
609 
610  // This code duplicates much of the code in the previous big simulation loop, and will
611  // produce the similar results. (It won't be identical, because of shaping and other
612  // factors; not every bit of charge in a channel ends up contributing to a hit.) The
613  // point is to show different methods of accessing information, not to produce
614  // profound physics results -- that part is up to you!
615 
616  // For the rest of this method, I'm going to assume you've read the comments in
617  // previous section.
618 
619  // Start by reading the Hits. We don't use art::ValidHandle here because there might
620  // be no hits in the input; e.g., we ran the simulation but did not run the
621  // reconstruction steps yet. If there are no hits we may as well skip the rest of this
622  // module.
623 
625  if (!event.getByLabel(fHitProducerLabel, hitHandle)) return;
626 
627  // Our goal is to accumulate the dE/dx of any particles associated with the hits that
628  // match our criteria: primary particles with the PDG code from the .fcl file. I don't
629  // know how many such particles there will be in a given event. We'll use a map, with
630  // track ID as the key, to hold the vectors of dE/dx information.
631  std::map<int, std::vector<double>> dEdxMap;
632 
633  auto const clock_data =
635 
636  // For every Hit:
637  for (auto const& hit : *hitHandle) {
638  // The channel associated with this hit.
639  auto hitChannelNumber = hit.Channel();
640 
641  // We have a hit. For this example let's just focus on the hits in the collection
642  // plane.
643  if (wireReadoutGeom.SignalType(hitChannelNumber) != geo::kCollection) continue;
644 
645  MF_LOG_DEBUG("AnalysisExample") << "Hit in collection plane" << std::endl;
646 
647  // In a few lines we're going to look for possible energy deposits that correspond
648  // to that hit. Determine a reasonable range of times that might correspond to those
649  // energy deposits.
650 
651  // In reconstruction, the channel waveforms are truncated. So we have to adjust the
652  // Hit TDC ticks to match those of the SimChannels, which were created during
653  // simulation.
654 
655  using TDC_t = sim::SimChannel::StoredTDC_t;
656  TDC_t start_tdc = clock_data.TPCTick2TDC(hit.StartTick());
657  TDC_t end_tdc = clock_data.TPCTick2TDC(hit.EndTick());
658  TDC_t hitStart_tdc = clock_data.TPCTick2TDC(hit.PeakTime() - 3. * hit.SigmaPeakTime());
659  TDC_t hitEnd_tdc = clock_data.TPCTick2TDC(hit.PeakTime() + 3. * hit.SigmaPeakTime());
660 
661  start_tdc = std::max(start_tdc, hitStart_tdc);
662  end_tdc = std::min(end_tdc, hitEnd_tdc);
663 
664  // In the simulation section, we started with particles to find channels with a
665  // matching track ID. Now we search in reverse: search the SimChannels for matching
666  // channel number, then look at the tracks inside the channel.
667 
668  for (auto const& channel : *simChannelHandle) {
669  auto simChannelNumber = channel.Channel();
670  if (simChannelNumber != hitChannelNumber) continue;
671 
672  MF_LOG_DEBUG("AnalysisExample") << "SimChannel number = " << simChannelNumber << std::endl;
673 
674  // The time slices in this channel.
675  auto const& timeSlices = channel.TDCIDEMap();
676 
677  // We want to look over the range of time slices in this channel that correspond
678  // to the range of hit times.
679 
680  // We have to create "dummy" time slices for the search.
681  sim::TDCIDE startTime;
682  sim::TDCIDE endTime;
683  startTime.first = start_tdc;
684  endTime.first = end_tdc;
685 
686  // Here are the fast searches:
687 
688  // Find a pointer to the first channel with time >= start_tdc.
689  auto const startPointer =
690  std::lower_bound(timeSlices.begin(), timeSlices.end(), startTime, TDCIDETimeCompare);
691 
692  // From that time slice, find the last channel with time < end_tdc.
693  auto const endPointer =
694  std::upper_bound(startPointer, timeSlices.end(), endTime, TDCIDETimeCompare);
695 
696  // Did we find anything? If not, skip.
697  if (startPointer == timeSlices.end() || startPointer == endPointer) continue;
698  MF_LOG_DEBUG("AnalysisExample")
699  << "Time slice start = " << (*startPointer).first << std::endl;
700 
701  // Loop over the channel times we found that match the hit times.
702  for (auto slicePointer = startPointer; slicePointer != endPointer; slicePointer++) {
703  auto const timeSlice = *slicePointer;
704  auto time = timeSlice.first;
705 
706  // How to debug a problem: Lots of print statements. There are debuggers such as
707  // gdb, but they can be tricky to use with shared libraries and don't work if
708  // you're using software that was compiled somewhere else (e.g., you're
709  // accessing LArSoft libraries via CVMFS).
710 
711  // The MF_LOG_DEBUG statements below are left over from when I was trying to
712  // solve a problem about hit timing. I could have deleted them, but decided to
713  // leave them to demonsrate what a typical debugging process looks like.
714 
715  MF_LOG_DEBUG("AnalysisExample")
716  << "Hit index = " << hit.LocalIndex() << " channel number = " << hitChannelNumber
717  << " start TDC tick = " << hit.StartTick() << " end TDC tick = " << hit.EndTick()
718  << " peak TDC tick = " << hit.PeakTime() << " sigma peak time = " << hit.SigmaPeakTime()
719  << " adjusted start TDC tick = " << clock_data.TPCTick2TDC(hit.StartTick())
720  << " adjusted end TDC tick = " << clock_data.TPCTick2TDC(hit.EndTick())
721  << " adjusted peak TDC tick = " << clock_data.TPCTick2TDC(hit.PeakTime())
722  << " adjusted start_tdc = " << start_tdc << " adjusted end_tdc = " << end_tdc
723  << " time = " << time << std::endl;
724 
725  // Loop over the energy deposits.
726  auto const& energyDeposits = timeSlice.second;
727  for (auto const& energyDeposit : energyDeposits) {
728  // Remember that map of MCParticles we created near the top of this method?
729  // Now we can use it. Search the map for the track ID associated with this
730  // energy deposit. Since a map is automatically sorted, we can use a fast
731  // binary search method, 'find()'.
732 
733  // By the way, the type of "search" is an iterator (to be specific, it's an
734  // std::map<int,simb::MCParticle*>::const_iterator, which makes you thankful
735  // for the "auto" keyword). If you're going to work with C++ containers,
736  // you'll have to learn about iterators eventually; do a web search on "STL
737  // iterator" to get started.
738  auto search = particleMap.find(energyDeposit.trackID);
739 
740  // Did we find this track ID in the particle map? It's possible for the
741  // answer to be "no"; some particles are too low in kinetic energy to be
742  // written by the simulation (see ${LARSIM_DIR}/job/simulationservices.fcl,
743  // parameter ParticleKineticEnergyCut).
744  if (search == particleMap.end()) continue;
745 
746  // "search" points to a pair in the map: <track ID, MCParticle*>
747  int trackID = (*search).first;
748  const simb::MCParticle& particle = *((*search).second);
749 
750  // Is this a primary particle, with a PDG code that matches the user input?
751  if (particle.Process() != "primary" || particle.PdgCode() != fSelectedPDG) continue;
752 
753  // Determine the dE/dx of this particle.
754  const TLorentzVector& positionStart = particle.Position(0);
755  TVector3 location(energyDeposit.x, energyDeposit.y, energyDeposit.z);
756  double distance = (location - positionStart.Vect()).Mag();
757  unsigned int bin = int(distance / fBinSize);
758  double energy = energyDeposit.numElectrons * fElectronsToGeV;
759 
760  // A feature of maps: if we refer to dEdxMap[trackID], and there's no such
761  // entry in the map yet, it will be automatically created with a zero-size
762  // vector. Test to see if the vector for this track ID is big enough.
763  //
764  // dEdxMap is a map, which is a slow container compared to a vector. If we are
765  // going to access the same element over and over, it is a good idea to find
766  // that element once, and then refer to that item directly. Since we don't
767  // care about the specific type of dEdxMap[trackID] (a vector, by the way), we
768  // can use "auto" to save some time. This must be a reference, since we want
769  // to change the original value in the map, and can't be constant.
770  auto& track_dEdX = dEdxMap[trackID];
771  if (track_dEdX.size() < bin + 1) {
772  // Increase the vector size, padding it with zeroes.
773  track_dEdX.resize(bin + 1, 0);
774  }
775 
776  // Add the energy to the dE/dx bins for this track.
777  track_dEdX[bin] += energy;
778 
779  } // loop over energy deposits
780  } // loop over time slices
781  } // for each SimChannel
782  } // for each Hit
783 
784  // We have a map of dE/dx vectors. Write each one of them to the reconstruction
785  // n-tuple.
786  for (const auto& dEdxEntry : dEdxMap) {
787  // Here, the map entries are <first,second>=<track ID, dE/dx vector>
788  fRecoTrackID = dEdxEntry.first;
789 
790  // This is an example of how we'd pick out the PDG code if there are multiple
791  // particle types or tracks in a single event allowed in the n-tuple.
792  fRecoPDG = particleMap[fRecoTrackID]->PdgCode();
793 
794  // Get the number of bins for this track.
795  const std::vector<double>& dEdx = dEdxEntry.second;
796  fRecoNdEdxBins = dEdx.size();
797 
798  // Copy this track's dE/dx information.
800 
801  // At this point, we've filled in all the reconstruction n-tuple's variables. Write
802  // it.
803  fReconstructionNtuple->Fill();
804  }
805 
806  // Think about the two big loops above, One starts from the particles then looks at
807  // the channels; the other starts with the hits and backtracks to the particles. What
808  // links the information in simb::MCParticle and sim::SimChannel is the track ID
809  // number assigned by the LArG4 simulation; what links sim::SimChannel and recob::Hit
810  // is the channel ID.
811 
812  // In general, that is not how objects in the LArSoft reconstruction chain are linked
813  // together. Most of them are linked using associations and the art::Assns class:
814  // <https://cdcvs.fnal.gov/redmine/projects/larsoft/wiki/Using_art_in_LArSoft#artAssns>
815 
816  // The web page may a bit difficult to understand (at least, it is for me), so let's
817  // try to simplify things:
818 
819  // - If you're doing an analysis task, you don't have to worry about creating
820  // art::Assns objects.
821 
822  // - You don't have to read the art:Assns objects on your own. There are helper
823  // classes (FindXXX) which will do that for you.
824 
825  // - There's only one helper you need: art::FindManyP. It will do what you want with a
826  // minimum of fuss.
827 
828  // - Associations are symmetric. If you see an art::Assns<ClassA,ClassB>, the helper
829  // classes will find all of ClassA associated with ClassB or all of ClassB
830  // associated with ClassA.
831 
832  // - To know what associations exist, you have to be a 'code detective'. The important
833  // clue is to look for a 'produces' line in the code that writes objects to an
834  // art::Event. For example, in ${LARSIM_DIR}/source/larsim/LArG4/LArG4_module.cc,
835  // you'll see this line:
836 
837  // produces< art::Assns<simb::MCTruth, simb::MCParticle> >();
838 
839  // That means a simulated event will have an association between simb::MCTruth (the
840  // primary particle produced by the event generator) and the simb::MCParticle
841  // objects (the secondary particles produced in the detector simulation).
842 
843  // Let's try it. The following statement will find the simb::MCTruth objects
844  // associated with the simb::MCParticle objects in the event (referenced by
845  // particleHandle):
846 
847  const art::FindManyP<simb::MCTruth> findManyTruth(
848  particleHandle, event, fSimulationProducerLabel);
849 
850  // Note that we still have to supply the module label of the step that created the
851  // association. Also note that we did not have to explicitly read in the simb::MCTruth
852  // objects from the art::Event object 'event'; FindManyP did that for us.
853 
854  // Also note that at this point art::FindManyP has already found all the simb::MCTruth
855  // associated with each of the particles in particleHandle. This is a slow process, so
856  // in general you want to do it only once. If we had a loop over the particles, we
857  // would still do this outside that loop.
858 
859  // Now we can query the 'findManyTruth' object to access the information. First, check
860  // that there wasn't some kind of error:
861 
862  if (!findManyTruth.isValid()) {
863  mf::LogError("AnalysisExample") << "findManyTruth simb::MCTruth for simb::MCParticle failed!";
864  }
865 
866  // I'm going to be lazy, and just look at the simb::MCTruth object associated with the
867  // first simb::MCParticle we read in. (The main reason I'm being lazy is that if I
868  // used the single-particle generator in prodsingle.fcl, every particle in the event
869  // is going to be associated with just the one primary particle from the event
870  // generator.)
871 
872  size_t particle_index = 0; // look at first particle in particleHandle's vector.
873 
874  // The result of FindManyP::at() is a vector of pointers, in this case
875  // simb::MCTruth*. In this case it will be a vector with just one entry; I could have
876  // used art::FindOneP instead. (This will be a vector of art::Ptr, which is a type of
877  // smart pointer; see
878  // <https://cdcvs.fnal.gov/redmine/projects/larsoft/wiki/Using_art_in_LArSoft#artPtrltTgt-and-artPtrVectorltTgt>
879  // To avoid unnecessary copying, and since art::FindManyP returns a constant
880  // reference, use "auto const&".
881 
882  auto const& truth = findManyTruth.at(particle_index);
883 
884  // Make sure there's no problem.
885  if (truth.empty()) {
886  mf::LogError("AnalysisExample")
887  << "Particle ID=" << particleHandle->at(particle_index).TrackId() << " has no primary!";
888  }
889 
890  // Use the message facility to write output. I don't have to write the event, run, or
891  // subrun numbers; the message facility takes care of that automatically. I'm "going
892  // at warp speed" with the vectors, pointers, and methods; see
893  // ${NUSIMDATA_INC}/nusimdata/SimulationBase/MCTruth.h and
894  // ${NUSIMDATA_INC}/nusimdata/SimulationBase/MCParticle.h for the nested calls I'm
895  // using.
896  mf::LogInfo("AnalysisExample") << "Particle ID=" << particleHandle->at(particle_index).TrackId()
897  << " primary PDG code=" << truth[0]->GetParticle(0).PdgCode();
898 
899  // Let's try a slightly more realistic example. Suppose I want to read in the
900  // clusters, and learn what hits are associated with them. Then I could backtrack from
901  // those hits to determine the dE/dx of the particles in the clusters. (Don't worry; I
902  // won't put you through the n-tuple creation procedure for a third time.)
903 
904  // First, read in the clusters (if there are any).
906  if (!event.getByLabel(fClusterProducerLabel, clusterHandle)) return;
907 
908  // Now use the associations to find which hits are associated with which
909  // clusters. Note that this is not as trivial a query as the one with MCTruth, since
910  // it's possible for a hit to be assigned to more than one cluster.
911 
912  // We have to include fClusterProducerLabel, since that's the step that created the
913  // art::Assns<recob::Hit,recob::Cluster> object; look at the modules in
914  // ${LARRECO_DIR}/source/larreco/ClusterFinder/ and search for the 'produces'
915  // lines. (I did not know this before I wrote these comments. I had to be a code
916  // detective and use UNIX tools like 'grep' and 'find' to locate those routines.)
917  const art::FindManyP<recob::Hit> findManyHits(clusterHandle, event, fClusterProducerLabel);
918 
919  if (!findManyHits.isValid()) {
920  mf::LogError("AnalysisExample") << "findManyHits recob::Hit for recob::Cluster failed;"
921  << " cluster label='" << fClusterProducerLabel << "'";
922  }
923 
924  // Now we'll loop over the clusters to see the hits associated with each one. Note
925  // that we're not using a range-based for loop. That's because FindManyP::at() expects
926  // a number as an argument, so we might as well go through the cluster objects via
927  // numerical index instead.
928  for (size_t cluster_index = 0; cluster_index != clusterHandle->size(); cluster_index++) {
929  // In this case, FindManyP::at() will return a vector of pointers to recob::Hit that
930  // corresponds to the "cluster_index"-th cluster.
931  auto const& hits = findManyHits.at(cluster_index);
932 
933  // We have a vector of pointers to the hits associated with the cluster, but for
934  // this example I'm not going to do anything fancy with them. I'll just print out
935  // how many there are.
936  mf::LogInfo("AnalysisExample") << "Cluster ID=" << clusterHandle->at(cluster_index).ID()
937  << " has " << hits.size() << " hits";
938  }
939 
940  } // AnalysisExample::analyze()
int fSimNdEdxBins
Number of dE/dx bins in a given track.
const TLorentzVector & Position(const int i=0) const
Definition: MCParticle.h:220
int PdgCode() const
Definition: MCParticle.h:213
TH1D * fPDGCodeHist
PDG code of all particles.
int fRecoTrackID
GEANT ID of the particle being processed.
std::pair< unsigned short, std::vector< sim::IDE > > TDCIDE
List of energy deposits at the same time (on this channel)
Definition: SimChannel.h:118
MaybeLogger_< ELseverityLevel::ELsev_info, false > LogInfo
art::InputTag fHitProducerLabel
The name of the producer that created hits.
double fEndXYZT[4]
(x,y,z,t) of the true end of the particle
std::vector< double > fSimdEdxBins
The vector that will be used to accumulate dE/dx values as a function of range.
double fStartPE[4]
(Px,Py,Pz,E) at the true start of the particle
int fSimTrackID
GEANT ID of the particle being processed.
cout<< "Opened file "<< fin<< " ixs= "<< ixs<< endl;if(ixs==0) hhh=(TH1F *) fff-> Get("h1")
Definition: AddMC.C:8
TTree * fReconstructionNtuple
tuple with reconstructed data
MaybeLogger_< ELseverityLevel::ELsev_error, false > LogError
std::string Process() const
Definition: MCParticle.h:216
int fRecoNdEdxBins
Number of dE/dx bins in a given track.
void hits()
Definition: readHits.C:15
double fStartXYZT[4]
(x,y,z,t) of the true start of the particle
double energy
Definition: plottest35.C:25
std::vector< double > fRecodEdxBins
The vector that will be used to accumulate dE/dx values as a function of range.
float bin[41]
Definition: plottest35.C:14
int fRecoPDG
PDG ID of the particle being processed.
float dEdx(detinfo::DetectorClocksData const &clockData, detinfo::DetectorPropertiesData const &detProp, const TCSlice &slc, TP3D &tp3d)
Definition: PFPUtils.cxx:2671
Detector simulation of raw signals on wires.
int fSimPDG
PDG ID of the particle being processed.
int fRun
number of the run being processed
TH1D * fMomentumHist
momentum [GeV] of all selected particles
bool getByLabel(std::string const &label, std::string const &instance, Handle< PROD > &result) const
int fEvent
number of the event being processed
#define MF_LOG_DEBUG(id)
double fEndPE[4]
(Px,Py,Pz,E) at the true end of the particle
TTree * fSimulationNtuple
tuple with simulated data
double fElectronsToGeV
conversion factor
double fBinSize
For dE/dx work: the value of dx.
cet::coded_exception< error, detail::translate > exception
Definition: exception.h:33
TH1D * fTrackLengthHist
true length [cm] of all selected particles
int fSelectedPDG
PDG code of particle we&#39;ll focus on.
TDCIDE::first_type StoredTDC_t
Type for TDC tick used in the internal representation.
Definition: SimChannel.h:139
Signal from collection planes.
Definition: geo_types.h:148
void lar::example::AnalysisExample::beginJob ( )
overridevirtual

Reimplemented from art::EDAnalyzer.

Definition at line 352 of file AnalysisExample_module.cc.

References fEndPE, fEndXYZT, fEvent, fGeometryService, fMomentumHist, fPDGCodeHist, fRecodEdxBins, fRecoNdEdxBins, fReconstructionNtuple, fRecoPDG, fRecoTrackID, fRun, fSimdEdxBins, fSimNdEdxBins, fSimPDG, fSimTrackID, fSimulationNtuple, fStartPE, fStartXYZT, fSubRun, and fTrackLengthHist.

353  {
354  // Get the detector length, to determine the maximum bin edge of one of the
355  // histograms.
356  const double detectorLength = DetectorDiagonal(*fGeometryService);
357 
358  // Access art's TFileService, which will handle creating and writing histograms and
359  // n-tuples for us.
361 
362  // For TFileService, the arguments to 'make<whatever>' are the same as those passed to
363  // the 'whatever' constructor, provided 'whatever' is a ROOT class that TFileService
364  // recognizes.
365 
366  // Define the histograms. Putting semi-colons around the title causes it to be
367  // displayed as the x-axis label if the histogram is drawn (the format is "title;label
368  // on abscissae;label on ordinates").
369  fPDGCodeHist = tfs->make<TH1D>("pdgcodes", ";PDG Code;", 5000, -2500, 2500);
370  fMomentumHist = tfs->make<TH1D>("mom", ";particle Momentum (GeV);", 100, 0., 10.);
372  tfs->make<TH1D>("length", ";particle track length (cm);", 200, 0, detectorLength);
373 
374  // Define our n-tuples, which are limited forms of ROOT TTrees. Start with the TTree
375  // itself.
376  fSimulationNtuple = tfs->make<TTree>("AnalysisExampleSimulation", "AnalysisExampleSimulation");
378  tfs->make<TTree>("AnalysisExampleReconstruction", "AnalysisExampleReconstruction");
379 
380  // Define the branches (columns) of our simulation n-tuple. To write a variable, we
381  // give the address of the variable to TTree::Branch.
382  fSimulationNtuple->Branch("Event", &fEvent, "Event/I");
383  fSimulationNtuple->Branch("SubRun", &fSubRun, "SubRun/I");
384  fSimulationNtuple->Branch("Run", &fRun, "Run/I");
385  fSimulationNtuple->Branch("TrackID", &fSimTrackID, "TrackID/I");
386  fSimulationNtuple->Branch("PDG", &fSimPDG, "PDG/I");
387  // When we write arrays, we give the address of the array to TTree::Branch; in C++
388  // this is simply the array name.
389  fSimulationNtuple->Branch("StartXYZT", fStartXYZT, "StartXYZT[4]/D");
390  fSimulationNtuple->Branch("EndXYZT", fEndXYZT, "EndXYZT[4]/D");
391  fSimulationNtuple->Branch("StartPE", fStartPE, "StartPE[4]/D");
392  fSimulationNtuple->Branch("EndPE", fEndPE, "EndPE[4]/D");
393  // For a variable-length array: include the number of bins.
394  fSimulationNtuple->Branch("NdEdx", &fSimNdEdxBins, "NdEdx/I");
395  // ROOT branches can contain std::vector objects.
396  fSimulationNtuple->Branch("dEdx", &fSimdEdxBins);
397 
398  // A similar definition for the reconstruction n-tuple. Note that we use some of the
399  // same variables in both n-tuples.
400  fReconstructionNtuple->Branch("Event", &fEvent, "Event/I");
401  fReconstructionNtuple->Branch("SubRun", &fSubRun, "SubRun/I");
402  fReconstructionNtuple->Branch("Run", &fRun, "Run/I");
403  fReconstructionNtuple->Branch("TrackID", &fRecoTrackID, "TrackID/I");
404  fReconstructionNtuple->Branch("PDG", &fRecoPDG, "PDG/I");
405  fReconstructionNtuple->Branch("NdEdx", &fRecoNdEdxBins, "NdEdx/I");
406  fReconstructionNtuple->Branch("dEdx", &fRecodEdxBins);
407  }
int fSimNdEdxBins
Number of dE/dx bins in a given track.
TH1D * fPDGCodeHist
PDG code of all particles.
int fRecoTrackID
GEANT ID of the particle being processed.
double fEndXYZT[4]
(x,y,z,t) of the true end of the particle
geo::GeometryCore const * fGeometryService
pointer to Geometry provider
std::vector< double > fSimdEdxBins
The vector that will be used to accumulate dE/dx values as a function of range.
double fStartPE[4]
(Px,Py,Pz,E) at the true start of the particle
int fSimTrackID
GEANT ID of the particle being processed.
TTree * fReconstructionNtuple
tuple with reconstructed data
int fRecoNdEdxBins
Number of dE/dx bins in a given track.
double fStartXYZT[4]
(x,y,z,t) of the true start of the particle
std::vector< double > fRecodEdxBins
The vector that will be used to accumulate dE/dx values as a function of range.
int fRecoPDG
PDG ID of the particle being processed.
int fSimPDG
PDG ID of the particle being processed.
int fRun
number of the run being processed
TH1D * fMomentumHist
momentum [GeV] of all selected particles
int fEvent
number of the event being processed
double fEndPE[4]
(Px,Py,Pz,E) at the true end of the particle
TTree * fSimulationNtuple
tuple with simulated data
TH1D * fTrackLengthHist
true length [cm] of all selected particles
void lar::example::AnalysisExample::beginRun ( const art::Run run)
override

Definition at line 410 of file AnalysisExample_module.cc.

References fElectronsToGeV, and sim::LArG4Parameters::GeVToElectrons().

411  {
412  // How to convert from number of electrons to GeV. The ultimate source of this
413  // conversion factor is
414  // ${LARCOREOBJ_INC}/larcoreobj/SimpleTypesAndConstants/PhysicalConstants.h. But
415  // sim::LArG4Parameters might in principle ask a database for it.
417  fElectronsToGeV = 1. / larParameters->GeVToElectrons();
418  }
double fElectronsToGeV
conversion factor
double GeVToElectrons() const
template<typename T , BranchType BT>
ProductToken< T > art::ModuleBase::consumes ( InputTag const &  tag)
protectedinherited

Definition at line 61 of file ModuleBase.h.

References art::ModuleBase::collector_, and art::ConsumesCollector::consumes().

62  {
63  return collector_.consumes<T, BT>(tag);
64  }
ConsumesCollector collector_
Definition: ModuleBase.h:56
ProductToken< T > consumes(InputTag const &)
ConsumesCollector & art::ModuleBase::consumesCollector ( )
protectedinherited

Definition at line 57 of file ModuleBase.cc.

References art::ModuleBase::collector_.

58  {
59  return collector_;
60  }
ConsumesCollector collector_
Definition: ModuleBase.h:56
template<typename T , BranchType BT>
void art::ModuleBase::consumesMany ( )
protectedinherited

Definition at line 75 of file ModuleBase.h.

References art::ModuleBase::collector_, and art::ConsumesCollector::consumesMany().

76  {
77  collector_.consumesMany<T, BT>();
78  }
ConsumesCollector collector_
Definition: ModuleBase.h:56
template<typename Element , BranchType = InEvent>
ViewToken<Element> art::ModuleBase::consumesView ( InputTag const &  )
protectedinherited
template<typename T , BranchType BT>
ViewToken<T> art::ModuleBase::consumesView ( InputTag const &  tag)
inherited

Definition at line 68 of file ModuleBase.h.

References art::ModuleBase::collector_, and art::ConsumesCollector::consumesView().

69  {
70  return collector_.consumesView<T, BT>(tag);
71  }
ConsumesCollector collector_
Definition: ModuleBase.h:56
ViewToken< Element > consumesView(InputTag const &)
void art::detail::Analyzer::doBeginJob ( SharedResources const &  resources)
inherited

Definition at line 25 of file Analyzer.cc.

Referenced by art::detail::Analyzer::Analyzer().

26  {
27  setupQueues(resources);
28  ProcessingFrame const frame{ScheduleID{}};
29  beginJobWithFrame(frame);
30  }
virtual void beginJobWithFrame(ProcessingFrame const &)=0
virtual void setupQueues(SharedResources const &)=0
bool art::detail::Analyzer::doBeginRun ( RunPrincipal rp,
ModuleContext const &  mc 
)
inherited

Definition at line 68 of file Analyzer.cc.

References art::ModuleContext::scheduleID().

Referenced by art::detail::Analyzer::Analyzer().

69  {
70  ProcessingFrame const frame{mc.scheduleID()};
71  beginRunWithFrame(std::as_const(rp).makeRun(mc), frame);
72  return true;
73  }
virtual void beginRunWithFrame(Run const &, ProcessingFrame const &)=0
bool art::detail::Analyzer::doBeginSubRun ( SubRunPrincipal srp,
ModuleContext const &  mc 
)
inherited

Definition at line 84 of file Analyzer.cc.

References art::ModuleContext::scheduleID().

Referenced by art::detail::Analyzer::Analyzer().

85  {
86  ProcessingFrame const frame{mc.scheduleID()};
87  beginSubRunWithFrame(std::as_const(srp).makeSubRun(mc), frame);
88  return true;
89  }
virtual void beginSubRunWithFrame(SubRun const &, ProcessingFrame const &)=0
void art::detail::Analyzer::doEndJob ( )
inherited

Definition at line 33 of file Analyzer.cc.

Referenced by art::detail::Analyzer::Analyzer().

34  {
35  ProcessingFrame const frame{ScheduleID{}};
36  endJobWithFrame(frame);
37  }
virtual void endJobWithFrame(ProcessingFrame const &)=0
bool art::detail::Analyzer::doEndRun ( RunPrincipal rp,
ModuleContext const &  mc 
)
inherited

Definition at line 76 of file Analyzer.cc.

References art::ModuleContext::scheduleID().

Referenced by art::detail::Analyzer::Analyzer().

77  {
78  ProcessingFrame const frame{mc.scheduleID()};
79  endRunWithFrame(std::as_const(rp).makeRun(mc), frame);
80  return true;
81  }
virtual void endRunWithFrame(Run const &, ProcessingFrame const &)=0
bool art::detail::Analyzer::doEndSubRun ( SubRunPrincipal srp,
ModuleContext const &  mc 
)
inherited

Definition at line 92 of file Analyzer.cc.

References art::ModuleContext::scheduleID().

Referenced by art::detail::Analyzer::Analyzer().

93  {
94  ProcessingFrame const frame{mc.scheduleID()};
95  endSubRunWithFrame(std::as_const(srp).makeSubRun(mc), frame);
96  return true;
97  }
virtual void endSubRunWithFrame(SubRun const &, ProcessingFrame const &)=0
bool art::detail::Analyzer::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 
)
inherited

Definition at line 100 of file Analyzer.cc.

References e, and art::ModuleContext::scheduleID().

Referenced by art::detail::Analyzer::Analyzer().

105  {
106  auto const e = std::as_const(ep).makeEvent(mc);
107  if (wantEvent(mc.scheduleID(), e)) {
108  ++counts_run;
109  ProcessingFrame const frame{mc.scheduleID()};
110  analyzeWithFrame(e, frame);
111  ++counts_passed;
112  }
113  return true;
114  }
bool wantEvent(ScheduleID id, Event const &e) const
Definition: Observer.cc:63
Float_t e
Definition: plot.C:35
virtual void analyzeWithFrame(Event const &, ProcessingFrame const &)=0
void art::detail::Analyzer::doRespondToCloseInputFile ( FileBlock const &  fb)
inherited

Definition at line 47 of file Analyzer.cc.

Referenced by art::detail::Analyzer::Analyzer().

48  {
49  ProcessingFrame const frame{ScheduleID{}};
51  }
TFile fb("Li6.root")
virtual void respondToCloseInputFileWithFrame(FileBlock const &, ProcessingFrame const &)=0
void art::detail::Analyzer::doRespondToCloseOutputFiles ( FileBlock const &  fb)
inherited

Definition at line 61 of file Analyzer.cc.

Referenced by art::detail::Analyzer::Analyzer().

62  {
63  ProcessingFrame const frame{ScheduleID{}};
65  }
virtual void respondToCloseOutputFilesWithFrame(FileBlock const &, ProcessingFrame const &)=0
TFile fb("Li6.root")
void art::detail::Analyzer::doRespondToOpenInputFile ( FileBlock const &  fb)
inherited

Definition at line 40 of file Analyzer.cc.

Referenced by art::detail::Analyzer::Analyzer().

41  {
42  ProcessingFrame const frame{ScheduleID{}};
44  }
TFile fb("Li6.root")
virtual void respondToOpenInputFileWithFrame(FileBlock const &, ProcessingFrame const &)=0
void art::detail::Analyzer::doRespondToOpenOutputFiles ( FileBlock const &  fb)
inherited

Definition at line 54 of file Analyzer.cc.

Referenced by art::detail::Analyzer::Analyzer().

55  {
56  ProcessingFrame const frame{ScheduleID{}};
58  }
virtual void respondToOpenOutputFilesWithFrame(FileBlock const &, ProcessingFrame const &)=0
TFile fb("Li6.root")
std::array< std::vector< ProductInfo >, NumBranchTypes > const & art::ModuleBase::getConsumables ( ) const
inherited

Definition at line 43 of file ModuleBase.cc.

References art::ModuleBase::collector_, and art::ConsumesCollector::getConsumables().

44  {
45  return collector_.getConsumables();
46  }
ConsumesCollector collector_
Definition: ModuleBase.h:56
std::array< std::vector< ProductInfo >, NumBranchTypes > const & getConsumables() const
Handle< TriggerResults > art::Observer::getTriggerResults ( Event const &  e) const
protectedinherited

Definition at line 75 of file Observer.cc.

References art::ProductRetriever::get(), and art::Observer::selectors_.

Referenced by art::OutputModule::doWriteEvent(), and art::Observer::wantAllEvents().

76  {
77  if (selectors_) {
78  return selectors_->getOneTriggerResults(e);
79  }
80 
81  // The following applies for cases where no SelectEvents entries
82  // exist.
83  Handle<TriggerResults> h;
84  if (e.get(empty_process_name, h)) {
85  return h;
86  }
87  return Handle<TriggerResults>{};
88  }
Float_t e
Definition: plot.C:35
std::optional< detail::ProcessAndEventSelectors > selectors_
Definition: Observer.h:79
std::unique_ptr< Worker > art::ModuleBase::makeWorker ( WorkerParams const &  wp)
inherited

Definition at line 37 of file ModuleBase.cc.

References art::ModuleBase::doMakeWorker(), and art::NumBranchTypes.

38  {
39  return doMakeWorker(wp);
40  }
virtual std::unique_ptr< Worker > doMakeWorker(WorkerParams const &wp)=0
template<typename T , BranchType BT>
ProductToken< T > art::ModuleBase::mayConsume ( InputTag const &  tag)
protectedinherited

Definition at line 82 of file ModuleBase.h.

References art::ModuleBase::collector_, and art::ConsumesCollector::mayConsume().

83  {
84  return collector_.mayConsume<T, BT>(tag);
85  }
ProductToken< T > mayConsume(InputTag const &)
ConsumesCollector collector_
Definition: ModuleBase.h:56
template<typename T , BranchType BT>
void art::ModuleBase::mayConsumeMany ( )
protectedinherited

Definition at line 96 of file ModuleBase.h.

References art::ModuleBase::collector_, and art::ConsumesCollector::mayConsumeMany().

97  {
98  collector_.mayConsumeMany<T, BT>();
99  }
ConsumesCollector collector_
Definition: ModuleBase.h:56
template<typename Element , BranchType = InEvent>
ViewToken<Element> art::ModuleBase::mayConsumeView ( InputTag const &  )
protectedinherited
template<typename T , BranchType BT>
ViewToken<T> art::ModuleBase::mayConsumeView ( InputTag const &  tag)
inherited

Definition at line 89 of file ModuleBase.h.

References art::ModuleBase::collector_, and art::ConsumesCollector::mayConsumeView().

90  {
91  return collector_.mayConsumeView<T, BT>(tag);
92  }
ConsumesCollector collector_
Definition: ModuleBase.h:56
ViewToken< Element > mayConsumeView(InputTag const &)
ModuleDescription const & art::ModuleBase::moduleDescription ( ) const
inherited

Definition at line 13 of file ModuleBase.cc.

References art::errors::LogicError.

Referenced by art::OutputModule::doRespondToOpenInputFile(), art::OutputModule::doWriteEvent(), art::Modifier::fillProductDescriptions(), art::OutputModule::makePlugins_(), art::OutputWorker::OutputWorker(), reco::shower::LArPandoraModularShowerCreation::produce(), art::Modifier::registerProducts(), and art::OutputModule::registerProducts().

14  {
15  if (md_.has_value()) {
16  return *md_;
17  }
18 
20  "There was an error while calling moduleDescription().\n"}
21  << "The moduleDescription() base-class member function cannot be called\n"
22  "during module construction. To determine which module is "
23  "responsible\n"
24  "for calling it, find the '<module type>:<module "
25  "label>@Construction'\n"
26  "tag in the message prefix above. Please contact artists@fnal.gov\n"
27  "for guidance.\n";
28  }
cet::coded_exception< errors::ErrorCodes, ExceptionDetail::translate > Exception
Definition: Exception.h:66
std::optional< ModuleDescription > md_
Definition: ModuleBase.h:55
string const & art::Observer::processName ( ) const
protectedinherited

Definition at line 57 of file Observer.cc.

References art::Observer::process_name_.

Referenced by art::FileDumperOutput::printPrincipal().

58  {
59  return process_name_;
60  }
std::string process_name_
Definition: Observer.h:76
void art::ModuleBase::setModuleDescription ( ModuleDescription const &  md)
inherited

Definition at line 31 of file ModuleBase.cc.

References art::ModuleBase::md_.

32  {
33  md_ = md;
34  }
std::optional< ModuleDescription > md_
Definition: ModuleBase.h:55
void art::ModuleBase::sortConsumables ( std::string const &  current_process_name)
inherited

Definition at line 49 of file ModuleBase.cc.

References art::ModuleBase::collector_, and art::ConsumesCollector::sortConsumables().

50  {
51  // Now that we know we have seen all the consumes declarations,
52  // sort the results for fast lookup later.
53  collector_.sortConsumables(current_process_name);
54  }
ConsumesCollector collector_
Definition: ModuleBase.h:56
void sortConsumables(std::string const &current_process_name)
bool art::Observer::wantAllEvents ( ) const
inlineprotectednoexceptinherited

Definition at line 31 of file Observer.h.

References e, art::Observer::getTriggerResults(), art::Observer::wantAllEvents_, and art::Observer::wantEvent().

32  {
33  return wantAllEvents_;
34  }
bool wantAllEvents_
Definition: Observer.h:75
bool art::Observer::wantEvent ( ScheduleID  id,
Event const &  e 
) const
protectedinherited

Definition at line 63 of file Observer.cc.

References art::Observer::rejectors_, art::Observer::selectors_, and art::Observer::wantAllEvents_.

Referenced by art::OutputModule::doEvent(), art::OutputModule::doWriteEvent(), and art::Observer::wantAllEvents().

64  {
65  if (wantAllEvents_) {
66  return true;
67  }
68  bool const select_event = selectors_ ? selectors_->matchEvent(id, e) : true;
69  bool const reject_event =
70  rejectors_ ? rejectors_->matchEvent(id, e) : false;
71  return select_event and not reject_event;
72  }
bool wantAllEvents_
Definition: Observer.h:75
std::optional< detail::ProcessAndEventSelectors > rejectors_
Definition: Observer.h:80
Float_t e
Definition: plot.C:35
std::optional< detail::ProcessAndEventSelectors > selectors_
Definition: Observer.h:79

Member Data Documentation

double lar::example::AnalysisExample::fBinSize
private

For dE/dx work: the value of dx.

Definition at line 250 of file AnalysisExample_module.cc.

Referenced by analyze().

art::InputTag lar::example::AnalysisExample::fClusterProducerLabel
private

The name of the producer that created clusters

Definition at line 247 of file AnalysisExample_module.cc.

Referenced by AnalysisExample(), and analyze().

double lar::example::AnalysisExample::fElectronsToGeV
private

conversion factor

Definition at line 304 of file AnalysisExample_module.cc.

Referenced by analyze(), and beginRun().

double lar::example::AnalysisExample::fEndPE[4]
private

(Px,Py,Pz,E) at the true end of the particle

Definition at line 280 of file AnalysisExample_module.cc.

Referenced by analyze(), and beginJob().

double lar::example::AnalysisExample::fEndXYZT[4]
private

(x,y,z,t) of the true end of the particle

Definition at line 278 of file AnalysisExample_module.cc.

Referenced by analyze(), and beginJob().

int lar::example::AnalysisExample::fEvent
private

number of the event being processed

Definition at line 264 of file AnalysisExample_module.cc.

Referenced by analyze(), and beginJob().

geo::GeometryCore const* lar::example::AnalysisExample::fGeometryService
private

pointer to Geometry provider

Definition at line 303 of file AnalysisExample_module.cc.

Referenced by AnalysisExample(), and beginJob().

art::InputTag lar::example::AnalysisExample::fHitProducerLabel
private

The name of the producer that created hits.

Definition at line 246 of file AnalysisExample_module.cc.

Referenced by AnalysisExample(), and analyze().

TH1D* lar::example::AnalysisExample::fMomentumHist
private

momentum [GeV] of all selected particles

Definition at line 254 of file AnalysisExample_module.cc.

Referenced by analyze(), and beginJob().

TH1D* lar::example::AnalysisExample::fPDGCodeHist
private

PDG code of all particles.

Definition at line 253 of file AnalysisExample_module.cc.

Referenced by analyze(), and beginJob().

std::vector<double> lar::example::AnalysisExample::fRecodEdxBins
private

The vector that will be used to accumulate dE/dx values as a function of range.

Definition at line 298 of file AnalysisExample_module.cc.

Referenced by analyze(), and beginJob().

int lar::example::AnalysisExample::fRecoNdEdxBins
private

Number of dE/dx bins in a given track.

Definition at line 295 of file AnalysisExample_module.cc.

Referenced by analyze(), and beginJob().

TTree* lar::example::AnalysisExample::fReconstructionNtuple
private

tuple with reconstructed data

Definition at line 259 of file AnalysisExample_module.cc.

Referenced by analyze(), and beginJob().

int lar::example::AnalysisExample::fRecoPDG
private

PDG ID of the particle being processed.

Definition at line 291 of file AnalysisExample_module.cc.

Referenced by analyze(), and beginJob().

int lar::example::AnalysisExample::fRecoTrackID
private

GEANT ID of the particle being processed.

Definition at line 292 of file AnalysisExample_module.cc.

Referenced by analyze(), and beginJob().

int lar::example::AnalysisExample::fRun
private

number of the run being processed

Definition at line 265 of file AnalysisExample_module.cc.

Referenced by analyze(), and beginJob().

int lar::example::AnalysisExample::fSelectedPDG
private

PDG code of particle we'll focus on.

Definition at line 249 of file AnalysisExample_module.cc.

Referenced by analyze().

std::vector<double> lar::example::AnalysisExample::fSimdEdxBins
private

The vector that will be used to accumulate dE/dx values as a function of range.

Definition at line 286 of file AnalysisExample_module.cc.

Referenced by analyze(), and beginJob().

int lar::example::AnalysisExample::fSimNdEdxBins
private

Number of dE/dx bins in a given track.

Definition at line 283 of file AnalysisExample_module.cc.

Referenced by analyze(), and beginJob().

int lar::example::AnalysisExample::fSimPDG
private

PDG ID of the particle being processed.

Definition at line 271 of file AnalysisExample_module.cc.

Referenced by analyze(), and beginJob().

int lar::example::AnalysisExample::fSimTrackID
private

GEANT ID of the particle being processed.

Definition at line 272 of file AnalysisExample_module.cc.

Referenced by analyze(), and beginJob().

TTree* lar::example::AnalysisExample::fSimulationNtuple
private

tuple with simulated data

Definition at line 258 of file AnalysisExample_module.cc.

Referenced by analyze(), and beginJob().

art::InputTag lar::example::AnalysisExample::fSimulationProducerLabel
private

The name of the producer that tracked simulated particles through the detector

Definition at line 244 of file AnalysisExample_module.cc.

Referenced by AnalysisExample(), and analyze().

double lar::example::AnalysisExample::fStartPE[4]
private

(Px,Py,Pz,E) at the true start of the particle

Definition at line 279 of file AnalysisExample_module.cc.

Referenced by analyze(), and beginJob().

double lar::example::AnalysisExample::fStartXYZT[4]
private

(x,y,z,t) of the true start of the particle

Definition at line 277 of file AnalysisExample_module.cc.

Referenced by analyze(), and beginJob().

int lar::example::AnalysisExample::fSubRun
private

number of the sub-run being processed

Definition at line 266 of file AnalysisExample_module.cc.

Referenced by analyze(), and beginJob().

TH1D* lar::example::AnalysisExample::fTrackLengthHist
private

true length [cm] of all selected particles

Definition at line 255 of file AnalysisExample_module.cc.

Referenced by analyze(), and beginJob().

int lar::example::AnalysisExample::fTriggerOffset
private

(units of ticks) time of expected neutrino event

Definition at line 305 of file AnalysisExample_module.cc.

Referenced by AnalysisExample().


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