Belle II Software development
NNWaveFitter Class Reference

The class uses a neural network to find a probability distribution of arrival times for a sextet of APX25 signal samples. More...

#include <NNWaveFitter.h>

Public Types

typedef std::pair< double, double > nnBoundsType
 Bounds type used to hold network parameter bounds used in training the network.
 

Public Member Functions

 NNWaveFitter (std::string xmlData="")
 Constructor constructs the wavefitter from data in xml file.
 
void setNetwrok (const std::string &xmlData)
 Set proper network definition file.
 
std::shared_ptr< nnFitterBinDatagetFit (const apvSamples &samples, double tau)
 Fitting method Send data and get rseult structure.
 
const NNWaveFitToolgetFitTool () const
 Get a handle to a NNWaveFit object.
 
const nnFitterBinDatagetBinCenters () const
 Get bin times of the network output.
 
const nnFitterBinsgetBins () const
 Get bins of netwrok output.
 
nnBoundsType getWaveWidthBounds () const
 Get width bounds Width bounds are minimum and maximum width used in training the network.
 
nnBoundsType getAmplitudeBounds () const
 Get amplitude bounds.
 
nnBoundsType getTimeShiftBounds () const
 Get time shift bounds.
 
bool isValid () const
 Is this fitter working? Return false if the fitter was not initialized properly.
 
bool checkCoefficients (const std::string &dumpname, double tol=1.0e-10)
 Check NN data against a dump from Python.
 

Private Types

typedef std::function< double(double)> activationType
 Activation functions type.
 
typedef std::map< size_t, std::pair< Eigen::MatrixXd, Eigen::VectorXd > > networkWeightsType
 We use map to store network layers since we don't know if we'll be reading them in correct order.
 
typedef std::vector< Eigen::VectorXd > layerStatesType
 Storage for internal states.
 

Private Member Functions

Eigen::VectorXd softmax (const Eigen::VectorXd &input)
 Softmax function, normalization for the network's output layer.
 
int readNetworkData (const std::string &xmlFileName)
 The method that actually reads the xml file.
 

Private Attributes

activationType relu
 Rectifier activation.
 
activationType sigmoid
 Sigmoid activation.
 
bool m_isValid
 true if fitter was properly initialized
 
std::size_t m_nLayers
 number of NN layers, read from xml
 
nnFitterBinData m_binCenters
 centers of NN time bins
 
nnFitterBins m_bins
 NN time bin boundaries.
 
layerStatesType m_layerStates
 vectors of layer states
 
std::vector< std::size_t > m_layerSizes
 NN layer sizes.
 
networkWeightsType m_networkCoefs
 NN weights and intercepts.
 
activationType m_activation
 Network activation function.
 
nnBoundsType m_amplitudeBounds
 Amplitude range of the network.
 
nnBoundsType m_waveWidthBounds
 Waveform width range of the network.
 
nnBoundsType m_timeShiftBounds
 Time shift range of the network.
 
TauEncoder m_tauCoder
 Tau encoder class instance to scale tau values.
 
WaveformShape m_wave
 Wave function used in training the network.
 
std::shared_ptr< NNWaveFitToolm_fitTool
 FitterTool object allowing calculations on network fits.
 

Detailed Description

The class uses a neural network to find a probability distribution of arrival times for a sextet of APX25 signal samples.

The input to the network are sample S/N ratios, i.e., we expect sigma^2 = 1, plus width of the waveform from calibration data. The network(s) is currently trained using Python's scikit-learn package and results are saved in a PMML-like xml file. The NNWaveFitter reads the contents of this xml and on request provides the following data to the calling routine:

  1. Signal's arrival time estimate and its error,
  2. Signal's amplitude and its error,
  3. chi2 of the fit
  4. (Binned) probability distribution for arrival time.

The network is a multiclass classifier and the time bin probabilities are its main output. Time and amplitude are calculated from this distribution. The probability distribution is useful by itself in background suppression and u/v hit pairing.

