Classes | Public Types | Public Member Functions | Private Member Functions | Private Attributes | List of all members
calo::GnocchiCalorimetry Class Reference
Inheritance diagram for calo::GnocchiCalorimetry:
art::EDProducer art::detail::Producer art::detail::LegacyModule art::Modifier art::ModuleBase art::ProductRegistryHelper

Classes

struct  Config
 

Public Types

using Parameters = art::EDProducer::Table< Config >
 
- Public Types inherited from art::EDProducer
using ModuleType = EDProducer
 
using WorkerType = WorkerT< EDProducer >
 
- Public Types inherited from art::detail::Producer
template<typename UserConfig , typename KeysToIgnore = void>
using Table = Modifier::Table< UserConfig, KeysToIgnore >
 
- Public Types inherited from art::Modifier
template<typename UserConfig , typename UserKeysToIgnore = void>
using Table = ProducerTable< UserConfig, detail::ModuleConfig, UserKeysToIgnore >
 

Public Member Functions

 GnocchiCalorimetry (Parameters const &config)
 
void produce (art::Event &evt) override
 
- Public Member Functions inherited from art::EDProducer
 EDProducer (fhicl::ParameterSet const &pset)
 
template<typename Config >
 EDProducer (Table< Config > const &config)
 
std::string workerType () const
 
- Public Member Functions inherited from art::detail::Producer
virtual ~Producer () noexcept
 
 Producer (fhicl::ParameterSet const &)
 
 Producer (Producer const &)=delete
 
 Producer (Producer &&)=delete
 
Produceroperator= (Producer const &)=delete
 
Produceroperator= (Producer &&)=delete
 
void doBeginJob (SharedResources const &resources)
 
void doEndJob ()
 
void doRespondToOpenInputFile (FileBlock const &fb)
 
void doRespondToCloseInputFile (FileBlock const &fb)
 
void doRespondToOpenOutputFiles (FileBlock const &fb)
 
void doRespondToCloseOutputFiles (FileBlock const &fb)
 
bool doBeginRun (RunPrincipal &rp, ModuleContext const &mc)
 
bool doEndRun (RunPrincipal &rp, ModuleContext const &mc)
 
bool doBeginSubRun (SubRunPrincipal &srp, ModuleContext const &mc)
 
bool doEndSubRun (SubRunPrincipal &srp, ModuleContext const &mc)
 
bool doEvent (EventPrincipal &ep, ModuleContext const &mc, std::atomic< std::size_t > &counts_run, std::atomic< std::size_t > &counts_passed, std::atomic< std::size_t > &counts_failed)
 
- Public Member Functions inherited from art::Modifier
 ~Modifier () noexcept
 
 Modifier ()
 
 Modifier (Modifier const &)=delete
 
 Modifier (Modifier &&)=delete
 
Modifieroperator= (Modifier const &)=delete
 
Modifieroperator= (Modifier &&)=delete
 
- Public Member Functions inherited from art::ModuleBase
virtual ~ModuleBase () noexcept
 
 ModuleBase ()
 
ModuleDescription const & moduleDescription () const
 
void setModuleDescription (ModuleDescription const &)
 
std::array< std::vector< ProductInfo >, NumBranchTypes > const & getConsumables () const
 
void sortConsumables (std::string const &current_process_name)
 
template<typename T , BranchType BT>
ViewToken< T > consumesView (InputTag const &tag)
 
template<typename T , BranchType BT>
ViewToken< T > mayConsumeView (InputTag const &tag)
 

Private Member Functions

std::vector< std::vector< unsigned > > OrganizeHits (const std::vector< art::Ptr< recob::Hit >> &hits, const std::vector< const recob::TrackHitMeta * > &thms, const recob::Track &track, unsigned nplanes)
 
std::vector< std::vector< unsigned > > OrganizeHitsIndividual (const std::vector< art::Ptr< recob::Hit >> &hits, const std::vector< const recob::TrackHitMeta * > &thms, const recob::Track &track, unsigned nplanes)
 
std::vector< std::vector< unsigned > > OrganizeHitsSnippets (const std::vector< art::Ptr< recob::Hit >> &hits, const std::vector< const recob::TrackHitMeta * > &thms, const recob::Track &track, unsigned nplanes)
 
