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

Public Member Functions

 ShSeg (fhicl::ParameterSet const &p)
 
 ShSeg (ShSeg const &)=delete
 
 ShSeg (ShSeg &&)=delete
 
ShSegoperator= (ShSeg const &)=delete
 
ShSegoperator= (ShSeg &&)=delete
 
void beginJob () override
 
void reconfigure (fhicl::ParameterSet const &p)
 
void produce (art::Event &e) 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

bool InsideFidVol (TLorentzVector const &pvtx) const
 
void FilterOutSmallParts (detinfo::DetectorPropertiesData const &detProp, double r2d, const std::vector< art::Ptr< recob::Hit > > &hits_in, std::vector< art::Ptr< recob::Hit > > &hits_out)
 
bool GetCloseHits (detinfo::DetectorPropertiesData const &detProp, double r2d, const std::vector< art::Ptr< recob::Hit > > &hits_in, std::vector< size_t > &used, std::vector< art::Ptr< recob::Hit > > &hits_out)
 
bool Has (const std::vector< size_t > &v, size_t idx)
 
void CorrOffset (detinfo::DetectorPropertiesData const &detProp, TVector3 &vec, const simb::MCParticle &particle)
 
recob::Track ConvertFrom (pma::Track3D const &src)
 
void Smearel ()
 
void Smearph ()
 
bool BuildSegMC (art::Event &e)
 

Private Attributes

std::vector< pma::Track3D * > fPmatracks
 
double fFidVolCut
 
double fR0
 
double fR1
 
short isdata
 
std::string fHitsModuleLabel
 
calo::CalorimetryAlg fCalorimetryAlg
 
pma::ProjectionMatchingAlg fProjectionMatchingAlg
 
int fTrkindex
 
int fEvNumber
 
double fDedxavg
 
double fDedx
 
TTree * fEvTree
 
TTree * fShTree
 

Additional Inherited Members

- 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 >
 
- 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 61 of file ShSeg_module.cc.

Constructor & Destructor Documentation

dunefd::ShSeg::ShSeg ( fhicl::ParameterSet const &  p)
explicit

Definition at line 128 of file ShSeg_module.cc.

128  : EDProducer{p},
129  fCalorimetryAlg(p.get<fhicl::ParameterSet>("CalorimetryAlg")),
130  fProjectionMatchingAlg(p.get< fhicl::ParameterSet >("ProjectionMatchingAlg"))
131  {
132  fR0 = 0.0; fR1 = 0.0;
133  fTrkindex = 0;
134  fEvNumber = -1; fDedxavg = -1.0; fDedx = -1.0;
135 
136  fPmatracks.clear();
137  this->reconfigure(p);
138 
139  produces< std::vector<recob::Track> >();
140  produces< std::vector<recob::SpacePoint> >();
141  produces< art::Assns<recob::Track, recob::Hit> >();
142  produces< art::Assns<recob::Track, recob::SpacePoint> >();
143  produces< art::Assns<recob::SpacePoint, recob::Hit> >();
144  }
EDProducer(fhicl::ParameterSet const &pset)
Definition: EDProducer.h:20
std::vector< pma::Track3D * > fPmatracks
p
Definition: test.py:223
void reconfigure(fhicl::ParameterSet const &p)
pma::ProjectionMatchingAlg fProjectionMatchingAlg
calo::CalorimetryAlg fCalorimetryAlg
dunefd::ShSeg::ShSeg ( ShSeg const &  )
delete
dunefd::ShSeg::ShSeg ( ShSeg &&  )
delete

Member Function Documentation

void dunefd::ShSeg::beginJob ( )
overridevirtual

Reimplemented from art::EDProducer.

Definition at line 148 of file ShSeg_module.cc.

149 {
151 
152  fEvTree = tfs->make<TTree>("Egtestev", "egtestev");
153  fEvTree->Branch("fEvNumber", &fEvNumber, "fEvNumber/I");
154  fEvTree->Branch("fDedxavg", &fDedxavg, "fDedxavg/D");
155 
156  fShTree = tfs->make<TTree>("Egtestsh", "egtestsh");
157  fShTree->Branch("fDedx", &fDedx, "fDedx/D");
158 }
bool dunefd::ShSeg::BuildSegMC ( art::Event e)
private

Definition at line 284 of file ShSeg_module.cc.

