LArSoft  v09_90_00
Liquid Argon Software toolkit - https://larsoft.org/
LArPandoraModularShowerCreation_module.cc
Go to the documentation of this file.
1 //################################################################################
2 //### Name: LArPandoraModularShowerCreation ###
3 //### Author: Dominic Barker and Ed Tyley (e.tyley@sheffield.ac.uk) ###
4 //### Date: 15.05.19 ###
5 //### Description: Generic Shower Charaterisation module which allows the ###
6 //### the user choose which tool to calculate shower metrics. ###
7 //### For a complete shower the tools must define (use the exact ###
8 //### name) the following metrics in the shower property holder: ###
9 //### ShowerStartPosition ###
10 //### ShowerDirection ###
11 //### ShowerEnergy ###
12 //### ShowerdEdx ###
13 //### ShowerLength ###
14 //### ShowerOpeningAngle ###
15 //################################################################################
16 
17 //Framework includes
21 
22 //LArSoft includes
30 
31 namespace reco::shower {
32  class LArPandoraModularShowerCreation;
33 }
34 
35 //Class
36 
38 public:
40 
41 private:
42  void produce(art::Event& evt);
43 
44  //This function returns the art::Ptr to the data object InstanceName.
45  //In the background it uses the PtrMaker which requires the element index of
46  //the unique ptr (iter).
47  template <class T>
48  art::Ptr<T> GetProducedElementPtr(const std::string& InstanceName,
49  const reco::shower::ShowerElementHolder& ShowerEleHolder,
50  const int& iter = -1);
51 
52  //fcl object names
53  unsigned int fNumPlanes;
56  const int fVerbose;
57  const bool fUseAllParticles;
58 
59  //tool tags which calculate the characteristics of the shower
60  const std::string fShowerStartPositionLabel;
61  const std::string fShowerDirectionLabel;
62  const std::string fShowerEnergyLabel;
63  const std::string fShowerLengthLabel;
64  const std::string fShowerOpeningAngleLabel;
65  const std::string fShowerdEdxLabel;
66  const std::string fShowerBestPlaneLabel;
67 
68  //fcl tools
69  std::vector<std::unique_ptr<ShowerRecoTools::IShowerTool>> fShowerTools;
70  std::vector<std::string> fShowerToolNames;
71 
72  //map to the unique ptrs to
74 
75  // Required services
77 };
78 
79 //This function returns the art::Ptr to the data object InstanceName.
80 //In the background it uses the PtrMaker which requires the element index of
81 //the unique ptr (iter).
82 template <class T>
84  const std::string& InstanceName,
85  const reco::shower::ShowerElementHolder& ShowerEleHolder,
86  const int& iter)
87 {
88 
89  bool check_element = ShowerEleHolder.CheckElement(InstanceName);
90  if (!check_element) {
91  throw cet::exception("LArPandoraModularShowerCreation")
92  << "To get a element that does not exist" << std::endl;
93  }
94 
95  bool check_ptr = uniqueproducerPtrs.CheckUniqueProduerPtr(InstanceName);
96  if (!check_ptr) {
97  throw cet::exception("LArPandoraModularShowerCreation")
98  << "Tried to get a ptr that does not exist" << std::endl;
99  }
100 
101  //Get the number of the shower we are on.
102  int index;
103  if (iter != -1) { index = iter; }
104  else {
105  index = ShowerEleHolder.GetShowerNumber();
106  }
107 
108  //Make the ptr
109  art::Ptr<T> artptr = uniqueproducerPtrs.GetArtPtr<T>(InstanceName, index);
110  return artptr;
111 }
112 
114  fhicl::ParameterSet const& pset)
115  : EDProducer{pset}
116  , fPFParticleLabel(pset.get<art::InputTag>("PFParticleLabel"))
117  , fAllowPartialShowers(pset.get<bool>("AllowPartialShowers"))
118  , fVerbose(pset.get<int>("Verbose", 0))
119  , fUseAllParticles(pset.get<bool>("UseAllParticles", false))
120  , fShowerStartPositionLabel(pset.get<std::string>("ShowerStartPositionLabel"))
121  , fShowerDirectionLabel(pset.get<std::string>("ShowerDirectionLabel"))
122  , fShowerEnergyLabel(pset.get<std::string>("ShowerEnergyLabel"))
123  , fShowerLengthLabel(pset.get<std::string>("ShowerLengthLabel"))
124  , fShowerOpeningAngleLabel(pset.get<std::string>("ShowerOpeningAngleLabel"))
125  , fShowerdEdxLabel(pset.get<std::string>("ShowerdEdxLabel"))
126  , fShowerBestPlaneLabel(pset.get<std::string>("ShowerBestPlaneLabel"))
127 {
128  //Intialise the tools
129  auto tool_psets = pset.get<std::vector<fhicl::ParameterSet>>("ShowerFinderTools");
130  for (auto& tool_pset : tool_psets) {
131 
132  const std::string tool_name(tool_pset.get<std::string>("tool_type"));
133  // If the PFPaticle label is not set for the tool, make it use the one for the module
134  // Note we also need to set the label in the Alg, via the base tool
135  if (!tool_pset.has_key("PFParticleLabel")) {
136 
137  // I cannot pass an art::InputTag as it is mangled, so lets make a string instead
138  const std::string PFParticleLabelString(fPFParticleLabel.label() + ":" +
139  fPFParticleLabel.instance() + ":" +
141 
142  tool_pset.put<std::string>("PFParticleLabel", PFParticleLabelString);
143  fhicl::ParameterSet base_pset = tool_pset.get<fhicl::ParameterSet>("BaseTools");
144  fhicl::ParameterSet alg_pset = base_pset.get<fhicl::ParameterSet>("LArPandoraShowerAlg");
145  alg_pset.put<std::string>("PFParticleLabel", PFParticleLabelString);
146  base_pset.put_or_replace<fhicl::ParameterSet>("LArPandoraShowerAlg", alg_pset);
147  tool_pset.put_or_replace<fhicl::ParameterSet>("BaseTools", base_pset);
148 
149  if (tool_pset.has_key("LArPandoraShowerCheatingAlg")) {
150  fhicl::ParameterSet cheat_alg_pset =
151  tool_pset.get<fhicl::ParameterSet>("LArPandoraShowerCheatingAlg");
152  cheat_alg_pset.put<std::string>("PFParticleLabel", PFParticleLabelString);
153  cheat_alg_pset.put_or_replace<fhicl::ParameterSet>("LArPandoraShowerAlg", alg_pset);
154  tool_pset.put_or_replace<fhicl::ParameterSet>("LArPandoraShowerCheatingAlg",
155  cheat_alg_pset);
156  }
157  }
158 
159  // If we have not explicitly set verboseness for a given tool, use global level
160  if (!tool_pset.has_key("Verbose")) { tool_pset.put<int>("Verbose", fVerbose); }
161 
162  fShowerTools.push_back(art::make_tool<ShowerRecoTools::IShowerTool>(tool_pset));
163  fShowerToolNames.push_back(tool_name);
164 
165  fNumPlanes = fGeom->Nplanes();
166  }
167 
168  // Initialise the EDProducer ptr in the tools
169  std::vector<std::string> SetupTools;
170  for (unsigned int i = 0; i < fShowerTools.size(); ++i) {
171  if (std::find(SetupTools.begin(), SetupTools.end(), fShowerToolNames[i]) != SetupTools.end()) {
172  continue;
173  }
174  fShowerTools[i]->SetPtr(&producesCollector());
175  fShowerTools[i]->InitaliseProducerPtr(uniqueproducerPtrs);
176  fShowerTools[i]->InitialiseProducers();
177  }
178 
179  //Initialise the other paramters.
180 
181  produces<std::vector<recob::Shower>>();
182  produces<art::Assns<recob::Shower, recob::Hit>>();
183  produces<art::Assns<recob::Shower, recob::Cluster>>();
184  produces<art::Assns<recob::Shower, recob::SpacePoint>>();
185  produces<art::Assns<recob::Shower, recob::PFParticle>>();
186 
187  // Output -- showers and associations with hits and clusters
188  uniqueproducerPtrs.SetShowerUniqueProduerPtr(type<std::vector<recob::Shower>>(), "shower");
190  "clusterAssociationsbase");
192  "hitAssociationsbase");
194  "spShowerAssociationsbase");
196  "pfShowerAssociationsbase");
197 
199 }
200 
202 {
203 
204  //Ptr makers for the products
206  reco::shower::ShowerElementHolder showerEleHolder;
207 
208  //Get the PFParticles
209  auto const pfpHandle = evt.getValidHandle<std::vector<recob::PFParticle>>(fPFParticleLabel);
210  std::vector<art::Ptr<recob::PFParticle>> pfps;
211  art::fill_ptr_vector(pfps, pfpHandle);
212 
213  //Handle to access the pandora hits assans
214  auto const clusterHandle = evt.getValidHandle<std::vector<recob::Cluster>>(fPFParticleLabel);
215 
216  //Get the assoications to hits, clusters and spacespoints
217  const art::FindManyP<recob::Hit>& fmh =
218  showerEleHolder.GetFindManyP<recob::Hit>(clusterHandle, evt, fPFParticleLabel);
219  const art::FindManyP<recob::Cluster>& fmcp =
220  showerEleHolder.GetFindManyP<recob::Cluster>(pfpHandle, evt, fPFParticleLabel);
221  const art::FindManyP<recob::SpacePoint>& fmspp =
222  showerEleHolder.GetFindManyP<recob::SpacePoint>(pfpHandle, evt, fPFParticleLabel);
223 
224  //Holder to pass to the functions, contains the 6 properties of the shower
225  // - Start Poistion
226  // - Direction
227  // - Initial Track
228  // - Initial Track Hits
229  // - Energy
230  // - dEdx
231  // - Length
232  // - Opening Angle
233 
234  int shower_iter = 0;
235  //Loop of the pf particles
236  for (auto const& pfp : pfps) {
237 
238  //Update the shower iterator
239  showerEleHolder.SetShowerNumber(shower_iter);
240 
241  //loop only over showers unless otherwise specified
242  if (!fUseAllParticles && pfp->PdgCode() != 11 && pfp->PdgCode() != 22) continue;
243 
244  //Get the associated hits,clusters and spacepoints
245  const std::vector<art::Ptr<recob::Cluster>> showerClusters = fmcp.at(pfp.key());
246  const std::vector<art::Ptr<recob::SpacePoint>> showerSpacePoints = fmspp.at(pfp.key());
247 
248  // Check the pfp has at least 1 cluster (i.e. not a pfp neutrino)
249  if (!showerClusters.size()) continue;
250 
251  if (fVerbose > 1)
252  mf::LogInfo("LArPandoraModularShowerCreation")
253  << "Running on shower: " << shower_iter << std::endl;
254 
255  //Calculate the shower properties
256  //Loop over the shower tools
257  int err = 0;
258  for (unsigned int i = 0; i < fShowerTools.size(); i++) {
259 
260  //Calculate the metric
261  if (fVerbose > 1)
262  mf::LogInfo("LArPandoraModularShowerCreation")
263  << "Running shower tool: " << fShowerToolNames[i] << std::endl;
264  std::string evd_disp_append = fShowerToolNames[i] + "_iteration" + std::to_string(0) + "_" +
265  this->moduleDescription().moduleLabel();
266 
267  err = fShowerTools[i]->RunShowerTool(pfp, evt, showerEleHolder, evd_disp_append);
268 
269  if (err && fVerbose) {
270  mf::LogError("LArPandoraModularShowerCreation")
271  << "Error in shower tool: " << fShowerToolNames[i] << " with code: " << err << std::endl;
272  }
273  }
274 
275  //If we are are not allowing partial shower check all of the things
276  if (!fAllowPartialShowers) {
277  // If we recieved an error call from a tool return;
278 
279  // Check everything we need is in the shower element holder
280  if (!showerEleHolder.CheckElement(fShowerStartPositionLabel)) {
281  if (fVerbose)
282  mf::LogError("LArPandoraModularShowerCreation")
283  << "The start position is not set in the element holder. bailing" << std::endl;
284  continue;
285  }
286  if (!showerEleHolder.CheckElement(fShowerDirectionLabel)) {
287  if (fVerbose)
288  mf::LogError("LArPandoraModularShowerCreation")
289  << "The direction is not set in the element holder. bailing" << std::endl;
290  continue;
291  }
292  if (!showerEleHolder.CheckElement(fShowerEnergyLabel)) {
293  if (fVerbose)
294  mf::LogError("LArPandoraModularShowerCreation")
295  << "The energy is not set in the element holder. bailing" << std::endl;
296  continue;
297  }
298  if (!showerEleHolder.CheckElement(fShowerdEdxLabel)) {
299  if (fVerbose)
300  mf::LogError("LArPandoraModularShowerCreation")
301  << "The dEdx is not set in the element holder. bailing" << std::endl;
302  continue;
303  }
304  if (!showerEleHolder.CheckElement(fShowerBestPlaneLabel)) {
305  if (fVerbose)
306  mf::LogError("LArPandoraModularShowerCreation")
307  << "The BestPlane is not set in the element holder. bailing" << std::endl;
308  continue;
309  }
310  if (!showerEleHolder.CheckElement(fShowerLengthLabel)) {
311  if (fVerbose)
312  mf::LogError("LArPandoraModularShowerCreation")
313  << "The length is not set in the element holder. bailing" << std::endl;
314  continue;
315  }
316  if (!showerEleHolder.CheckElement(fShowerOpeningAngleLabel)) {
317  if (fVerbose)
318  mf::LogError("LArPandoraModularShowerCreation")
319  << "The opening angle is not set in the element holder. bailing" << std::endl;
320  continue;
321  }
322 
323  //Check All of the products that have been asked to be checked.
324  bool elements_are_set = showerEleHolder.CheckAllElementTags();
325  if (!elements_are_set) {
326  if (fVerbose)
327  mf::LogError("LArPandoraModularShowerCreation")
328  << "Not all the elements in the property holder which should be set are not. Bailing. "
329  << std::endl;
330  continue;
331  }
332 
334  bool producers_are_set = uniqueproducerPtrs.CheckAllProducedElements(showerEleHolder);
335  if (!producers_are_set) {
336  if (fVerbose)
337  mf::LogError("LArPandoraModularShowerCreation")
338  << "Not all the elements in the property holder which are produced are not set. "
339  "Bailing. "
340  << std::endl;
341  continue;
342  }
343  }
344 
345  //Get the properties
346  geo::Point_t ShowerStartPosition(-999, -999, -999);
347  geo::Vector_t ShowerDirection(-999, -999, -999);
348  std::vector<double> ShowerEnergy(fNumPlanes, -999);
349  std::vector<double> ShowerdEdx(fNumPlanes, -999);
350  int BestPlane(-999);
351  double ShowerLength(-999);
352  double ShowerOpeningAngle(-999);
353 
354  geo::Point_t ShowerStartPositionErr(-999, -999, -999);
355  geo::Vector_t ShowerDirectionErr(-999, -999, -999);
356  std::vector<double> ShowerEnergyErr(fNumPlanes, -999);
357  std::vector<double> ShowerdEdxErr(fNumPlanes, -999);
358 
359  err = 0;
360  if (showerEleHolder.CheckElement(fShowerStartPositionLabel))
361  err += showerEleHolder.GetElementAndError(
362  fShowerStartPositionLabel, ShowerStartPosition, ShowerStartPositionErr);
363  if (showerEleHolder.CheckElement(fShowerDirectionLabel))
364  err += showerEleHolder.GetElementAndError(
365  fShowerDirectionLabel, ShowerDirection, ShowerDirectionErr);
366  if (showerEleHolder.CheckElement(fShowerEnergyLabel))
367  err += showerEleHolder.GetElementAndError(fShowerEnergyLabel, ShowerEnergy, ShowerEnergyErr);
368  if (showerEleHolder.CheckElement(fShowerdEdxLabel))
369  err += showerEleHolder.GetElementAndError(fShowerdEdxLabel, ShowerdEdx, ShowerdEdxErr);
370  if (showerEleHolder.CheckElement(fShowerBestPlaneLabel))
371  err += showerEleHolder.GetElement(fShowerBestPlaneLabel, BestPlane);
372  if (showerEleHolder.CheckElement(fShowerLengthLabel))
373  err += showerEleHolder.GetElement(fShowerLengthLabel, ShowerLength);
374  if (showerEleHolder.CheckElement(fShowerOpeningAngleLabel))
375  err += showerEleHolder.GetElement(fShowerOpeningAngleLabel, ShowerOpeningAngle);
376 
377  if (err) {
378  throw cet::exception("LArPandoraModularShowerCreation")
379  << "Error in LArPandoraModularShowerCreation Module. A Check on a shower property failed "
380  << std::endl;
381  }
382 
383  if (fVerbose > 1) {
384  //Check the shower
385  std::cout << "Shower Vertex: X:" << ShowerStartPosition.X()
386  << " Y: " << ShowerStartPosition.Y() << " Z: " << ShowerStartPosition.Z()
387  << std::endl;
388  std::cout << "Shower Direction: X:" << ShowerDirection.X() << " Y: " << ShowerDirection.Y()
389  << " Z: " << ShowerDirection.Z() << std::endl;
390  std::cout << "Shower dEdx:";
391  for (unsigned int i = 0; i < fNumPlanes; i++) {
392  std::cout << " Plane " << i << ": " << ShowerdEdx.at(i);
393  }
394  std::cout << std::endl;
395  std::cout << "Shower Energy:";
396  for (unsigned int i = 0; i < fNumPlanes; i++) {
397  std::cout << " Plane " << i << ": " << ShowerEnergy.at(i);
398  }
399  std::cout << std::endl;
400  std::cout << "Shower Best Plane: " << BestPlane << std::endl;
401  std::cout << "Shower Length: " << ShowerLength << std::endl;
402  std::cout << "Shower Opening Angle: " << ShowerOpeningAngle << std::endl;
403 
404  //Print what has been created in the shower
405  showerEleHolder.PrintElements();
406  }
407 
408  if (ShowerdEdx.size() != fNumPlanes) {
409  throw cet::exception("LArPandoraModularShowerCreation")
410  << "dEdx vector is wrong size: " << ShowerdEdx.size()
411  << " compared to Nplanes: " << fNumPlanes << std::endl;
412  }
413  if (ShowerEnergy.size() != fNumPlanes) {
414  throw cet::exception("LArPandoraModularShowerCreation")
415  << "Energy vector is wrong size: " << ShowerEnergy.size()
416  << " compared to Nplanes: " << fNumPlanes << std::endl;
417  }
418 
419  //Make the shower
420  using namespace geo::vect;
421  recob::Shower shower(convertTo<TVector3>(ShowerDirection),
422  convertTo<TVector3>(ShowerDirectionErr),
423  convertTo<TVector3>(ShowerStartPosition),
424  convertTo<TVector3>(ShowerDirectionErr),
425  ShowerEnergy,
426  ShowerEnergyErr,
427  ShowerdEdx,
428  ShowerdEdxErr,
429  BestPlane,
431  ShowerLength,
432  ShowerOpeningAngle);
433  showerEleHolder.SetElement(shower, "shower");
434  ++shower_iter;
435  art::Ptr<recob::Shower> ShowerPtr =
436  this->GetProducedElementPtr<recob::Shower>("shower", showerEleHolder);
437 
438  //Associate the pfparticle
440  ShowerPtr, pfp, "pfShowerAssociationsbase");
441 
442  //Add the hits for each "cluster"
443  for (auto const& cluster : showerClusters) {
444 
445  //Associate the clusters
446  std::vector<art::Ptr<recob::Hit>> ClusterHits = fmh.at(cluster.key());
448  ShowerPtr, cluster, "clusterAssociationsbase");
449 
450  //Associate the hits
451  for (auto const& hit : ClusterHits) {
453  ShowerPtr, hit, "hitAssociationsbase");
454  }
455  }
456 
457  //Associate the spacepoints
458  for (auto const& sp : showerSpacePoints) {
460  ShowerPtr, sp, "spShowerAssociationsbase");
461  }
462 
463  //Loop over the tool data products and add them.
464  uniqueproducerPtrs.AddDataProducts(showerEleHolder);
465 
466  //AddAssociations
467  int assn_err = 0;
468  for (auto const& fShowerTool : fShowerTools) {
469  //AddAssociations
470  assn_err += fShowerTool->AddAssociations(pfp, evt, showerEleHolder);
471  }
472  if (!fAllowPartialShowers && assn_err > 0) {
473  if (fVerbose)
474  mf::LogError("LArPandoraModularShowerCreation")
475  << "A association failed and not allowing partial showers. The association will not be "
476  "added to the event "
477  << std::endl;
478  }
479 
480  //Reset the showerproperty holder.
481  showerEleHolder.ClearShower();
482  }
483 
484  //Put everything in the event.
486 
487  //Reset the ptrs to the data products
489 }
490 
ROOT::Math::DisplacementVector3D< ROOT::Math::Cartesian3D< double >, ROOT::Math::GlobalCoordinateSystemTag > Vector_t
Type for representation of momenta in 3D space.
Definition: geo_vectors.h:160
MaybeLogger_< ELseverityLevel::ELsev_info, false > LogInfo
void AddDataProducts(const reco::shower::ShowerElementHolder &selement_holder)
std::string const & moduleLabel() const
Declaration of signal hit object.
EDProducer(fhicl::ParameterSet const &pset)
Definition: EDProducer.cc:6
void SetElement(T &dataproduct, const std::string &Name, bool checktag=false)
std::string const & instance() const noexcept
Definition: InputTag.cc:85
constexpr int kBogusI
obviously bogus integer value
Set of hits with a 2D structure.
Definition: Cluster.h:69
double ShowerEnergy(const ShowerStruct3D &ss3)
Definition: TCShower.cxx:3899
MaybeLogger_< ELseverityLevel::ELsev_error, false > LogError
Cluster finding and building.
std::string const & process() const noexcept
Definition: InputTag.cc:91
bool CheckUniqueProduerPtr(const std::string &Name) const
art::Ptr< T > GetProducedElementPtr(const std::string &InstanceName, const reco::shower::ShowerElementHolder &ShowerEleHolder, const int &iter=-1)
std::string const & label() const noexcept
Definition: InputTag.cc:79
Utilities to manipulate geometry vectors.The utilities include generic vector interface facilities al...
#define DEFINE_ART_MODULE(klass)
Definition: ModuleMacros.h:65
void AddSingle(A &a, B &b, const std::string &Name)
decltype(auto) constexpr to_string(T &&obj)
ADL-aware version of std::to_string.
T get(std::string const &key) const
Definition: ParameterSet.h:314
bool CheckElement(const std::string &Name) const
int GetElement(const std::string &Name, T &Element) const
Declaration of cluster object.
Detector simulation of raw signals on wires.
ProducesCollector & producesCollector() noexcept
bool CheckAllProducedElements(reco::shower::ShowerElementHolder &selement_holder) const
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
ValidHandle< PROD > getValidHandle(InputTag const &tag) const
Utility object to perform functions of association.
unsigned int Nplanes(TPCID const &tpcid=tpc_zero) const
Returns the total number of planes in the specified TPC.
Definition: GeometryCore.h:977
std::vector< std::unique_ptr< ShowerRecoTools::IShowerTool > > fShowerTools
int GetElementAndError(const std::string &Name, T &Element, T2 &ElementErr) const
2D representation of charge deposited in the TDC/wire plane
Definition: Hit.h:46
art::Ptr< T > GetArtPtr(const std::string &Name, const int &iter) const
TCEvent evt
Definition: DataStructs.cxx:8
void fill_ptr_vector(std::vector< Ptr< T >> &ptrs, H const &h)
Definition: Ptr.h:306
void put_or_replace(std::string const &key)
int SetShowerUniqueProduerPtr(type< T >, const std::string &Name, const std::string &Instance="")
const art::FindManyP< T1 > & GetFindManyP(const art::ValidHandle< std::vector< T2 >> &handle, const art::Event &evt, const art::InputTag &moduleTag)
Definition: fwd.h:26
void put(std::string const &key)
cet::coded_exception< error, detail::translate > exception
Definition: exception.h:33
ModuleDescription const & moduleDescription() const
Definition: ModuleBase.cc:13