bool HitIsValid (const art::Ptr< recob::Hit > hit, const recob::TrackHitMeta *thm, const recob::Track &track)
 
geo::Point_t GetLocation (const recob::Track &track, const art::Ptr< recob::Hit > hit, const recob::TrackHitMeta *meta)
 
geo::Point_t GetLocationAtWires (const recob::Track &track, const art::Ptr< recob::Hit > hit, const recob::TrackHitMeta *meta)
 
geo::Point_t WireToTrajectoryPosition (const geo::Point_t &loc, const geo::TPCID &tpc)
 
geo::Point_t TrajectoryToWirePosition (const geo::Point_t &loc, const geo::TPCID &tpc)
 
double GetPitch (const recob::Track &track, const art::Ptr< recob::Hit > hit, const recob::TrackHitMeta *meta)
 
double GetCharge (const art::Ptr< recob::Hit > hit)
 
double GetEfield (const detinfo::DetectorPropertiesData &dprop, const recob::Track &track, const art::Ptr< recob::Hit > hit, const recob::TrackHitMeta *meta)
 

Private Attributes

Config fConfig
 
CalorimetryAlg fCaloAlg
 

Additional Inherited Members

- Static Public Member Functions inherited from art::EDProducer
static void commitEvent (EventPrincipal &ep, Event &e)
 
- Protected Member Functions inherited from art::ModuleBase
ConsumesCollectorconsumesCollector ()
 
template<typename T , BranchType = InEvent>
ProductToken< T > consumes (InputTag const &)
 
template<typename Element , BranchType = InEvent>
ViewToken< Element > consumesView (InputTag const &)
 
template<typename T , BranchType = InEvent>
void consumesMany ()
 
template<typename T , BranchType = InEvent>
ProductToken< T > mayConsume (InputTag const &)
 
template<typename Element , BranchType = InEvent>
ViewToken< Element > mayConsumeView (InputTag const &)
 
template<typename T , BranchType = InEvent>
void mayConsumeMany ()
 

Detailed Description

Definition at line 57 of file GnocchiCalorimetry_module.cc.

Member Typedef Documentation

Definition at line 126 of file GnocchiCalorimetry_module.cc.

Constructor & Destructor Documentation

calo::GnocchiCalorimetry::GnocchiCalorimetry ( Parameters const &  config)
explicit

Definition at line 160 of file GnocchiCalorimetry_module.cc.

160  :
161  EDProducer{param},
162  fConfig(param()),
164 {
165  produces< std::vector<anab::Calorimetry> >();
166  produces< art::Assns<recob::Track, anab::Calorimetry> >();
167 
168 }
EDProducer(fhicl::ParameterSet const &pset)
Definition: EDProducer.h:20
fhicl::Table< calo::CalorimetryAlg::Config > CalorimetryAlgConfig

Member Function Documentation

double calo::GnocchiCalorimetry::GetCharge ( const art::Ptr< recob::Hit hit)
private

Definition at line 540 of file GnocchiCalorimetry_module.cc.

540  {
541  switch (fConfig.ChargeMethod()) {
543  return hit->Integral();
545  return hit->PeakAmplitude();
547  return hit->SummedADC();
548  default:
549  return 0.;
550  }
551  return 0.;
552 }
float Integral() const
Integral under the calibrated signal waveform of the hit, in tick x ADC units.
Definition: Hit.h:224
float PeakAmplitude() const
The estimated amplitude of the hit at its peak, in ADC units.
Definition: Hit.h:221
float SummedADC() const
The sum of calibrated ADC counts of the hit (0. by default)
Definition: Hit.h:223
double calo::GnocchiCalorimetry::GetEfield ( const detinfo::DetectorPropertiesData dprop,
const recob::Track track,
const art::Ptr< recob::Hit hit,
const recob::TrackHitMeta meta 
)
private

Definition at line 554 of file GnocchiCalorimetry_module.cc.

