peakingduck::core namespace

Defines all core parts of the library (algorithms, data structures).

Copyright

UK Atomic Energy Authority (UKAEA) - 2019-20

using peakingduck::core::Array1D = Eigen::Array<Scalar, Size, 1>
using peakingduck::core::Array1Di = Array1D<int>
using peakingduck::core::Array1Df = Array1D<float>
using peakingduck::core::Array1Dd = Array1D<double>
using peakingduck::core::DefaultType = double
using peakingduck::core::PeakList = std::vector<PeakInfo<ValueType>>
const int peakingduck::core::ArrayTypeDynamic = Eigen::Dynamic
template<typename T = DefaultType, int Size = ArrayTypeDynamic>
NumericalData<T, Size> peakingduck::core::combine(const NumericalData<T, ArrayTypeDynamic> &one, const NumericalData<T, ArrayTypeDynamic> &two)

Combine (concatenate) arrays into another.

template<typename T = DefaultType, int InputSize = ArrayTypeDynamic, int WindowSize = ArrayTypeDynamic>
NumericalData<T, WindowSize> peakingduck::core::window(const NumericalData<T, InputSize> &data, int centerindex, int nouter = 5, int ninner = 0, bool includeindex = true)

Given a list of values take nouter points either side of the index given and ignore ninner points.

Examples:

  1. values = [8, 2, 5, 2, 6, 6, 9, 23, 12] index = 4 nouter = 3 ninner = 0 includeindex = True

    => [2, 5, 2, 6, 6, 9, 23]

  2. values = [8, 2, 5, 2, 6, 6, 9, 23, 12] index = 4 nouter = 3 ninner = 0 includeindex = False

    => [2, 5, 2, 6, 9, 23]

  3. values = [8, 2, 5, 2, 6, 6, 9, 23, 12] index = 4 nouter = 3 ninner = 1 includeindex = True

    => [2, 5, 6, 9, 23]

  4. values = [8, 2, 5, 2, 6, 6, 9, 23, 12] index = 4 nouter = 3 ninner = 1 includeindex = False

    => [2, 5, 9, 23]

Therefore:

  • ninner >= 0

  • ninner <= nouter

  • index >= nouter

  • index < values.size()

It will clip at (0, len(values))

template<typename T = DefaultType, int Size = ArrayTypeDynamic>
struct peakingduck::core::peakingduck::core::ChunkedThresholdPeakFilter : public peakingduck::core::IProcess<T, Size>
#include <peaking.hpp>

Simple threshold local/chunked peak filter.

Public Functions

ChunkedThresholdPeakFilter(T percentThreshold, size_t chunkSize = 10)
NumericalData<T, Size> go(const NumericalData<T, Size> &data) const final override
template<typename T, template<typename> class crtpType>
class peakingduck::core::peakingduck::core::crtp
#include <crtp.hpp>

Represents the base CRTP class to allow objects to have certain abilities in their interface.

Public Functions

T &underlying()
T const &underlying() const
template<typename T = DefaultType, int Size = ArrayTypeDynamic>
struct peakingduck::core::peakingduck::core::GlobalThresholdPeakFilter : public peakingduck::core::IProcess<T, Size>
#include <peaking.hpp>

Simple threshold global peak filter.

Public Functions

GlobalThresholdPeakFilter(T percentThreshold)
NumericalData<T, Size> go(const NumericalData<T, Size> &data) const final override
template<typename XScalar, typename YScalar>
class peakingduck::core::peakingduck::core::Histogram
#include <spectral.hpp>

Represents a basic 1D histogram Energies vs values or Channel vs values.

Subclassed by peakingduck::core::Spectrum< XScalar, YScalar >

Public Functions

Histogram()
Histogram(const NumericalData<XScalar> &X, const NumericalData<YScalar> &Y)
Histogram(const Histogram &other)
Histogram(Histogram &&other) = default
Histogram &operator=(const Histogram &other) = default
Histogram &operator=(Histogram &&other) noexcept = default
~Histogram()
NumericalData<XScalar> X() const
NumericalData<YScalar> Y() const

Protected Attributes

NumericalData<XScalar> _X
NumericalData<YScalar> _Y
template<typename ValueType = DefaultType, int Size = ArrayTypeDynamic>
struct peakingduck::core::peakingduck::core::IPeakFinder
#include <peaking.hpp>

Interface for peak finding algorithms.

