LArSoft  v09_90_00
Liquid Argon Software toolkit - https://larsoft.org/
shower::TCShowerAlg Class Reference

#include "TCShowerAlg.h"

Public Member Functions

 TCShowerAlg (fhicl::ParameterSet const &pset)
 
int makeShowers (detinfo::DetectorClocksData const &dataClock, detinfo::DetectorPropertiesData const &detProp, std::vector< art::Ptr< recob::PFParticle >> const &pfplist, std::vector< art::Ptr< recob::Cluster >> const &clusterlist, std::vector< art::Ptr< recob::Hit >> const &hitlist, art::FindManyP< recob::Hit > const &cls_fm, art::FindManyP< recob::Cluster > const &clspfp_fm, art::FindManyP< recob::Vertex > const &vtxpfp_fm, art::FindManyP< recob::Cluster > const &hitcls_fm, art::FindManyP< recob::Track > const &trkpfp_fm)
 

Public Attributes

geo::Vector_t shwDir
 
TVector3 dcosVtxErr
 
geo::Point_t shwvtx
 
TVector3 xyzErr
 
std::vector< double > totalEnergy
 
std::vector< double > totalEnergyErr
 
std::vector< double > dEdx
 
std::vector< double > dEdxErr
 
int bestplane
 
std::vector< art::Ptr< recob::Hit > > showerHits
 

Private Member Functions

int goodHit (detinfo::DetectorClocksData const &dataClock, detinfo::DetectorPropertiesData const &detProp, art::Ptr< recob::Hit > const &, double maxDist, double minDistVert, std::map< geo::PlaneID, double > const &trk_wire1, std::map< geo::PlaneID, double > const &trk_tick1, std::map< geo::PlaneID, double > const &trk_wire2, std::map< geo::PlaneID, double > const &trk_tick2) const
 
int goodHit (detinfo::DetectorClocksData const &dataClock, detinfo::DetectorPropertiesData const &detProp, art::Ptr< recob::Hit > const &, double maxDist, double minDistVert, std::map< geo::PlaneID, double > const &trk_wire1, std::map< geo::PlaneID, double > const &trk_tick1, std::map< geo::PlaneID, double > const &trk_wire2, std::map< geo::PlaneID, double > const &trk_tick2, int &pull) const
 
bool addShowerHit (art::Ptr< recob::Hit > hit, std::vector< art::Ptr< recob::Hit >> showerhits) const
 

Private Attributes

calo::CalorimetryAlg fCalorimetryAlg
 
pma::ProjectionMatchingAlg fProjectionMatchingAlg
 

Detailed Description

Definition at line 38 of file TCShowerAlg.h.

Constructor & Destructor Documentation

shower::TCShowerAlg::TCShowerAlg ( fhicl::ParameterSet const &  pset)

Definition at line 44 of file TCShowerAlg.cxx.

45  : fCalorimetryAlg(pset.get<fhicl::ParameterSet>("CalorimetryAlg"))
46  , fProjectionMatchingAlg(pset.get<fhicl::ParameterSet>("ProjectionMatchingAlg"))
47  {}
pma::ProjectionMatchingAlg fProjectionMatchingAlg
Definition: TCShowerAlg.h:67
calo::CalorimetryAlg fCalorimetryAlg
Definition: TCShowerAlg.h:66

Member Function Documentation

bool shower::TCShowerAlg::addShowerHit ( art::Ptr< recob::Hit hit,
std::vector< art::Ptr< recob::Hit >>  showerhits 
) const
private

Definition at line 446 of file TCShowerAlg.cxx.

References art::Ptr< T >::key().

Referenced by makeShowers().

