diff --git a/include/BilinearTransform.h b/include/BilinearTransform.h new file mode 100644 index 0000000..daf689f --- /dev/null +++ b/include/BilinearTransform.h @@ -0,0 +1,48 @@ +#pragma once + +#include + +namespace fratio { + +template +struct BilinearTransform { + using SubType = internal::complex_sub_type_t; + static_assert(std::is_floating_point::value, "This struct can only accept floating point types (real and complex)."); + + static void SToZ(SubType fs, const T& sPlanePole, T& zPlanePole); + static void SToZ(SubType fs, const std::vector& sPlanePoles, std::vector& zPlanePoles); + static void ZToS(SubType fs, const T& zPlanePole, T& sPlanePole); + static void ZToS(SubType fs, const std::vector& zPlanePoles, std::vector& sPlanePoles); +}; + +template +void BilinearTransform::SToZ(SubType fs, const T& sPlanePole, T& zPlanePole) +{ + T scalePole = sPlanePole / (2 * fs); + zPlanePole = (T(1) + scalePole) / (T(1) - scalePole); +} + +template +void BilinearTransform::SToZ(SubType fs, const std::vector& sPlanePoles, std::vector& zPlanePoles) +{ + assert(sPlanePoles.size() == zPlanePoles.size()); + for (size_t k = 0; k < sPlanePoles.size(); ++k) + SToZ(fs, sPlanePoles[k], zPlanePoles[k]); +} + +template +void BilinearTransform::ZToS(SubType fs, const T& zPlanePole, T& sPlanePole) +{ + T invPole = T(1) / zPlanePole; + sPlanePole = 2 * fs * (T(1) - invPole) / (T(1) + invPole); +} + +template +void BilinearTransform::ZToS(SubType fs, const std::vector& zPlanePoles, std::vector& sPlanePoles) +{ + assert(zPlanePoles.size() == sPlanePoles.size()); + for (size_t k = 0; k < sPlanePoles.size(); ++k) + ZToS(fs, zPlanePoles[k], sPlanePoles[k]); +} + +} // namespace fratio \ No newline at end of file diff --git a/include/Butterworth.h b/include/Butterworth.h index 0594a41..1ee1361 100644 --- a/include/Butterworth.h +++ b/include/Butterworth.h @@ -11,7 +11,8 @@ template class Butterworth : public GenericFilter { public: enum class Type { - LowPass + LowPass, + HighPass }; // public: @@ -26,6 +27,7 @@ private: void initialize(size_t order, T fc, T fs); void computeDigitalRep(); void updateCoeffSize(); + void transformFilter(); private: Type m_type; diff --git a/include/Butterworth.inl b/include/Butterworth.inl index 9697b2c..d80bc71 100644 --- a/include/Butterworth.inl +++ b/include/Butterworth.inl @@ -1,3 +1,4 @@ +#include "BilinearTransform.h" #include "polynome_functions.h" #include #include @@ -58,8 +59,7 @@ void Butterworth::computeDigitalRep() std::complex scalePole; for (size_t k = 1; k <= m_order; ++k) { scalePole = scaleFactor * std::complex(-std::sin(thetaK(k)), std::cos(thetaK(k))); - scalePole /= 2 * m_fs; - m_poles[k - 1] = (T(1) + scalePole) / (T(1) - scalePole); + BilinearTransform>::SToZ(m_fs, scalePole, m_poles[k - 1]); } std::vector> numPoles(m_order, std::complex(-1)); @@ -89,4 +89,17 @@ void Butterworth::updateCoeffSize() resetFilter(); } +template +void Butterworth::transformFilter() +{ + switch (m_type) { + case Type::HighPass: + break; + + case Type::LowPass: + default: + break; + } +} + } // namespace fratio \ No newline at end of file diff --git a/include/GenericFilter.h b/include/GenericFilter.h index d744153..106704b 100644 --- a/include/GenericFilter.h +++ b/include/GenericFilter.h @@ -6,11 +6,10 @@ namespace fratio { -template ::value && !std::is_const::value>> -class GenericFilter; - template -class GenericFilter { +class GenericFilter { + static_assert(std::is_floating_point::value && !std::is_const::value, "Only accept non-complex floating point types."); + public: T stepFilter(T data); std::vector filter(const std::vector& data); diff --git a/include/fratio.h b/include/fratio.h index 0fb47d3..8a98c23 100644 --- a/include/fratio.h +++ b/include/fratio.h @@ -7,6 +7,7 @@ namespace fratio { +// Filters using DigitalFilterf = DigitalFilter; using DigitalFilterd = DigitalFilter; using MovingAveragef = MovingAverage; @@ -14,6 +15,7 @@ using MovingAveraged = MovingAverage; using Butterworthf = Butterworth; using Butterworthd = Butterworth; +// Polynome helper functions using VietaAlgof = VietaAlgo; using VietaAlgod = VietaAlgo; using VietaAlgoi = VietaAlgo; @@ -21,4 +23,10 @@ using VietaAlgocf = VietaAlgo>; using VietaAlgocd = VietaAlgo>; using VietaAlgoci = VietaAlgo>; +// Bilinear transformation functions +using BilinearTransformf = BilinearTransform; +using BilinearTransformd = BilinearTransform; +using BilinearTransformcf = BilinearTransform>; +using BilinearTransformcd = BilinearTransform>; + } // namespace fratio \ No newline at end of file diff --git a/include/polynome_functions.h b/include/polynome_functions.h index 72c4e3a..8b62ac2 100644 --- a/include/polynome_functions.h +++ b/include/polynome_functions.h @@ -6,11 +6,10 @@ namespace fratio { -template ::value || is_complex_t::value>> -struct VietaAlgo; - template -struct VietaAlgo { +struct VietaAlgo { + static_assert(std::is_arithmetic>::value, "This struct can only accept arithmetic types or complex."); + // Vieta's computation: https://en.wikipedia.org/wiki/Vieta%27s_formulas static std::vector polyCoeffFromRoot(const std::vector& poles); }; diff --git a/include/type_checks.h b/include/type_checks.h index 2d71353..0a4cb8c 100644 --- a/include/type_checks.h +++ b/include/type_checks.h @@ -4,12 +4,29 @@ namespace fratio { -template -struct is_complex_t : public std::false_type { -}; +namespace internal { -template -struct is_complex_t> : public std::true_type { -}; + // Based on https://stackoverflow.com/a/41449096/9774052 + template + std::true_type is_tt_impl(std::complex); + std::false_type is_tt_impl(...); + + template + using is_complex = decltype(is_tt_impl(std::declval())); + + template + struct sub_type { + using type = T; + }; + + template + struct sub_type { + using type = typename T::value_type; + }; + + template + using complex_sub_type_t = typename sub_type::value>::type; + +} // namespace internal } // namespace fratio \ No newline at end of file diff --git a/tests/ButterworthFilterTests.cpp b/tests/ButterworthFilterTests.cpp index 0f455b9..c141238 100644 --- a/tests/ButterworthFilterTests.cpp +++ b/tests/ButterworthFilterTests.cpp @@ -8,7 +8,7 @@ DISABLE_CONVERSION_WARNING_BEGIN template struct System { - std::vector data = { 1., 2., 3., 4., 5., 6., 7., 8. }; + std::vector data = { 1, 2, 3, 4, 5, 6, 7, 8 }; size_t order = 5; T fc = 10; T fs = 100; @@ -88,7 +88,7 @@ BOOST_FIXTURE_TEST_CASE(BUTTERWORTH_FILTER_DOUBLE, System) bf2.setFilterParameters(order, fc, fs); filteredData.resize(0); - for (float d : data) + for (double d : data) filteredData.push_back(bf2.stepFilter(d)); for (size_t i = 0; i < filteredData.size(); ++i) diff --git a/tests/DigitalFilterTests.cpp b/tests/DigitalFilterTests.cpp index a1e6f78..71dd2d3 100644 --- a/tests/DigitalFilterTests.cpp +++ b/tests/DigitalFilterTests.cpp @@ -8,8 +8,8 @@ DISABLE_CONVERSION_WARNING_BEGIN template struct System { - std::vector data = { 1., 2., 3., 4. }; - std::vector aCoeff = { 1., -0.99993717 }; + std::vector data = { 1, 2, 3, 4 }; + std::vector aCoeff = { 1, -0.99993717 }; std::vector bCoeff = { 0.99996859, -0.99996859 }; std::vector results = { 0.99996859, 1.999874351973491, 2.999717289867956, 3.999497407630634 }; }; diff --git a/tests/GenericFilterTests.cpp b/tests/GenericFilterTests.cpp index 1ba3906..63f7af0 100644 --- a/tests/GenericFilterTests.cpp +++ b/tests/GenericFilterTests.cpp @@ -5,12 +5,12 @@ BOOST_AUTO_TEST_CASE(FilterThrows) { - BOOST_REQUIRE_THROW(fratio::DigitalFilterd({}, { 1., 2. }), std::runtime_error); - BOOST_REQUIRE_THROW(fratio::DigitalFilterd({ 1., 2. }, {}), std::runtime_error); - BOOST_REQUIRE_THROW(fratio::DigitalFilterd({ 0. }, { 1. }), std::invalid_argument); + BOOST_REQUIRE_THROW(fratio::DigitalFilterd({}, { 1, 2 }), std::runtime_error); + BOOST_REQUIRE_THROW(fratio::DigitalFilterd({ 1, 2 }, {}), std::runtime_error); + BOOST_REQUIRE_THROW(fratio::DigitalFilterd({ 0 }, { 1 }), std::invalid_argument); auto df = fratio::DigitalFilterd(); - BOOST_REQUIRE_THROW(df.setCoeff({}, { 1., 2. }), std::runtime_error); - BOOST_REQUIRE_THROW(df.setCoeff({ 1., 2. }, {}), std::runtime_error); - BOOST_REQUIRE_THROW(df.setCoeff({ 0. }, { 1. }), std::invalid_argument); + BOOST_REQUIRE_THROW(df.setCoeff({}, { 1, 2 }), std::runtime_error); + BOOST_REQUIRE_THROW(df.setCoeff({ 1, 2 }, {}), std::runtime_error); + BOOST_REQUIRE_THROW(df.setCoeff({ 0 }, { 1 }), std::invalid_argument); } diff --git a/tests/MovingAverageFilterTests.cpp b/tests/MovingAverageFilterTests.cpp index fea7437..3c55388 100644 --- a/tests/MovingAverageFilterTests.cpp +++ b/tests/MovingAverageFilterTests.cpp @@ -5,7 +5,7 @@ template struct System { - std::vector data = { 1., 2., 3., 4., 5., 6. }; + std::vector data = { 1, 2, 3, 4, 5, 6 }; size_t windowSize = 4; std::vector results = { 0.25, 0.75, 1.5, 2.5, 3.5, 4.5 }; }; diff --git a/tests/polynome_functions_tests.cpp b/tests/polynome_functions_tests.cpp index bdfac26..fe4ec3c 100644 --- a/tests/polynome_functions_tests.cpp +++ b/tests/polynome_functions_tests.cpp @@ -19,13 +19,13 @@ DISABLE_CONVERSION_WARNING_BEGIN template struct SystemFloat { std::vector data = { 0.32, -0.0518, 41.4, 0.89 }; - std::vector results = { 1., -42.558199999999999, 48.171601999999993, -9.181098159999999, -0.610759296 }; + std::vector results = { 1, -42.558199999999999, 48.171601999999993, -9.181098159999999, -0.610759296 }; }; template struct SystemCFloat { std::vector> data = { { 1.35, 0.2 }, { -1.5, 4.45 }, { 12.8, -3.36 }, { 5.156, 2.12 } }; - std::vector> results = { { 1., 0. }, { -17.806, -3.41 }, { 73.2776, 99.20074 }, { 101.857496, -444.634694 }, { -269.1458768, 388.7308864 } }; + std::vector> results = { { 1, 0 }, { -17.806, -3.41 }, { 73.2776, 99.20074 }, { 101.857496, -444.634694 }, { -269.1458768, 388.7308864 } }; }; DISABLE_CONVERSION_WARNING_END