LArSoft  v09_90_00
Liquid Argon Software toolkit - https://larsoft.org/
simb::MCTrajectory Class Reference

#include "MCTrajectory.h"

Public Types

typedef std::vector< std::pair< TLorentzVector, TLorentzVector > > list_type
 
typedef list_type::value_type value_type
 
typedef list_type::iterator iterator
 
typedef list_type::const_iterator const_iterator
 
typedef list_type::reverse_iterator reverse_iterator
 
typedef list_type::const_reverse_iterator const_reverse_iterator
 
typedef list_type::size_type size_type
 
typedef list_type::difference_type difference_type
 
typedef std::vector< std::pair< size_t, unsigned char > > ProcessMap
 

Public Member Functions

 MCTrajectory ()
 
 MCTrajectory (const TLorentzVector &vertex, const TLorentzVector &momentum)
 
const TLorentzVector & Position (const size_type) const
 The accessor methods described above. More...
 
const TLorentzVector & Momentum (const size_type) const
 
double X (const size_type i) const
 
double Y (const size_type i) const
 
double Z (const size_type i) const
 
double T (const size_type i) const
 
double Px (const size_type i) const
 
double Py (const size_type i) const
 
double Pz (const size_type i) const
 
double E (const size_type i) const
 
double TotalLength () const
 
iterator begin ()
 
const_iterator begin () const
 
iterator end ()
 
const_iterator end () const
 
reverse_iterator rbegin ()
 
const_reverse_iterator rbegin () const
 
reverse_iterator rend ()
 
const_reverse_iterator rend () const
 
size_type size () const
 
bool empty () const
 
void swap (simb::MCTrajectory &other)
 
void clear ()
 
const value_typeoperator[] (const size_type i) const
 
const value_typeat (const size_type i) const
 
void push_back (value_type const &v)
 
void push_back (TLorentzVector const &p, TLorentzVector const &m)
 
void Add (TLorentzVector const &p, TLorentzVector const &m)
 
void Add (TLorentzVector const &p, TLorentzVector const &m, std::string const &process, bool keepTransportation=false)
 
unsigned char ProcessToKey (std::string const &process) const
 
std::string KeyToProcess (unsigned char const &key) const
 
ProcessMap const & TrajectoryProcesses () const
 
void Sparsify (double margin=.1, bool keep_second_to_last=false)
 

Private Attributes

list_type ftrajectory
 The list of trajectory points. More...
 
ProcessMap fTrajectoryProcess
 

Friends

std::ostream & operator<< (std::ostream &output, const MCTrajectory &)
 

Detailed Description

Definition at line 58 of file MCTrajectory.h.

Member Typedef Documentation

typedef list_type::const_iterator simb::MCTrajectory::const_iterator

Definition at line 66 of file MCTrajectory.h.

typedef list_type::const_reverse_iterator simb::MCTrajectory::const_reverse_iterator

Definition at line 68 of file MCTrajectory.h.

typedef list_type::difference_type simb::MCTrajectory::difference_type

Definition at line 70 of file MCTrajectory.h.

typedef list_type::iterator simb::MCTrajectory::iterator

Definition at line 65 of file MCTrajectory.h.

typedef std::vector< std::pair<TLorentzVector, TLorentzVector> > simb::MCTrajectory::list_type