554  {
555  auto const* sce = lar::providerFrom<spacecharge::SpaceChargeService>();
556 
557  double EField = dprop.Efield();
558  if (sce->EnableSimEfieldSCE() && fConfig.FieldDistortionEfield()) {
559  // Gets relative E field Distortions
560  geo::Vector_t EFieldOffsets = sce->GetEfieldOffsets(GetLocation(track, hit, meta));
561  // Add 1 in X direction as this is the direction of the drift field
562  EFieldOffsets = EFieldOffsets + geo::Vector_t{1, 0, 0};
563  // Convert to Absolute E Field from relative
564  EFieldOffsets = EField * EFieldOffsets;
565  // We only care about the magnitude for recombination
566  EField = EFieldOffsets.r();
567  }
568  return EField;
569 
570 }
ROOT::Math::DisplacementVector3D< ROOT::Math::Cartesian3D< double >, ROOT::Math::GlobalCoordinateSystemTag > Vector_t
Type for representation of momenta in 3D space.
Definition: geo_vectors.h:164
double Efield(unsigned int planegap=0) const
kV/cm
geo::Point_t GetLocation(const recob::Track &track, const art::Ptr< recob::Hit > hit, const recob::TrackHitMeta *meta)
geo::Point_t calo::GnocchiCalorimetry::GetLocation ( const recob::Track track,
const art::Ptr< recob::Hit hit,
const recob::TrackHitMeta meta 
)
private

Definition at line 442 of file GnocchiCalorimetry_module.cc.

442  {
443  geo::Point_t loc = track.LocationAtPoint(meta->Index());
445 }
geo::WireID WireID() const
Definition: Hit.h:233
Point_t const & LocationAtPoint(size_t i) const
Definition: Track.h:126
geo::Point_t WireToTrajectoryPosition(const geo::Point_t &loc, const geo::TPCID &tpc)
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:184
unsigned int Index() const
Hit index along the track trajectory.
Definition: TrackHitMeta.h:55
geo::Point_t calo::GnocchiCalorimetry::GetLocationAtWires ( const recob::Track track,
const art::Ptr< recob::Hit hit,
const recob::TrackHitMeta meta 
)
private

Definition at line 464 of file GnocchiCalorimetry_module.cc.

464  {
465  geo::Point_t loc = track.LocationAtPoint(meta->Index());
467 }
geo::WireID WireID() const
Definition: Hit.h:233
Point_t const & LocationAtPoint(size_t i) const
Definition: Track.h:126
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:184
geo::Point_t TrajectoryToWirePosition(const geo::Point_t &loc, const geo::TPCID &tpc)
unsigned int Index() const
Hit index along the track trajectory.
Definition: TrackHitMeta.h:55
double calo::GnocchiCalorimetry::GetPitch ( const recob::Track track,
const art::Ptr< recob::Hit hit,
const recob::TrackHitMeta meta 
)
private

Definition at line 489 of file GnocchiCalorimetry_module.cc.