285 {
286  bool result = true;
287 
288  std::vector<art::Ptr<recob::Hit> > hitlist;
289  auto hitListHandle = e.getHandle< std::vector<recob::Hit> >(fHitsModuleLabel);
290  if (hitListHandle)
291  art::fill_ptr_vector(hitlist, hitListHandle);
292 
295  const sim::ParticleList& plist = pi_serv->ParticleList();
296 
297  std::vector<const simb::MCParticle * > primaries = plist.GetPrimaries();
298  if (primaries.size() != 1) return false;
299 
300  const simb::MCParticle* firstel = primaries[0];
301  if ((firstel->PdgCode() != 11) && (firstel->PdgCode() != -11) && (firstel->PdgCode() != 22)) return false;
302 
303  TLorentzVector startingp = primaries[0]->Position();
304  TVector3 primaryvtx(startingp.X(), startingp.Y(), startingp.Z());
305 
306  // pretend that there is vertex which can hide situation in it.
307  // cut hits randomly wth gauss dist (both: el, gammas)
308  Smearel();
309  if ((firstel->PdgCode() == 22) && (firstel->EndProcess() == "conv"))
310  {
311  startingp = primaries[0]->EndPosition();
312  TVector3 startp(startingp.X(), startingp.Y(), startingp.Z());
313  // smear with normal dist (only for gammas)
314  Smearph();
315  }
316 
317  // check fiducial volume
318  if (!InsideFidVol(startingp)) return false;
319 
320  TLorentzVector mom = primaries[0]->Momentum();
321  TVector3 momvec3(mom.Px(), mom.Py(), mom.Pz());
322  TVector3 dir = momvec3 * (1 / momvec3.Mag());
323 
324  TVector3 firstpoint(startingp.X(), startingp.Y(), startingp.Z());
325  TVector3 secondpoint = firstpoint + dir;
326 
327  // offset corrections
328  auto const clockData = art::ServiceHandle<detinfo::DetectorClocksService const>()->DataFor(e);
329  auto const detProp = art::ServiceHandle<detinfo::DetectorPropertiesService const>()->DataFor(e, clockData);
330  CorrOffset(detProp, primaryvtx, *firstel);
331  CorrOffset(detProp, firstpoint, *firstel);
332  CorrOffset(detProp, secondpoint, *firstel);
333 
334  double startvtx[3] = {firstpoint.X(), firstpoint.Y(), firstpoint.Z()};
335 
336  // check if it is inside
337  geo::TPCID tpcid = geom->FindTPCAtPosition(startvtx);
338  if (!tpcid.isValid) return false;
339 
340  // try to build seg, it is based on MC truth
341  size_t tpc = tpcid.TPC;
342  size_t cryo = geom->FindCryostatAtPosition(startvtx);
343  pma::Track3D* iniseg = new pma::Track3D();
344  iniseg->AddNode(detProp, firstpoint, tpc, cryo);
345  iniseg->AddNode(detProp, secondpoint, tpc, cryo);
346 
347  for (size_t h = 0; h < hitlist.size(); h++)
348  if (hitlist[h]->WireID().TPC == tpc)
349  iniseg->push_back(detProp, hitlist[h]);
350 
351  iniseg->MakeProjection();
352  iniseg->SortHits();
353 
354  // only fraction of hits are interesting
355  size_t hi = 0;
356  while (hi < iniseg->size())
357  {
358  pma::Hit3D* hit3d = (*iniseg)[hi];
359 
360  if ((hit3d->GetSegFraction() > 20) || (hit3d->GetSegFraction() < 0))
361  {
362  iniseg->release_at(hi);
363  continue;
364  }
365  hi++;
366  }
367 
368  // check hits size in the track
369  if (iniseg->size() < 5) return false;
370 
371  // 0: geo::kU, 1: geo::kV, 2: geo::kZ
372  const geo::TPCGeo& tpcgeo = geom->GetElement(tpcid);
373  double maxdist = 0.0; size_t bestview = 0;
374 
375  for (size_t view = 0; view < tpcgeo.Nplanes(); ++view)
376  {
377  std::map< size_t, std::vector< double > > ex;
378  iniseg->GetRawdEdxSequence(ex, view);
379 
380  TVector2 proj_i = pma::GetProjectionToPlane(firstpoint, view, tpc, cryo);
381  TVector2 proj_f = pma::GetProjectionToPlane(secondpoint, view, tpc, cryo);
382  double dist = std::sqrt(pma::Dist2(proj_i, proj_f));
383  if ((dist > maxdist) && (ex.size() > 0))
384  {
385  maxdist = dist;
386  bestview = view;
387  }
388  }
389 
390 
391  fPmatracks.push_back(iniseg);
392  /************************************/
393 
394  iniseg->CompleteMissingWires(detProp, bestview);
395  std::map< size_t, std::vector< double > > dedx;
396  iniseg->GetRawdEdxSequence(dedx, bestview);
397  double sumdx = 0.0; fDedxavg = 0.0;
398 
399  double rmin = 1.0e9;
400  for (auto & v: dedx)
401  if (v.second[7] < rmin) rmin = v.second[7];
402 
403  double rmax = rmin + 3.0;
404 
405  for (auto & v: dedx)
406  if (((v.second[7] + fR1) > fR0) && (v.second[7] < 15))
407  {
408  double dE = fCalorimetryAlg.dEdx_AREA(clockData, detProp, v.second[5], v.second[1], bestview);
409  double range = v.second[7] + (fR0 - fR1);
410 
411  v.second[5] = dE;
412  v.second[7] = range;
413 
414  if ((v.second[5] > 0) && (v.second[6] > 0))
415  {
416  fDedx = v.second[5] / v.second[6];
417  fShTree->Fill();
418 
419  if (v.second[7] < rmax)
420  {
421  fDedxavg += v.second[5];
422  sumdx += v.second[6];
423  }
424  }
425  }
426 
427  if (sumdx > 0) fDedxavg = fDedxavg / sumdx;
428  fEvTree->Fill();
429 
430  return result;
431 }
void MakeProjection()
int PdgCode() const
Definition: MCParticle.h:212
CryostatGeo const & GetElement(geo::CryostatID const &cryoid) const
static QCString result
double Dist2(const TVector2 &v1, const TVector2 &v2)
Definition: Utilities.cxx:37
Handle< PROD > getHandle(SelectorBase const &) const
Definition: DataViewImpl.h:382
unsigned int Nplanes() const
Number of planes in this tpc.
Definition: TPCGeo.h:165
bool isValid
Whether this ID points to a valid element.
Definition: geo_types.h:211
Geometry information for a single TPC.
Definition: TPCGeo.h:38
string dir
geo::TPCID FindTPCAtPosition(double const worldLoc[3]) const
Returns the ID of the TPC at specified location.
decltype(auto) constexpr size(T &&obj)
ADL-aware version of std::size.
Definition: StdUtils.h:92
void AddNode(pma::Node3D *node)
pma::Hit3D * release_at(size_t index)
Definition: PmaTrack3D.cxx:343
IDparameter< geo::WireID > WireID
Member type of validated geo::WireID parameter.
std::string EndProcess() const
Definition: MCParticle.h:216
void CorrOffset(detinfo::DetectorPropertiesData const &detProp, TVector3 &vec, const simb::MCParticle &particle)
std::vector< pma::Track3D * > fPmatracks
The data type to uniquely identify a TPC.
Definition: geo_types.h:386
float GetSegFraction() const noexcept
Definition: PmaHit3D.h:143
const sim::ParticleList & ParticleList() const
size_t CompleteMissingWires(detinfo::DetectorPropertiesData const &detProp, unsigned int view)
void push_back(pma::Hit3D *hit)
Definition: PmaTrack3D.h:87
double GetRawdEdxSequence(std::map< size_t, std::vector< double >> &dedx, unsigned int view=geo::kZ, unsigned int skip=0, bool inclDisabled=false) const
constexpr double dist(const TReal *x, const TReal *y, const unsigned int dimension)
TVector2 GetProjectionToPlane(const TVector3 &p, unsigned int plane, unsigned int tpc, unsigned int cryo)
Definition: Utilities.cxx:270
std::string fHitsModuleLabel
double dEdx_AREA(detinfo::DetectorClocksData const &clock_data, detinfo::DetectorPropertiesData const &det_prop, recob::Hit const &hit, double pitch, double T0=0) const
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
size_t size() const
Definition: PmaTrack3D.h:108
geo::CryostatID::CryostatID_t FindCryostatAtPosition(geo::Point_t const &worldLoc) const
Returns the index of the cryostat at specified location.
bool InsideFidVol(TLorentzVector const &pvtx) const
calo::CalorimetryAlg fCalorimetryAlg
recob::Track dunefd::ShSeg::ConvertFrom ( pma::Track3D const &  src)
private

