LArSoft  v10_04_05
Liquid Argon Software toolkit - https://larsoft.org/
LineCluster_module.cc
Go to the documentation of this file.
1 
10 // Framework libraries
15 #include "fhiclcpp/ParameterSet.h"
16 
17 // LArSoft includes
22 
23 // ... more includes in the implementation section
24 
25 namespace cluster {
38  class LineCluster : public art::EDProducer {
39  public:
40  explicit LineCluster(fhicl::ParameterSet const& pset);
41 
42  private:
43  void produce(art::Event& evt) override;
44 
45  ClusterCrawlerAlg fCCAlg; // define ClusterCrawlerAlg object
46 
48 
51 
52  }; // class LineCluster
53 
54 } // namespace cluster
55 
56 //******************************************************************************
57 //*** implementation
58 //***
59 
60 // C/C++ standard libraries
61 #include <memory> // std::move()
62 
63 // Framework libraries
67 
68 //LArSoft includes
70 #include "lardata/ArtDataHelper/HitCreator.h" // recob::HitCollectionAssociator
76 
77 namespace cluster {
78 
79  //----------------------------------------------------------------------------
81  : EDProducer{pset}, fCCAlg{pset.get<fhicl::ParameterSet>("ClusterCrawlerAlg")}
82  {
83  fHitFinderLabel = pset.get<art::InputTag>("HitFinderModuleLabel");
84  fDoWireAssns = pset.get<bool>("DoWireAssns", true);
85  fDoRawDigitAssns = pset.get<bool>("DoRawDigitAssns", false);
86 
87  // let HitCollectionAssociator declare that we are going to produce
88  // hits and associations with wires and raw digits
89  // (with no particular product label)
91  producesCollector(), "", fDoWireAssns, fDoRawDigitAssns);
92 
93  produces<std::vector<recob::Cluster>>();
94  produces<std::vector<recob::Vertex>>();
95  produces<std::vector<recob::EndPoint2D>>();
96  produces<art::Assns<recob::Cluster, recob::Hit>>();
97  produces<art::Assns<recob::Cluster, recob::Vertex, unsigned short>>();
98  produces<art::Assns<recob::Cluster, recob::EndPoint2D, unsigned short>>();
99  } // LineCluster::LineCluster()
100 
101  //----------------------------------------------------------------------------
103  {
104  // fetch the wires needed by CCHitFinder
105 
106  // make this accessible to ClusterCrawler_module
107  auto hitVecHandle = evt.getValidHandle<std::vector<recob::Hit>>(fHitFinderLabel);
108 
109  // look for clusters in all planes
110  auto const clock_data =
112  auto const det_prop =
114  fCCAlg.RunCrawler(clock_data, det_prop, *hitVecHandle);
115 
116  auto FinalHits = std::make_unique<std::vector<recob::Hit>>(std::move(fCCAlg.YieldHits()));
117 
118  // shcol contains the hit collection
119  // and its associations to wires and raw digits;
120  // we get the association to raw digits through wire associations
122  std::vector<recob::Cluster> sccol;
123  std::vector<recob::Vertex> sv3col;
124  std::vector<recob::EndPoint2D> sv2col;
125 
126  std::unique_ptr<art::Assns<recob::Cluster, recob::Hit>> hc_assn(
128  std::unique_ptr<art::Assns<recob::Cluster, recob::Vertex, unsigned short>> cv_assn(
130  std::unique_ptr<art::Assns<recob::Cluster, recob::EndPoint2D, unsigned short>> cep_assn(
132 
133  std::vector<ClusterCrawlerAlg::ClusterStore> const& Clusters = fCCAlg.GetClusters();
134 
135  // make EndPoints (aka 2D vertices)
136  std::vector<ClusterCrawlerAlg::VtxStore> const& EndPts = fCCAlg.GetEndPoints();
137  std::vector<unsigned int> indxToIndx(EndPts.size());
138  auto const& wireReadoutGeom = art::ServiceHandle<geo::WireReadout>()->Get();
139  unsigned short vtxID = 0, end, wire, ivx;
140  for (ivx = 0; ivx < EndPts.size(); ++ivx) {
141  if (EndPts[ivx].NClusters == 0) continue;
142  indxToIndx[ivx] = vtxID;
143  ++vtxID;
144  wire = (0.5 + EndPts[ivx].Wire);
145  geo::PlaneID plID = ClusterCrawlerAlg::DecodeCTP(EndPts[ivx].CTP);
146  geo::WireID wID = geo::WireID(plID.Cryostat, plID.TPC, plID.Plane, wire);
147  geo::View_t view = wireReadoutGeom.Plane(wID).View();
148  sv2col.emplace_back((double)EndPts[ivx].Time, // Time
149  wID, // WireID
150  0, // strength - not relevant
151  vtxID, // ID
152  view, // View
153  0); // total charge - not relevant
154  } // iv
155  // convert 2D Vertex vector to unique_ptrs
156  std::unique_ptr<std::vector<recob::EndPoint2D>> v2col(
157  new std::vector<recob::EndPoint2D>(std::move(sv2col)));
158 
159  // make 3D vertices
160  std::vector<ClusterCrawlerAlg::Vtx3Store> const& Vertices = fCCAlg.GetVertices();
161  double xyz[3] = {0, 0, 0};
162  vtxID = 0;
163  for (ClusterCrawlerAlg::Vtx3Store const& vtx3 : Vertices) {
164  // ignore incomplete vertices
165  if (vtx3.Ptr2D[0] < 0) continue;
166  if (vtx3.Ptr2D[1] < 0) continue;
167  if (vtx3.Ptr2D[2] < 0) continue;
168  ++vtxID;
169  xyz[0] = vtx3.X;
170  xyz[1] = vtx3.Y;
171  xyz[2] = vtx3.Z;
172  sv3col.emplace_back(xyz, vtxID);
173  } // 3D vertices
174  // convert Vertex vector to unique_ptrs
175  std::unique_ptr<std::vector<recob::Vertex>> v3col(
176  new std::vector<recob::Vertex>(std::move(sv3col)));
177 
178  // make the clusters and associations
179  float sumChg, sumADC;
180  unsigned int clsID = 0, nclhits;
181  for (unsigned int icl = 0; icl < Clusters.size(); ++icl) {
182  ClusterCrawlerAlg::ClusterStore const& clstr = Clusters[icl];
183  if (clstr.ID < 0) continue;
184  ++clsID;
185  sumChg = 0;
186  sumADC = 0;
188  unsigned short plane = planeID.Plane;
189  nclhits = clstr.tclhits.size();
190  std::vector<unsigned int> clsHitIndices;
191  // correct the hit indices to refer to the valid hits that were just added
192  for (unsigned int itt = 0; itt < nclhits; ++itt) {
193  unsigned int iht = clstr.tclhits[itt];
194  recob::Hit const& hit = FinalHits->at(iht);
195  sumChg += hit.Integral();
196  sumADC += hit.ROISummedADC();
197  } // itt
198  // get the wire, plane from a hit
199  unsigned int iht = clstr.tclhits[0];
200 
201  geo::View_t view = FinalHits->at(iht).View();
202  sccol.emplace_back((float)clstr.BeginWir, // Start wire
203  0, // sigma start wire
204  clstr.BeginTim, // start tick
205  0, // sigma start tick
206  clstr.BeginChg, // start charge
207  clstr.BeginAng, // start angle
208  0, // start opening angle (0 for line-like clusters)
209  (float)clstr.EndWir, // end wire
210  0, // sigma end wire
211  clstr.EndTim, // end tick
212  0, // sigma end tick
213  clstr.EndChg, // end charge
214  clstr.EndAng, // end angle
215  0, // end opening angle (0 for line-like clusters)
216  sumChg, // integral
217  0, // sigma integral
218  sumADC, // summed ADC
219  0, // sigma summed ADC
220  nclhits, // n hits
221  0, // wires over hits
222  0, // width (0 for line-like clusters)
223  clsID, // ID
224  view, // view
225  planeID, // plane
226  recob::Cluster::Sentry // sentry
227  );
228  // make the cluster - hit association
229  if (!util::CreateAssn(
230  evt, *hc_assn, sccol.size() - 1, clstr.tclhits.begin(), clstr.tclhits.end())) {
232  << "Failed to associate hit " << iht << " with cluster " << icl;
233  } // exception
234  // make the cluster - EndPoint2D and Vertex associations
235  if (clstr.BeginVtx >= 0) {
236  end = 0;
237  if (!util::CreateAssnD(evt, *cep_assn, clsID - 1, indxToIndx[clstr.BeginVtx], end)) {
239  << "Failed to associate cluster " << clsID << " with EndPoint2D " << clstr.BeginVtx;
240  } // exception
241  // See if this endpoint is associated with a 3D vertex
242  unsigned short vtxIndex = 0;
243  for (ClusterCrawlerAlg::Vtx3Store const& vtx3 : Vertices) {
244  // ignore incomplete vertices
245  if (vtx3.Ptr2D[0] < 0) continue;
246  if (vtx3.Ptr2D[1] < 0) continue;
247  if (vtx3.Ptr2D[2] < 0) continue;
248  if (vtx3.Ptr2D[plane] == clstr.BeginVtx) {
249  if (!util::CreateAssnD(evt, *cv_assn, clsID - 1, vtxIndex, end)) {
251  << "Failed to associate cluster " << icl << " with vertex";
252  } // exception
253  break;
254  } // vertex match
255  ++vtxIndex;
256  } // 3D vertices
257  } // clstr.BeginVtx >= 0
258  if (clstr.EndVtx >= 0) {
259  end = 1;
260  if (!util::CreateAssnD(evt, *cep_assn, clsID - 1, indxToIndx[clstr.EndVtx], end)) {
262  << "Failed to associate cluster " << clsID << " with EndPoint2D " << clstr.BeginVtx;
263  } // exception
264  // See if this endpoint is associated with a 3D vertex
265  unsigned short vtxIndex = 0;
266  for (ClusterCrawlerAlg::Vtx3Store const& vtx3 : Vertices) {
267  // ignore incomplete vertices
268  if (vtx3.Ptr2D[0] < 0) continue;
269  if (vtx3.Ptr2D[1] < 0) continue;
270  if (vtx3.Ptr2D[2] < 0) continue;
271  if (vtx3.Ptr2D[plane] == clstr.EndVtx) {
272  if (!util::CreateAssnD(evt, *cv_assn, clsID - 1, vtxIndex, end)) {
274  << "Failed to associate cluster " << icl << " with endpoint";
275  } // exception
276  break;
277  } // vertex match
278  ++vtxIndex;
279  } // 3D vertices
280  } // clstr.BeginVtx >= 0
281  } // icl
282 
283  // convert cluster vector to unique_ptrs
284  std::unique_ptr<std::vector<recob::Cluster>> ccol(
285  new std::vector<recob::Cluster>(std::move(sccol)));
286 
287  shcol.use_hits(std::move(FinalHits));
288 
289  // clean up
291 
292  // move the hit collection and the associations into the event:
293  shcol.put_into(evt);
294  evt.put(std::move(ccol));
295  evt.put(std::move(hc_assn));
296  evt.put(std::move(v2col));
297  evt.put(std::move(v3col));
298  evt.put(std::move(cv_assn));
299  evt.put(std::move(cep_assn));
300 
301  } // LineCluster::produce()
302 
303  //----------------------------------------------------------------------------
305 
306 } // namespace cluster
bool CreateAssnD(art::Event &evt, art::Assns< T, U, D > &assn, size_t first_index, size_t second_index, typename art::Assns< T, U, D >::data_t &&data)
Creates a single one-to-one association with associated data.
void produce(art::Event &evt) override
enum geo::_plane_proj View_t
Enumerate the possible plane projections.
Declaration of signal hit object.
EDProducer(fhicl::ParameterSet const &pset)
Definition: EDProducer.cc:6
The data type to uniquely identify a Plane.
Definition: geo_types.h:364
std::vector< ClusterStore > const & GetClusters() const
Returns a constant reference to the clusters found.
CryostatID_t Cryostat
Index of cryostat.
Definition: geo_types.h:195
cout<< "Opened file "<< fin<< " ixs= "<< ixs<< endl;if(ixs==0) hhh=(TH1F *) fff-> Get("h1")
Definition: AddMC.C:8
float Integral() const
Integral under the calibrated signal waveform of the hit, in tick x ADC units.
Definition: Hit.h:254
Cluster finding and building.
struct of temporary 3D vertices
static void declare_products(art::ProducesCollector &collector, std::string instance_name="", bool doWireAssns=true, bool doRawDigitAssns=true)
Declares the hit products we are going to fill.
Definition: HitCreator.cxx:261
static const SentryArgument_t Sentry
An instance of the sentry object.
Definition: Cluster.h:174
static geo::PlaneID DecodeCTP(CTP_t CTP)
PutHandle< PROD > put(std::unique_ptr< PROD > &&edp, std::string const &instance={})
Definition: Event.h:77
decltype(auto) constexpr end(T &&obj)
ADL-aware version of std::end.
Definition: StdUtils.h:77
Helper functions to create a hit.
std::vector< recob::Hit > && YieldHits()
Returns (and loses) the collection of reconstructed hits.
std::vector< VtxStore > const & GetEndPoints() const
Returns a constant reference to the 2D end points found.
void use_hits(std::unique_ptr< std::vector< recob::Hit >> &&srchits)
Uses the specified collection as data product.
Definition: HitCreator.cxx:521
#define DEFINE_ART_MODULE(klass)
Definition: ModuleMacros.h:65
IDparameter< geo::WireID > WireID
Member type of validated geo::WireID parameter.
void put_into(art::Event &)
Moves the data into the event.
Definition: HitCreator.h:858
art::InputTag fHitFinderLabel
label of module producing input hits
A class handling a collection of hits and its associations.
Definition: HitCreator.h:795
PlaneID_t Plane
Index of the plane within its TPC.
Definition: geo_types.h:373
Declaration of cluster object.
Definition of data types for geometry description.
ClusterCrawlerAlg fCCAlg
Detector simulation of raw signals on wires.
ProducesCollector & producesCollector() noexcept
Produces clusters by ClusterCrawler algorithm.
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.
cet::coded_exception< errors::ErrorCodes, ExceptionDetail::translate > Exception
Definition: Exception.h:66
float ROISummedADC() const
The sum of calibrated ADC counts of the ROI (0. by default)
Definition: Hit.h:246
ValidHandle< PROD > getValidHandle(InputTag const &tag) const
Utility object to perform functions of association.
std::vector< Vtx3Store > const & GetVertices() const
Returns a constant reference to the 3D vertices found.
void RunCrawler(detinfo::DetectorClocksData const &clock_data, detinfo::DetectorPropertiesData const &det_prop, std::vector< recob::Hit > const &srchits)
LineCluster(fhicl::ParameterSet const &pset)
2D representation of charge deposited in the TDC/wire plane
Definition: Hit.h:46
TCEvent evt
Definition: DataStructs.cxx:8
TPCID_t TPC
Index of the TPC within its cryostat.
Definition: geo_types.h:315