LArSoft  v10_04_05
Liquid Argon Software toolkit - https://larsoft.org/
hit::FFTHitFinder Class Reference
Inheritance diagram for hit::FFTHitFinder:
art::EDProducer art::detail::Producer art::detail::LegacyModule art::Modifier art::ModuleBase art::ProductRegistryHelper

Public Types

using ModuleType = EDProducer
 
template<typename UserConfig , typename KeysToIgnore = void>
using Table = Modifier::Table< UserConfig, KeysToIgnore >
 

Public Member Functions

 FFTHitFinder (fhicl::ParameterSet const &pset)
 
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)
 
void fillProductDescriptions ()
 
void registerProducts (ProductDescriptions &productsToRegister)
 
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)
 
std::unique_ptr< Worker > makeWorker (WorkerParams const &wp)
 
template<typename T , BranchType BT>
ViewToken< T > consumesView (InputTag const &tag)
 
template<typename T , BranchType BT>
ViewToken< T > mayConsumeView (InputTag const &tag)
 

Protected Member Functions

ConsumesCollector & consumesCollector ()
 
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 ()
 

Private Member Functions

void produce (art::Event &evt) override
 

Private Attributes

std::string fCalDataModuleLabel
 
double fMinSigInd
 Induction signal height threshold. More...
 
double fMinSigCol
 Collection signal height threshold. More...
 
double fIndWidth
 Initial width for induction fit. More...
 
double fColWidth
 Initial width for collection fit. More...
 
double fIndMinWidth
 Minimum induction hit width. More...
 
double fColMinWidth
 Minimum collection hit width. More...
 
int fMaxMultiHit
 maximum hits for multi fit More...
 
int fAreaMethod
 Type of area calculation. More...
 
std::vector< double > fAreaNorms
 factors for converting area to same units as peak height More...
 

Detailed Description

Definition at line 38 of file FFTHitFinder_module.cc.

Member Typedef Documentation

Definition at line 17 of file EDProducer.h.

template<typename UserConfig , typename KeysToIgnore = void>
using art::detail::Producer::Table = Modifier::Table<UserConfig, KeysToIgnore>
inherited

Definition at line 26 of file Producer.h.

Constructor & Destructor Documentation

hit::FFTHitFinder::FFTHitFinder ( fhicl::ParameterSet const &  pset)
explicit

Definition at line 60 of file FFTHitFinder_module.cc.

References recob::HitAndAssociationsWriterBase::declare_products(), fAreaMethod, fAreaNorms, fCalDataModuleLabel, fColMinWidth, fColWidth, fIndMinWidth, fIndWidth, fMaxMultiHit, fMinSigCol, fMinSigInd, and art::ProductRegistryHelper::producesCollector().

60  : EDProducer{pset}
61  {
62  fCalDataModuleLabel = pset.get<std::string>("CalDataModuleLabel");
63  fMinSigInd = pset.get<double>("MinSigInd");
64  fMinSigCol = pset.get<double>("MinSigCol");
65  fIndWidth = pset.get<double>("IndWidth");
66  fColWidth = pset.get<double>("ColWidth");
67  fIndMinWidth = pset.get<double>("IndMinWidth");
68  fColMinWidth = pset.get<double>("ColMinWidth");
69  fMaxMultiHit = pset.get<int>("MaxMultiHit");
70  fAreaMethod = pset.get<int>("AreaMethod");
71  fAreaNorms = pset.get<std::vector<double>>("AreaNorms");
72 
73  // let HitCollectionCreator declare that we are going to produce
74  // hits and associations with wires and raw digits
75  // (with no particular product label)
77  }
double fIndMinWidth
Minimum induction hit width.
double fMinSigInd
Induction signal height threshold.
EDProducer(fhicl::ParameterSet const &pset)
Definition: EDProducer.cc:6
double fIndWidth
Initial width for induction fit.
int fMaxMultiHit
maximum hits for multi fit
double fColWidth
Initial width for collection fit.
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
std::vector< double > fAreaNorms
factors for converting area to same units as peak height
std::string fCalDataModuleLabel
double fColMinWidth
Minimum collection hit width.
ProducesCollector & producesCollector() noexcept
double fMinSigCol
Collection signal height threshold.
int fAreaMethod
Type of area calculation.

