LArSoft  v10_04_05
Liquid Argon Software toolkit - https://larsoft.org/
GraphCluster_module.cc
Go to the documentation of this file.
1 
11 #include <vector>
12 
13 #ifdef __ROOTCLING__
14 namespace art {
15  class EDProducer;
16  class Event;
17  class PtrVector;
18  class Ptr;
19  class ServiceHandle;
20 }
21 
22 #include "fhiclcpp/fwd.h"
23 
24 namespace recob {
25  class Hit;
26 }
27 #else
34 #endif
35 
37 //
38 // GraphCluster class
39 //
40 // andrzej.szelc@yale.edu
41 // ellen.klein@yale.edu
42 //
43 // This dummy producer is designed to create a hitlist and maybe cluster from EVD input
45 
46 // Framework includes
48 #include "fhiclcpp/ParameterSet.h"
49 
50 // LArSoft Includes
58 
59 namespace util {
60  class PxPoint;
61 }
62 
63 namespace evd {
64 
65  class InfoTransfer;
66 
67  class GraphCluster : public art::EDProducer {
68 
69  public:
70  explicit GraphCluster(fhicl::ParameterSet const&);
71  void produce(art::Event& evt);
72 
73  private:
75 
76  void GetStartEndHits(unsigned int plane, recob::Hit* starthit, recob::Hit* endhit);
77  void GetStartEndHits(unsigned int plane);
78 
79  void GetHitList(unsigned int plane, art::PtrVector<recob::Hit>& ptrhitlist);
80 
81  std::vector<util::PxLine> GetSeedLines();
82 
83  unsigned int fNPlanes;
84 
85  int TestFlag;
86  int fRun;
87  int fSubRun;
88  int fEvent;
89 
90  std::vector<recob::Hit*> starthit;
91  std::vector<recob::Hit*> endhit;
92 
93  std::vector<util::PxLine> startendpoints;
94  }; // class GraphCluster
95 
96  //-------------------------------------------------
97  GraphCluster::GraphCluster(fhicl::ParameterSet const& pset)
98  : EDProducer{pset}, fGClAlg(pset.get<fhicl::ParameterSet>("GraphClusterAlg"))
99  {
100  produces<std::vector<recob::Cluster>>();
101  produces<art::Assns<recob::Cluster, recob::Hit>>();
102  produces<std::vector<art::PtrVector<recob::Cluster>>>();
103 
104  auto const& wireReadoutGeom = art::ServiceHandle<geo::WireReadout>()->Get();
105  fNPlanes = wireReadoutGeom.Nplanes();
106  starthit.resize(fNPlanes);
107  endhit.resize(fNPlanes);
108 
109  startendpoints.resize(fNPlanes);
110  }
111 
112  //
113  //-------------------------------------------------
117  {
118 
119  std::unique_ptr<std::vector<recob::Cluster>> Graphcol(new std::vector<recob::Cluster>);
120  std::unique_ptr<art::Assns<recob::Cluster, recob::Hit>> hassn(
122  std::unique_ptr<std::vector<art::PtrVector<recob::Cluster>>> classn(
124 
125  auto const& wireReadoutGeom = art::ServiceHandle<geo::WireReadout>()->Get();
126 
127  // check if evt and run numbers check out, etc...
128  if (fGClAlg.CheckValidity(evt) == -1) { return; }
129 
130  for (unsigned int ip = 0; ip < fNPlanes; ip++) {
131  startendpoints[ip] = util::PxLine(); //assign empty PxLine
132  }
133 
134  std::vector<art::PtrVector<recob::Hit>> hitlist;
135  hitlist.resize(fNPlanes);
136 
137  for (unsigned int ip = 0; ip < fNPlanes; ip++) {
138 
139  fGClAlg.GetHitListAndEndPoints(ip, hitlist[ip], startendpoints[ip]);
140 
141  if (hitlist[ip].size() == 0) continue;
142 
143  if (hitlist[ip].size() > 0 && !(TestFlag == -1)) //old event or transfer not ready
144  {
145  double swterror = 0., ewterror = 0.;
146 
147  if (startendpoints[ip].w0 == 0) swterror = 999;
148 
149  if (startendpoints[ip].t1 == 0) ewterror = 999;
150 
151  std::cout << " clustering @ " << startendpoints[ip].w0 << " +/- " << swterror << " "
152  << startendpoints[ip].t0 << " +/- " << swterror << " " << startendpoints[ip].w1
153  << " +/- " << ewterror << " " << startendpoints[ip].t1 << " +/- " << ewterror
154  << std::endl;
155 
156  // this is all we can do easily without getting the full-blown
157  // ClusterParamsAlg (that means lareventdisplay has to depend on larreco!)
158  lar::util::StatCollector<float> integral, summedADC;
159  for (art::Ptr<recob::Hit> const& hit : hitlist[ip]) {
160  integral.add(hit->Integral());
161  summedADC.add(hit->ROISummedADC());
162  } // for
163 
164  // get the plane ID from the first hit
165  geo::PlaneID planeID = hitlist[ip].front()->WireID().planeID();
166  Graphcol->emplace_back(startendpoints[ip].w0,
167  swterror,
168  startendpoints[ip].t0,
169  swterror,
170  0., // start_charge
171  0., // start_angle
172  0., // start_opening
173  startendpoints[ip].w1,
174  ewterror,
175  startendpoints[ip].t1,
176  ewterror,
177  0., // end_charge
178  0., // end_angle
179  0., // end_opening
180  integral.Sum(), // integral
181  integral.RMS(), // integral_stddev
182  summedADC.Sum(), // summedADC
183  summedADC.RMS(), // summedADC_stddev
184  hitlist[ip].size(), // n_hits
185  0., // multiple_hit_density
186  0., // width
187  ip,
188  wireReadoutGeom.Plane({planeID.asTPCID(), ip}).View(),
189  planeID,
191 
192  // associate the hits to this cluster
193  util::CreateAssn(evt, *Graphcol, hitlist[ip], *hassn);
194  }
195 
196  } // end of loop on planes
197 
199  cvec.reserve(fNPlanes);
200 
201  for (unsigned int ip = 0; ip < fNPlanes; ip++) {
202  art::ProductID aid = evt.getProductID<std::vector<recob::Cluster>>();
203  art::Ptr<recob::Cluster> aptr(aid, ip, evt.productGetter(aid));
204  cvec.push_back(aptr);
205  }
206 
207  classn->push_back(cvec);
208 
209  evt.put(std::move(Graphcol));
210  evt.put(std::move(hassn));
211  evt.put(std::move(classn));
212 
213  return;
214  } // end of produce
215 
217 
218 } //end of evd namespace
code to link reconstructed objects back to the MC truth information
void reserve(size_type n)
Definition: PtrVector.h:337
Namespace for general, non-LArSoft-specific utilities.
Definition: PIDAAlg.h:26
TTree * t1
Definition: plottest35.C:26
Reconstruction base classes.
ProductID getProductID(std::string const &instance_name="") const
Declaration of signal hit object.
std::vector< recob::Hit * > starthit
The data type to uniquely identify a Plane.
Definition: geo_types.h:364
GraphClusterAlg fGClAlg
cout<< "Opened file "<< fin<< " ixs= "<< ixs<< endl;if(ixs==0) hhh=(TH1F *) fff-> Get("h1")
Definition: AddMC.C:8
Classes gathering simple statistics.
Weight_t RMS() const
Returns the root mean square.
static const SentryArgument_t Sentry
An instance of the sentry object.
Definition: Cluster.h:174
PutHandle< PROD > put(std::unique_ptr< PROD > &&edp, std::string const &instance={})
Definition: Event.h:77
decltype(auto) constexpr size(T &&obj)
ADL-aware version of std::size.
Definition: StdUtils.h:101
auto vector(Vector const &v)
Returns a manipulator which will print the specified array.
Definition: DumpUtils.h:289
Ptr(H, T) -> Ptr< detail::not_map_vector_t< typename H::element_type >>
EDProductGetter const * productGetter(ProductID const pid) const
std::vector< recob::Hit * > endhit
LArSoft includes.
#define DEFINE_ART_MODULE(klass)
Definition: ModuleMacros.h:65
void push_back(Ptr< U > const &p)
Definition: PtrVector.h:435
Weight_t Sum() const
Returns the weighted sum of the values.
int CheckValidity(art::Event &evt)
std::vector< util::PxLine > startendpoints
Declaration of cluster object.
Definition of data types for geometry description.
Detector simulation of raw signals on wires.
void GetHitListAndEndPoints(unsigned int plane, art::PtrVector< recob::Hit > &ptrhitlist, util::PxLine &startendpoints)
bool CreateAssn(art::Event &evt, std::vector< T > const &a, art::Ptr< U > const &b, art::Assns< U, T > &assn, std::string a_instance, size_t index=UINT_MAX)
Creates a single one-to-one association.
void produce(art::Event &evt)
Utility object to perform functions of association.
Encapsulate the construction of a single detector plane .
Definition: MVAAlg.h:12
2D representation of charge deposited in the TDC/wire plane
Definition: Hit.h:46
TCEvent evt
Definition: DataStructs.cxx:8
Collects statistics on a single quantity (weighted)
constexpr TPCID const & asTPCID() const
Conversion to PlaneID (for convenience of notation).
Definition: geo_types.h:409
void add(Data_t value, Weight_t weight=Weight_t(1.0))
Adds one entry with specified value and weight.