Working generic filter with tests.

topic/diffentiators
Vincent Samy 2018-10-23 19:22:57 +09:00
rodzic 322b4e10db
commit 3c710edc1b
5 zmienionych plików z 198 dodań i 143 usunięć

Wyświetl plik

@ -1,12 +1,14 @@
# Version minimum
cmake_minimum_required(VERSION 3)
cmake_minimum_required(VERSION 2.8)
include(cmake/base.cmake)
include(cmake/eigen.cmake)
include(cmake/msvc-specific.cmake)
set(PROJECT_NAME fratio)
set(PROJECT_DESCRIPTION "Filter using rational transfer function")
set(PROJECT_URL "...")
SET(PROJECT_DEBUG_POSTFIX "_d")
#SET(CXX_DISABLE_WERROR True)
set(DOXYGEN_USE_MATHJAX "NO")
@ -16,9 +18,12 @@ set(CMAKE_CXX_STANDARD 14)
setup_project()
option(DISABLE_TESTS "Disable unit tests." OFF)
# for MSVC
if(MSVC)
set(CMAKE_CXX_FLAGS_RELEASE "${CMAKE_CXX_FLAGS_RELEASE} /W1 /MP")
set(CMAKE_MSVCIDE_RUN_PATH "\$(SolutionDir)/src/\$(Configuration)")
endif()
# Eigen
@ -27,4 +32,9 @@ search_for_eigen()
add_subdirectory(src)
if(NOT ${DISABLE_TESTS})
add_subdirectory(tests)
endif()
pkg_config_append_libs(${PROJECT_NAME})
setup_project_finalize()

Wyświetl plik

@ -1,46 +1,44 @@
#pragma once
#include <Eigen/Core>
#include <stddef.h>
#include <vector>
namespace msfc {
namespace fratio {
namespace filt {
class GenericFilter {
public:
GenericFilter() = default;
GenericFilter(size_t nData);
GenericFilter(size_t nData, const std::vector<double>& aCoeff, const std::vector<double>& bCoeff);
class GenericFilter {
public:
GenericFilter() = default;
GenericFilter(size_t nData);
GenericFilter(size_t nData, const std::vector<double>& aCoeff, const std::vector<double>& bCoeff);
void setNData(size_t nData);
void setCoeff(const std::vector<double>& aCoeff, const std::vector<double>& bCoeff);
void getCoeff(std::vector<double>& aCoeff, std::vector<double>& bCoeff) const noexcept;
size_t aOrder() const noexcept { return m_aCoeff.size(); }
size_t bOrder() const noexcept { return m_bCoeff.size(); }
void setNData(size_t nData);
void setCoeff(const std::vector<double>& aCoeff, const std::vector<double>& bCoeff);
void getCoeff(std::vector<double>& aCoeff, std::vector<double>& bCoeff) const noexcept;
size_t filterOrder() const noexcept { return m_order; }
// bool stepFilter(const Eigen::VectorXd& data);
// https://stackoverflow.com/questions/50511549/meaning-of-rational-transfer-function-underlying-matlab-filter-or-scipy-signal-f
double stepFilter(double data);
std::vector<double> filter(const std::vector<double>& data);
void resetFilter();
// bool stepFilter(const Eigen::VectorXd& data);
bool stepFilter(double data);
// Eigen::VectorXd results() const noexcept;
// std::vector<double> results() const noexcept { return m_filteredData; }
// Eigen::VectorXd results() const noexcept;
double results() const noexcept;
private:
void checkCoeffs(const std::vector<double>& aCoeff, const std::vector<double>& bCoeff);
void normalize();
void shiftData();
private:
void normalize();
void shiftData();
protected:
std::vector<double> m_aCoeff;
std::vector<double> m_bCoeff;
protected:
size_t m_nACoeffFilteredData;
size_t m_nBCoeffFilteredData;
std::vector<double> m_aCoeff;
std::vector<double> m_bCoeff;
std::vector<double> m_filteredData;
std::vector<double> m_rawData;
// Eigen::MatrixXd m_filteredData;
// Eigen::MatrixXd m_rawData;
};
std::vector<double> m_filteredData;
std::vector<double> m_rawData;
// Eigen::MatrixXd m_filteredData;
// Eigen::MatrixXd m_rawData;
};
} // namespace filt
} // namespace msfc
} // namespace fratio

