Simple bezier curves

by Erik Möller, November 2009

When I wrote the editor for Emberwind it didn't take long until we realized placing lines of pickups as individual objects quickly became a tedious maintenance task as the levels were playtested and refined. We also needed to simply define mortar trajectories and so I started to look at adding some curve editing support to the editor.
editor curve

My first thought was bezier curves as that's pretty much standard nowadays. I wanted something simple, I wanted to easily calculate the arc length, get a point at a certain length along the curve and also be able to get the normal at a specified point so I could make "wide curves" of pickups.
I stumbled upon a neat approach [1] to generate quadratic uniform b-spline curves presented by George Chaikin in 1974. The idea is exremely simple, start with your control polygon, for each non endpoint CV create two new CVs at 1/3 the distance to the neighboring CV and discard the original CV. Rinse and repeat until it's smooth enough for your purpose. Sure, bezier curves represented in some analytical form are cool and elegant, but the simplicity of this approach outweighs it all for my purposes.

This little animation shows how a curve is getting progressively smoother (using canvas if you're using any of the standardized browsers.. ie. not Internet Explorer). If you haven't seen the light yet I've supplied a gif ;)


My implementation is a namespace containing a few template functions to do the required operations on a CV list. The actual CV is a template parameter to easily work with any 2d vector class. It's been compiled on VC2005 and GCC (as well as being ported to javascript for the illustration above).

[1] CHAIKIN, G.
An algorithm for high speed curve generation.
Computer Graphics and Image Processing 3 (1974), 346-349.


Source listing chaikin_curve.h

			      
#ifndef CHAIKIN_CURVE_H
#define CHAIKIN_CURVE_H

#include <vector>

namespace ChaikinCurve
{

  template <typename CV>
  void Subdivide(std::vector<CV> &handles, int subdivs)
  {
    if (handles.empty()) return;
    do
    {
      size_t numHandles = handles.size();

      // keep the first point
      handles.push_back(handles.front());
      for (unsigned int i = 0; i < numHandles - 1; ++i) 
      {
        // get 2 original points
        const CV& p0 = handles[i];
        const CV& p1 = handles[i + 1];
        // average the 2 original points to create 2 new points. For each
        // CV, another 2 verts are created.
        CV Q(0.75f * p0.x + 0.25f * p1.x, 0.75f * p0.y + 0.25f * p1.y);
        CV R(0.25f * p0.x + 0.75f * p1.x, 0.25f * p0.y + 0.75f * p1.y);

        handles.push_back(Q);
        handles.push_back(R);
      }
      // keep the last point
      handles.push_back(handles[numHandles - 1]);

      // update the points array
      handles.erase(handles.begin(), handles.begin() + numHandles);
    } while (--subdivs > 0);
  }

  template <typename CV>
  float GetLength(const std::vector<CV> &points)
  {
    float len = 0.0f;
    // gcc and msvc doesn't agree on the syntax, gcc needs the extra "typename" to 
    // parse it properly.
    if (points.size() >= 2)
    for (typename std::vector<CV>::const_iterator i = points.begin() + 1; i != points.end(); ++i)
    {
      CV diff(*i - *(i - 1)); 
      len += sqrtf((float)(diff.x * diff.x + diff.y * diff.y));
    }
    return len;
  }

  template <typename CV>
  CV GetPointAtLength(const std::vector<CV> &points, float len)
  {
    if (points.empty()) return CV(0,0);
    if (points.size() == 1) return points.front();
    for (int index = 0; index != (int)points.size() - 1; ++index)
    {
      CV diff(points[index + 1] - points[index]);
      float segLen = sqrtf((float)(diff.x * diff.x + diff.y * diff.y));
      if (segLen > len)
        return CV(points[index].x + diff.x * len / segLen,  points[index].y + diff.y * len / segLen);
      else
        len -= segLen;
    }
    return points.back();
  }

  template <typename CV>
  CV GetDirAtParam(const std::vector<CV> &points, float param)
  {
    if (points.size() < 2) return CV(0,0);
    float totalLen = GetLength(points);
    float tgtLen = param * totalLen;
    for (int index = 0; index != (int)points.size() - 1; ++index)
    {
      CV diff(points[index + 1] - points[index]);
      float segLen = sqrtf((float)(diff.x * diff.x + diff.y * diff.y));
      if (segLen > tgtLen)
        return diff;
      else
        tgtLen -= segLen;
    }
    return points[points.size()-1] - points[points.size()-2];
  }

  template <typename CV>
  float GetDXAtLength(const std::vector<CV> &points, float len)
  {
    if (points.size() < 2) return 0;
    for (int index = 0; index != (int)points.size() - 1; ++index)
    {
      CV diff(points[index + 1] - points[index]);
      float segLen = sqrtf((float)(diff.x * diff.x + diff.y * diff.y));
      if (segLen > len)
        return diff.x / fabsf((float)diff.y);
      else
        len -= segLen;
    }
    return 0;
  }

  template <typename CV>
  void GetEvenlySpacedPoints(const std::vector<CV> &handles, int count, std::vector<CV> &points, std::vector<CV> *normals)
  {
    std::vector<CV> tmp;
    tmp = handles;
    Subdivide(tmp, 3);
    float len = GetLength(tmp);
    float spacing = len / (count - 1);
    points.push_back(tmp.front());
    if (normals)
    {
      CV dir = GetDirAtParam(tmp, 0);
      normals->push_back(CV(-dir.y, dir.x));
    }
    for (int i = 1; i < count - 1; ++i)
    {
      points.push_back(GetPointAtLength(tmp, i * spacing));
      if (normals)
      {
        CV dir = GetDirAtParam(tmp, i * spacing / len);
        normals->push_back(CV(-dir.y, dir.x));
      }
    }
    points.push_back(tmp.back());
    if (normals)
    {
      CV dir = GetDirAtParam(tmp, 1);
      normals->push_back(CV(-dir.y, dir.x));
    }
  }

}

#endif //CHAIKIN_CURVE_H