Definition at line 253 of file ShSeg_module.cc.

254 {
255  std::vector< TVector3 > xyz, dircos;
256 
257  for (size_t i = 0; i < src.size(); ++i)
258  if (src[i]->IsEnabled())
259  {
260  xyz.push_back(src[i]->Point3D());
261 
262  if (i < src.size() - 1)
263  {
264  TVector3 dc(src[i + 1]->Point3D());
265  dc -= src[i]->Point3D();
266  dc *= 1.0 / dc.Mag();
267  dircos.push_back(dc);
268  }
269  else dircos.push_back(dircos.back());
270  }
271 
272  if (xyz.size() != dircos.size())
273  {
274  mf::LogError("IniSegReco") << "pma::Track3D to recob::Track conversion problem.";
275  }
278  recob::Track::Flags_t(xyz.size()), false),
280 }
ROOT::Math::SMatrix< Double32_t, 5, 5, ROOT::Math::MatRepSym< Double32_t, 5 > > SMatrixSym55
Definition: TrackingTypes.h:85
std::vector< Vector_t > convertCollToVector(std::vector< Vector > const &coll)
Definition: TrackingTypes.h:58
::fhicl::TupleAs< Point(::geo::Length_t,::geo::Length_t,::geo::Length_t)> Point3D
Atom object for reading a 3D point or vector (centimeters).
TrackTrajectory::Flags_t Flags_t
Definition: Track.h:68
MaybeLogger_< ELseverityLevel::ELsev_error, false > LogError
A trajectory in space reconstructed from hits.
std::vector< Point_t > convertCollToPoint(std::vector< Point > const &coll)
Definition: TrackingTypes.h:55
TrackCollectionProxyElement< TrackCollProxy > Track
Proxy to an element of a proxy collection of recob::Track objects.
Definition: Track.h:1036
void dunefd::ShSeg::CorrOffset ( detinfo::DetectorPropertiesData const &  detProp,
TVector3 &  vec,
const simb::MCParticle particle 
)
private