Member Function Documentation

template<typename T , BranchType BT>
ProductToken< T > art::ModuleBase::consumes ( InputTag const &  tag)
protectedinherited

Definition at line 61 of file ModuleBase.h.

References art::ModuleBase::collector_, and art::ConsumesCollector::consumes().

62  {
63  return collector_.consumes<T, BT>(tag);
64  }
ConsumesCollector collector_
Definition: ModuleBase.h:56
ProductToken< T > consumes(InputTag const &)
ConsumesCollector & art::ModuleBase::consumesCollector ( )
protectedinherited

Definition at line 57 of file ModuleBase.cc.

References art::ModuleBase::collector_.

58  {
59  return collector_;
60  }
ConsumesCollector collector_
Definition: ModuleBase.h:56
template<typename T , BranchType BT>
void art::ModuleBase::consumesMany ( )
protectedinherited

Definition at line 75 of file ModuleBase.h.

References art::ModuleBase::collector_, and art::ConsumesCollector::consumesMany().

76  {
77  collector_.consumesMany<T, BT>();
78  }
ConsumesCollector collector_
Definition: ModuleBase.h:56
template<typename Element , BranchType = InEvent>
ViewToken<Element> art::ModuleBase::consumesView ( InputTag const &  )
protectedinherited
template<typename T , BranchType BT>
ViewToken<T> art::ModuleBase::consumesView ( InputTag const &  tag)
inherited

Definition at line 68 of file ModuleBase.h.

References art::ModuleBase::collector_, and art::ConsumesCollector::consumesView().

69  {
70  return collector_.consumesView<T, BT>(tag);
71  }
ConsumesCollector collector_
Definition: ModuleBase.h:56
ViewToken< Element > consumesView(InputTag const &)
void art::detail::Producer::doBeginJob ( SharedResources const &  resources)
inherited

Definition at line 22 of file Producer.cc.

References art::detail::Producer::beginJobWithFrame(), and art::detail::Producer::setupQueues().

23  {
24  setupQueues(resources);
25  ProcessingFrame const frame{ScheduleID{}};
26  beginJobWithFrame(frame);
27  }
virtual void setupQueues(SharedResources const &)=0
virtual void beginJobWithFrame(ProcessingFrame const &)=0
bool art::detail::Producer::doBeginRun ( RunPrincipal rp,
ModuleContext const &  mc 
)
inherited

Definition at line 65 of file Producer.cc.

References art::detail::Producer::beginRunWithFrame(), art::RangeSet::forRun(), art::RunPrincipal::makeRun(), r, art::RunPrincipal::runID(), and art::ModuleContext::scheduleID().

66  {
67  auto r = rp.makeRun(mc, RangeSet::forRun(rp.runID()));
68  ProcessingFrame const frame{mc.scheduleID()};
69  beginRunWithFrame(r, frame);
70  r.commitProducts();
71  return true;
72  }
TRandom r
Definition: spectrum.C:23
virtual void beginRunWithFrame(Run &, ProcessingFrame const &)=0
static RangeSet forRun(RunID)
Definition: RangeSet.cc:51
bool art::detail::Producer::doBeginSubRun ( SubRunPrincipal srp,
ModuleContext const &  mc 
)
inherited

Definition at line 85 of file Producer.cc.

References art::detail::Producer::beginSubRunWithFrame(), art::RangeSet::forSubRun(), art::SubRunPrincipal::makeSubRun(), art::ModuleContext::scheduleID(), and art::SubRunPrincipal::subRunID().

86  {
87  auto sr = srp.makeSubRun(mc, RangeSet::forSubRun(srp.subRunID()));
88  ProcessingFrame const frame{mc.scheduleID()};
89  beginSubRunWithFrame(sr, frame);
90  sr.commitProducts();
91  return true;
92  }
virtual void beginSubRunWithFrame(SubRun &, ProcessingFrame const &)=0
static RangeSet forSubRun(SubRunID)
Definition: RangeSet.cc:57
void art::detail::Producer::doEndJob ( )
inherited

Definition at line 30 of file Producer.cc.

References art::detail::Producer::endJobWithFrame().

