views:

97

answers:

3

I've got a vector of samples that form a curve. Let's imagine there are 1000 points in it. If I want to stretch it to fill 1500 points, what is the simplest algorithm that gives decent results? I'm looking for something that is just a few lines of C/C++.

I'll always want to increase the size of the vector, and the new vector can be anywhere from 1.1x to 50x the size of the current vector.

Thanks!

+2  A: 

what is the simplest algorithm that gives decent results?

Catmull-Rom splines. (if you want a smooth curve)

http://www.mvps.org/directx/articles/catmull/
http://en.wikipedia.org/wiki/Cubic_Hermite_spline

For each new item calculate fractional position in old array, use use fractional part (f - floor(f)) as interpolation factor, and "integer" (i.e. floor(f)) part to find nearest elements.

That is assuming that you're operating on data that can be mathematically interpolated (floats). If data cannot be interpolated (strings), then the only solution is to use nearest available element of old array.

You'll need some tweaking if points in array aren't evenly distributed.

SigTerm
A: 

Simplest option I can think of is just a fn that expands the array based on mean averages, so:

x,y,z

becomes

x, avg(x,y), y, avg (y,z), z

If you need more data points, just run it multiple times on the vector.

mikurski
..And if new array isn't 2x times larger than old one, results will be incorrect. It isn't much better than linear interpolation - you won't get smooth curve this way.
SigTerm
+1  A: 

Here's C++ for linear and quadratic interpolation.
interp1( 5.3, a, n ) is a[5] + .3 * (a[6] - a[5]), .3 of the way from a[5] to a[6];
interp1array( a, 1000, b, 1500 ) would stretch a to b.
interp2( 5.3, a, n ) draws a parabola through the 3 nearest points a[4] a[5] a[6]: smoother than interp1 but still fast.
(Splines use 4 nearest points, smoother yet; if you read python, see basic-spline-interpolation-in-a-few-lines-of-numpy.

// linear, quadratic interpolation in arrays
// from interpol.py denis 2010-07-23 July

#include <stdio.h>
#include <stdlib.h>

    // linear interpolate x in an array
// inline
float interp1( float x, float a[], int n )
{
    if( x <= 0 )  return a[0];
    if( x >= n - 1 )  return a[n-1];
    int j = int(x);
    return a[j] + (x - j) * (a[j+1] - a[j]);
}

    // linear interpolate array a[] -> array b[]
void inter1parray( float a[], int n, float b[], int m )
{
    float step = float( n - 1 ) / (m - 1);
    for( int j = 0; j < m; j ++ ){
        b[j] = interp1( j*step, a, n );
    }
}

//..............................................................................
    // parabola through 3 points, -1 < x < 1
float parabola( float x, float f_1, float f0, float f1 )
{
    if( x <= -1 )  return f_1; 
    if( x >= 1 )  return f1; 
    float l = f0 - x * (f_1 - f0);
    float r = f0 + x * (f1 - f0);
    return (l + r + x * (r - l)) / 2;
}

    // quadratic interpolate x in an array
float interp2( float x, float a[], int n )
{
    if( x <= .5  ||  x >= n - 1.5 )
        return interp1( x, a, n );
    int j = int( x + .5 );
    float t = 2 * (x - j);  // -1 .. 1
    return parabola( t, (a[j-1] + a[j]) / 2, a[j], (a[j] + a[j+1]) / 2 );
}

    // quadratic interpolate array a[] -> array b[]
void interp2array( float a[], int n, float b[], int m )
{
    float step = float( n - 1 ) / (m - 1);
    for( int j = 0; j < m; j ++ ){
        b[j] = interp2( j*step, a, n );
    }
}

int main( int argc, char* argv[] )
{
        // a.out [n m] --
    int n = 10, m = 100;
    int *ns[] = { &n, &m, 0 },
        **np = ns;
    char* arg;
    for( argv ++;  (arg = *argv) && *np;  argv ++, np ++ )
        **np = atoi( arg );
    printf( "n: %d  m: %d\n", n, m );

    float a[n], b[m];
    for( int j = 0; j < n; j ++ ){
        a[j] = j * j;
    }
    interp2array( a, n, b, m );  // a[] -> b[]

    for( int j = 0; j < m; j ++ ){
        printf( "%.1f ", b[j] );
    }
    printf( "\n" );
}
Denis
I haven't thoroughly tested this code, so I can't vouch for its correctness, but this is exactly the kind of answer I was looking for. Thanks!
twk
You're welcome. If it works, tell the boss; if not, tell me
Denis