All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Public Member Functions | Private Attributes | List of all members
reco_tool::PeakFitterGaussian Class Reference
Inheritance diagram for reco_tool::PeakFitterGaussian:
reco_tool::IPeakFitter

Public Member Functions

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

Private Attributes

const double fMinWidth
 minimum initial width for gaussian fit More...
 
const double fMaxWidthMult
 multiplier for max width for gaussian fit More...
 
const double fPeakRange
 set range limits for peak center More...
 
const double fAmpRange
 set range limit for peak amplitude More...
 
const bool fFloatBaseline
 Allow baseline to "float" away from zero. More...
 
const bool fOutputHistograms
 If true will generate summary style histograms. More...
 
TH1F * fNumCandHitsHist
 
TH1F * fROISizeHist
 
TH1F * fCandPeakPositionHist
 
TH1F * fCandPeakWidHist
 
TH1F * fCandPeakAmpitudeHist
 
TH1F * fCandBaselineHist
 
TH1F * fFitPeakPositionHist
 
TH1F * fFitPeakWidHist
 
TH1F * fFitPeakAmpitudeHist
 
TH1F * fFitBaselineHist
 
BaselinedGausFitCache fFitCache
 Preallocated ROOT functions for the fits. More...
 
TH1F fHistogram
 
const geo::GeometryCorefGeometry = lar::providerFrom<geo::Geometry>()
 
- Private Attributes inherited from reco_tool::IPeakFitter
float peakCenterError
 
float peakSigma
 
float peakSigmaError
 
float peakAmplitude
 
float peakAmplitudeError
 
float peakTauLeft
 
float peakTauLeftError
 
float peakTauRight
 
float peakTauRightError
 
float peakBaseline
 
float peakBaselineError
 

Additional Inherited Members