31  {
32  ProcessingFrame const frame{ScheduleID{}};
33  endJobWithFrame(frame);
34  }
virtual void endJobWithFrame(ProcessingFrame const &)=0
bool art::detail::Producer::doEndRun ( RunPrincipal rp,
ModuleContext const &  mc 
)
inherited

Definition at line 75 of file Producer.cc.

References art::detail::Producer::endRunWithFrame(), art::RunPrincipal::makeRun(), r, art::ModuleContext::scheduleID(), and art::Principal::seenRanges().

76  {
77  auto r = rp.makeRun(mc, rp.seenRanges());
78  ProcessingFrame const frame{mc.scheduleID()};
79  endRunWithFrame(r, frame);
80  r.commitProducts();
81  return true;
82  }
TRandom r
Definition: spectrum.C:23
virtual void endRunWithFrame(Run &, ProcessingFrame const &)=0
bool art::detail::Producer::doEndSubRun ( SubRunPrincipal srp,
ModuleContext const &  mc 
)
inherited

Definition at line 95 of file Producer.cc.

References art::detail::Producer::endSubRunWithFrame(), art::SubRunPrincipal::makeSubRun(), art::ModuleContext::scheduleID(), and art::Principal::seenRanges().

96  {
97  auto sr = srp.makeSubRun(mc, srp.seenRanges());
98  ProcessingFrame const frame{mc.scheduleID()};
99  endSubRunWithFrame(sr, frame);
100  sr.commitProducts();
101  return true;
102  }
virtual void endSubRunWithFrame(SubRun &, ProcessingFrame const &)=0
bool art::detail::Producer::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 
)
inherited

Definition at line 105 of file Producer.cc.

References art::detail::Producer::checkPutProducts_, e, art::EventPrincipal::makeEvent(), art::detail::Producer::produceWithFrame(), and art::ModuleContext::scheduleID().

110  {
111  auto e = ep.makeEvent(mc);
112  ++counts_run;
113  ProcessingFrame const frame{mc.scheduleID()};
114  produceWithFrame(e, frame);
115  e.commitProducts(checkPutProducts_, &expectedProducts<InEvent>());
116  ++counts_passed;
117  return true;
118  }
bool const checkPutProducts_
Definition: Producer.h:70
Float_t e
Definition: plot.C:35
virtual void produceWithFrame(Event &, ProcessingFrame const &)=0
void art::detail::Producer::doRespondToCloseInputFile ( FileBlock const &  fb)
inherited

Definition at line 44 of file Producer.cc.

References art::detail::Producer::respondToCloseInputFileWithFrame().

45  {
46  ProcessingFrame const frame{ScheduleID{}};
48  }
virtual void respondToCloseInputFileWithFrame(FileBlock const &, ProcessingFrame const &)=0
TFile fb("Li6.root")
void art::detail::Producer::doRespondToCloseOutputFiles ( FileBlock const &  fb)
inherited

Definition at line 58 of file Producer.cc.

References art::detail::Producer::respondToCloseOutputFilesWithFrame().

59  {
60  ProcessingFrame const frame{ScheduleID{}};
62  }
virtual void respondToCloseOutputFilesWithFrame(FileBlock const &, ProcessingFrame const &)=0
TFile fb("Li6.root")
void art::detail::Producer::doRespondToOpenInputFile ( FileBlock const &  fb)
inherited

Definition at line 37 of file Producer.cc.

References art::detail::Producer::respondToOpenInputFileWithFrame().

38  {
39  ProcessingFrame const frame{ScheduleID{}};
41  }
virtual void respondToOpenInputFileWithFrame(FileBlock const &, ProcessingFrame const &)=0
TFile fb("Li6.root")
void art::detail::Producer::doRespondToOpenOutputFiles ( FileBlock const &  fb)
inherited

Definition at line 51 of file Producer.cc.

References art::detail::Producer::respondToOpenOutputFilesWithFrame().

52  {
53  ProcessingFrame const frame{ScheduleID{}};
55  }
virtual void respondToOpenOutputFilesWithFrame(FileBlock const &, ProcessingFrame const &)=0
TFile fb("Li6.root")
void art::Modifier::fillProductDescriptions ( )
inherited

Definition at line 10 of file Modifier.cc.

References art::ProductRegistryHelper::fillDescriptions(), and art::ModuleBase::moduleDescription().

