Belle II Software  release-05-02-19
NNWaveFitTool.h
1 /**************************************************************************
2  * BASF2 (Belle Analysis Framework 2) *
3  * Copyright(C) 2010 - Belle II Collaboration *
4  * *
5  * Author: The Belle II Collaboration *
6  * Contributors: Peter Kvasnicka *
7  * *
8  * This software is provided "as is" without any warranty. *
9  **************************************************************************/
10 
11 #pragma once
12 #ifndef _SVD_RECONSTRUCTION_NNWAVEFITTOOL_H
13 #define _SVD_RECONSTRUCTION_NNWAVEFITTOOL_H
14 
15 #include <vector>
16 #include <memory>
17 #include <tuple>
18 #include <algorithm>
19 #include <numeric>
20 #include <functional>
21 #include <svd/simulation/SVDSimulationTools.h>
22 
23 namespace Belle2 {
28  namespace SVD {
29 
32  typedef std::vector<double> nnFitterBinData;
33 
35  typedef std::vector<double> nnFitterBins;
36 
42  public:
48  m_bins(bins)
49  {
50  // Create table of probabiilities for EmpiricalDistributionFunction calculation
51  m_pTable.resize(m_bins.size());
52  auto ibin = m_bins.begin();
53  auto ipt = m_pTable.begin();
54  double cp = 0.0;
55  *ipt++ = std::make_pair(*ibin++, cp);
56  for (auto prob : p) {
57  cp += prob;
58  *ipt++ = std::make_pair(*ibin++, cp);
59  }
60  }
61 
66  double operator()(double t)
67  {
68  // pTable is sorted by .first
69  if (t > m_pTable.back().first) return 1.0;
70  if (t < m_pTable[0].first) return 0.0;
71  auto it = lower_bound(m_pTable.begin(), m_pTable.end(), std::make_pair(t, 0.0));
72  if (it == m_pTable.begin()) return it->second;
73  auto it2 = it;
74  --it2;
75  double result = it2->second + (it->second - it2->second) * (t - it2->first) / (it->first - it2->first);
76  return result;
77  }
78 
79  private:
81  std::vector<std::pair<double, double> > m_pTable;
83  }; // class EmpiricalDistributionFunction
84 
85  // ========================================================================
86  // NNWaveFitTool - object to handle computations with neural network fit.
87  // ------------------------------------------------------------------------
93  class NNWaveFitTool {
94  public:
100  NNWaveFitTool(const nnFitterBins& bins, const nnFitterBinData& binCenters,
101  const WaveGenerator& wave):
102  m_bins(bins), m_binCenters(binCenters), m_waveGenerator(wave)
103  {
104  m_altBinData.resize(m_binCenters.size());
105  }
106 
107  //------- Operations on probability arrays -------//
108 
112  const nnFitterBins& getBins() const { return m_bins; }
113 
117  const nnFitterBinData& getBinCenters() const { return m_binCenters; }
118 
125  void shiftInTime(nnFitterBinData& p, double timeShift);
126 
132  {
133  std::transform(p.begin(), p.end(),
134  p1.begin(), p.begin(), std::multiplies<double>());
135  normalize(p);
136  }
137 
144  std::shared_ptr<nnFitterBinData> pFromInterval(double left, double right);
145 
146  //------- Waveform parameters ------------------------------------//
147 
152  std::tuple<double, double> getTimeShift(const nnFitterBinData& p);
153 
162  std::tuple<double, double, double>
163  getAmplitudeChi2(const apvSamples& samples, double timeShift, double tau);
164 
171 
172  protected:
173 
181  {
182  double pnorm = std::accumulate(p.begin(), p.end(), 0.0);
183  // If the norm is too small, return default distribution.
184  if (pnorm < 1.0e-10) {
185  double uniformP = 1.0 / std::distance(p.begin(), p.end());
186  std::fill(p.begin(), p.end(), uniformP);
187  } else {
188  // Normalize if the norm makes sense.
189  std::transform(p.begin(), p.end(), p.begin(),
190  std::bind2nd(std::divides<double>(), pnorm));
191  }
192  }
193 
194  private:
200  }; // class NNWaveFitTool
201 
202  } // SVD namespace
204 } // Belle2 namespace
205 
206 #endif
Belle2::SVD::NNWaveFitTool::normalize
void normalize(nnFitterBinData &p)
Normalize the probability distribution.
Definition: NNWaveFitTool.h:180
Belle2::SVD::NNWaveFitTool::pFromInterval
std::shared_ptr< nnFitterBinData > pFromInterval(double left, double right)
Convert a uniform distribution to time-shift-like pdf.
Definition: NNWaveFitTool.cc:68
Belle2::SVD::EmpiricalDistributionFunction::m_pTable
std::vector< std::pair< double, double > > m_pTable
(bin, prob) pairs
Definition: NNWaveFitTool.h:81
Belle2::SVD::NNWaveFitTool::m_binCenters
const nnFitterBinData & m_binCenters
Centers of bins.
Definition: NNWaveFitTool.h:196
Belle2::SVD::EmpiricalDistributionFunction
Empirical distribution function object is basic for mainpulation of probabilities.
Definition: NNWaveFitTool.h:41
Belle2::SVD::WaveGenerator
Waveform generator This is a functor to calculate signal values.
Definition: SVDSimulationTools.h:102
Belle2::SVD::NNWaveFitTool::shiftInTime
void shiftInTime(nnFitterBinData &p, double timeShift)
Shift the probability array in time.
Definition: NNWaveFitTool.cc:54
Belle2::SVD::NNWaveFitTool::pLessThan
double pLessThan(nnFitterBinData p1, nnFitterBinData p2)
Return the probability X < Y, where X and Y are random variables with binned pdfs p1 and p2.
Definition: NNWaveFitTool.cc:80
Belle2::SVD::EmpiricalDistributionFunction::operator()
double operator()(double t)
EmpiricalDistributionFunction(t) gives edf value at time t, linearly interpolated from cummulative bi...
Definition: NNWaveFitTool.h:66
Belle2::SVD::EmpiricalDistributionFunction::m_bins
const nnFitterBins & m_bins
array of bin boundaries.
Definition: NNWaveFitTool.h:80
Belle2::SVD::NNWaveFitTool::NNWaveFitTool
NNWaveFitTool(const nnFitterBins &bins, const nnFitterBinData &binCenters, const WaveGenerator &wave)
Constructor takes arrays of signals and bin probabilities.
Definition: NNWaveFitTool.h:100
Belle2::SVD::nnFitterBins
std::vector< double > nnFitterBins
Vector of bin edges, nnFitterBinData.size() + 1.
Definition: NNWaveFitTool.h:35
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::SVD::apvSamples
std::array< apvSampleBaseType, nAPVSamples > apvSamples
vector od apvSample BaseType objects
Definition: SVDSimulationTools.h:41
Belle2::SVD::NNWaveFitTool::m_altBinData
nnFitterBinData m_altBinData
Bin data array for re-use.
Definition: NNWaveFitTool.h:197
Belle2::SVD::NNWaveFitTool::multiply
void multiply(nnFitterBinData &p, const nnFitterBinData &p1)
Multiply probabilities.
Definition: NNWaveFitTool.h:131
Belle2::SVD::NNWaveFitTool::getAmplitudeChi2
std::tuple< double, double, double > getAmplitudeChi2(const apvSamples &samples, double timeShift, double tau)
Return std::tuple with signal amplitude, its error, and chi2 of the fit.
Definition: NNWaveFitTool.cc:29
Belle2::SVD::EmpiricalDistributionFunction::EmpiricalDistributionFunction
EmpiricalDistributionFunction(const nnFitterBinData &p, const nnFitterBins &bins)
Constructor constructs edf object from a probability distribution.
Definition: NNWaveFitTool.h:47
Belle2::SVD::NNWaveFitTool::m_waveGenerator
WaveGenerator m_waveGenerator
APV signal generator.
Definition: NNWaveFitTool.h:199
Belle2::SVD::NNWaveFitTool::m_bins
const nnFitterBins & m_bins
Edges of bins.
Definition: NNWaveFitTool.h:195
Belle2::SVD::NNWaveFitTool::getBins
const nnFitterBins & getBins() const
Get edges of time shift bins.
Definition: NNWaveFitTool.h:112
Belle2::SVD::NNWaveFitTool::getBinCenters
const nnFitterBinData & getBinCenters() const
Get mean bin time shifts.
Definition: NNWaveFitTool.h:117
Belle2::SVD::NNWaveFitTool
The class holds arrays of bins and bin centers, and a wave generator object containing information on...
Definition: NNWaveFitTool.h:93
Belle2::SVD::nnFitterBinData
std::vector< double > nnFitterBinData
Vector of values defined for bins, such as bin times or bin probabilities.
Definition: NNWaveFitTool.h:32
Belle2::SVD::NNWaveFitTool::m_altSamples
apvSamples m_altSamples
Array of 6 apv samples for re-use.
Definition: NNWaveFitTool.h:198
Belle2::SVD::NNWaveFitTool::getTimeShift
std::tuple< double, double > getTimeShift(const nnFitterBinData &p)
Return std::tuple containing time shift and its error.
Definition: NNWaveFitTool.cc:18