LArSoft  v07_13_02
Liquid Argon Software toolkit - http://larsoft.org/
TrackShowerSeparationAlg.h
Go to the documentation of this file.
1 // Class: TrackShowerSeparationAlg
3 // File: TrackShowerSeparationAlg.h
4 // Author: Mike Wallbank (m.wallbank@sheffield.ac.uk), November 2015
5 //
6 // Track/shower separation class.
7 // Provides methods for removing hits associated with track-like
8 // objects.
9 // To be run after track reconstruction, before shower reconstruction.
11 
12 #ifndef TrackShowerSeparationAlg_hxx
13 #define TrackShowerSeparationAlg_hxx
14 
15 // Framework
16 #include "fhiclcpp/ParameterSet.h"
24 
25 // larsoft
33 
34 // ROOT
35 #include "TTree.h"
36 #include "TMath.h"
37 #include "TPrincipal.h"
38 #include "TGraph2D.h"
39 #include "TPolyLine3D.h"
40 #include "TCanvas.h"
41 
42 namespace shower {
43  class TrackShowerSeparationAlg;
44  class ReconTrack;
45 }
46 
48  public:
49 
50  ReconTrack(int id) {
51  fID = id;
52  fTrack = false;
53  fShower = false;
54  fShowerTrack = false;
55  fShowerCone = false;
56  }
57 
58  // Setters
59  void SetVertex(TVector3 vertex) { fVertex = vertex; }
60  void SetEnd(TVector3 end) { fEnd = end; }
61  void SetLength(double length) { fLength = length; }
62  void SetVertexDir(TVector3 vertexDir) { fVertexDir = vertexDir; }
63  void SetDirection(TVector3 direction) { fDirection = direction; }
65  void SetSpacePoints(std::vector<art::Ptr<recob::SpacePoint> > spacePoints) { fSpacePoints = spacePoints; }
66 
67  void AddForwardTrack(int track) { if (std::find(fForwardConeTracks.begin(), fForwardConeTracks.end(), track) == fForwardConeTracks.end()) fForwardConeTracks.push_back(track); }
68  void AddBackwardTrack(int track) { if (std::find(fBackwardConeTracks.begin(), fBackwardConeTracks.end(), track) == fBackwardConeTracks.end()) fBackwardConeTracks.push_back(track); }
69  void AddShowerTrack(int track) { fShowerTracks.push_back(track); }
70 
71  void AddForwardSpacePoint(int spacePoint) { fForwardSpacePoints.push_back(spacePoint); }
72  void AddBackwardSpacePoint(int spacePoint) { fBackwardSpacePoints.push_back(spacePoint); }
73  void AddCylinderSpacePoint(int spacePoint) { fCylinderSpacePoints.push_back(spacePoint); }
74  void AddSphereSpacePoint(int spacePoint) { fSphereSpacePoints.push_back(spacePoint); }
75  void AddIsolationSpacePoint(int spacePoint, double distance) { fIsolationSpacePoints[spacePoint] = distance; }
76 
77  // Getters
78  int ID() const { return fID; }
79  TVector3 Vertex() const { return fVertex; }
80  TVector3 End() const { return fEnd; }
81  double Length() const { return fLength; }
82  TVector3 VertexDirection() const { return fVertexDir; }
83  TVector3 Direction() const { return fDirection; }
84  const std::vector<art::Ptr<recob::Hit> >& Hits() const { return fHits; }
85  const std::vector<art::Ptr<recob::SpacePoint> >& SpacePoints() const { return fSpacePoints; }
86 
87  void FlipTrack() {
88  TVector3 tmp = fEnd;
89  fEnd = fVertex;
90  fVertex = tmp;
91  fDirection *= -1;
92  }
93 
94  void MakeShower() {
95  if (fTrack)
96  this->MakeShowerTrack();
97  else
98  this->MakeShowerCone();
99  }
101  fShower = true;
102  fShowerTrack = true;
103  fShowerCone = false;
104  fTrack = false;
105  }
106  void MakeShowerCone() {
107  fShower = true;
108  fShowerCone = true;
109  fShowerTrack = false;
110  fTrack = false;
111  }
112  void MakeTrack() {
113  fTrack = true;
114  fShower = false;
115  fShowerTrack = false;
116  fShowerCone = false;
117  }
118 
119  bool IsShower() const { return fShower; }
120  bool IsShowerTrack() const { return fShowerTrack; }
121  bool IsShowerCone() const { return fShowerCone; }
122  bool IsTrack() const { return fTrack; }
123  bool IsUndetermined() const { return !fTrack and !fShower; }
124 
125  int TrackConeSize() const { return (int)fForwardConeTracks.size() - (int)fBackwardConeTracks.size(); }
126  bool ShowerTrackCandidate() const { return TrackConeSize() > 5; }
127  const std::vector<int>& ShowerTracks() const { return fShowerTracks; }
128  const std::vector<int>& ForwardConeTracks() const { return fForwardConeTracks; }
129 
130  int ConeSize() const { return (int)fForwardSpacePoints.size() - (int)fBackwardSpacePoints.size(); }
131  int ForwardSpacePoints() const { return fForwardSpacePoints.size(); }
132  int NumCylinderSpacePoints() const { return fCylinderSpacePoints.size(); }
133  double CylinderSpacePointRatio() const { return (double)fCylinderSpacePoints.size()/(double)fSpacePoints.size(); }
134  int NumSphereSpacePoints() const { return fSphereSpacePoints.size(); }
135  //double SphereSpacePointDensity() const { return (double)fSphereSpacePoints.size()/((double)fSpacePoints.size()); }
136  double SphereSpacePointDensity(double scale) const { return (double)fSphereSpacePoints.size()/(4*TMath::Pi()*TMath::Power((scale*fLength/2.),3)/3.); }
137  double IsolationSpacePointDistance() const { std::vector<double> distances;
138  std::transform(fIsolationSpacePoints.begin(), fIsolationSpacePoints.end(), std::back_inserter(distances), [](const std::pair<int,double>& p){return p.second;});
139  return TMath::Mean(distances.begin(), distances.end()); }
140 
141  private:
142 
143  int fID;
144  TVector3 fVertex;
145  TVector3 fEnd;
146  double fLength;
147  TVector3 fVertexDir;
148  TVector3 fDirection;
149  std::vector<art::Ptr<recob::Hit> > fHits;
150  std::vector<art::Ptr<recob::SpacePoint> > fSpacePoints;
151 
152  std::vector<int> fForwardConeTracks;
153  std::vector<int> fBackwardConeTracks;
154  std::vector<int> fShowerTracks;
155 
156  std::vector<int> fForwardSpacePoints;
157  std::vector<int> fBackwardSpacePoints;
158  std::vector<int> fCylinderSpacePoints;
159  std::vector<int> fSphereSpacePoints;
160  std::map<int,double> fIsolationSpacePoints;
161 
162  bool fShower;
165  bool fTrack;
166 
167 };
168 
170  public:
171 
173 
175  void reconfigure(fhicl::ParameterSet const& pset);
176 
177  std::vector<art::Ptr<recob::Hit> > SelectShowerHits(int event,
179  const std::vector<art::Ptr<recob::Track> >& tracks,
180  const std::vector<art::Ptr<recob::SpacePoint> >& spacePoints,
181  const art::FindManyP<recob::Hit>& fmht,
182  const art::FindManyP<recob::Track>& fmth,
184  const art::FindManyP<recob::Track>& fmtsp);
185 
186  private:
187 
189  std::vector<int> InitialTrackLikeSegment(std::map<int,std::unique_ptr<ReconTrack> >& reconTracks);
190 
192  TVector3 Gradient(const std::vector<TVector3>& points, const std::unique_ptr<TVector3>& dir);
193 
195  TVector3 Gradient(const art::Ptr<recob::Track>& track);
196 
198  TVector3 Gradient(const std::vector<art::Ptr<recob::SpacePoint> >& spacePoints);
199 
202  TVector3 ProjPoint(const TVector3& point, const TVector3& direction, const TVector3& origin = TVector3(0,0,0));
203 
205  TVector3 SpacePointPos(const art::Ptr<recob::SpacePoint>& spacePoint);
206 
208  double SpacePointsRMS(const std::vector<art::Ptr<recob::SpacePoint> >& spacePoints);
209 
210  // Parameters
211  int fDebug;
212 
213  // Properties
214  double fConeAngle;
216 
217  // Cuts
219  double fCylinderCut;
221 
222 };
223 
224 #endif
225 
226 
227 
228 
229 
230 
231 
232 
233 
234 
235 
236 
237 
238  /* // --------------------------- OLD (late 2015) ------------------------------- */
239 
240  /* public: */
241 
242  /* /// Takes specific previously reconstructed quantites and removes hits which are considered track-like */
243  /* /// Returns a vector of hits which are determined to be shower-like */
244  /* std::vector<art::Ptr<recob::Hit> > RemoveTrackHits(const std::vector<art::Ptr<recob::Hit> >& initialHits, */
245  /* const std::vector<art::Ptr<recob::Track> >& tracks, */
246  /* const std::vector<art::Ptr<recob::SpacePoint> >& spacePoints, */
247  /* const std::vector<art::Ptr<recob::Vertex> >& vertices, */
248  /* const art::FindManyP<recob::Track>& fmth, */
249  /* const art::FindManyP<recob::Track>& fmtsp, */
250  /* const art::FindManyP<recob::Hit>& fmh, */
251  /* int event, */
252  /* int run); */
253 
254  /* /// Uses information from Pandora reconstruction to separate track-like and shower-like hits */
255  /* /// Returns a vector of hits which are determined to be shower-like */
256  /* std::vector<art::Ptr<recob::Hit> > RemoveTrackHits(const std::vector<art::Ptr<recob::Hit> >& hits, */
257  /* const std::vector<art::Ptr<recob::PFParticle> > pfParticles, */
258  /* const art::FindManyP<recob::Cluster>& fmc, */
259  /* const art::FindManyP<recob::Hit>& fmh); */
260 
261 
262  /* private: */
263 
264  /* /// Fill the output container with all the hits not associated with track-like objects */
265  /* std::vector<art::Ptr<recob::Hit> > FillHitsToCluster(const std::vector<art::Ptr<recob::Hit> >& initialHits, */
266  /* const art::FindManyP<recob::Track>& fmt); */
267 
268  /* /// Find the true track most likely associated with this hit */
269  /* int FindTrackID(const art::Ptr<recob::Hit>& hit); */
270 
271  /* /// Find the true track most likely associated with this set of hits */
272  /* int FindTrueTrack(const std::vector<art::Ptr<recob::Hit> >& trackHits); */
273 
274  /* /// Finds the space points surrounding the track but not part of it */
275  /* std::vector<art::Ptr<recob::SpacePoint> > GetSurroundingSpacePoints(const std::vector<art::Ptr<recob::SpacePoint> >& spacePoints, */
276  /* const art::FindManyP<recob::Track>& fmt, */
277  /* unsigned int trackID); */
278 
279  /* /// Look for space points near the track and within a narrow cone */
280  /* std::vector<art::Ptr<recob::SpacePoint> > GetSpacePointsInCone(const std::vector<art::Ptr<recob::SpacePoint> >& spacePoints, */
281  /* const TVector3& trackEnd, */
282  /* const TVector3& trackDirection); */
283 
284  /* /// Determines whether or not a track is actually the start of a shower */
285  /* bool IdentifyShowerLikeTrack(const TVector3& end, */
286  /* const TVector3& direction, */
287  /* const std::vector<art::Ptr<recob::SpacePoint> >& spacePoints); */
288 
289  /* /// Attempt to identify tracks in the event by using the centre of the event */
290  /* void IdentifyTracksFromEventCentre(const std::vector<art::Ptr<recob::Track> >& tracks, */
291  /* const std::vector<art::Ptr<recob::SpacePoint> >& spacePoints, */
292  /* const art::FindManyP<recob::Track>& fmtsp); */
293 
294  /* /// Identifies tracks which start just after previously identified tracks end */
295  /* void IdentifyTracksNearTracks(const std::vector<art::Ptr<recob::Track> >& tracks); */
296 
297  /* /// Identifies hadron-like tracks which originate from near the interaction vertex */
298  /* void IdentifyTracksNearVertex(const art::Ptr<recob::Vertex>& vertex, */
299  /* const std::vector<art::Ptr<recob::Track> >& tracks, */
300  /* const std::vector<art::Ptr<recob::SpacePoint> >& spacePoints, */
301  /* const art::FindManyP<recob::Track>& fmtsp); */
302 
303  /* /// Finds the spread of a set of space points about their central axis */
304  /* double SpacePointSpread(const std::vector<art::Ptr<recob::SpacePoint> >& spacePoints); */
305 
306  /* std::vector<int> fTrackLikeIDs, fShowerLikeIDs; */
307 
308  /* // Configurable parameters */
309  /* double fAngleCut, fDistanceCut, fVertexProximityCut, fTrackProximityCut, fAvTrackHitDistance; */
310 
311  /* art::ServiceHandle<cheat::BackTracker> backtracker; */
312  /* art::ServiceHandle<art::TFileService> tfs; */
313 
314  /* TTree* ftree; */
315  /* double Distance, Angle, Length, AvDistance; */
316  /* int Event, Run, TrackID, pdg, NSpacePoints; */
void SetDirection(TVector3 direction)
const std::vector< art::Ptr< recob::SpacePoint > > & SpacePoints() const
void SetEnd(TVector3 end)
const std::vector< int > & ForwardConeTracks() const
double CylinderSpacePointRatio() const
const std::vector< int > & ShowerTracks() const
void AddForwardSpacePoint(int spacePoint)
Declaration of signal hit object.
void SetLength(double length)
void AddBackwardTrack(int track)
double IsolationSpacePointDistance() const
std::map< int, double > fIsolationSpacePoints
Float_t tmp
Definition: plot.C:37
void SetSpacePoints(std::vector< art::Ptr< recob::SpacePoint > > spacePoints)
void AddForwardTrack(int track)
std::vector< int > fBackwardSpacePoints
auto vector(Vector const &v)
Returns a manipulator which will print the specified array.
Definition: DumpUtils.h:265
Double_t scale
Definition: plot.C:25
void SetVertex(TVector3 vertex)
void hits()
Definition: readHits.C:15
Provides recob::Track data product.
std::vector< int > fForwardConeTracks
std::vector< art::Ptr< recob::SpacePoint > > fSpacePoints
const std::vector< art::Ptr< recob::Hit > > & Hits() const
std::vector< int > fSphereSpacePoints
Declaration of cluster object.
void SetVertexDir(TVector3 vertexDir)
void AddCylinderSpacePoint(int spacePoint)
std::vector< int > fShowerTracks
std::vector< int > fForwardSpacePoints
TVector3 VertexDirection() const
TDirectory * dir
Definition: macro.C:5
double SphereSpacePointDensity(double scale) const
std::vector< art::Ptr< recob::Hit > > fHits
void AddSphereSpacePoint(int spacePoint)
std::vector< evd::details::RawDigitInfo_t >::const_iterator end(RawDigitCacheDataClass const &cache)
void AddIsolationSpacePoint(int spacePoint, double distance)
void SetHits(std::vector< art::Ptr< recob::Hit > > hits)
TVector3 Direction() const
void AddShowerTrack(int track)
Float_t track
Definition: plot.C:34
void AddBackwardSpacePoint(int spacePoint)
std::vector< int > fCylinderSpacePoints
constexpr Point origin()
Returns a origin position with a point of the specified type.
Definition: geo_vectors.h:230
Event finding and building.
std::vector< int > fBackwardConeTracks
vertex reconstruction