- Private Types inherited from reco_tool::IPeakFitter
using PeakFitParams_t = struct PeakFitParams{float peakCenter
 
using PeakParamsVec = std::vector< PeakFitParams_t >
 
using PeakParamsVec = std::vector< PeakFitParams_t >
 
- Private Member Functions inherited from reco_tool::IPeakFitter
virtual ~IPeakFitter () noexcept=default
 
virtual void configure (const fhicl::ParameterSet &pset)=0
 
virtual ~IPeakFitter ()=default
 

Detailed Description

Definition at line 58 of file PeakFitterGaussian_tool.cc.

Constructor & Destructor Documentation

reco_tool::PeakFitterGaussian::PeakFitterGaussian ( const fhicl::ParameterSet &  pset)
explicit

Definition at line 98 of file PeakFitterGaussian_tool.cc.

98  :
99  fMinWidth(pset.get<double>("MinWidth", 0.5)),
100  fMaxWidthMult (pset.get<double>("MaxWidthMult", 3.)),
101  fPeakRange(pset.get<double>("PeakRangeFact", 2.)),
102  fAmpRange(pset.get<double>("PeakAmpRange", 2.)),
103  fFloatBaseline(pset.get< bool >("FloatBaseline", false)),
104  fOutputHistograms(pset.get< bool >("OutputHistograms", false))
105 {
106  fHistogram = TH1F("PeakFitterHitSignal","",500,0.,500.);
107 
108  fHistogram.Sumw2();
109 
110  std::string function = "Gaus(0)";
111 
112  // If asked, define the global histograms
113  if (fOutputHistograms)
114  {
115  // Access ART's TFileService, which will handle creating and writing
116  // histograms and n-tuples for us.
117  art::ServiceHandle<art::TFileService> tfs;
118 
119  // Make a directory for these histograms
120  art::TFileDirectory dir = tfs->mkdir("PeakFit");
121 
122  fNumCandHitsHist = dir.make<TH1F>("NumCandHits", "# Candidate Hits", 100, 0., 100.);
123  fROISizeHist = dir.make<TH1F>("ROISize", "ROI Size", 400, 0., 400.);
124  fCandPeakPositionHist = dir.make<TH1F>("CPeakPosition", "Peak Position", 200, 0., 400.);
125  fCandPeakWidHist = dir.make<TH1F>("CPeadWidth", "Peak Width", 100, 0., 25.);
126  fCandPeakAmpitudeHist = dir.make<TH1F>("CPeakAmplitude", "Peak Amplitude", 100, 0., 200.);
127  fCandBaselineHist = dir.make<TH1F>("CBaseline", "Baseline", 200, -25., 25.);
128  fFitPeakPositionHist = dir.make<TH1F>("FPeakPosition", "Peak Position", 200, 0., 400.);
129  fFitPeakWidHist = dir.make<TH1F>("FPeadWidth", "Peak Width", 100, 0., 25.);
130  fFitPeakAmpitudeHist = dir.make<TH1F>("FPeakAmplitude", "Peak Amplitude", 100, 0., 200.);
131  fFitBaselineHist = dir.make<TH1F>("FBaseline", "Baseline", 200, -25., 25.);
132  }
133 
134  return;
135 }
const double fPeakRange
set range limits for peak center
const double fMaxWidthMult
multiplier for max width for gaussian fit
const bool fOutputHistograms
If true will generate summary style histograms.
const double fAmpRange
set range limit for peak amplitude
tuple dir
Definition: dropbox.py:28
const bool fFloatBaseline
Allow baseline to &quot;float&quot; away from zero.
art::ServiceHandle< art::TFileService > tfs
const double fMinWidth
minimum initial width for gaussian fit

Member Function Documentation

void reco_tool::PeakFitterGaussian::findPeakParameters ( const std::vector< float > &  roiSignalVec,
const ICandidateHitFinder::HitCandidateVec hitCandidateVec,
PeakParamsVec peakParamsVec,
double &  chi2PerNDF,
int &  NDF 
) const
overridevirtual

Implements reco_tool::IPeakFitter.

Definition at line 138 of file PeakFitterGaussian_tool.cc.

143 {
144  // The following is a translation of the original FitGaussians function in the original
145  // GausHitFinder module originally authored by Jonathan Asaadi
146  //
147  // *** NOTE: this algorithm assumes the reference time for input hit candidates is to
148  // the first tick of the input waveform (ie 0)
149  //
150  if (hitCandidateVec.empty()) return;
151 
152  // in case of a fit failure, set the chi-square to infinity
153  chi2PerNDF = std::numeric_limits<double>::infinity();
154 
155  int startTime = hitCandidateVec.front().startTick;
156  int endTime = hitCandidateVec.back().stopTick;
157  int roiSize = endTime - startTime;
158 
159  // Check to see if we need a bigger histogram for fitting
160  if (roiSize > fHistogram.GetNbinsX())
161  {
162  std::string histName = "PeakFitterHitSignal_" + std::to_string(roiSize);
163  fHistogram = TH1F(histName.c_str(),"",roiSize,0.,roiSize);
164  fHistogram.Sumw2();
165  }
166 
167  for(int idx = 0; idx < roiSize; idx++) fHistogram.SetBinContent(idx+1,roiSignalVec[startTime+idx]);
168 
169  // Build the string to describe the fit formula
170 #if 0
171  std::string equation = "gaus(0)";
172 
173  for(size_t idx = 1; idx < hitCandidateVec.size(); idx++) equation += "+gaus(" + std::to_string(3*idx) + ")";
174 
175  // Set the baseline if so desired
176  float baseline(0.);
177 
178  if (fFloatBaseline)
179  {
180  baseline = roiSignalVec[startTime];
181 
182  equation += "+" + std::to_string(baseline);
183  }
184 
185  // Now define the complete function to fit
186  TF1 Gaus("Gaus",equation.c_str(),0,roiSize,TF1::EAddToList::kNo);
187 #else
188  unsigned int const nGaus = hitCandidateVec.size();
189  assert(fFitCache.Get(nGaus));
190  TF1& Gaus = *(fFitCache.Get(nGaus));
191 
192  // Set the baseline if so desired
193  float baseline(0.);
194 
195  if (fFloatBaseline)
196  {
197  baseline = roiSignalVec[startTime];
198  Gaus.SetParameter(nGaus * 3, baseline);
199  Gaus.SetParLimits( nGaus * 3, baseline - 12., baseline + 12.);
200  }
201  else Gaus.FixParameter(nGaus * 3, baseline); // last parameter is the baseline
202 
203 #endif // 0
204 
205  if (fOutputHistograms)
206  {
207  fNumCandHitsHist->Fill(hitCandidateVec.size(), 1.);
208  fROISizeHist->Fill(roiSize, 1.);
209  fCandBaselineHist->Fill(baseline, 1.);
210  }
211 
212  // ### Setting the parameters for the Gaussian Fit ###
213  int parIdx{0};
214  for(auto const& candidateHit : hitCandidateVec)
215  {
216  double const peakMean = candidateHit.hitCenter - float(startTime);
217  double const peakWidth = candidateHit.hitSigma;
218  double const amplitude = candidateHit.hitHeight - baseline;
219  double const meanLowLim = std::max(peakMean - fPeakRange * peakWidth, 0.);
220  double const meanHiLim = std::min(peakMean + fPeakRange * peakWidth, double(roiSize));
221 
222  if (fOutputHistograms)
223  {
224  fCandPeakPositionHist->Fill(peakMean, 1.);
225  fCandPeakWidHist->Fill(peakWidth, 1.);
226  fCandPeakAmpitudeHist->Fill(amplitude, 1.);
227  }
228 
229  Gaus.SetParameter( parIdx, amplitude);
230  Gaus.SetParameter(1+parIdx, peakMean);
231  Gaus.SetParameter(2+parIdx, peakWidth);
232  Gaus.SetParLimits( parIdx, 0.1 * amplitude, fAmpRange * amplitude);
233  Gaus.SetParLimits(1+parIdx, meanLowLim, meanHiLim);
234  Gaus.SetParLimits(2+parIdx, std::max(fMinWidth, 0.1 * peakWidth), fMaxWidthMult * peakWidth);
235 
236  parIdx += 3;
237  }
238 
239  int fitResult{-1};
240 
241  try
242  { fitResult = fHistogram.Fit(&Gaus,"QNWB","", 0., roiSize);}
243  catch(...)
244  {mf::LogWarning("GausHitFinder") << "Fitter failed finding a hit";}
245 
246  // If the fit result is not zero there was an error
247  if (!fitResult)
248  {
249  // ##################################################
250  // ### Getting the fitted parameters from the fit ###
251  // ##################################################
252  chi2PerNDF = (Gaus.GetChisquare() / Gaus.GetNDF());
253  NDF = Gaus.GetNDF();
254 
255  int parIdx { 0 };
256  for(size_t idx = 0; idx < hitCandidateVec.size(); idx++)
257  {
258  PeakFitParams_t peakParams;
259 
260  peakParams.peakAmplitude = Gaus.GetParameter(parIdx);
261  peakParams.peakAmplitudeError = Gaus.GetParError( parIdx);
262  peakParams.peakCenter = Gaus.GetParameter(parIdx + 1) + float(startTime);
263  peakParams.peakCenterError = Gaus.GetParError( parIdx + 1);
264  peakParams.peakSigma = Gaus.GetParameter(parIdx + 2);
265  peakParams.peakSigmaError = Gaus.GetParError( parIdx + 2);
266 
267  if (fOutputHistograms)
268  {
269  fFitPeakPositionHist->Fill(peakParams.peakCenter, 1.);
270  fFitPeakWidHist->Fill(peakParams.peakSigma, 1.);
271  fFitPeakAmpitudeHist->Fill(peakParams.peakAmplitude, 1.);
272  }
273 
274  peakParamsVec.emplace_back(peakParams);
275 
276  parIdx += 3;
277  }
278 
279  if (fOutputHistograms) fFitBaselineHist->Fill(Gaus.GetParameter(3*nGaus), 1.);
280  }
281 #if 0
282  Gaus.Delete();
283 #endif // 0
284  return;
285 }
const double fPeakRange
set range limits for peak center
const double fMaxWidthMult
multiplier for max width for gaussian fit
const bool fOutputHistograms
If true will generate summary style histograms.
const double fAmpRange
set range limit for peak amplitude
BaselinedGausFitCache fFitCache
Preallocated ROOT functions for the fits.
BEGIN_PROLOG baseline
virtual TF1 * Get(size_t nFunc)
Returns a function sum of nFunc base functions.
const bool fFloatBaseline
Allow baseline to &quot;float&quot; away from zero.
std::string to_string(WindowPattern const &pattern)
const double fMinWidth
minimum initial width for gaussian fit

Member Data Documentation

const double reco_tool::PeakFitterGaussian::fAmpRange
private

set range limit for peak amplitude

Definition at line 74 of file PeakFitterGaussian_tool.cc.

TH1F* reco_tool::PeakFitterGaussian::fCandBaselineHist
private

Definition at line 83 of file PeakFitterGaussian_tool.cc.

TH1F* reco_tool::PeakFitterGaussian::fCandPeakAmpitudeHist
private

Definition at line 82 of file PeakFitterGaussian_tool.cc.

TH1F* reco_tool::PeakFitterGaussian::fCandPeakPositionHist
private

Definition at line 80 of file PeakFitterGaussian_tool.cc.

TH1F* reco_tool::PeakFitterGaussian::fCandPeakWidHist
private

Definition at line 81 of file PeakFitterGaussian_tool.cc.

TH1F* reco_tool::PeakFitterGaussian::fFitBaselineHist
private

Definition at line 87 of file PeakFitterGaussian_tool.cc.

BaselinedGausFitCache reco_tool::PeakFitterGaussian::fFitCache
mutableprivate

Preallocated ROOT functions for the fits.

Definition at line 89 of file PeakFitterGaussian_tool.cc.

TH1F* reco_tool::PeakFitterGaussian::fFitPeakAmpitudeHist
private

Definition at line 86 of file PeakFitterGaussian_tool.cc.

TH1F* reco_tool::PeakFitterGaussian::fFitPeakPositionHist
private

Definition at line 84 of file PeakFitterGaussian_tool.cc.

TH1F* reco_tool::PeakFitterGaussian::fFitPeakWidHist
private

Definition at line 85 of file PeakFitterGaussian_tool.cc.

const bool reco_tool::PeakFitterGaussian::fFloatBaseline
private

Allow baseline to "float" away from zero.

Definition at line 75 of file PeakFitterGaussian_tool.cc.

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

Definition at line 93 of file PeakFitterGaussian_tool.cc.

TH1F reco_tool::PeakFitterGaussian::fHistogram
mutableprivate

Definition at line 91 of file PeakFitterGaussian_tool.cc.

const double reco_tool::PeakFitterGaussian::fMaxWidthMult
private

multiplier for max width for gaussian fit

Definition at line 72 of file PeakFitterGaussian_tool.cc.

const double reco_tool::PeakFitterGaussian::fMinWidth
private

minimum initial width for gaussian fit

Definition at line 71 of file PeakFitterGaussian_tool.cc.

TH1F* reco_tool::PeakFitterGaussian::fNumCandHitsHist
private

Definition at line 78 of file PeakFitterGaussian_tool.cc.

const bool reco_tool::PeakFitterGaussian::fOutputHistograms
private

If true will generate summary style histograms.

Definition at line 76 of file PeakFitterGaussian_tool.cc.

const double reco_tool::PeakFitterGaussian::fPeakRange
private

set range limits for peak center

Definition at line 73 of file PeakFitterGaussian_tool.cc.

TH1F* reco_tool::PeakFitterGaussian::fROISizeHist
private

Definition at line 79 of file PeakFitterGaussian_tool.cc.


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