Definition at line 435 of file ShSeg_module.cc.

437 {
438  float corrt0x = particle.T() * 1.e-3 * detProp.DriftVelocity();
439 
440  float px = vec.X();
441  if (px > 0) px -= corrt0x;
442  else px += corrt0x;
443  vec.SetX(px);
444 }
double T(const int i=0) const
Definition: MCParticle.h:224
void dunefd::ShSeg::FilterOutSmallParts ( detinfo::DetectorPropertiesData const &  detProp,
double  r2d,
const std::vector< art::Ptr< recob::Hit > > &  hits_in,
std::vector< art::Ptr< recob::Hit > > &  hits_out 
)
private

Definition at line 526 of file ShSeg_module.cc.

531 {
532  size_t min_size = hits_in.size() / 5;
533  if (min_size < 3) min_size = 3;
534 
535  std::vector<size_t> used;
536  std::vector< art::Ptr<recob::Hit> > close_hits;
537 
538  while (GetCloseHits(detProp, r2d, hits_in, used, close_hits))
539  {
540  if (close_hits.size() > min_size)
541  for (auto h : close_hits) hits_out.push_back(h);
542  }
543 }
bool GetCloseHits(detinfo::DetectorPropertiesData const &detProp, double r2d, const std::vector< art::Ptr< recob::Hit > > &hits_in, std::vector< size_t > &used, std::vector< art::Ptr< recob::Hit > > &hits_out)
bool dunefd::ShSeg::GetCloseHits ( detinfo::DetectorPropertiesData const &  detProp,
double  r2d,
const std::vector< art::Ptr< recob::Hit > > &  hits_in,
std::vector< size_t > &  used,
std::vector< art::Ptr< recob::Hit > > &  hits_out 
)
private

Definition at line 547 of file ShSeg_module.cc.