The input xml file has, with a few extensions, the official PMML format (http://dmg.org), which can be created and read by most machine learning libraries.

The fitter generates a NNWaveFit object to aid in computations with bin probabilities, such as shifting, multiplication, calculations of time shift, amplitude and their errors.

Definition at line 61 of file NNWaveFitter.h.

Member Typedef Documentation

◆ activationType

typedef std::function<double(double)> activationType
private

Activation functions type.

Definition at line 143 of file NNWaveFitter.h.

◆ layerStatesType

typedef std::vector<Eigen::VectorXd> layerStatesType
private

Storage for internal states.

Definition at line 177 of file NNWaveFitter.h.

◆ networkWeightsType

typedef std::map< size_t, std::pair< Eigen::MatrixXd, Eigen::VectorXd > > networkWeightsType
private

We use map to store network layers since we don't know if we'll be reading them in correct order.

Definition at line 174 of file NNWaveFitter.h.

◆ nnBoundsType

typedef std::pair<double, double> nnBoundsType

Bounds type used to hold network parameter bounds used in training the network.

Also, this is the range guaranteed network applicability.

Definition at line 70 of file NNWaveFitter.h.

Constructor & Destructor Documentation

◆ NNWaveFitter()

NNWaveFitter ( std::string  xmlData = "")

Constructor constructs the wavefitter from data in xml file.

Parameters
xmlDatastring containing the network definition xml tree

Definition at line 260 of file NNWaveFitter.cc.

261{
263 if (xmlData != "") setNetwrok(xmlData);
264}
void setNetwrok(const std::string &xmlData)
Set proper network definition file.
WaveformShape m_wave
Wave function used in training the network.
Definition: NNWaveFitter.h:192
double w_betaprime(double t)
Beta-prime waveform shape, x^alpha/(1+x)^beta.

Member Function Documentation

◆ checkCoefficients()

bool checkCoefficients ( const std::string &  dumpname,
double  tol = 1.0e-10 
)

Check NN data against a dump from Python.

Parameters
dumpnameFilename of a text dump of network coefficients.
tolTolerance for float comparisons
Returns
true if all comparisons pass.

Definition at line 192 of file NNWaveFitter.cc.

193{
194 bool result = true;
195 ifstream dump(dumpname);
196
197 // Read the dumpfile line by line and compare coefficients.
198 string line;
199 getline(dump, line); // Header
200 // The dimension of the matrix is
201 // (m_layerSizes[iLayer], m_layerSizes[iLayer-1])
202 // and the dimension of the intercepts is
203 // m_layerSizes[iLayer]
204 for (size_t iLayer = 1; iLayer < m_nLayers; ++iLayer) {
205 getline(dump, line); // Layer
206 getline(dump, line); // Matrix
207 // the matrix is m_networkCoefs[iLayer].first, and has
208 // m_layerSizes[iLayer] rows with m_layerSizes[iLayer-1] entries
209 // in a row (columns)
210 for (size_t iRow = 0; iRow < m_layerSizes[iLayer]; ++iRow) {
211 getline(dump, line);
212 istringstream iline(line);
213 for (size_t iCol = 0; iCol < m_layerSizes[iLayer - 1]; ++iCol) {
214 double value = 0;
215 iline >> value;
216 if (fabs(value - m_networkCoefs[iLayer].first(iRow, iCol)) > tol) {
217 B2DEBUG(90, "Mismatch in weights in layer: " << iLayer
218 << " row: " << iRow << " col: " << iCol
219 << " C++: " << m_networkCoefs[iLayer].first(iRow, iCol)
220 << " Python: " << value
221 << " Diff: " << m_networkCoefs[iLayer].first(iRow, iCol) - value
222 );
223 result = false;
224 }
225 }
226 }
227
228 getline(dump, line); // Intercepts
229 // this is m_networkCoefs[iLayer].second, and has m_layerSizes[iLayer]
230 // numbers in a single row.
231 getline(dump, line);
232 istringstream iline(line);
233 for (size_t iRow = 0; iRow < m_layerSizes[iLayer]; ++iRow) {
234 double value = 0;
235 iline >> value;
236 if (fabs(value - m_networkCoefs[iLayer].second(iRow)) > tol) {
237 B2DEBUG(90, "Mismatch in intercepts in layer: " << iLayer
238 << " row: " << iRow
239 << " C++: " << m_networkCoefs[iLayer].second(iRow)
240 << " Python: " << value
241 << " Diff: " << m_networkCoefs[iLayer].second(iRow) - value
242 );
243
244 result = false;
245 }
246 }
247 }
248 return result;
249}
std::size_t m_nLayers
number of NN layers, read from xml
Definition: NNWaveFitter.h:181
std::vector< std::size_t > m_layerSizes
NN layer sizes.
Definition: NNWaveFitter.h:185
networkWeightsType m_networkCoefs
NN weights and intercepts.
Definition: NNWaveFitter.h:186

◆ getAmplitudeBounds()

nnBoundsType getAmplitudeBounds ( ) const
inline

Get amplitude bounds.

Returns
std::pair of amplitude bounds used in simulation

Definition at line 121 of file NNWaveFitter.h.

121{ return m_amplitudeBounds; }
nnBoundsType m_amplitudeBounds
Amplitude range of the network.
Definition: NNWaveFitter.h:188

◆ getBinCenters()

const nnFitterBinData & getBinCenters ( ) const
inline

Get bin times of the network output.

Returns
array of mean bin times (averaged over training data), properly sorted.

Definition at line 104 of file NNWaveFitter.h.

104{ return m_binCenters; }
nnFitterBinData m_binCenters
centers of NN time bins
Definition: NNWaveFitter.h:182

◆ getBins()

const nnFitterBins & getBins ( ) const
inline

Get bins of netwrok output.

Returns
array of bin edges (length of timebase + 1)

Definition at line 109 of file NNWaveFitter.h.

109{ return m_bins; }
nnFitterBins m_bins
NN time bin boundaries.
Definition: NNWaveFitter.h:183

◆ getFit()

std::shared_ptr< nnFitterBinData > getFit ( const apvSamples samples,
double  tau 
)

Fitting method Send data and get rseult structure.

If the fitter is not properly initialized, empty structure with valid = false will be returned, no warning - warning is only given for unsuccessful initialization.

Parameters
samplesArray of 6 apv samples
tauThe wave width for the current strip (unencoded ns)
Returns
Pointer to array of bin probabilities.

Definition at line 267 of file NNWaveFitter.cc.

268{
269 // Create the fit result object
270 std::shared_ptr<nnFitterBinData> result(new nnFitterBinData);
271 if (!m_isValid) return result;
272
273 copy(samples.begin(), samples.end(), &m_layerStates[0](0));
274 // Add tau in the last input cell
275 m_layerStates[0][samples.size()] = m_tauCoder.encodeTau(tau);
276 ostringstream os;
277 os << "Fitting with tau = " << tau << " encoded to " << m_tauCoder.encodeTau(tau) << endl
278 << "Samples: ";
279 copy(samples.begin(), samples.end(), ostream_iterator<double>(os, " "));
280 os << endl;
281 m_tauCoder.print(os);
282 B2DEBUG(100, os.str());
283 // Calculate network response
284 os.str() = "";
285 os << "Layer states: " << endl;
286 for (size_t iLayer = 1; iLayer < m_nLayers; ++iLayer) {
287 m_layerStates[iLayer] = (
288 m_networkCoefs[iLayer].first * m_layerStates[iLayer - 1]
289 + m_networkCoefs[iLayer].second
290 ).unaryExpr(relu);
291 os << m_layerStates[iLayer] << endl;
292 }
293 B2DEBUG(100, os.str());
294
295 // Output transformation 1: softmax normalization
296
298
299 result->resize(m_layerStates[m_nLayers - 1].size());
300 for (size_t i = 0; i < result->size(); ++i)
301 (*result)[i] = m_layerStates[m_nLayers - 1][i];
302
303 os.str() = "";
304 os << "Result: " << endl;
305 copy(result->begin(), result->end(), ostream_iterator<double>(os, " "));
306 os << endl;
307 B2DEBUG(100, os.str());
308
309 return result;
310}
TauEncoder m_tauCoder
Tau encoder class instance to scale tau values.
Definition: NNWaveFitter.h:191
layerStatesType m_layerStates
vectors of layer states
Definition: NNWaveFitter.h:184
activationType relu
Rectifier activation.
Definition: NNWaveFitter.h:145
bool m_isValid
true if fitter was properly initialized
Definition: NNWaveFitter.h:180
Eigen::VectorXd softmax(const Eigen::VectorXd &input)
Softmax function, normalization for the network's output layer.
Definition: NNWaveFitter.h:155
void print(std::ostringstream &os) const
print tool
double encodeTau(double tau) const
Return encoded value for a waveform width tau value.
std::vector< double > nnFitterBinData
Vector of values defined for bins, such as bin times or bin probabilities.
Definition: NNWaveFitTool.h:30

◆ getFitTool()

const NNWaveFitTool & getFitTool ( ) const
inline

Get a handle to a NNWaveFit object.

The object provides a set of tools to maniulate probabilities and calculate things from probability distributions.

Definition at line 98 of file NNWaveFitter.h.

98{ return *m_fitTool; }
std::shared_ptr< NNWaveFitTool > m_fitTool
FitterTool object allowing calculations on network fits.
Definition: NNWaveFitter.h:193

◆ getTimeShiftBounds()

nnBoundsType getTimeShiftBounds ( ) const
inline

Get time shift bounds.

Returns
std::pair of time shift bounds used in simulation

Definition at line 126 of file NNWaveFitter.h.

126{ return m_timeShiftBounds; }
nnBoundsType m_timeShiftBounds
Time shift range of the network.
Definition: NNWaveFitter.h:190

◆ getWaveWidthBounds()

nnBoundsType getWaveWidthBounds ( ) const
inline

Get width bounds Width bounds are minimum and maximum width used in training the network.

The bounds are used to encode width on input to the fitter.

Definition at line 116 of file NNWaveFitter.h.

116{ return m_waveWidthBounds; }
nnBoundsType m_waveWidthBounds
Waveform width range of the network.
Definition: NNWaveFitter.h:189

◆ isValid()

bool isValid ( ) const
inline

Is this fitter working? Return false if the fitter was not initialized properly.

Definition at line 131 of file NNWaveFitter.h.

131{ return m_isValid; }

◆ readNetworkData()

int readNetworkData ( const std::string &  xmlFileName)
private

The method that actually reads the xml file.

Counting things in NN:

Parameters
xmlFileNameName of the source xml file with network data
Returns
0: OK, -1: something went wrong.

Layers: input (hidden layers) output = n_hidden_lyaers + 2 = m_nLayers 0,1,2,... m_nLayers - 1. Layer sizes: input is 6, output is given by time binning in training, so we know terms 0 and m_nLayers - 1, and we can even establish hidden layer sizes beforehand from the repr of the classifier! Weights are about connections between layers, so we have w01, w12,.. w(m_nLayers-2)(m_nLayers-1) so we denote w1, w2, ... w(m_nLayers - 1) Intercepts are related to neurons where connections converge and activation is applied. So we have i1,...i(m_nLayers - 1).

Definition at line 36 of file NNWaveFitter.cc.

37{
38 using namespace boost::property_tree;
39 using boost::property_tree::ptree;
40
41 ptree propertyTree;
42
43 try {
44 B2DEBUG(400, "Reading xml");
45 stringstream ss;
46 ss << xmlData;
47 read_xml(ss, propertyTree);
48 } catch (const ptree_error& e) {
49 B2ERROR("Failed to parse xml data: " << e.what());
50 } catch (std::exception const& ex) {
51 B2ERROR("STD excpetion " << ex.what() << " in parsing.");
52 return -1;
53 }
54 // 1. Check that we have a classifier with relu activation.
55 try {
56 // 1a. Classifier
57 string function_name = propertyTree.get<string>("PMML.NeuralNetwork.<xmlattr>.functionName");
58 if (function_name != "classification") {
59 B2ERROR("This is incorrect network, expected multiclass classifier.");
60 return -1;
61 }
62 B2DEBUG(400, "Claasifier confirmed.");
63 // 1b. Activation
64 string activationFunction = propertyTree.get<string>("PMML.NeuralNetwork.<xmlattr>.activationFunction");
65 if (activationFunction != "rectifier") {
66 B2ERROR("Expected rectifier (relu) activation, found " << activationFunction << "instead.");
67 return -1;
68 }
69 B2DEBUG(400, "Activation confirmed.");
70 // 1c.Waveform type FIXME: no logic or checking here:
71 string waveFormType = propertyTree.get<string>("PMML.MiningBuildTask.Training.Waveform");
72 B2DEBUG(400, "Waveform set to " << waveFormType);
73 } catch (const ptree_error& e) {
74 B2ERROR("PropertyTree excpetion : " << e.what() << " in network check.");
75 return -1;
76 }
77
78 // 2. Read network topology
79 try {
80 string pathString = "PMML.MiningBuildTask.NetworkParameters";
81 map<size_t, size_t> m_layer_sizes_map;
82 for (const auto& layer_tag : propertyTree.get_child(pathString)) {
83 if (layer_tag.first == "<xmlattr>") continue;
84 size_t layer_no = layer_tag.second.get<size_t>("<xmlattr>.number");
85 size_t layer_size = layer_tag.second.get<size_t>("<xmlattr>.size");
86 m_layer_sizes_map.insert(make_pair(layer_no, layer_size));
87 }
88 // Once we have layer sizes, we can set weight matrices and intercept
89 // vectors to proper sizes.
90 m_nLayers = m_layer_sizes_map.size();
91 m_layerSizes.resize(m_nLayers);
92 for (size_t iLayer = 0; iLayer < m_nLayers; ++iLayer)
93 m_layerSizes[iLayer] = m_layer_sizes_map[iLayer];
94 for (size_t iLayer = 1; iLayer < m_nLayers; ++iLayer) {
95 m_networkCoefs.insert(make_pair(iLayer, make_pair(
96 MatrixXd(m_layerSizes[iLayer], m_layerSizes[iLayer - 1]),
97 VectorXd(m_layerSizes[iLayer])
98 )
99 ));
100 }
101 for (size_t iLayer = 0; iLayer < m_nLayers; ++iLayer)
102 m_layerStates.push_back(VectorXd(m_layerSizes[iLayer]));
103 B2DEBUG(400, "Network topology read.");
104 } catch (const ptree_error& e) {
105 B2ERROR("PropertyTree excpetion : " << e.what() << " when reading network topology.");
106 return -1;
107 }
108
109 // Read training parameter bounds
110 try {
111 string pathString("PMML.MiningBuildTask.Training");
112 for (const auto& param_tag : propertyTree.get_child(pathString)) {
113 if (param_tag.first == "Parameter") {
114 string valueString(param_tag.second.get<string>("<xmlattr>.value"));
115 double low = param_tag.second.get<double>("<xmlattr>.low");
116 double high = param_tag.second.get<double>("<xmlattr>.high");
117 if (valueString == "amplitude") {
118 m_amplitudeBounds.first = low;
119 m_amplitudeBounds.second = high;
120 } else if (valueString == "t0") {
121 m_timeShiftBounds.first = low * Unit::ns;
122 m_timeShiftBounds.second = high * Unit::ns;
123 } else if (valueString == "tau") {
124 m_waveWidthBounds.first = low * Unit::ns;
125 m_waveWidthBounds.second = high * Unit::ns;
126 }
127 }
128 }
129 // Set bounds in the tau encoder.
130 // Tau bounds are raw values!
133 B2DEBUG(400, "Read parameter bounds.");
134 } catch (const ptree_error& e) {
135 B2ERROR("PropertyTree excpetion: " << e.what() << " when reading parameter bounds.");
136 return -1;
137 }
138
139 // Reading output binning.
140 try {
141 string pathString("PMML.DataDictionary.DataField");
143 m_bins.resize(m_layerSizes[m_nLayers - 1] + 1);
144 for (const auto& value_tag : propertyTree.get_child(pathString)) {
145 if (value_tag.first == "Value") {
146 size_t i = value_tag.second.get<size_t>("<xmlattr>.value");
147 if (i == 1) // numbering is from 1!
148 m_bins[0] = value_tag.second.get<double>("<xmlattr>.lower") * Unit::ns;
149 m_bins[i] = value_tag.second.get<double>("<xmlattr>.upper") * Unit::ns;
150 m_binCenters[i - 1] = value_tag.second.get<double>("<xmlattr>.midpoint") * Unit::ns;
151 }
152 }
153 m_fitTool = std::shared_ptr<NNWaveFitTool>(new NNWaveFitTool(m_bins, m_binCenters, WaveGenerator(m_wave)));
154 B2DEBUG(400, "Outputs done.");
155 } catch (const ptree_error& e) {
156 B2ERROR("PropertyTree excpetion: " << e.what() << " when reading bin data.");
157 return -1;
158 }
159
160 // Reading neurons.
161 try {
162 string pathString("PMML.NeuralNetwork");
163 for (const auto& nl_tag : propertyTree.get_child(pathString)) {
164 if (nl_tag.first != "NeuralLayer") continue;
165 B2DEBUG(400, "Reading neural layers " << nl_tag.first << " " << nl_tag.second.size());
166 for (const auto& neuron_tag : nl_tag.second) {
167 if (neuron_tag.first != "Neuron") continue;
168 double bias = neuron_tag.second.get<double>("<xmlattr>.bias");
169 string sid = neuron_tag.second.get<string>("<xmlattr>.id");
170 size_t layer = stoi(sid.substr(0, 1)); // !! Won't work with > 9 !!
171 size_t pos = stoi(sid.substr(2, sid.size()));
172 B2DEBUG(400, "Reading neurons " << layer << "/" << pos << " bias: " << bias);
173 m_networkCoefs[layer].second(pos - 1) = bias;
174 for (const auto& con_tag : neuron_tag.second) {
175 if (con_tag.first != "Con") continue;
176 double weight = con_tag.second.get<double>("<xmlattr>.weight");
177 string sid2 = con_tag.second.get<string>("<xmlattr>.from");
178 size_t pos2 = stoi(sid2.substr(sid2.find('/') + 1, sid2.size()));
179 m_networkCoefs[layer].first(pos - 1, pos2 - 1) = weight;
180 B2DEBUG(400, "Reading connections " << sid2 << " weight: " << weight);
181 }
182 }
183 }
184 B2DEBUG(400, "Neurons done.");
185 } catch (const ptree_error& e) {
186 B2ERROR("PropertyTree excpetion: " << e.what() << " when reading neurons.");
187 return -1;
188 }
189 return 0;
190}
The class holds arrays of bins and bin centers, and a wave generator object containing information on...
Definition: NNWaveFitTool.h:91
void setBounds(double min_amplitude, double max_amplitude, double min_tau, double max_tau)
Set encoder baounds (e.g.
Waveform generator This is a functor to calculate APV samples from waveform.
static const double ns
Standard of [time].
Definition: Unit.h:48

◆ setNetwrok()

void setNetwrok ( const std::string &  xmlData)

Set proper network definition file.

Parameters
xmlDatastring containing the network definition xml tree

Definition at line 251 of file NNWaveFitter.cc.

252{
253 if ((xmlData != "") && (readNetworkData(xmlData) == 0)) {
254 m_isValid = true;
255 } else m_isValid = false;
256 // We don't issue any additional warnings here.
257}
int readNetworkData(const std::string &xmlFileName)
The method that actually reads the xml file.
Definition: NNWaveFitter.cc:36

◆ softmax()

Eigen::VectorXd softmax ( const Eigen::VectorXd &  input)
inlineprivate

Softmax function, normalization for the network's output layer.

This is a strange softmax that maps zero values back to zero, but this is what scikit-learn returns.

Definition at line 155 of file NNWaveFitter.h.

156 {
157 auto output = input.array().unaryExpr(
158 [](double x)->double { return ((x > 0.0) ? exp(x) : 0.0); }
159 );
160 double norm = output.sum();
161 return output / norm;
162 }

Member Data Documentation

◆ m_activation

activationType m_activation
private

Network activation function.

Definition at line 187 of file NNWaveFitter.h.

◆ m_amplitudeBounds

nnBoundsType m_amplitudeBounds
private

Amplitude range of the network.

Definition at line 188 of file NNWaveFitter.h.

◆ m_binCenters

nnFitterBinData m_binCenters
private

centers of NN time bins

Definition at line 182 of file NNWaveFitter.h.

◆ m_bins

nnFitterBins m_bins
private

NN time bin boundaries.

Definition at line 183 of file NNWaveFitter.h.

◆ m_fitTool

std::shared_ptr<NNWaveFitTool> m_fitTool
private

FitterTool object allowing calculations on network fits.

Definition at line 193 of file NNWaveFitter.h.

◆ m_isValid

bool m_isValid
private

true if fitter was properly initialized

Definition at line 180 of file NNWaveFitter.h.

◆ m_layerSizes

std::vector<std::size_t> m_layerSizes
private

NN layer sizes.

Definition at line 185 of file NNWaveFitter.h.

◆ m_layerStates

layerStatesType m_layerStates
private

vectors of layer states

Definition at line 184 of file NNWaveFitter.h.

◆ m_networkCoefs

networkWeightsType m_networkCoefs
private

NN weights and intercepts.

Definition at line 186 of file NNWaveFitter.h.

◆ m_nLayers

std::size_t m_nLayers
private

number of NN layers, read from xml

Definition at line 181 of file NNWaveFitter.h.

◆ m_tauCoder

TauEncoder m_tauCoder
private

Tau encoder class instance to scale tau values.

Definition at line 191 of file NNWaveFitter.h.

◆ m_timeShiftBounds

nnBoundsType m_timeShiftBounds
private

Time shift range of the network.

Definition at line 190 of file NNWaveFitter.h.

◆ m_wave

WaveformShape m_wave
private

Wave function used in training the network.

Definition at line 192 of file NNWaveFitter.h.

◆ m_waveWidthBounds

nnBoundsType m_waveWidthBounds
private

Waveform width range of the network.

Definition at line 189 of file NNWaveFitter.h.

◆ relu

activationType relu
private
Initial value:
=
[](double x) -> double { return std::max(double(0.0), x); }

Rectifier activation.

Definition at line 145 of file NNWaveFitter.h.

◆ sigmoid

activationType sigmoid
private
Initial value:
=
[](double x) -> double { double e = std::exp(x); return e / (1.0 + e); }

Sigmoid activation.

Definition at line 148 of file NNWaveFitter.h.


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