489  {
491  auto const* sce = lar::providerFrom<spacecharge::SpaceChargeService>();
492 
493  double angleToVert = geom->WireAngleToVertical(hit->View(), hit->WireID().TPC, hit->WireID().Cryostat) - 0.5*::util::pi<>();
494 
496 
497  // "dir" should be the direction that the wires see. If the track already has the field
498  // distortion corrections applied, then we need to de-apply them to get the direction as
499  // seen by the wire planes
500  if (sce->EnableCalSpatialSCE() && fConfig.FieldDistortion() && fConfig.TrackIsFieldDistortionCorrected()) {
501  geo::Point_t loc = track.LocationAtPoint(meta->Index());
502 
503  // compute the dir of the track trajectory
504  geo::Vector_t track_dir = track.DirectionAtPoint(meta->Index());
505  geo::Point_t loc_mdx = loc - track_dir * (geom->WirePitch(hit->View()) / 2.);
506  geo::Point_t loc_pdx = loc + track_dir * (geom->WirePitch(hit->View()) / 2.);
507 
508  loc_mdx = TrajectoryToWirePosition(loc_mdx, hit->WireID());
509  loc_pdx = TrajectoryToWirePosition(loc_pdx, hit->WireID());
510 
511  // Direction at wires
512  dir = (loc_pdx - loc_mdx) / (loc_mdx - loc_pdx).r();
513  }
514  // If there is no space charge or the track is not yet corrected, then the dir
515  // is the track is what we want
516  else {
517  dir = track.DirectionAtPoint(meta->Index());
518  }
519 
520  double cosgamma = std::abs(std::sin(angleToVert)*dir.Y() + std::cos(angleToVert)*dir.Z());
521  double pitch;
522  if (cosgamma) {
523  pitch = geom->WirePitch(hit->View())/cosgamma;
524  }
525  else {
526  pitch = 0.;
527  }
528 
529  // now take the pitch computed on the wires and correct it back to the particle trajectory
530  geo::Point_t loc_w = GetLocationAtWires(track, hit, meta);
531 
532  geo::Point_t locw_pdx_traj = WireToTrajectoryPosition(loc_w + pitch*dir, hit->WireID());
533  geo::Point_t loc = WireToTrajectoryPosition(loc_w, hit->WireID());
534 
535  pitch = (locw_pdx_traj - loc).R();
536 
537  return pitch;
538 }
geo::WireID WireID() const
Definition: Hit.h:233
Point_t const & LocationAtPoint(size_t i) const
Definition: Track.h:126
CryostatID_t Cryostat
Index of cryostat.
Definition: geo_types.h:212
geo::View_t View() const
View for the plane of the hit.
Definition: Hit.h:232
ROOT::Math::DisplacementVector3D< ROOT::Math::Cartesian3D< double >, ROOT::Math::GlobalCoordinateSystemTag > Vector_t
Type for representation of momenta in 3D space.
Definition: geo_vectors.h:164
string dir
geo::Length_t WirePitch(geo::PlaneID const &planeid) const
Returns the distance between two consecutive wires.
geo::Point_t GetLocationAtWires(const recob::Track &track, const art::Ptr< recob::Hit > hit, const recob::TrackHitMeta *meta)
T abs(T value)
geo::Point_t WireToTrajectoryPosition(const geo::Point_t &loc, const geo::TPCID &tpc)
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:184
geo::Point_t TrajectoryToWirePosition(const geo::Point_t &loc, const geo::TPCID &tpc)
unsigned int Index() const
Hit index along the track trajectory.
Definition: TrackHitMeta.h:55
Vector_t DirectionAtPoint(size_t i) const
Definition: Track.h:134
TPCID_t TPC
Index of the TPC within its cryostat.
Definition: geo_types.h:406
double WireAngleToVertical(geo::View_t view, geo::TPCID const &tpcid) const
Returns the angle of the wires in the specified view from vertical.
bool calo::GnocchiCalorimetry::HitIsValid ( const art::Ptr< recob::Hit hit,
const recob::TrackHitMeta thm,
const recob::Track track 
)
private

Definition at line 436 of file GnocchiCalorimetry_module.cc.

436  {
437  if (thm->Index() == int_max_as_unsigned_int) return false;
438  if (!track.HasValidPoint(thm->Index())) return false;
439  return true;
440 }
bool HasValidPoint(size_t i) const
Definition: Track.h:111
unsigned int Index() const
Hit index along the track trajectory.
Definition: TrackHitMeta.h:55
std::vector< std::vector< unsigned > > calo::GnocchiCalorimetry::OrganizeHits ( const std::vector< art::Ptr< recob::Hit >> &  hits,
const std::vector< const recob::TrackHitMeta * > &  thms,
const recob::Track track,
unsigned  nplanes 
)
private

Definition at line 350 of file GnocchiCalorimetry_module.cc.

352  {
353  // charge is computed per hit -- we organize hits indivudally
355  return OrganizeHitsIndividual(hits, thms, track, nplanes);
356  }
357  // charge is computed per snippet -- we organize hits by snippet
358  else {
359  return OrganizeHitsSnippets(hits, thms, track, nplanes);
360  }
361 }
std::vector< std::vector< unsigned > > OrganizeHitsSnippets(const std::vector< art::Ptr< recob::Hit >> &hits, const std::vector< const recob::TrackHitMeta * > &thms, const recob::Track &track, unsigned nplanes)
std::vector< std::vector< unsigned > > OrganizeHitsIndividual(const std::vector< art::Ptr< recob::Hit >> &hits, const std::vector< const recob::TrackHitMeta * > &thms, const recob::Track &track, unsigned nplanes)
std::vector< std::vector< unsigned > > calo::GnocchiCalorimetry::OrganizeHitsIndividual ( const std::vector< art::Ptr< recob::Hit >> &  hits,
const std::vector< const recob::TrackHitMeta * > &  thms,
const recob::Track track,
unsigned  nplanes 
)
private