Operates on numerical data (filtered or unfiltered) Never mutates the input (always const process) returns a list of peaks - PeakList

Subclassed by peakingduck::core::SimplePeakFinder< ValueType, Size >

Public Functions

~IPeakFinder()
PeakList<ValueType> find(const NumericalData<ValueType, Size> &data) const = 0

Identifies potential peaks in the data.

template<typename T = DefaultType, int Size = ArrayTypeDynamic>
struct peakingduck::core::peakingduck::core::IProcess
#include <process.hpp>

Interface for all process algorithms.

Operates on numerical data Never mutates the input (always const process) returns a new numerical array

Subclassed by peakingduck::core::ChunkedThresholdPeakFilter< T, Size >, peakingduck::core::GlobalThresholdPeakFilter< T, Size >, peakingduck::core::MovingAveragePeakFilter< T, Size >, peakingduck::core::MovingAverageSmoother< T, Size >, peakingduck::core::WeightedMovingAverageSmoother< T, Size >

Public Functions

~IProcess()
NumericalData<T, Size> go(const NumericalData<T, Size> &data) const = 0
template<typename T = DefaultType, int Size = ArrayTypeDynamic>
struct peakingduck::core::peakingduck::core::IProcessManager
#include <process.hpp>

A general process manager interface.

Subclassed by peakingduck::core::SimpleProcessManager< T, Size >

Public Functions

~IProcessManager()
IProcessManager &append(const std::shared_ptr<IProcess<T, Size>> &process) = 0
NumericalData<T, Size> run(const NumericalData<T, Size> &data) const = 0
size_t size() const = 0
void reset() = 0
template<typename T = DefaultType, int Size = ArrayTypeDynamic>
struct peakingduck::core::peakingduck::core::MovingAveragePeakFilter : public peakingduck::core::IProcess<T, Size>
#include <peaking.hpp>

Simple moving average peak filter.

Public Functions

MovingAveragePeakFilter(int windowsize)
NumericalData<T, Size> go(const NumericalData<T, Size> &data) const final override
template<typename T = DefaultType, int Size = ArrayTypeDynamic>
struct peakingduck::core::peakingduck::core::MovingAverageSmoother : public peakingduck::core::IProcess<T, Size>
#include <smoothing.hpp>

Simple moving average smoother.

Can we support windowsize given at compile time too? It would certainly help with unit tests.

Public Functions

MovingAverageSmoother(int windowsize)
NumericalData<T, Size> go(const NumericalData<T, Size> &data) const final override
template<typename T = DefaultType, int Size = ArrayTypeDynamic>
struct peakingduck::core::peakingduck::core::NumericalData : private Array1D<T, Size>, public peakingduck::core::NumericalFunctions<NumericalData<T, Size>>
#include <numerical.hpp>

Represents a 1-dimensional data structure (basically a 1D Eigen array) Dynamic array - most use cases will be determined at runtime (I am assuming). We don’t want anyone to know we are using Eigen beyond this file, since (in theory) it should make it easier to change library if need be. We only really need the array datastructure from Eigen and not much else and instead of reinventing the wheel, we wrap Eigen array.

We wrap this with private inheritance on the Eigen type but there are a lot of methods to expose, easy to add when/if we need them.

Eigen array is pretty good, it has things like sqrt, exp on array coefficients, but we need to extend this to other functions, so we use CRTP to do this.

For all of this, you may ask why not just use Eigen and use an alias? Well for one, we don’t need all of Eigen just the array, and not all of the array type (we require a simpler interface). Additionally, at some point we may wish to use another data structure as std::array for example. In this case we just change the NumericalData class to wrap that instead. If we change the alias this could break existing interfaces and APIs, causing big changes later on. Since this datastructure is fundamental to everything we need to make sure that we have this sorted properly first!

Public Types

using value_type = T
using BaseEigenArray = Array1D<value_type, Size>

Public Functions

