Synopsis

#include <boost/math/interpolators/catmull_rom.hpp>

namespace boost{ namespace math{

    template<class Point>
    class catmull_rom
    {
    public:
        catmull_rom(const Point* const points, size_t num_points, bool closed = false, Real alpha = (Real) 1/ (Real) 2)

        catmull_rom(std::vector<Point>&& points, bool closed = false, Real alpha = (Real) 1/ (Real) 2)

        Real operator()(Real s) const;

        Real max_parameter() const;

        Real parameter_at_point(size_t i) const;

        Point prime(Real s) const;
    };

}}

Description

Catmull-Rom splines are a family of interpolating curves which are commonly used in computer graphics and animation. Catmull-Rom splines enjoy the following properties:

The catmull_rom class provided by boost creates a cubic Catmull-Rom spline from an array of points in any dimension. Since there are numerous ways to represent a point in n-dimensional space, the class attempts to be flexible by templating on the point type. The requirements on the point type are discussing in more detail below, but roughly, it must have a dereference operator defined (i.e., p[0] is not a syntax error), it must be able to be dereferenced up to dimension -1, and p[i] is of type Real, define value_type, and the free function size(). These requirements are met by std::vector and std::array. The basic usage is shown here:

std::vector<std::array<Real, 3>> points(4);
points[0] = {0,0,0};
points[1] = {1,0,0};
points[2] = {0,1,0};
points[3] = {0,0,1};
catmull_rom<std::array<Real, 3>> cr(points.data(), points.size());
// Interpolate at s = 0.1:
auto point = cr(0.1);

Memory conscious programmers may enjoy using the move constructor instead:

catmull_rom<std::array<Real, 3>> cr(std::move(points));

The spline can be either open or closed, closed meaning that there is some t such that P(t) = P(0). The default is open, but this can be easily changed:

// closed = true
catmull_rom<std::array<Real, 3>> cr(points.data(), points.size(), true);

Inside catmull_rom, the Catmull-Rom curve is represented as closed. This is because an open Catmull-Rom curve is implicitly closed, but the closing point is the zero vector. There is no reason to suppose that the zero vector is a better closing point than the endpoint (or any other point, for that matter), so traditionally Catmull-Rom splines leave the segment between the first and second point undefined, as well as the segment between the second-to-last and last point. We find this property of the traditional implementation of Catmull-Rom splines annoying and confusing to the user. Hence internally, we close the curve so that the first and last segments are defined. Of course, this causes the tangent vectors to the first and last points to be bizarre. This is a "pick your poison" design decision-either the curve cannot interpolate in its first and last segments, or the tangents along the first and last segments are meaningless.

Since the routine internally represents the curve as closed, a question arises: Why does the user have to specify if the curve is open or closed? The answer is that the parameterization is chosen by the routine, so it is of interest to the user to understand the values where a meaningful result is returned.

Real max_s = cr.max_parameter();

If you attempt to interpolate for s > max_s, an exception is thrown. If the curve is closed, then cr(max_s) = p0, where p0 is the first point on the curve. If the curve is open, then cr(max_s) = pf, where pf is the final point on the curve.

The Catmull-Rom curve admits an infinite number of parameterizations. The default parameterization of the catmull_rom class is the so-called centripedal parameterization. This parameterization has been shown to be the only parameterization that does not form cusps or self-intersections within segments. However, for advanced users, other parameterizations can be chosen using the alpha parameter:

// alpha = 1 is the "chordal" parameterization.
catmull_rom<std::array<double, 3>> cr(points.data(), points.size(), false, 1.0);

Finally, the tangent vector to any point of the curve can be computed via

double s = 0.1;
Point tangent = cr.prime(s);

Since the magnitude of the tangent vector is dependent on the parameterization, it is not as meaningful as (say) arc-length parameterization. However, its direction is meaningful, so the user may wish to normalize this result.

Examples

Performance

The following performance numbers were generated for a call to the Catmull-Rom interpolation method. The number that follows the slash is the number of points passed to the interpolant.

Run on 2700 MHz CPU
CPU Caches:
  L1 Data 32K (x2)
  L1 Instruction 32K (x2)
  L2 Unified 262K (x2)
  L3 Unified 3145K (x1)
---------------------------------------------------------
Benchmark                              Time           CPU
---------------------------------------------------------
BM_CatmullRom<double>/4               20 ns         20 ns
BM_CatmullRom<double>/8               21 ns         21 ns
BM_CatmullRom<double>/16              23 ns         23 ns
BM_CatmullRom<double>/32              24 ns         24 ns
BM_CatmullRom<double>/64              27 ns         27 ns
BM_CatmullRom<double>/128             27 ns         27 ns
BM_CatmullRom<double>/256             30 ns         30 ns
BM_CatmullRom<double>/512             32 ns         31 ns
BM_CatmullRom<double>/1024            33 ns         33 ns
BM_CatmullRom<double>/2048            34 ns         34 ns
BM_CatmullRom<double>/4096            36 ns         36 ns
BM_CatmullRom<double>/8192            38 ns         38 ns
BM_CatmullRom<double>/16384           39 ns         39 ns
BM_CatmullRom<double>/32768           40 ns         40 ns
BM_CatmullRom<double>/65536           45 ns         44 ns
BM_CatmullRom<double>/131072          46 ns         46 ns
BM_CatmullRom<double>/262144          50 ns         50 ns
BM_CatmullRom<double>/524288          53 ns         52 ns
BM_CatmullRom<double>/1048576         58 ns         57 ns
BM_CatmullRom<double>_BigO          2.97 lgN       2.97 lgN
BM_CatmullRom<double>_RMS             19 %         19 %