Definition at line 363 of file GnocchiCalorimetry_module.cc.

365  {
366  std::vector<std::vector<unsigned>> ret(nplanes);
367  for (unsigned i = 0; i < hits.size(); i++) {
368  if (HitIsValid(hits[i], thms[i], track)) {
369  ret[hits[i]->WireID().Plane].push_back(i);
370  }
371  }
372 
373  return ret;
374 }
bool HitIsValid(const art::Ptr< recob::Hit > hit, const recob::TrackHitMeta *thm, const recob::Track &track)
std::vector< std::vector< unsigned > > calo::GnocchiCalorimetry::OrganizeHitsSnippets ( const std::vector< art::Ptr< recob::Hit >> &  hits,
const std::vector< const recob::TrackHitMeta * > &  thms,
const recob::Track track,
unsigned  nplanes 
)
private

Definition at line 376 of file GnocchiCalorimetry_module.cc.

378  {
379  // In this case, we need to only accept one hit in each snippet
380  // Snippets are counted by the Start, End, and Wire. If all these are the same for a hit, then they are on the same snippet.
381  //
382  // If there are multiple valid hits on the same snippet, we need a way to pick the best one.
383  // (TODO: find a good way). The current method is to take the one with the highest charge integral.
384  struct HitIdentifier {
385  int startTick;
386  int endTick;
387  int wire;
388  float integral;
389 
390  // construct
391  explicit HitIdentifier(const recob::Hit &hit):
392  startTick(hit.StartTick()),
393  endTick(hit.EndTick()),
394  wire(hit.WireID().Wire),
395  integral(hit.Integral())
396  {}
397 
398  // Defines whether two hits are on the same snippet
399  inline bool operator==(const HitIdentifier& rhs) const {
400  return startTick == rhs.startTick && endTick == rhs.endTick && wire == rhs.wire;
401  }
402 
403  // Defines which hit to pick between two both on the same snippet
404  inline bool operator>(const HitIdentifier& rhs) const {
405  return integral > rhs.integral;
406  }
407  };
408 
409  std::vector<std::vector<unsigned>> ret(nplanes);
410  std::vector<std::vector<HitIdentifier>> hit_idents(nplanes);
411  for (unsigned i = 0; i < hits.size(); i++) {
412  if (HitIsValid(hits[i], thms[i], track)) {
413  HitIdentifier this_ident(*hits[i]);
414 
415  // check if we have found a hit on this snippet before
416  bool found_snippet = false;
417  for (unsigned j = 0; j < ret[hits[i]->WireID().Plane].size(); j++) {
418  if (this_ident == hit_idents[hits[i]->WireID().Plane][j]) {
419  found_snippet = true;
420  if (this_ident > hit_idents[hits[i]->WireID().Plane][j]) {
421  ret[hits[i]->WireID().Plane][j] = i;
422  hit_idents[hits[i]->WireID().Plane][j] = this_ident;
423  }
424  break;
425  }
426  }
427  if (!found_snippet) {
428  ret[hits[i]->WireID().Plane].push_back(i);
429  hit_idents[hits[i]->WireID().Plane].push_back(this_ident);
430  }
431  }
432  }
433  return ret;
434 }
bool operator>(ScheduleID const left, ScheduleID const right) noexcept
Definition: ScheduleID.cc:53
bool HitIsValid(const art::Ptr< recob::Hit > hit, const recob::TrackHitMeta *thm, const recob::Track &track)
IDparameter< geo::WireID > WireID
Member type of validated geo::WireID parameter.
Detector simulation of raw signals on wires.
2D representation of charge deposited in the TDC/wire plane
Definition: Hit.h:48
bool operator==(ModuleKeyAndType const &a, ModuleKeyAndType const &b) noexcept
void calo::GnocchiCalorimetry::produce ( art::Event evt)
overridevirtual

Implements art::EDProducer.

Definition at line 170 of file GnocchiCalorimetry_module.cc.

