LArSoft  v09_90_00
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...
 
virtual void beginJob () override
 
virtual void beginRun (const art::Run &run) override
 
virtual 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
 
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
 

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 178 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 374 of file AnalysisExample_module.cc.

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

375  : EDAnalyzer(config)
376  , fSimulationProducerLabel(config().SimulationLabel())
377  , fHitProducerLabel(config().HitLabel())
378  , fClusterProducerLabel(config().ClusterLabel())
379  , fSelectedPDG(config().PDGcode())
380  , fBinSize(config().BinSize())
381  {
382  // Get a pointer to the geometry service provider.
383  fGeometryService = lar::providerFrom<geo::Geometry>();
384 
385  auto const clock_data =
387  fTriggerOffset = trigger_offset(clock_data);
388 
389  // Since art 2.8, you can and should tell beforehand, here in the
390  // constructor, all the data the module is going to read ("consumes") or
391  // might read
392  // ("may_consume"). Diligence here will in the future help the framework
393  // execute modules in parallel, making sure the order is correct.
394  consumes<std::vector<simb::MCParticle>>(fSimulationProducerLabel);
395  consumes<std::vector<sim::SimChannel>>(fSimulationProducerLabel);
396  consumes<art::Assns<simb::MCTruth, simb::MCParticle>>(fSimulationProducerLabel);
397  consumes<std::vector<recob::Hit>>(fHitProducerLabel);
398  consumes<std::vector<recob::Cluster>>(fClusterProducerLabel);
399  consumes<art::Assns<recob::Cluster, recob::Hit>>(fHitProducerLabel);
400  }
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)
overridevirtual

Definition at line 481 of file AnalysisExample_module.cc.

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