448  {
449 
450  for (size_t i = 0; i < showerhits.size(); ++i) {
451  if (hit.key() == showerhits[i].key()) return false;
452  }
453 
454  return true;
455 
456  } // addShowerHit
key_type key() const noexcept
Definition: Ptr.h:166
int shower::TCShowerAlg::goodHit ( detinfo::DetectorClocksData const &  dataClock,
detinfo::DetectorPropertiesData const &  detProp,
art::Ptr< recob::Hit > const &  hit,
double  maxDist,
double  minDistVert,
std::map< geo::PlaneID, double > const &  trk_wire1,
std::map< geo::PlaneID, double > const &  trk_tick1,
std::map< geo::PlaneID, double > const &  trk_wire2,
std::map< geo::PlaneID, double > const &  trk_tick2 
) const
private

Definition at line 356 of file TCShowerAlg.cxx.

Referenced by makeShowers().

365  {
366  int pull = 0;
367  return goodHit(clockData,
368  detProp,
369  hit,
370  maxDist,
371  minDistVert,
372  trk_wire1,
373  trk_tick1,
374  trk_wire2,
375  trk_tick2,
376  pull);
377 
378  } // goodHit
int goodHit(detinfo::DetectorClocksData const &dataClock, detinfo::DetectorPropertiesData const &detProp, art::Ptr< recob::Hit > const &, double maxDist, double minDistVert, std::map< geo::PlaneID, double > const &trk_wire1, std::map< geo::PlaneID, double > const &trk_tick1, std::map< geo::PlaneID, double > const &trk_wire2, std::map< geo::PlaneID, double > const &trk_tick2) const
int shower::TCShowerAlg::goodHit ( detinfo::DetectorClocksData const &  dataClock,
detinfo::DetectorPropertiesData const &  detProp,
art::Ptr< recob::Hit > const &  hit,
double  maxDist,
double  minDistVert,
std::map< geo::PlaneID, double > const &  trk_wire1,
std::map< geo::PlaneID, double > const &  trk_tick1,
std::map< geo::PlaneID, double > const &  trk_wire2,
std::map< geo::PlaneID, double > const &  trk_tick2,
int &  pull 
) const
private

Definition at line 385 of file TCShowerAlg.cxx.

References util::abs(), larg4::dist(), detinfo::DetectorPropertiesData::DriftVelocity(), detinfo::DetectorPropertiesData::Efield(), recob::Hit::PeakTime(), detinfo::sampling_rate(), detinfo::DetectorPropertiesData::Temperature(), geo::WireID::Wire, recob::Hit::WireID(), geo::GeometryCore::WirePitch(), x1, x2, y1, and y2.

