From 3057163e8e22bfdd91dcde01cdc426c4c191a0e5 Mon Sep 17 00:00:00 2001 From: kunitoki Date: Thu, 21 May 2026 20:57:09 +0200 Subject: [PATCH 1/3] Oversampler and Resampler --- .../yup_dsp/resampling/yup_CircularBuffer.h | 95 +++++ modules/yup_dsp/resampling/yup_Oversampler.h | 365 ++++++++++++++++++ modules/yup_dsp/resampling/yup_Resampler.h | 242 ++++++++++++ modules/yup_dsp/resampling/yup_SincTable.h | 178 +++++++++ modules/yup_dsp/yup_dsp.h | 6 + tests/yup_dsp.cpp | 5 + tests/yup_dsp/yup_CircularBuffer.cpp | 164 ++++++++ tests/yup_dsp/yup_Oversampler.cpp | 241 ++++++++++++ tests/yup_dsp/yup_Resampler.cpp | 251 ++++++++++++ tests/yup_dsp/yup_SincTable.cpp | 176 +++++++++ 10 files changed, 1723 insertions(+) create mode 100644 modules/yup_dsp/resampling/yup_CircularBuffer.h create mode 100644 modules/yup_dsp/resampling/yup_Oversampler.h create mode 100644 modules/yup_dsp/resampling/yup_Resampler.h create mode 100644 modules/yup_dsp/resampling/yup_SincTable.h create mode 100644 tests/yup_dsp/yup_CircularBuffer.cpp create mode 100644 tests/yup_dsp/yup_Oversampler.cpp create mode 100644 tests/yup_dsp/yup_Resampler.cpp create mode 100644 tests/yup_dsp/yup_SincTable.cpp diff --git a/modules/yup_dsp/resampling/yup_CircularBuffer.h b/modules/yup_dsp/resampling/yup_CircularBuffer.h new file mode 100644 index 000000000..f5687b5e6 --- /dev/null +++ b/modules/yup_dsp/resampling/yup_CircularBuffer.h @@ -0,0 +1,95 @@ +/* + ============================================================================== + + This file is part of the YUP library. + Copyright (c) 2026 - kunitoki@gmail.com + + YUP is an open source library subject to open-source licensing. + + The code included in this file is provided under the terms of the ISC license + http://www.isc.org/downloads/software-support-policy/isc-license. Permission + to use, copy, modify, and/or distribute this software for any purpose with or + without fee is hereby granted provided that the above copyright notice and + this permission notice appear in all copies. + + YUP IS PROVIDED "AS IS" WITHOUT ANY WARRANTY, AND ALL WARRANTIES, WHETHER + EXPRESSED OR IMPLIED, INCLUDING MERCHANTABILITY AND FITNESS FOR PURPOSE, ARE + DISCLAIMED. + + ============================================================================== +*/ + +#pragma once + +namespace yup +{ + +//============================================================================== +/** + Fixed-size circular (ring) buffer for efficient sample history access. + + Provides O(1) push and random-read operations. Logically, index 0 refers to + the oldest element and index BufferSize - 1 to the most recently pushed one. + + Intended as a building block for block-based DSP algorithms that require a + sliding window of sample history, such as polyphase resampling and FIR + convolution. + + @tparam SampleType Element type stored in the buffer. + @tparam BufferSize Compile-time capacity (must be > 0). +*/ +template +class CircularBuffer +{ +public: + static_assert (BufferSize > 0, "CircularBuffer BufferSize must be greater than zero"); + + //============================================================================== + /** Default constructor. */ + CircularBuffer() noexcept + { + buffer.fill (SampleType {}); + } + + /** Constructs a buffer with all entries initialised to initValue. */ + explicit CircularBuffer (SampleType initValue) noexcept + { + buffer.fill (initValue); + } + + //============================================================================== + /** Inserts value, overwriting the oldest entry and advancing the write pointer. */ + forcedinline void push (SampleType value) noexcept + { + buffer[static_cast (oldestIndex)] = value; + + if (++oldestIndex == BufferSize) + oldestIndex = 0; + } + + /** Returns the element at logical index (0 = oldest, BufferSize-1 = newest). */ + forcedinline SampleType& operator[] (int index) noexcept + { + return buffer[static_cast ((oldestIndex + index) % BufferSize)]; + } + + /** Returns the element at logical index (0 = oldest, BufferSize-1 = newest). */ + const forcedinline SampleType& operator[] (int index) const noexcept + { + return buffer[static_cast ((oldestIndex + index) % BufferSize)]; + } + + //============================================================================== + /** Resets all entries to zero and rewinds the write pointer. */ + void clear() noexcept + { + buffer.fill (SampleType {}); + oldestIndex = 0; + } + +private: + std::array (BufferSize)> buffer {}; + int oldestIndex = 0; +}; + +} // namespace yup diff --git a/modules/yup_dsp/resampling/yup_Oversampler.h b/modules/yup_dsp/resampling/yup_Oversampler.h new file mode 100644 index 000000000..824d8f178 --- /dev/null +++ b/modules/yup_dsp/resampling/yup_Oversampler.h @@ -0,0 +1,365 @@ +/* + ============================================================================== + + This file is part of the YUP library. + Copyright (c) 2026 - kunitoki@gmail.com + + YUP is an open source library subject to open-source licensing. + + The code included in this file is provided under the terms of the ISC license + http://www.isc.org/downloads/software-support-policy/isc-license. Permission + to use, copy, modify, and/or distribute this software for any purpose with or + without fee is hereby granted provided that the above copyright notice and + this permission notice appear in all copies. + + YUP IS PROVIDED "AS IS" WITHOUT ANY WARRANTY, AND ALL WARRANTIES, WHETHER + EXPRESSED OR IMPLIED, INCLUDING MERCHANTABILITY AND FITNESS FOR PURPOSE, ARE + DISCLAIMED. + + ============================================================================== +*/ + +#pragma once + +namespace yup +{ + +//============================================================================== +/** + Multi-channel integer-factor oversampler using windowed sinc interpolation. + + Oversampler up- and downsamples audio by an integer factor with + bandlimited interpolation and anti-aliasing. Internal per-channel history + buffers allow seamless multi-block (real-time) operation. + + Typical usage in an audio effect: + @code + yup::Oversampler os; + os.prepare (44100.0, 2, 512); + + // Inside your audio callback: + os.upsample (inputPtrs, numChannels, numSamples); + os.processOversampledBlock ([&] (auto& buf) + { + applyDistortion (buf); // buf is std::vector>& + }); + os.downsample (outputPtrs, numChannels, numSamples); + @endcode + + @tparam SampleType Audio sample type (float or double). + @tparam OversampleFactor Integer upsample ratio (2, 4, 8, …). + @tparam SincRadius Half-width of the sinc kernel in original-rate samples. + Higher values give better stopband rejection at the + cost of more computation. + @tparam CoeffType Precision for internal filter coefficients (default double). +*/ +template +class Oversampler +{ +public: + static_assert (OversampleFactor >= 2, "OversampleFactor must be at least 2"); + static_assert (SincRadius >= 1, "SincRadius must be at least 1"); + + //============================================================================== + /** Default constructor. Call prepare() before any processing. */ + Oversampler() = default; + + /** Destructor. */ + ~Oversampler() = default; + + //============================================================================== + /** + Prepares the oversampler for processing. + + Configures the internal windowed sinc table and allocates per-channel + history and staging buffers. Must be called before upsample() or downsample(). + + @param sampleRate Input sample rate in Hz. + @param maxChannels Maximum number of audio channels. + @param maxBlockSize Maximum input block size in samples. + */ + void prepare (double sampleRate, int maxChannels, int maxBlockSize) + { + jassert (sampleRate > 0.0 && maxChannels > 0 && maxBlockSize > 0); + + sincTable.configure (static_cast (sampleRate)); + sincTable.applyKaiserWindow (CoeffType (5)); + + const int maxInterpolated = maxBlockSize * OversampleFactor; + + interpolBeginBufs.assign (maxChannels, CircularBuffer {}); + interpolEndBufs.assign (maxChannels, CircularBuffer {}); + decimBeginBufs.assign (maxChannels, CircularBuffer {}); + decimEndBufs.assign (maxChannels, CircularBuffer {}); + + xInterp.assign (maxChannels, std::vector (static_cast (maxBlockSize + SincRadius), SampleType {})); + xDecim.assign (maxChannels, std::vector (static_cast (maxInterpolated + SincRadius * OversampleFactor), SampleType {})); + oversampledBuffer.assign (maxChannels, std::vector (static_cast (maxInterpolated), SampleType {})); + currentOversampledSize = 0; + } + + /** + Resets all internal processing state. + + Clears all history buffers so that a fresh processing session can begin + without artifacts from a previous session. Filter coefficients are + preserved; there is no need to call prepare() again. + */ + void reset() noexcept + { + for (auto& b : interpolBeginBufs) + b.clear(); + for (auto& b : interpolEndBufs) + b.clear(); + for (auto& b : decimBeginBufs) + b.clear(); + for (auto& b : decimEndBufs) + b.clear(); + + for (auto& ch : xInterp) + std::fill (ch.begin(), ch.end(), SampleType {}); + for (auto& ch : xDecim) + std::fill (ch.begin(), ch.end(), SampleType {}); + for (auto& ch : oversampledBuffer) + std::fill (ch.begin(), ch.end(), SampleType {}); + currentOversampledSize = 0; + } + + //============================================================================== + /** + Upsample an input block into the internal oversampled buffer. + + After this call the internal buffer holds numSamples * OversampleFactor + bandlimited interpolated samples per channel, accessible via + getOversampledChannelData() or processOversampledBlock(). + + @param input Array of read pointers, one per channel (channel-major). + @param numChannels Number of channels to process (must be <= maxChannels from prepare()). + @param numSamples Number of input samples per channel (must be <= maxBlockSize). + */ + void upsample (const SampleType* const* input, int numChannels, int numSamples) noexcept + { + jassert (numChannels > 0 && numSamples > 0); + jassert (numChannels <= static_cast (xInterp.size())); + jassert (numSamples + SincRadius <= static_cast (xInterp[0].size())); + + for (int ch = 0; ch < numChannels; ++ch) + { + auto& xBuf = xInterp[static_cast (ch)]; + auto& endBuf = interpolEndBufs[static_cast (ch)]; + + for (int i = 0; i < numSamples + SincRadius; ++i) + { + xBuf[static_cast (i)] = (i >= SincRadius) + ? input[ch][i - SincRadius] + : endBuf[i]; + } + } + + currentOversampledSize = numSamples * OversampleFactor; + + for (int ch = 0; ch < numChannels; ++ch) + { + auto& xBuf = xInterp[static_cast (ch)]; + auto& outBuf = oversampledBuffer[static_cast (ch)]; + auto& beginBuf = interpolBeginBufs[static_cast (ch)]; + auto& endBuf = interpolEndBufs[static_cast (ch)]; + + outBuf[0] = xBuf[0]; + + for (int k = 1; k < currentOversampledSize; ++k) + { + const int delta = k % OversampleFactor; + const int index = k / OversampleFactor; + + if (delta != 0) + { + CoeffType acc = CoeffType (0); + + for (int n = -SincRadius; n <= 0; ++n) + acc += sincTable (n, delta) + * static_cast (xBuf[static_cast (index - n)]); + + for (int n = 1; n <= SincRadius; ++n) + acc += sincTable (n, delta) + * static_cast (beginBuf[SincRadius - n]); + + outBuf[static_cast (k)] = static_cast (acc); + } + else + { + outBuf[static_cast (k)] = xBuf[static_cast (index)]; + beginBuf.push (xBuf[static_cast (index - 1)]); + } + } + + beginBuf.push (xBuf[static_cast (numSamples - 1)]); + + for (int i = 0; i < SincRadius; ++i) + endBuf.push (xBuf[static_cast (numSamples + i)]); + } + } + + /** + Downsample the internal oversampled buffer into an output block. + + Applies a lowpass anti-aliasing FIR to the oversampled data and decimates + by OversampleFactor. Must be called after the oversampled buffer has been + processed (e.g. via processOversampledBlock()). + + @param output Array of write pointers, one per channel. + @param numChannels Number of channels to write (must be <= maxChannels from prepare()). + @param numSamples Number of output samples per channel (must match the numSamples + passed to the preceding upsample() call). + */ + void downsample (SampleType* const* output, int numChannels, int numSamples) noexcept + { + jassert (numChannels > 0 && numSamples > 0); + jassert (numChannels <= static_cast (xDecim.size())); + jassert (currentOversampledSize > 0); + + const int interpolatedSize = currentOversampledSize; + + for (int ch = 0; ch < numChannels; ++ch) + { + auto& xBuf = xDecim[static_cast (ch)]; + auto& inBuf = oversampledBuffer[static_cast (ch)]; + auto& dEndBuf = decimEndBufs[static_cast (ch)]; + + for (int i = 0; i < interpolatedSize + SincRadius * OversampleFactor; ++i) + { + xBuf[static_cast (i)] = (i >= SincRadius * OversampleFactor) + ? inBuf[static_cast (i - SincRadius * OversampleFactor)] + : dEndBuf[i]; + } + } + + for (int ch = 0; ch < numChannels; ++ch) + { + auto& xBuf = xDecim[static_cast (ch)]; + auto& beginBuf = decimBeginBufs[static_cast (ch)]; + auto& dEndBuf = decimEndBufs[static_cast (ch)]; + + for (int k = 0; k < numSamples; ++k) + { + const int index = OversampleFactor * k; + CoeffType acc = CoeffType (0); + + for (int n = 1; n <= SincRadius * OversampleFactor; ++n) + acc += sincTable[n] + * static_cast (beginBuf[SincRadius * OversampleFactor - n]); + + for (int n = 0; n >= -(SincRadius * OversampleFactor); --n) + acc += sincTable[n] + * static_cast (xBuf[static_cast (index - n)]); + + for (int i = 0; i < OversampleFactor; ++i) + beginBuf.push (xBuf[static_cast (index + i)]); + + output[ch][k] = static_cast (acc / static_cast (OversampleFactor)); + } + + for (int i = 0; i < SincRadius * OversampleFactor; ++i) + dEndBuf.push (xBuf[static_cast (interpolatedSize + i)]); + } + } + + //============================================================================== + /** + Invokes a callback with the internal oversampled multi-channel buffer. + + The callback receives a reference to the internal + `std::vector>`, where each inner vector has + getOversampledNumSamples() elements. Use this to apply processing at the + elevated sample rate. + + @param callback Callable with signature + `void(std::vector>&)`. + */ + template + void processOversampledBlock (Callable&& callback) + { + callback (oversampledBuffer); + } + + //============================================================================== + /** + Returns a writable pointer to the data for a single oversampled channel. + + @param channel Zero-based channel index. + @return Pointer to getOversampledNumSamples() contiguous samples, + or nullptr if the channel index is out of range or + prepare() has not been called. + */ + forcedinline SampleType* getOversampledChannelData (int channel) noexcept + { + if (channel < 0 || channel >= static_cast (oversampledBuffer.size())) + return nullptr; + + return oversampledBuffer[static_cast (channel)].data(); + } + + /** + Returns a read-only pointer to the data for a single oversampled channel. + + @param channel Zero-based channel index. + @return Pointer to getOversampledNumSamples() contiguous samples, + or nullptr if the channel index is out of range. + */ + const forcedinline SampleType* getOversampledChannelData (int channel) const noexcept + { + if (channel < 0 || channel >= static_cast (oversampledBuffer.size())) + return nullptr; + + return oversampledBuffer[static_cast (channel)].data(); + } + + /** + Returns the number of samples currently in each oversampled channel. + + Equal to the numSamples argument of the most recent upsample() call + multiplied by OversampleFactor. Returns 0 before the first upsample() + call or after reset(). + */ + forcedinline int getOversampledNumSamples() const noexcept + { + return currentOversampledSize; + } + + /** + Returns the processing latency introduced by the oversampler. + + @return Latency in input-rate samples (= 2 * SincRadius). + */ + forcedinline int getLatencyInSamples() const noexcept + { + return 2 * SincRadius; + } + +private: + //============================================================================== + SincTable sincTable; + + std::vector> interpolBeginBufs; + std::vector> interpolEndBufs; + std::vector> decimBeginBufs; + std::vector> decimEndBufs; + + std::vector> xInterp; + std::vector> xDecim; + std::vector> oversampledBuffer; + int currentOversampledSize = 0; + + YUP_DECLARE_NON_COPYABLE_WITH_LEAK_DETECTOR (Oversampler) +}; + +//============================================================================== +/** @name Convenience type aliases for common oversampling configurations */ +using Oversampler2xFloat = Oversampler; /**< 2x oversampler, float, 8-tap radius */ +using Oversampler4xFloat = Oversampler; /**< 4x oversampler, float, 8-tap radius */ +using Oversampler8xFloat = Oversampler; /**< 8x oversampler, float, 8-tap radius */ +using Oversampler2xDouble = Oversampler; /**< 2x oversampler, double, 8-tap radius */ +using Oversampler4xDouble = Oversampler; /**< 4x oversampler, double, 8-tap radius */ +using Oversampler8xDouble = Oversampler; /**< 8x oversampler, double, 8-tap radius */ + +} // namespace yup diff --git a/modules/yup_dsp/resampling/yup_Resampler.h b/modules/yup_dsp/resampling/yup_Resampler.h new file mode 100644 index 000000000..cc1ac76e7 --- /dev/null +++ b/modules/yup_dsp/resampling/yup_Resampler.h @@ -0,0 +1,242 @@ +/* + ============================================================================== + + This file is part of the YUP library. + Copyright (c) 2026 - kunitoki@gmail.com + + YUP is an open source library subject to open-source licensing. + + The code included in this file is provided under the terms of the ISC license + http://www.isc.org/downloads/software-support-policy/isc-license. Permission + to use, copy, modify, and/or distribute this software for any purpose with or + without fee is hereby granted provided that the above copyright notice and + this permission notice appear in all copies. + + YUP IS PROVIDED "AS IS" WITHOUT ANY WARRANTY, AND ALL WARRANTIES, WHETHER + EXPRESSED OR IMPLIED, INCLUDING MERCHANTABILITY AND FITNESS FOR PURPOSE, ARE + DISCLAIMED. + + ============================================================================== +*/ + +#pragma once + +namespace yup +{ + +//============================================================================== +/** + Multi-channel asynchronous resampler for arbitrary sample-rate conversion. + + Resampler converts between any two sample rates (including non-integer + ratios such as 44100 Hz → 48000 Hz) using a polyphase windowed sinc + filter with high-resolution phase lookup. Per-channel history buffers + ensure seamless multi-block (real-time) processing. + + Key behaviours: + - Automatically scales gain to compensate when downsampling. + - Phase state persists across blocks, so consecutive calls produce a + continuous, gapless output stream. + - Call reset() to restart the phase accumulator (e.g. after a transport loop). + + @code + yup::Resampler resampler; + resampler.prepare (44100.0, 48000.0, 2, 512); + + // Inside your audio callback: + int produced = resampler.resample (inputPtrs, outputPtrs, 2, numSamples); + @endcode + + @tparam SampleType Audio sample type (float or double). + @tparam SincRadius Half-width of the sinc kernel in source-rate samples. + @tparam Resolution Number of fractional phase sub-steps (default 256). + Higher values give more accurate interpolation at the + cost of a larger lookup table. + @tparam CoeffType Precision for internal filter coefficients (default double). +*/ +template +class Resampler +{ +public: + static_assert (SincRadius >= 1, "SincRadius must be at least 1"); + static_assert (Resolution >= 2, "Resolution must be at least 2"); + + //============================================================================== + /** Default constructor. Call prepare() before any processing. */ + Resampler() = default; + + /** Destructor. */ + ~Resampler() = default; + + //============================================================================== + /** + Prepares the resampler for processing. + + Configures the internal windowed sinc table, allocates per-channel + history buffers, and initialises the phase accumulator. Must be called + before resample(). + + @param sourceSampleRate Sample rate of the input signal in Hz. + @param targetSampleRate Desired output sample rate in Hz. + @param maxChannels Maximum number of audio channels. + @param maxBlockSize Maximum number of input samples per call to resample(). + */ + void prepare (double sourceSampleRate, double targetSampleRate, int maxChannels, int maxBlockSize) + { + jassert (sourceSampleRate > 0.0 && targetSampleRate > 0.0); + jassert (maxChannels > 0 && maxBlockSize > 0); + + const CoeffType cutoff = static_cast ( + std::min (sourceSampleRate, targetSampleRate) / 2.0); + + sincTable.configureWithCutoff (cutoff, static_cast (sourceSampleRate)); + sincTable.applyKaiserWindow (CoeffType (5)); + + currentPhase = 0.0; + oversampleFactor = targetSampleRate / sourceSampleRate; + maxOutputSamples = static_cast (maxBlockSize * oversampleFactor) + 1; + + beginBufs.assign (maxChannels, CircularBuffer {}); + endBufs.assign (maxChannels, CircularBuffer {}); + + xBufs.assign (maxChannels, std::vector (static_cast (maxBlockSize + SincRadius + 1), SampleType {})); + } + + /** + Resets the phase accumulator and all history buffers. + + Call this when restarting processing after a discontinuity such as a + transport loop. The filter coefficients remain valid; there is no need + to call prepare() again. + */ + void reset() noexcept + { + currentPhase = 0.0; + + for (auto& b : beginBufs) + b.clear(); + + for (auto& b : endBufs) + b.clear(); + + for (auto& ch : xBufs) + std::fill (ch.begin(), ch.end(), SampleType {}); + } + + //============================================================================== + /** + Converts numSamples input samples into output samples at the target rate. + + @param input Array of read pointers, one per channel (channel-major). + @param output Array of write pointers, one per channel. The caller + must allocate at least ceil(numSamples * targetRate / sourceRate) + 1 + samples per channel in the output buffers. + @param numChannels Number of channels to process. + @param numSamples Number of input samples per channel. + @return Number of output samples written per channel. + */ + int resample (const SampleType* const* input, SampleType* const* output, int numChannels, int numSamples) noexcept + { + jassert (numChannels > 0 && numSamples > 0); + jassert (numChannels <= static_cast (xBufs.size())); + jassert (numSamples + SincRadius + 1 <= static_cast (xBufs[0].size())); + + for (int ch = 0; ch < numChannels; ++ch) + { + auto& xBuf = xBufs[static_cast (ch)]; + auto& eBuf = endBufs[static_cast (ch)]; + + for (int i = 0; i < numSamples + SincRadius; ++i) + { + xBuf[static_cast (i)] = (i >= SincRadius) + ? input[ch][i - SincRadius] + : eBuf[i]; + } + } + + const int outputCount = static_cast ((numSamples - currentPhase) * oversampleFactor); + const CoeffType gainScale = (oversampleFactor < 1.0) + ? static_cast (oversampleFactor) + : CoeffType (1); + + for (int ch = 0; ch < numChannels; ++ch) + { + auto& xBuf = xBufs[static_cast (ch)]; + auto& bBuf = beginBufs[static_cast (ch)]; + + for (int k = 0; k < outputCount; ++k) + { + const double virtualIndex = static_cast (k) / oversampleFactor + currentPhase; + const int index = static_cast (virtualIndex); + const int delta = static_cast ((virtualIndex - index) * Resolution); + + if (delta != 0) + { + CoeffType acc = CoeffType (0); + + for (int n = -SincRadius; n <= 0; ++n) + acc += sincTable (n, delta) + * static_cast (xBuf[static_cast (index - n)]); + + for (int n = 1; n <= SincRadius; ++n) + acc += sincTable (n, delta) + * static_cast (bBuf[SincRadius - n]); + + output[ch][k] = static_cast (acc * gainScale); + } + else + { + output[ch][k] = static_cast ( + static_cast (xBuf[static_cast (index)]) * gainScale); + + if (index > 0) + bBuf.push (xBuf[static_cast (index - 1)]); + else + bBuf.push (bBuf[SincRadius - 1]); + } + } + + beginBufs[static_cast (ch)].push (xBuf[static_cast (numSamples - 1)]); + + for (int i = 0; i < SincRadius; ++i) + endBufs[static_cast (ch)].push (xBuf[static_cast (numSamples + i)]); + } + + currentPhase = (currentPhase + static_cast (outputCount) / oversampleFactor) - numSamples; + return outputCount; + } + + //============================================================================== + /** + Returns the processing latency introduced by the resampler. + + @return Latency in input-rate samples (= SincRadius). + */ + int getLatencyInSamples() const noexcept + { + return SincRadius; + } + +private: + //============================================================================== + SincTable sincTable; + + std::vector> beginBufs; + std::vector> endBufs; + std::vector> xBufs; + + double currentPhase = 0.0; + double oversampleFactor = 1.0; + int maxOutputSamples = 0; + + YUP_DECLARE_NON_COPYABLE_WITH_LEAK_DETECTOR (Resampler) +}; + +//============================================================================== +/** @name Convenience type aliases for common resampling configurations */ +///@{ +using ResamplerFloat = Resampler; /**< Resampler for float samples, 8-tap radius */ +using ResamplerDouble = Resampler; /**< Resampler for double samples, 8-tap radius */ +///@} + +} // namespace yup diff --git a/modules/yup_dsp/resampling/yup_SincTable.h b/modules/yup_dsp/resampling/yup_SincTable.h new file mode 100644 index 000000000..b2f06dafb --- /dev/null +++ b/modules/yup_dsp/resampling/yup_SincTable.h @@ -0,0 +1,178 @@ +/* + ============================================================================== + + This file is part of the YUP library. + Copyright (c) 2026 - kunitoki@gmail.com + + YUP is an open source library subject to open-source licensing. + + The code included in this file is provided under the terms of the ISC license + http://www.isc.org/downloads/software-support-policy/isc-license. Permission + to use, copy, modify, and/or distribute this software for any purpose with or + without fee is hereby granted provided that the above copyright notice and + this permission notice appear in all copies. + + YUP IS PROVIDED "AS IS" WITHOUT ANY WARRANTY, AND ALL WARRANTIES, WHETHER + EXPRESSED OR IMPLIED, INCLUDING MERCHANTABILITY AND FITNESS FOR PURPOSE, ARE + DISCLAIMED. + + ============================================================================== +*/ + +#pragma once + +namespace yup +{ + +//============================================================================== +/** + Pre-computed, windowed sinc lookup table for polyphase sample-rate conversion. + + Stores the positive half of a symmetric windowed sinc kernel. The + OversampleFactor parameter controls the fractional-phase resolution: table + entry (tap, delta) corresponds to the sinc value at + t = tap + delta / OversampleFactor. + + @tparam CoeffType Precision for stored values (float or double). + @tparam OversampleFactor Number of fractional phase sub-steps per tap. + @tparam SincRadius Number of taps on each side of the kernel centre. +*/ +template +class SincTable +{ +public: + static_assert (OversampleFactor > 0, "SincTable OversampleFactor must be > 0"); + static_assert (SincRadius > 0, "SincTable SincRadius must be > 0"); + + /** Number of entries in the half-kernel storage. */ + static constexpr int tableSize = (SincRadius + 1) * OversampleFactor; + + //============================================================================== + /** Default constructor. */ + SincTable() noexcept = default; + + //============================================================================== + /** + Computes the sinc kernel with cutoff at sampleRate / 2. + + Suitable for integer-factor upsampling where the anti-image cutoff + equals the input Nyquist frequency. + + @param sampleRate Input sample rate in Hz. + */ + void configure (CoeffType sampleRate) noexcept + { + const CoeffType fc = sampleRate / CoeffType (2); + const CoeffType T = CoeffType (1) / sampleRate; + + for (int i = 0; i <= SincRadius; ++i) + { + for (int delta = 0; delta < OversampleFactor; ++delta) + { + const CoeffType t = static_cast (i) + + static_cast (delta) / static_cast (OversampleFactor); + + if (t == CoeffType (0)) + (*this) (i, delta) = CoeffType (1); + else + (*this) (i, delta) = std::sin (MathConstants::twoPi * fc * t * T) + / (MathConstants::twoPi * fc * t * T); + } + } + } + + /** + Computes the sinc kernel with an explicit cutoff frequency. + + Use this when the cutoff must be lower than sampleRate / 2, e.g. for + downsampling where the target Nyquist is the cutoff. + + @param cutoff Anti-aliasing cutoff in Hz (must be in (0, sampleRate/2]). + @param sampleRate Input sample rate in Hz. + */ + void configureWithCutoff (CoeffType cutoff, CoeffType sampleRate) noexcept + { + jassert (cutoff > CoeffType (0) && cutoff <= sampleRate / CoeffType (2)); + + const CoeffType T = CoeffType (1) / sampleRate; + + for (int i = 0; i <= SincRadius; ++i) + { + for (int delta = 0; delta < OversampleFactor; ++delta) + { + const CoeffType t = static_cast (i) + + static_cast (delta) / static_cast (OversampleFactor); + + if (t == CoeffType (0)) + (*this) (i, delta) = CoeffType (1); + else + (*this) (i, delta) = std::sin (MathConstants::twoPi * cutoff * t * T) + / (MathConstants::twoPi * cutoff * t * T); + } + } + } + + /** + Multiplies the stored half-kernel by the second half of a Kaiser window. + + The full symmetric window has 2 * tableSize samples; this method applies + sample indices [tableSize .. 2*tableSize-1] (the descending half) to the + stored positive-half sinc values. + + @param beta Kaiser window shape parameter (higher = more side-lobe suppression). + */ + void applyKaiserWindow (CoeffType beta = CoeffType (5)) noexcept + { + constexpr int N = tableSize * 2; + + for (int i = 0; i < tableSize; ++i) + table[static_cast (i)] *= + WindowFunctions::kaiser (i + tableSize, N, beta); + } + + //============================================================================== + /** Returns the kernel value at absolute (mirrored) index i. */ + forcedinline CoeffType& operator[] (int i) noexcept + { + const int idx = (i < 0) ? -i : i; + return table[static_cast (idx)]; + } + + /** Returns the kernel value at absolute (mirrored) index i. */ + const forcedinline CoeffType& operator[] (int i) const noexcept + { + const int idx = (i < 0) ? -i : i; + return table[static_cast (idx)]; + } + + /** + Returns the kernel value at (tap, delta), where the fractional phase is + t = tap + delta / OversampleFactor. + + Negative tap indices are resolved by mirroring (symmetric kernel). + Negative delta with tap == 0 is also resolved by symmetry. + */ + forcedinline CoeffType& operator() (int tap, int delta) noexcept + { + if (tap < 0) + return table[static_cast ((-tap) * OversampleFactor - delta)]; + if (tap == 0 && delta < 0) + return table[static_cast (-delta)]; // sinc(-t) == sinc(t) + return table[static_cast (tap * OversampleFactor + delta)]; + } + + /** @copydoc operator()(int, int) */ + const forcedinline CoeffType& operator() (int tap, int delta) const noexcept + { + if (tap < 0) + return table[static_cast ((-tap) * OversampleFactor - delta)]; + if (tap == 0 && delta < 0) + return table[static_cast (-delta)]; // sinc(-t) == sinc(t) + return table[static_cast (tap * OversampleFactor + delta)]; + } + +private: + std::array (tableSize)> table {}; +}; + +} // namespace yup diff --git a/modules/yup_dsp/yup_dsp.h b/modules/yup_dsp/yup_dsp.h index 34083bc5b..1ccd00033 100644 --- a/modules/yup_dsp/yup_dsp.h +++ b/modules/yup_dsp/yup_dsp.h @@ -176,3 +176,9 @@ // Time-stretching and pitch-shifting #include "stretching/yup_TimeStretchProcessor.h" + +// Oversampling and sample-rate conversion +#include "resampling/yup_CircularBuffer.h" +#include "resampling/yup_SincTable.h" +#include "resampling/yup_Oversampler.h" +#include "resampling/yup_Resampler.h" diff --git a/tests/yup_dsp.cpp b/tests/yup_dsp.cpp index 0c4e2adc8..8fdfc9627 100644 --- a/tests/yup_dsp.cpp +++ b/tests/yup_dsp.cpp @@ -22,6 +22,7 @@ #include "yup_dsp/yup_BiquadCascade.cpp" #include "yup_dsp/yup_BiquadFilter.cpp" #include "yup_dsp/yup_ButterworthFilter.cpp" +#include "yup_dsp/yup_CircularBuffer.cpp" #include "yup_dsp/yup_DirectFIR.cpp" #include "yup_dsp/yup_FFTProcessor.cpp" #include "yup_dsp/yup_FilterDesigner.cpp" @@ -31,12 +32,16 @@ #include "yup_dsp/yup_LinkwitzRileyFilter.cpp" #include "yup_dsp/yup_LoudnessFilter.cpp" #include "yup_dsp/yup_NoiseGenerators.cpp" +#include "yup_dsp/yup_Oversampler.cpp" #include "yup_dsp/yup_PartitionedConvolver.cpp" #include "yup_dsp/yup_RbjFilter.cpp" +#include "yup_dsp/yup_Resampler.cpp" +#include "yup_dsp/yup_SincTable.cpp" #include "yup_dsp/yup_SoftClipper.cpp" #include "yup_dsp/yup_SpectrumAnalyzerState.cpp" #include "yup_dsp/yup_StateVariableFilter.cpp" #include "yup_dsp/yup_WindowFunctions.cpp" + #if YUP_MODULE_AVAILABLE_bungee_library #include "yup_dsp/yup_TimeStretchProcessor.cpp" #endif diff --git a/tests/yup_dsp/yup_CircularBuffer.cpp b/tests/yup_dsp/yup_CircularBuffer.cpp new file mode 100644 index 000000000..602edd069 --- /dev/null +++ b/tests/yup_dsp/yup_CircularBuffer.cpp @@ -0,0 +1,164 @@ +/* + ============================================================================== + + This file is part of the YUP library. + Copyright (c) 2026 - kunitoki@gmail.com + + YUP is an open source library subject to open-source licensing. + + The code included in this file is provided under the terms of the ISC license + http://www.isc.org/downloads/software-support-policy/isc-license. Permission + to use, copy, modify, and/or distribute this software for any purpose with or + without fee is hereby granted provided that the above copyright notice and + this permission notice appear in all copies. + + YUP IS PROVIDED "AS IS" WITHOUT ANY WARRANTY, AND ALL WARRANTIES, WHETHER + EXPRESSED OR IMPLIED, INCLUDING MERCHANTABILITY AND FITNESS FOR PURPOSE, ARE + DISCLAIMED. + + ============================================================================== +*/ + +#include + +#include + +namespace yup::test +{ + +//============================================================================== +class CircularBufferTest : public ::testing::Test +{ +}; + +//============================================================================== +TEST_F (CircularBufferTest, DefaultConstructionZeroInitialises) +{ + CircularBuffer buf; + for (int i = 0; i < 4; ++i) + EXPECT_EQ (buf[i], 0.0f); +} + +TEST_F (CircularBufferTest, CustomInitValueFillsBuffer) +{ + CircularBuffer buf (42); + for (int i = 0; i < 3; ++i) + EXPECT_EQ (buf[i], 42); +} + +TEST_F (CircularBufferTest, PushAdvancesWritePointerAndOldestIsAtZero) +{ + CircularBuffer buf; + buf.push (1); + buf.push (2); + buf.push (3); + buf.push (4); + + // After 4 pushes into a size-4 buffer, index 0 is the oldest (1) and 3 is newest (4) + EXPECT_EQ (buf[0], 1); + EXPECT_EQ (buf[1], 2); + EXPECT_EQ (buf[2], 3); + EXPECT_EQ (buf[3], 4); +} + +TEST_F (CircularBufferTest, PushWrapsAroundCorrectly) +{ + CircularBuffer buf; + buf.push (10); + buf.push (20); + buf.push (30); + // Buffer full: [10, 20, 30], oldest=10 at index 0 + + buf.push (40); + // Oldest entry (10) overwritten: [20, 30, 40] + EXPECT_EQ (buf[0], 20); + EXPECT_EQ (buf[1], 30); + EXPECT_EQ (buf[2], 40); +} + +TEST_F (CircularBufferTest, PushMaintainsNewestAtHighestIndex) +{ + CircularBuffer buf; + for (int i = 1; i <= 10; ++i) + buf.push (i); + + // After 10 pushes into size-5: newest=10, oldest=6 + EXPECT_EQ (buf[4], 10); + EXPECT_EQ (buf[3], 9); + EXPECT_EQ (buf[2], 8); + EXPECT_EQ (buf[1], 7); + EXPECT_EQ (buf[0], 6); +} + +TEST_F (CircularBufferTest, ClearResetsToZeroAndRewindsPointer) +{ + CircularBuffer buf; + buf.push (1.0f); + buf.push (2.0f); + buf.push (3.0f); + + buf.clear(); + + for (int i = 0; i < 4; ++i) + EXPECT_EQ (buf[i], 0.0f); + + // After clear, push should start from logical index 0 again + buf.push (7.0f); + EXPECT_EQ (buf[0], 7.0f); + for (int i = 1; i < 4; ++i) + EXPECT_EQ (buf[i], 0.0f); +} + +TEST_F (CircularBufferTest, WriteViaSubscriptOperator) +{ + CircularBuffer buf; + buf.push (1); + buf.push (2); + buf.push (3); + buf.push (4); + + buf[2] = 99; + EXPECT_EQ (buf[2], 99); + // Surrounding entries unchanged + EXPECT_EQ (buf[1], 2); + EXPECT_EQ (buf[3], 4); +} + +TEST_F (CircularBufferTest, ConstAccessMatchesNonConst) +{ + CircularBuffer buf; + buf.push (1.5f); + buf.push (2.5f); + buf.push (3.5f); + + const auto& constBuf = buf; + for (int i = 0; i < 3; ++i) + EXPECT_EQ (constBuf[i], buf[i]); +} + +TEST_F (CircularBufferTest, SizeOneFunctionsCorrectly) +{ + CircularBuffer buf; + buf.push (5); + EXPECT_EQ (buf[0], 5); + + buf.push (9); + EXPECT_EQ (buf[0], 9); + + buf.clear(); + EXPECT_EQ (buf[0], 0); +} + +TEST_F (CircularBufferTest, DoublePrecisionType) +{ + CircularBuffer buf; + buf.push (1.1); + buf.push (2.2); + buf.push (3.3); + buf.push (4.4); + + EXPECT_DOUBLE_EQ (buf[0], 1.1); + EXPECT_DOUBLE_EQ (buf[3], 4.4); +} + +} // namespace yup::test diff --git a/tests/yup_dsp/yup_Oversampler.cpp b/tests/yup_dsp/yup_Oversampler.cpp new file mode 100644 index 000000000..dccfddf0e --- /dev/null +++ b/tests/yup_dsp/yup_Oversampler.cpp @@ -0,0 +1,241 @@ +/* + ============================================================================== + + This file is part of the YUP library. + Copyright (c) 2026 - kunitoki@gmail.com + + YUP is an open source library subject to open-source licensing. + + The code included in this file is provided under the terms of the ISC license + http://www.isc.org/downloads/software-support-policy/isc-license. Permission + to use, copy, modify, and/or distribute this software for any purpose with or + without fee is hereby granted provided that the above copyright notice and + this permission notice appear in all copies. + + YUP IS PROVIDED "AS IS" WITHOUT ANY WARRANTY, AND ALL WARRANTIES, WHETHER + EXPRESSED OR IMPLIED, INCLUDING MERCHANTABILITY AND FITNESS FOR PURPOSE, ARE + DISCLAIMED. + + ============================================================================== +*/ + +#include + +#include + +#include +#include + +namespace yup::test +{ + +//============================================================================== +class OversamplerTest : public ::testing::Test +{ +protected: + static constexpr double sampleRate = 44100.0; + static constexpr int maxChannels = 2; + static constexpr int blockSize = 256; + + void SetUp() override + { + os2x.prepare (sampleRate, maxChannels, blockSize); + os4x.prepare (sampleRate, maxChannels, blockSize); + } + + float calculateRMS (const float* data, int numSamples) const + { + float sum = 0.0f; + for (int i = 0; i < numSamples; ++i) + sum += data[i] * data[i]; + return std::sqrt (sum / static_cast (numSamples)); + } + + void fillSine (std::vector& buf, float frequency, float amplitude = 1.0f) const + { + for (std::size_t i = 0; i < buf.size(); ++i) + { + buf[i] = amplitude * std::sin (MathConstants::twoPi * frequency * static_cast (i) / static_cast (sampleRate)); + } + } + + void fillDC (std::vector& buf, float value) const + { + std::fill (buf.begin(), buf.end(), value); + } + + Oversampler os2x; + Oversampler os4x; +}; + +//============================================================================== +TEST_F (OversamplerTest, DefaultConstructionDoesNotCrash) +{ + Oversampler os; + EXPECT_EQ (os.getOversampledNumSamples(), 0); + EXPECT_EQ (os.getLatencyInSamples(), 16); + EXPECT_EQ (os.getOversampledChannelData (0), nullptr); +} + +TEST_F (OversamplerTest, PrepareAllocatesOversampledBuffer) +{ + EXPECT_EQ (os2x.getOversampledNumSamples(), 0); + + std::vector ch0 (blockSize, 0.0f); + const float* inputPtrs[] = { ch0.data() }; + os2x.upsample (inputPtrs, 1, blockSize); + + EXPECT_EQ (os2x.getOversampledNumSamples(), blockSize * 2); + EXPECT_NE (os2x.getOversampledChannelData (0), nullptr); +} + +TEST_F (OversamplerTest, LatencyReturnsCorrectValue) +{ + EXPECT_EQ (os2x.getLatencyInSamples(), 16); // 2 * SincRadius = 2 * 8 + EXPECT_EQ (os4x.getLatencyInSamples(), 16); +} + +TEST_F (OversamplerTest, ResetClearsOversampledSize) +{ + std::vector ch0 (blockSize, 1.0f); + const float* inputPtrs[] = { ch0.data() }; + os2x.upsample (inputPtrs, 1, blockSize); + + ASSERT_EQ (os2x.getOversampledNumSamples(), blockSize * 2); + + os2x.reset(); + EXPECT_EQ (os2x.getOversampledNumSamples(), 0); +} + +TEST_F (OversamplerTest, UpsampleDCSignalHasCorrectMagnitude) +{ + constexpr float dcValue = 0.5f; + std::vector ch0 (blockSize, dcValue); + const float* inputPtrs[] = { ch0.data() }; + + // Warm up the filter (several blocks to flush transient) + for (int b = 0; b < 5; ++b) + os2x.upsample (inputPtrs, 1, blockSize); + + const float* outData = os2x.getOversampledChannelData (0); + ASSERT_NE (outData, nullptr); + + // Second half of the block should be at steady state + const int halfSize = blockSize; // = blockSize * 2 / 2 + float rms = calculateRMS (outData + halfSize, halfSize); + EXPECT_NEAR (rms, dcValue, 0.05f); +} + +TEST_F (OversamplerTest, ProcessOversampledBlockCallbackReceivesCorrectSize) +{ + std::vector ch0 (blockSize, 0.0f), ch1 (blockSize, 0.0f); + const float* inputPtrs[] = { ch0.data(), ch1.data() }; + os2x.upsample (inputPtrs, 2, blockSize); + + int callbackChannels = 0; + int callbackSamples = 0; + os2x.processOversampledBlock ([&] (auto& buf) + { + callbackChannels = static_cast (buf.size()); + if (! buf.empty()) + callbackSamples = static_cast (buf[0].size()); + }); + + EXPECT_EQ (callbackChannels, maxChannels); + EXPECT_EQ (callbackSamples, blockSize * 2); +} + +TEST_F (OversamplerTest, UpsampleThenDownsamplePreservesLowFrequencySine) +{ + constexpr float frequency = 440.0f; // A4 - well below Nyquist/4 + std::vector input (blockSize); + fillSine (input, frequency); + + std::vector output (blockSize, 0.0f); + const float* inPtrs[] = { input.data() }; + float* outPtrs[] = { output.data() }; + + // Warm up to flush latency + for (int b = 0; b < 10; ++b) + { + os2x.upsample (inPtrs, 1, blockSize); + os2x.downsample (outPtrs, 1, blockSize); + } + + // Compare RMS energy - should be preserved + float rmsIn = calculateRMS (input.data(), blockSize); + float rmsOut = calculateRMS (output.data(), blockSize); + EXPECT_NEAR (rmsOut, rmsIn, rmsIn * 0.1f); // within 10% +} + +TEST_F (OversamplerTest, DecimationFiltersOversampledDomainHighFrequency) +{ + // Upsample silence to get a clean oversampled buffer + std::vector silence (blockSize, 0.0f); + const float* silPtrs[] = { silence.data() }; + + for (int b = 0; b < 5; ++b) + os2x.upsample (silPtrs, 1, blockSize); + + // Inject a tone ABOVE original Nyquist directly into the oversampled buffer. + // This simulates distortion harmonics created at the elevated sample rate. + // The decimation filter should attenuate this before downsampling. + const double oversampledRate = sampleRate * 2.0; + const double highFreq = sampleRate * 0.6; // above original Nyquist (22050 Hz) + float* oversampledData = os2x.getOversampledChannelData (0); + const int oversampledLen = os2x.getOversampledNumSamples(); + constexpr float injectedAmplitude = 0.5f; + + for (int i = 0; i < oversampledLen; ++i) + oversampledData[i] += injectedAmplitude * static_cast (std::sin (MathConstants::twoPi * highFreq * static_cast (i) / oversampledRate)); + + std::vector output (blockSize, 0.0f); + float* outPtrs[] = { output.data() }; + os2x.downsample (outPtrs, 1, blockSize); + + float rmsOut = calculateRMS (output.data(), blockSize); + + // The anti-aliasing filter should substantially reduce the injected tone + EXPECT_LT (rmsOut, injectedAmplitude * 0.5f); +} + +TEST_F (OversamplerTest, OversampledChannelDataNotNullAfterUpsample) +{ + std::vector ch0 (blockSize, 0.0f); + const float* inputPtrs[] = { ch0.data() }; + os2x.upsample (inputPtrs, 1, blockSize); + + EXPECT_NE (os2x.getOversampledChannelData (0), nullptr); + EXPECT_EQ (os2x.getOversampledChannelData (1), nullptr); // channel 1 not prepared + EXPECT_EQ (os2x.getOversampledChannelData (-1), nullptr); // invalid index +} + +TEST_F (OversamplerTest, FourXOversamplerHasCorrectOutputSize) +{ + std::vector ch0 (blockSize, 0.0f); + const float* inputPtrs[] = { ch0.data() }; + os4x.upsample (inputPtrs, 1, blockSize); + + EXPECT_EQ (os4x.getOversampledNumSamples(), blockSize * 4); +} + +//============================================================================== +TEST (OversamplerTypeAliasTest, TypeAliasesCompile) +{ + Oversampler2xFloat a; + Oversampler4xFloat b; + Oversampler8xFloat c; + Oversampler2xDouble d; + Oversampler4xDouble e; + + // Prepare briefly to confirm the types are usable + a.prepare (44100.0, 1, 64); + b.prepare (44100.0, 1, 64); + c.prepare (44100.0, 1, 64); + d.prepare (44100.0, 1, 64); + e.prepare (44100.0, 1, 64); + + SUCCEED(); +} + +} // namespace yup::test diff --git a/tests/yup_dsp/yup_Resampler.cpp b/tests/yup_dsp/yup_Resampler.cpp new file mode 100644 index 000000000..b314aa9b0 --- /dev/null +++ b/tests/yup_dsp/yup_Resampler.cpp @@ -0,0 +1,251 @@ +/* + ============================================================================== + + This file is part of the YUP library. + Copyright (c) 2026 - kunitoki@gmail.com + + YUP is an open source library subject to open-source licensing. + + The code included in this file is provided under the terms of the ISC license + http://www.isc.org/downloads/software-support-policy/isc-license. Permission + to use, copy, modify, and/or distribute this software for any purpose with or + without fee is hereby granted provided that the above copyright notice and + this permission notice appear in all copies. + + YUP IS PROVIDED "AS IS" WITHOUT ANY WARRANTY, AND ALL WARRANTIES, WHETHER + EXPRESSED OR IMPLIED, INCLUDING MERCHANTABILITY AND FITNESS FOR PURPOSE, ARE + DISCLAIMED. + + ============================================================================== +*/ + +#include + +#include + +#include +#include + +namespace yup::test +{ + +//============================================================================== +class ResamplerTest : public ::testing::Test +{ +protected: + static constexpr double sourceSampleRate = 44100.0; + static constexpr double targetSampleRate = 48000.0; + static constexpr int maxChannels = 2; + static constexpr int blockSize = 256; + + void SetUp() override + { + resamplerUp.prepare (sourceSampleRate, targetSampleRate, maxChannels, blockSize); + resamplerDown.prepare (targetSampleRate, sourceSampleRate, maxChannels, blockSize); + } + + float calculateRMS (const float* data, int numSamples) const + { + float sum = 0.0f; + for (int i = 0; i < numSamples; ++i) + sum += data[i] * data[i]; + return std::sqrt (sum / static_cast (numSamples)); + } + + void fillSine (std::vector& buf, double frequency, double sr, float amplitude = 1.0f) const + { + for (std::size_t i = 0; i < buf.size(); ++i) + { + buf[i] = amplitude * static_cast (std::sin (MathConstants::twoPi * frequency * static_cast (i) / sr)); + } + } + + void fillDC (std::vector& buf, float value) const + { + std::fill (buf.begin(), buf.end(), value); + } + + // Max output buffer size for a 44100 -> 48000 conversion + static constexpr int maxOutputSize = static_cast (blockSize * 48000.0 / 44100.0) + 2; + + Resampler resamplerUp; // 44100 -> 48000 + Resampler resamplerDown; // 48000 -> 44100 +}; + +//============================================================================== +TEST_F (ResamplerTest, DefaultConstructionDoesNotCrash) +{ + Resampler r; + EXPECT_EQ (r.getLatencyInSamples(), 8); +} + +TEST_F (ResamplerTest, LatencyReturnsCorrectValue) +{ + EXPECT_EQ (resamplerUp.getLatencyInSamples(), 8); // SincRadius = 8 + EXPECT_EQ (resamplerDown.getLatencyInSamples(), 8); +} + +TEST_F (ResamplerTest, UpsampleProducesMoreOutputThanInput) +{ + std::vector input (blockSize, 0.0f); + std::vector output (maxOutputSize, 0.0f); + const float* inPtrs[] = { input.data() }; + float* outPtrs[] = { output.data() }; + + int produced = resamplerUp.resample (inPtrs, outPtrs, 1, blockSize); + + // 44100 -> 48000: expect approx blockSize * 48000/44100 output samples + const int expectedApprox = static_cast (blockSize * targetSampleRate / sourceSampleRate); + EXPECT_GT (produced, 0); + EXPECT_NEAR (produced, expectedApprox, 2); +} + +TEST_F (ResamplerTest, DownsampleProducesFewerOutputThanInput) +{ + std::vector input (blockSize, 0.0f); + std::vector output (blockSize, 0.0f); + const float* inPtrs[] = { input.data() }; + float* outPtrs[] = { output.data() }; + + int produced = resamplerDown.resample (inPtrs, outPtrs, 1, blockSize); + + // 48000 -> 44100: expect fewer output samples than input + const int expectedApprox = static_cast (blockSize * sourceSampleRate / targetSampleRate); + EXPECT_GT (produced, 0); + EXPECT_NEAR (produced, expectedApprox, 2); +} + +TEST_F (ResamplerTest, ResampleDCPreservesMagnitude) +{ + constexpr float dcValue = 0.5f; + std::vector input (blockSize); + fillDC (input, dcValue); + + std::vector output (maxOutputSize, 0.0f); + const float* inPtrs[] = { input.data() }; + float* outPtrs[] = { output.data() }; + + // Warm up to flush transient + for (int b = 0; b < 10; ++b) + resamplerUp.resample (inPtrs, outPtrs, 1, blockSize); + + int produced = resamplerUp.resample (inPtrs, outPtrs, 1, blockSize); + ASSERT_GT (produced, 0); + + float rms = calculateRMS (output.data(), produced); + EXPECT_NEAR (rms, dcValue, 0.05f); +} + +TEST_F (ResamplerTest, ResampleSinePreservesFrequencyContent) +{ + constexpr double frequency = 440.0; // A4 + std::vector input (blockSize); + fillSine (input, frequency, sourceSampleRate); + + std::vector output (maxOutputSize, 0.0f); + const float* inPtrs[] = { input.data() }; + float* outPtrs[] = { output.data() }; + + // Warm up + for (int b = 0; b < 10; ++b) + resamplerUp.resample (inPtrs, outPtrs, 1, blockSize); + + int produced = resamplerUp.resample (inPtrs, outPtrs, 1, blockSize); + ASSERT_GT (produced, 0); + + float rmsIn = calculateRMS (input.data(), blockSize); + float rmsOut = calculateRMS (output.data(), produced); + + EXPECT_NEAR (rmsOut, rmsIn, rmsIn * 0.1f); // within 10% +} + +TEST_F (ResamplerTest, ResetRestoresCleanState) +{ + std::vector input (blockSize); + fillSine (input, 440.0, sourceSampleRate); + + std::vector output1 (maxOutputSize, 0.0f); + std::vector output2 (maxOutputSize, 0.0f); + const float* inPtrs[] = { input.data() }; + float* outPtrs1[] = { output1.data() }; + float* outPtrs2[] = { output2.data() }; + + // Run a few blocks, then reset, then run again from beginning + int produced1 = resamplerUp.resample (inPtrs, outPtrs1, 1, blockSize); + + resamplerUp.reset(); + + int produced2 = resamplerUp.resample (inPtrs, outPtrs2, 1, blockSize); + + // After reset, first block output should match the very first run + EXPECT_EQ (produced1, produced2); + EXPECT_NEAR (output1[0], output2[0], 1e-4f); +} + +TEST_F (ResamplerTest, TwoChannelResamplingProducesConsistentOutput) +{ + std::vector ch0 (blockSize), ch1 (blockSize); + fillSine (ch0, 440.0, sourceSampleRate, 0.7f); + fillSine (ch1, 880.0, sourceSampleRate, 0.5f); + + std::vector out0 (maxOutputSize, 0.0f); + std::vector out1 (maxOutputSize, 0.0f); + const float* inPtrs[] = { ch0.data(), ch1.data() }; + float* outPtrs[] = { out0.data(), out1.data() }; + + // Warm up + for (int b = 0; b < 5; ++b) + resamplerUp.resample (inPtrs, outPtrs, 2, blockSize); + + int produced = resamplerUp.resample (inPtrs, outPtrs, 2, blockSize); + ASSERT_GT (produced, 0); + + float rms0 = calculateRMS (out0.data(), produced); + float rms1 = calculateRMS (out1.data(), produced); + + // Both channels should have non-trivial output + EXPECT_GT (rms0, 0.1f); + EXPECT_GT (rms1, 0.1f); + + // RMS should be roughly proportional to input amplitude + EXPECT_GT (rms0, rms1); +} + +TEST_F (ResamplerTest, IdentityResamplePreservesSignal) +{ + // Same input and output rate should be a near-identity operation + Resampler identity; + identity.prepare (44100.0, 44100.0, 1, blockSize); + + std::vector input (blockSize); + fillSine (input, 440.0, 44100.0); + + std::vector output (blockSize + 4, 0.0f); + const float* inPtrs[] = { input.data() }; + float* outPtrs[] = { output.data() }; + + // Warm up + for (int b = 0; b < 10; ++b) + identity.resample (inPtrs, outPtrs, 1, blockSize); + + int produced = identity.resample (inPtrs, outPtrs, 1, blockSize); + ASSERT_GT (produced, 0); + + float rmsIn = calculateRMS (input.data(), blockSize); + float rmsOut = calculateRMS (output.data(), produced); + EXPECT_NEAR (rmsOut, rmsIn, rmsIn * 0.05f); // within 5% +} + +//============================================================================== +TEST (ResamplerTypeAliasTest, TypeAliasesCompile) +{ + ResamplerFloat a; + ResamplerDouble b; + + a.prepare (44100.0, 48000.0, 1, 64); + b.prepare (44100.0, 48000.0, 1, 64); + + SUCCEED(); +} + +} // namespace yup::test diff --git a/tests/yup_dsp/yup_SincTable.cpp b/tests/yup_dsp/yup_SincTable.cpp new file mode 100644 index 000000000..a11ec8a09 --- /dev/null +++ b/tests/yup_dsp/yup_SincTable.cpp @@ -0,0 +1,176 @@ +/* + ============================================================================== + + This file is part of the YUP library. + Copyright (c) 2026 - kunitoki@gmail.com + + YUP is an open source library subject to open-source licensing. + + The code included in this file is provided under the terms of the ISC license + http://www.isc.org/downloads/software-support-policy/isc-license. Permission + to use, copy, modify, and/or distribute this software for any purpose with or + without fee is hereby granted provided that the above copyright notice and + this permission notice appear in all copies. + + YUP IS PROVIDED "AS IS" WITHOUT ANY WARRANTY, AND ALL WARRANTIES, WHETHER + EXPRESSED OR IMPLIED, INCLUDING MERCHANTABILITY AND FITNESS FOR PURPOSE, ARE + DISCLAIMED. + + ============================================================================== +*/ + +#include + +#include + +#include + +namespace yup::test +{ + +//============================================================================== +class SincTableTest : public ::testing::Test +{ +protected: + static constexpr double sampleRate = 44100.0; + static constexpr int radius = 8; + static constexpr int factor = 4; // oversampling factor / resolution + static constexpr double tolerance = 1e-9; +}; + +//============================================================================== +TEST_F (SincTableTest, TableSizeIsCorrect) +{ + EXPECT_EQ ((SincTable::tableSize), (radius + 1) * factor); +} + +TEST_F (SincTableTest, ConfigureSetsUnityAtZero) +{ + SincTable table; + table.configure (static_cast (sampleRate)); + + // (tap=0, delta=0) is the kernel centre — sinc(0) = 1 + EXPECT_NEAR (table (0, 0), 1.0, tolerance); +} + +TEST_F (SincTableTest, ConfigureWithCutoffSetsUnityAtZero) +{ + SincTable table; + table.configureWithCutoff (sampleRate / 2.0, sampleRate); + + EXPECT_NEAR (table (0, 0), 1.0, tolerance); +} + +TEST_F (SincTableTest, KernelDecaysAwayFromCentre) +{ + SincTable table; + table.configure (static_cast (sampleRate)); + + // Magnitude at tap=1, delta=0 must be strictly less than the centre value + EXPECT_LT (std::abs (table (1, 0)), std::abs (table (0, 0))); + EXPECT_LT (std::abs (table (2, 0)), std::abs (table (1, 0))); +} + +TEST_F (SincTableTest, SymmetryViaOperatorBracket) +{ + SincTable table; + table.configure (static_cast (sampleRate)); + + // operator[] mirrors negative indices: table[-i] == table[i] + for (int i = 1; i <= radius * factor; ++i) + EXPECT_DOUBLE_EQ (table[-i], table[i]); +} + +TEST_F (SincTableTest, SymmetryViaTapDelta) +{ + SincTable table; + table.configure (static_cast (sampleRate)); + + // Negative tap mirrors positive tap: table(-n, d) == table(n, -d) (by sinc symmetry) + for (int n = 1; n <= radius; ++n) + { + for (int d = 0; d < factor; ++d) + { + // table(-n, d) should equal table(n, -d... handled via tap=0, delta<0 branch at boundary) + // Direct: negative tap accesses (-tap)*factor - delta = n*factor - d + // That equals table[n*factor - d] which == table(n, -d) only for d==0 or d mirrored + // Simplest check: table(-n, 0) == table(n, 0) + if (d == 0) + EXPECT_DOUBLE_EQ (table (-n, 0), table (n, 0)); + } + } +} + +TEST_F (SincTableTest, NegativeDeltaAtTapZeroUsesSymmetry) +{ + SincTable table; + table.configure (static_cast (sampleRate)); + + // table(0, -d) must equal table(0, d) — sinc(-t) == sinc(t) + for (int d = 1; d < factor; ++d) + EXPECT_DOUBLE_EQ (table (0, -d), table (0, d)); +} + +TEST_F (SincTableTest, ApplyKaiserWindowReducesOuterTaps) +{ + SincTable table; + table.configure (static_cast (sampleRate)); + + // Record value at a tap far from centre before windowing + const double beforeWindow = std::abs (table (radius, 0)); + + table.applyKaiserWindow (5.0); + + const double afterWindow = std::abs (table (radius, 0)); + + // Kaiser window tapers to zero at the edges, so outer taps must shrink + EXPECT_LT (afterWindow, beforeWindow); +} + +TEST_F (SincTableTest, ApplyKaiserWindowPreservesCentre) +{ + SincTable table; + table.configure (static_cast (sampleRate)); + + // Centre before windowing + const double before = table (0, 0); + + table.applyKaiserWindow (5.0); + + // The Kaiser window at the centre (index tableSize in the full 2*tableSize window) + // is close to 1 but not exactly 1 — just verify it is still close + EXPECT_NEAR (table (0, 0), before, 0.01); +} + +TEST_F (SincTableTest, ConfigureWithLowerCutoffReducesBandwidth) +{ + SincTable highCut; + highCut.configure (static_cast (sampleRate)); // cutoff = sr/2 + + SincTable lowCut; + lowCut.configureWithCutoff (sampleRate / 4.0, sampleRate); // cutoff = sr/4 + + // Lower cutoff → longer sinc zeros → tap-1 magnitude should be smaller for low-cut + EXPECT_LT (std::abs (lowCut (1, 0)), std::abs (highCut (1, 0))); +} + +TEST_F (SincTableTest, FloatPrecisionInstantiates) +{ + SincTable table; + table.configure (44100.0f); + table.applyKaiserWindow (5.0f); + + EXPECT_NEAR (table (0, 0), 1.0f, 1e-4f); +} + +TEST_F (SincTableTest, HighResolutionInstantiates) +{ + SincTable table; + table.configure (static_cast (sampleRate)); + table.applyKaiserWindow (5.0); + + EXPECT_EQ ((SincTable::tableSize), 9 * 256); + EXPECT_NEAR (table (0, 0), 1.0, tolerance); +} + +} // namespace yup::test From 403d257a734a20e530a89406897e999a223020dd Mon Sep 17 00:00:00 2001 From: kunitoki Date: Thu, 21 May 2026 23:31:05 +0200 Subject: [PATCH 2/3] Fixes --- modules/yup_dsp/resampling/yup_Oversampler.h | 17 +++++++++---- modules/yup_dsp/resampling/yup_SincTable.h | 11 +++++---- tests/yup_dsp/yup_CircularBuffer.cpp | 8 +++---- tests/yup_dsp/yup_SincTable.cpp | 25 ++++++++++---------- 4 files changed, 35 insertions(+), 26 deletions(-) diff --git a/modules/yup_dsp/resampling/yup_Oversampler.h b/modules/yup_dsp/resampling/yup_Oversampler.h index 824d8f178..59fd89f90 100644 --- a/modules/yup_dsp/resampling/yup_Oversampler.h +++ b/modules/yup_dsp/resampling/yup_Oversampler.h @@ -96,6 +96,7 @@ class Oversampler xDecim.assign (maxChannels, std::vector (static_cast (maxInterpolated + SincRadius * OversampleFactor), SampleType {})); oversampledBuffer.assign (maxChannels, std::vector (static_cast (maxInterpolated), SampleType {})); currentOversampledSize = 0; + currentNumChannels = 0; } /** @@ -123,6 +124,7 @@ class Oversampler for (auto& ch : oversampledBuffer) std::fill (ch.begin(), ch.end(), SampleType {}); currentOversampledSize = 0; + currentNumChannels = 0; } //============================================================================== @@ -157,6 +159,7 @@ class Oversampler } currentOversampledSize = numSamples * OversampleFactor; + currentNumChannels = numChannels; for (int ch = 0; ch < numChannels; ++ch) { @@ -288,12 +291,13 @@ class Oversampler @param channel Zero-based channel index. @return Pointer to getOversampledNumSamples() contiguous samples, - or nullptr if the channel index is out of range or - prepare() has not been called. + or nullptr if the channel index is out of range, prepare() + has not been called, or the channel was not processed by + the most recent upsample() call. */ forcedinline SampleType* getOversampledChannelData (int channel) noexcept { - if (channel < 0 || channel >= static_cast (oversampledBuffer.size())) + if (channel < 0 || channel >= currentNumChannels) return nullptr; return oversampledBuffer[static_cast (channel)].data(); @@ -304,11 +308,13 @@ class Oversampler @param channel Zero-based channel index. @return Pointer to getOversampledNumSamples() contiguous samples, - or nullptr if the channel index is out of range. + or nullptr if the channel index is out of range or the + channel was not processed by the most recent upsample() + call. */ const forcedinline SampleType* getOversampledChannelData (int channel) const noexcept { - if (channel < 0 || channel >= static_cast (oversampledBuffer.size())) + if (channel < 0 || channel >= currentNumChannels) return nullptr; return oversampledBuffer[static_cast (channel)].data(); @@ -349,6 +355,7 @@ class Oversampler std::vector> xDecim; std::vector> oversampledBuffer; int currentOversampledSize = 0; + int currentNumChannels = 0; YUP_DECLARE_NON_COPYABLE_WITH_LEAK_DETECTOR (Oversampler) }; diff --git a/modules/yup_dsp/resampling/yup_SincTable.h b/modules/yup_dsp/resampling/yup_SincTable.h index b2f06dafb..dc4046c98 100644 --- a/modules/yup_dsp/resampling/yup_SincTable.h +++ b/modules/yup_dsp/resampling/yup_SincTable.h @@ -115,19 +115,20 @@ class SincTable /** Multiplies the stored half-kernel by the second half of a Kaiser window. - The full symmetric window has 2 * tableSize samples; this method applies - sample indices [tableSize .. 2*tableSize-1] (the descending half) to the - stored positive-half sinc values. + The full symmetric window has 2 * tableSize - 1 samples, so the stored + center coefficient is exactly aligned with the window center and remains + unchanged. @param beta Kaiser window shape parameter (higher = more side-lobe suppression). */ void applyKaiserWindow (CoeffType beta = CoeffType (5)) noexcept { - constexpr int N = tableSize * 2; + constexpr int N = tableSize * 2 - 1; + constexpr int center = tableSize - 1; for (int i = 0; i < tableSize; ++i) table[static_cast (i)] *= - WindowFunctions::kaiser (i + tableSize, N, beta); + WindowFunctions::kaiser (center + i, N, beta); } //============================================================================== diff --git a/tests/yup_dsp/yup_CircularBuffer.cpp b/tests/yup_dsp/yup_CircularBuffer.cpp index 602edd069..b4a817778 100644 --- a/tests/yup_dsp/yup_CircularBuffer.cpp +++ b/tests/yup_dsp/yup_CircularBuffer.cpp @@ -90,7 +90,7 @@ TEST_F (CircularBufferTest, PushMaintainsNewestAtHighestIndex) EXPECT_EQ (buf[0], 6); } -TEST_F (CircularBufferTest, ClearResetsToZeroAndRewindsPointer) +TEST_F (CircularBufferTest, ClearResetsToZeroAndPreservesHistoryOrdering) { CircularBuffer buf; buf.push (1.0f); @@ -102,11 +102,11 @@ TEST_F (CircularBufferTest, ClearResetsToZeroAndRewindsPointer) for (int i = 0; i < 4; ++i) EXPECT_EQ (buf[i], 0.0f); - // After clear, push should start from logical index 0 again + // After clear, a new push becomes the newest logical entry. buf.push (7.0f); - EXPECT_EQ (buf[0], 7.0f); - for (int i = 1; i < 4; ++i) + for (int i = 0; i < 3; ++i) EXPECT_EQ (buf[i], 0.0f); + EXPECT_EQ (buf[3], 7.0f); } TEST_F (CircularBufferTest, WriteViaSubscriptOperator) diff --git a/tests/yup_dsp/yup_SincTable.cpp b/tests/yup_dsp/yup_SincTable.cpp index a11ec8a09..550bc6563 100644 --- a/tests/yup_dsp/yup_SincTable.cpp +++ b/tests/yup_dsp/yup_SincTable.cpp @@ -61,14 +61,15 @@ TEST_F (SincTableTest, ConfigureWithCutoffSetsUnityAtZero) EXPECT_NEAR (table (0, 0), 1.0, tolerance); } -TEST_F (SincTableTest, KernelDecaysAwayFromCentre) +TEST_F (SincTableTest, KernelDecaysAwayFromCenter) { SincTable table; table.configure (static_cast (sampleRate)); - // Magnitude at tap=1, delta=0 must be strictly less than the centre value - EXPECT_LT (std::abs (table (1, 0)), std::abs (table (0, 0))); - EXPECT_LT (std::abs (table (2, 0)), std::abs (table (1, 0))); + // Use a fractional phase so the full-band sinc is not sampled exactly at + // its integer zero crossings. + EXPECT_LT (std::abs (table (1, 1)), std::abs (table (0, 1))); + EXPECT_LT (std::abs (table (2, 1)), std::abs (table (1, 1))); } TEST_F (SincTableTest, SymmetryViaOperatorBracket) @@ -127,22 +128,21 @@ TEST_F (SincTableTest, ApplyKaiserWindowReducesOuterTaps) EXPECT_LT (afterWindow, beforeWindow); } -TEST_F (SincTableTest, ApplyKaiserWindowPreservesCentre) +TEST_F (SincTableTest, ApplyKaiserWindowPreservesCenter) { SincTable table; table.configure (static_cast (sampleRate)); - // Centre before windowing + // Center before windowing const double before = table (0, 0); table.applyKaiserWindow (5.0); - // The Kaiser window at the centre (index tableSize in the full 2*tableSize window) - // is close to 1 but not exactly 1 — just verify it is still close - EXPECT_NEAR (table (0, 0), before, 0.01); + // The Kaiser window is exactly 1 at the center. + EXPECT_DOUBLE_EQ (table (0, 0), before); } -TEST_F (SincTableTest, ConfigureWithLowerCutoffReducesBandwidth) +TEST_F (SincTableTest, ConfigureWithLowerCutoffWidensMainLobe) { SincTable highCut; highCut.configure (static_cast (sampleRate)); // cutoff = sr/2 @@ -150,8 +150,9 @@ TEST_F (SincTableTest, ConfigureWithLowerCutoffReducesBandwidth) SincTable lowCut; lowCut.configureWithCutoff (sampleRate / 4.0, sampleRate); // cutoff = sr/4 - // Lower cutoff → longer sinc zeros → tap-1 magnitude should be smaller for low-cut - EXPECT_LT (std::abs (lowCut (1, 0)), std::abs (highCut (1, 0))); + // Lower cutoff -> wider time-domain main lobe, so near-center fractional + // taps have a larger magnitude. + EXPECT_GT (std::abs (lowCut (1, 1)), std::abs (highCut (1, 1))); } TEST_F (SincTableTest, FloatPrecisionInstantiates) From 799286e21257abe23ace900ab3cc78447f98412d Mon Sep 17 00:00:00 2001 From: kunitoki Date: Fri, 22 May 2026 11:05:12 +0200 Subject: [PATCH 3/3] Added Oversampling and Resampling to the Spectrum analyzer --- .../source/examples/SpectrumAnalyzer.h | 311 ++++++++++++++++-- 1 file changed, 278 insertions(+), 33 deletions(-) diff --git a/examples/graphics/source/examples/SpectrumAnalyzer.h b/examples/graphics/source/examples/SpectrumAnalyzer.h index 36e663f0d..aa53dbc89 100644 --- a/examples/graphics/source/examples/SpectrumAnalyzer.h +++ b/examples/graphics/source/examples/SpectrumAnalyzer.h @@ -25,6 +25,10 @@ #include #include +#include +#include +#include + //============================================================================== class SignalGenerator @@ -39,12 +43,23 @@ class SignalGenerator brownNoise }; + enum class SweepPlaybackMode + { + direct, + resampled, + resampledOversampled2x, + resampledOversampled4x, + resampledOversampled8x + }; + SignalGenerator() : sampleRate (44100.0) + , sourceSampleRate (sampleRate * resamplerSourceRateMultiplier) , frequency (440.0) , phase (0.0) , amplitude (0.5f) - , signalType (SignalType::singleTone) + , signalType (SignalType::frequencySweep) + , sweepPlaybackMode (SweepPlaybackMode::resampled) , sweepStartFreq (20.0) , sweepEndFreq (22000.0) , sweepDurationSeconds (10.0) @@ -54,6 +69,8 @@ class SignalGenerator , smoothedFrequency (440.0) , smoothedAmplitude (0.5f) { + initialiseSineTable(); + // Initialize pink noise filter state for (int i = 0; i < 7; ++i) pinkFilters[i] = 0.0; @@ -61,11 +78,20 @@ class SignalGenerator // Set default smoothing time (50ms) smoothedFrequency.reset (sampleRate, 0.05); smoothedAmplitude.reset (sampleRate, 0.05); + + prepareResampling (defaultMaxBlockSize); + } + + void prepare (double newSampleRate, int maxBlockSize) + { + setSampleRate (newSampleRate); + prepareResampling (maxBlockSize); } void setSampleRate (double newSampleRate) { - sampleRate = newSampleRate; + sampleRate = newSampleRate > 0.0 ? newSampleRate : 44100.0; + sourceSampleRate = sampleRate * resamplerSourceRateMultiplier; updatePhaseIncrement(); // Update smoothing times with new sample rate @@ -90,15 +116,24 @@ class SignalGenerator { signalType = type; if (type == SignalType::frequencySweep) - sweepProgress = 0.0; + resetSweepPlaybackState(); + } + + void setSweepPlaybackMode (SweepPlaybackMode mode) + { + if (sweepPlaybackMode == mode) + return; + + sweepPlaybackMode = mode; + resetSweepPlaybackState(); } void setSweepParameters (double startFreq, double endFreq, double durationSeconds) { sweepStartFreq = startFreq; sweepEndFreq = endFreq; - sweepDurationSeconds = durationSeconds; - sweepProgress = 0.0; + sweepDurationSeconds = yup::jmax (0.001, durationSeconds); + resetSweepPlaybackState(); } void setSmoothingTime (float timeInSeconds) @@ -107,6 +142,21 @@ class SignalGenerator smoothedAmplitude.reset (sampleRate, timeInSeconds); } + void renderNextBlock (float* output, int numSamples) + { + if (output == nullptr || numSamples <= 0) + return; + + if (signalType == SignalType::frequencySweep && sweepPlaybackMode != SweepPlaybackMode::direct) + { + renderResampledSweepBlock (output, numSamples); + return; + } + + for (int sample = 0; sample < numSamples; ++sample) + output[sample] = getNextSample(); + } + float getNextSample() { // Get smoothed parameter values for this sample @@ -138,34 +188,181 @@ class SignalGenerator } private: - float generateSine (double freq) + static constexpr int wavetableSize = 2048; + static constexpr int defaultMaxBlockSize = 512; + static constexpr double resamplerSourceRateMultiplier = 2.0; + + void initialiseSineTable() + { + for (int i = 0; i <= wavetableSize; ++i) + { + const double phaseInRadians = yup::MathConstants::twoPi * static_cast (i) / static_cast (wavetableSize); + sineTable[static_cast (i)] = static_cast (std::sin (phaseInRadians)); + } + } + + void prepareResampling (int maxBlockSize) + { + maxOutputBlockSize = yup::jmax (1, maxBlockSize); + sourceBlockCapacity = maxOutputBlockSize * static_cast (resamplerSourceRateMultiplier) + 64; + resampledBlockCapacity = maxOutputBlockSize + 64; + + sourceBuffer.assign (static_cast (sourceBlockCapacity), 0.0f); + silenceBuffer.assign (static_cast (sourceBlockCapacity), 0.0f); + resampledBuffer.assign (static_cast (resampledBlockCapacity), 0.0f); + + resampler.prepare (sourceSampleRate, sampleRate, 1, sourceBlockCapacity); + oversampler2x.prepare (sourceSampleRate, 1, sourceBlockCapacity); + oversampler4x.prepare (sourceSampleRate, 1, sourceBlockCapacity); + oversampler8x.prepare (sourceSampleRate, 1, sourceBlockCapacity); + + resetSweepPlaybackState(); + } + + void ensurePreparedForBlock (int numSamples) + { + if (numSamples > maxOutputBlockSize) + prepareResampling (numSamples); + } + + void resetSweepPlaybackState() + { + phase = 0.0; + sweepProgress = 0.0; + pendingOutputSamples = 0; + pendingReadPosition = 0; + + resampler.reset(); + oversampler2x.reset(); + oversampler4x.reset(); + oversampler8x.reset(); + } + + void renderResampledSweepBlock (float* output, int numSamples) { - // Calculate phase increment for the smoothed frequency - double currentPhaseIncrement = yup::MathConstants::twoPi * freq / sampleRate; + ensurePreparedForBlock (numSamples); + + int outputPosition = 0; + + while (outputPosition < numSamples) + { + if (pendingReadPosition >= pendingOutputSamples) + { + const int remainingOutputSamples = numSamples - outputPosition; + const int sourceSamplesNeeded = yup::jmin (sourceBlockCapacity, + yup::jmax (1, remainingOutputSamples * static_cast (resamplerSourceRateMultiplier) + 32)); + + renderSourceSweepBlock (sourceSamplesNeeded); + + const float* inputPtrs[] = { sourceBuffer.data() }; + float* outputPtrs[] = { resampledBuffer.data() }; + + pendingOutputSamples = resampler.resample (inputPtrs, outputPtrs, 1, sourceSamplesNeeded); + pendingReadPosition = 0; - float sample = std::sin (phase); - phase += currentPhaseIncrement; + if (pendingOutputSamples <= 0) + { + while (outputPosition < numSamples) + output[outputPosition++] = 0.0f; - if (phase >= yup::MathConstants::twoPi) - phase -= yup::MathConstants::twoPi; + return; + } + } + + const int samplesToCopy = yup::jmin (numSamples - outputPosition, pendingOutputSamples - pendingReadPosition); + + for (int i = 0; i < samplesToCopy; ++i) + { + const float currentAmp = smoothedAmplitude.getNextValue(); + output[outputPosition++] = resampledBuffer[static_cast (pendingReadPosition++)] * currentAmp; + } + } + } + + void renderSourceSweepBlock (int numSamples) + { + switch (sweepPlaybackMode) + { + case SweepPlaybackMode::direct: + case SweepPlaybackMode::resampled: + renderWavetableSweepBlock (sourceBuffer.data(), numSamples, sourceSampleRate); + break; + + case SweepPlaybackMode::resampledOversampled2x: + renderOversampledSourceBlock (oversampler2x, 2, numSamples); + break; + + case SweepPlaybackMode::resampledOversampled4x: + renderOversampledSourceBlock (oversampler4x, 4, numSamples); + break; + + case SweepPlaybackMode::resampledOversampled8x: + renderOversampledSourceBlock (oversampler8x, 8, numSamples); + break; + } + } + + template + void renderOversampledSourceBlock (OversamplerType& oversampler, int oversampleFactor, int numSamples) + { + const float* inputPtrs[] = { silenceBuffer.data() }; + oversampler.upsample (inputPtrs, 1, numSamples); + + auto* oversampledData = oversampler.getOversampledChannelData (0); + renderWavetableSweepBlock (oversampledData, oversampler.getOversampledNumSamples(), sourceSampleRate * static_cast (oversampleFactor)); + + float* outputPtrs[] = { sourceBuffer.data() }; + oversampler.downsample (outputPtrs, 1, numSamples); + } + + void renderWavetableSweepBlock (float* output, int numSamples, double generationSampleRate) + { + for (int i = 0; i < numSamples; ++i) + output[i] = generateSweepAtRate (generationSampleRate); + } + + float readSineTable() const + { + const double tablePosition = phase * static_cast (wavetableSize); + const int index = static_cast (tablePosition); + const float fraction = static_cast (tablePosition - static_cast (index)); + + const float sample0 = sineTable[static_cast (index)]; + const float sample1 = sineTable[static_cast (index + 1)]; + + return sample0 + (sample1 - sample0) * fraction; + } + + void advancePhase (double freq, double generationSampleRate) + { + phase += freq / generationSampleRate; + + while (phase >= 1.0) + phase -= 1.0; + } + + float generateSine (double freq) + { + float sample = readSineTable(); + advancePhase (freq, sampleRate); return sample; } float generateSweep() + { + return generateSweepAtRate (sampleRate); + } + + float generateSweepAtRate (double generationSampleRate) { // Linear frequency sweep double currentFreq = sweepStartFreq + (sweepEndFreq - sweepStartFreq) * sweepProgress; - double currentPhaseIncrement = yup::MathConstants::twoPi * currentFreq / sampleRate; - - float sample = std::sin (phase); - phase += currentPhaseIncrement; - - if (phase >= yup::MathConstants::twoPi) - phase -= yup::MathConstants::twoPi; + float sample = readSineTable(); + advancePhase (currentFreq, generationSampleRate); // Update sweep progress - sweepProgress += 1.0 / (sweepDurationSeconds * sampleRate); + sweepProgress += 1.0 / (sweepDurationSeconds * generationSampleRate); if (sweepProgress >= 1.0) sweepProgress = 0.0; // Loop the sweep @@ -205,16 +402,20 @@ class SignalGenerator void updatePhaseIncrement() { - phaseIncrement = yup::MathConstants::twoPi * frequency / sampleRate; + phaseIncrement = frequency / sampleRate; } + std::array sineTable; + double sampleRate; + double sourceSampleRate; double frequency; double phase; double phaseIncrement = 0.0; float amplitude; SignalType signalType; + SweepPlaybackMode sweepPlaybackMode; // Sweep parameters double sweepStartFreq, sweepEndFreq, sweepDurationSeconds; @@ -228,6 +429,20 @@ class SignalGenerator // Smoothed parameter values yup::SmoothedValue smoothedFrequency; yup::SmoothedValue smoothedAmplitude; + + // Resampling state + yup::ResamplerFloat resampler; + yup::Oversampler2xFloat oversampler2x; + yup::Oversampler4xFloat oversampler4x; + yup::Oversampler8xFloat oversampler8x; + std::vector sourceBuffer; + std::vector silenceBuffer; + std::vector resampledBuffer; + int maxOutputBlockSize = 0; + int sourceBlockCapacity = 0; + int resampledBlockCapacity = 0; + int pendingOutputSamples = 0; + int pendingReadPosition = 0; }; //============================================================================== @@ -326,11 +541,14 @@ class SpectrumAnalyzerDemo int numSamples, const yup::AudioIODeviceCallbackContext& context) override { - // Generate test audio samples + if (monoOutputBuffer.size() < static_cast (numSamples)) + monoOutputBuffer.resize (static_cast (numSamples), 0.0f); + + signalGenerator.renderNextBlock (monoOutputBuffer.data(), numSamples); + for (int sample = 0; sample < numSamples; ++sample) { - // Generate audio sample using signal generator - const float audioSample = signalGenerator.getNextSample(); + const float audioSample = monoOutputBuffer[static_cast (sample)]; // Output to all channels for (int channel = 0; channel < numOutputChannels; ++channel) @@ -344,12 +562,14 @@ class SpectrumAnalyzerDemo void audioDeviceAboutToStart (yup::AudioIODevice* device) override { double sampleRate = device->getCurrentSampleRate(); + const int maxBlockSize = yup::jmax (1, device->getCurrentBufferSizeSamples()); // Setup signal generator - signalGenerator.setSampleRate (sampleRate); + signalGenerator.prepare (sampleRate, maxBlockSize); signalGenerator.setFrequency (currentFrequency); signalGenerator.setAmplitude (currentAmplitude); signalGenerator.setSweepParameters (20.0, 22000.0, sweepDurationSeconds); + monoOutputBuffer.assign (static_cast (maxBlockSize), 0.0f); // Configure spectrum analyzer analyzerComponent.setSampleRate (sampleRate); @@ -374,11 +594,15 @@ class SpectrumAnalyzerDemo // Signal type selector signalTypeCombo = std::make_unique ("SignalType"); signalTypeCombo->addItem ("Single Tone", 1); - signalTypeCombo->addItem ("Sweep", 2); - signalTypeCombo->addItem ("White Noise", 3); - signalTypeCombo->addItem ("Pink Noise", 4); - signalTypeCombo->addItem ("Brown Noise", 5); - signalTypeCombo->setSelectedId (1); + signalTypeCombo->addItem ("Sweep Direct", 2); + signalTypeCombo->addItem ("Sweep Resampler", 3); + signalTypeCombo->addItem ("Sweep 2x", 4); + signalTypeCombo->addItem ("Sweep 4x", 5); + signalTypeCombo->addItem ("Sweep 8x", 6); + signalTypeCombo->addItem ("White Noise", 7); + signalTypeCombo->addItem ("Pink Noise", 8); + signalTypeCombo->addItem ("Brown Noise", 9); + signalTypeCombo->setSelectedId (3); signalTypeCombo->onSelectedItemChanged = [this] { updateSignalType(); @@ -533,6 +757,8 @@ class SpectrumAnalyzerDemo label->setFont (labelFont); addAndMakeVisible (*label); } + + updateSignalType(); } void setupAudio() @@ -609,6 +835,7 @@ class SpectrumAnalyzerDemo void updateSignalType() { SignalGenerator::SignalType signalType = SignalGenerator::SignalType::singleTone; + auto sweepPlaybackMode = SignalGenerator::SweepPlaybackMode::resampled; switch (signalTypeCombo->getSelectedId()) { @@ -617,22 +844,39 @@ class SpectrumAnalyzerDemo break; case 2: signalType = SignalGenerator::SignalType::frequencySweep; + sweepPlaybackMode = SignalGenerator::SweepPlaybackMode::direct; break; case 3: - signalType = SignalGenerator::SignalType::whiteNoise; + signalType = SignalGenerator::SignalType::frequencySweep; + sweepPlaybackMode = SignalGenerator::SweepPlaybackMode::resampled; break; case 4: - signalType = SignalGenerator::SignalType::pinkNoise; + signalType = SignalGenerator::SignalType::frequencySweep; + sweepPlaybackMode = SignalGenerator::SweepPlaybackMode::resampledOversampled2x; break; case 5: + signalType = SignalGenerator::SignalType::frequencySweep; + sweepPlaybackMode = SignalGenerator::SweepPlaybackMode::resampledOversampled4x; + break; + case 6: + signalType = SignalGenerator::SignalType::frequencySweep; + sweepPlaybackMode = SignalGenerator::SweepPlaybackMode::resampledOversampled8x; + break; + case 7: + signalType = SignalGenerator::SignalType::whiteNoise; + break; + case 8: + signalType = SignalGenerator::SignalType::pinkNoise; + break; + case 9: signalType = SignalGenerator::SignalType::brownNoise; break; } signalGenerator.setSignalType (signalType); + signalGenerator.setSweepPlaybackMode (sweepPlaybackMode); // Enable/disable frequency and sweep controls based on signal type - bool isToneOrSweep = (signalType == SignalGenerator::SignalType::singleTone || signalType == SignalGenerator::SignalType::frequencySweep); frequencySlider->setEnabled (signalType == SignalGenerator::SignalType::singleTone); sweepDurationSlider->setEnabled (signalType == SignalGenerator::SignalType::frequencySweep); } @@ -743,6 +987,7 @@ class SpectrumAnalyzerDemo std::unique_ptr fftInfoLabel; yup::OwnedArray parameterLabels; + std::vector monoOutputBuffer; // Parameters double currentFrequency = 440.0;