482  {
483  // Start by fetching some basic event information for our n-tuple.
484  fEvent = event.id().event();
485  fRun = event.run();
486  fSubRun = event.subRun();
487 
488  // This is the standard method of reading multiple objects
489  // associated with the same event; see
490  // <https://cdcvs.fnal.gov/redmine/projects/larsoft/wiki/Using_art_in_LArSoft>
491  // for more information.
492  //
493  // Define a "handle" to point to a vector of the objects.
495 
496  // Then tell the event to fill the vector with all the objects of
497  // that type produced by a particular producer.
498  //
499  // Note that if I don't test whether this is successful, and there
500  // aren't any simb::MCParticle objects, then the first time we
501  // access particleHandle, art will display a "ProductNotFound"
502  // exception message and, depending on art settings, it may skip
503  // all processing for the rest of this event (including any
504  // subsequent analysis steps) or stop the execution.
505  if (!event.getByLabel(fSimulationProducerLabel, particleHandle)) {
506  // If we have no MCParticles at all in an event, then we're in
507  // big trouble. Throw an exception to force this module to
508  // stop. Try to include enough information for the user to
509  // figure out what's going on. Note that when we throw a
510  // cet::exception, the run and event number will automatically
511  // be displayed.
512  //
513  // __LINE__ and __FILE__ are values computed by the compiler.
514 
515  throw cet::exception("AnalysisExample")
516  << " No simb::MCParticle objects in this event - "
517  << " Line " << __LINE__ << " in file " << __FILE__ << std::endl;
518  }
519 
520  // Get all the simulated channels for the event. These channels
521  // include the energy deposited for each simulated track.
522  //
523  // Here we use a different method to access objects:
524  // art::ValidHandle. Using this method, if there aren't any
525  // objects of the given type (sim::SimChannel in this case) in the
526  // input file for this art::Event, art will throw a
527  // ProductNotFound exception.
528  //
529  // The "auto" type means that the C++ compiler will determine the
530  // appropriate type for us, based on the return type of
531  // art::Event::getValidHandle<T>(). The "auto" keyword is a great
532  // timesaver, especially with frameworks like LArSoft which often
533  // have complicated data product structures.
534 
535  auto simChannelHandle =
536  event.getValidHandle<std::vector<sim::SimChannel>>(fSimulationProducerLabel);
537 
538  //
539  // Let's compute the variables for the simulation n-tuple first.
540  //
541 
542  // The MCParticle objects are not necessarily in any particular
543  // order. Since we may have to search the list of particles, let's
544  // put them into a map, a sorted container that will make
545  // searching fast and easy. To save both space and time, the map
546  // will not contain a copy of the MCParticle, but a pointer to it.
547  std::map<int, const simb::MCParticle*> particleMap;
548 
549  // This is a "range-based for loop" in the 2011 version of C++; do
550  // a web search on "c++ range based for loop" for more
551  // information. Here's how it breaks down:
552 
553  // - A range-based for loop operates on a container.
554  // "particleHandle" is not a container; it's a pointer to a
555  // container. If we want C++ to "see" a container, we have to
556  // dereference the pointer, like this: *particleHandle.
557 
558  // - The loop variable that is set to each object in the container
559  // is named "particle". As for the loop variable's type:
560 
561  // - To save a little bit of typing, and to make the code easier
562  // to maintain, we're going to let the C++ compiler deduce the
563  // type of what's in the container (simb::MCParticle objects
564  // in this case), so we use "auto".
565 
566  // - We do _not_ want to change the contents of the container,
567  // so we use the "const" keyword to make sure.
568 
569  // - We don't need to copy each object from the container into
570  // the variable "particle". It's sufficient to refer to the
571  // object by its address. So we use the reference operator "&"
572  // to tell the compiler to just copy the address, not the
573  // entire object.
574 
575  // It sounds complicated, but the end result is that we loop over
576  // the list of particles in the art::Event in the most efficient
577  // way possible.
578 
579  for (auto const& particle : (*particleHandle)) {
580  // For the methods you can call to get particle information,
581  // see ${NUSIMDATA_INC}/nusimdata/SimulationBase/MCParticle.h.
582  fSimTrackID = particle.TrackId();
583 
584  // Add the address of the MCParticle to the map, with the
585  // track ID as the key.
586  particleMap[fSimTrackID] = &particle;
587 
588  // Histogram the PDG code of every particle in the event.
589  fSimPDG = particle.PdgCode();
590  fPDGCodeHist->Fill(fSimPDG);
591 
592  // For this example, we want to fill the n-tuples and histograms
593  // only with information from the primary particles in the
594  // event, whose PDG codes match a value supplied in the .fcl file.
595  if (particle.Process() != "primary" || fSimPDG != fSelectedPDG) continue;
596 
597  // A particle has a trajectory, consisting of a set of
598  // 4-positions and 4-mommenta.
599  const size_t numberTrajectoryPoints = particle.NumberTrajectoryPoints();
600 
601  // For trajectories, as for vectors and arrays, the first
602  // point is #0, not #1.
603  const int last = numberTrajectoryPoints - 1;
604  const TLorentzVector& positionStart = particle.Position(0);
605  const TLorentzVector& positionEnd = particle.Position(last);
606  const TLorentzVector& momentumStart = particle.Momentum(0);
607  const TLorentzVector& momentumEnd = particle.Momentum(last);
608 
609  // Make a histogram of the starting momentum.
610  fMomentumHist->Fill(momentumStart.P());
611 
612  // Fill arrays with the 4-values. (Don't be fooled by
613  // the name of the method; it just puts the numbers from
614  // the 4-vector into the array.)
615  positionStart.GetXYZT(fStartXYZT);
616  positionEnd.GetXYZT(fEndXYZT);
617  momentumStart.GetXYZT(fStartPE);
618  momentumEnd.GetXYZT(fEndPE);
619 
620  // Use a polar-coordinate view of the 4-vectors to
621  // get the track length.
622  const double trackLength = (positionEnd - positionStart).Rho();
623 
624  // Let's print some debug information in the job output to see
625  // that everything is fine. MF_LOG_DEBUG() is a messagefacility
626  // macro that prints stuff when the message level is set to
627  // standard_debug in the .fcl file.
628  MF_LOG_DEBUG("AnalysisExample") << "Track length: " << trackLength << " cm";
629 
630  // Fill a histogram of the track length.
631  fTrackLengthHist->Fill(trackLength);
632 
633  MF_LOG_DEBUG("AnalysisExample")
634  << "track ID=" << fSimTrackID << " (PDG ID: " << fSimPDG << ") " << trackLength
635  << " cm long, momentum " << momentumStart.P() << " GeV/c, has " << numberTrajectoryPoints
636  << " trajectory points";
637 
638  // Determine the number of dE/dx bins for the n-tuple.
639  fSimNdEdxBins = int(trackLength / fBinSize) + 1;
640  // Initialize the vector of dE/dx bins to be empty.
641  fSimdEdxBins.clear();
642 
643  // To look at the energy deposited by this particle's track,
644  // we loop over the SimChannel objects in the event.
645  for (auto const& channel : (*simChannelHandle)) {
646  // Get the numeric ID associated with this channel. (The
647  // channel number is a 32-bit unsigned int, which normally
648  // requires a special data type. Let's use "auto" so we
649  // don't have to remember "raw::ChannelID_t".
650  auto const channelNumber = channel.Channel();
651 
652  // A little care: There is more than one plane that reacts
653  // to charge in the TPC. We only want to include the
654  // energy from the collection plane. Note:
655  // geo::kCollection is defined in
656  // ${LARCOREOBJ_INC}/larcoreobj/SimpleTypesAndConstants/geo_types.h
657  if (fGeometryService->SignalType(channelNumber) != geo::kCollection) continue;
658 
659  // Each channel has a map inside it that connects a time
660  // slice to energy deposits in the detector. We'll use
661  // "auto", but it's worth noting that the full type of
662  // this map is
663  // std::map<unsigned short, std::vector<sim::IDE>>
664  auto const& timeSlices = channel.TDCIDEMap();
665 
666  // For every time slice in this channel:
667  for (auto const& timeSlice : timeSlices) {
668  // Each entry in a map is a pair<first,second>. For
669  // the timeSlices map, the 'first' is a time slice
670  // number. The 'second' is a vector of IDE objects.
671  auto const& energyDeposits = timeSlice.second;
672 
673  // Loop over the energy deposits. An "energy deposit"
674  // object is something that knows how much
675  // charge/energy was deposited in a small volume, by
676  // which particle, and where. The type of
677  // 'energyDeposit' will be sim::IDE, which is defined
678  // in
679  // ${LARDATAOBJ_INC}/lardataobj/Simulation/SimChannel.h.
680  for (auto const& energyDeposit : energyDeposits) {
681  // Check if the track that deposited the
682  // energy matches the track of the particle.
683  if (energyDeposit.trackID != fSimTrackID) continue;
684  // Get the (x,y,z) of the energy deposit.
685  TVector3 location(energyDeposit.x, energyDeposit.y, energyDeposit.z);
686 
687  // Distance from the start of the track.
688  const double distance = (location - positionStart.Vect()).Mag();
689 
690  // Into which bin of the dE/dx array do we add the energy?
691  const unsigned int bin = (unsigned int)(distance / fBinSize);
692 
693  // Is the dE/dx array big enough to include this bin?
694  if (fSimdEdxBins.size() < bin + 1) {
695  // Increase the array size to accomodate
696  // the new bin, padding it with zeros.
697  fSimdEdxBins.resize(bin + 1, 0.);
698  }
699 
700  // Add the energy deposit to that bin. (If you look at the
701  // definition of sim::IDE, you'll see that there's another
702  // way to get the energy. Are the two methods equivalent?
703  // Compare the results and see!)
704  fSimdEdxBins[bin] += energyDeposit.numElectrons * fElectronsToGeV;
705 
706  } // For each energy deposit
707  } // For each time slice
708  } // For each SimChannel
709 
710  // At this point we've filled in the values of all the
711  // variables and arrays we want to write to the n-tuple
712  // for this particle. The following command actually
713  // writes those values.
714  fSimulationNtuple->Fill();
715 
716  } // loop over all particles in the event.
717 
718  //
719  // Reconstruction n-tuple
720  //
721 
722  // All of the above is based on objects entirely produced by the
723  // simulation. Let's try to do something based on reconstructed
724  // objects. A Hit (see ${LARDATAOBJ_INC}/lardataobj/RecoBase/Hit.h)
725  // is a 2D object in a plane.
726 
727  // This code duplicates much of the code in the previous big
728  // simulation loop, and will produce the similar results. (It
729  // won't be identical, because of shaping and other factors; not
730  // every bit of charge in a channel ends up contributing to a
731  // hit.) The point is to show different methods of accessing
732  // information, not to produce profound physics results -- that
733  // part is up to you!
734 
735  // For the rest of this method, I'm going to assume you've read
736  // the comments in previous section; I won't repeat all the C++
737  // coding tricks and whatnot.
738 
739  // Start by reading the Hits. We don't use art::ValidHandle here
740  // because there might be no hits in the input; e.g., we ran the
741  // simulation but did not run the reconstruction steps yet. If
742  // there are no hits we may as well skip the rest of this module.
743 
745  if (!event.getByLabel(fHitProducerLabel, hitHandle)) return;
746 
747  // Our goal is to accumulate the dE/dx of any particles associated
748  // with the hits that match our criteria: primary particles with
749  // the PDG code from the .fcl file. I don't know how many such
750  // particles there will be in a given event. We'll use a map, with
751  // track ID as the key, to hold the vectors of dE/dx information.
752  std::map<int, std::vector<double>> dEdxMap;
753 
754  auto const clock_data =
756 
757  // For every Hit:
758  for (auto const& hit : (*hitHandle)) {
759  // The channel associated with this hit.
760  auto hitChannelNumber = hit.Channel();
761 
762  // We have a hit. For this example let's just focus on the
763  // hits in the collection plane.
764  if (fGeometryService->SignalType(hitChannelNumber) != geo::kCollection) continue;
765 
766  MF_LOG_DEBUG("AnalysisExample") << "Hit in collection plane" << std::endl;
767 
768  // In a few lines we're going to look for possible energy
769  // deposits that correspond to that hit. Determine a
770  // reasonable range of times that might correspond to those
771  // energy deposits.
772 
773  // In reconstruction, the channel waveforms are truncated. So
774  // we have to adjust the Hit TDC ticks to match those of the
775  // SimChannels, which were created during simulation.
776 
777  // Save a bit of typing, while still allowing for potential
778  // changes in the definitions of types in
779  // $LARDATAOBJ_DIR/source/lardataobj/Simulation/SimChannel.h
780 
781  typedef sim::SimChannel::StoredTDC_t TDC_t;
782  TDC_t start_tdc = clock_data.TPCTick2TDC(hit.StartTick());
783  TDC_t end_tdc = clock_data.TPCTick2TDC(hit.EndTick());
784  TDC_t hitStart_tdc = clock_data.TPCTick2TDC(hit.PeakTime() - 3. * hit.SigmaPeakTime());
785  TDC_t hitEnd_tdc = clock_data.TPCTick2TDC(hit.PeakTime() + 3. * hit.SigmaPeakTime());
786 
787  start_tdc = std::max(start_tdc, hitStart_tdc);
788  end_tdc = std::min(end_tdc, hitEnd_tdc);
789 
790  // In the simulation section, we started with particles to find
791  // channels with a matching track ID. Now we search in reverse:
792  // search the SimChannels for matching channel number, then look
793  // at the tracks inside the channel.
794 
795  for (auto const& channel : (*simChannelHandle)) {
796  auto simChannelNumber = channel.Channel();
797  if (simChannelNumber != hitChannelNumber) continue;
798 
799  MF_LOG_DEBUG("AnalysisExample")
800  << "SimChannel number = " << simChannelNumber << std::endl;
801 
802  // The time slices in this channel.
803  auto const& timeSlices = channel.TDCIDEMap();
804 
805  // We want to look over the range of time slices in this
806  // channel that correspond to the range of hit times.
807 
808  // To do this, we're going to use some fast STL search
809  // methods; STL algorithms are usually faster than the
810  // ones you might write yourself. The price for this speed
811  // is a bit of code complexity: in particular, we need to
812  // a custom comparison function, TDCIDETimeCompare, to
813  // define a "less-than" function for the searches.
814 
815  // For an introduction to STL algorithms, see
816  // <http://www.learncpp.com/cpp-tutorial/16-4-stl-algorithms-overview/>.
817  // For a list of available STL algorithms, see
818  // <http://en.cppreference.com/w/cpp/algorithm>
819 
820  // We have to create "dummy" time slices for the search.
821  sim::TDCIDE startTime;
822  sim::TDCIDE endTime;
823  startTime.first = start_tdc;
824  endTime.first = end_tdc;
825 
826  // Here are the fast searches:
827 
828  // Find a pointer to the first channel with time >= start_tdc.
829  auto const startPointer =
830  std::lower_bound(timeSlices.begin(), timeSlices.end(), startTime, TDCIDETimeCompare);
831 
832  // From that time slice, find the last channel with time < end_tdc.
833  auto const endPointer =
834  std::upper_bound(startPointer, timeSlices.end(), endTime, TDCIDETimeCompare);
835 
836  // Did we find anything? If not, skip.
837  if (startPointer == timeSlices.end() || startPointer == endPointer) continue;
838  MF_LOG_DEBUG("AnalysisExample")
839  << "Time slice start = " << (*startPointer).first << std::endl;
840 
841  // Loop over the channel times we found that match the hit
842  // times.
843  for (auto slicePointer = startPointer; slicePointer != endPointer; slicePointer++) {
844  auto const timeSlice = *slicePointer;
845  auto time = timeSlice.first;
846 
847  // How to debug a problem: Lots of print statements. There are
848  // debuggers such as gdb, but they can be tricky to use with
849  // shared libraries and don't work if you're using software
850  // that was compiled somewhere else (e.g., you're accessing
851  // LArSoft libraries via CVMFS).
852 
853  // The MF_LOG_DEBUG statements below are left over from when I
854  // was trying to solve a problem about hit timing. I could
855  // have deleted them, but decided to leave them to demonsrate
856  // what a typical debugging process looks like.
857 
858  MF_LOG_DEBUG("AnalysisExample")
859  << "Hit index = " << hit.LocalIndex() << " channel number = " << hitChannelNumber
860  << " start TDC tick = " << hit.StartTick() << " end TDC tick = " << hit.EndTick()
861  << " peak TDC tick = " << hit.PeakTime()
862  << " sigma peak time = " << hit.SigmaPeakTime()
863  << " adjusted start TDC tick = " << clock_data.TPCTick2TDC(hit.StartTick())
864  << " adjusted end TDC tick = " << clock_data.TPCTick2TDC(hit.EndTick())
865  << " adjusted peak TDC tick = " << clock_data.TPCTick2TDC(hit.PeakTime())
866  << " adjusted start_tdc = " << start_tdc << " adjusted end_tdc = " << end_tdc
867  << " time = " << time << std::endl;
868 
869  // Loop over the energy deposits.
870  auto const& energyDeposits = timeSlice.second;
871  for (auto const& energyDeposit : energyDeposits) {
872  // Remember that map of MCParticles we created
873  // near the top of this method? Now we can use
874  // it. Search the map for the track ID associated
875  // with this energy deposit. Since a map is
876  // automatically sorted, we can use a fast binary
877  // search method, 'find()'.
878 
879  // By the way, the type of "search" is an iterator
880  // (to be specific, it's an
881  // std::map<int,simb::MCParticle*>::const_iterator,
882  // which makes you thankful for the "auto"
883  // keyword). If you're going to work with C++
884  // containers, you'll have to learn about
885  // iterators eventually; do a web search on "STL
886  // iterator" to get started.
887  auto search = particleMap.find(energyDeposit.trackID);
888 
889  // Did we find this track ID in the particle map?
890  // It's possible for the answer to be "no"; some
891  // particles are too low in kinetic energy to be
892  // written by the simulation (see
893  // ${LARSIM_DIR}/job/simulationservices.fcl,
894  // parameter ParticleKineticEnergyCut).
895  if (search == particleMap.end()) continue;
896 
897  // "search" points to a pair in the map: <track ID, MCParticle*>
898  int trackID = (*search).first;
899  const simb::MCParticle& particle = *((*search).second);
900 
901  // Is this a primary particle, with a PDG code that
902  // matches the user input?
903  if (particle.Process() != "primary" || particle.PdgCode() != fSelectedPDG) continue;
904 
905  // Determine the dE/dx of this particle.
906  const TLorentzVector& positionStart = particle.Position(0);
907  TVector3 location(energyDeposit.x, energyDeposit.y, energyDeposit.z);
908  double distance = (location - positionStart.Vect()).Mag();
909  unsigned int bin = int(distance / fBinSize);
910  double energy = energyDeposit.numElectrons * fElectronsToGeV;
911 
912  // A feature of maps: if we refer to
913  // dEdxMap[trackID], and there's no such entry in
914  // the map yet, it will be automatically created
915  // with a zero-size vector. Test to see if the
916  // vector for this track ID is big enough.
917  //
918  // dEdxMap is a map, which is a slow container
919  // compared to a vector. If we are going to access
920  // the same element over and over, it is a good
921  // idea to find that element once, and then refer
922  // to that item directly. Since we don't care
923  // about the specific type of dEdxMap[trackID] (a
924  // vector, by the way), we can use "auto" to save
925  // some time. This must be a reference, since we
926  // want to change the original value in the map,
927  // and can't be constant.
928  auto& track_dEdX = dEdxMap[trackID];
929  if (track_dEdX.size() < bin + 1) {
930  // Increase the vector size, padding it with
931  // zeroes.
932  track_dEdX.resize(bin + 1, 0);
933  }
934 
935  // Add the energy to the dE/dx bins for this track.
936  track_dEdX[bin] += energy;
937 
938  } // loop over energy deposits
939  } // loop over time slices
940  } // for each SimChannel
941  } // for each Hit
942 
943  // We have a map of dE/dx vectors. Write each one of them to the
944  // reconstruction n-tuple.
945  for (const auto& dEdxEntry : dEdxMap) {
946  // Here, the map entries are <first,second>=<track ID, dE/dx vector>
947  fRecoTrackID = dEdxEntry.first;
948 
949  // This is an example of how we'd pick out the PDG code if
950  // there are multiple particle types or tracks in a single
951  // event allowed in the n-tuple.
952  fRecoPDG = particleMap[fRecoTrackID]->PdgCode();
953 
954  // Get the number of bins for this track.
955  const std::vector<double>& dEdx = dEdxEntry.second;
956  fRecoNdEdxBins = dEdx.size();
957 
958  // Copy this track's dE/dx information.
960 
961  // At this point, we've filled in all the reconstruction
962  // n-tuple's variables. Write it.
963  fReconstructionNtuple->Fill();
964  }
965 
966  // Think about the two big loops above, One starts from the
967  // particles then looks at the channels; the other starts with the
968  // hits and backtracks to the particles. What links the
969  // information in simb::MCParticle and sim::SimChannel is the
970  // track ID number assigned by the LArG4 simulation; what links
971  // sim::SimChannel and recob::Hit is the channel ID.
972 
973  // In general, that is not how objects in the LArSoft
974  // reconstruction chain are linked together. Most of them are
975  // linked using associations and the art::Assns class:
976  // <https://cdcvs.fnal.gov/redmine/projects/larsoft/wiki/Using_art_in_LArSoft#artAssns>
977 
978  // The web page may a bit difficult to understand (at least, it is
979  // for me), so let's try to simplify things:
980 
981  // - If you're doing an analysis task, you don't have to worry
982  // about creating art::Assns objects.
983 
984  // - You don't have to read the art:Assns objects on your
985  // own. There are helper classes (FindXXX) which will do that for
986  // you.
987 
988  // - There's only one helper you need: art::FindManyP. It will do
989  // what you want with a minimum of fuss.
990 
991  // - Associations are symmetric. If you see an
992  // art::Assns<ClassA,ClassB>, the helper classes will find all of
993  // ClassA associated with ClassB or all of ClassB associated with
994  // ClassA.
995 
996  // - To know what associations exist, you have to be a 'code
997  // detective'. The important clue is to look for a 'produces' line
998  // in the code that writes objects to an art::Event. For example,
999  // in ${LARSIM_DIR}/source/larsim/LArG4/LArG4_module.cc, you'll see this
1000  // line:
1001 
1002  // produces< art::Assns<simb::MCTruth, simb::MCParticle> >();
1003 
1004  // That means a simulated event will have an association between
1005  // simb::MCTruth (the primary particle produced by the event
1006  // generator) and the simb::MCParticle objects (the secondary
1007  // particles produced in the detector simulation).
1008 
1009  // Let's try it. The following statement will find the
1010  // simb::MCTruth objects associated with the simb::MCParticle
1011  // objects in the event (referenced by particleHandle):
1012 
1013  const art::FindManyP<simb::MCTruth> findManyTruth(
1014  particleHandle, event, fSimulationProducerLabel);
1015 
1016  // Note that we still have to supply the module label of the step
1017  // that created the association. Also note that we did not have to
1018  // explicitly read in the simb::MCTruth objects from the
1019  // art::Event object 'event'; FindManyP did that for us.
1020 
1021  // Also note that at this point art::FindManyP has already found
1022  // all the simb::MCTruth associated with each of the particles in
1023  // particleHandle. This is a slow process, so in general you want
1024  // to do it only once. If we had a loop over the particles, we
1025  // would still do this outside that loop.
1026 
1027  // Now we can query the 'findManyTruth' object to access the
1028  // information. First, check that there wasn't some kind of error:
1029 
1030  if (!findManyTruth.isValid()) {
1031  mf::LogError("AnalysisExample")
1032  << "findManyTruth simb::MCTruth for simb::MCParticle failed!";
1033  }
1034 
1035  // I'm going to be lazy, and just look at the simb::MCTruth object
1036  // associated with the first simb::MCParticle we read in. (The
1037  // main reason I'm being lazy is that if I used the
1038  // single-particle generator in prodsingle.fcl, every particle in
1039  // the event is going to be associated with just the one primary
1040  // particle from the event generator.)
1041 
1042  size_t particle_index = 0; // look at first particle in
1043  // particleHandle's vector.
1044 
1045  // I'm using "auto" to save on typing. The result of
1046  // FindManyP::at() is a vector of pointers, in this case
1047  // simb::MCTruth*. In this case it will be a vector with just one
1048  // entry; I could have used art::FindOneP instead. (This will be a
1049  // vector of art::Ptr, which is a type of smart pointer; see
1050  // <https://cdcvs.fnal.gov/redmine/projects/larsoft/wiki/Using_art_in_LArSoft#artPtrltTgt-and-artPtrVectorltTgt>
1051  // To avoid unnecessary copying, and since art::FindManyP returns
1052  // a constant reference, use "auto const&".
1053 
1054  auto const& truth = findManyTruth.at(particle_index);
1055 
1056  // Make sure there's no problem.
1057  if (truth.empty()) {
1058  mf::LogError("AnalysisExample")
1059  << "Particle ID=" << particleHandle->at(particle_index).TrackId() << " has no primary!";
1060  }
1061 
1062  // Use the message facility to write output. I don't have to write
1063  // the event, run, or subrun numbers; the message facility takes
1064  // care of that automatically. I'm "going at warp speed" with the
1065  // vectors, pointers, and methods; see
1066  // ${NUSIMDATA_INC}/nusimdata/SimulationBase/MCTruth.h and
1067  // ${NUSIMDATA_INC}/nusimdata/SimulationBase/MCParticle.h for the
1068  // nested calls I'm using.
1069  mf::LogInfo("AnalysisExample")
1070  << "Particle ID=" << particleHandle->at(particle_index).TrackId()
1071  << " primary PDG code=" << truth[0]->GetParticle(0).PdgCode();
1072 
1073  // Let's try a slightly more realistic example. Suppose I want to
1074  // read in the clusters, and learn what hits are associated with
1075  // them. Then I could backtrack from those hits to determine the
1076  // dE/dx of the particles in the clusters. (Don't worry; I won't
1077  // put you through the n-tuple creation procedure for a third
1078  // time.)
1079 
1080  // First, read in the clusters (if there are any).
1082  if (!event.getByLabel(fClusterProducerLabel, clusterHandle)) return;
1083 
1084  // Now use the associations to find which hits are associated with
1085  // which clusters. Note that this is not as trivial a query as the
1086  // one with MCTruth, since it's possible for a hit to be assigned
1087  // to more than one cluster.
1088 
1089  // We have to include fClusterProducerLabel, since that's the step
1090  // that created the art::Assns<recob::Hit,recob::Cluster> object;
1091  // look at the modules in ${LARRECO_DIR}/source/larreco/ClusterFinder/
1092  // and search for the 'produces' lines. (I did not know this before I
1093  // wrote these comments. I had to be a code detective and use UNIX
1094  // tools like 'grep' and 'find' to locate those routines.)
1095  const art::FindManyP<recob::Hit> findManyHits(clusterHandle, event, fClusterProducerLabel);
1096 
1097  if (!findManyHits.isValid()) {
1098  mf::LogError("AnalysisExample") << "findManyHits recob::Hit for recob::Cluster failed;"
1099  << " cluster label='" << fClusterProducerLabel << "'";
1100  }
1101 
1102  // Now we'll loop over the clusters to see the hits associated
1103  // with each one. Note that we're not using a range-based for
1104  // loop. That's because FindManyP::at() expects a number as an
1105  // argument, so we might as well go through the cluster objects
1106  // via numerical index instead.
1107  for (size_t cluster_index = 0; cluster_index != clusterHandle->size(); cluster_index++) {
1108  // In this case, FindManyP::at() will return a vector of
1109  // pointers to recob::Hit that corresponds to the
1110  // "cluster_index"-th cluster.
1111  auto const& hits = findManyHits.at(cluster_index);
1112 
1113  // We have a vector of pointers to the hits associated
1114  // with the cluster, but for this example I'm not going to
1115  // do anything fancy with them. I'll just print out how
1116  // many there are.
1117  mf::LogInfo("AnalysisExample") << "Cluster ID=" << clusterHandle->at(cluster_index).ID()
1118  << " has " << hits.size() << " hits";
1119  }
1120 
1121  } // 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
geo::GeometryCore const * fGeometryService
pointer to Geometry provider
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
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
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:2675
Detector simulation of raw signals on wires.
int fSimPDG
PDG ID of the particle being processed.
int fRun
number of the run being processed
SigType_t SignalType(PlaneID const &pid) const
Returns the type of signal on the channels of specified TPC plane.
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:152
void lar::example::AnalysisExample::beginJob ( )
overridevirtual