Some type definitions to make life easier, and to help "hide" the implementation details. (If you're not familiar with STL, you can ignore these definitions.)

Definition at line 63 of file MCTrajectory.h.

typedef std::vector< std::pair<size_t, unsigned char> > simb::MCTrajectory::ProcessMap

Definition at line 71 of file MCTrajectory.h.

typedef list_type::reverse_iterator simb::MCTrajectory::reverse_iterator

Definition at line 67 of file MCTrajectory.h.

typedef list_type::size_type simb::MCTrajectory::size_type

Definition at line 69 of file MCTrajectory.h.

typedef list_type::value_type simb::MCTrajectory::value_type

Definition at line 64 of file MCTrajectory.h.

Constructor & Destructor Documentation

simb::MCTrajectory::MCTrajectory ( )

Standard constructor: Start with initial position and momentum of the particle.

Definition at line 24 of file MCTrajectory.cxx.

25  : ftrajectory()
26  {}
list_type ftrajectory
The list of trajectory points.
Definition: MCTrajectory.h:77
simb::MCTrajectory::MCTrajectory ( const TLorentzVector &  vertex,
const TLorentzVector &  momentum 
)

Definition at line 29 of file MCTrajectory.cxx.

References ftrajectory.

31  {
32  ftrajectory.push_back( value_type( position, momentum ) );
33  }
list_type::value_type value_type
Definition: MCTrajectory.h:64
list_type ftrajectory
The list of trajectory points.
Definition: MCTrajectory.h:77

Member Function Documentation

void simb::MCTrajectory::Add ( TLorentzVector const &  p,
TLorentzVector const &  m 
)
void simb::MCTrajectory::Add ( TLorentzVector const &  p,
TLorentzVector const &  m,
std::string const &  process,
bool  keepTransportation = false 
)

Definition at line 136 of file MCTrajectory.cxx.

References ftrajectory, fTrajectoryProcess, ProcessToKey(), and push_back().

140  {
141  // add the the momentum and position, then get the location of the added
142  // bits to store the process
143  this->push_back(p, m);
144 
145  size_t insertLoc = ftrajectory.size() - 1;
146 
147  auto key = this->ProcessToKey(process);
148 
149  // only add a process to the list if it is defined, ie one of the values
150  // allowed in the ProcessToKey() method
151  //
152  // Also, keep 10 (transportation) if the flag allows
153  if(key > 0 && (key != 10 || keepTransportation))
154  fTrajectoryProcess.push_back(std::make_pair(insertLoc, key));
155 
156  return;
157  }
void push_back(value_type const &v)
unsigned char ProcessToKey(std::string const &process) const
ProcessMap fTrajectoryProcess
Definition: MCTrajectory.h:78
list_type ftrajectory
The list of trajectory points.
Definition: MCTrajectory.h:77
const simb::MCTrajectory::value_type & simb::MCTrajectory::at ( const size_type  i) const
inline

Definition at line 175 of file MCTrajectory.h.

References Add(), ftrajectory, and push_back().

Referenced by nnet::TrainingDataAlg::collectVtxFlags(), and Sparsify().

176  { return ftrajectory.at(i); }
list_type ftrajectory
The list of trajectory points.
Definition: MCTrajectory.h:77
simb::MCTrajectory::iterator simb::MCTrajectory::begin ( )
inline

Standard STL methods, to make this class look like an STL map. Again, if you don't know STL, you can just ignore these methods.

Definition at line 158 of file MCTrajectory.h.

References ftrajectory.

Referenced by simb::operator<<().

158 { return ftrajectory.begin(); }
list_type ftrajectory
The list of trajectory points.
Definition: MCTrajectory.h:77
simb::MCTrajectory::const_iterator simb::MCTrajectory::begin ( ) const
inline

Definition at line 159 of file MCTrajectory.h.

References ftrajectory.

159 { return ftrajectory.begin(); }
list_type ftrajectory
The list of trajectory points.
Definition: MCTrajectory.h:77
void simb::MCTrajectory::clear ( )
inline

Definition at line 168 of file MCTrajectory.h.

References ftrajectory.

168 { ftrajectory.clear(); }
list_type ftrajectory
The list of trajectory points.
Definition: MCTrajectory.h:77
double simb::MCTrajectory::E ( const size_type  i) const
inline

Definition at line 156 of file MCTrajectory.h.

References Momentum().

156 { return Momentum(i).E(); }
const TLorentzVector & Momentum(const size_type) const
bool simb::MCTrajectory::empty ( ) const
inline
simb::MCTrajectory::iterator simb::MCTrajectory::end ( void  )
inline

Definition at line 160 of file MCTrajectory.h.

References ftrajectory.

Referenced by simb::operator<<().

160 { return ftrajectory.end(); }
list_type ftrajectory
The list of trajectory points.
Definition: MCTrajectory.h:77
simb::MCTrajectory::const_iterator simb::MCTrajectory::end ( void  ) const
inline

Definition at line 161 of file MCTrajectory.h.

References ftrajectory.

161 { return ftrajectory.end(); }
list_type ftrajectory
The list of trajectory points.
Definition: MCTrajectory.h:77
std::string simb::MCTrajectory::KeyToProcess ( unsigned char const &  key) const

Definition at line 117 of file MCTrajectory.cxx.

Referenced by nnet::TrainingDataAlg::collectVtxFlags().

118  {
119  std::string process("Unknown");
120 
121  if (key == 1) process = "hadElastic";
122  else if(key == 2) process = "pi-Inelastic";
123  else if(key == 3) process = "pi+Inelastic";
124  else if(key == 4) process = "kaon-Inelastic";
125  else if(key == 5) process = "kaon+Inelastic";
126  else if(key == 6) process = "protonInelastic";
127  else if(key == 7) process = "neutronInelastic";
128  else if(key == 8) process = "CoulombScat";
129  else if(key == 9) process = "nCapture";
130  else if(key == 10) process = "Transportation";
131 
132  return process;
133  }
const TLorentzVector & simb::MCTrajectory::Momentum ( const size_type  index) const

Definition at line 44 of file MCTrajectory.cxx.

References ftrajectory.

Referenced by E(), simb::MCParticle::Momentum(), Px(), Py(), and Pz().

45  {
46  const_iterator i = ftrajectory.begin();
47  std::advance(i,index);
48  return (*i).second;
49  }
list_type ftrajectory
The list of trajectory points.
Definition: MCTrajectory.h:77
const simb::MCTrajectory::value_type & simb::MCTrajectory::operator[] ( const size_type  i) const
inline

Definition at line 172 of file MCTrajectory.h.

References ftrajectory.

173  { return ftrajectory[i]; }
list_type ftrajectory
The list of trajectory points.
Definition: MCTrajectory.h:77
const TLorentzVector & simb::MCTrajectory::Position ( const size_type  index) const

The accessor methods described above.

Definition at line 36 of file MCTrajectory.cxx.

References ftrajectory.

Referenced by simb::MCParticle::Position(), T(), TotalLength(), X(), Y(), and Z().

37  {
38  const_iterator i = ftrajectory.begin();
39  std::advance(i,index);
40  return (*i).first;
41  }
list_type ftrajectory
The list of trajectory points.
Definition: MCTrajectory.h:77
unsigned char simb::MCTrajectory::ProcessToKey ( std::string const &  process) const

Definition at line 97 of file MCTrajectory.cxx.

Referenced by Add().

98  {
99  unsigned char key = 0;
100 
101  if (process.compare("hadElastic") == 0) key = 1;
102  else if(process.compare("pi-Inelastic") == 0) key = 2;
103  else if(process.compare("pi+Inelastic") == 0) key = 3;
104  else if(process.compare("kaon-Inelastic") == 0) key = 4;
105  else if(process.compare("kaon+Inelastic") == 0) key = 5;
106  else if(process.compare("protonInelastic") == 0) key = 6;
107  else if(process.compare("neutronInelastic") == 0) key = 7;
108  else if(process.compare("CoulombScat") == 0) key = 8;
109  else if(process.compare("nCapture") == 0) key = 9;
110  else if(process.compare("Transportation") == 0)
111  key = 10;
112 
113  return key;
114  }
void simb::MCTrajectory::push_back ( value_type const &  v)

The only "set" methods for this class; once you've added a trajectory point, you can't take it back.

Referenced by Add(), and at().

void simb::MCTrajectory::push_back ( TLorentzVector const &  p,
TLorentzVector const &  m 
)
double simb::MCTrajectory::Px ( const size_type  i) const
inline

Definition at line 153 of file MCTrajectory.h.

References Momentum().

153 { return Momentum(i).Px(); }
const TLorentzVector & Momentum(const size_type) const
double simb::MCTrajectory::Py ( const size_type  i) const
inline

Definition at line 154 of file MCTrajectory.h.

References Momentum().

154 { return Momentum(i).Py(); }
const TLorentzVector & Momentum(const size_type) const
double simb::MCTrajectory::Pz ( const size_type  i) const
inline

Definition at line 155 of file MCTrajectory.h.

References Momentum().

155 { return Momentum(i).Pz(); }
const TLorentzVector & Momentum(const size_type) const
simb::MCTrajectory::reverse_iterator simb::MCTrajectory::rbegin ( )
inline

Definition at line 162 of file MCTrajectory.h.

References ftrajectory.

Referenced by Sparsify().

162 { return ftrajectory.rbegin(); }
list_type ftrajectory
The list of trajectory points.
Definition: MCTrajectory.h:77
simb::MCTrajectory::const_reverse_iterator simb::MCTrajectory::rbegin ( ) const
inline

Definition at line 163 of file MCTrajectory.h.

References ftrajectory.

163 { return ftrajectory.rbegin(); }
list_type ftrajectory
The list of trajectory points.
Definition: MCTrajectory.h:77
simb::MCTrajectory::reverse_iterator simb::MCTrajectory::rend ( )
inline

Definition at line 164 of file MCTrajectory.h.

References ftrajectory.

164 { return ftrajectory.rend(); }
list_type ftrajectory
The list of trajectory points.
Definition: MCTrajectory.h:77
simb::MCTrajectory::const_reverse_iterator simb::MCTrajectory::rend ( ) const
inline

Definition at line 165 of file MCTrajectory.h.

References ftrajectory.

165 { return ftrajectory.rend(); }
list_type ftrajectory
The list of trajectory points.
Definition: MCTrajectory.h:77
void simb::MCTrajectory::Sparsify ( double  margin = .1,
bool  keep_second_to_last = false 
)

Remove points from trajectory. Straight line interpolation between the remaining points will pass no further than margin from removed points.

Definition at line 160 of file MCTrajectory.cxx.

References at(), dir, ftrajectory, fTrajectoryProcess, rbegin(), and size().

Referenced by simb::MCParticle::SparsifyTrajectory().

161  {
162  // This is a divide-and-conquer algorithm. If the straight line between two
163  // points is close enough to all the intermediate points, then just keep
164  // the endpoints. Otherwise, divide the range in two and try again.
165 
166  // We keep the ranges that need checking in "toCheck". If a range is good
167  // as-is, we put just the starting point in "done". The end-point will be
168  // taken care of by the next range.
169 
170  // Need at least three points to think of removing one
171  // D.R. -- let's up this to four points before we start removing
172  // -- this is helpful when retrieving the energy of the particle prior
173  // -- to a final interaction : (Start, p1, ..., p_(n-1), End)
174  if(size() <= 3 && keep_second_to_last) return;
175  else if(size() <= 2) return;
176 
177  // Deal in terms of distance-squared to save some sqrts
178  margin *= margin;
179 
180  // Deque because we add things still to check on the end, and pop things
181  // we've checked from the front.
182  std::deque<std::pair<int, int> > toCheck;
183  // Start off by trying to replace the whole trajectory with just the
184  // endpoints.
185  toCheck.push_back(std::make_pair(0, size()-1));
186 
187  // use a std::set to keep track of which indices of points we want to
188  // keep because the set does not allow duplicates and it keeps items in
189  // order
190  std::set<int> done;
191  if (keep_second_to_last)
192  done.insert(size()-2); // -- D.R. store the penultimate point
193 
194  while(!toCheck.empty()){
195  const int loIdx = toCheck.front().first;
196  const int hiIdx = toCheck.front().second;
197  toCheck.pop_front();
198 
199  // Should never have been given a degenerate range
200  if(hiIdx < loIdx+2)
201  throw cet::exception("MCTrajectory") << "Degnerate range in Sparsify method";
202 
203  const TVector3 loVec = at(loIdx).first.Vect();
204  const TVector3 hiVec = at(hiIdx).first.Vect();
205 
206  const TVector3 dir = (hiVec-loVec).Unit();
207 
208  // Are all the points in between close enough?
209  bool ok = true;
210  for(int i = loIdx+1; i < hiIdx; ++i){
211  const TVector3 toHere = at(i).first.Vect()-loVec;
212  // Perpendicular distance^2 from the line joining loVec to hiVec
213  const double impact = (toHere-dir.Dot(toHere)*dir).Mag2();
214  if(impact > margin){ok = false; break;}
215  }
216 
217  if(ok){
218  // These points adequately represent this range
219  done.insert(loIdx);
220  }
221  else{
222  // Split in half
223  const int midIdx = (loIdx+hiIdx)/2;
224  // Should never have a range this small
225  if(midIdx == loIdx)
226  throw cet::exception("MCTrajectory") << "Midpoint in sparsification is same as lowpoint";
227  if(midIdx == hiIdx)
228  throw cet::exception("MCTrajectory") << "Midpoint in sparsification is same as hipoint";
229 
230  // The range can be small enough that upon splitting, the new ranges
231  // are degenerate, and should just be written out straight away. Check
232  // for those cases.
233 
234  if(midIdx == loIdx+1){
235  done.insert(loIdx);
236  }
237  else{
238  toCheck.push_back(std::make_pair(loIdx, midIdx));
239  }
240 
241  if(midIdx == hiIdx-1){
242  done.insert(midIdx);
243  }
244  else{
245  toCheck.push_back(std::make_pair(midIdx, hiIdx));
246  }
247  }
248  } // end while
249 
250  // now make sure we have not left out any of the indices where interesting
251  // processes were recorded
252  std::map< size_t, unsigned char> processMap;
253  for(auto itr : fTrajectoryProcess){
254  done.insert(itr.first);
255  processMap[itr.first] = itr.second;
256  }
257 
258  // Look up the trajectory points at the stored indices, write them to a new
259  // trajectory
260  const unsigned int I = done.size();
261  list_type newTraj;
262  newTraj.reserve(I+1);
263 
264  // make a new process map as well based on these points
265  ProcessMap newProcMap;
266 
267  for(auto idx : done){
268  newTraj.push_back(at(idx));
269  if(processMap.count(idx) > 0){
270  newProcMap.push_back(std::make_pair(newTraj.size() - 1,
271  processMap.find(idx)->second)
272  );
273  }
274  }
275 
276  // Remember to add the very last point in if it hasn't already been added
277  if(done.count(ftrajectory.size() - 1) < 1) newTraj.push_back(*rbegin());
278 
279  // Replace trajectory and fTrajectoryProcess with new versions
280  std::swap(ftrajectory, newTraj );
281  std::swap(fTrajectoryProcess, newProcMap);
282 
283  return;
284  }
reverse_iterator rbegin()
Definition: MCTrajectory.h:162
const value_type & at(const size_type i) const
Definition: MCTrajectory.h:175
std::vector< std::pair< TLorentzVector, TLorentzVector > > list_type
Definition: MCTrajectory.h:63
ProcessMap fTrajectoryProcess
Definition: MCTrajectory.h:78
std::vector< std::pair< size_t, unsigned char > > ProcessMap
Definition: MCTrajectory.h:71
size_type size() const
Definition: MCTrajectory.h:166
TDirectory * dir
Definition: macro.C:5
void swap(lar::deep_const_fwd_iterator_nested< CITER, INNERCONTEXTRACT > &a, lar::deep_const_fwd_iterator_nested< CITER, INNERCONTEXTRACT > &b)
cet::coded_exception< error, detail::translate > exception
Definition: exception.h:33
list_type ftrajectory
The list of trajectory points.
Definition: MCTrajectory.h:77
void simb::MCTrajectory::swap ( simb::MCTrajectory other)
inline

Definition at line 169 of file MCTrajectory.h.

References ftrajectory.

170  { ftrajectory.swap( other.ftrajectory ); }
list_type ftrajectory
The list of trajectory points.
Definition: MCTrajectory.h:77
double simb::MCTrajectory::T ( const size_type  i) const
inline

Definition at line 152 of file MCTrajectory.h.

References Position().

152 { return Position(i).T(); }
const TLorentzVector & Position(const size_type) const
The accessor methods described above.
double simb::MCTrajectory::TotalLength ( ) const

Definition at line 52 of file MCTrajectory.cxx.

References larg4::dist(), n, Position(), and size().

Referenced by sim::dump::DumpMCParticle(), and cosmic::BeamFlashTrackMatchTaggerAlg::RunHypothesisComparison().

53  {
54  const int N = size();
55  if(N < 2) return 0;
56 
57  // We take the sum of the straight lines between the trajectory points
58  double dist = 0;
59  for(int n = 0; n < N-1; ++n){
60  dist += (Position(n+1)-Position(n)).Vect().Mag();
61  }
62 
63  return dist;
64  }
const TLorentzVector & Position(const size_type) const
The accessor methods described above.
size_type size() const
Definition: MCTrajectory.h:166
constexpr double dist(const TReal *x, const TReal *y, const unsigned int dimension)
Char_t n[5]
simb::MCTrajectory::ProcessMap const & simb::MCTrajectory::TrajectoryProcesses ( ) const
inline

Definition at line 188 of file MCTrajectory.h.

References fTrajectoryProcess.

Referenced by nnet::TrainingDataAlg::collectVtxFlags().

188 { return fTrajectoryProcess; }
ProcessMap fTrajectoryProcess
Definition: MCTrajectory.h:78
double simb::MCTrajectory::X ( const size_type  i) const
inline

Definition at line 149 of file MCTrajectory.h.

References Position().

Referenced by evdb_tool::DrawLArVoxel3D::Draw(), evd::SimulationDrawer::MCTruth3D(), and evd::SimulationDrawer::MCTruthOrtho().

149 { return Position(i).X(); }
const TLorentzVector & Position(const size_type) const
The accessor methods described above.
double simb::MCTrajectory::Y ( const size_type  i) const
inline

Definition at line 150 of file MCTrajectory.h.

References Position().

Referenced by evdb_tool::DrawLArVoxel3D::Draw(), evd::SimulationDrawer::MCTruth3D(), and evd::SimulationDrawer::MCTruthOrtho().

150 { return Position(i).Y(); }
const TLorentzVector & Position(const size_type) const
The accessor methods described above.
double simb::MCTrajectory::Z ( const size_type  i) const
inline

Definition at line 151 of file MCTrajectory.h.

References Position().

Referenced by evdb_tool::DrawLArVoxel3D::Draw(), evd::SimulationDrawer::MCTruth3D(), and evd::SimulationDrawer::MCTruthOrtho().

151 { return Position(i).Z(); }
const TLorentzVector & Position(const size_type) const
The accessor methods described above.

Friends And Related Function Documentation

std::ostream& operator<< ( std::ostream &  output,
const MCTrajectory list 
)
friend

Definition at line 67 of file MCTrajectory.cxx.

68  {
69  // Determine a field width for the voxel number.
70  MCTrajectory::size_type numberOfTrajectories = list.size();
71  int numberOfDigits = (int) std::log10( (double) numberOfTrajectories ) + 1;
72 
73  // A simple header.
74  output.width( numberOfDigits );
75  output << "#" << ": < position (x,y,z,t), momentum (Px,Py,Pz,E) >" << std::endl;
76 
77  // Write each trajectory point on a separate line.
78  MCTrajectory::size_type nTrajectory = 0;
79  for ( MCTrajectory::const_iterator trajectory = list.begin(); trajectory != list.end(); ++trajectory, ++nTrajectory ){
80  output.width( numberOfDigits );
81  output << nTrajectory << ": "
82  << "< (" << (*trajectory).first.X()
83  << "," << (*trajectory).first.Y()
84  << "," << (*trajectory).first.Z()
85  << "," << (*trajectory).first.T()
86  << ") , (" << (*trajectory).second.Px()
87  << "," << (*trajectory).second.Py()
88  << "," << (*trajectory).second.Pz()
89  << "," << (*trajectory).second.E()
90  << ") >" << std::endl;
91  }
92 
93  return output;
94  }
list_type::size_type size_type
Definition: MCTrajectory.h:69
list_type::const_iterator const_iterator
Definition: MCTrajectory.h:66

Member Data Documentation

list_type simb::MCTrajectory::ftrajectory
private

The list of trajectory points.

Definition at line 77 of file MCTrajectory.h.

Referenced by Add(), at(), begin(), clear(), empty(), end(), MCTrajectory(), Momentum(), operator[](), Position(), rbegin(), rend(), size(), Sparsify(), and swap().

ProcessMap simb::MCTrajectory::fTrajectoryProcess
private

map of the scattering process to index in ftrajectory for a given point

Definition at line 78 of file MCTrajectory.h.

Referenced by Add(), Sparsify(), and TrajectoryProcesses().


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