170  {
171  // Get services
173  auto const clock_data = art::ServiceHandle<detinfo::DetectorClocksService const>()->DataFor(evt);
174  auto const det_prop =
176 
177  size_t nplanes = geom->Nplanes();
178 
179  // Define output collections
180  std::unique_ptr< std::vector<anab::Calorimetry> > outputCalo(new std::vector<anab::Calorimetry>);
181  std::unique_ptr< art::Assns<recob::Track, anab::Calorimetry> > outputCaloAssn(new art::Assns<recob::Track, anab::Calorimetry>);
182 
183  // collect input
184  art::Handle< std::vector<recob::Track> > trackListHandle;
185  std::vector<art::Ptr<recob::Track> > tracklist;
186  if (evt.getByLabel(fConfig.TrackModuleLabel(), trackListHandle)) {
187  art::fill_ptr_vector(tracklist, trackListHandle);
188  }
189 
190  // get the label to collect this hits
192  art::FindManyP<recob::Hit, recob::TrackHitMeta> fmHits(trackListHandle, evt, hitLabel);
193 
194  // must be valid if the T0 module label is non-empty
195  art::FindManyP<anab::T0> fmT0s(trackListHandle, evt, fConfig.T0ModuleLabel());
196 
197  // iterate over all the tracks
198  for (unsigned trk_i = 0; trk_i < tracklist.size(); trk_i++) {
199  const recob::Track &track = *tracklist[trk_i];
200 
201  // collect input for this track
202  const std::vector<art::Ptr<recob::Hit>> &hits = fmHits.at(trk_i);
203  const std::vector<const recob::TrackHitMeta *> &thms = fmHits.data(trk_i);
204 
205  double T0 = 0;
206  if (fConfig.T0ModuleLabel().size()) {
207  const std::vector<art::Ptr<anab::T0>> &this_t0s = fmT0s.at(trk_i);
208  if (this_t0s.size()) T0 = this_t0s.at(0)->Time();
209  }
210 
211  // organize the hits by plane
212  std::vector<std::vector<unsigned>> hit_indices = OrganizeHits(hits, thms, track, nplanes);
213 
214  for (unsigned plane_i = 0; plane_i < nplanes; plane_i++) {
215 
216  float kinetic_energy = 0.;
217  std::vector<float> dEdxs;
218  std::vector<float> dQdxs;
219  std::vector<float> resranges;
220  std::vector<float> deadwireresranges;
221  float range = 0.;
222  std::vector<float> pitches;
223  std::vector<geo::Point_t> xyzs;
224  std::vector<size_t> tp_indices;
225  geo::PlaneID plane;
226 
227  // setup the plane ID
228  plane.Plane = plane_i;
229  plane.TPC = 0; // arbitrary -- tracks can cross TPC boundaries
230  plane.Cryostat = fConfig.Cryostat();
231  plane.isValid = true;
232 
233  std::vector<float> lengths;
234  for (unsigned hit_i = 0; hit_i < hit_indices[plane_i].size(); hit_i++) {
235  unsigned hit_index = hit_indices[plane_i][hit_i];
236 
237  // Get the location of this point
238  geo::Point_t location = GetLocation(track, hits[hit_index], thms[hit_index]);
239 
240  // Get the pitch
241  double pitch = GetPitch(track, hits[hit_index], thms[hit_index]);
242 
243  // And the charge
244  double charge = GetCharge(hits[hit_index]);
245 
246  // Get the EField
247  double EField = GetEfield(det_prop, track, hits[hit_index], thms[hit_index]);
248 
249  double dQdx = charge / pitch;
250 
251  // turn into dEdx
253  fCaloAlg.dEdx_AMP(clock_data, det_prop, dQdx, hits[hit_index]->PeakTime(), hits[hit_index]->WireID().Plane, T0, EField) : \
254  fCaloAlg.dEdx_AREA(clock_data, det_prop, dQdx, hits[hit_index]->PeakTime(), hits[hit_index]->WireID().Plane, T0, EField);
255 
256  // save the length between each pair of hits
257  if (xyzs.size() == 0) {
258  lengths.push_back(0.);
259  }
260  else {
261  lengths.push_back((location - xyzs.back()).r());
262  }
263 
264  // save stuff
265  dEdxs.push_back(dEdx);
266  dQdxs.push_back(dQdx);
267  pitches.push_back(pitch);
268  xyzs.push_back(location);
269  kinetic_energy += dEdx * pitch;
270 
271  // TODO: FIXME
272  // It seems weird that the "trajectory-point-index" actually is the
273  // index of the hit... is this a bug in the documentation
274  // of anab::Calorimetry?
275  //
276  // i.e. -- I think this piece of code should actually be:
277  // tp_indices.push_back(thms[hit_index]->Index());
278  tp_indices.push_back(hits[hit_index].key());
279 
280 
281  } // end iterate over hits
282 
283 
284  // turn the lengths vector into a residual-range vector and total length
285  if (lengths.size() > 1) {
286  range = std::accumulate(lengths.begin(), lengths.end(), 0.);
287 
288  // check the direction that the hits are going in the track:
289  // upstream (end-start) or downstream (start-end)
290  bool is_downstream = \
291  (track.Trajectory().Start() - xyzs[0]).r() + (track.Trajectory().End() - xyzs.back()).r() <
292  (track.Trajectory().End() - xyzs[0]).r() + (track.Trajectory().Start() - xyzs.back()).r();
293 
294  resranges.resize(lengths.size());
295  if (is_downstream) {
296  resranges[lengths.size() - 1] = lengths.back() / 2.;
297  for (int i_len = lengths.size() - 2; i_len >= 0; i_len --) {
298  resranges[i_len] = resranges[i_len+1] + lengths[i_len+1];
299  }
300  }
301  else {
302  resranges[0] = lengths[1] / 2.;
303  for (unsigned i_len = 1; i_len < lengths.size(); i_len ++) {
304  resranges[i_len] = resranges[i_len-1] + lengths[i_len];
305  }
306  }
307  }
308 
309  // save the Calorimetry output
310  //
311  // Bogus if less than two hits on this plane
312  if (lengths.size() > 1) {
313  outputCalo->push_back(anab::Calorimetry(kinetic_energy,
314  dEdxs,
315  dQdxs,
316  resranges,
317  deadwireresranges,
318  range,
319  pitches,
320  xyzs,
321  tp_indices,
322  plane));
323  }
324  else {
325  outputCalo->push_back(anab::Calorimetry(util::kBogusD,
326  {},
327  {},
328  {},
329  {},
331  {},
332  {},
333  {},
334  plane));
335  }
336 
337  util::CreateAssn(*this, evt, *outputCalo, tracklist[trk_i], *outputCaloAssn);
338 
339  } // end iterate over planes
340 
341  } // end iterate over tracks
342 
343  evt.put(std::move(outputCalo));
344  evt.put(std::move(outputCaloAssn));
345 
346  return;
347 
348 }
std::string string
Definition: nybbler.cc:12
const recob::TrackTrajectory & Trajectory() const
Access to the stored recob::TrackTrajectory.
Definition: Track.h:98
The data type to uniquely identify a Plane.
Definition: geo_types.h:472
bool isValid
Whether this ID points to a valid element.
Definition: geo_types.h:211
CryostatID_t Cryostat
Index of cryostat.
Definition: geo_types.h:212
fhicl::Atom< std::string > TrackModuleLabel
fhicl::Atom< std::string > T0ModuleLabel
bool getByLabel(std::string const &label, std::string const &instance, Handle< PROD > &result) const
Definition: DataViewImpl.h:633
double dEdx(float dqdx, float Efield)
Definition: doAna.cpp:21
unsigned int Nplanes(unsigned int tpc=0, unsigned int cstat=0) const
Returns the total number of wire planes in the specified TPC.
double dEdx_AMP(detinfo::DetectorClocksData const &clock_data, detinfo::DetectorPropertiesData const &det_prop, recob::Hit const &hit, double pitch, double T0=0) const
def key(type, name=None)
Definition: graph.py:13
def move(depos, offset)
Definition: depos.py:107
double GetEfield(const detinfo::DetectorPropertiesData &dprop, const recob::Track &track, const art::Ptr< recob::Hit > hit, const recob::TrackHitMeta *meta)
ProductID put(std::unique_ptr< PROD > &&edp, std::string const &instance={})
Definition: DataViewImpl.h:686
bool CreateAssn(PRODUCER const &prod, art::Event &evt, std::vector< T > const &a, art::Ptr< U > const &b, art::Assns< U, T > &assn, std::string a_instance, size_t indx=UINT_MAX)
Creates a single one-to-one association.
PlaneID_t Plane
Index of the plane within its TPC.
Definition: geo_types.h:493
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:184
std::vector< std::vector< unsigned > > OrganizeHits(const std::vector< art::Ptr< recob::Hit >> &hits, const std::vector< const recob::TrackHitMeta * > &thms, const recob::Track &track, unsigned nplanes)
Point_t const & End() const
Returns the position of the last valid point of the trajectory [cm].
double dEdx_AREA(detinfo::DetectorClocksData const &clock_data, detinfo::DetectorPropertiesData const &det_prop, recob::Hit const &hit, double pitch, double T0=0) const
double GetPitch(const recob::Track &track, const art::Ptr< recob::Hit > hit, const recob::TrackHitMeta *meta)
double GetCharge(const art::Ptr< recob::Hit > hit)
constexpr double kBogusD
obviously bogus double value
geo::Point_t GetLocation(const recob::Track &track, const art::Ptr< recob::Hit > hit, const recob::TrackHitMeta *meta)
void fill_ptr_vector(std::vector< Ptr< T >> &ptrs, H const &h)
Definition: Ptr.h:297
TPCID_t TPC
Index of the TPC within its cryostat.
Definition: geo_types.h:406
fhicl::Atom< std::string > AssocHitModuleLabel
Point_t const & Start() const
Returns the position of the first valid point of the trajectory [cm].
Track from a non-cascading particle.A recob::Track consists of a recob::TrackTrajectory, plus additional members relevant for a "fitted" track:
Definition: Track.h:49
geo::Point_t calo::GnocchiCalorimetry::TrajectoryToWirePosition ( const geo::Point_t loc,
const geo::TPCID tpc 
)
private

