From 69d55d7e03996484cc899de1e21172a7a4532d6b Mon Sep 17 00:00:00 2001 From: Chris Robinson Date: Fri, 4 Dec 2020 09:42:13 -0800 Subject: Move the filters to core --- CMakeLists.txt | 12 +- al/source.cpp | 4 +- alc/alc.cpp | 4 +- alc/alcmain.h | 2 +- alc/alu.cpp | 6 +- alc/bformatdec.cpp | 12 +- alc/bformatdec.h | 14 +- alc/effects/autowah.cpp | 2 +- alc/effects/convolution.cpp | 4 +- alc/effects/distortion.cpp | 5 +- alc/effects/echo.cpp | 12 +- alc/effects/equalizer.cpp | 5 +- alc/effects/modulator.cpp | 5 +- alc/effects/reverb.cpp | 7 +- alc/filters/biquad.cpp | 163 ------------------- alc/filters/biquad.h | 144 ----------------- alc/filters/nfc.cpp | 383 -------------------------------------------- alc/filters/nfc.h | 63 -------- alc/filters/splitter.cpp | 113 ------------- alc/filters/splitter.h | 36 ----- alc/front_stablizer.h | 2 +- alc/hrtf.cpp | 2 +- alc/hrtf.h | 2 +- alc/voice.cpp | 6 +- alc/voice.h | 6 +- core/filters/biquad.cpp | 163 +++++++++++++++++++ core/filters/biquad.h | 144 +++++++++++++++++ core/filters/nfc.cpp | 383 ++++++++++++++++++++++++++++++++++++++++++++ core/filters/nfc.h | 63 ++++++++ core/filters/splitter.cpp | 113 +++++++++++++ core/filters/splitter.h | 36 +++++ 31 files changed, 955 insertions(+), 961 deletions(-) delete mode 100644 alc/filters/biquad.cpp delete mode 100644 alc/filters/biquad.h delete mode 100644 alc/filters/nfc.cpp delete mode 100644 alc/filters/nfc.h delete mode 100644 alc/filters/splitter.cpp delete mode 100644 alc/filters/splitter.h create mode 100644 core/filters/biquad.cpp create mode 100644 core/filters/biquad.h create mode 100644 core/filters/nfc.cpp create mode 100644 core/filters/nfc.h create mode 100644 core/filters/splitter.cpp create mode 100644 core/filters/splitter.h diff --git a/CMakeLists.txt b/CMakeLists.txt index 97004962..8c08af6b 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -633,6 +633,12 @@ set(ALC_OBJS core/bufferline.h core/devformat.cpp core/devformat.h + core/filters/biquad.h + core/filters/biquad.cpp + core/filters/nfc.cpp + core/filters/nfc.h + core/filters/splitter.cpp + core/filters/splitter.h core/mastering.cpp core/mastering.h @@ -674,12 +680,6 @@ set(ALC_OBJS alc/effects/pshifter.cpp alc/effects/reverb.cpp alc/effects/vmorpher.cpp - alc/filters/biquad.h - alc/filters/biquad.cpp - alc/filters/nfc.cpp - alc/filters/nfc.h - alc/filters/splitter.cpp - alc/filters/splitter.h alc/fmt_traits.cpp alc/fmt_traits.h alc/fpu_ctrl.cpp diff --git a/al/source.cpp b/al/source.cpp index 1af39da2..5dcd1cca 100644 --- a/al/source.cpp +++ b/al/source.cpp @@ -59,10 +59,10 @@ #include "backends/base.h" #include "bformatdec.h" #include "buffer.h" +#include "core/filters/nfc.h" +#include "core/filters/splitter.h" #include "event.h" #include "filter.h" -#include "filters/nfc.h" -#include "filters/splitter.h" #include "inprogext.h" #include "logging.h" #include "math_defs.h" diff --git a/alc/alc.cpp b/alc/alc.cpp index faf6392e..87c0578e 100644 --- a/alc/alc.cpp +++ b/alc/alc.cpp @@ -82,10 +82,10 @@ #include "compat.h" #include "core/devformat.h" #include "core/mastering.h" +#include "core/filters/nfc.h" +#include "core/filters/splitter.h" #include "cpu_caps.h" #include "effects/base.h" -#include "filters/nfc.h" -#include "filters/splitter.h" #include "fpu_ctrl.h" #include "front_stablizer.h" #include "hrtf.h" diff --git a/alc/alcmain.h b/alc/alcmain.h index 7949baf2..42808f26 100644 --- a/alc/alcmain.h +++ b/alc/alcmain.h @@ -25,7 +25,7 @@ #include "atomic.h" #include "core/bufferline.h" #include "core/devformat.h" -#include "filters/splitter.h" +#include "core/filters/splitter.h" #include "hrtf.h" #include "inprogext.h" #include "intrusive_ptr.h" diff --git a/alc/alu.cpp b/alc/alu.cpp index 18a26a38..c13365ec 100644 --- a/alc/alu.cpp +++ b/alc/alu.cpp @@ -61,12 +61,12 @@ #include "bs2b.h" #include "core/bsinc_tables.h" #include "core/devformat.h" +#include "core/filters/biquad.h" +#include "core/filters/nfc.h" +#include "core/filters/splitter.h" #include "core/mastering.h" #include "cpu_caps.h" #include "effects/base.h" -#include "filters/biquad.h" -#include "filters/nfc.h" -#include "filters/splitter.h" #include "fpu_ctrl.h" #include "front_stablizer.h" #include "hrtf.h" diff --git a/alc/bformatdec.cpp b/alc/bformatdec.cpp index b2a2aec9..6f229947 100644 --- a/alc/bformatdec.cpp +++ b/alc/bformatdec.cpp @@ -10,12 +10,10 @@ #include #include -#include "AL/al.h" - #include "almalloc.h" #include "alu.h" #include "ambdec.h" -#include "filters/splitter.h" +#include "core/filters/splitter.h" #include "front_stablizer.h" #include "math_defs.h" #include "opthelpers.h" @@ -33,7 +31,7 @@ constexpr std::array Ambi3DDecoderHFScale3O{{ 5.89792205e-01f, 8.79693856e-01f, 1.00000000e+00f, 1.00000000e+00f }}; -inline auto GetDecoderHFScales(ALuint order) noexcept -> const std::array& +inline auto GetDecoderHFScales(uint order) noexcept -> const std::array& { if(order >= 3) return Ambi3DDecoderHFScale3O; if(order == 2) return Ambi3DDecoderHFScale2O; @@ -52,7 +50,7 @@ inline auto GetAmbiScales(AmbDecScale scaletype) noexcept BFormatDec::BFormatDec(const AmbDecConf *conf, const bool allow_2band, const size_t inchans, - const ALuint srate, const ALuint (&chanmap)[MAX_OUTPUT_CHANNELS], + const uint srate, const uint (&chanmap)[MAX_OUTPUT_CHANNELS], std::unique_ptr stablizer) : mStablizer{std::move(stablizer)}, mDualBand{allow_2band && (conf->FreqBands == 2)} , mChannelDec{inchans} @@ -269,7 +267,7 @@ void BFormatDec::processStablize(const al::span OutBuffer, } -auto BFormatDec::GetHFOrderScales(const ALuint in_order, const ALuint out_order) noexcept +auto BFormatDec::GetHFOrderScales(const uint in_order, const uint out_order) noexcept -> std::array { std::array ret{}; @@ -286,7 +284,7 @@ auto BFormatDec::GetHFOrderScales(const ALuint in_order, const ALuint out_order) } std::unique_ptr BFormatDec::Create(const AmbDecConf *conf, const bool allow_2band, - const size_t inchans, const ALuint srate, const ALuint (&chanmap)[MAX_OUTPUT_CHANNELS], + const size_t inchans, const uint srate, const uint (&chanmap)[MAX_OUTPUT_CHANNELS], std::unique_ptr stablizer) { return std::unique_ptr{new(FamCount(inchans)) diff --git a/alc/bformatdec.h b/alc/bformatdec.h index 280ca2fc..adc69dc3 100644 --- a/alc/bformatdec.h +++ b/alc/bformatdec.h @@ -5,14 +5,12 @@ #include #include -#include "AL/al.h" - -#include "alcmain.h" #include "almalloc.h" #include "alspan.h" #include "ambidefs.h" +#include "core/bufferline.h" #include "core/devformat.h" -#include "filters/splitter.h" +#include "core/filters/splitter.h" struct AmbDecConf; struct FrontStablizer; @@ -44,7 +42,7 @@ class BFormatDec { public: BFormatDec(const AmbDecConf *conf, const bool allow_2band, const size_t inchans, - const ALuint srate, const ALuint (&chanmap)[MAX_OUTPUT_CHANNELS], + const uint srate, const uint (&chanmap)[MAX_OUTPUT_CHANNELS], std::unique_ptr stablizer); BFormatDec(const size_t inchans, const al::span coeffs, const al::span coeffslf, std::unique_ptr stablizer); @@ -61,11 +59,11 @@ public: const size_t SamplesToDo); /* Retrieves per-order HF scaling factors for "upsampling" ambisonic data. */ - static std::array GetHFOrderScales(const ALuint in_order, - const ALuint out_order) noexcept; + static std::array GetHFOrderScales(const uint in_order, + const uint out_order) noexcept; static std::unique_ptr Create(const AmbDecConf *conf, const bool allow_2band, - const size_t inchans, const ALuint srate, const ALuint (&chanmap)[MAX_OUTPUT_CHANNELS], + const size_t inchans, const uint srate, const uint (&chanmap)[MAX_OUTPUT_CHANNELS], std::unique_ptr stablizer); static std::unique_ptr Create(const size_t inchans, const al::span coeffs, const al::span coeffslf, diff --git a/alc/effects/autowah.cpp b/alc/effects/autowah.cpp index 5ac51e32..f2ffab44 100644 --- a/alc/effects/autowah.cpp +++ b/alc/effects/autowah.cpp @@ -29,7 +29,7 @@ #include "alcmain.h" #include "alcontext.h" #include "alu.h" -#include "filters/biquad.h" +#include "core/filters/biquad.h" #include "vecmat.h" namespace { diff --git a/alc/effects/convolution.cpp b/alc/effects/convolution.cpp index e191e7bc..2442f97e 100644 --- a/alc/effects/convolution.cpp +++ b/alc/effects/convolution.cpp @@ -8,7 +8,6 @@ #include "AL/al.h" #include "AL/alc.h" -#include "al/auxeffectslot.h" #include "alcmain.h" #include "alcomplex.h" #include "alcontext.h" @@ -17,8 +16,9 @@ #include "ambidefs.h" #include "bformatdec.h" #include "buffer_storage.h" +#include "core/filters/splitter.h" #include "effects/base.h" -#include "filters/splitter.h" +#include "effectslot.h" #include "fmt_traits.h" #include "logging.h" #include "polyphase_resampler.h" diff --git a/alc/effects/distortion.cpp b/alc/effects/distortion.cpp index 757244c5..65f8977b 100644 --- a/alc/effects/distortion.cpp +++ b/alc/effects/distortion.cpp @@ -24,11 +24,10 @@ #include #include -#include "al/auxeffectslot.h" #include "alcmain.h" #include "alcontext.h" -#include "alu.h" -#include "filters/biquad.h" +#include "core/filters/biquad.h" +#include "effectslot.h" namespace { diff --git a/alc/effects/echo.cpp b/alc/effects/echo.cpp index a50d3c61..c030ac5b 100644 --- a/alc/effects/echo.cpp +++ b/alc/effects/echo.cpp @@ -25,17 +25,19 @@ #include -#include "al/auxeffectslot.h" -#include "al/filter.h" +#include "AL/efx.h" + #include "alcmain.h" #include "alcontext.h" -#include "alu.h" -#include "filters/biquad.h" +#include "core/filters/biquad.h" +#include "effectslot.h" #include "vector.h" namespace { +constexpr float LowpassFreqRef{5000.0f}; + struct EchoState final : public EffectState { al::vector mSampleBuffer; @@ -95,7 +97,7 @@ void EchoState::update(const ALCcontext *context, const EffectSlot *slot, mTap[1].delay = float2uint(props->Echo.LRDelay*frequency + 0.5f) + mTap[0].delay; const float gainhf{maxf(1.0f - props->Echo.Damping, 0.0625f)}; /* Limit -24dB */ - mFilter.setParamsFromSlope(BiquadType::HighShelf, LOWPASSFREQREF/frequency, gainhf, 1.0f); + mFilter.setParamsFromSlope(BiquadType::HighShelf, LowpassFreqRef/frequency, gainhf, 1.0f); mFeedGain = props->Echo.Feedback; diff --git a/alc/effects/equalizer.cpp b/alc/effects/equalizer.cpp index 19d38498..c311a941 100644 --- a/alc/effects/equalizer.cpp +++ b/alc/effects/equalizer.cpp @@ -26,11 +26,10 @@ #include #include -#include "al/auxeffectslot.h" #include "alcmain.h" #include "alcontext.h" -#include "alu.h" -#include "filters/biquad.h" +#include "core/filters/biquad.h" +#include "effectslot.h" #include "vecmat.h" diff --git a/alc/effects/modulator.cpp b/alc/effects/modulator.cpp index 7e4f9fc0..a0af9890 100644 --- a/alc/effects/modulator.cpp +++ b/alc/effects/modulator.cpp @@ -26,11 +26,10 @@ #include #include -#include "al/auxeffectslot.h" #include "alcmain.h" #include "alcontext.h" -#include "alu.h" -#include "filters/biquad.h" +#include "core/filters/biquad.h" +#include "effectslot.h" #include "vecmat.h" diff --git a/alc/effects/reverb.cpp b/alc/effects/reverb.cpp index 6471b210..a4b423c7 100644 --- a/alc/effects/reverb.cpp +++ b/alc/effects/reverb.cpp @@ -29,13 +29,12 @@ #include #include -#include "al/auxeffectslot.h" -#include "al/listener.h" #include "alcmain.h" #include "alcontext.h" -#include "alu.h" +#include "alnumeric.h" #include "bformatdec.h" -#include "filters/biquad.h" +#include "core/filters/biquad.h" +#include "effectslot.h" #include "vector.h" #include "vecmat.h" diff --git a/alc/filters/biquad.cpp b/alc/filters/biquad.cpp deleted file mode 100644 index fefdc8e1..00000000 --- a/alc/filters/biquad.cpp +++ /dev/null @@ -1,163 +0,0 @@ - -#include "config.h" - -#include "biquad.h" - -#include -#include -#include - -#include "opthelpers.h" - - -template -void BiquadFilterR::setParams(BiquadType type, Real f0norm, Real gain, Real rcpQ) -{ - // Limit gain to -100dB - assert(gain > 0.00001f); - - const Real w0{al::MathDefs::Tau() * f0norm}; - const Real sin_w0{std::sin(w0)}; - const Real cos_w0{std::cos(w0)}; - const Real alpha{sin_w0/2.0f * rcpQ}; - - Real sqrtgain_alpha_2; - Real a[3]{ 1.0f, 0.0f, 0.0f }; - Real b[3]{ 1.0f, 0.0f, 0.0f }; - - /* Calculate filter coefficients depending on filter type */ - switch(type) - { - case BiquadType::HighShelf: - sqrtgain_alpha_2 = 2.0f * std::sqrt(gain) * alpha; - b[0] = gain*((gain+1.0f) + (gain-1.0f)*cos_w0 + sqrtgain_alpha_2); - b[1] = -2.0f*gain*((gain-1.0f) + (gain+1.0f)*cos_w0 ); - b[2] = gain*((gain+1.0f) + (gain-1.0f)*cos_w0 - sqrtgain_alpha_2); - a[0] = (gain+1.0f) - (gain-1.0f)*cos_w0 + sqrtgain_alpha_2; - a[1] = 2.0f* ((gain-1.0f) - (gain+1.0f)*cos_w0 ); - a[2] = (gain+1.0f) - (gain-1.0f)*cos_w0 - sqrtgain_alpha_2; - break; - case BiquadType::LowShelf: - sqrtgain_alpha_2 = 2.0f * std::sqrt(gain) * alpha; - b[0] = gain*((gain+1.0f) - (gain-1.0f)*cos_w0 + sqrtgain_alpha_2); - b[1] = 2.0f*gain*((gain-1.0f) - (gain+1.0f)*cos_w0 ); - b[2] = gain*((gain+1.0f) - (gain-1.0f)*cos_w0 - sqrtgain_alpha_2); - a[0] = (gain+1.0f) + (gain-1.0f)*cos_w0 + sqrtgain_alpha_2; - a[1] = -2.0f* ((gain-1.0f) + (gain+1.0f)*cos_w0 ); - a[2] = (gain+1.0f) + (gain-1.0f)*cos_w0 - sqrtgain_alpha_2; - break; - case BiquadType::Peaking: - b[0] = 1.0f + alpha * gain; - b[1] = -2.0f * cos_w0; - b[2] = 1.0f - alpha * gain; - a[0] = 1.0f + alpha / gain; - a[1] = -2.0f * cos_w0; - a[2] = 1.0f - alpha / gain; - break; - - case BiquadType::LowPass: - b[0] = (1.0f - cos_w0) / 2.0f; - b[1] = 1.0f - cos_w0; - b[2] = (1.0f - cos_w0) / 2.0f; - a[0] = 1.0f + alpha; - a[1] = -2.0f * cos_w0; - a[2] = 1.0f - alpha; - break; - case BiquadType::HighPass: - b[0] = (1.0f + cos_w0) / 2.0f; - b[1] = -(1.0f + cos_w0); - b[2] = (1.0f + cos_w0) / 2.0f; - a[0] = 1.0f + alpha; - a[1] = -2.0f * cos_w0; - a[2] = 1.0f - alpha; - break; - case BiquadType::BandPass: - b[0] = alpha; - b[1] = 0.0f; - b[2] = -alpha; - a[0] = 1.0f + alpha; - a[1] = -2.0f * cos_w0; - a[2] = 1.0f - alpha; - break; - } - - mA1 = a[1] / a[0]; - mA2 = a[2] / a[0]; - mB0 = b[0] / a[0]; - mB1 = b[1] / a[0]; - mB2 = b[2] / a[0]; -} - -template -void BiquadFilterR::process(const al::span src, Real *dst) -{ - const Real b0{mB0}; - const Real b1{mB1}; - const Real b2{mB2}; - const Real a1{mA1}; - const Real a2{mA2}; - Real z1{mZ1}; - Real z2{mZ2}; - - /* Processing loop is Transposed Direct Form II. This requires less storage - * compared to Direct Form I (only two delay components, instead of a four- - * sample history; the last two inputs and outputs), and works better for - * floating-point which favors summing similarly-sized values while being - * less bothered by overflow. - * - * See: http://www.earlevel.com/main/2003/02/28/biquads/ - */ - auto proc_sample = [b0,b1,b2,a1,a2,&z1,&z2](Real input) noexcept -> Real - { - const Real output{input*b0 + z1}; - z1 = input*b1 - output*a1 + z2; - z2 = input*b2 - output*a2; - return output; - }; - std::transform(src.cbegin(), src.cend(), dst, proc_sample); - - mZ1 = z1; - mZ2 = z2; -} - -template -void BiquadFilterR::dualProcess(BiquadFilterR &other, const al::span src, - Real *dst) -{ - const Real b00{mB0}; - const Real b01{mB1}; - const Real b02{mB2}; - const Real a01{mA1}; - const Real a02{mA2}; - const Real b10{other.mB0}; - const Real b11{other.mB1}; - const Real b12{other.mB2}; - const Real a11{other.mA1}; - const Real a12{other.mA2}; - Real z01{mZ1}; - Real z02{mZ2}; - Real z11{other.mZ1}; - Real z12{other.mZ2}; - - auto proc_sample = [b00,b01,b02,a01,a02,b10,b11,b12,a11,a12,&z01,&z02,&z11,&z12](Real input) noexcept -> Real - { - const Real tmpout{input*b00 + z01}; - z01 = input*b01 - tmpout*a01 + z02; - z02 = input*b02 - tmpout*a02; - input = tmpout; - - const Real output{input*b10 + z11}; - z11 = input*b11 - output*a11 + z12; - z12 = input*b12 - output*a12; - return output; - }; - std::transform(src.cbegin(), src.cend(), dst, proc_sample); - - mZ1 = z01; - mZ2 = z02; - other.mZ1 = z11; - other.mZ2 = z12; -} - -template class BiquadFilterR; -template class BiquadFilterR; diff --git a/alc/filters/biquad.h b/alc/filters/biquad.h deleted file mode 100644 index 3ce70cb3..00000000 --- a/alc/filters/biquad.h +++ /dev/null @@ -1,144 +0,0 @@ -#ifndef FILTERS_BIQUAD_H -#define FILTERS_BIQUAD_H - -#include -#include -#include -#include - -#include "alspan.h" -#include "math_defs.h" - - -/* Filters implementation is based on the "Cookbook formulae for audio - * EQ biquad filter coefficients" by Robert Bristow-Johnson - * http://www.musicdsp.org/files/Audio-EQ-Cookbook.txt - */ -/* Implementation note: For the shelf and peaking filters, the specified gain - * is for the centerpoint of the transition band. This better fits EFX filter - * behavior, which expects the shelf's reference frequency to reach the given - * gain. To set the gain for the shelf or peak itself, use the square root of - * the desired linear gain (or halve the dB gain). - */ - -enum class BiquadType { - /** EFX-style low-pass filter, specifying a gain and reference frequency. */ - HighShelf, - /** EFX-style high-pass filter, specifying a gain and reference frequency. */ - LowShelf, - /** Peaking filter, specifying a gain and reference frequency. */ - Peaking, - - /** Low-pass cut-off filter, specifying a cut-off frequency. */ - LowPass, - /** High-pass cut-off filter, specifying a cut-off frequency. */ - HighPass, - /** Band-pass filter, specifying a center frequency. */ - BandPass, -}; - -template -class BiquadFilterR { - /* Last two delayed components for direct form II. */ - Real mZ1{0.0f}, mZ2{0.0f}; - /* Transfer function coefficients "b" (numerator) */ - Real mB0{1.0f}, mB1{0.0f}, mB2{0.0f}; - /* Transfer function coefficients "a" (denominator; a0 is pre-applied). */ - Real mA1{0.0f}, mA2{0.0f}; - - void setParams(BiquadType type, Real f0norm, Real gain, Real rcpQ); - - /** - * Calculates the rcpQ (i.e. 1/Q) coefficient for shelving filters, using - * the reference gain and shelf slope parameter. - * \param gain 0 < gain - * \param slope 0 < slope <= 1 - */ - static Real rcpQFromSlope(Real gain, Real slope) - { return std::sqrt((gain + 1.0f/gain)*(1.0f/slope - 1.0f) + 2.0f); } - - /** - * Calculates the rcpQ (i.e. 1/Q) coefficient for filters, using the - * normalized reference frequency and bandwidth. - * \param f0norm 0 < f0norm < 0.5. - * \param bandwidth 0 < bandwidth - */ - static Real rcpQFromBandwidth(Real f0norm, Real bandwidth) - { - const Real w0{al::MathDefs::Tau() * f0norm}; - return 2.0f*std::sinh(std::log(Real{2.0f})/2.0f*bandwidth*w0/std::sin(w0)); - } - -public: - void clear() noexcept { mZ1 = mZ2 = 0.0f; } - - /** - * Sets the filter state for the specified filter type and its parameters. - * - * \param type The type of filter to apply. - * \param f0norm The normalized reference frequency (ref / sample_rate). - * This is the center point for the Shelf, Peaking, and BandPass filter - * types, or the cutoff frequency for the LowPass and HighPass filter - * types. - * \param gain The gain for the reference frequency response. Only used by - * the Shelf and Peaking filter types. - * \param slope Slope steepness of the transition band. - */ - void setParamsFromSlope(BiquadType type, Real f0norm, Real gain, Real slope) - { - gain = std::max(gain, 0.001f); /* Limit -60dB */ - setParams(type, f0norm, gain, rcpQFromSlope(gain, slope)); - } - - /** - * Sets the filter state for the specified filter type and its parameters. - * - * \param type The type of filter to apply. - * \param f0norm The normalized reference frequency (ref / sample_rate). - * This is the center point for the Shelf, Peaking, and BandPass filter - * types, or the cutoff frequency for the LowPass and HighPass filter - * types. - * \param gain The gain for the reference frequency response. Only used by - * the Shelf and Peaking filter types. - * \param bandwidth Normalized bandwidth of the transition band. - */ - void setParamsFromBandwidth(BiquadType type, Real f0norm, Real gain, Real bandwidth) - { setParams(type, f0norm, gain, rcpQFromBandwidth(f0norm, bandwidth)); } - - void copyParamsFrom(const BiquadFilterR &other) - { - mB0 = other.mB0; - mB1 = other.mB1; - mB2 = other.mB2; - mA1 = other.mA1; - mA2 = other.mA2; - } - - void process(const al::span src, Real *dst); - /** Processes this filter and the other at the same time. */ - void dualProcess(BiquadFilterR &other, const al::span src, Real *dst); - - /* Rather hacky. It's just here to support "manual" processing. */ - std::pair getComponents() const noexcept { return {mZ1, mZ2}; } - void setComponents(Real z1, Real z2) noexcept { mZ1 = z1; mZ2 = z2; } - Real processOne(const Real in, Real &z1, Real &z2) const noexcept - { - const Real out{in*mB0 + z1}; - z1 = in*mB1 - out*mA1 + z2; - z2 = in*mB2 - out*mA2; - return out; - } -}; - -template -struct DualBiquadR { - BiquadFilterR &f0, &f1; - - void process(const al::span src, Real *dst) - { f0.dualProcess(f1, src, dst); } -}; - -using BiquadFilter = BiquadFilterR; -using DualBiquad = DualBiquadR; - -#endif /* FILTERS_BIQUAD_H */ diff --git a/alc/filters/nfc.cpp b/alc/filters/nfc.cpp deleted file mode 100644 index 9a28517c..00000000 --- a/alc/filters/nfc.cpp +++ /dev/null @@ -1,383 +0,0 @@ - -#include "config.h" - -#include "nfc.h" - -#include - -#include "opthelpers.h" - - -/* Near-field control filters are the basis for handling the near-field effect. - * The near-field effect is a bass-boost present in the directional components - * of a recorded signal, created as a result of the wavefront curvature (itself - * a function of sound distance). Proper reproduction dictates this be - * compensated for using a bass-cut given the playback speaker distance, to - * avoid excessive bass in the playback. - * - * For real-time rendered audio, emulating the near-field effect based on the - * sound source's distance, and subsequently compensating for it at output - * based on the speaker distances, can create a more realistic perception of - * sound distance beyond a simple 1/r attenuation. - * - * These filters do just that. Each one applies a low-shelf filter, created as - * the combination of a bass-boost for a given sound source distance (near- - * field emulation) along with a bass-cut for a given control/speaker distance - * (near-field compensation). - * - * Note that it is necessary to apply a cut along with the boost, since the - * boost alone is unstable in higher-order ambisonics as it causes an infinite - * DC gain (even first-order ambisonics requires there to be no DC offset for - * the boost to work). Consequently, ambisonics requires a control parameter to - * be used to avoid an unstable boost-only filter. NFC-HOA defines this control - * as a reference delay, calculated with: - * - * reference_delay = control_distance / speed_of_sound - * - * This means w0 (for input) or w1 (for output) should be set to: - * - * wN = 1 / (reference_delay * sample_rate) - * - * when dealing with NFC-HOA content. For FOA input content, which does not - * specify a reference_delay variable, w0 should be set to 0 to apply only - * near-field compensation for output. It's important that w1 be a finite, - * positive, non-0 value or else the bass-boost will become unstable again. - * Also, w0 should not be too large compared to w1, to avoid excessively loud - * low frequencies. - */ - -namespace { - -constexpr float B[5][4] = { - { 0.0f }, - { 1.0f }, - { 3.0f, 3.0f }, - { 3.6778f, 6.4595f, 2.3222f }, - { 4.2076f, 11.4877f, 5.7924f, 9.1401f } -}; - -NfcFilter1 NfcFilterCreate1(const float w0, const float w1) noexcept -{ - NfcFilter1 nfc{}; - float b_00, g_0; - float r; - - nfc.base_gain = 1.0f; - nfc.gain = 1.0f; - - /* Calculate bass-boost coefficients. */ - r = 0.5f * w0; - b_00 = B[1][0] * r; - g_0 = 1.0f + b_00; - - nfc.gain *= g_0; - nfc.b1 = 2.0f * b_00 / g_0; - - /* Calculate bass-cut coefficients. */ - r = 0.5f * w1; - b_00 = B[1][0] * r; - g_0 = 1.0f + b_00; - - nfc.base_gain /= g_0; - nfc.gain /= g_0; - nfc.a1 = 2.0f * b_00 / g_0; - - return nfc; -} - -void NfcFilterAdjust1(NfcFilter1 *nfc, const float w0) noexcept -{ - const float r{0.5f * w0}; - const float b_00{B[1][0] * r}; - const float g_0{1.0f + b_00}; - - nfc->gain = nfc->base_gain * g_0; - nfc->b1 = 2.0f * b_00 / g_0; -} - - -NfcFilter2 NfcFilterCreate2(const float w0, const float w1) noexcept -{ - NfcFilter2 nfc{}; - float b_10, b_11, g_1; - float r; - - nfc.base_gain = 1.0f; - nfc.gain = 1.0f; - - /* Calculate bass-boost coefficients. */ - r = 0.5f * w0; - b_10 = B[2][0] * r; - b_11 = B[2][1] * r * r; - g_1 = 1.0f + b_10 + b_11; - - nfc.gain *= g_1; - nfc.b1 = (2.0f*b_10 + 4.0f*b_11) / g_1; - nfc.b2 = 4.0f * b_11 / g_1; - - /* Calculate bass-cut coefficients. */ - r = 0.5f * w1; - b_10 = B[2][0] * r; - b_11 = B[2][1] * r * r; - g_1 = 1.0f + b_10 + b_11; - - nfc.base_gain /= g_1; - nfc.gain /= g_1; - nfc.a1 = (2.0f*b_10 + 4.0f*b_11) / g_1; - nfc.a2 = 4.0f * b_11 / g_1; - - return nfc; -} - -void NfcFilterAdjust2(NfcFilter2 *nfc, const float w0) noexcept -{ - const float r{0.5f * w0}; - const float b_10{B[2][0] * r}; - const float b_11{B[2][1] * r * r}; - const float g_1{1.0f + b_10 + b_11}; - - nfc->gain = nfc->base_gain * g_1; - nfc->b1 = (2.0f*b_10 + 4.0f*b_11) / g_1; - nfc->b2 = 4.0f * b_11 / g_1; -} - - -NfcFilter3 NfcFilterCreate3(const float w0, const float w1) noexcept -{ - NfcFilter3 nfc{}; - float b_10, b_11, g_1; - float b_00, g_0; - float r; - - nfc.base_gain = 1.0f; - nfc.gain = 1.0f; - - /* Calculate bass-boost coefficients. */ - r = 0.5f * w0; - b_10 = B[3][0] * r; - b_11 = B[3][1] * r * r; - b_00 = B[3][2] * r; - g_1 = 1.0f + b_10 + b_11; - g_0 = 1.0f + b_00; - - nfc.gain *= g_1 * g_0; - nfc.b1 = (2.0f*b_10 + 4.0f*b_11) / g_1; - nfc.b2 = 4.0f * b_11 / g_1; - nfc.b3 = 2.0f * b_00 / g_0; - - /* Calculate bass-cut coefficients. */ - r = 0.5f * w1; - b_10 = B[3][0] * r; - b_11 = B[3][1] * r * r; - b_00 = B[3][2] * r; - g_1 = 1.0f + b_10 + b_11; - g_0 = 1.0f + b_00; - - nfc.base_gain /= g_1 * g_0; - nfc.gain /= g_1 * g_0; - nfc.a1 = (2.0f*b_10 + 4.0f*b_11) / g_1; - nfc.a2 = 4.0f * b_11 / g_1; - nfc.a3 = 2.0f * b_00 / g_0; - - return nfc; -} - -void NfcFilterAdjust3(NfcFilter3 *nfc, const float w0) noexcept -{ - const float r{0.5f * w0}; - const float b_10{B[3][0] * r}; - const float b_11{B[3][1] * r * r}; - const float b_00{B[3][2] * r}; - const float g_1{1.0f + b_10 + b_11}; - const float g_0{1.0f + b_00}; - - nfc->gain = nfc->base_gain * g_1 * g_0; - nfc->b1 = (2.0f*b_10 + 4.0f*b_11) / g_1; - nfc->b2 = 4.0f * b_11 / g_1; - nfc->b3 = 2.0f * b_00 / g_0; -} - - -NfcFilter4 NfcFilterCreate4(const float w0, const float w1) noexcept -{ - NfcFilter4 nfc{}; - float b_10, b_11, g_1; - float b_00, b_01, g_0; - float r; - - nfc.base_gain = 1.0f; - nfc.gain = 1.0f; - - /* Calculate bass-boost coefficients. */ - r = 0.5f * w0; - b_10 = B[4][0] * r; - b_11 = B[4][1] * r * r; - b_00 = B[4][2] * r; - b_01 = B[4][3] * r * r; - g_1 = 1.0f + b_10 + b_11; - g_0 = 1.0f + b_00 + b_01; - - nfc.gain *= g_1 * g_0; - nfc.b1 = (2.0f*b_10 + 4.0f*b_11) / g_1; - nfc.b2 = 4.0f * b_11 / g_1; - nfc.b3 = (2.0f*b_00 + 4.0f*b_01) / g_0; - nfc.b4 = 4.0f * b_01 / g_0; - - /* Calculate bass-cut coefficients. */ - r = 0.5f * w1; - b_10 = B[4][0] * r; - b_11 = B[4][1] * r * r; - b_00 = B[4][2] * r; - b_01 = B[4][3] * r * r; - g_1 = 1.0f + b_10 + b_11; - g_0 = 1.0f + b_00 + b_01; - - nfc.base_gain /= g_1 * g_0; - nfc.gain /= g_1 * g_0; - nfc.a1 = (2.0f*b_10 + 4.0f*b_11) / g_1; - nfc.a2 = 4.0f * b_11 / g_1; - nfc.a3 = (2.0f*b_00 + 4.0f*b_01) / g_0; - nfc.a4 = 4.0f * b_01 / g_0; - - return nfc; -} - -void NfcFilterAdjust4(NfcFilter4 *nfc, const float w0) noexcept -{ - const float r{0.5f * w0}; - const float b_10{B[4][0] * r}; - const float b_11{B[4][1] * r * r}; - const float b_00{B[4][2] * r}; - const float b_01{B[4][3] * r * r}; - const float g_1{1.0f + b_10 + b_11}; - const float g_0{1.0f + b_00 + b_01}; - - nfc->gain = nfc->base_gain * g_1 * g_0; - nfc->b1 = (2.0f*b_10 + 4.0f*b_11) / g_1; - nfc->b2 = 4.0f * b_11 / g_1; - nfc->b3 = (2.0f*b_00 + 4.0f*b_01) / g_0; - nfc->b4 = 4.0f * b_01 / g_0; -} - -} // namespace - -void NfcFilter::init(const float w1) noexcept -{ - first = NfcFilterCreate1(0.0f, w1); - second = NfcFilterCreate2(0.0f, w1); - third = NfcFilterCreate3(0.0f, w1); - fourth = NfcFilterCreate4(0.0f, w1); -} - -void NfcFilter::adjust(const float w0) noexcept -{ - NfcFilterAdjust1(&first, w0); - NfcFilterAdjust2(&second, w0); - NfcFilterAdjust3(&third, w0); - NfcFilterAdjust4(&fourth, w0); -} - - -void NfcFilter::process1(const al::span src, float *RESTRICT dst) -{ - const float gain{first.gain}; - const float b1{first.b1}; - const float a1{first.a1}; - float z1{first.z[0]}; - auto proc_sample = [gain,b1,a1,&z1](const float in) noexcept -> float - { - const float y{in*gain - a1*z1}; - const float out{y + b1*z1}; - z1 += y; - return out; - }; - std::transform(src.cbegin(), src.cend(), dst, proc_sample); - first.z[0] = z1; -} - -void NfcFilter::process2(const al::span src, float *RESTRICT dst) -{ - const float gain{second.gain}; - const float b1{second.b1}; - const float b2{second.b2}; - const float a1{second.a1}; - const float a2{second.a2}; - float z1{second.z[0]}; - float z2{second.z[1]}; - auto proc_sample = [gain,b1,b2,a1,a2,&z1,&z2](const float in) noexcept -> float - { - const float y{in*gain - a1*z1 - a2*z2}; - const float out{y + b1*z1 + b2*z2}; - z2 += z1; - z1 += y; - return out; - }; - std::transform(src.cbegin(), src.cend(), dst, proc_sample); - second.z[0] = z1; - second.z[1] = z2; -} - -void NfcFilter::process3(const al::span src, float *RESTRICT dst) -{ - const float gain{third.gain}; - const float b1{third.b1}; - const float b2{third.b2}; - const float b3{third.b3}; - const float a1{third.a1}; - const float a2{third.a2}; - const float a3{third.a3}; - float z1{third.z[0]}; - float z2{third.z[1]}; - float z3{third.z[2]}; - auto proc_sample = [gain,b1,b2,b3,a1,a2,a3,&z1,&z2,&z3](const float in) noexcept -> float - { - float y{in*gain - a1*z1 - a2*z2}; - float out{y + b1*z1 + b2*z2}; - z2 += z1; - z1 += y; - - y = out - a3*z3; - out = y + b3*z3; - z3 += y; - return out; - }; - std::transform(src.cbegin(), src.cend(), dst, proc_sample); - third.z[0] = z1; - third.z[1] = z2; - third.z[2] = z3; -} - -void NfcFilter::process4(const al::span src, float *RESTRICT dst) -{ - const float gain{fourth.gain}; - const float b1{fourth.b1}; - const float b2{fourth.b2}; - const float b3{fourth.b3}; - const float b4{fourth.b4}; - const float a1{fourth.a1}; - const float a2{fourth.a2}; - const float a3{fourth.a3}; - const float a4{fourth.a4}; - float z1{fourth.z[0]}; - float z2{fourth.z[1]}; - float z3{fourth.z[2]}; - float z4{fourth.z[3]}; - auto proc_sample = [gain,b1,b2,b3,b4,a1,a2,a3,a4,&z1,&z2,&z3,&z4](const float in) noexcept -> float - { - float y{in*gain - a1*z1 - a2*z2}; - float out{y + b1*z1 + b2*z2}; - z2 += z1; - z1 += y; - - y = out - a3*z3 - a4*z4; - out = y + b3*z3 + b4*z4; - z4 += z3; - z3 += y; - return out; - }; - std::transform(src.cbegin(), src.cend(), dst, proc_sample); - fourth.z[0] = z1; - fourth.z[1] = z2; - fourth.z[2] = z3; - fourth.z[3] = z4; -} diff --git a/alc/filters/nfc.h b/alc/filters/nfc.h deleted file mode 100644 index 2327c966..00000000 --- a/alc/filters/nfc.h +++ /dev/null @@ -1,63 +0,0 @@ -#ifndef FILTER_NFC_H -#define FILTER_NFC_H - -#include - -#include "alspan.h" - - -struct NfcFilter1 { - float base_gain, gain; - float b1, a1; - float z[1]; -}; -struct NfcFilter2 { - float base_gain, gain; - float b1, b2, a1, a2; - float z[2]; -}; -struct NfcFilter3 { - float base_gain, gain; - float b1, b2, b3, a1, a2, a3; - float z[3]; -}; -struct NfcFilter4 { - float base_gain, gain; - float b1, b2, b3, b4, a1, a2, a3, a4; - float z[4]; -}; - -class NfcFilter { - NfcFilter1 first; - NfcFilter2 second; - NfcFilter3 third; - NfcFilter4 fourth; - -public: - /* NOTE: - * w0 = speed_of_sound / (source_distance * sample_rate); - * w1 = speed_of_sound / (control_distance * sample_rate); - * - * Generally speaking, the control distance should be approximately the - * average speaker distance, or based on the reference delay if outputing - * NFC-HOA. It must not be negative, 0, or infinite. The source distance - * should not be too small relative to the control distance. - */ - - void init(const float w1) noexcept; - void adjust(const float w0) noexcept; - - /* Near-field control filter for first-order ambisonic channels (1-3). */ - void process1(const al::span src, float *RESTRICT dst); - - /* Near-field control filter for second-order ambisonic channels (4-8). */ - void process2(const al::span src, float *RESTRICT dst); - - /* Near-field control filter for third-order ambisonic channels (9-15). */ - void process3(const al::span src, float *RESTRICT dst); - - /* Near-field control filter for fourth-order ambisonic channels (16-24). */ - void process4(const al::span src, float *RESTRICT dst); -}; - -#endif /* FILTER_NFC_H */ diff --git a/alc/filters/splitter.cpp b/alc/filters/splitter.cpp deleted file mode 100644 index 5cc670b7..00000000 --- a/alc/filters/splitter.cpp +++ /dev/null @@ -1,113 +0,0 @@ - -#include "config.h" - -#include "splitter.h" - -#include -#include -#include - -#include "math_defs.h" -#include "opthelpers.h" - - -template -void BandSplitterR::init(Real f0norm) -{ - const Real w{f0norm * al::MathDefs::Tau()}; - const Real cw{std::cos(w)}; - if(cw > std::numeric_limits::epsilon()) - mCoeff = (std::sin(w) - 1.0f) / cw; - else - mCoeff = cw * -0.5f; - - mLpZ1 = 0.0f; - mLpZ2 = 0.0f; - mApZ1 = 0.0f; -} - -template -void BandSplitterR::process(const al::span input, Real *hpout, Real *lpout) -{ - const Real ap_coeff{mCoeff}; - const Real lp_coeff{mCoeff*0.5f + 0.5f}; - Real lp_z1{mLpZ1}; - Real lp_z2{mLpZ2}; - Real ap_z1{mApZ1}; - auto proc_sample = [ap_coeff,lp_coeff,&lp_z1,&lp_z2,&ap_z1,&lpout](const Real in) noexcept -> Real - { - /* Low-pass sample processing. */ - Real d{(in - lp_z1) * lp_coeff}; - Real lp_y{lp_z1 + d}; - lp_z1 = lp_y + d; - - d = (lp_y - lp_z2) * lp_coeff; - lp_y = lp_z2 + d; - lp_z2 = lp_y + d; - - *(lpout++) = lp_y; - - /* All-pass sample processing. */ - Real ap_y{in*ap_coeff + ap_z1}; - ap_z1 = in - ap_y*ap_coeff; - - /* High-pass generated from removing low-passed output. */ - return ap_y - lp_y; - }; - std::transform(input.cbegin(), input.cend(), hpout, proc_sample); - mLpZ1 = lp_z1; - mLpZ2 = lp_z2; - mApZ1 = ap_z1; -} - -template -void BandSplitterR::processHfScale(const al::span samples, const Real hfscale) -{ - const Real ap_coeff{mCoeff}; - const Real lp_coeff{mCoeff*0.5f + 0.5f}; - Real lp_z1{mLpZ1}; - Real lp_z2{mLpZ2}; - Real ap_z1{mApZ1}; - auto proc_sample = [hfscale,ap_coeff,lp_coeff,&lp_z1,&lp_z2,&ap_z1](const Real in) noexcept -> Real - { - /* Low-pass sample processing. */ - Real d{(in - lp_z1) * lp_coeff}; - Real lp_y{lp_z1 + d}; - lp_z1 = lp_y + d; - - d = (lp_y - lp_z2) * lp_coeff; - lp_y = lp_z2 + d; - lp_z2 = lp_y + d; - - /* All-pass sample processing. */ - Real ap_y{in*ap_coeff + ap_z1}; - ap_z1 = in - ap_y*ap_coeff; - - /* High-pass generated by removing the low-passed signal, which is then - * scaled and added back to the low-passed signal. - */ - return (ap_y-lp_y)*hfscale + lp_y; - }; - std::transform(samples.begin(), samples.end(), samples.begin(), proc_sample); - mLpZ1 = lp_z1; - mLpZ2 = lp_z2; - mApZ1 = ap_z1; -} - -template -void BandSplitterR::applyAllpass(const al::span samples) const -{ - const Real coeff{mCoeff}; - Real z1{0.0f}; - auto proc_sample = [coeff,&z1](const Real in) noexcept -> Real - { - const Real out{in*coeff + z1}; - z1 = in - out*coeff; - return out; - }; - std::transform(samples.begin(), samples.end(), samples.begin(), proc_sample); -} - - -template class BandSplitterR; -template class BandSplitterR; diff --git a/alc/filters/splitter.h b/alc/filters/splitter.h deleted file mode 100644 index 18ab4998..00000000 --- a/alc/filters/splitter.h +++ /dev/null @@ -1,36 +0,0 @@ -#ifndef FILTER_SPLITTER_H -#define FILTER_SPLITTER_H - -#include - -#include "alspan.h" - - -/* Band splitter. Splits a signal into two phase-matching frequency bands. */ -template -class BandSplitterR { - Real mCoeff{0.0f}; - Real mLpZ1{0.0f}; - Real mLpZ2{0.0f}; - Real mApZ1{0.0f}; - -public: - BandSplitterR() = default; - BandSplitterR(const BandSplitterR&) = default; - BandSplitterR(Real f0norm) { init(f0norm); } - - void init(Real f0norm); - void clear() noexcept { mLpZ1 = mLpZ2 = mApZ1 = 0.0f; } - void process(const al::span input, Real *hpout, Real *lpout); - - void processHfScale(const al::span samples, const Real hfscale); - - /* The all-pass portion of the band splitter. Applies the same phase shift - * without splitting the signal. Note that each use of this method is - * indepedent, it does not track history between calls. - */ - void applyAllpass(const al::span samples) const; -}; -using BandSplitter = BandSplitterR; - -#endif /* FILTER_SPLITTER_H */ diff --git a/alc/front_stablizer.h b/alc/front_stablizer.h index 94882835..0fedeb50 100644 --- a/alc/front_stablizer.h +++ b/alc/front_stablizer.h @@ -6,7 +6,7 @@ #include "almalloc.h" #include "core/bufferline.h" -#include "filters/splitter.h" +#include "core/filters/splitter.h" struct FrontStablizer { diff --git a/alc/hrtf.cpp b/alc/hrtf.cpp index 12cea416..01dc342f 100644 --- a/alc/hrtf.cpp +++ b/alc/hrtf.cpp @@ -47,7 +47,7 @@ #include "alnumeric.h" #include "aloptional.h" #include "alspan.h" -#include "filters/splitter.h" +#include "core/filters/splitter.h" #include "logging.h" #include "math_defs.h" #include "opthelpers.h" diff --git a/alc/hrtf.h b/alc/hrtf.h index 3b9a272c..c26947de 100644 --- a/alc/hrtf.h +++ b/alc/hrtf.h @@ -11,7 +11,7 @@ #include "ambidefs.h" #include "atomic.h" #include "core/bufferline.h" -#include "filters/splitter.h" +#include "core/filters/splitter.h" #include "intrusive_ptr.h" #include "vector.h" diff --git a/alc/voice.cpp b/alc/voice.cpp index 54775a72..8b5cd02f 100644 --- a/alc/voice.cpp +++ b/alc/voice.cpp @@ -47,10 +47,10 @@ #include "alstring.h" #include "alu.h" #include "core/devformat.h" +#include "core/filters/biquad.h" +#include "core/filters/nfc.h" +#include "core/filters/splitter.h" #include "cpu_caps.h" -#include "filters/biquad.h" -#include "filters/nfc.h" -#include "filters/splitter.h" #include "fmt_traits.h" #include "hrtf.h" #include "inprogext.h" diff --git a/alc/voice.h b/alc/voice.h index 8a7d12d1..00012de3 100644 --- a/alc/voice.h +++ b/alc/voice.h @@ -8,9 +8,9 @@ #include "alu.h" #include "buffer_storage.h" #include "core/devformat.h" -#include "filters/biquad.h" -#include "filters/nfc.h" -#include "filters/splitter.h" +#include "core/filters/biquad.h" +#include "core/filters/nfc.h" +#include "core/filters/splitter.h" #include "hrtf.h" #include "mixer/defs.h" diff --git a/core/filters/biquad.cpp b/core/filters/biquad.cpp new file mode 100644 index 00000000..fefdc8e1 --- /dev/null +++ b/core/filters/biquad.cpp @@ -0,0 +1,163 @@ + +#include "config.h" + +#include "biquad.h" + +#include +#include +#include + +#include "opthelpers.h" + + +template +void BiquadFilterR::setParams(BiquadType type, Real f0norm, Real gain, Real rcpQ) +{ + // Limit gain to -100dB + assert(gain > 0.00001f); + + const Real w0{al::MathDefs::Tau() * f0norm}; + const Real sin_w0{std::sin(w0)}; + const Real cos_w0{std::cos(w0)}; + const Real alpha{sin_w0/2.0f * rcpQ}; + + Real sqrtgain_alpha_2; + Real a[3]{ 1.0f, 0.0f, 0.0f }; + Real b[3]{ 1.0f, 0.0f, 0.0f }; + + /* Calculate filter coefficients depending on filter type */ + switch(type) + { + case BiquadType::HighShelf: + sqrtgain_alpha_2 = 2.0f * std::sqrt(gain) * alpha; + b[0] = gain*((gain+1.0f) + (gain-1.0f)*cos_w0 + sqrtgain_alpha_2); + b[1] = -2.0f*gain*((gain-1.0f) + (gain+1.0f)*cos_w0 ); + b[2] = gain*((gain+1.0f) + (gain-1.0f)*cos_w0 - sqrtgain_alpha_2); + a[0] = (gain+1.0f) - (gain-1.0f)*cos_w0 + sqrtgain_alpha_2; + a[1] = 2.0f* ((gain-1.0f) - (gain+1.0f)*cos_w0 ); + a[2] = (gain+1.0f) - (gain-1.0f)*cos_w0 - sqrtgain_alpha_2; + break; + case BiquadType::LowShelf: + sqrtgain_alpha_2 = 2.0f * std::sqrt(gain) * alpha; + b[0] = gain*((gain+1.0f) - (gain-1.0f)*cos_w0 + sqrtgain_alpha_2); + b[1] = 2.0f*gain*((gain-1.0f) - (gain+1.0f)*cos_w0 ); + b[2] = gain*((gain+1.0f) - (gain-1.0f)*cos_w0 - sqrtgain_alpha_2); + a[0] = (gain+1.0f) + (gain-1.0f)*cos_w0 + sqrtgain_alpha_2; + a[1] = -2.0f* ((gain-1.0f) + (gain+1.0f)*cos_w0 ); + a[2] = (gain+1.0f) + (gain-1.0f)*cos_w0 - sqrtgain_alpha_2; + break; + case BiquadType::Peaking: + b[0] = 1.0f + alpha * gain; + b[1] = -2.0f * cos_w0; + b[2] = 1.0f - alpha * gain; + a[0] = 1.0f + alpha / gain; + a[1] = -2.0f * cos_w0; + a[2] = 1.0f - alpha / gain; + break; + + case BiquadType::LowPass: + b[0] = (1.0f - cos_w0) / 2.0f; + b[1] = 1.0f - cos_w0; + b[2] = (1.0f - cos_w0) / 2.0f; + a[0] = 1.0f + alpha; + a[1] = -2.0f * cos_w0; + a[2] = 1.0f - alpha; + break; + case BiquadType::HighPass: + b[0] = (1.0f + cos_w0) / 2.0f; + b[1] = -(1.0f + cos_w0); + b[2] = (1.0f + cos_w0) / 2.0f; + a[0] = 1.0f + alpha; + a[1] = -2.0f * cos_w0; + a[2] = 1.0f - alpha; + break; + case BiquadType::BandPass: + b[0] = alpha; + b[1] = 0.0f; + b[2] = -alpha; + a[0] = 1.0f + alpha; + a[1] = -2.0f * cos_w0; + a[2] = 1.0f - alpha; + break; + } + + mA1 = a[1] / a[0]; + mA2 = a[2] / a[0]; + mB0 = b[0] / a[0]; + mB1 = b[1] / a[0]; + mB2 = b[2] / a[0]; +} + +template +void BiquadFilterR::process(const al::span src, Real *dst) +{ + const Real b0{mB0}; + const Real b1{mB1}; + const Real b2{mB2}; + const Real a1{mA1}; + const Real a2{mA2}; + Real z1{mZ1}; + Real z2{mZ2}; + + /* Processing loop is Transposed Direct Form II. This requires less storage + * compared to Direct Form I (only two delay components, instead of a four- + * sample history; the last two inputs and outputs), and works better for + * floating-point which favors summing similarly-sized values while being + * less bothered by overflow. + * + * See: http://www.earlevel.com/main/2003/02/28/biquads/ + */ + auto proc_sample = [b0,b1,b2,a1,a2,&z1,&z2](Real input) noexcept -> Real + { + const Real output{input*b0 + z1}; + z1 = input*b1 - output*a1 + z2; + z2 = input*b2 - output*a2; + return output; + }; + std::transform(src.cbegin(), src.cend(), dst, proc_sample); + + mZ1 = z1; + mZ2 = z2; +} + +template +void BiquadFilterR::dualProcess(BiquadFilterR &other, const al::span src, + Real *dst) +{ + const Real b00{mB0}; + const Real b01{mB1}; + const Real b02{mB2}; + const Real a01{mA1}; + const Real a02{mA2}; + const Real b10{other.mB0}; + const Real b11{other.mB1}; + const Real b12{other.mB2}; + const Real a11{other.mA1}; + const Real a12{other.mA2}; + Real z01{mZ1}; + Real z02{mZ2}; + Real z11{other.mZ1}; + Real z12{other.mZ2}; + + auto proc_sample = [b00,b01,b02,a01,a02,b10,b11,b12,a11,a12,&z01,&z02,&z11,&z12](Real input) noexcept -> Real + { + const Real tmpout{input*b00 + z01}; + z01 = input*b01 - tmpout*a01 + z02; + z02 = input*b02 - tmpout*a02; + input = tmpout; + + const Real output{input*b10 + z11}; + z11 = input*b11 - output*a11 + z12; + z12 = input*b12 - output*a12; + return output; + }; + std::transform(src.cbegin(), src.cend(), dst, proc_sample); + + mZ1 = z01; + mZ2 = z02; + other.mZ1 = z11; + other.mZ2 = z12; +} + +template class BiquadFilterR; +template class BiquadFilterR; diff --git a/core/filters/biquad.h b/core/filters/biquad.h new file mode 100644 index 00000000..b2e2cfdb --- /dev/null +++ b/core/filters/biquad.h @@ -0,0 +1,144 @@ +#ifndef CORE_FILTERS_BIQUAD_H +#define CORE_FILTERS_BIQUAD_H + +#include +#include +#include +#include + +#include "alspan.h" +#include "math_defs.h" + + +/* Filters implementation is based on the "Cookbook formulae for audio + * EQ biquad filter coefficients" by Robert Bristow-Johnson + * http://www.musicdsp.org/files/Audio-EQ-Cookbook.txt + */ +/* Implementation note: For the shelf and peaking filters, the specified gain + * is for the centerpoint of the transition band. This better fits EFX filter + * behavior, which expects the shelf's reference frequency to reach the given + * gain. To set the gain for the shelf or peak itself, use the square root of + * the desired linear gain (or halve the dB gain). + */ + +enum class BiquadType { + /** EFX-style low-pass filter, specifying a gain and reference frequency. */ + HighShelf, + /** EFX-style high-pass filter, specifying a gain and reference frequency. */ + LowShelf, + /** Peaking filter, specifying a gain and reference frequency. */ + Peaking, + + /** Low-pass cut-off filter, specifying a cut-off frequency. */ + LowPass, + /** High-pass cut-off filter, specifying a cut-off frequency. */ + HighPass, + /** Band-pass filter, specifying a center frequency. */ + BandPass, +}; + +template +class BiquadFilterR { + /* Last two delayed components for direct form II. */ + Real mZ1{0.0f}, mZ2{0.0f}; + /* Transfer function coefficients "b" (numerator) */ + Real mB0{1.0f}, mB1{0.0f}, mB2{0.0f}; + /* Transfer function coefficients "a" (denominator; a0 is pre-applied). */ + Real mA1{0.0f}, mA2{0.0f}; + + void setParams(BiquadType type, Real f0norm, Real gain, Real rcpQ); + + /** + * Calculates the rcpQ (i.e. 1/Q) coefficient for shelving filters, using + * the reference gain and shelf slope parameter. + * \param gain 0 < gain + * \param slope 0 < slope <= 1 + */ + static Real rcpQFromSlope(Real gain, Real slope) + { return std::sqrt((gain + 1.0f/gain)*(1.0f/slope - 1.0f) + 2.0f); } + + /** + * Calculates the rcpQ (i.e. 1/Q) coefficient for filters, using the + * normalized reference frequency and bandwidth. + * \param f0norm 0 < f0norm < 0.5. + * \param bandwidth 0 < bandwidth + */ + static Real rcpQFromBandwidth(Real f0norm, Real bandwidth) + { + const Real w0{al::MathDefs::Tau() * f0norm}; + return 2.0f*std::sinh(std::log(Real{2.0f})/2.0f*bandwidth*w0/std::sin(w0)); + } + +public: + void clear() noexcept { mZ1 = mZ2 = 0.0f; } + + /** + * Sets the filter state for the specified filter type and its parameters. + * + * \param type The type of filter to apply. + * \param f0norm The normalized reference frequency (ref / sample_rate). + * This is the center point for the Shelf, Peaking, and BandPass filter + * types, or the cutoff frequency for the LowPass and HighPass filter + * types. + * \param gain The gain for the reference frequency response. Only used by + * the Shelf and Peaking filter types. + * \param slope Slope steepness of the transition band. + */ + void setParamsFromSlope(BiquadType type, Real f0norm, Real gain, Real slope) + { + gain = std::max(gain, 0.001f); /* Limit -60dB */ + setParams(type, f0norm, gain, rcpQFromSlope(gain, slope)); + } + + /** + * Sets the filter state for the specified filter type and its parameters. + * + * \param type The type of filter to apply. + * \param f0norm The normalized reference frequency (ref / sample_rate). + * This is the center point for the Shelf, Peaking, and BandPass filter + * types, or the cutoff frequency for the LowPass and HighPass filter + * types. + * \param gain The gain for the reference frequency response. Only used by + * the Shelf and Peaking filter types. + * \param bandwidth Normalized bandwidth of the transition band. + */ + void setParamsFromBandwidth(BiquadType type, Real f0norm, Real gain, Real bandwidth) + { setParams(type, f0norm, gain, rcpQFromBandwidth(f0norm, bandwidth)); } + + void copyParamsFrom(const BiquadFilterR &other) + { + mB0 = other.mB0; + mB1 = other.mB1; + mB2 = other.mB2; + mA1 = other.mA1; + mA2 = other.mA2; + } + + void process(const al::span src, Real *dst); + /** Processes this filter and the other at the same time. */ + void dualProcess(BiquadFilterR &other, const al::span src, Real *dst); + + /* Rather hacky. It's just here to support "manual" processing. */ + std::pair getComponents() const noexcept { return {mZ1, mZ2}; } + void setComponents(Real z1, Real z2) noexcept { mZ1 = z1; mZ2 = z2; } + Real processOne(const Real in, Real &z1, Real &z2) const noexcept + { + const Real out{in*mB0 + z1}; + z1 = in*mB1 - out*mA1 + z2; + z2 = in*mB2 - out*mA2; + return out; + } +}; + +template +struct DualBiquadR { + BiquadFilterR &f0, &f1; + + void process(const al::span src, Real *dst) + { f0.dualProcess(f1, src, dst); } +}; + +using BiquadFilter = BiquadFilterR; +using DualBiquad = DualBiquadR; + +#endif /* CORE_FILTERS_BIQUAD_H */ diff --git a/core/filters/nfc.cpp b/core/filters/nfc.cpp new file mode 100644 index 00000000..9a28517c --- /dev/null +++ b/core/filters/nfc.cpp @@ -0,0 +1,383 @@ + +#include "config.h" + +#include "nfc.h" + +#include + +#include "opthelpers.h" + + +/* Near-field control filters are the basis for handling the near-field effect. + * The near-field effect is a bass-boost present in the directional components + * of a recorded signal, created as a result of the wavefront curvature (itself + * a function of sound distance). Proper reproduction dictates this be + * compensated for using a bass-cut given the playback speaker distance, to + * avoid excessive bass in the playback. + * + * For real-time rendered audio, emulating the near-field effect based on the + * sound source's distance, and subsequently compensating for it at output + * based on the speaker distances, can create a more realistic perception of + * sound distance beyond a simple 1/r attenuation. + * + * These filters do just that. Each one applies a low-shelf filter, created as + * the combination of a bass-boost for a given sound source distance (near- + * field emulation) along with a bass-cut for a given control/speaker distance + * (near-field compensation). + * + * Note that it is necessary to apply a cut along with the boost, since the + * boost alone is unstable in higher-order ambisonics as it causes an infinite + * DC gain (even first-order ambisonics requires there to be no DC offset for + * the boost to work). Consequently, ambisonics requires a control parameter to + * be used to avoid an unstable boost-only filter. NFC-HOA defines this control + * as a reference delay, calculated with: + * + * reference_delay = control_distance / speed_of_sound + * + * This means w0 (for input) or w1 (for output) should be set to: + * + * wN = 1 / (reference_delay * sample_rate) + * + * when dealing with NFC-HOA content. For FOA input content, which does not + * specify a reference_delay variable, w0 should be set to 0 to apply only + * near-field compensation for output. It's important that w1 be a finite, + * positive, non-0 value or else the bass-boost will become unstable again. + * Also, w0 should not be too large compared to w1, to avoid excessively loud + * low frequencies. + */ + +namespace { + +constexpr float B[5][4] = { + { 0.0f }, + { 1.0f }, + { 3.0f, 3.0f }, + { 3.6778f, 6.4595f, 2.3222f }, + { 4.2076f, 11.4877f, 5.7924f, 9.1401f } +}; + +NfcFilter1 NfcFilterCreate1(const float w0, const float w1) noexcept +{ + NfcFilter1 nfc{}; + float b_00, g_0; + float r; + + nfc.base_gain = 1.0f; + nfc.gain = 1.0f; + + /* Calculate bass-boost coefficients. */ + r = 0.5f * w0; + b_00 = B[1][0] * r; + g_0 = 1.0f + b_00; + + nfc.gain *= g_0; + nfc.b1 = 2.0f * b_00 / g_0; + + /* Calculate bass-cut coefficients. */ + r = 0.5f * w1; + b_00 = B[1][0] * r; + g_0 = 1.0f + b_00; + + nfc.base_gain /= g_0; + nfc.gain /= g_0; + nfc.a1 = 2.0f * b_00 / g_0; + + return nfc; +} + +void NfcFilterAdjust1(NfcFilter1 *nfc, const float w0) noexcept +{ + const float r{0.5f * w0}; + const float b_00{B[1][0] * r}; + const float g_0{1.0f + b_00}; + + nfc->gain = nfc->base_gain * g_0; + nfc->b1 = 2.0f * b_00 / g_0; +} + + +NfcFilter2 NfcFilterCreate2(const float w0, const float w1) noexcept +{ + NfcFilter2 nfc{}; + float b_10, b_11, g_1; + float r; + + nfc.base_gain = 1.0f; + nfc.gain = 1.0f; + + /* Calculate bass-boost coefficients. */ + r = 0.5f * w0; + b_10 = B[2][0] * r; + b_11 = B[2][1] * r * r; + g_1 = 1.0f + b_10 + b_11; + + nfc.gain *= g_1; + nfc.b1 = (2.0f*b_10 + 4.0f*b_11) / g_1; + nfc.b2 = 4.0f * b_11 / g_1; + + /* Calculate bass-cut coefficients. */ + r = 0.5f * w1; + b_10 = B[2][0] * r; + b_11 = B[2][1] * r * r; + g_1 = 1.0f + b_10 + b_11; + + nfc.base_gain /= g_1; + nfc.gain /= g_1; + nfc.a1 = (2.0f*b_10 + 4.0f*b_11) / g_1; + nfc.a2 = 4.0f * b_11 / g_1; + + return nfc; +} + +void NfcFilterAdjust2(NfcFilter2 *nfc, const float w0) noexcept +{ + const float r{0.5f * w0}; + const float b_10{B[2][0] * r}; + const float b_11{B[2][1] * r * r}; + const float g_1{1.0f + b_10 + b_11}; + + nfc->gain = nfc->base_gain * g_1; + nfc->b1 = (2.0f*b_10 + 4.0f*b_11) / g_1; + nfc->b2 = 4.0f * b_11 / g_1; +} + + +NfcFilter3 NfcFilterCreate3(const float w0, const float w1) noexcept +{ + NfcFilter3 nfc{}; + float b_10, b_11, g_1; + float b_00, g_0; + float r; + + nfc.base_gain = 1.0f; + nfc.gain = 1.0f; + + /* Calculate bass-boost coefficients. */ + r = 0.5f * w0; + b_10 = B[3][0] * r; + b_11 = B[3][1] * r * r; + b_00 = B[3][2] * r; + g_1 = 1.0f + b_10 + b_11; + g_0 = 1.0f + b_00; + + nfc.gain *= g_1 * g_0; + nfc.b1 = (2.0f*b_10 + 4.0f*b_11) / g_1; + nfc.b2 = 4.0f * b_11 / g_1; + nfc.b3 = 2.0f * b_00 / g_0; + + /* Calculate bass-cut coefficients. */ + r = 0.5f * w1; + b_10 = B[3][0] * r; + b_11 = B[3][1] * r * r; + b_00 = B[3][2] * r; + g_1 = 1.0f + b_10 + b_11; + g_0 = 1.0f + b_00; + + nfc.base_gain /= g_1 * g_0; + nfc.gain /= g_1 * g_0; + nfc.a1 = (2.0f*b_10 + 4.0f*b_11) / g_1; + nfc.a2 = 4.0f * b_11 / g_1; + nfc.a3 = 2.0f * b_00 / g_0; + + return nfc; +} + +void NfcFilterAdjust3(NfcFilter3 *nfc, const float w0) noexcept +{ + const float r{0.5f * w0}; + const float b_10{B[3][0] * r}; + const float b_11{B[3][1] * r * r}; + const float b_00{B[3][2] * r}; + const float g_1{1.0f + b_10 + b_11}; + const float g_0{1.0f + b_00}; + + nfc->gain = nfc->base_gain * g_1 * g_0; + nfc->b1 = (2.0f*b_10 + 4.0f*b_11) / g_1; + nfc->b2 = 4.0f * b_11 / g_1; + nfc->b3 = 2.0f * b_00 / g_0; +} + + +NfcFilter4 NfcFilterCreate4(const float w0, const float w1) noexcept +{ + NfcFilter4 nfc{}; + float b_10, b_11, g_1; + float b_00, b_01, g_0; + float r; + + nfc.base_gain = 1.0f; + nfc.gain = 1.0f; + + /* Calculate bass-boost coefficients. */ + r = 0.5f * w0; + b_10 = B[4][0] * r; + b_11 = B[4][1] * r * r; + b_00 = B[4][2] * r; + b_01 = B[4][3] * r * r; + g_1 = 1.0f + b_10 + b_11; + g_0 = 1.0f + b_00 + b_01; + + nfc.gain *= g_1 * g_0; + nfc.b1 = (2.0f*b_10 + 4.0f*b_11) / g_1; + nfc.b2 = 4.0f * b_11 / g_1; + nfc.b3 = (2.0f*b_00 + 4.0f*b_01) / g_0; + nfc.b4 = 4.0f * b_01 / g_0; + + /* Calculate bass-cut coefficients. */ + r = 0.5f * w1; + b_10 = B[4][0] * r; + b_11 = B[4][1] * r * r; + b_00 = B[4][2] * r; + b_01 = B[4][3] * r * r; + g_1 = 1.0f + b_10 + b_11; + g_0 = 1.0f + b_00 + b_01; + + nfc.base_gain /= g_1 * g_0; + nfc.gain /= g_1 * g_0; + nfc.a1 = (2.0f*b_10 + 4.0f*b_11) / g_1; + nfc.a2 = 4.0f * b_11 / g_1; + nfc.a3 = (2.0f*b_00 + 4.0f*b_01) / g_0; + nfc.a4 = 4.0f * b_01 / g_0; + + return nfc; +} + +void NfcFilterAdjust4(NfcFilter4 *nfc, const float w0) noexcept +{ + const float r{0.5f * w0}; + const float b_10{B[4][0] * r}; + const float b_11{B[4][1] * r * r}; + const float b_00{B[4][2] * r}; + const float b_01{B[4][3] * r * r}; + const float g_1{1.0f + b_10 + b_11}; + const float g_0{1.0f + b_00 + b_01}; + + nfc->gain = nfc->base_gain * g_1 * g_0; + nfc->b1 = (2.0f*b_10 + 4.0f*b_11) / g_1; + nfc->b2 = 4.0f * b_11 / g_1; + nfc->b3 = (2.0f*b_00 + 4.0f*b_01) / g_0; + nfc->b4 = 4.0f * b_01 / g_0; +} + +} // namespace + +void NfcFilter::init(const float w1) noexcept +{ + first = NfcFilterCreate1(0.0f, w1); + second = NfcFilterCreate2(0.0f, w1); + third = NfcFilterCreate3(0.0f, w1); + fourth = NfcFilterCreate4(0.0f, w1); +} + +void NfcFilter::adjust(const float w0) noexcept +{ + NfcFilterAdjust1(&first, w0); + NfcFilterAdjust2(&second, w0); + NfcFilterAdjust3(&third, w0); + NfcFilterAdjust4(&fourth, w0); +} + + +void NfcFilter::process1(const al::span src, float *RESTRICT dst) +{ + const float gain{first.gain}; + const float b1{first.b1}; + const float a1{first.a1}; + float z1{first.z[0]}; + auto proc_sample = [gain,b1,a1,&z1](const float in) noexcept -> float + { + const float y{in*gain - a1*z1}; + const float out{y + b1*z1}; + z1 += y; + return out; + }; + std::transform(src.cbegin(), src.cend(), dst, proc_sample); + first.z[0] = z1; +} + +void NfcFilter::process2(const al::span src, float *RESTRICT dst) +{ + const float gain{second.gain}; + const float b1{second.b1}; + const float b2{second.b2}; + const float a1{second.a1}; + const float a2{second.a2}; + float z1{second.z[0]}; + float z2{second.z[1]}; + auto proc_sample = [gain,b1,b2,a1,a2,&z1,&z2](const float in) noexcept -> float + { + const float y{in*gain - a1*z1 - a2*z2}; + const float out{y + b1*z1 + b2*z2}; + z2 += z1; + z1 += y; + return out; + }; + std::transform(src.cbegin(), src.cend(), dst, proc_sample); + second.z[0] = z1; + second.z[1] = z2; +} + +void NfcFilter::process3(const al::span src, float *RESTRICT dst) +{ + const float gain{third.gain}; + const float b1{third.b1}; + const float b2{third.b2}; + const float b3{third.b3}; + const float a1{third.a1}; + const float a2{third.a2}; + const float a3{third.a3}; + float z1{third.z[0]}; + float z2{third.z[1]}; + float z3{third.z[2]}; + auto proc_sample = [gain,b1,b2,b3,a1,a2,a3,&z1,&z2,&z3](const float in) noexcept -> float + { + float y{in*gain - a1*z1 - a2*z2}; + float out{y + b1*z1 + b2*z2}; + z2 += z1; + z1 += y; + + y = out - a3*z3; + out = y + b3*z3; + z3 += y; + return out; + }; + std::transform(src.cbegin(), src.cend(), dst, proc_sample); + third.z[0] = z1; + third.z[1] = z2; + third.z[2] = z3; +} + +void NfcFilter::process4(const al::span src, float *RESTRICT dst) +{ + const float gain{fourth.gain}; + const float b1{fourth.b1}; + const float b2{fourth.b2}; + const float b3{fourth.b3}; + const float b4{fourth.b4}; + const float a1{fourth.a1}; + const float a2{fourth.a2}; + const float a3{fourth.a3}; + const float a4{fourth.a4}; + float z1{fourth.z[0]}; + float z2{fourth.z[1]}; + float z3{fourth.z[2]}; + float z4{fourth.z[3]}; + auto proc_sample = [gain,b1,b2,b3,b4,a1,a2,a3,a4,&z1,&z2,&z3,&z4](const float in) noexcept -> float + { + float y{in*gain - a1*z1 - a2*z2}; + float out{y + b1*z1 + b2*z2}; + z2 += z1; + z1 += y; + + y = out - a3*z3 - a4*z4; + out = y + b3*z3 + b4*z4; + z4 += z3; + z3 += y; + return out; + }; + std::transform(src.cbegin(), src.cend(), dst, proc_sample); + fourth.z[0] = z1; + fourth.z[1] = z2; + fourth.z[2] = z3; + fourth.z[3] = z4; +} diff --git a/core/filters/nfc.h b/core/filters/nfc.h new file mode 100644 index 00000000..33f67a5f --- /dev/null +++ b/core/filters/nfc.h @@ -0,0 +1,63 @@ +#ifndef CORE_FILTERS_NFC_H +#define CORE_FILTERS_NFC_H + +#include + +#include "alspan.h" + + +struct NfcFilter1 { + float base_gain, gain; + float b1, a1; + float z[1]; +}; +struct NfcFilter2 { + float base_gain, gain; + float b1, b2, a1, a2; + float z[2]; +}; +struct NfcFilter3 { + float base_gain, gain; + float b1, b2, b3, a1, a2, a3; + float z[3]; +}; +struct NfcFilter4 { + float base_gain, gain; + float b1, b2, b3, b4, a1, a2, a3, a4; + float z[4]; +}; + +class NfcFilter { + NfcFilter1 first; + NfcFilter2 second; + NfcFilter3 third; + NfcFilter4 fourth; + +public: + /* NOTE: + * w0 = speed_of_sound / (source_distance * sample_rate); + * w1 = speed_of_sound / (control_distance * sample_rate); + * + * Generally speaking, the control distance should be approximately the + * average speaker distance, or based on the reference delay if outputing + * NFC-HOA. It must not be negative, 0, or infinite. The source distance + * should not be too small relative to the control distance. + */ + + void init(const float w1) noexcept; + void adjust(const float w0) noexcept; + + /* Near-field control filter for first-order ambisonic channels (1-3). */ + void process1(const al::span src, float *RESTRICT dst); + + /* Near-field control filter for second-order ambisonic channels (4-8). */ + void process2(const al::span src, float *RESTRICT dst); + + /* Near-field control filter for third-order ambisonic channels (9-15). */ + void process3(const al::span src, float *RESTRICT dst); + + /* Near-field control filter for fourth-order ambisonic channels (16-24). */ + void process4(const al::span src, float *RESTRICT dst); +}; + +#endif /* CORE_FILTERS_NFC_H */ diff --git a/core/filters/splitter.cpp b/core/filters/splitter.cpp new file mode 100644 index 00000000..5cc670b7 --- /dev/null +++ b/core/filters/splitter.cpp @@ -0,0 +1,113 @@ + +#include "config.h" + +#include "splitter.h" + +#include +#include +#include + +#include "math_defs.h" +#include "opthelpers.h" + + +template +void BandSplitterR::init(Real f0norm) +{ + const Real w{f0norm * al::MathDefs::Tau()}; + const Real cw{std::cos(w)}; + if(cw > std::numeric_limits::epsilon()) + mCoeff = (std::sin(w) - 1.0f) / cw; + else + mCoeff = cw * -0.5f; + + mLpZ1 = 0.0f; + mLpZ2 = 0.0f; + mApZ1 = 0.0f; +} + +template +void BandSplitterR::process(const al::span input, Real *hpout, Real *lpout) +{ + const Real ap_coeff{mCoeff}; + const Real lp_coeff{mCoeff*0.5f + 0.5f}; + Real lp_z1{mLpZ1}; + Real lp_z2{mLpZ2}; + Real ap_z1{mApZ1}; + auto proc_sample = [ap_coeff,lp_coeff,&lp_z1,&lp_z2,&ap_z1,&lpout](const Real in) noexcept -> Real + { + /* Low-pass sample processing. */ + Real d{(in - lp_z1) * lp_coeff}; + Real lp_y{lp_z1 + d}; + lp_z1 = lp_y + d; + + d = (lp_y - lp_z2) * lp_coeff; + lp_y = lp_z2 + d; + lp_z2 = lp_y + d; + + *(lpout++) = lp_y; + + /* All-pass sample processing. */ + Real ap_y{in*ap_coeff + ap_z1}; + ap_z1 = in - ap_y*ap_coeff; + + /* High-pass generated from removing low-passed output. */ + return ap_y - lp_y; + }; + std::transform(input.cbegin(), input.cend(), hpout, proc_sample); + mLpZ1 = lp_z1; + mLpZ2 = lp_z2; + mApZ1 = ap_z1; +} + +template +void BandSplitterR::processHfScale(const al::span samples, const Real hfscale) +{ + const Real ap_coeff{mCoeff}; + const Real lp_coeff{mCoeff*0.5f + 0.5f}; + Real lp_z1{mLpZ1}; + Real lp_z2{mLpZ2}; + Real ap_z1{mApZ1}; + auto proc_sample = [hfscale,ap_coeff,lp_coeff,&lp_z1,&lp_z2,&ap_z1](const Real in) noexcept -> Real + { + /* Low-pass sample processing. */ + Real d{(in - lp_z1) * lp_coeff}; + Real lp_y{lp_z1 + d}; + lp_z1 = lp_y + d; + + d = (lp_y - lp_z2) * lp_coeff; + lp_y = lp_z2 + d; + lp_z2 = lp_y + d; + + /* All-pass sample processing. */ + Real ap_y{in*ap_coeff + ap_z1}; + ap_z1 = in - ap_y*ap_coeff; + + /* High-pass generated by removing the low-passed signal, which is then + * scaled and added back to the low-passed signal. + */ + return (ap_y-lp_y)*hfscale + lp_y; + }; + std::transform(samples.begin(), samples.end(), samples.begin(), proc_sample); + mLpZ1 = lp_z1; + mLpZ2 = lp_z2; + mApZ1 = ap_z1; +} + +template +void BandSplitterR::applyAllpass(const al::span samples) const +{ + const Real coeff{mCoeff}; + Real z1{0.0f}; + auto proc_sample = [coeff,&z1](const Real in) noexcept -> Real + { + const Real out{in*coeff + z1}; + z1 = in - out*coeff; + return out; + }; + std::transform(samples.begin(), samples.end(), samples.begin(), proc_sample); +} + + +template class BandSplitterR; +template class BandSplitterR; diff --git a/core/filters/splitter.h b/core/filters/splitter.h new file mode 100644 index 00000000..ba548c10 --- /dev/null +++ b/core/filters/splitter.h @@ -0,0 +1,36 @@ +#ifndef CORE_FILTERS_SPLITTER_H +#define CORE_FILTERS_SPLITTER_H + +#include + +#include "alspan.h" + + +/* Band splitter. Splits a signal into two phase-matching frequency bands. */ +template +class BandSplitterR { + Real mCoeff{0.0f}; + Real mLpZ1{0.0f}; + Real mLpZ2{0.0f}; + Real mApZ1{0.0f}; + +public: + BandSplitterR() = default; + BandSplitterR(const BandSplitterR&) = default; + BandSplitterR(Real f0norm) { init(f0norm); } + + void init(Real f0norm); + void clear() noexcept { mLpZ1 = mLpZ2 = mApZ1 = 0.0f; } + void process(const al::span input, Real *hpout, Real *lpout); + + void processHfScale(const al::span samples, const Real hfscale); + + /* The all-pass portion of the band splitter. Applies the same phase shift + * without splitting the signal. Note that each use of this method is + * indepedent, it does not track history between calls. + */ + void applyAllpass(const al::span samples) const; +}; +using BandSplitter = BandSplitterR; + +#endif /* CORE_FILTERS_SPLITTER_H */ -- cgit v1.2.3