Reimplemented from art::EDAnalyzer.

Definition at line 403 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.

404  {
405  // Get the detector length, to determine the maximum bin edge of one
406  // of the histograms.
407  const double detectorLength = DetectorDiagonal(*fGeometryService);
408 
409  // Access ART's TFileService, which will handle creating and writing
410  // histograms and n-tuples for us.
412 
413  // For TFileService, the arguments to 'make<whatever>' are the
414  // same as those passed to the 'whatever' constructor, provided
415  // 'whatever' is a ROOT class that TFileService recognizes.
416 
417  // Define the histograms. Putting semi-colons around the title
418  // causes it to be displayed as the x-axis label if the histogram
419  // is drawn (the format is "title;label on abscissae;label on ordinates").
420  fPDGCodeHist = tfs->make<TH1D>("pdgcodes", ";PDG Code;", 5000, -2500, 2500);
421  fMomentumHist = tfs->make<TH1D>("mom", ";particle Momentum (GeV);", 100, 0., 10.);
423  tfs->make<TH1D>("length", ";particle track length (cm);", 200, 0, detectorLength);
424 
425  // Define our n-tuples, which are limited forms of ROOT
426  // TTrees. Start with the TTree itself.
428  tfs->make<TTree>("AnalysisExampleSimulation", "AnalysisExampleSimulation");
430  tfs->make<TTree>("AnalysisExampleReconstruction", "AnalysisExampleReconstruction");
431 
432  // Define the branches (columns) of our simulation n-tuple. To
433  // write a variable, we give the address of the variable to
434  // TTree::Branch.
435  fSimulationNtuple->Branch("Event", &fEvent, "Event/I");
436  fSimulationNtuple->Branch("SubRun", &fSubRun, "SubRun/I");
437  fSimulationNtuple->Branch("Run", &fRun, "Run/I");
438  fSimulationNtuple->Branch("TrackID", &fSimTrackID, "TrackID/I");
439  fSimulationNtuple->Branch("PDG", &fSimPDG, "PDG/I");
440  // When we write arrays, we give the address of the array to
441  // TTree::Branch; in C++ this is simply the array name.
442  fSimulationNtuple->Branch("StartXYZT", fStartXYZT, "StartXYZT[4]/D");
443  fSimulationNtuple->Branch("EndXYZT", fEndXYZT, "EndXYZT[4]/D");
444  fSimulationNtuple->Branch("StartPE", fStartPE, "StartPE[4]/D");
445  fSimulationNtuple->Branch("EndPE", fEndPE, "EndPE[4]/D");
446  // For a variable-length array: include the number of bins.
447  fSimulationNtuple->Branch("NdEdx", &fSimNdEdxBins, "NdEdx/I");
448  // ROOT branches can contain std::vector objects.
449  fSimulationNtuple->Branch("dEdx", &fSimdEdxBins);
450 
451  // A similar definition for the reconstruction n-tuple. Note that we
452  // use some of the same variables in both n-tuples.
453  fReconstructionNtuple->Branch("Event", &fEvent, "Event/I");
454  fReconstructionNtuple->Branch("SubRun", &fSubRun, "SubRun/I");
455  fReconstructionNtuple->Branch("Run", &fRun, "Run/I");
456  fReconstructionNtuple->Branch("TrackID", &fRecoTrackID, "TrackID/I");
457  fReconstructionNtuple->Branch("PDG", &fRecoPDG, "PDG/I");
458  fReconstructionNtuple->Branch("NdEdx", &fRecoNdEdxBins, "NdEdx/I");
459  fReconstructionNtuple->Branch("dEdx", &fRecodEdxBins);
460  }
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
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
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)
overridevirtual