Definition at line 469 of file GnocchiCalorimetry_module.cc.

469  {
470  auto const* sce = lar::providerFrom<spacecharge::SpaceChargeService>();
472 
473  geo::Point_t ret = loc;
474 
475  if (sce->EnableCalSpatialSCE() && fConfig.FieldDistortion()) {
476  // Returned X is the drift -- multiply by the drift direction to undo this
477  int corr = geom->TPC(tpc.TPC).DriftDir()[0];
478 
479  geo::Vector_t offset = sce->GetPosOffsets(ret);
480 
481  ret.SetX(ret.X() + corr * fConfig.FieldDistortionCorrectionXSign() * offset.X());
482  ret.SetY(ret.Y() + offset.Y());
483  ret.SetZ(ret.Z() + offset.Z());
484  }
485 
486  return ret;
487 }
ROOT::Math::DisplacementVector3D< ROOT::Math::Cartesian3D< double >, ROOT::Math::GlobalCoordinateSystemTag > Vector_t
Type for representation of momenta in 3D space.
Definition: geo_vectors.h:164
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:184
TPCID_t TPC
Index of the TPC within its cryostat.
Definition: geo_types.h:406
geo::Point_t calo::GnocchiCalorimetry::WireToTrajectoryPosition ( const geo::Point_t loc,
const geo::TPCID tpc 
)
private