11  {
13  }
void fillDescriptions(ModuleDescription const &md)
ModuleDescription const & moduleDescription() const
Definition: ModuleBase.cc:13
std::array< std::vector< ProductInfo >, NumBranchTypes > const & art::ModuleBase::getConsumables ( ) const
inherited

Definition at line 43 of file ModuleBase.cc.

References art::ModuleBase::collector_, and art::ConsumesCollector::getConsumables().

44  {
45  return collector_.getConsumables();
46  }
ConsumesCollector collector_
Definition: ModuleBase.h:56
std::array< std::vector< ProductInfo >, NumBranchTypes > const & getConsumables() const
std::unique_ptr< Worker > art::ModuleBase::makeWorker ( WorkerParams const &  wp)
inherited

Definition at line 37 of file ModuleBase.cc.

References art::ModuleBase::doMakeWorker(), and art::NumBranchTypes.

38  {
39  return doMakeWorker(wp);
40  }
virtual std::unique_ptr< Worker > doMakeWorker(WorkerParams const &wp)=0
template<typename T , BranchType BT>
ProductToken< T > art::ModuleBase::mayConsume ( InputTag const &  tag)
protectedinherited

Definition at line 82 of file ModuleBase.h.

References art::ModuleBase::collector_, and art::ConsumesCollector::mayConsume().

83  {
84  return collector_.mayConsume<T, BT>(tag);
85  }
ProductToken< T > mayConsume(InputTag const &)
ConsumesCollector collector_
Definition: ModuleBase.h:56
template<typename T , BranchType BT>
void art::ModuleBase::mayConsumeMany ( )
protectedinherited

Definition at line 96 of file ModuleBase.h.

References art::ModuleBase::collector_, and art::ConsumesCollector::mayConsumeMany().

97  {
98  collector_.mayConsumeMany<T, BT>();
99  }
ConsumesCollector collector_
Definition: ModuleBase.h:56
template<typename Element , BranchType = InEvent>
ViewToken<Element> art::ModuleBase::mayConsumeView ( InputTag const &  )
protectedinherited
template<typename T , BranchType BT>
ViewToken<T> art::ModuleBase::mayConsumeView ( InputTag const &  tag)
inherited

Definition at line 89 of file ModuleBase.h.

References art::ModuleBase::collector_, and art::ConsumesCollector::mayConsumeView().

90  {
91  return collector_.mayConsumeView<T, BT>(tag);
92  }
ConsumesCollector collector_
Definition: ModuleBase.h:56
ViewToken< Element > mayConsumeView(InputTag const &)
ModuleDescription const & art::ModuleBase::moduleDescription ( ) const
inherited

Definition at line 13 of file ModuleBase.cc.

References art::errors::LogicError.

Referenced by art::OutputModule::doRespondToOpenInputFile(), art::OutputModule::doWriteEvent(), art::Modifier::fillProductDescriptions(), art::OutputModule::makePlugins_(), art::OutputWorker::OutputWorker(), reco::shower::LArPandoraModularShowerCreation::produce(), art::Modifier::registerProducts(), and art::OutputModule::registerProducts().

14  {
15  if (md_.has_value()) {
16  return *md_;
17  }
18 
20  "There was an error while calling moduleDescription().\n"}
21  << "The moduleDescription() base-class member function cannot be called\n"
22  "during module construction. To determine which module is "
23  "responsible\n"
24  "for calling it, find the '<module type>:<module "
25  "label>@Construction'\n"
26  "tag in the message prefix above. Please contact artists@fnal.gov\n"
27  "for guidance.\n";
28  }
cet::coded_exception< errors::ErrorCodes, ExceptionDetail::translate > Exception
Definition: Exception.h:66
std::optional< ModuleDescription > md_
Definition: ModuleBase.h:55
void hit::FFTHitFinder::produce ( art::Event evt)
overrideprivatevirtual
Todo:
  • just get the integral from the fit for totSig
Todo:
need to have a disambiguation algorithm somewhere in here
Todo:
  • multiplicity and local_index have to be determined

Implements art::EDProducer.

Definition at line 83 of file FFTHitFinder_module.cc.

