kopia lustrzana https://github.com/vsamy/DiFipp
Change Base filter class.
rodzic
3b927b8174
commit
09a697bb5c
|
@ -0,0 +1,132 @@
|
||||||
|
// Copyright (c) 2019, Vincent SAMY
|
||||||
|
// All rights reserved.
|
||||||
|
|
||||||
|
// Redistribution and use in source and binary forms, with or without
|
||||||
|
// modification, are permitted provided that the following conditions are met:
|
||||||
|
|
||||||
|
// 1. Redistributions of source code must retain the above copyright notice,
|
||||||
|
// this list of conditions and the following disclaimer.
|
||||||
|
// 2. Redistributions in binary form must reproduce the above copyright notice,
|
||||||
|
// this list of conditions and the following disclaimer in the documentation
|
||||||
|
// and/or other materials provided with the distribution.
|
||||||
|
|
||||||
|
// THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND
|
||||||
|
// ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED
|
||||||
|
// WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
|
||||||
|
// DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR
|
||||||
|
// ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES
|
||||||
|
// (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
|
||||||
|
// LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND
|
||||||
|
// ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
|
||||||
|
// (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
|
||||||
|
// SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
|
||||||
|
|
||||||
|
// The views and conclusions contained in the software and documentation are those
|
||||||
|
// of the authors and should not be interpreted as representing official policies,
|
||||||
|
// either expressed or implied, of the FreeBSD Project.
|
||||||
|
|
||||||
|
#pragma once
|
||||||
|
|
||||||
|
#include "gsl/gsl_assert.h"
|
||||||
|
#include "type_checks.h"
|
||||||
|
#include "typedefs.h"
|
||||||
|
#include <stddef.h>
|
||||||
|
#include <string>
|
||||||
|
|
||||||
|
namespace difi {
|
||||||
|
|
||||||
|
// TODO: noexcept(Function of gsl variable)
|
||||||
|
// TODO: constructor with universal refs
|
||||||
|
|
||||||
|
/*! \brief Low-level filter.
|
||||||
|
*
|
||||||
|
* It creates the basic and common functions of all linear filter that can written as a digital filter.
|
||||||
|
* This class can not be instantiated directly.
|
||||||
|
*
|
||||||
|
* \warning In Debug mode, all functions may throw if a filter is badly initialized.
|
||||||
|
* This not the case in Realese mode.
|
||||||
|
*
|
||||||
|
* \tparam T Floating type.
|
||||||
|
*/
|
||||||
|
template <typename T, typename Derived>
|
||||||
|
class BaseFilter {
|
||||||
|
static_assert(std::is_floating_point<T>::value && !std::is_const<T>::value, "Only accept non-complex floating point types.");
|
||||||
|
friend Derived;
|
||||||
|
|
||||||
|
public:
|
||||||
|
/*! \brief Reset the data and filtered data. */
|
||||||
|
void resetFilter() noexcept { derived().resetFilter(); };
|
||||||
|
|
||||||
|
/*!< \brief Return the filter type */
|
||||||
|
FilterType type() const noexcept { return (m_center == 0 ? Type::Forward : Type::Centered); }
|
||||||
|
/*! \brief Get digital filter coefficients.
|
||||||
|
*
|
||||||
|
* It will automatically resize the given vectors.
|
||||||
|
* \param[out] aCoeff Denominator coefficients of the filter in decreasing order.
|
||||||
|
* \param[out] bCoeff Numerator coefficients of the filter in decreasing order.
|
||||||
|
*/
|
||||||
|
void getCoeffs(vectX_t<T>& aCoeff, vectX_t<T>& bCoeff) const noexcept;
|
||||||
|
/*! \brief Return coefficients of the denominator polynome. */
|
||||||
|
const vectX_t<T>& aCoeff() const noexcept { return m_aCoeff; }
|
||||||
|
/*! \brief Return coefficients of the numerator polynome. */
|
||||||
|
const vectX_t<T>& bCoeff() const noexcept { return m_bCoeff; }
|
||||||
|
/*! \brief Return the order the denominator polynome order of the filter. */
|
||||||
|
Eigen::Index aOrder() const noexcept { return m_aCoeff.size(); }
|
||||||
|
/*! \brief Return the order the numerator polynome order of the filter. */
|
||||||
|
Eigen::Index bOrder() const noexcept { return m_bCoeff.size(); }
|
||||||
|
/*! \brief Return the initialization state of the filter0 */
|
||||||
|
bool isInitialized() const noexcept { return m_isInitialized; }
|
||||||
|
|
||||||
|
protected:
|
||||||
|
/*! \brief Set type of filter (one-sided or centered)
|
||||||
|
*
|
||||||
|
* \param type The filter type.
|
||||||
|
* \warning bCoeff must be set before.
|
||||||
|
*/
|
||||||
|
void setType(FilterType type);
|
||||||
|
/*! \brief Set the new coefficients of the filters.
|
||||||
|
*
|
||||||
|
* It awaits a universal reference.
|
||||||
|
* \param aCoeff Denominator coefficients of the filter in decreasing order.
|
||||||
|
* \param bCoeff Numerator coefficients of the filter in decreasing order.
|
||||||
|
*/
|
||||||
|
template <typename T2>
|
||||||
|
void setCoeffs(T2&& aCoeff, T2&& bCoeff);
|
||||||
|
/*! \brief Normalized the filter coefficients such that aCoeff(0) = 1. */
|
||||||
|
void normalizeCoeffs();
|
||||||
|
/*! \brief Check for bad coefficients.
|
||||||
|
*
|
||||||
|
* Set the filter status to ready is everything is fine.
|
||||||
|
* \param aCoeff Denominator coefficients of the filter.
|
||||||
|
* \param bCoeff Numerator coefficients of the filter.
|
||||||
|
* \return True if the filter status is set on READY.
|
||||||
|
*/
|
||||||
|
bool checkCoeffs(const vectX_t<T>& aCoeff, const vectX_t<T>& bCoeff, Type type);
|
||||||
|
|
||||||
|
private:
|
||||||
|
/*! \brief Default uninitialized constructor. */
|
||||||
|
BaseFilter() = default;
|
||||||
|
/*! \brief Constructor.
|
||||||
|
* \param aCoeff Denominator coefficients of the filter in decreasing order.
|
||||||
|
* \param bCoeff Numerator coefficients of the filter in decreasing order.
|
||||||
|
* \param center
|
||||||
|
*/
|
||||||
|
BaseFilter(const vectX_t<T>& aCoeff, const vectX_t<T>& bCoeff, FilterType type = FilterType::Forward);
|
||||||
|
/*! \brief Default destructor. */
|
||||||
|
virtual ~BaseFilter() = default;
|
||||||
|
|
||||||
|
Derived& derived() noexcept { return *static_cast<Derived*>(this); }
|
||||||
|
const Derived& derived() const noexcept { return *static_cast<const Derived*>(this); }
|
||||||
|
|
||||||
|
private:
|
||||||
|
Eigen::Index m_center = 0; /*!< Center of the filter. 0 is a one-sided filter. Default is 0. */
|
||||||
|
bool m_isInitialized = false; /*!< Initialization state of the filter. Default is false */
|
||||||
|
vectX_t<T> m_aCoeff; /*!< Denominator coefficients of the filter */
|
||||||
|
vectX_t<T> m_bCoeff; /*!< Numerator coefficients of the filter */
|
||||||
|
vectX_t<T> m_filteredData; /*!< Last set of filtered data */
|
||||||
|
vectX_t<T> m_rawData; /*!< Last set of non-filtered data */
|
||||||
|
};
|
||||||
|
|
||||||
|
} // namespace difi
|
||||||
|
|
||||||
|
#include "BaseFilter.tpp"
|
|
@ -0,0 +1,96 @@
|
||||||
|
// Copyright (c) 2019, Vincent SAMY
|
||||||
|
// All rights reserved.
|
||||||
|
|
||||||
|
// Redistribution and use in source and binary forms, with or without
|
||||||
|
// modification, are permitted provided that the following conditions are met:
|
||||||
|
|
||||||
|
// 1. Redistributions of source code must retain the above copyright notice,
|
||||||
|
// this list of conditions and the following disclaimer.
|
||||||
|
// 2. Redistributions in binary form must reproduce the above copyright notice,
|
||||||
|
// this list of conditions and the following disclaimer in the documentation
|
||||||
|
// and/or other materials provided with the distribution.
|
||||||
|
|
||||||
|
// THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND
|
||||||
|
// ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED
|
||||||
|
// WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
|
||||||
|
// DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR
|
||||||
|
// ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES
|
||||||
|
// (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
|
||||||
|
// LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND
|
||||||
|
// ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
|
||||||
|
// (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
|
||||||
|
// SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
|
||||||
|
|
||||||
|
// The views and conclusions contained in the software and documentation are those
|
||||||
|
// of the authors and should not be interpreted as representing official policies,
|
||||||
|
// either expressed or implied, of the FreeBSD Project.
|
||||||
|
|
||||||
|
#include <limits>
|
||||||
|
|
||||||
|
namespace difi {
|
||||||
|
|
||||||
|
// Public functions
|
||||||
|
|
||||||
|
template <typename T>
|
||||||
|
void BaseFilter<T>::setType(Type type)
|
||||||
|
{
|
||||||
|
Expects(type == Type::Centered ? m_bCoeff.size() > 2 && m_bCoeff.size() % 2 == 1 : true);
|
||||||
|
m_center = (type == Type::Forward ? 0 : (m_bCoeff.size() - 1) / 2);
|
||||||
|
}
|
||||||
|
|
||||||
|
template <typename T>
|
||||||
|
template <typename T2>
|
||||||
|
void BaseFilter<T>::setCoeffs(T2&& aCoeff, T2&& bCoeff)
|
||||||
|
{
|
||||||
|
static_assert(std::is_convertible_v<T2, vectX_t<T>>, "The coefficients types should be convertible to vectX_t<T>");
|
||||||
|
|
||||||
|
Expects(checkCoeffs(aCoeff, bCoeff, (m_center == 0 ? Type::Forward : Type::Centered)));
|
||||||
|
m_aCoeff = aCoeff;
|
||||||
|
m_bCoeff = bCoeff;
|
||||||
|
normalizeCoeffs();
|
||||||
|
resetFilter();
|
||||||
|
m_isInitialized = true;
|
||||||
|
}
|
||||||
|
|
||||||
|
template <typename T>
|
||||||
|
void BaseFilter<T>::getCoeffs(vectX_t<T>& aCoeff, vectX_t<T>& bCoeff) const noexcept
|
||||||
|
{
|
||||||
|
aCoeff = m_aCoeff;
|
||||||
|
bCoeff = m_bCoeff;
|
||||||
|
}
|
||||||
|
|
||||||
|
// Protected functions
|
||||||
|
|
||||||
|
template <typename T>
|
||||||
|
BaseFilter<T>::BaseFilter(const vectX_t<T>& aCoeff, const vectX_t<T>& bCoeff, FilterType type)
|
||||||
|
: m_aCoeff(aCoeff)
|
||||||
|
, m_bCoeff(bCoeff)
|
||||||
|
, m_filteredData(aCoeff.size())
|
||||||
|
, m_rawData(bCoeff.size())
|
||||||
|
{
|
||||||
|
Expects(checkCoeffs(aCoeff, bCoeff, type));
|
||||||
|
m_center = (type == Type::Forward ? 0 : (bCoeff.size() - 1) / 2);
|
||||||
|
normalizeCoeffs();
|
||||||
|
resetFilter();
|
||||||
|
m_isInitialized = true;
|
||||||
|
}
|
||||||
|
|
||||||
|
template <typename T>
|
||||||
|
void BaseFilter<T>::normalizeCoeffs()
|
||||||
|
{
|
||||||
|
T a0 = m_aCoeff(0);
|
||||||
|
if (std::abs(a0 - T(1)) < std::numeric_limits<T>::epsilon())
|
||||||
|
return;
|
||||||
|
|
||||||
|
m_aCoeff /= a0;
|
||||||
|
m_bCoeff /= a0;
|
||||||
|
}
|
||||||
|
|
||||||
|
template <typename T>
|
||||||
|
bool BaseFilter<T>::checkCoeffs(const vectX_t<T>& aCoeff, const vectX_t<T>& bCoeff, FilterType type)
|
||||||
|
{
|
||||||
|
bool centering = (type == Type::Centered ? (bCoeff.size() % 2 == 1) : true);
|
||||||
|
return aCoeff.size() > 0 && std::abs(aCoeff[0]) > std::numeric_limits<T>::epsilon() && bCoeff.size() > 0 && centering;
|
||||||
|
}
|
||||||
|
|
||||||
|
} // namespace difi
|
|
@ -60,7 +60,7 @@ public:
|
||||||
/*! \brief Function to help you design a Butterworth filter.
|
/*! \brief Function to help you design a Butterworth filter.
|
||||||
*
|
*
|
||||||
* It finds optimal values of the order and cut-off frequency.
|
* It finds optimal values of the order and cut-off frequency.
|
||||||
* \warning Works only for low-pass and high-pass filters.
|
* \note Works only for low-pass and high-pass filters.
|
||||||
* \see http://www.matheonics.com/Tutorials/Butterworth.html#Paragraph_3.2
|
* \see http://www.matheonics.com/Tutorials/Butterworth.html#Paragraph_3.2
|
||||||
* \see https://www.mathworks.com/help/signal/ref/buttord.html#d120e11079
|
* \see https://www.mathworks.com/help/signal/ref/buttord.html#d120e11079
|
||||||
* \param wPass Pass band edge.
|
* \param wPass Pass band edge.
|
||||||
|
|
|
@ -50,6 +50,7 @@ public:
|
||||||
: GenericFilter<T>(aCoeff, bCoeff)
|
: GenericFilter<T>(aCoeff, bCoeff)
|
||||||
{
|
{
|
||||||
}
|
}
|
||||||
|
|
||||||
void setCoefficients(vectX_t<T>&& aCoeff, vectX_t<T>&& bCoeff)
|
void setCoefficients(vectX_t<T>&& aCoeff, vectX_t<T>&& bCoeff)
|
||||||
{
|
{
|
||||||
setCoeffs(std::forward(aCoeff), std::forward(bCoeff));
|
setCoeffs(std::forward(aCoeff), std::forward(bCoeff));
|
||||||
|
|
|
@ -27,37 +27,12 @@
|
||||||
|
|
||||||
#pragma once
|
#pragma once
|
||||||
|
|
||||||
#include "gsl/gsl_assert.h"
|
#include "BaseFilter.h"
|
||||||
#include "type_checks.h"
|
|
||||||
#include "typedefs.h"
|
|
||||||
#include <stddef.h>
|
|
||||||
#include <string>
|
|
||||||
|
|
||||||
namespace difi {
|
namespace difi {
|
||||||
|
|
||||||
// TODO: noexcept(Function of gsl variable)
|
|
||||||
// TODO: constructor with universal refs
|
|
||||||
|
|
||||||
/*! \brief Low-level filter.
|
|
||||||
*
|
|
||||||
* It creates the basic and common functions of all linear filter that can written as a digital filter.
|
|
||||||
* This class can not be instantiated directly.
|
|
||||||
*
|
|
||||||
* \warning In Debug mode, all functions may throw if a filter is badly initialized.
|
|
||||||
* This not the case in Realese mode.
|
|
||||||
*
|
|
||||||
* \tparam T Floating type.
|
|
||||||
*/
|
|
||||||
template <typename T>
|
template <typename T>
|
||||||
class GenericFilter {
|
class GenericFilter : BaseFilter<T, GenericFilter<T>> {
|
||||||
static_assert(std::is_floating_point<T>::value && !std::is_const<T>::value, "Only accept non-complex floating point types.");
|
|
||||||
|
|
||||||
public:
|
|
||||||
enum class Type {
|
|
||||||
OneSided,
|
|
||||||
Centered
|
|
||||||
};
|
|
||||||
|
|
||||||
public:
|
public:
|
||||||
/*! \brief Filter a new data.
|
/*! \brief Filter a new data.
|
||||||
*
|
*
|
||||||
|
@ -73,75 +48,48 @@ public:
|
||||||
* \return Filtered signal.
|
* \return Filtered signal.
|
||||||
*/
|
*/
|
||||||
vectX_t<T> filter(const vectX_t<T>& data);
|
vectX_t<T> filter(const vectX_t<T>& data);
|
||||||
/*! \brief Reset the data and filtered data. */
|
|
||||||
void resetFilter() noexcept;
|
|
||||||
|
|
||||||
/*!< \brief Return the filter type */
|
void resetFilter() noexcept
|
||||||
Type type() const noexcept { return (m_center == 0 ? Type::OneSided : Type::Centered); }
|
|
||||||
/*! \brief Get digital filter coefficients.
|
|
||||||
*
|
|
||||||
* It will automatically resize the given vectors.
|
|
||||||
* \param[out] aCoeff Denominator coefficients of the filter in decreasing order.
|
|
||||||
* \param[out] bCoeff Numerator coefficients of the filter in decreasing order.
|
|
||||||
*/
|
|
||||||
void getCoeffs(vectX_t<T>& aCoeff, vectX_t<T>& bCoeff) const noexcept;
|
|
||||||
/*! \brief Return coefficients of the denominator polynome. */
|
|
||||||
const vectX_t<T>& aCoeff() const noexcept { return m_aCoeff; }
|
|
||||||
/*! \brief Return coefficients of the numerator polynome. */
|
|
||||||
const vectX_t<T>& bCoeff() const noexcept { return m_bCoeff; }
|
|
||||||
/*! \brief Return the order the denominator polynome order of the filter. */
|
|
||||||
Eigen::Index aOrder() const noexcept { return m_aCoeff.size(); }
|
|
||||||
/*! \brief Return the order the numerator polynome order of the filter. */
|
|
||||||
Eigen::Index bOrder() const noexcept { return m_bCoeff.size(); }
|
|
||||||
/*! \brief Return the initialization state of the filter0 */
|
|
||||||
bool isInitialized() const noexcept { return m_isInitialized; }
|
|
||||||
|
|
||||||
protected:
|
protected:
|
||||||
/*! \brief Default uninitialized constructor. */
|
|
||||||
GenericFilter() = default;
|
GenericFilter() = default;
|
||||||
/*! \brief Constructor.
|
GenericFilter(const vectX_t<T>& aCoeff, const vectX_t<T>& bCoeff, Type type = Type::Forward)
|
||||||
* \param aCoeff Denominator coefficients of the filter in decreasing order.
|
: BaseFilter(aCoeff, bCoeff, type)
|
||||||
* \param bCoeff Numerator coefficients of the filter in decreasing order.
|
{}
|
||||||
* \param center
|
};
|
||||||
*/
|
|
||||||
GenericFilter(const vectX_t<T>& aCoeff, const vectX_t<T>& bCoeff, Type type = Type::OneSided);
|
|
||||||
/*! \brief Default destructor. */
|
|
||||||
virtual ~GenericFilter() = default;
|
|
||||||
|
|
||||||
/*! \brief Set type of filter (one-sided or centered)
|
template <typename T>
|
||||||
|
class TVGenericFilter : BaseFilter<T, GenericFilter<T>> {
|
||||||
|
public:
|
||||||
|
/*! \brief Filter a new data.
|
||||||
*
|
*
|
||||||
* \param type The filter type.
|
* This function is practical for online application that does not know the whole signal in advance.
|
||||||
* \warning bCoeff must be set before.
|
* \param data New data to filter.
|
||||||
|
* \return Filtered data.
|
||||||
*/
|
*/
|
||||||
void setType(Type type);
|
T stepFilter(const T& data, const T& time);
|
||||||
/*! \brief Set the new coefficients of the filters.
|
/*! \brief Filter a signal.
|
||||||
*
|
|
||||||
* It awaits a universal reference.
|
|
||||||
* \param aCoeff Denominator coefficients of the filter in decreasing order.
|
|
||||||
* \param bCoeff Numerator coefficients of the filter in decreasing order.
|
|
||||||
*/
|
|
||||||
template <typename T2>
|
|
||||||
void setCoeffs(T2&& aCoeff, T2&& bCoeff);
|
|
||||||
/*! \brief Normalized the filter coefficients such that aCoeff(0) = 1. */
|
|
||||||
void normalizeCoeffs();
|
|
||||||
/*! \brief Check for bad coefficients.
|
|
||||||
*
|
*
|
||||||
* Set the filter status to ready is everything is fine.
|
* Filter all data given by the signal.
|
||||||
* \param aCoeff Denominator coefficients of the filter.
|
* \param data Signal.
|
||||||
* \param bCoeff Numerator coefficients of the filter.
|
* \return Filtered signal.
|
||||||
* \return True if the filter status is set on READY.
|
|
||||||
*/
|
*/
|
||||||
bool checkCoeffs(const vectX_t<T>& aCoeff, const vectX_t<T>& bCoeff, Type type);
|
vectX_t<T> filter(const vectX_t<T>& data, const vectX_t<T>& time);
|
||||||
|
|
||||||
|
protected:
|
||||||
|
GenericFilter() = default;
|
||||||
|
GenericFilter(const vectX_t<T>& aCoeff, const vectX_t<T>& bCoeff, Type type = Type::Forward)
|
||||||
|
: BaseFilter(aCoeff, bCoeff, type)
|
||||||
|
{
|
||||||
|
m_timers.resize(m_bCoeffs.size());
|
||||||
|
m_timeDiffs.resize(m_bCoeffs.size());
|
||||||
|
}
|
||||||
|
|
||||||
private:
|
private:
|
||||||
Eigen::Index m_center = 0; /*!< Center of the filter. 0 is a one-sided filter. Default is 0. */
|
vectX_t<T> m_timers;
|
||||||
bool m_isInitialized = false; /*!< Initialization state of the filter. Default is false */
|
vectX_t<T> m_timeDiffs;
|
||||||
vectX_t<T> m_aCoeff; /*!< Denominator coefficients of the filter */
|
|
||||||
vectX_t<T> m_bCoeff; /*!< Numerator coefficients of the filter */
|
|
||||||
vectX_t<T> m_filteredData; /*!< Last set of filtered data */
|
|
||||||
vectX_t<T> m_rawData; /*!< Last set of non-filtered data */
|
|
||||||
};
|
};
|
||||||
|
|
||||||
} // namespace difi
|
} // namespace difi
|
||||||
|
|
||||||
#include "GenericFilter.tpp"
|
#include "GenericFilter.tpp"
|
|
@ -25,12 +25,8 @@
|
||||||
// of the authors and should not be interpreted as representing official policies,
|
// of the authors and should not be interpreted as representing official policies,
|
||||||
// either expressed or implied, of the FreeBSD Project.
|
// either expressed or implied, of the FreeBSD Project.
|
||||||
|
|
||||||
#include <limits>
|
|
||||||
|
|
||||||
namespace difi {
|
namespace difi {
|
||||||
|
|
||||||
// Public functions
|
|
||||||
|
|
||||||
template <typename T>
|
template <typename T>
|
||||||
T GenericFilter<T>::stepFilter(const T& data)
|
T GenericFilter<T>::stepFilter(const T& data)
|
||||||
{
|
{
|
||||||
|
@ -42,10 +38,10 @@ T GenericFilter<T>::stepFilter(const T& data)
|
||||||
for (Eigen::Index i = m_filteredData.size() - 1; i > 0; --i)
|
for (Eigen::Index i = m_filteredData.size() - 1; i > 0; --i)
|
||||||
m_filteredData(i) = m_filteredData(i - 1);
|
m_filteredData(i) = m_filteredData(i - 1);
|
||||||
|
|
||||||
m_rawData[0] = data;
|
m_rawData(0) = data;
|
||||||
m_filteredData[0] = 0;
|
m_filteredData(0) = 0;
|
||||||
m_filteredData[m_center] = m_bCoeff.dot(m_rawData) - m_aCoeff.dot(m_filteredData);
|
m_filteredData(m_center) = m_bCoeff.dot(m_rawData) - m_aCoeff.dot(m_filteredData);
|
||||||
return m_filteredData[m_center];
|
return m_filteredData(m_center);
|
||||||
}
|
}
|
||||||
|
|
||||||
template <typename T>
|
template <typename T>
|
||||||
|
@ -66,65 +62,54 @@ void GenericFilter<T>::resetFilter() noexcept
|
||||||
}
|
}
|
||||||
|
|
||||||
template <typename T>
|
template <typename T>
|
||||||
void GenericFilter<T>::setType(Type type)
|
T TVGenericFilter<T>::stepFilter(const T& data, const T& time)
|
||||||
{
|
{
|
||||||
Expects(type == Type::Centered ? m_bCoeff.size() > 2 && m_bCoeff.size() % 2 == 1 : true);
|
Expects(m_isInitialized);
|
||||||
m_center = (type == Type::OneSided ? 0 : (m_bCoeff.size() - 1) / 2);
|
|
||||||
|
// Slide data (can't use SIMD, but should be small)
|
||||||
|
for (Eigen::Index i = m_rawData.size() - 1; i > 0; --i) {
|
||||||
|
m_rawData(i) = m_rawData(i - 1);
|
||||||
|
m_timers(i) = m_timers(i - 1);
|
||||||
|
}
|
||||||
|
for (Eigen::Index i = m_filteredData.size() - 1; i > 0; --i)
|
||||||
|
m_filteredData(i) = m_filteredData(i - 1);
|
||||||
|
|
||||||
|
m_timers(0) = time;
|
||||||
|
if (m_center == 0) {
|
||||||
|
for (Eigen::Index i = m_rawData.size() - 1; i > 0; --i)
|
||||||
|
m_timeDiffs(i) = m_timeDiffs(i - 1);
|
||||||
|
m_timeDiffs(0) = time - m_timers(1);
|
||||||
|
} else {
|
||||||
|
const Eigen::Index S = data.size() - 1;
|
||||||
|
const Eigen::Index M = S / 2;
|
||||||
|
m_timeDiffs(M) = T(1);
|
||||||
|
for (Eigen::Index i = 0; i < M; ++i)
|
||||||
|
m_timeDiffs(i) = m_timers(i) - m_timers(S - i);
|
||||||
|
}
|
||||||
|
m_rawData(0) = data;
|
||||||
|
m_filteredData(0) = 0;
|
||||||
|
m_filteredData(m_center) = (m_bCoeff.cwiseQuotient(m_timeDiffs)).dot(m_rawData) - m_aCoeff.dot(m_filteredData);
|
||||||
|
return m_filteredData(m_center);
|
||||||
}
|
}
|
||||||
|
|
||||||
template <typename T>
|
template <typename T>
|
||||||
template <typename T2>
|
vectX_t<T> TVGenericFilter<T>::filter(const vectX_t<T>& data, const vectX_t<T>& time)
|
||||||
void GenericFilter<T>::setCoeffs(T2&& aCoeff, T2&& bCoeff)
|
|
||||||
{
|
{
|
||||||
static_assert(std::is_convertible_v<T2, vectX_t<T>>, "The coefficients types should be convertible to vectX_t<T>");
|
Expects(m_isInitialized);
|
||||||
|
Expects(data.size() == time.size());
|
||||||
Expects(checkCoeffs(aCoeff, bCoeff, (m_center == 0 ? Type::OneSided : Type::Centered)));
|
vectX_t<T> results(data.size());
|
||||||
m_aCoeff = aCoeff;
|
for (Eigen::Index i = 0; i < data.size(); ++i)
|
||||||
m_bCoeff = bCoeff;
|
results(i) = stepFilter(data(i), time(i));
|
||||||
normalizeCoeffs();
|
return results;
|
||||||
resetFilter();
|
|
||||||
m_isInitialized = true;
|
|
||||||
}
|
}
|
||||||
|
|
||||||
template <typename T>
|
template <typename T>
|
||||||
void GenericFilter<T>::getCoeffs(vectX_t<T>& aCoeff, vectX_t<T>& bCoeff) const noexcept
|
void GenericFilter<T>::resetFilter() noexcept
|
||||||
{
|
{
|
||||||
aCoeff = m_aCoeff;
|
m_filteredData.setZero(m_aCoeff.size());
|
||||||
bCoeff = m_bCoeff;
|
m_rawData.setZero(m_bCoeff.size());
|
||||||
}
|
m_timers.setZero(m_bCoeff.size());
|
||||||
|
m_timerDiffs.setZero(m_bCoeff.size());
|
||||||
// Protected functions
|
|
||||||
|
|
||||||
template <typename T>
|
|
||||||
GenericFilter<T>::GenericFilter(const vectX_t<T>& aCoeff, const vectX_t<T>& bCoeff, Type type)
|
|
||||||
: m_aCoeff(aCoeff)
|
|
||||||
, m_bCoeff(bCoeff)
|
|
||||||
, m_filteredData(aCoeff.size())
|
|
||||||
, m_rawData(bCoeff.size())
|
|
||||||
{
|
|
||||||
Expects(checkCoeffs(aCoeff, bCoeff, type));
|
|
||||||
m_center = (type == Type::OneSided ? 0 : (bCoeff.size() - 1) / 2);
|
|
||||||
normalizeCoeffs();
|
|
||||||
resetFilter();
|
|
||||||
m_isInitialized = true;
|
|
||||||
}
|
|
||||||
|
|
||||||
template <typename T>
|
|
||||||
void GenericFilter<T>::normalizeCoeffs()
|
|
||||||
{
|
|
||||||
T a0 = m_aCoeff(0);
|
|
||||||
if (std::abs(a0 - T(1)) < std::numeric_limits<T>::epsilon())
|
|
||||||
return;
|
|
||||||
|
|
||||||
m_aCoeff /= a0;
|
|
||||||
m_bCoeff /= a0;
|
|
||||||
}
|
|
||||||
|
|
||||||
template <typename T>
|
|
||||||
bool GenericFilter<T>::checkCoeffs(const vectX_t<T>& aCoeff, const vectX_t<T>& bCoeff, Type type)
|
|
||||||
{
|
|
||||||
bool centering = (type == Type::Centered ? (bCoeff.size() % 2 == 1) : true);
|
|
||||||
return aCoeff.size() > 0 && std::abs(aCoeff[0]) > std::numeric_limits<T>::epsilon() && bCoeff.size() > 0 && centering;
|
|
||||||
}
|
}
|
||||||
|
|
||||||
} // namespace difi
|
} // namespace difi
|
|
@ -226,7 +226,7 @@ public:
|
||||||
T timestep() const noexcept { return std::pow(bCoeff()(0) / CoeffGetter<T, N>()(0), T(1) / Order); }
|
T timestep() const noexcept { return std::pow(bCoeff()(0) / CoeffGetter<T, N>()(0), T(1) / Order); }
|
||||||
};
|
};
|
||||||
|
|
||||||
template <typename T, size_t N, template<class, class> typename CoeffGetter>
|
template <typename T, size_t N, int Order, template<class, class> typename CoeffGetter>
|
||||||
class CentralDifferentiator : public GenericFilter<T> {
|
class CentralDifferentiator : public GenericFilter<T> {
|
||||||
public:
|
public:
|
||||||
CentralDifferentiator()
|
CentralDifferentiator()
|
||||||
|
|
|
@ -43,4 +43,9 @@ using vectN_t = Eigen::Matrix<T, N, 1>; /*!< Fixed Eigen column-vector */
|
||||||
template <typename T>
|
template <typename T>
|
||||||
using vectXc_t = vectX_t<std::complex<T>>; /*!< Eigen complex column-vector */
|
using vectXc_t = vectX_t<std::complex<T>>; /*!< Eigen complex column-vector */
|
||||||
|
|
||||||
|
enum class FilterType {
|
||||||
|
Forward,
|
||||||
|
Centered
|
||||||
|
};
|
||||||
|
|
||||||
} // namespace difi
|
} // namespace difi
|
Ładowanie…
Reference in New Issue