553 {
554 
555  hits_out.clear();
556 
557  const double gapMargin = 5.0; // can be changed to f(id_tpc1, id_tpc2)
558  size_t idx = 0;
559 
560  while ((idx < hits_in.size()) && Has(used, idx)) idx++;
561 
562  if (idx < hits_in.size())
563  {
564  hits_out.push_back(hits_in[idx]);
565  used.push_back(idx);
566 
567  double r2d2 = r2d*r2d;
568  double gapMargin2 = sqrt(2 * gapMargin*gapMargin);
569  gapMargin2 = (gapMargin2 + r2d)*(gapMargin2 + r2d);
570 
571  bool collect = true;
572  while (collect)
573  {
574  collect = false;
575  for (size_t i = 0; i < hits_in.size(); i++)
576  if (!Has(used, i))
577  {
578  art::Ptr<recob::Hit> hi = hits_in[i];
579  TVector2 hi_cm = pma::WireDriftToCm(detProp, hi->WireID().Wire, hi->PeakTime(), hi->WireID().Plane, hi->WireID().TPC, hi->WireID().Cryostat);
580 
581  bool accept = false;
582  //for (auto const& ho : hits_out)
583  for (size_t idx_o = 0; idx_o < hits_out.size(); idx_o++)
584  {
585  art::Ptr<recob::Hit> ho = hits_out[idx_o];
586 
587  double d2 = pma::Dist2(
588  hi_cm, pma::WireDriftToCm(detProp, ho->WireID().Wire, ho->PeakTime(), ho->WireID().Plane, ho->WireID().TPC, ho->WireID().Cryostat));
589 
590  if (hi->WireID().TPC == ho->WireID().TPC)
591  {
592  if (d2 < r2d2) { accept = true; break; }
593  }
594  else
595  {
596  if (d2 < gapMargin2) { accept = true; break; }
597  }
598  }
599  if (accept)
600  {
601  collect = true;
602  hits_out.push_back(hi);
603  used.push_back(i);
604  }
605  }
606  }
607  return true;
608  }
609  else return false;
610 }
bool Has(const std::vector< size_t > &v, size_t idx)
double Dist2(const TVector2 &v1, const TVector2 &v2)
Definition: Utilities.cxx:37
geo::WireID WireID() const
Definition: Hit.h:233
CryostatID_t Cryostat
Index of cryostat.
Definition: geo_types.h:212
TVector2 WireDriftToCm(detinfo::DetectorPropertiesData const &detProp, unsigned int wire, float drift, unsigned int plane, unsigned int tpc, unsigned int cryo)
Definition: Utilities.cxx:294
WireID_t Wire
Index of the wire within its plane.
Definition: geo_types.h:580
PlaneID_t Plane
Index of the plane within its TPC.
Definition: geo_types.h:493
float PeakTime() const
Time of the signal peak, in tick units.
Definition: Hit.h:218
TPCID_t TPC
Index of the TPC within its cryostat.
Definition: geo_types.h:406
bool dunefd::ShSeg::Has ( const std::vector< size_t > &  v,
size_t  idx 
)
private

Definition at line 518 of file ShSeg_module.cc.

519 {
520  for (auto c : v) if (c == idx) return true;
521  return false;
522 }
bool dunefd::ShSeg::InsideFidVol ( TLorentzVector const &  pvtx) const
private

Definition at line 472 of file ShSeg_module.cc.

473 {
475  double vtx[3] = {pvtx.X(), pvtx.Y(), pvtx.Z()};
476  bool inside = false;
477 
478  geo::TPCID idtpc = geom->FindTPCAtPosition(vtx);
479 
480  if (geom->HasTPC(idtpc))
481  {
482  const geo::TPCGeo& tpcgeo = geom->GetElement(idtpc);
483  double minx = tpcgeo.MinX(); double maxx = tpcgeo.MaxX();
484  double miny = tpcgeo.MinY(); double maxy = tpcgeo.MaxY();
485  double minz = tpcgeo.MinZ(); double maxz = tpcgeo.MaxZ();
486 
487  //x
488  double dista = fabs(minx - pvtx.X());
489  double distb = fabs(pvtx.X() - maxx);
490 
491  if ((pvtx.X() > minx) && (pvtx.X() < maxx) &&
492  (dista > fFidVolCut) && (distb > fFidVolCut))
493  {
494  inside = true;
495  }
496  else { inside = false; }
497 
498  //y
499  dista = fabs(maxy - pvtx.Y());
500  distb = fabs(pvtx.Y() - miny);
501  if (inside && (pvtx.Y() > miny) && (pvtx.Y() < maxy) &&
502  (dista > fFidVolCut) && (distb > fFidVolCut)) inside = true;
503  else inside = false;
504 
505  //z
506  dista = fabs(maxz - pvtx.Z());
507  distb = fabs(pvtx.Z() - minz);
508  if (inside && (pvtx.Z() > minz) && (pvtx.Z() < maxz) &&
509  (dista > fFidVolCut) && (distb > fFidVolCut)) inside = true;
510  else inside = false;
511  }
512 
513  return inside;
514 }
CryostatGeo const & GetElement(geo::CryostatID const &cryoid) const
double MinX() const
Returns the world x coordinate of the start of the box.
Definition: BoxBoundedGeo.h:88
Geometry information for a single TPC.
Definition: TPCGeo.h:38
double MaxX() const
Returns the world x coordinate of the end of the box.
Definition: BoxBoundedGeo.h:91
geo::TPCID FindTPCAtPosition(double const worldLoc[3]) const
Returns the ID of the TPC at specified location.
double MinZ() const
Returns the world z coordinate of the start of the box.
The data type to uniquely identify a TPC.
Definition: geo_types.h:386
double MaxY() const
Returns the world y coordinate of the end of the box.
bool HasTPC(geo::TPCID const &tpcid) const
Returns whether we have the specified TPC.
double MaxZ() const
Returns the world z coordinate of the end of the box.
double MinY() const
Returns the world y coordinate of the start of the box.
ShSeg& dunefd::ShSeg::operator= ( ShSeg const &  )
delete
ShSeg& dunefd::ShSeg::operator= ( ShSeg &&  )
delete
void dunefd::ShSeg::produce ( art::Event e)
overridevirtual