Wyświetl plik

@ -1,126 +1,106 @@
#include "msfc/filter/GenericFilter.h"
#include "msfc/logging.h"
#include "GenericFilter.h"
namespace msfc {
#include <sstream>
namespace filt {
namespace fratio {
GenericFilter::GenericFilter(size_t nData)
// : m_filteredData(nData, order)
// , m_rawData(nData, order)
{
}
GenericFilter::GenericFilter(size_t nData)
{
}
GenericFilter::GenericFilter(size_t nData, const std::vector<double>& aCoeff, const std::vector<double>& bCoeff)
: m_aCoeff(aCoeff)
, m_bCoeff(bCoeff)
, m_filteredData(aCoeff.size())
, m_rawData(bCoeff.size())
{
assert(aCoeff.size() > 0);
assert(bCoeff.size() > 0);
normalize();
}
GenericFilter::GenericFilter(size_t nData, const std::vector<double>& aCoeff, const std::vector<double>& bCoeff)
: m_aCoeff(aCoeff)
, m_bCoeff(bCoeff)
, m_filteredData(aCoeff.size(), 0)
, m_rawData(bCoeff.size(), 0)
{
checkCoeffs(aCoeff, bCoeff);
normalize();
}
void GenericFilter::setNData(size_t nData)
{
// m_filteredData.resize(nData, m_filteredData.cols());
// m_rawData.resize(nData, nData.cols());
}
void GenericFilter::setNData(size_t nData)
{
// m_filteredData.resize(nData, m_filteredData.cols());
// m_rawData.resize(nData, nData.cols());
}
void GenericFilter::setCoeff(const std::vector<double>& aCoeff, const std::vector<double>& bCoeff)
{
assert(aCoeff.size() > 0);
assert(bCoeff.size() > 0);
void GenericFilter::setCoeff(const std::vector<double>& aCoeff, const std::vector<double>& bCoeff)
{
checkCoeffs(aCoeff, bCoeff);
m_nACoeffFilteredData = 0;
m_nBCoeffFilteredData = 0;
m_aCoeff = aCoeff;
m_bCoeff = bCoeff;
m_filteredData.resize(aCoeff.size());
m_rawData.resize(bCoeff.size());
m_aCoeff = aCoeff;
m_bCoeff = bCoeff;
resetFilter();
normalize();
}
normalize();
}
void GenericFilter::getCoeff(std::vector<double>& aCoeff, std::vector<double>& bCoeff) const noexcept
{
aCoeff = m_aCoeff;
bCoeff = m_bCoeff;
}
void GenericFilter::getCoeff(std::vector<double>& aCoeff, std::vector<double>& bCoeff) const noexcept
{
aCoeff = m_aCoeff;
bCoeff = m_bCoeff;
}
// bool GenericFilter::stepFilter(const Eigen::VectorXd& data)
// {
// if (m_filteredData.rows() != data.size()) {
// LOG_ERROR("Bad data size. Expected vector of size " << m_filteredData.rows() << ", got a vector of size " << data.size());
// return false;
// }
double GenericFilter::stepFilter(double data)
{
// Slide data
for (auto rit1 = m_rawData.rbegin(), rit2 = m_rawData.rbegin() + 1; rit2 != m_rawData.rend(); ++rit1, ++rit2)
*rit1 = *rit2;
for (auto rit1 = m_filteredData.rbegin(), rit2 = m_filteredData.rbegin() + 1; rit2 != m_filteredData.rend(); ++rit1, ++rit2)
*rit1 = *rit2;
// if (m_nFilteredData == 0) {
// m_rawData.row(0) = data;
// m_filteredData.row(0).noalias() = m_bCoeff(0) * data;
// ++m_nFilteredData;
// } else if (m_nFilteredData < m_order) {
// m_filteredData.row(m_nFilteredData).noalias() = m_bCoeff(0) * data;
// for (size_t i = 1; i <= m_nFilteredData; ++i) {
// m_filteredData.row(m_nFilteredData).noalias() += m_bCoeff(i) * m_rawData.row(i - 1);
// }
// ++m_nFilteredData;
// }
double filtData = 0.;
m_rawData[0] = data;
// if (m_nFilteredData >= m_order) {
for (size_t k = 0; k < m_bCoeff.size(); ++k)
filtData += m_bCoeff[k] * m_rawData[k];
for (size_t k = 1; k < m_aCoeff.size(); ++k)
filtData -= m_aCoeff[k] * m_filteredData[k];
// } else {
// for (size_t i = 0; i < m_nFilteredData; ++i) {
// }
m_filteredData[0] = filtData;
// ++m_nFilteredData;
// }
// }
return filtData;
}
// https://stackoverflow.com/questions/50511549/meaning-of-rational-transfer-function-underlying-matlab-filter-or-scipy-signal-f
bool GenericFilter::stepFilter(double data)
{
double filtData;
for (size_t i = 0; i < m_nBCoeffFilteredData; ++i)
filtData += m_bCoeff[i] * m_rawData[m_nBCoeffFilteredData - i];
for (size_t i = 1; i < m_nACoeffFilteredData; ++i)
filtData -= m_aCoeff[i] * m_filteredData[m_nACoeffFilteredData - i];
std::vector<double> GenericFilter::filter(const std::vector<double>& data)
{
std::vector<double> results;
results.reserve(data.size());
m_filteredData[m_nACoeffFilteredData] = filtData;
m_rawData[m_nBCoeffFilteredData] = data;
++m_nACoeffFilteredData;
++m_nBCoeffFilteredData;
for (double d : data)
results.emplace_back(stepFilter(d));
if (m_nACoeffFilteredData == m_filteredData.size()) {
double* fd = m_filteredData.data();
for (size_t i = 0; i < m_filteredData.size() - 1; ++i)
*(fd) = *(++fd);
return results;
}
--m_nACoeffFilteredData;
}
void GenericFilter::resetFilter()
{
m_filteredData.assign(m_aCoeff.size(), 0);
m_rawData.assign(m_bCoeff.size(), 0);
}
if (m_nBCoeffFilteredData == m_rawData.size()) {
double* rd = m_rawData.data();
for (size_t i = 0; i < m_rawData.size() - 1; ++i)
*(rd) = *(++rd);
void GenericFilter::normalize()
{
double a0 = m_aCoeff.front();
if (std::abs(a0) < 1e-6) // Divide by zero
throw std::invalid_argument("By filtering value for coefficient a0. Should be superior to 1e-6");
--m_nBCoeffFilteredData;
}
}
for (double& a : m_aCoeff)
a /= a0;
for (double& b : m_bCoeff)
b /= a0;
}
void GenericFilter::normalize()
{
double a0 = m_aCoeff.front();
if (std::abs(a0) < 1e-6) // Divide by zero
LOG_ERROR_AND_THROW(std::invalid_argument, "By filtering value for coefficient a0. Should be superior to 1e-6");
void GenericFilter::checkCoeffs(const std::vector<double>& aCoeff, const std::vector<double>& bCoeff)
{
std::stringstream err;
if (aCoeff.size() == 0)
err << "The size of coefficient 'a' should greater than 0\n";
if (bCoeff.size() == 0)
err << "The size of coefficient 'b' should greater than 0\n";
for (double& a : m_aCoeff)
a /= a0;
for (double& b : m_bCoeff)
b /= a0;
}
if (err.str().size() > 0)
throw std::runtime_error(err.str());
}
} // namespace filt
} // namespace msfc
} // namespace fratio