References recob::Wire::Channel(), DEFINE_ART_MODULE, recob::HitCollectionCreator::emplace_back(), fAreaMethod, fAreaNorms, fCalDataModuleLabel, fColMinWidth, fColWidth, fIndMinWidth, fIndWidth, fMaxMultiHit, fMinSigCol, fMinSigInd, Get, art::ProductRetriever::getByLabel(), raw::InvalidChannelID, geo::kCollection, geo::kInduction, recob::HitCollectionCreator::put_into(), recob::Wire::Signal(), and util::size().

84  {
85  // this object contains the hit collection
86  // and its associations to wires and raw digits:
88 
89  // Read in the wire List object(s).
91  evt.getByLabel(fCalDataModuleLabel, wireVecHandle);
92  auto const& wireReadoutGeom = art::ServiceHandle<geo::WireReadout const>()->Get();
93 
94  // also get the raw digits associated with wires
95  art::FindOneP<raw::RawDigit> WireToRawDigits(wireVecHandle, evt, fCalDataModuleLabel);
96 
97  std::vector<int> startTimes; // stores time of 1st local minimum
98  std::vector<int> maxTimes; // stores time of local maximum
99  std::vector<int> endTimes; // stores time of 2nd local minimum
100  int time = 0; // current time bin
101  int minTimeHolder = 0; // current start time
102  raw::ChannelID_t channel = raw::InvalidChannelID; // channel number
103  bool maxFound = false; // Flag for whether a peak > threshold has been found
104  double threshold = 0.; // minimum signal size for id'ing a hit
105  double fitWidth = 0.; // hit fit width initial value
106  double minWidth = 0.; // minimum hit width
107  std::string eqn = "gaus(0)"; // string for equation to fit
108  geo::SigType_t sigType = geo::kInduction; // type of plane we are looking at
109  std::stringstream numConv;
110 
111  //loop over wires
112  for (unsigned int wireIter = 0; wireIter < wireVecHandle->size(); wireIter++) {
113 
114  art::Ptr<recob::Wire> wire(wireVecHandle, wireIter);
115  startTimes.clear();
116  maxTimes.clear();
117  endTimes.clear();
118  std::vector<float> signal(wire->Signal());
119  std::vector<float>::iterator timeIter; // iterator for time bins
120  time = 0;
121  minTimeHolder = 0;
122  maxFound = false;
123  channel = wire->Channel();
124  sigType = wireReadoutGeom.SignalType(channel);
125 
126  //Set the appropriate signal widths and thresholds
127  if (sigType == geo::kInduction) {
128  threshold = fMinSigInd;
129  fitWidth = fIndWidth;
130  minWidth = fIndMinWidth;
131  }
132  else if (sigType == geo::kCollection) {
133  threshold = fMinSigCol;
134  fitWidth = fColWidth;
135  minWidth = fColMinWidth;
136  }
137  // loop over signal
138  for (timeIter = signal.begin(); timeIter + 2 < signal.end(); timeIter++) {
139  //test if timeIter+1 is a local minimum
140  if (*timeIter > *(timeIter + 1) && *(timeIter + 1) < *(timeIter + 2)) {
141  //only add points if already found a local max above threshold.
142  if (maxFound) {
143  endTimes.push_back(time + 1);
144  maxFound = false;
145  //keep these in case new hit starts right away
146  minTimeHolder = time + 2;
147  }
148  else
149  minTimeHolder = time + 1;
150  }
151  //if not a minimum, test if we are at a local maximum
152  //if so, and the max value is above threshold, add it and proceed.
153  else if (*timeIter < *(timeIter + 1) && *(timeIter + 1) > *(timeIter + 2) &&
154  *(timeIter + 1) > threshold) {
155  maxFound = true;
156  maxTimes.push_back(time + 1);
157  startTimes.push_back(minTimeHolder);
158  }
159  time++;
160  } //end loop over signal vec
161 
162  //if no inflection found before end, but peak found add end point
163  while (maxTimes.size() > endTimes.size())
164  endTimes.push_back(signal.size() - 1);
165  if (startTimes.size() == 0) continue;
166 
167  //All code below does the fitting, adding of hits
168  //to the hit vector and when all wires are complete
169  //saving them
170  double totSig(0); // stores the total hit signal
171  double startT(0); // stores the start time
172  double endT(0); // stores the end time
173  int numHits(0); // number of consecutive hits being fitted
174  int size(0); // size of data vector for fit
175  int hitIndex(0); // index of current hit in sequence
176  double amplitude(0), position(0), width(0); //fit parameters
177  double amplitudeErr(0), positionErr(0), widthErr(0); //fit errors
178  double goodnessOfFit(0), chargeErr(0); //Chi2/NDF and error on charge
179  double minPeakHeight(0); //lowest peak height in multi-hit fit
180 
181  //stores gaussian paramters first index is the hit number
182  //the second refers to height, position, and width respectively
183  std::vector<double> hitSig;
184 
185  //add found hits to hit vector
186  while (hitIndex < (signed)startTimes.size()) {
187 
188  startT = endT = 0;
189  numHits = 1;
190  minPeakHeight = signal[maxTimes[hitIndex]];
191 
192  //consider adding pulse to group of consecutive hits if:
193  //1 less than max consecutive hits
194  //2 we are not at the last point in the signal vector
195  //3 the height of the dip between the two is greater than threshold/2
196  //4 and there is no gap between them
197  while (numHits < fMaxMultiHit && numHits + hitIndex < (signed)endTimes.size() &&
198  signal[endTimes[hitIndex + numHits - 1]] > threshold / 2.0 &&
199  startTimes[hitIndex + numHits] - endTimes[hitIndex + numHits - 1] < 2) {
200 
201  if (signal[maxTimes[hitIndex + numHits]] < minPeakHeight)
202  minPeakHeight = signal[maxTimes[hitIndex + numHits]];
203 
204  ++numHits;
205  }
206 
207  //finds the first point > 1/2 the smallest peak
208  startT = startTimes[hitIndex];
209 
210  while (signal[(int)startT] < minPeakHeight / 2.0)
211  ++startT;
212 
213  //finds the first point from the end > 1/2 the smallest peak
214  endT = endTimes[hitIndex + numHits - 1];
215 
216  while (signal[(int)endT] < minPeakHeight / 2.0)
217  --endT;
218  size = (int)(endT - startT);
219  TH1D hitSignal("hitSignal", "", size, startT, endT);
220  for (int i = (int)startT; i < (int)endT; ++i)
221  hitSignal.Fill(i, signal[i]);
222 
223  //build the TFormula
224  eqn = "gaus(0)";
225 
226  for (int i = 3; i < numHits * 3; i += 3) {
227  eqn.append("+gaus(");
228  numConv.str("");
229  numConv << i;
230  eqn.append(numConv.str());
231  eqn.append(")");
232  }
233 
234  TF1 gSum("gSum", eqn.c_str(), 0, size);
235 
236  if (numHits > 1) {
237  TArrayD data(numHits * numHits);
238  TVectorD amps(numHits);
239  for (int i = 0; i < numHits; ++i) {
240  amps[i] = signal[maxTimes[hitIndex + i]];
241  for (int j = 0; j < numHits; j++)
242  data[i + numHits * j] =
243  TMath::Gaus(maxTimes[hitIndex + j], maxTimes[hitIndex + i], fitWidth);
244  } //end loop over hits
245 
246  //This section uses a linear approximation in order to get an
247  //initial value of the individual hit amplitudes
248  try {
249  TMatrixD h(numHits, numHits);
250  h.Use(numHits, numHits, data.GetArray());
251  TDecompSVD a(h);
252  a.Solve(amps);
253  }
254  catch (...) {
255  mf::LogInfo("FFTHitFinder") << "TDcompSVD failed";
256  hitIndex += numHits;
257  continue;
258  }
259 
260  for (int i = 0; i < numHits; ++i) {
261  //if the approximation makes a peak vanish
262  //set initial height as average of threshold and
263  //raw peak height
264  if (amps[i] > 0)
265  amplitude = amps[i];
266  else
267  amplitude = 0.5 * (threshold + signal[maxTimes[hitIndex + i]]);
268  gSum.SetParameter(3 * i, amplitude);
269  gSum.SetParameter(1 + 3 * i, maxTimes[hitIndex + i]);
270  gSum.SetParameter(2 + 3 * i, fitWidth);
271  gSum.SetParLimits(3 * i, 0.0, 3.0 * amplitude);
272  gSum.SetParLimits(1 + 3 * i, startT, endT);
273  gSum.SetParLimits(2 + 3 * i, 0.0, 10.0 * fitWidth);
274  } //end loop over hits
275  } //end if numHits > 1
276  else {
277  gSum.SetParameters(signal[maxTimes[hitIndex]], maxTimes[hitIndex], fitWidth);
278  gSum.SetParLimits(0, 0.0, 1.5 * signal[maxTimes[hitIndex]]);
279  gSum.SetParLimits(1, startT, endT);
280  gSum.SetParLimits(2, 0.0, 10.0 * fitWidth);
281  }
282 
284  hitSignal.Fit(&gSum, "QNRW", "", startT, endT);
285  for (int hitNumber = 0; hitNumber < numHits; ++hitNumber) {
286  totSig = 0;
287  if (gSum.GetParameter(3 * hitNumber) > threshold / 2.0 &&
288  gSum.GetParameter(3 * hitNumber + 2) > minWidth) {
289  amplitude = gSum.GetParameter(3 * hitNumber);
290  position = gSum.GetParameter(3 * hitNumber + 1);
291  width = gSum.GetParameter(3 * hitNumber + 2);
292  amplitudeErr = gSum.GetParError(3 * hitNumber);
293  positionErr = gSum.GetParError(3 * hitNumber + 1);
294  widthErr = gSum.GetParError(3 * hitNumber + 2);
295  goodnessOfFit = gSum.GetChisquare() / (double)gSum.GetNDF();
296  int DoF = gSum.GetNDF();
297 
298  //estimate error from area of Gaussian
299  chargeErr = std::sqrt(TMath::Pi()) * (amplitudeErr * width + widthErr * amplitude);
300 
301  hitSig.resize(size);
302 
303  for (int sigPos = 0; sigPos < size; ++sigPos) {
304  hitSig[sigPos] = amplitude * TMath::Gaus(sigPos + startT, position, width);
305  totSig += hitSig[(int)sigPos];
306  }
307 
308  if (fAreaMethod)
309  totSig = std::sqrt(2 * TMath::Pi()) * amplitude * width / fAreaNorms[(size_t)sigType];
310 
311  // get the WireID for this hit
312  std::vector<geo::WireID> wids = wireReadoutGeom.ChannelToWire(channel);
314  // for now, just take the first option returned from ChannelToWire
315  geo::WireID wid = wids[0];
316 
317  // make the hit
318  recob::HitCreator hit(*wire, // wire
319  wid, // wireID
320  (int)startT, // start_tick
321  (int)endT, // end_tick
322  width, // rms
323  position, // peak_time
324  positionErr, // sigma_peak_time
325  amplitude, // peak_amplitude
326  amplitudeErr, // sigma_peak_amplitude
327  totSig, // hit_integral
328  chargeErr, // hit_sigma_integral
329  std::accumulate // summedADC
330  (signal.begin() + (int)startT, signal.begin() + (int)endT, 0.),
331  1, // multiplicity
332  -1, // local_index
334  goodnessOfFit, // goodness_of_fit
335  DoF // dof
336  );
337 
338  // get the object associated with the original hit
339  art::Ptr<raw::RawDigit> rawdigits = WireToRawDigits.at(wireIter);
340 
341  hcol.emplace_back(hit.move(), wire, rawdigits);
342  } //end if over threshold
343  } //end loop over hits
344  hitIndex += numHits;
345  } // end while on hitIndex<(signed)startTimes.size()
346 
347  } // while on Wires
348 
349  // put the hit collection and associations into the event
350  hcol.put_into(evt);
351 
352  } // End of produce()
double fIndMinWidth
Minimum induction hit width.
intermediate_table::iterator iterator
MaybeLogger_< ELseverityLevel::ELsev_info, false > LogInfo
double fMinSigInd
Induction signal height threshold.
double fIndWidth
Initial width for induction fit.
cout<< "Opened file "<< fin<< " ixs= "<< ixs<< endl;if(ixs==0) hhh=(TH1F *) fff-> Get("h1")
Definition: AddMC.C:8
int fMaxMultiHit
maximum hits for multi fit
double fColWidth
Initial width for collection fit.
std::vector< double > fAreaNorms
factors for converting area to same units as peak height
decltype(auto) constexpr size(T &&obj)
ADL-aware version of std::size.
Definition: StdUtils.h:101
Class managing the creation of a new recob::Hit object.
Definition: HitCreator.h:87
std::string fCalDataModuleLabel
constexpr ChannelID_t InvalidChannelID
ID of an invalid channel.
Definition: RawTypes.h:31
A class handling a collection of hits and its associations.
Definition: HitCreator.h:489
Signal from induction planes.
Definition: geo_types.h:147
enum geo::_plane_sigtype SigType_t
Enumerate the possible plane projections.
double fColMinWidth
Minimum collection hit width.
Detector simulation of raw signals on wires.
bool getByLabel(std::string const &label, std::string const &instance, Handle< PROD > &result) const
double fMinSigCol
Collection signal height threshold.
int fAreaMethod
Type of area calculation.
unsigned int ChannelID_t
Type representing the ID of a readout channel.
Definition: RawTypes.h:28
Definition: fwd.h:26
Signal from collection planes.
Definition: geo_types.h:148
void art::Modifier::registerProducts ( ProductDescriptions productsToRegister)
inherited