395  {
397 
398  double wirePitch = geom->WirePitch(hit->WireID());
399  double tickToDist = detProp.DriftVelocity(detProp.Efield(), detProp.Temperature());
400  tickToDist *= 1.e-3 * sampling_rate(clockData); // 1e-3 is conversion of 1/us to 1/ns
401  double UnitsPerTick = tickToDist / wirePitch;
402 
403  double x0 = hit->WireID().Wire;
404  double y0 = hit->PeakTime() * UnitsPerTick;
405 
406  double x1 = trk_wire1.at(hit->WireID());
407  double y1 = trk_tick1.at(hit->WireID()) * UnitsPerTick;
408 
409  double x2 = trk_wire2.at(hit->WireID());
410  double y2 = trk_tick2.at(hit->WireID()) * UnitsPerTick;
411 
412  double distToVert = std::hypot(x0 - x1, y0 - y1);
413  if (distToVert < minDistVert) return -1;
414 
415  if (x2 > x1) {
416  if (distToVert < 50)
417  maxDist /= 4;
418  else if (distToVert < 100)
419  maxDist /= 2; // trying to exclude photon showers
420  }
421 
422  // exclude cluster if it's "behind" the vertex
423  double a = std::hypot(x2 - x1, y2 - y1);
424  double b = std::hypot(x0 - x1, y0 - y1);
425  double c = std::hypot(x0 - x2, y0 - y2);
426  double costheta = -(pow(c, 2) - pow(a, 2) - pow(b, 2)) / (2 * a * b);
427  if (costheta < 0) return -1;
428 
429  double dist =
430  std::abs((y2 - y1) * x0 - (x2 - x1) * y0 + x2 * y1 - y2 * x1) / std::hypot(y2 - y1, x2 - x1);
431 
432  if (dist < maxDist) {
433  if (((y2 - y1) * x0 - (x2 - x1) * y0 + x2 * y1 - y2 * x1) > 0)
434  pull = 1;
435  else
436  pull = -1;
437  return 1;
438  }
439 
440  return 0;
441 
442  } // goodHit
Float_t y1[n_points_granero]
Definition: compare.C:5
Float_t x1[n_points_granero]
Definition: compare.C:5
constexpr auto abs(T v)
Returns the absolute value of the argument.
WireID_t Wire
Index of the wire within its plane.
Definition: geo_types.h:563
geo::WireID const & WireID() const
Initial tdc tick for hit.
Definition: Hit.h:280
Float_t y2[n_points_geant4]
Definition: compare.C:26
float PeakTime() const
Time of the signal peak, in tick units.
Definition: Hit.h:220
constexpr double dist(const TReal *x, const TReal *y, const unsigned int dimension)
Float_t x2[n_points_geant4]
Definition: compare.C:26
double sampling_rate(DetectorClocksData const &data)
Returns the period of the TPC readout electronics clock.
Length_t WirePitch(PlaneID const &planeid=plane_zero) const
Returns the distance between two consecutive wires.
int shower::TCShowerAlg::makeShowers ( detinfo::DetectorClocksData const &  dataClock,
detinfo::DetectorPropertiesData const &  detProp,
std::vector< art::Ptr< recob::PFParticle >> const &  pfplist,
std::vector< art::Ptr< recob::Cluster >> const &  clusterlist,
std::vector< art::Ptr< recob::Hit >> const &  hitlist,
art::FindManyP< recob::Hit > const &  cls_fm,
art::FindManyP< recob::Cluster > const &  clspfp_fm,
art::FindManyP< recob::Vertex > const &  vtxpfp_fm,
art::FindManyP< recob::Cluster > const &  hitcls_fm,
art::FindManyP< recob::Track > const &  trkpfp_fm 
)

Definition at line 49 of file TCShowerAlg.cxx.

References util::abs(), addShowerHit(), detinfo::DetectorPropertiesData::ConvertXToTicks(), dEdx, dEdxErr, recob::Track::End(), recob::Track::EndDirection(), goodHit(), geo::GeometryCore::Iterate(), recob::Vertex::position(), showerHits, shwDir, shwvtx, totalEnergy, totalEnergyErr, lar::dump::vector(), recob::Track::Vertex(), recob::Track::VertexDirection(), and geo::GeometryCore::WireCoordinate().

Referenced by shower::TCShower::getShowersWithoutSlices(), and shower::TCShower::getShowersWithSlices().

