From 2e123702513ab59efb5ac56564bf7178492fb399 Mon Sep 17 00:00:00 2001 From: Thomas Harte Date: Sat, 11 Nov 2017 12:20:11 -0500 Subject: [PATCH] Resolves some of the dangling C-isms remaining in my FIR filter, and introduces filter composition. --- .../Internals/Shaders/IntermediateShader.cpp | 4 +- SignalProcessing/FIRFilter.cpp | 98 ++++++++++++------- SignalProcessing/FIRFilter.hpp | 68 ++++++------- 3 files changed, 97 insertions(+), 73 deletions(-) diff --git a/Outputs/CRT/Internals/Shaders/IntermediateShader.cpp b/Outputs/CRT/Internals/Shaders/IntermediateShader.cpp index 0185afb98..5032a3d1b 100644 --- a/Outputs/CRT/Internals/Shaders/IntermediateShader.cpp +++ b/Outputs/CRT/Internals/Shaders/IntermediateShader.cpp @@ -334,9 +334,9 @@ void IntermediateShader::set_filter_coefficients(float sampling_rate, float cuto unsigned int taps = 11; // unsigned int taps = 21; // while(1) { - float coefficients[21]; +// float coefficients[21]; SignalProcessing::FIRFilter luminance_filter(taps, sampling_rate, 0.0f, cutoff_frequency, SignalProcessing::FIRFilter::DefaultAttenuation); - luminance_filter.get_coefficients(coefficients); + std::vector coefficients = luminance_filter.get_coefficients(); // int sample = 0; // int c = 0; diff --git a/SignalProcessing/FIRFilter.cpp b/SignalProcessing/FIRFilter.cpp index 8531ab3cc..1ce6a0ce2 100644 --- a/SignalProcessing/FIRFilter.cpp +++ b/SignalProcessing/FIRFilter.cpp @@ -46,7 +46,7 @@ float FIRFilter::ino(float a) { return s; } -void FIRFilter::coefficients_for_idealised_filter_response(short *filterCoefficients, float *A, float attenuation, unsigned int numberOfTaps) { +void FIRFilter::coefficients_for_idealised_filter_response(short *filter_coefficients, float *A, float attenuation, size_t number_of_taps) { /* calculate alpha, which is the Kaiser-Bessel window shape factor */ float a; // to take the place of alpha in the normal derivation @@ -59,46 +59,46 @@ void FIRFilter::coefficients_for_idealised_filter_response(short *filterCoeffici a = 0.5842f * powf(attenuation - 21.0f, 0.4f) + 0.7886f * (attenuation - 21.0f); } - float *filterCoefficientsFloat = new float[numberOfTaps]; + std::vector filter_coefficients_float(number_of_taps); /* work out the right hand side of the filter coefficients */ - unsigned int Np = (numberOfTaps - 1) / 2; + size_t Np = (number_of_taps - 1) / 2; float I0 = ino(a); - float NpSquared = static_cast(Np * Np); - for(unsigned int i = 0; i <= Np; i++) { - filterCoefficientsFloat[Np + i] = + float Np_squared = static_cast(Np * Np); + for(unsigned int i = 0; i <= Np; ++i) { + filter_coefficients_float[Np + i] = A[i] * - ino(a * sqrtf(1.0f - (static_cast(i * i) / NpSquared) )) / + ino(a * sqrtf(1.0f - (static_cast(i * i) / Np_squared) )) / I0; } /* coefficients are symmetrical, so copy from right hand side to left side */ - for(unsigned int i = 0; i < Np; i++) { - filterCoefficientsFloat[i] = filterCoefficientsFloat[numberOfTaps - 1 - i]; + for(size_t i = 0; i < Np; ++i) { + filter_coefficients_float[i] = filter_coefficients_float[number_of_taps - 1 - i]; } /* scale back up so that we retain 100% of input volume */ float coefficientTotal = 0.0f; - for(unsigned int i = 0; i < numberOfTaps; i++) { - coefficientTotal += filterCoefficientsFloat[i]; + for(size_t i = 0; i < number_of_taps; ++i) { + coefficientTotal += filter_coefficients_float[i]; } /* we'll also need integer versions, potentially */ float coefficientMultiplier = 1.0f / coefficientTotal; - for(unsigned int i = 0; i < numberOfTaps; i++) { - filterCoefficients[i] = (short)(filterCoefficientsFloat[i] * kCSKaiserBesselFilterFixedMultiplier * coefficientMultiplier); - } - - delete[] filterCoefficientsFloat; -} - -void FIRFilter::get_coefficients(float *coefficients) { - for(unsigned int i = 0; i < number_of_taps_; i++) { - coefficients[i] = static_cast(filter_coefficients_[i]) / kCSKaiserBesselFilterFixedMultiplier; + for(size_t i = 0; i < number_of_taps; ++i) { + filter_coefficients[i] = (short)(filter_coefficients_float[i] * FixedMultiplier * coefficientMultiplier); } } -FIRFilter::FIRFilter(unsigned int number_of_taps, float input_sample_rate, float low_frequency, float high_frequency, float attenuation) { +std::vector FIRFilter::get_coefficients() const { + std::vector coefficients; + for(auto short_coefficient: filter_coefficients_) { + coefficients.push_back(static_cast(short_coefficient) / FixedMultiplier); + } + return coefficients; +} + +FIRFilter::FIRFilter(size_t number_of_taps, float input_sample_rate, float low_frequency, float high_frequency, float attenuation) { // we must be asked to filter based on an odd number of // taps, and at least three if(number_of_taps < 3) number_of_taps = 3; @@ -108,30 +108,52 @@ FIRFilter::FIRFilter(unsigned int number_of_taps, float input_sample_rate, float number_of_taps |= 1; // store instance variables - number_of_taps_ = number_of_taps; - filter_coefficients_ = new short[number_of_taps_]; + filter_coefficients_.resize(number_of_taps); /* calculate idealised filter response */ - unsigned int Np = (number_of_taps - 1) / 2; - float twoOverSampleRate = 2.0f / input_sample_rate; + size_t Np = (number_of_taps - 1) / 2; + float two_over_sample_rate = 2.0f / input_sample_rate; - float *A = new float[Np+1]; + std::vector A(Np+1); A[0] = 2.0f * (high_frequency - low_frequency) / input_sample_rate; - for(unsigned int i = 1; i <= Np; i++) { - float iPi = static_cast(i) * static_cast(M_PI); + for(unsigned int i = 1; i <= Np; ++i) { + float i_pi = static_cast(i) * static_cast(M_PI); A[i] = ( - sinf(twoOverSampleRate * iPi * high_frequency) - - sinf(twoOverSampleRate * iPi * low_frequency) - ) / iPi; + sinf(two_over_sample_rate * i_pi * high_frequency) - + sinf(two_over_sample_rate * i_pi * low_frequency) + ) / i_pi; } - FIRFilter::coefficients_for_idealised_filter_response(filter_coefficients_, A, attenuation, number_of_taps_); - - /* clean up */ - delete[] A; + FIRFilter::coefficients_for_idealised_filter_response(filter_coefficients_.data(), A.data(), attenuation, number_of_taps); } -FIRFilter::~FIRFilter() { - delete[] filter_coefficients_; +FIRFilter::FIRFilter(const std::vector &coefficients) { + for(auto coefficient: coefficients) { + filter_coefficients_.push_back(static_cast(coefficient * FixedMultiplier)); + } +} + +FIRFilter FIRFilter::operator+(const FIRFilter &rhs) const { + std::vector coefficients = get_coefficients(); + std::vector rhs_coefficients = rhs.get_coefficients(); + + std::vector sum; + for(size_t i = 0; i < coefficients.size(); ++i) { + sum.push_back((coefficients[i] + rhs_coefficients[i]) / 2.0f); + } + + return FIRFilter(sum); +} + +FIRFilter FIRFilter::operator*(const FIRFilter &rhs) const { + std::vector coefficients = get_coefficients(); + std::vector rhs_coefficients = rhs.get_coefficients(); + + std::vector sum; + for(size_t i = 0; i < coefficients.size(); ++i) { + sum.push_back(coefficients[i] * rhs_coefficients[i]); + } + + return FIRFilter(sum); } diff --git a/SignalProcessing/FIRFilter.hpp b/SignalProcessing/FIRFilter.hpp index 8f7afaf39..d29fe652f 100644 --- a/SignalProcessing/FIRFilter.hpp +++ b/SignalProcessing/FIRFilter.hpp @@ -9,34 +9,24 @@ #ifndef FIRFilter_hpp #define FIRFilter_hpp -/* - - The FIR filter takes a 1d PCM signal with - a given sample rate and filters it according - to a specified filter (band pass only at - present, more to come if required). The number - of taps (ie, samples considered simultaneously - to make an output sample) is configurable; - smaller numbers permit a filter that operates - more quickly and with less lag but less - effectively. - - FIR filters are window functions; expected use is - to point sample an input that has been subject to - a filter. - -*/ - #ifdef __APPLE__ #include #endif +#include + namespace SignalProcessing { +/*! + The FIR filter takes a 1d PCM signal with a given sample rate and applies a band-pass filter to it. + + The number of taps (ie, samples considered simultaneously to make an output sample) is configurable; + smaller numbers permit a filter that operates more quickly and with less lag but less effectively. +*/ class FIRFilter { private: - static constexpr float kCSKaiserBesselFilterFixedMultiplier = 32767.0f; - static constexpr int kCSKaiserBesselFilterFixedShift = 15; + static constexpr float FixedMultiplier = 32767.0f; + static constexpr int FixedShift = 15; public: /*! @@ -48,9 +38,8 @@ class FIRFilter { @param high_frequency The highest frequency of signal to retain in the output. @param attenuation The attenuation of the discarded frequencies. */ - FIRFilter(unsigned int number_of_taps, float input_sample_rate, float low_frequency, float high_frequency, float attenuation); - - ~FIRFilter(); + FIRFilter(size_t number_of_taps, float input_sample_rate, float low_frequency, float high_frequency, float attenuation); + FIRFilter(const std::vector &coefficients); /*! A suggested default attenuation value. */ constexpr static float DefaultAttenuation = 60.0f; @@ -61,31 +50,44 @@ class FIRFilter { @param src The source buffer to apply the filter to. @returns The result of applying the filter. */ - inline short apply(const short *src) { + inline short apply(const short *src) const { #ifdef __APPLE__ short result; - vDSP_dotpr_s1_15(filter_coefficients_, 1, src, 1, &result, number_of_taps_); + vDSP_dotpr_s1_15(filter_coefficients_.data(), 1, src, 1, &result, filter_coefficients_.size()); return result; #else int outputValue = 0; - for(unsigned int c = 0; c < number_of_taps_; c++) { + for(size_t c = 0; c < filter_coefficients_.size(); ++c) { outputValue += filter_coefficients_[c] * src[c]; } - return (short)(outputValue >> kCSKaiserBesselFilterFixedShift); + return (short)(outputValue >> FixedShift); #endif } - inline unsigned int get_number_of_taps() { - return number_of_taps_; + /*! @returns The number of taps used by this filter. */ + inline size_t get_number_of_taps() const { + return filter_coefficients_.size(); } - void get_coefficients(float *coefficients); + /*! @returns The weighted coefficients that describe this filter. */ + std::vector get_coefficients() const; + + /*! + @returns A filter that would have the effect of adding (and scaling) the outputs of the two filters. + Defined only if both have the same number of taps. + */ + FIRFilter operator+(const FIRFilter &) const; + + /*! + @returns A filter that would have the effect of applying the two filters in succession. + Defined only if both have the same number of taps. + */ + FIRFilter operator*(const FIRFilter &) const; private: - short *filter_coefficients_; - unsigned int number_of_taps_; + std::vector filter_coefficients_; - static void coefficients_for_idealised_filter_response(short *filterCoefficients, float *A, float attenuation, unsigned int numberOfTaps); + static void coefficients_for_idealised_filter_response(short *filterCoefficients, float *A, float attenuation, size_t numberOfTaps); static float ino(float a); };