LArSoft  v09_90_00
Liquid Argon Software toolkit - https://larsoft.org/
reco_tool::PeakFitterMrqdt Class Reference
Inheritance diagram for reco_tool::PeakFitterMrqdt:
reco_tool::IPeakFitter

Public Member Functions

 PeakFitterMrqdt (const fhicl::ParameterSet &pset)
 
void findPeakParameters (const std::vector< float > &, const ICandidateHitFinder::HitCandidateVec &, PeakParamsVec &, double &, int &) const override
 

Private Types

using PeakParamsVec = std::vector< PeakFitParams_t >
 

Private Attributes

const double fMinWidth
 
const double fMaxWidthMult
 
const double fPeakRange
 
const double fAmpRange
 
std::unique_ptr< gshf::MarqFitAlgfMarqFitAlg
 
const geo::GeometryCorefGeometry = lar::providerFrom<geo::Geometry>()
 

Detailed Description

Definition at line 19 of file PeakFitterMrqdt_tool.cc.

Member Typedef Documentation

Definition at line 34 of file IPeakFitter.h.

Constructor & Destructor Documentation

reco_tool::PeakFitterMrqdt::PeakFitterMrqdt ( const fhicl::ParameterSet pset)
explicit

Definition at line 44 of file PeakFitterMrqdt_tool.cc.

45  : fMinWidth(pset.get<double>("MinWidth", 0.5))
46  , fMaxWidthMult(pset.get<double>("MaxWidthMult", 3.))
47  , fPeakRange(pset.get<double>("PeakRangeFact", 2.))
48  , fAmpRange(pset.get<double>("PeakAmpRange", 2.))
49  {}
T get(std::string const &key) const
Definition: ParameterSet.h:314

Member Function Documentation

void reco_tool::PeakFitterMrqdt::findPeakParameters ( const std::vector< float > &  signal,
const ICandidateHitFinder::HitCandidateVec fhc_vec,
PeakParamsVec mhpp_vec,
double &  chi2PerNDF,
int &  NDF 
) const
overridevirtual

Implements reco_tool::IPeakFitter.

Definition at line 53 of file PeakFitterMrqdt_tool.cc.

References util::abs(), DEFINE_ART_CLASS_TOOL, e, fAmpRange, fMarqFitAlg, fMaxWidthMult, fMinWidth, fPeakRange, reco_tool::IPeakFitter::PeakFitParams_t::peakAmplitude, reco_tool::IPeakFitter::PeakFitParams_t::peakAmplitudeError, reco_tool::IPeakFitter::PeakFitParams_t::peakCenter, reco_tool::IPeakFitter::PeakFitParams_t::peakCenterError, reco_tool::IPeakFitter::PeakFitParams_t::peakSigma, reco_tool::IPeakFitter::PeakFitParams_t::peakSigmaError, and y.