59  {
60  totalEnergy.resize(2);
61  totalEnergyErr.resize(2);
62  dEdx.resize(2);
63  dEdxErr.resize(2);
64 
66 
67  std::vector<pfpStuff> allpfps;
68 
69  // put together pfparticle information
70  for (size_t i = 0; i < pfplist.size(); ++i) {
71  pfpStuff thispfp;
72  thispfp.hits.clear();
73  thispfp.clsIDs.clear();
74  thispfp.pfp = pfplist[i];
75 
76  std::vector<art::Ptr<recob::Vertex>> thisvtxlist = vtxpfp_fm.at(pfplist[i].key());
77  if (thisvtxlist.size()) thispfp.vtx = thisvtxlist[0];
78 
79  std::vector<art::Ptr<recob::Track>> thistrklist = trkpfp_fm.at(pfplist[i].key());
80  if (thistrklist.size()) thispfp.trk = thistrklist[0];
81 
82  std::vector<art::Ptr<recob::Cluster>> thisclusterlist = clspfp_fm.at(pfplist[i].key());
83  std::vector<int> clustersize;
84 
85  for (size_t j = 0; j < thisclusterlist.size(); ++j) {
86 
87  thispfp.clsIDs.push_back(thisclusterlist[j]->ID());
88 
89  std::vector<art::Ptr<recob::Hit>> thishitlist = cls_fm.at(thisclusterlist[j].key());
90  clustersize.push_back((int)thishitlist.size());
91 
92  for (size_t k = 0; k < thishitlist.size(); ++k) {
93  thispfp.hits.push_back(thishitlist[k]);
94  } // loop through hits
95 
96  } // loop through clusters
97 
98  if (clustersize.size() == 3) {
99  if (!thispfp.vtx) continue;
100  if (!thispfp.trk) continue;
101 
102  allpfps.push_back(thispfp);
103 
104  double tick = detProp.ConvertXToTicks(thispfp.vtx->position().X(), geo::PlaneID(0, 0, 2));
105  int wire = geom->WireCoordinate(thispfp.vtx->position(), geo::PlaneID(0, 0, 2));
106 
107  std::cout << "pfp " << thispfp.pfp->Self() + 1 << " cluster sizes " << clustersize[0] << ":"
108  << clustersize[1] << ":" << clustersize[2] << " vertex " << thispfp.vtx->ID()
109  << " " << tick << ":" << wire << " z " << thispfp.vtx->position().Z()
110  << std::endl;
111 
112  } // add pfp to list
113 
114  } // loop through pfparticles
115 
116  std::sort(allpfps.begin(), allpfps.end(), comparePFP);
117  std::reverse(allpfps.begin(), allpfps.end());
118 
119  std::cout << "sorted pfps: ";
120  for (size_t i = 0; i < allpfps.size(); ++i)
121  std::cout << allpfps[i].pfp->Self() + 1 << " ";
122  std::cout << std::endl;
123 
124  bool showerCandidate = false;
125 
126  for (size_t i = 0; i < allpfps.size(); ++i) {
127 
128  showerHits.clear();
129 
130  art::Ptr<recob::Vertex> pfpvtx = allpfps[i].vtx;
131  art::Ptr<recob::Track> pfptrk = allpfps[i].trk;
132  std::vector<art::Ptr<recob::Hit>> pfphits = allpfps[i].hits;
133  std::vector<art::Ptr<recob::Cluster>> pfpcls = clspfp_fm.at(allpfps[i].pfp.key());
134 
135  std::cout << "pfp " << allpfps[i].pfp->Self() + 1 << " hits " << pfphits.size() << std::endl;
136 
137  TVector3 vtx;
138  vtx[0] = pfpvtx->position().X();
139  vtx[1] = pfpvtx->position().Y();
140  vtx[2] = pfpvtx->position().Z();
141 
142  if (pfptrk->Vertex().Z() < pfptrk->End().Z()) {
143  shwvtx = pfptrk->Vertex();
144  shwDir = pfptrk->VertexDirection();
145  }
146  else {
147  shwvtx = pfptrk->End();
148  shwDir = -pfptrk->EndDirection();
149  }
150 
151  int tolerance = 60; // how many shower like cluster you need to define a shower
152  double pullTolerance = 0.7; // hits should be evenly distributed around the track
153  double maxDist = 20; // how far a shower like cluster can be from the track
154  double minDistVert = 15; // exclude tracks near the vertex
155 
156  // TODO: count number of shower-like hits hits "behind" vertex
157  // if high, pick an earlier cluster and refit hits
158  // this is going to take some restructuring
159 
160  if (pfphits.size() < 30) continue;
161  if (pfphits.size() > 500) continue;
162  // adjust tolerances for short tracks
163  if (pfphits.size() < 90) {
164  tolerance = 50;
165  pullTolerance = 0.9;
166  }
167  if (pfphits.size() > 400)
168  tolerance = 200;
169  else if (pfphits.size() > 100)
170  tolerance = 100; // RSF added used to be 100, 100
171 
172  // add pfp hits to shower
173  for (size_t ii = 0; ii < pfphits.size(); ++ii) {
174  if (addShowerHit(pfphits[ii], showerHits)) showerHits.push_back(pfphits[ii]);
175  } // loop over pfphits
176 
177  int nShowerHits = 0;
178  double showerHitPull = 0;
179 
180  geo::Point_t pfpStart = shwvtx;
181  geo::Point_t pfpPt2 = shwvtx + shwDir; // a second point along the track
182 
183  // track vertex
184  std::map<geo::PlaneID, double> trk_tick1;
185  std::map<geo::PlaneID, double> trk_wire1;
186 
187  // second track point
188  std::map<geo::PlaneID, double> trk_tick2;
189  std::map<geo::PlaneID, double> trk_wire2;
190 
191  for (auto const& planeid : geom->Iterate<geo::PlaneID>()) {
192  trk_tick1[planeid] = detProp.ConvertXToTicks(pfpStart.X(), planeid);
193  trk_wire1[planeid] = geom->WireCoordinate(pfpStart, planeid);
194  trk_tick2[planeid] = detProp.ConvertXToTicks(pfpPt2.X(), planeid);
195  trk_wire2[planeid] = geom->WireCoordinate(pfpPt2, planeid);
196  }
197 
198  for (size_t j = 0; j < clusterlist.size(); ++j) {
199  std::vector<art::Ptr<recob::Hit>> cls_hitlist = cls_fm.at(clusterlist[j].key());
200 
201  if (clusterlist[j]->ID() > 0 && cls_hitlist.size() > 10) continue;
202  if (cls_hitlist.size() > 50) continue;
203 
204  bool isGoodCluster = false; // true if the hit belongs to a cluster that
205  // should be added to the shower
206 
207  bool skipit = false; // skip clusters already in the pfp
208  for (size_t k = 0; k < allpfps[i].clsIDs.size(); ++k) {
209  if (allpfps[i].clsIDs[k] == clusterlist[j]->ID()) skipit = true;
210  }
211  if (skipit) continue;
212 
213  for (size_t jj = 0; jj < cls_hitlist.size(); ++jj) {
214  int isGoodHit = goodHit(clockData,
215  detProp,
216  cls_hitlist[jj],
217  maxDist,
218  minDistVert,
219  trk_wire1,
220  trk_tick1,
221  trk_wire2,
222  trk_tick2);
223 
224  if (isGoodHit == -1) {
225  isGoodCluster = false;
226  break;
227  }
228  else if (isGoodHit == 1) {
229  isGoodCluster = true;
230  }
231 
232  } // loop over hits in cluster
233 
234  // add hits to shower
235  if (isGoodCluster) {
236  for (size_t jj = 0; jj < cls_hitlist.size(); ++jj) {
237  nShowerHits++;
238 
239  int showerHitPullAdd = 0;
240  goodHit(clockData,
241  detProp,
242  cls_hitlist[jj],
243  maxDist,
244  minDistVert,
245  trk_wire1,
246  trk_tick1,
247  trk_wire2,
248  trk_tick2,
249  showerHitPullAdd);
250  showerHitPull += showerHitPullAdd;
251 
252  if (addShowerHit(cls_hitlist[jj], showerHits)) showerHits.push_back(cls_hitlist[jj]);
253 
254  } // loop over hits in cluster
255  } // cluster contains hit close to track
256 
257  } // loop through cluserlist
258 
259  showerHitPull /= nShowerHits;
260 
261  std::cout << "shower hits " << showerHits.size() << " " << nShowerHits << " shower pull "
262  << showerHitPull << std::endl;
263 
264  if (nShowerHits > tolerance && std::abs(showerHitPull) < pullTolerance) {
265  showerCandidate = true;
266  std::cout << "SHOWER CANDIDATE" << std::endl;
267  // loop over hits to find those that aren't associated with any clusters
268  if (nShowerHits > 400) maxDist *= 2; // TODO: optimize this threshold
269  for (size_t k = 0; k < hitlist.size(); ++k) {
270  std::vector<art::Ptr<recob::Cluster>> hit_clslist = hitcls_fm.at(hitlist[k].key());
271  if (hit_clslist.size()) continue;
272 
273  int isGoodHit = goodHit(clockData,
274  detProp,
275  hitlist[k],
276  maxDist * 2,
277  minDistVert * 2,
278  trk_wire1,
279  trk_tick1,
280  trk_wire2,
281  trk_tick2);
282  if (isGoodHit == 1 && addShowerHit(hitlist[k], showerHits))
283  showerHits.push_back(hitlist[k]);
284  } // loop over hits
285 
286  // loop over clusters to see if any fall within the shower
287  for (size_t k = 0; k < clusterlist.size(); ++k) {
288  std::vector<art::Ptr<recob::Hit>> cls_hitlist = cls_fm.at(clusterlist[k].key());
289  if (clusterlist[k]->ID() > 0 && cls_hitlist.size() > 50) continue;
290 
291  double thisDist = maxDist;
292  double thisMin = minDistVert;
293 
294  if (cls_hitlist.size() < 10) {
295  thisDist *= 4;
296  thisMin *= 4;
297  }
298  else if (cls_hitlist.size() < 30)
299  thisDist *= 2;
300 
301  int nhits = 0;
302  int ngoodhits = 0;
303 
304  // are the cluster hits close?
305  for (size_t kk = 0; kk < cls_hitlist.size(); ++kk) {
306  nhits++;
307  int isGoodHit = goodHit(clockData,
308  detProp,
309  cls_hitlist[kk],
310  thisDist,
311  thisMin,
312  trk_wire1,
313  trk_tick1,
314  trk_wire2,
315  trk_tick2);
316  if (isGoodHit == -1) {
317  ngoodhits = 0;
318  break;
319  }
320  else if (isGoodHit == 1) {
321  ngoodhits++;
322  }
323  } // loop over cluster hits
324 
325  double fracGood = (double)ngoodhits / nhits;
326 
327  bool isGoodTrack = fracGood > 0.4;
328 
329  if (isGoodTrack) {
330  for (size_t kk = 0; kk < cls_hitlist.size(); ++kk) {
331  if (addShowerHit(cls_hitlist[kk], showerHits)) showerHits.push_back(cls_hitlist[kk]);
332  } // loop over hits to add them to showe
333  }
334 
335  } // loop through clusterlist
336 
337  } // decide if shower
338 
339  if (showerCandidate) {
340  std::cout << "THIS IS THE SHOWER PFP: " << allpfps[i].pfp->Self() + 1 << std::endl;
341  break;
342  }
343 
344  } // loop through allpfps
345 
346  if (showerCandidate) return 1;
347 
348  return 0;
349  } // makeShowers
details::range_type< T > Iterate() const
Initializes the specified ID with the ID of the first cryostat.
Definition: GeometryCore.h:541
Length_t WireCoordinate(Point_t const &pos, PlaneID const &planeid) const
Returns the index of the nearest wire to the specified position.
std::vector< double > totalEnergyErr
Definition: TCShowerAlg.h:46
bool addShowerHit(art::Ptr< recob::Hit > hit, std::vector< art::Ptr< recob::Hit >> showerhits) const
std::vector< double > dEdx
Definition: TCShowerAlg.h:47
The data type to uniquely identify a Plane.
Definition: geo_types.h:463
constexpr auto abs(T v)
Returns the absolute value of the argument.
Vector_t VertexDirection() const
Access to track direction at different points.
Definition: Track.h:166
geo::Point_t shwvtx
Definition: TCShowerAlg.h:43
IDparameter< geo::PlaneID > PlaneID
Member type of validated geo::PlaneID parameter.
auto vector(Vector const &v)
Returns a manipulator which will print the specified array.
Definition: DumpUtils.h:289
geo::Vector_t shwDir
Definition: TCShowerAlg.h:41
Point_t const & Vertex() const
Access to track position at different points.
Definition: Track.h:158
tick_as<> tick
Tick number, represented by std::ptrdiff_t.
Definition: electronics.h:73
ROOT::Math::PositionVector3D< ROOT::Math::Cartesian3D< double >, ROOT::Math::GlobalCoordinateSystemTag > Point_t
Type for representation of position in physical 3D space.
Definition: geo_vectors.h:180
Vector_t EndDirection() const
Access to track direction at different points.
Definition: Track.h:167
std::vector< double > totalEnergy
Definition: TCShowerAlg.h:45
Point_t const & End() const
Access to track position at different points.
Definition: Track.h:159
int goodHit(detinfo::DetectorClocksData const &dataClock, detinfo::DetectorPropertiesData const &detProp, art::Ptr< recob::Hit > const &, double maxDist, double minDistVert, std::map< geo::PlaneID, double > const &trk_wire1, std::map< geo::PlaneID, double > const &trk_tick1, std::map< geo::PlaneID, double > const &trk_wire2, std::map< geo::PlaneID, double > const &trk_tick2) const
std::vector< art::Ptr< recob::Hit > > showerHits
Definition: TCShowerAlg.h:50
std::vector< double > dEdxErr
Definition: TCShowerAlg.h:48
const Point_t & position() const
Return vertex 3D position.
Definition: Vertex.h:64