template<typename OtherDerived>
NumericalData(const Eigen::ArrayBase<OtherDerived> &other)
NumericalData()
NumericalData(const std::vector<value_type> &other)
template<typename OtherDerived>
NumericalData &operator=(const Eigen::ArrayBase<OtherDerived> &other)
template<typename OtherDerived>
NumericalData &operator=(const Eigen::EigenBase<OtherDerived> &other)
template<typename OtherDerived>
NumericalData &operator=(const Eigen::ReturnByValue<OtherDerived> &other)
NumericalData &operator=(const std::vector<T> &other)
NumericalData operator+(const T &scalar) const
NumericalData operator+(const NumericalData &rhs) const
const NumericalData &operator+=(const NumericalData &rhs)
NumericalData operator-(const T &scalar) const
NumericalData operator-(const NumericalData &rhs) const
const NumericalData &operator-=(const NumericalData &rhs)
NumericalData operator*(const T &scalar) const
NumericalData operator*(const NumericalData &rhs) const
const NumericalData &operator*=(const NumericalData &rhs)
NumericalData operator/(const T &scalar) const
NumericalData operator/(const NumericalData &rhs) const
const NumericalData &operator/=(const NumericalData &rhs)
void from_vector(const std::vector<value_type> &raw)
std::vector<value_type> to_vector() const
NumericalData<value_type> slice(int sindex, int eindex) const
NumericalData<value_type> operator()(int sindex, int eindex) const
NumericalData ramp(const value_type &threshold) const

A simple function for filtering values above a certain threshold (>=). This is useful to remove entries that are negative for example.

Returns a new array

NumericalData &rampInPlace(const value_type &threshold)

A simple function for filtering values above a certain threshold (>=). This is useful to remove entries that are negative for example.

Mutates underlying array

Friends

friend NumericalData operator+(const T &scalar, const NumericalData &rhs)
friend NumericalData operator-(const T &scalar, const NumericalData &rhs)
friend NumericalData operator*(const T &scalar, const NumericalData &rhs)
friend NumericalData operator/(const T &scalar, const NumericalData &rhs)
template<class Derived>
struct peakingduck::core::peakingduck::core::NumericalFunctions : public peakingduck::core::crtp<Derived, NumericalFunctions>
#include <numericalfunctions.hpp>

To extend the NumericalData type with certain numerical abilities, we add it to this using CRTP to keep it in the interface but it doesn’t require us to change the underlying data structure.

Public Functions

decltype(auto) stddev(int ddof = 0) const
Derived LLS() const

log(log(sqrt(value + 1) + 1) + 1) Returns a new array

Derived &LLSInPlace()

log(log(sqrt(value + 1) + 1) + 1) Changes the underlying array

Derived inverseLLS() const

exp(exp(sqrt(value + 1) + 1) + 1) Returns a new array

Derived &inverseLLSInPlace()

exp(exp(sqrt(value + 1) + 1) + 1) Changes the underlying array

Derived symmetricNeighbourOp(const std::function<void(int, int, const Derived&, Derived&)> &operation, int order = 1, ) const

For each element calculate a new value from the symmetric neigbour values at a given order. Take the i-order point and the i+order point and apply a function to that value End points are not counted (stay as original) - max(0, i-j) and min(i+j, len(array))

Returns a new array

Derived gradient(int order = 1) const

For each element calculate the numerical gradient value from the adjacent elements at a given order. Take the i-1 point and the i+1 point and determine the grad = (array[i+1]-array[i-1])/2.0. End points are handled differently as first point grad = array[1]-array[0] End points are handled differently as last point grad = array[-1]-array[-2].

For example, given an array: [ 1. , 2.0, 4.0, 7.0, 11.0, 16.0 ] The 1-st order gradient would be: [ 1. , 1.5, 2.5, 3.5, 4.5, 5. ] The 2-nd order gradient would be: [ 0.5, 0.75, 1.0, 1.0, 0.75, 0.5 ] The 3-rd order gradient would be: [ 0.25, 0.25, 0.125, -0.125, -0.25, -0.25 ]

Returns a new array

Derived &gradientInPlace(int order = 1)

Computes the numerical gradient in place.

See: NumericalFunctions::gradient

Mutates underlying array

Derived midpoint(int order = 1) const

For each element calculate the midpoint value from the adjacent elements at a given order. Take the i-order point and the i+order point and determine the average = (array[i-j]+array[i+j])/2.0. End points are not counted (stay as original) - max(0, i-j) and min(i+j, len(array))

For example, given an array: [1, 4, 6, 2, 4, 2, 5] we have the midpoints for order 0: [1, 4, 6, 2, 4, 2, 5] we have the midpoints for order 1: [1, 3.5, 3, 5, 2, 4.5, 5] we have the midpoints for order 2: [1, 4, 2.5, 3, 5.5, 2, 5] we have the midpoints for order 3: [1, 4, 6, 3, 4, 2, 5] we have the midpoints for order 4+: [1, 4, 6, 2, 4, 2, 5]