Implements art::EDProducer.

Definition at line 171 of file ShSeg_module.cc.

172 {
173  isdata = e.isRealData();
174  fEvNumber = e.id().event();
175 
176  std::unique_ptr< std::vector< recob::Track > > tracks(new std::vector< recob::Track >);
177  std::unique_ptr< std::vector< recob::SpacePoint > > allsp(new std::vector< recob::SpacePoint >);
178 
179  std::unique_ptr< art::Assns< recob::Track, recob::Hit > > trk2hit(new art::Assns< recob::Track, recob::Hit >);
180  std::unique_ptr< art::Assns< recob::Track, recob::SpacePoint > > trk2sp(new art::Assns< recob::Track, recob::SpacePoint >);
181  std::unique_ptr< art::Assns< recob::SpacePoint, recob::Hit > > sp2hit(new art::Assns< recob::SpacePoint, recob::Hit >);
182 
183  if (!isdata) BuildSegMC(e);
184 
185  // if (!isdata) BuildSegReco(e);
186 
187  if (fPmatracks.size())
188  {
189  size_t spStart = 0, spEnd = 0;
190  double sp_pos[3], sp_err[6];
191  for (size_t i = 0; i < 6; i++) sp_err[i] = 1.0;
192 
193  fTrkindex = 0;
194  for (auto trk : fPmatracks)
195  {
196  tracks->push_back(ConvertFrom(*trk));
197  fTrkindex++;
198 
199  std::vector< art::Ptr< recob::Hit > > hits2d;
201  spStart = allsp->size();
202 
203  for (int h = trk->size() - 1; h >= 0; h--)
204  {
205  pma::Hit3D* h3d = (*trk)[h];
206  if (!h3d->IsEnabled()) continue;
207  hits2d.push_back(h3d->Hit2DPtr());
208 
209  if ((h == 0) ||
210  (sp_pos[0] != h3d->Point3D().X()) ||
211  (sp_pos[1] != h3d->Point3D().Y()) ||
212  (sp_pos[2] != h3d->Point3D().Z()))
213  {
214  if (sp_hits.size()) // hits assigned to the previous sp
215  {
216  util::CreateAssn(*this, e, *allsp, sp_hits, *sp2hit);
217  sp_hits.clear();
218  }
219  sp_pos[0] = h3d->Point3D().X();
220  sp_pos[1] = h3d->Point3D().Y();
221  sp_pos[2] = h3d->Point3D().Z();
222  allsp->push_back(recob::SpacePoint(sp_pos, sp_err, 1.0));
223  }
224  sp_hits.push_back(h3d->Hit2DPtr());
225  }
226  if (sp_hits.size()) // hits assigned to the last sp
227  {
228  util::CreateAssn(*this, e, *allsp, sp_hits, *sp2hit);
229  }
230  spEnd = allsp->size();
231 
232  if (hits2d.size())
233  {
234  util::CreateAssn(*this, e, *tracks, *allsp, *trk2sp, spStart, spEnd);
235  util::CreateAssn(*this, e, *tracks, hits2d, *trk2hit);
236  }
237  }
238 
239  // data prods done, delete all pma::Track3D's
240  for (size_t t = 0; t < fPmatracks.size(); t++) delete fPmatracks[t];
241  fPmatracks.clear();
242  }
243 
244  e.put(std::move(tracks));
245  e.put(std::move(allsp));
246  e.put(std::move(trk2hit));
247  e.put(std::move(trk2sp));
248  e.put(std::move(sp2hit));
249 }
bool BuildSegMC(art::Event &e)
bool IsEnabled() const noexcept
Definition: PmaHit3D.h:161
TVector3 const & Point3D() const
Definition: PmaHit3D.h:55
bool isRealData() const
void push_back(Ptr< U > const &p)
Definition: PtrVector.h:435
recob::Track ConvertFrom(pma::Track3D const &src)
std::vector< pma::Track3D * > fPmatracks
def move(depos, offset)
Definition: depos.py:107
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.
size_type size() const
Definition: PtrVector.h:302
Definition: tracks.py:1
art::Ptr< recob::Hit > const & Hit2DPtr() const
Definition: PmaHit3D.h:49
EventNumber_t event() const
Definition: EventID.h:116
void clear()
Definition: PtrVector.h:533
EventID id() const
Definition: Event.cc:34
void dunefd::ShSeg::reconfigure ( fhicl::ParameterSet const &  p)