Member Data Documentation

int shower::TCShowerAlg::bestplane

Definition at line 49 of file TCShowerAlg.h.

Referenced by shower::TCShower::produce().

TVector3 shower::TCShowerAlg::dcosVtxErr

Definition at line 42 of file TCShowerAlg.h.

Referenced by shower::TCShower::produce().

std::vector<double> shower::TCShowerAlg::dEdx

Definition at line 47 of file TCShowerAlg.h.

Referenced by makeShowers(), and shower::TCShower::produce().

std::vector<double> shower::TCShowerAlg::dEdxErr

Definition at line 48 of file TCShowerAlg.h.

Referenced by makeShowers(), and shower::TCShower::produce().

calo::CalorimetryAlg shower::TCShowerAlg::fCalorimetryAlg
private

Definition at line 66 of file TCShowerAlg.h.

pma::ProjectionMatchingAlg shower::TCShowerAlg::fProjectionMatchingAlg
private

Definition at line 67 of file TCShowerAlg.h.

std::vector<art::Ptr<recob::Hit> > shower::TCShowerAlg::showerHits

Definition at line 50 of file TCShowerAlg.h.

Referenced by makeShowers(), and shower::TCShower::produce().

geo::Vector_t shower::TCShowerAlg::shwDir

Definition at line 41 of file TCShowerAlg.h.

Referenced by makeShowers(), and shower::TCShower::produce().

geo::Point_t shower::TCShowerAlg::shwvtx

Definition at line 43 of file TCShowerAlg.h.

Referenced by makeShowers(), and shower::TCShower::produce().

std::vector<double> shower::TCShowerAlg::totalEnergy

Definition at line 45 of file TCShowerAlg.h.

Referenced by makeShowers(), and shower::TCShower::produce().

std::vector<double> shower::TCShowerAlg::totalEnergyErr

Definition at line 46 of file TCShowerAlg.h.

Referenced by makeShowers(), and shower::TCShower::produce().

TVector3 shower::TCShowerAlg::xyzErr

Definition at line 44 of file TCShowerAlg.h.

Referenced by shower::TCShower::produce().


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