class Spline
{
  float[] m_x;
  float[] m_y;

  float[] m_a;

  float[] m_u;
  float[] m_d;
  float[] m_l;

  Cubic[] m_funcs;

  Spline( float[] x, float[] y )
  {
    m_x = x;
    m_y = y;

    m_a = new float[ m_x.length ];
    m_a[0] = m_y[1] - m_y[0];
    for ( int i = 1; i < ( m_x.length - 1 ); ++i )
    {
      m_a[i] = 3.0 * (m_y[i+1] - m_y[i-1]);
    }
    m_a[m_x.length - 1] = 3.0 * (m_y[ m_x.length-1 ] - m_y[ m_x.length-2 ] );

    m_u = new float[ m_x.length ];
    m_d = new float[ m_x.length ];
    m_l = new float[ m_x.length ];

    m_l[0] = 0; 
    m_d[0] = 2; 
    m_u[0] = 1;
    for ( int i = 0; i < m_x.length-1; ++i )
    {
      m_l[i] = 1; 
      m_d[i] = 4; 
      m_u[i] = 1;
    }
    m_l[ m_x.length-1 ] = 1; 
    m_d[ m_x.length-1] = 4; 
    m_u[ m_x.length-1] = 0;

    float beta;
    for ( int i = 1; i < m_x.length; ++i )
    {
      beta = ( -1.0 / m_d[i-1] ) * m_l[i];
      m_l[i] = 0;
      m_d[i] += m_u[i-1] * beta;
      m_a[i] += m_a[i-1] * beta;
    }

    m_a[ m_x.length-1 ] /= m_d[ m_x.length-1 ];
    m_d[ m_x.length-1 ] = 1;
    for ( int i = m_x.length-2; i >= 0; --i )
    {
      beta = ( -1.0 * m_u[i] );
      m_u[i] = 0;
      m_a[i] += (m_a[i+1] * beta);
      m_a[i] /= m_d[i];
      m_d[i] = 1.0;
    }

    m_funcs = new Cubic[ m_x.length - 1 ];
    float a, b, c, d;
    for ( int i = 0; i < m_funcs.length; ++i )
    {
      a = m_y[i];
      b = m_a[i];
      c = 3 * ( m_y[i+1] - m_y[i] ) - 2 * m_a[i] - m_a[i+1];
      d = 2 * ( m_y[i] - m_y[ i+1 ] ) + m_a[i] + m_a[ i+1 ];
      m_funcs[i] = new Cubic( a, b, c, d, m_x[i], m_x[i+1] );
    }
  }

  float evaluate( float x )
  {
    if ( x < m_x[0] ) return m_y[0];
    if ( x > m_x[m_x.length-1] ) return m_y[ m_x.length-1 ];
    int i = 0;
    while( x > m_x[i+1] ) ++i;
    return m_funcs[i].evaluate( x );
  }
}

class Cubic
{
  float m_a;
  float m_b;
  float m_c;
  float m_d;

  float m_x0;
  float m_x1;
  float m_dist;

  Cubic( float a, float b, float c, float d, float x0, float x1 )
  {
    m_a = a;
    m_b = b;
    m_c = c;
    m_d = d;
    m_x0 = x0;
    m_x1 = x1;
    m_dist = 1.0/(m_x1 - m_x0);
  }

  boolean in( float x )
  {
    return ( (x >= m_x0) && (x <= m_x1) );
  }

  float evaluate( float x )
  {
    float t = (x - m_x0) * m_dist;
    float y = m_a + t * ( m_b + t * ( m_c + t * m_d ) );
    return y;
  }

}
  

float[] t = { 0, 1, 2, 3, 4, 5, 6, 7, 8, 9 };
float[] x;
float[] y;

Spline xSpline;
Spline ySpline;

void setup()
{
  size( 400, 400 );
  x = new float[ t.length ];
  y = new float[ t.length ];
  for ( int i = 0; i < t.length; ++i )
  {
    x[i] = random( 0, width );
    y[i] = random( 0, height );
  }
 
  xSpline = new Spline( t, x );
  ySpline = new Spline( t, y ); 
}

void draw()
{
  background(0);
  stroke( 255 );
  strokeWeight( 1 );
  for ( float i = t[0]; i < t[t.length-1]; i += 0.01 )
  {
    line( xSpline.evaluate(i), ySpline.evaluate(i), xSpline.evaluate(i+0.01), ySpline.evaluate(i+0.01) );
  }
  
  stroke( 100, 240, 230 );
  strokeWeight( 5 );
  for ( int i = 0; i < t.length; ++i )
  {
    point( x[i], y[i] );
  }

}