Definition at line 470 of file AnalysisExample_module.cc.

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

471  {
472  // How to convert from number of electrons to GeV. The ultimate
473  // source of this conversion factor is
474  // ${LARCOREOBJ_INC}/larcoreobj/SimpleTypesAndConstants/PhysicalConstants.h.
475  // But sim::LArG4Parameters might in principle ask a database for it.
477  fElectronsToGeV = 1. / larParameters->GeVToElectrons();
478  }
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 297 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 294 of file AnalysisExample_module.cc.

Referenced by AnalysisExample(), and analyze().

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

conversion factor

Definition at line 353 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 327 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 325 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 311 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 352 of file AnalysisExample_module.cc.

Referenced by AnalysisExample(), analyze(), and beginJob().

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

The name of the producer that created hits.

Definition at line 293 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 301 of file AnalysisExample_module.cc.

Referenced by analyze(), and beginJob().

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

PDG code of all particles.

Definition at line 300 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 347 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 343 of file AnalysisExample_module.cc.

Referenced by analyze(), and beginJob().

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

tuple with reconstructed data

Definition at line 306 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 339 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 340 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 312 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 296 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 334 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 330 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 318 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 319 of file AnalysisExample_module.cc.

Referenced by analyze(), and beginJob().

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

tuple with simulated data

Definition at line 305 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 291 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 326 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 324 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 313 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 302 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 354 of file AnalysisExample_module.cc.

Referenced by AnalysisExample().


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