Definition at line 162 of file ShSeg_module.cc.

163 {
164  fHitsModuleLabel = p.get< std::string >("HitsModuleLabel");
165  fFidVolCut = p.get< double >("FidVolCut");
166  return;
167 }
std::string string
Definition: nybbler.cc:12
p
Definition: test.py:223
std::string fHitsModuleLabel
void dunefd::ShSeg::Smearel ( )
private

Definition at line 448 of file ShSeg_module.cc.

449 {
450  std::random_device rd;
451  std::mt19937 gen(rd());
452 
453  double mean = 0.0; double sigma = 0.7;
454  std::normal_distribution<double> dist (mean, sigma);
455 
456  fR0 = fabs(dist(gen));
457 }
gen
Definition: demo.py:24
constexpr double dist(const TReal *x, const TReal *y, const unsigned int dimension)
double mean(sqlite3 *db, std::string const &table_name, std::string const &column_name)
Definition: statistics.cc:16
void dunefd::ShSeg::Smearph ( )
private

Definition at line 461 of file ShSeg_module.cc.

462 {
463  std::random_device rd;
464  std::mt19937 gen(rd());
465 
466  std::uniform_real_distribution<double> distribution(0.0, fR0);
467  fR1 = distribution(gen);
468 }
gen
Definition: demo.py:24

Member Data Documentation

calo::CalorimetryAlg dunefd::ShSeg::fCalorimetryAlg
private

Definition at line 113 of file ShSeg_module.cc.

double dunefd::ShSeg::fDedx
private

Definition at line 119 of file ShSeg_module.cc.

double dunefd::ShSeg::fDedxavg
private

Definition at line 118 of file ShSeg_module.cc.

int dunefd::ShSeg::fEvNumber
private

Definition at line 117 of file ShSeg_module.cc.

TTree* dunefd::ShSeg::fEvTree
private

Definition at line 121 of file ShSeg_module.cc.

double dunefd::ShSeg::fFidVolCut
private

Definition at line 108 of file ShSeg_module.cc.

std::string dunefd::ShSeg::fHitsModuleLabel
private

Definition at line 112 of file ShSeg_module.cc.

std::vector< pma::Track3D* > dunefd::ShSeg::fPmatracks
private

Definition at line 106 of file ShSeg_module.cc.

pma::ProjectionMatchingAlg dunefd::ShSeg::fProjectionMatchingAlg
private

Definition at line 114 of file ShSeg_module.cc.

double dunefd::ShSeg::fR0
private

Definition at line 109 of file ShSeg_module.cc.

double dunefd::ShSeg::fR1
private

Definition at line 109 of file ShSeg_module.cc.

TTree* dunefd::ShSeg::fShTree
private

Definition at line 122 of file ShSeg_module.cc.

int dunefd::ShSeg::fTrkindex
private

Definition at line 116 of file ShSeg_module.cc.

short dunefd::ShSeg::isdata
private

Definition at line 110 of file ShSeg_module.cc.


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