Wyświetl plik

@ -1,19 +1,25 @@
include(${CMAKE_SOURCE_DIR}/cmake/boost.cmake)
enable_testing()
set(BOOST_COMPONENTS system unit_test_framework)
set(BOOST_COMPONENTS unit_test_framework)
search_for_boost()
add_definitions(-DBOOST_TEST_DYN_LINK)
include_directories(${Boost_INCLUDE_DIRS})
include_directories(${CMAKE_SOURCE_DIR}/include)
link_directories(${Boost_LIBRARY_DIR} ${Boost_LIBRARY_DIRS})
macro(addTest testName)
add_executable(${testName} ${testName}.cpp systems.h tools.h)
if(NOT MSVC)
target_link_libraries(${name} PUBLIC ${Boost_LIBRARIES})
endif(MSVC)
add_executable(${testName} ${testName}.cpp)
if(MSVC)
target_link_libraries(${testName} PUBLIC ${PROJECT_NAME})
else()
target_link_libraries(${testName} PUBLIC ${Boost_LIBRARIES} ${PROJECT_NAME})
endif()
add_test(${testName}Unit ${testName})
# Adding a project configuration file (for MSVC only)
GENERATE_MSVC_DOT_USER_FILE(${name})
GENERATE_MSVC_DOT_USER_FILE(${testName})
endmacro(addTest)
addTest(FilterTest)