Definition at line 447 of file GnocchiCalorimetry_module.cc.

447  {
448  auto const* sce = lar::providerFrom<spacecharge::SpaceChargeService>();
449 
450  geo::Point_t ret = loc;
451 
452  if (sce->EnableCalSpatialSCE() && fConfig.FieldDistortion()) {
453  geo::Vector_t offset = sce->GetCalPosOffsets(ret, tpc.TPC);
454 
455  ret.SetX(ret.X() + fConfig.FieldDistortionCorrectionXSign() * offset.X());
456  ret.SetY(ret.Y() + offset.Y());
457  ret.SetZ(ret.Z() + offset.Z());
458  }
459 
460  return ret;
461 
462 }
ROOT::Math::DisplacementVector3D< ROOT::Math::Cartesian3D< double >, ROOT::Math::GlobalCoordinateSystemTag > Vector_t
Type for representation of momenta in 3D space.
Definition: geo_vectors.h:164
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:184
TPCID_t TPC
Index of the TPC within its cryostat.
Definition: geo_types.h:406

Member Data Documentation

CalorimetryAlg calo::GnocchiCalorimetry::fCaloAlg
private

Definition at line 135 of file GnocchiCalorimetry_module.cc.

Config calo::GnocchiCalorimetry::fConfig
private

Definition at line 134 of file GnocchiCalorimetry_module.cc.


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