From 16b964020fdf9deda6262e7dd9048e36acc0912e Mon Sep 17 00:00:00 2001 From: Hans Baier Date: Fri, 24 Jul 2009 05:27:49 +0000 Subject: interpolation.cc/h: Remove all failed and obsolete attempts, leave linear and cubic git-svn-id: svn://localhost/ardour2/branches/3.0@5424 d708f5d6-7413-0410-9779-e7cbd77b26cf --- libs/ardour/ardour/interpolation.h | 104 +------------- libs/ardour/interpolation.cc | 239 -------------------------------- libs/ardour/tests/interpolation-test.cc | 192 +++++++++---------------- libs/ardour/tests/interpolation-test.h | 15 +- 4 files changed, 71 insertions(+), 479 deletions(-) (limited to 'libs/ardour') diff --git a/libs/ardour/ardour/interpolation.h b/libs/ardour/ardour/interpolation.h index 4ff3163cc6..98b6e7e775 100644 --- a/libs/ardour/ardour/interpolation.h +++ b/libs/ardour/ardour/interpolation.h @@ -38,51 +38,6 @@ class Interpolation { } }; -// 40.24 fixpoint math -#define FIXPOINT_ONE 0x1000000 - -class FixedPointLinearInterpolation : public Interpolation { - protected: - /// speed in fixed point math - uint64_t phi; - - /// target speed in fixed point math - uint64_t target_phi; - - std::vector last_phase; - - // Fixed point is just an integer with an implied scaling factor. - // In 40.24 the scaling factor is 2^24 = 16777216, - // so a value of 10*2^24 (in integer space) is equivalent to 10.0. - // - // The advantage is that addition and modulus [like x = (x + y) % 2^40] - // have no rounding errors and no drift, and just require a single integer add. - // (swh) - - static const int64_t fractional_part_mask = 0xFFFFFF; - static const Sample binary_scaling_factor = 16777216.0f; - - public: - - FixedPointLinearInterpolation () : phi (FIXPOINT_ONE), target_phi (FIXPOINT_ONE) {} - - void set_speed (double new_speed) { - target_phi = (uint64_t) (FIXPOINT_ONE * fabs(new_speed)); - phi = target_phi; - } - - uint64_t get_phi() { return phi; } - uint64_t get_target_phi() { return target_phi; } - uint64_t get_last_phase() { assert(last_phase.size()); return last_phase[0]; } - void set_last_phase(uint64_t phase) { assert(last_phase.size()); last_phase[0] = phase; } - - void add_channel_to (int input_buffer_size, int output_buffer_size); - void remove_channel_from (); - - nframes_t interpolate (int channel, nframes_t nframes, Sample* input, Sample* output); - void reset (); -}; - class LinearInterpolation : public Interpolation { protected: @@ -92,7 +47,7 @@ class LinearInterpolation : public Interpolation { class CubicInterpolation : public Interpolation { protected: - // shamelessly ripped from Steve Harris' swh-plugins + // shamelessly ripped from Steve Harris' swh-plugins (ladspa-util.h) static inline float cube_interp(const float fr, const float inm1, const float in, const float inp1, const float inp2) { @@ -105,63 +60,6 @@ class CubicInterpolation : public Interpolation { nframes_t interpolate (int channel, nframes_t nframes, Sample* input, Sample* output); }; - -/** - * @class SplineInterpolation - * - * @brief interpolates using cubic spline interpolation over an input period - * - * Splines are piecewise cubic functions between each samples, - * where the cubic polynomials' values, first and second derivatives are equal - * on each sample point. - * - * The interpolation polynomial in the i-th interval then has the form - * p_i(x) = a3 (x - i)^3 + a2 (x - i)^2 + a1 (x - i) + a0 - * = ((a3 * (x - i) + a2) * (x - i) + a1) * (x - i) + a0 - * where - * a3 = (M[i+1] - M[i]) / 6 - * a2 = M[i] / 2 - * a1 = y[i+1] - y[i] - M[i+1]/6 - M[i]/3 - * a0 = y[i] - * - * The M's are calculated recursively: - * M[i+2] = 6.0 * (y[i] - 2y[i+1] + y[i+2]) - 4M[i+1] - M[i] - * - */ -class SplineInterpolation : public Interpolation { - protected: - double M[2]; - - public: - void reset (); - SplineInterpolation(); - nframes_t interpolate (int channel, nframes_t nframes, Sample* input, Sample* output); -}; - -class LibSamplerateInterpolation : public Interpolation { - protected: - std::vector state; - std::vector data; - - int error; - - void reset_state (); - - public: - LibSamplerateInterpolation (); - ~LibSamplerateInterpolation (); - - void set_speed (double new_speed); - void set_target_speed (double new_speed) {} - double speed () const { return _speed; } - - void add_channel_to (int input_buffer_size, int output_buffer_size); - void remove_channel_from (); - - nframes_t interpolate (int channel, nframes_t nframes, Sample* input, Sample* output); - void reset() { reset_state (); } -}; - } // namespace ARDOUR #endif diff --git a/libs/ardour/interpolation.cc b/libs/ardour/interpolation.cc index c3a45a0401..9a45d560c0 100644 --- a/libs/ardour/interpolation.cc +++ b/libs/ardour/interpolation.cc @@ -5,69 +5,6 @@ using namespace ARDOUR; -nframes_t -FixedPointLinearInterpolation::interpolate (int channel, nframes_t nframes, Sample *input, Sample *output) -{ - // the idea behind phase is that when the speed is not 1.0, we have to - // interpolate between samples and then we have to store where we thought we were. - // rather than being at sample N or N+1, we were at N+0.8792922 - // so the "phase" element, if you want to think about this way, - // varies from 0 to 1, representing the "offset" between samples - uint64_t the_phase = last_phase[channel]; - - // acceleration - int64_t phi_delta; - - // phi = fixed point speed - if (phi != target_phi) { - phi_delta = ((int64_t)(target_phi - phi)) / nframes; - } else { - phi_delta = 0; - } - - // index in the input buffers - nframes_t i = 0; - - for (nframes_t outsample = 0; outsample < nframes; ++outsample) { - i = the_phase >> 24; - Sample fractional_phase_part = (the_phase & fractional_part_mask) / binary_scaling_factor; - - if (input && output) { - // Linearly interpolate into the output buffer - output[outsample] = - input[i] * (1.0f - fractional_phase_part) + - input[i+1] * fractional_phase_part; - } - - the_phase += phi + phi_delta; - } - - last_phase[channel] = (the_phase & fractional_part_mask); - - // playback distance - return i; -} - -void -FixedPointLinearInterpolation::add_channel_to (int /*input_buffer_size*/, int /*output_buffer_size*/) -{ - last_phase.push_back (0); -} - -void -FixedPointLinearInterpolation::remove_channel_from () -{ - last_phase.pop_back (); -} - -void -FixedPointLinearInterpolation::reset() -{ - for (size_t i = 0; i <= last_phase.size(); i++) { - last_phase[i] = 0; - } -} - nframes_t LinearInterpolation::interpolate (int channel, nframes_t nframes, Sample *input, Sample *output) @@ -144,179 +81,3 @@ CubicInterpolation::interpolate (int channel, nframes_t nframes, Sample *input, return i; } - -SplineInterpolation::SplineInterpolation() -{ - reset (); -} - -void SplineInterpolation::reset() -{ - Interpolation::reset(); - M[0] = 0.0; - M[1] = 0.0; - M[2] = 0.0; -} - -nframes_t -SplineInterpolation::interpolate (int channel, nframes_t nframes, Sample *input, Sample *output) -{ - - // now interpolate - // index in the input buffers - nframes_t i = 0, delta_i = 0; - - double acceleration; - double distance = 0.0; - - if (_speed != _target_speed) { - acceleration = _target_speed - _speed; - } else { - acceleration = 0.0; - } - - distance = phase[channel]; - assert(distance >= 0.0 && distance < 1.0); - - for (nframes_t outsample = 0; outsample < nframes; outsample++) { - i = floor(distance); - - double x = double(distance) - double(i); - - // if distance is something like 0.999999999999 - // it will get rounded to 1 in the conversion to float above - while (x >= 1.0) { - x -= 1.0; - i++; - } - - assert(x >= 0.0 && x < 1.0); - - if (input && output) { - // if i changed, recalculate coefficients - if (delta_i == 1) { - // if i changed, rotate the M's - M[0] = M[1]; - M[1] = M[2]; - M[2] = 6.0 * (input[i] - 2.0*input[i+1] + input[i+2]) - 4.0*M[1] - M[0]; - printf("\ny[%d] = %lf\n", i, input[i]); - printf("y[%d] = %lf\n", i+1, input[i+1]); - printf("y[%d] = %lf\n\n", i+2, input[i+2]); - printf("M[2] = %lf M[1] = %lf M[0] = %lf y-term: %lf M-term: %lf\n", - M[2], M[1], M[0], 6.0 * (input[i] - 2.0*input[i+1] + input[i+2]), - - 4.0*M[1] - M[0]); - } - double a3 = (M[1] - M[0]) / 6.0; - double a2 = M[0] / 2.0; - double a1 = input[i+1] - input[i] - (M[1] + 2.0*M[0]) / 6.0; - double a0 = input[i]; - // interpolate into the output buffer - output[outsample] = ((a3*x + a2)*x + a1)*x + a0; - //printf( "input[%d/%d] = %lf/%lf distance: %lf output[%d] = %lf\n", i, i+1, input[i], input[i+1], distance, outsample, output[outsample]); - - } - distance += _speed + acceleration; - - delta_i = floor(distance) - i; - } - - i = floor(distance); - phase[channel] = distance - floor(distance); - assert (phase[channel] >= 0.0 && phase[channel] < 1.0); - - return i; -} - -LibSamplerateInterpolation::LibSamplerateInterpolation() : state (0) -{ - _speed = 1.0; -} - -LibSamplerateInterpolation::~LibSamplerateInterpolation() -{ - for (size_t i = 0; i < state.size(); i++) { - state[i] = src_delete (state[i]); - } -} - -void -LibSamplerateInterpolation::set_speed (double new_speed) -{ - _speed = new_speed; - for (size_t i = 0; i < state.size(); i++) { - src_set_ratio (state[i], 1.0/_speed); - } -} - -void -LibSamplerateInterpolation::reset_state () -{ - printf("INTERPOLATION: reset_state()\n"); - for (size_t i = 0; i < state.size(); i++) { - if (state[i]) { - src_reset (state[i]); - } else { - state[i] = src_new (SRC_SINC_FASTEST, 1, &error); - } - } -} - -void -LibSamplerateInterpolation::add_channel_to (int input_buffer_size, int output_buffer_size) -{ - SRC_DATA* newdata = new SRC_DATA; - - /* Set up sample rate converter info. */ - newdata->end_of_input = 0 ; - - newdata->input_frames = input_buffer_size; - newdata->output_frames = output_buffer_size; - - newdata->input_frames_used = 0 ; - newdata->output_frames_gen = 0 ; - - newdata->src_ratio = 1.0/_speed; - - data.push_back (newdata); - state.push_back (0); - - reset_state (); -} - -void -LibSamplerateInterpolation::remove_channel_from () -{ - SRC_DATA* d = data.back (); - delete d; - data.pop_back (); - if (state.back ()) { - src_delete (state.back ()); - } - state.pop_back (); - reset_state (); -} - -nframes_t -LibSamplerateInterpolation::interpolate (int channel, nframes_t nframes, Sample *input, Sample *output) -{ - if (!data.size ()) { - printf ("ERROR: trying to interpolate with no channels\n"); - return 0; - } - - data[channel]->data_in = input; - data[channel]->data_out = output; - - data[channel]->input_frames = nframes * _speed; - data[channel]->output_frames = nframes; - data[channel]->src_ratio = 1.0/_speed; - - if ((error = src_process (state[channel], data[channel]))) { - printf ("\nError : %s\n\n", src_strerror (error)); - exit (1); - } - - //printf("INTERPOLATION: channel %d input_frames_used: %d\n", channel, data[channel]->input_frames_used); - - return data[channel]->input_frames_used; -} diff --git a/libs/ardour/tests/interpolation-test.cc b/libs/ardour/tests/interpolation-test.cc index 6a5bcd5ed4..45e350d12b 100644 --- a/libs/ardour/tests/interpolation-test.cc +++ b/libs/ardour/tests/interpolation-test.cc @@ -16,10 +16,7 @@ InterpolationTest::linearInterpolationTest () for (int i = 0; 3*i < NUM_SAMPLES - 1024;) { linear.set_speed (double(1.0)/double(3.0)); linear.set_target_speed (double(1.0)/double(3.0)); - //printf ("Interpolate: input: %d, output: %d, i: %d\n", input + i, output + i, i); result = linear.interpolate (0, 1024, input + i, output + i*3); - //printf ("Result: %d\n", result); - //CPPUNIT_ASSERT_EQUAL ((uint32_t)((NUM_SAMPLES - 100) * interpolation.speed()), result); i += result; } @@ -87,146 +84,87 @@ InterpolationTest::linearInterpolationTest () } /* for (int i=0; i < NUM_SAMPLES; ++i) { - cout << "input[" << i << "] = " << input[i] << " output[" << i << "] = " << output[i] << endl; - } - */ + cout << i << " " << output[i] << endl; + } + */ } void -InterpolationTest::splineInterpolationTest () +InterpolationTest::cubicInterpolationTest () { nframes_t result = 0; - cout << "\nspline Interpolation Test\n"; - - cout << "\nSpeed: 1/2" << endl; - spline.reset(); - spline.set_speed (0.5); - int one_period = 1024; - - /* - - for (int i = 0; 2 * i < NUM_SAMPLES - one_period;) { - result = spline.interpolate (0, one_period, input + i, output + 2*i); + cout << "\nCubic Interpolation Test\n"; + + cout << "\nSpeed: 1/3"; + for (int i = 0; 3*i < NUM_SAMPLES - 1024;) { + cubic.set_speed (double(1.0)/double(3.0)); + cubic.set_target_speed (double(1.0)/double(3.0)); + result = cubic.interpolate (0, 1024, input + i, output + i*3); i += result; } - for (int i=0; i < NUM_SAMPLES - one_period; ++i) { - //cout << "input[" << i << "] = " << input[i] << " output[" << i << "] = " << output[i] << endl; - if (i % 200 == 0) { CPPUNIT_ASSERT_EQUAL (double(1.0), double(output[i])); } - else if (i % 2 == 0) { CPPUNIT_ASSERT_EQUAL (double(0.0), double(output[i])); } - } - */ + cout << "\nSpeed: 1.0"; + cubic.reset(); + cubic.set_speed (1.0); + cubic.set_target_speed (cubic.speed()); + result = cubic.interpolate (0, NUM_SAMPLES, input, output); + CPPUNIT_ASSERT_EQUAL ((uint32_t)(NUM_SAMPLES * cubic.speed()), result); + for (int i = 0; i < NUM_SAMPLES; i += INTERVAL) { + CPPUNIT_ASSERT_EQUAL (1.0f, output[i]); + } - // square wave - for (int i = 0; i < NUM_SAMPLES; ++i) { - if (i % (INTERVAL/2) < INTERVAL/4 ) { - input[i] = 1.0f; - } else { - input[i] = 0.0f; - } - output[i] = 0.0f; + cout << "\nSpeed: 0.5"; + cubic.reset(); + cubic.set_speed (0.5); + cubic.set_target_speed (cubic.speed()); + result = cubic.interpolate (0, NUM_SAMPLES, input, output); + CPPUNIT_ASSERT_EQUAL ((uint32_t)(NUM_SAMPLES * cubic.speed()), result); + for (int i = 0; i < NUM_SAMPLES; i += (INTERVAL / cubic.speed() +0.5)) { + CPPUNIT_ASSERT_EQUAL (1.0f, output[i]); } + cout << "\nSpeed: 0.2"; + cubic.reset(); + cubic.set_speed (0.2); + cubic.set_target_speed (cubic.speed()); + result = cubic.interpolate (0, NUM_SAMPLES, input, output); + CPPUNIT_ASSERT_EQUAL ((uint32_t)(NUM_SAMPLES * cubic.speed()), result); + + cout << "\nSpeed: 0.02"; + cubic.reset(); + cubic.set_speed (0.02); + cubic.set_target_speed (cubic.speed()); + result = cubic.interpolate (0, NUM_SAMPLES, input, output); + CPPUNIT_ASSERT_EQUAL ((uint32_t)(NUM_SAMPLES * cubic.speed()), result); - /* - //sine wave - for (int i = 0; i < NUM_SAMPLES; ++i) { - input[i] = sin(double(i) * M_2_PI / INTERVAL * 10.0); - } + /* This one fails due too error accumulation + cout << "\nSpeed: 0.002"; + cubic.reset(); + cubic.set_speed (0.002); + cubic.set_target_speed (cubic.speed()); + result = cubic.interpolate (0, NUM_SAMPLES, input, output); + cubic.speed(); + CPPUNIT_ASSERT_EQUAL ((nframes_t)(NUM_SAMPLES * cubic.speed()), result); */ - one_period = 512; - - cout << "\nSpeed: 1/60" << endl; - spline.reset(); - spline.set_speed (1.0/90.0); - - - for (int i = 0, o = 0; 90 * i < NUM_SAMPLES - one_period; o++) { - result = spline.interpolate (0, one_period, input + i, output + o * one_period); - //printf ("Result: %d\n", result); - i += result; - } - - for (int i=0; i < NUM_SAMPLES - one_period; ++i) { - cout << i << " " << output[i] << endl; - //if (i % 333 == 0) { CPPUNIT_ASSERT_EQUAL (double(1.0), double(output[i])); } - //else if (i % 2 == 0) { CPPUNIT_ASSERT_EQUAL (double(0.0), double(output[i])); } + cout << "\nSpeed: 2.0"; + cubic.reset(); + cubic.set_speed (2.0); + cubic.set_target_speed (cubic.speed()); + result = cubic.interpolate (0, NUM_SAMPLES / 2, input, output); + CPPUNIT_ASSERT_EQUAL ((uint32_t)(NUM_SAMPLES / 2 * cubic.speed()), result); + for (int i = 0; i < NUM_SAMPLES / 2; i += (INTERVAL / cubic.speed() +0.5)) { + CPPUNIT_ASSERT_EQUAL (1.0f, output[i]); } -} - -void -InterpolationTest::libSamplerateInterpolationTest () -{ - nframes_t result; - - cout << "\nLibSamplerate Interpolation Test\n"; -/* - cout << "\nSpeed: 1.0"; - interpolation.set_speed (1.0); - for (int i = 0; i < NUM_SAMPLES;) { - interpolation.set_speed (1.0); - result = interpolation.interpolate (0, INTERVAL/10, input + i, output + i); - CPPUNIT_ASSERT_EQUAL ((uint32_t)(INTERVAL/10 * interpolation.speed()), result); - i += result; - } - - for (int i = 0; i < NUM_SAMPLES; i += INTERVAL) { - CPPUNIT_ASSERT_EQUAL (1.0f, output[i+1]); - } -*/ - - cout << "\nSpeed: 0.5"; - for (int i = 0; i < NUM_SAMPLES;) { - interpolation.set_speed (0.5); - //printf ("Interpolate: input: %d, output: %d, i: %d\n", input + i, output + i, i); - result = interpolation.interpolate (0, NUM_SAMPLES - 100, input + i, output + i); - printf ("Result: %d\n", result); - //CPPUNIT_ASSERT_EQUAL ((uint32_t)((NUM_SAMPLES - 100) * interpolation.speed()), result); - //i += result; - break; - } - - for (int i=0; i < NUM_SAMPLES; ++i) { - cout << "input[" << i << "] = " << input[i] << " output[" << i << "] = " << output[i] << endl; - } - - cout << "\nSpeed: 0.2"; - interpolation.set_speed (0.2); - result = interpolation.interpolate (0, NUM_SAMPLES, input, output); - CPPUNIT_ASSERT_EQUAL ((uint32_t)(NUM_SAMPLES * interpolation.speed()), result); - - cout << "\nSpeed: 0.02"; - interpolation.set_speed (0.02); - result = interpolation.interpolate (0, NUM_SAMPLES, input, output); - CPPUNIT_ASSERT_EQUAL ((uint32_t)(NUM_SAMPLES * interpolation.speed()), result); - - cout << "\nSpeed: 0.002"; - interpolation.set_speed (0.002); - result = interpolation.interpolate (0, NUM_SAMPLES, input, output); - CPPUNIT_ASSERT_EQUAL ((uint32_t)(NUM_SAMPLES * interpolation.speed()), result); - - cout << "\nSpeed: 2.0"; - interpolation.set_speed (2.0); - result = interpolation.interpolate (0, NUM_SAMPLES / 2, input, output); - CPPUNIT_ASSERT_EQUAL ((uint32_t)(NUM_SAMPLES / 2 * interpolation.speed()), result); - for (int i = 0; i < NUM_SAMPLES / 2; i += (INTERVAL / interpolation.speed() +0.5)) { - CPPUNIT_ASSERT_EQUAL (1.0f, output[i]); - } - - cout << "\nSpeed: 10.0"; - interpolation.set_speed (10.0); - result = interpolation.interpolate (0, NUM_SAMPLES / 10, input, output); - CPPUNIT_ASSERT_EQUAL ((uint32_t)(NUM_SAMPLES / 10 * interpolation.speed()), result); - for (int i = 0; i < NUM_SAMPLES / 10; i += (INTERVAL / interpolation.speed() +0.5)) { - CPPUNIT_ASSERT_EQUAL (1.0f, output[i]); - } - /* - for (int i=0; i < NUM_SAMPLES; ++i) { - cout << "input[" << i << "] = " << input[i] << " output[" << i << "] = " << output[i] << endl; - } - */ + cout << "\nSpeed: 10.0"; + cubic.set_speed (10.0); + cubic.set_target_speed (cubic.speed()); + result = cubic.interpolate (0, NUM_SAMPLES / 10, input, output); + CPPUNIT_ASSERT_EQUAL ((uint32_t)(NUM_SAMPLES / 10 * cubic.speed()), result); + for (int i = 0; i < NUM_SAMPLES / 10; i += (INTERVAL / cubic.speed() +0.5)) { + CPPUNIT_ASSERT_EQUAL (1.0f, output[i]); + } } diff --git a/libs/ardour/tests/interpolation-test.h b/libs/ardour/tests/interpolation-test.h index 07cc3ab4f7..639c8fc956 100644 --- a/libs/ardour/tests/interpolation-test.h +++ b/libs/ardour/tests/interpolation-test.h @@ -26,9 +26,8 @@ class InterpolationTest : public CppUnit::TestFixture { CPPUNIT_TEST_SUITE(InterpolationTest); - CPPUNIT_TEST(splineInterpolationTest); - //CPPUNIT_TEST(linearInterpolationTest); - //CPPUNIT_TEST(libSamplerateInterpolationTest); + CPPUNIT_TEST(cubicInterpolationTest); + CPPUNIT_TEST(linearInterpolationTest); CPPUNIT_TEST_SUITE_END(); #define NUM_SAMPLES 1000000 @@ -38,8 +37,7 @@ class InterpolationTest : public CppUnit::TestFixture ARDOUR::Sample output[NUM_SAMPLES]; ARDOUR::LinearInterpolation linear; - ARDOUR::SplineInterpolation spline; - ARDOUR::LibSamplerateInterpolation interpolation; + ARDOUR::CubicInterpolation cubic; public: @@ -53,15 +51,12 @@ class InterpolationTest : public CppUnit::TestFixture output[i] = 0.0f; } linear.add_channel_to (NUM_SAMPLES, NUM_SAMPLES); - spline.add_channel_to (NUM_SAMPLES, NUM_SAMPLES); - interpolation.add_channel_to (NUM_SAMPLES, NUM_SAMPLES); + cubic.add_channel_to (NUM_SAMPLES, NUM_SAMPLES); } void tearDown() { } void linearInterpolationTest(); - void splineInterpolationTest(); - void libSamplerateInterpolationTest(); - + void cubicInterpolationTest(); }; -- cgit v1.2.3