resampler: fix SincResampler
diff --git a/src/flowgraph/resampler/ContinuousResampler.h b/src/flowgraph/resampler/ContinuousResampler.h
index f17342e..de7a973 100644
--- a/src/flowgraph/resampler/ContinuousResampler.h
+++ b/src/flowgraph/resampler/ContinuousResampler.h
@@ -44,7 +44,7 @@
mPhase += mPhaseIncrement;
}
- float getPhase() {
+ double getPhase() {
return (float) mPhase;
}
diff --git a/src/flowgraph/resampler/MultiChannelResampler.cpp b/src/flowgraph/resampler/MultiChannelResampler.cpp
index 01e60d8..42c756b 100644
--- a/src/flowgraph/resampler/MultiChannelResampler.cpp
+++ b/src/flowgraph/resampler/MultiChannelResampler.cpp
@@ -91,9 +91,9 @@
}
void MultiChannelResampler::writeFrame(const float *frame) {
- // Advance cursor before write so that cursor points to last written frame in read.
- if (++mCursor >= getNumTaps()) {
- mCursor = 0;
+ // Move cursor before write so that cursor points to last written frame in read.
+ if (--mCursor < 0) {
+ mCursor = getNumTaps() - 1;
}
float *dest = &mX[mCursor * getChannelCount()];
int offset = getNumTaps() * getChannelCount();
@@ -118,3 +118,43 @@
float MultiChannelResampler::calculateWindowedSinc(float radians, int spread) {
return sinc(radians) * hammingWindow(radians, spread); // TODO try Kaiser window
}
+
+
+// Generate coefficients in the order they will be used by readFrame().
+// This is more complicated but readFrame() is called repeatedly and should be optimized.
+void MultiChannelResampler::generateCoefficients(int32_t inputRate,
+ int32_t outputRate,
+ int32_t numRows,
+ double phaseIncrement,
+ float normalizedCutoff) {
+ mCoefficients.resize(getNumTaps() * numRows);
+ int cursor = 0;
+ double phase = 0.0;
+ const float cutoffScaler = std::min(1.0f, normalizedCutoff * outputRate / inputRate);
+ const int numTapsHalf = getNumTaps() / 2; // numTaps must be even.
+ for (int i = 0; i < numRows; i++) {
+ float tapPhase = phase - numTapsHalf;
+ float gain = 0.0;
+ int gainCursor = cursor;
+ for (int tap = 0; tap < getNumTaps(); tap++) {
+ float radians = tapPhase * M_PI;
+ float coefficient = sinc(cutoffScaler * radians) * hammingWindow(radians, numTapsHalf);
+ mCoefficients.at(cursor++) = coefficient;
+ gain += coefficient;
+ tapPhase += 1.0;
+ }
+ phase += phaseIncrement;
+ while (phase >= 1.0) {
+ phase -= 1.0;
+ }
+
+ // Correct for gain variations. // TODO review
+ //printf("gain at %d was %f\n", i, gain);
+ float gainCorrection = 1.0 / gain; // normalize the gain
+ for (int tap = 0; tap < getNumTaps(); tap++) {
+ float scaledCoefficient = mCoefficients.at(gainCursor) * gainCorrection;
+ //printf("scaledCoefficient[%2d] = %10.6f\n", tap, scaledCoefficient);
+ mCoefficients.at(gainCursor++) = scaledCoefficient;
+ }
+ }
+}
diff --git a/src/flowgraph/resampler/MultiChannelResampler.h b/src/flowgraph/resampler/MultiChannelResampler.h
index dd1653e..a18ed6a 100644
--- a/src/flowgraph/resampler/MultiChannelResampler.h
+++ b/src/flowgraph/resampler/MultiChannelResampler.h
@@ -215,7 +215,22 @@
std::vector<float> mX;
std::vector<float> mSingleFrame;
+ std::vector<float> mCoefficients;
static constexpr int kMaxCoefficients = 8 * 1024;
+
+
+ /**
+ * Generate the filter coefficients in optimal order.
+ * @param inputRate
+ * @param outputRate
+ * @param normalizedCutoff filter cutoff frequency normalized to Nyquist rate of output
+ */
+ void generateCoefficients(int32_t inputRate,
+ int32_t outputRate,
+ int32_t numRows,
+ double phaseIncrement,
+ float normalizedCutoff);
+
private:
// max coefficients for polyphase filter
diff --git a/src/flowgraph/resampler/PolyphaseResampler.cpp b/src/flowgraph/resampler/PolyphaseResampler.cpp
index 8e3ce83..18ebe4c 100644
--- a/src/flowgraph/resampler/PolyphaseResampler.cpp
+++ b/src/flowgraph/resampler/PolyphaseResampler.cpp
@@ -24,65 +24,36 @@
: MultiChannelResampler(builder)
{
assert((getNumTaps() % 4) == 0); // Required for loop unrolling.
- generateCoefficients(builder.getInputRate(), builder.getOutputRate(), builder.getNormalizedCutoff());
-}
-// Generate coefficients in the order they will be used by readFrame().
-// This is more complicated but readFrame() is called repeatedly and should be optimized.
-void PolyphaseResampler::generateCoefficients(int32_t inputRate,
- int32_t outputRate,
- float normalizedCutoff) {
- {
- // Reduce sample rates to the smallest ratio.
- // For example 44100/48000 would become 147/160.
- IntegerRatio ratio(inputRate, outputRate);
- ratio.reduce();
- mNumerator = ratio.getNumerator();
- mDenominator = ratio.getDenominator();
- mIntegerPhase = mDenominator;
- mCoefficients.resize(getNumTaps() * ratio.getDenominator());
- }
- int cursor = 0;
- double phase = 0.0;
+ int32_t inputRate = builder.getInputRate();
+ int32_t outputRate = builder.getOutputRate();
+
+ // Reduce sample rates to the smallest ratio.
+ // For example 44100/48000 would become 147/160.
+ IntegerRatio ratio(inputRate, outputRate);
+ ratio.reduce();
+ mNumerator = ratio.getNumerator();
+ mDenominator = ratio.getDenominator();
+ mIntegerPhase = mDenominator;
+
+ int32_t numRows = mDenominator;
double phaseIncrement = (double) inputRate / (double) outputRate;
- const float cutoffScaler = std::min(1.0f, normalizedCutoff * outputRate / inputRate);
- const int spread = getNumTaps() / 2; // numTaps must be even.
- for (int i = 0; i < mDenominator; i++) {
- float tapPhase = phase - spread;
- float gain = 0.0;
- int gainCursor = cursor;
- for (int tap = 0; tap < getNumTaps(); tap++) {
- float radians = tapPhase * M_PI;
- float coefficient = sinc(cutoffScaler * radians) * hammingWindow(radians, spread);
- mCoefficients.at(cursor++) = coefficient;
- gain += coefficient;
- tapPhase += 1.0;
- }
- phase += phaseIncrement;
- while (phase >= 1.0) {
- phase -= 1.0;
- }
-
- // Correct for gain variations. // TODO review
- //printf("gain at %d was %f\n", i, gain);
- float gainCorrection = 1.0 / gain; // normalize the gain
- for (int tap = 0; tap < getNumTaps(); tap++) {
- float scaledCoefficient = mCoefficients.at(gainCursor) * gainCorrection;
- //printf("scaledCoefficient[%2d] = %10.6f\n", tap, scaledCoefficient);
- mCoefficients.at(gainCursor++) = scaledCoefficient;
- }
- }
+ generateCoefficients(inputRate, outputRate,
+ numRows, phaseIncrement,
+ builder.getNormalizedCutoff());
}
void PolyphaseResampler::readFrame(float *frame) {
// Clear accumulator for mixing.
std::fill(mSingleFrame.begin(), mSingleFrame.end(), 0.0);
+// printf("PolyphaseResampler: mCoefficientCursor = %4d\n", mCoefficientCursor);
// Multiply input times windowed sinc function.
float *coefficients = &mCoefficients[mCoefficientCursor];
float *xFrame = &mX[mCursor * getChannelCount()];
for (int i = 0; i < mNumTaps; i++) {
float coefficient = *coefficients++;
+// printf("PolyphaseResampler: coeff = %10.6f, xFrame[0] = %10.6f\n", coefficient, xFrame[0]);
for (int channel = 0; channel < getChannelCount(); channel++) {
mSingleFrame[channel] += *xFrame++ * coefficient;
}
diff --git a/src/flowgraph/resampler/PolyphaseResampler.h b/src/flowgraph/resampler/PolyphaseResampler.h
index a4cd5bf..fa267a1 100644
--- a/src/flowgraph/resampler/PolyphaseResampler.h
+++ b/src/flowgraph/resampler/PolyphaseResampler.h
@@ -54,23 +54,11 @@
protected:
- std::vector<float> mCoefficients;
int32_t mCoefficientCursor = 0;
int32_t mIntegerPhase = 0;
int32_t mNumerator = 0;
int32_t mDenominator = 0;
-private:
-
- /**
- * Generate the filter coefficients in optimal order.
- * @param inputRate
- * @param outputRate
- * @param normalizedCutoff filter cutoff frequency normalized to Nyquist rate of output
- */
- void generateCoefficients(int32_t inputRate,
- int32_t outputRate,
- float normalizedCutoff);
};
}
diff --git a/src/flowgraph/resampler/PolyphaseResamplerStereo.cpp b/src/flowgraph/resampler/PolyphaseResamplerStereo.cpp
index 4618f57..5c73a8e 100644
--- a/src/flowgraph/resampler/PolyphaseResamplerStereo.cpp
+++ b/src/flowgraph/resampler/PolyphaseResamplerStereo.cpp
@@ -25,7 +25,6 @@
assert(builder.getChannelCount() == STEREO);
}
-
void PolyphaseResamplerStereo::writeFrame(const float *frame) {
// Move cursor before write so that cursor points to last written frame in read.
if (--mCursor < 0) {
diff --git a/src/flowgraph/resampler/SincResampler.cpp b/src/flowgraph/resampler/SincResampler.cpp
index 444b5d1..9acd192 100644
--- a/src/flowgraph/resampler/SincResampler.cpp
+++ b/src/flowgraph/resampler/SincResampler.cpp
@@ -23,54 +23,33 @@
: ContinuousResampler(builder)
, mSingleFrame2(builder.getChannelCount()) {
assert((getNumTaps() % 4) == 0); // Required for loop unrolling.
- generateCoefficients(builder.getInputRate(), builder.getOutputRate(), builder.getNormalizedCutoff());
+ mNumRows = kMaxCoefficients / getNumTaps(); // no guard row needed
+// printf("SincResampler: numRows = %d\n", mNumRows);
+ double phaseIncrement = 1.0 / mNumRows;
+ generateCoefficients(builder.getInputRate(),
+ builder.getOutputRate(),
+ mNumRows,
+ phaseIncrement,
+ builder.getNormalizedCutoff());
}
-// Generate coefficients from min to max.
-void SincResampler::generateCoefficients(int32_t inputRate,
- int32_t outputRate,
- float normalizedCutoff) {
- mNumSeries = kMaxCoefficients / getNumTaps();
- mTablePhaseScaler = mNumSeries / 2.0f;
- mCoefficients.resize(getNumTaps() * mNumSeries);
-
- int cursor = 0;
- const float cutoffScaler = std::min(1.0f, normalizedCutoff * outputRate / inputRate);
- const int numTapsHalf = getNumTaps() / 2; // numTaps must be even.
- for (int i = 0; i < mNumSeries; i++) {
- float tapPhase = ((float)i / mNumSeries) - numTapsHalf;
- float gain = 0.0;
- int gainCursor = cursor;
- for (int tap = 0; tap < getNumTaps(); tap++) {
- float radians = tapPhase * M_PI;
- float coefficient = sinc(cutoffScaler * radians) * hammingWindow(radians, numTapsHalf);
- mCoefficients.at(cursor++) = coefficient;
- gain += coefficient;
- tapPhase += 1.0;
- }
-
- // Correct for gain variations. // TODO review
- //printf("gain at %d was %f\n", i, gain);
- float gainCorrection = 1.0 / gain; // normalize the gain
- for (int tap = 0; tap < getNumTaps(); tap++) {
- float scaledCoefficient = mCoefficients.at(gainCursor) * gainCorrection;
- //printf("scaledCoefficient[%2d] = %10.6f\n", tap, scaledCoefficient);
- mCoefficients.at(gainCursor++) = scaledCoefficient;
- }
- }
-}
-
+// Multiply input times windowed sinc function.
void SincResampler::readFrame(float *frame) {
// Clear accumulator for mixing.
std::fill(mSingleFrame.begin(), mSingleFrame.end(), 0.0);
+ std::fill(mSingleFrame2.begin(), mSingleFrame2.end(), 0.0);
- // Multiply input times windowed sinc function.
- float tablePhase = getPhase() * mTablePhaseScaler;
- long index1 = lrintf(tablePhase);
+ // Determine indices into coefficients table.
+ double tablePhase = getPhase() * mNumRows;
+ int index1 = static_cast<int>(floor(tablePhase));
+ if (index1 >= mNumRows) { // no guard row needed because we wrap the indices
+ tablePhase -= mNumRows;
+ index1 -= mNumRows;
+ }
float *coefficients1 = &mCoefficients[index1 * getNumTaps()];
- int index2 = (index1 + 1) * getNumTaps();
- if (index2 >= mCoefficients.size()) {
- index2 = 0;
+ int index2 = index1 + 1;
+ if (index2 >= mNumRows) { // no guard row needed because we wrap the indices
+ index2 -= mNumRows;
}
float *coefficients2 = &mCoefficients[index2 * getNumTaps()];
float *xFrame = &mX[mCursor * getChannelCount()];
@@ -79,7 +58,7 @@
float coefficient2 = *coefficients2++;
for (int channel = 0; channel < getChannelCount(); channel++) {
float sample = *xFrame++;
- mSingleFrame[channel] += sample* coefficient1;
+ mSingleFrame[channel] += sample * coefficient1;
mSingleFrame2[channel] += sample * coefficient2;
}
}
@@ -92,51 +71,3 @@
frame[channel] = low + (fraction * (high - low));
}
}
-
-//void SincResampler::readFrame(float *frame) {
-// // Clear accumulator for mix.
-// for (int channel = 0; channel < getChannelCount(); channel++) {
-// mSingleFrame[channel] = 0.0;
-// }
-// float phase = getPhase();
-// // Multiply input times windowed sinc function.
-// int xIndex = (mCursor + mNumTaps) * getChannelCount();
-// for (int i = 0; i < mNumTaps; i++) {
-// // TODO Use consecutive coefficient precomputed and then interpolate the result.
-// float coefficient = interpolateWindowedSinc(phase);
-// float *xFrame = &mX[xIndex];
-// for (int channel = 0; channel < getChannelCount(); channel++) {
-// mSingleFrame[channel] += coefficient * xFrame[channel];
-// }
-// xIndex -= getChannelCount();
-// phase += 1.0;
-// }
-// // Copy accumulator to output.
-// for (int channel = 0; channel < getChannelCount(); channel++) {
-// frame[channel] = mSingleFrame[channel];
-// }
-//}
-
-//float SincResampler::interpolateWindowedSinc(float phase) {
-// // convert from 0 to 2*kSpread range to 0 to table size range
-// float tablePhase = phase * mTablePhaseScaler;
-// //uint16_t tableIndex = fast_int_from_float(tablePhase);
-// long tableIndex = lrintf(tablePhase);
-// //int tableIndex = int(tablePhase);
-// float low = mWindowedSinc[tableIndex];
-// float high = mWindowedSinc[tableIndex + 1]; // OK because of guard point
-// float fraction = tablePhase - tableIndex;
-// return low + (fraction * (high - low));
-//}
-
-//void SincResampler::generateLookupTable() {
- // Place related coefficients together for faster convolution, like for Polyphase
- // but divide each zero crossing into 32 or 64 steps.
- // TODO store only half of function and use symmetry
- // By iterating over the table size we also set the guard point.
-// for (int i = 0; i < mWindowedSinc.size(); i++) {
-// float phase = (i * 2.0 * kSpread) / kNumPoints;
-// float radians = (phase - kSpread) * M_PI;
-// mWindowedSinc[i] = calculateWindowedSinc(radians, kSpread);
-// }
-//}
diff --git a/src/flowgraph/resampler/SincResampler.h b/src/flowgraph/resampler/SincResampler.h
index a195c84..cbfcc0b 100644
--- a/src/flowgraph/resampler/SincResampler.h
+++ b/src/flowgraph/resampler/SincResampler.h
@@ -32,34 +32,10 @@
void readFrame(float *frame) override;
- /**
- * @param phase between 0.0 and 2*kSpread
- * @return windowedSinc
- */
- float interpolateWindowedSinc(float phase);
-
protected:
- std::vector<float> mWindowedSinc;
-
-private:
-
- /**
- * Generate the filter coefficients in optimal order.
- * @param inputRate
- * @param outputRate
- * @param normalizedCutoff filter cutoff frequency normalized to Nyquist rate of output
- */
- void generateCoefficients(int32_t inputRate,
- int32_t outputRate,
- float normalizedCutoff);
-
- std::vector<float> mCoefficients;
- int32_t mNumSeries = 0;
- std::vector<float> mSingleFrame2; // for interpolation
-
-
- float mTablePhaseScaler = 0.0f;
+ std::vector<float> mSingleFrame2; // for interpolation
+ int32_t mNumRows = 0;
};
diff --git a/src/flowgraph/resampler/SincResamplerStereo.cpp b/src/flowgraph/resampler/SincResamplerStereo.cpp
index 1a8d1fc..5af7ff7 100644
--- a/src/flowgraph/resampler/SincResamplerStereo.cpp
+++ b/src/flowgraph/resampler/SincResamplerStereo.cpp
@@ -14,6 +14,8 @@
* limitations under the License.
*/
+#include <math.h>
+
#include "SincResamplerStereo.h"
using namespace resampler;
@@ -26,32 +28,53 @@
}
void SincResamplerStereo::writeFrame(const float *frame) {
- int xIndex = mCursor * STEREO;
- int offset = mNumTaps * STEREO;
- float *dest = &mX[xIndex];
- // Write each channel twice so we avoid having to wrap when running the FIR.
- dest[0] = dest[offset] = frame[0];
- dest[1] = dest[1 + offset] = frame[1];
- if (++mCursor >= mNumTaps) {
- mCursor = 0;
+ // Move cursor before write so that cursor points to last written frame in read.
+ if (--mCursor < 0) {
+ mCursor = getNumTaps() - 1;
}
+ float *dest = &mX[mCursor * STEREO];
+ const int offset = mNumTaps * STEREO;
+ // Write each channel twice so we avoid having to wrap when running the FIR.
+ const float left = frame[0];
+ const float right = frame[1];
+ // Put ordered writes together.
+ dest[0] = left;
+ dest[1] = right;
+ dest[offset] = left;
+ dest[1 + offset] = right;
}
+// Multiply input times windowed sinc function.
void SincResamplerStereo::readFrame(float *frame) {
-// float left = 0.0;
-// float right = 0.0;
-// float phase = getPhase();
-// // Multiply input times windowed sinc function.
-// int xIndex = (mCursor + mNumTaps) * STEREO;
-// for (int i = 0; i < mNumTaps; i++) {
-// float coefficient = interpolateWindowedSinc(phase);
-// float *xFrame = &mX[xIndex];
-// left += coefficient * xFrame[0];
-// right += coefficient * xFrame[1];
-// xIndex -= STEREO;
-// phase += 1.0;
-// }
-// // Copy accumulator to output.
-// frame[0] = left;
-// frame[1] = right;
+ // Clear accumulator for mixing.
+ std::fill(mSingleFrame.begin(), mSingleFrame.end(), 0.0);
+ std::fill(mSingleFrame2.begin(), mSingleFrame2.end(), 0.0);
+
+ // Determine indices into coefficients table.
+ double tablePhase = getPhase() * mNumRows;
+ int index1 = static_cast<int>(floor(tablePhase));
+ float *coefficients1 = &mCoefficients[index1 * getNumTaps()];
+ int index2 = (index1 + 1);
+ if (index2 >= mNumRows) { // no guard row needed because we wrap the indices
+ index2 = 0;
+ }
+ float *coefficients2 = &mCoefficients[index2 * getNumTaps()];
+ float *xFrame = &mX[mCursor * getChannelCount()];
+ for (int i = 0; i < mNumTaps; i++) {
+ float coefficient1 = *coefficients1++;
+ float coefficient2 = *coefficients2++;
+ for (int channel = 0; channel < getChannelCount(); channel++) {
+ float sample = *xFrame++;
+ mSingleFrame[channel] += sample * coefficient1;
+ mSingleFrame2[channel] += sample * coefficient2;
+ }
+ }
+
+ // Interpolate and copy to output.
+ float fraction = tablePhase - index1;
+ for (int channel = 0; channel < getChannelCount(); channel++) {
+ float low = mSingleFrame[channel];
+ float high = mSingleFrame2[channel];
+ frame[channel] = low + (fraction * (high - low));
+ }
}