58  {
59  if (fhc_vec.empty()) return;
60 
61  const float chiCut = 1e-3;
62  float lambda = 0.001; /* Marquardt damping parameter */
63  float chiSqr = std::numeric_limits<float>::max(), dchiSqr = std::numeric_limits<float>::max();
64  int nParams = 0;
65 
66  int startTime = fhc_vec.front().startTick;
67  int endTime = fhc_vec.back().stopTick;
68 
69  int roiSize = endTime - startTime;
70 
71  // init to a large number in case the fit fails
72  chi2PerNDF = chiSqr;
73 
74  std::vector<float> y(roiSize);
75  std::vector<float> p(3 * fhc_vec.size());
76  std::vector<float> plimmin(3 * fhc_vec.size());
77  std::vector<float> plimmax(3 * fhc_vec.size());
78  std::vector<float> perr(3 * fhc_vec.size());
79 
80  /* choose the fit function and set the parameters */
81  nParams = 0;
82 
83  for (size_t ih = 0; ih < fhc_vec.size(); ih++) {
84  float const peakMean =
85  fhc_vec[ih].hitCenter - 0.5 - (float)startTime; //shift by 0.5 to account for bin width
86  float const peakWidth = fhc_vec[ih].hitSigma;
87  float const amplitude = fhc_vec[ih].hitHeight;
88  float const meanLowLim = fmax(peakMean - fPeakRange * peakWidth, 0.);
89  float const meanHiLim = fmin(peakMean + fPeakRange * peakWidth, (float)roiSize);
90  p[0 + nParams] = amplitude;
91  p[1 + nParams] = peakMean;
92  p[2 + nParams] = peakWidth;
93 
94  plimmin[0 + nParams] = amplitude * 0.1;
95  plimmin[1 + nParams] = meanLowLim;
96  plimmin[2 + nParams] = std::max(fMinWidth, 0.1 * peakWidth);
97 
98  plimmax[0 + nParams] = amplitude * fAmpRange;
99  plimmax[1 + nParams] = meanHiLim;
100  plimmax[2 + nParams] = fMaxWidthMult * peakWidth;
101 
102  nParams += 3;
103  }
104  int fitResult = -1;
105 
106  for (size_t idx = 0; idx < size_t(roiSize); idx++) {
107  y[idx] = signal[startTime + idx];
108  }
109 
110  int trial = 0;
111  lambda = -1.; /* initialize lambda on first call */
112  do {
113  fitResult = fMarqFitAlg->gshf::MarqFitAlg::mrqdtfit(
114  lambda, &p[0], &plimmin[0], &plimmax[0], &y[0], nParams, roiSize, chiSqr, dchiSqr);
115  trial++;
116  if (fitResult || (trial > 100)) break;
117  } while (fabs(dchiSqr) >= chiCut);
118 
119  if (!fitResult) {
120  int fitResult2 =
121  fMarqFitAlg->gshf::MarqFitAlg::cal_perr(&p[0], &y[0], nParams, roiSize, &perr[0]);
122  if (!fitResult2) {
123  NDF = roiSize - nParams;
124  chi2PerNDF = chiSqr / NDF;
125  int parIdx = 0;
126  for (size_t i = 0; i < fhc_vec.size(); i++) {
127 
128  /* stand alone method
129  mhpp_vec[i].peakAmplitude = p[parIdx + 0];
130  mhpp_vec[i].peakAmplitudeError = perr[parIdx + 0];
131  mhpp_vec[i].peakCenter = p[parIdx + 1] + 0.5 + float(startTime);
132  mhpp_vec[i].peakCenterError = perr[parIdx + 1];
133  mhpp_vec[i].peakSigma = p[parIdx + 2];
134  mhpp_vec[i].peakSigmaError = perr[parIdx + 2];
135  */
136 
137  PeakFitParams_t mhpp;
138  mhpp.peakAmplitude = p[parIdx + 0];
139  mhpp.peakAmplitudeError = perr[parIdx + 0];
140  mhpp.peakCenter =
141  p[parIdx + 1] + 0.5 + float(startTime); //shift by 0.5 to account for bin width
142  mhpp.peakCenterError = perr[parIdx + 1];
143  mhpp.peakSigma = std::abs(p[parIdx + 2]);
144  mhpp.peakSigmaError = perr[parIdx + 2];
145 
146  mhpp_vec.emplace_back(mhpp);
147 
148  parIdx += 3;
149  }
150 
151  // fitStat=0;
152  }
153  }
154  return;
155  }
std::unique_ptr< gshf::MarqFitAlg > fMarqFitAlg
Float_t y
Definition: compare.C:6
constexpr auto abs(T v)
Returns the absolute value of the argument.
Float_t e
Definition: plot.C:35

Member Data Documentation

const double reco_tool::PeakFitterMrqdt::fAmpRange
private

Definition at line 34 of file PeakFitterMrqdt_tool.cc.

Referenced by findPeakParameters().

const geo::GeometryCore* reco_tool::PeakFitterMrqdt::fGeometry = lar::providerFrom<geo::Geometry>()
private

Definition at line 38 of file PeakFitterMrqdt_tool.cc.

std::unique_ptr<gshf::MarqFitAlg> reco_tool::PeakFitterMrqdt::fMarqFitAlg
private

Definition at line 36 of file PeakFitterMrqdt_tool.cc.

Referenced by findPeakParameters().

const double reco_tool::PeakFitterMrqdt::fMaxWidthMult
private

Definition at line 32 of file PeakFitterMrqdt_tool.cc.

Referenced by findPeakParameters().

const double reco_tool::PeakFitterMrqdt::fMinWidth
private

Definition at line 31 of file PeakFitterMrqdt_tool.cc.

Referenced by findPeakParameters().

const double reco_tool::PeakFitterMrqdt::fPeakRange
private

Definition at line 33 of file PeakFitterMrqdt_tool.cc.

Referenced by findPeakParameters().


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