Returns a new array

Derived &midpointInPlace(int order = 1)

For each element calculate the midpoint value from the adjacent elements at a given order.

See: NumericalFunctions::midpoint

Mutates underlying array

template<class Iterator>
Derived snip(Iterator first, Iterator last) const

Sensitive Nonlinear Iterative Peak (SNIP) algorithm for estimating backgrounds ref needed here:

Allows any form of iterations given an iterator

Returns a new array

Derived snip(int niterations) const

Sensitive Nonlinear Iterative Peak (SNIP) algorithm for estimating backgrounds ref needed here:

does via increasing window only (ToDo: need to allow decreasing window)

Deprecate this!

Returns a new array

Derived &snipInPlace(int niterations)

Sensitive Nonlinear Iterative Peak (SNIP) algorithm for estimating backgrounds.

See: NumericalFunctions::snip

Mutates underlying array

template<typename ValueType = DefaultType>
struct peakingduck::core::peakingduck::core::PeakInfo
#include <peaking.hpp>

Simple struct for holding peak info.

Stores the: value = value i.e. count, flux index = the corresponding index/channel in the data

Subclassed by peakingduck::core::PeakInfoWithProp< ValueType, PropType >

Public Functions

PeakInfo(size_t pindex, ValueType pvalue)

Public Members

const size_t index
const ValueType value
template<typename ValueType = DefaultType, typename PropType = DefaultType>
struct peakingduck::core::peakingduck::core::PeakInfoWithProp : public peakingduck::core::PeakInfo<ValueType>
#include <peaking.hpp>

Struct extends basic PeakInfo with a property value.

Stores the: property = property i.e. energy, time value = value i.e. count, flux index = the corresponding index/channel in the data

Public Functions

PeakInfoWithProp(size_t pindex, ValueType pvalue, PropType pprop)

Public Members

const PropType property
template<typename ValueType = DefaultType, int Size = ArrayTypeDynamic>
struct peakingduck::core::peakingduck::core::SimplePeakFinder : public peakingduck::core::IPeakFinder<ValueType, Size>
#include <peaking.hpp>

Interface for peak finding algorithms.

Operates on numerical data (filtered or unfiltered) Never mutates the input (always const process) returns a list of peaks - PeakList

Public Functions

SimplePeakFinder(ValueType percentThreshold)
~SimplePeakFinder()
PeakList<ValueType> find(const NumericalData<ValueType, Size> &data) const override

Identifies potential peaks in the data based on a simple global threshold of max coefficient.

template<typename T = DefaultType, int Size = ArrayTypeDynamic>
struct peakingduck::core::peakingduck::core::SimpleProcessManager : public peakingduck::core::IProcessManager<T, Size>
#include <process.hpp>

A simple process manager.

Public Functions

~SimpleProcessManager()
IProcessManager<T, Size> &append(const std::shared_ptr<IProcess<T, Size>> &process) override
NumericalData<T, Size> run(const NumericalData<T, Size> &data) const override
size_t size() const override
void reset() override
template<typename XScalar, typename YScalar>
struct peakingduck::core::peakingduck::core::Spectrum : public peakingduck::core::Histogram<XScalar, YScalar>
#include <spectral.hpp>

Represents a basic 1D histogram Energies vs values or Channel vs values.

Public Functions

template<class Iterator>
void removeBackground(Iterator first, Iterator last)
template<class Iterator>
NumericalData<YScalar> estimateBackground(Iterator first, Iterator last) const
template<typename T = DefaultType, int Size = ArrayTypeDynamic>
struct peakingduck::core::peakingduck::core::WeightedMovingAverageSmoother : public peakingduck::core::IProcess<T, Size>
#include <smoothing.hpp>

Weighted moving average smoother.

Uses weights determined, with windowsize = N if N=1, weights=[1] -> [1/N] if N=2, weights=[1,1] -> [1/2, 1/2] if N=3, weights=[1,2,1] -> [1/4, 2/4, 1/4] if N=4, weights=[1,2,2,1] -> [1/6, 2/6, 2/6, 1/6] if N=5, weights=[1,2,3,2,1] -> [1/9, 2/9, 3/9, 2/9, 1/9] ….

Public Functions

WeightedMovingAverageSmoother(int windowsize)
NumericalData<T, Size> go(const NumericalData<T, Size> &data) const final override