Definition at line 16 of file Modifier.cc.

References art::ModuleBase::moduleDescription(), and art::ProductRegistryHelper::registerProducts().

17  {
18  ProductRegistryHelper::registerProducts(productsToRegister,
20  }
void registerProducts(ProductDescriptions &productsToRegister, ModuleDescription const &md)
ModuleDescription const & moduleDescription() const
Definition: ModuleBase.cc:13
void art::ModuleBase::setModuleDescription ( ModuleDescription const &  md)
inherited

Definition at line 31 of file ModuleBase.cc.

References art::ModuleBase::md_.

32  {
33  md_ = md;
34  }
std::optional< ModuleDescription > md_
Definition: ModuleBase.h:55
void art::ModuleBase::sortConsumables ( std::string const &  current_process_name)
inherited

Definition at line 49 of file ModuleBase.cc.

References art::ModuleBase::collector_, and art::ConsumesCollector::sortConsumables().

50  {
51  // Now that we know we have seen all the consumes declarations,
52  // sort the results for fast lookup later.
53  collector_.sortConsumables(current_process_name);
54  }
ConsumesCollector collector_
Definition: ModuleBase.h:56
void sortConsumables(std::string const &current_process_name)

Member Data Documentation

int hit::FFTHitFinder::fAreaMethod
private

Type of area calculation.

Definition at line 54 of file FFTHitFinder_module.cc.

Referenced by FFTHitFinder(), and produce().

std::vector<double> hit::FFTHitFinder::fAreaNorms
private

factors for converting area to same units as peak height

Definition at line 55 of file FFTHitFinder_module.cc.

Referenced by FFTHitFinder(), and produce().

std::string hit::FFTHitFinder::fCalDataModuleLabel
private

Definition at line 46 of file FFTHitFinder_module.cc.

Referenced by FFTHitFinder(), and produce().

double hit::FFTHitFinder::fColMinWidth
private

Minimum collection hit width.

Definition at line 52 of file FFTHitFinder_module.cc.

Referenced by FFTHitFinder(), and produce().

double hit::FFTHitFinder::fColWidth
private

Initial width for collection fit.

Definition at line 50 of file FFTHitFinder_module.cc.

Referenced by FFTHitFinder(), and produce().

double hit::FFTHitFinder::fIndMinWidth
private

Minimum induction hit width.

Definition at line 51 of file FFTHitFinder_module.cc.

Referenced by FFTHitFinder(), and produce().

double hit::FFTHitFinder::fIndWidth
private

Initial width for induction fit.

Definition at line 49 of file FFTHitFinder_module.cc.

Referenced by FFTHitFinder(), and produce().

int hit::FFTHitFinder::fMaxMultiHit
private

maximum hits for multi fit

Definition at line 53 of file FFTHitFinder_module.cc.

Referenced by FFTHitFinder(), and produce().

double hit::FFTHitFinder::fMinSigCol
private

Collection signal height threshold.

Definition at line 48 of file FFTHitFinder_module.cc.

Referenced by FFTHitFinder(), and produce().

double hit::FFTHitFinder::fMinSigInd
private

Induction signal height threshold.

Definition at line 47 of file FFTHitFinder_module.cc.

Referenced by FFTHitFinder(), and produce().


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