Wyświetl plik

@ -0,0 +1,61 @@
// This file is part of copra.
// copra is free software: you can redistribute it and/or
// modify
// it under the terms of the GNU Lesser General Public License as published by
// the Free Software Foundation, either version 3 of the License, or
// (at your option) any later version.
//
// copra is distributed in the hope that it will be useful,
// but WITHOUT ANY WARRANTY; without even the implied warranty of
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
// GNU Lesser General Public License for more details.
//
// You should have received a copy of the GNU Lesser General Public License
// along with copra. If not, see
// <http://www.gnu.org/licenses/>.
#define BOOST_TEST_MODULE FilterTest
#include "GenericFilter.h"
#include <Eigen/Core>
#include <boost/test/unit_test.hpp>
struct System1 {
std::vector<double> data = { 1., 2., 3., 4. };
std::vector<double> aCoeff = { 1., -0.99993717 };
std::vector<double> bCoeff = { 0.99996859, -0.99996859 };
std::vector<double> results = { 0.99996859, 1.999874351973491, 2.999717289867956, 3.999497407630634 };
};
BOOST_AUTO_TEST_CASE(FilterThrows)
{
BOOST_REQUIRE_THROW(fratio::GenericFilter(1, {}, { 1., 2. }), std::runtime_error);
BOOST_REQUIRE_THROW(fratio::GenericFilter(1, { 1., 2. }, {}), std::runtime_error);
BOOST_REQUIRE_THROW(fratio::GenericFilter(1, { 0. }, { 1. }), std::invalid_argument);
auto gf = fratio::GenericFilter();
BOOST_REQUIRE_THROW(gf.setCoeff({}, { 1., 2. }), std::runtime_error);
BOOST_REQUIRE_THROW(gf.setCoeff({ 1., 2. }, {}), std::runtime_error);
BOOST_REQUIRE_THROW(gf.setCoeff({ 0. }, { 1. }), std::invalid_argument);
}
BOOST_FIXTURE_TEST_CASE(SimpleSystem, System1)
{
auto gf = fratio::GenericFilter(1, aCoeff, bCoeff);
BOOST_REQUIRE_EQUAL(aCoeff.size(), gf.aOrder());
BOOST_REQUIRE_EQUAL(bCoeff.size(), gf.bOrder());
std::vector<double> filteredData;
for (double d : data)
filteredData.push_back(gf.stepFilter(d));
for (size_t i = 0; i < filteredData.size(); ++i)
BOOST_CHECK_SMALL(std::abs(filteredData[i] - results[i]), 1e-14);
gf.resetFilter();
filteredData = gf.filter(data);
for (size_t i = 0; i < filteredData.size(); ++i)
BOOST_CHECK_SMALL(std::abs(filteredData[i] - results[i]), 1e-14);
}