Biblioteca de código abierto de C++ para el ajuste de curvas

Se requieren cuatro puntos para describir de manera única una curva cúbica (el primer artículo que ha vinculado cubre ese caso). Tiene más de cuatro puntos, por lo que es poco probable que obtenga un ajuste perfecto: se requerirá algún tipo de compromiso o compensación. ¡Bienvenido al arte negro de la optimización numérica!

El segundo artículo que ha vinculado presenta este tema usando C# y Math.NET.

Si está feliz de tratar el optimizador de propósito general como una caja negra, entonces prefiero Excel Solver (y un buen tutorial) para aprender sobre lo que puede hacer un optimizador. Básicamente definimos una función, luego el optimizador la llama repetidamente con diferentes parámetros intentando obtener el valor de retorno más bajo posible.

Recomiendo dlib . Aquí uso dlib para ajustar algunos puntos a una curva Bezier. Logré compilarlo con bastante facilidad en Visual Studio 2013, pero probablemente necesitará ajustarlo para que funcione con un compilador diferente. Las partes de dlib que se usan aquí están todas en los encabezados, no hay ningún archivo .lib al que deba vincularse.

#include <stdio.h>
#include "../dlib/optimization.h"

typedef dlib::matrix<double,0,1> column_vector;

struct Point
{
    double x;
    double y;
};

static Point points[] =
{
    { 254 , 102 },
    { 226 , 63  },
    { 185 , 49  },
    { 146 , 74  },
    { 142 , 119 },
    { 117 , 169 },
    { 86  , 214 },
    { 40  , 200 },
};

//
// cubicBezier
//
// p1 - start point
// c1 - first control point
// c2 - second control point
// p2 - end point
//
double cubicBezier(double p1, double c1, double c2, double p2, double t)
{
    double s = (1 - t);

    double v = 0;
    v = v + (1 * p1 * s * s * s);
    v = v + (3 * c1 * s * s * t);
    v = v + (3 * c2 * s * t * t);
    v = v + (1 * p2 * t * t * t);

    return v;        
};

Point cubicBezier(double p1x, double p1y, double c1x, double c1y,
                  double c2x, double c2y, double p2x, double p2y,
                  double t)
{
    Point pt;
    pt.x = cubicBezier(p1x, c1x, c2x, p2x, t);
    pt.y = cubicBezier(p1y, c1y, c2y, p2y, t);
    return pt;
}

// Any distance function can be used for optimisation.  This one, where we want
// to find the least-squares is most common. 
double dist(double x1, double y1, double x2, double y2)
{
    double x = x2 - x1;
    double y = y2 - y1;

    return (x * x) + (y * y);
}

// This is the function that the optimiser calls repeatedly with different
// parameters, attempting to get the lowest return value it can.
double rateCurve(const column_vector& params)
{
    double p1x = points[0].x;
    double p1y = points[0].y;
    double c1x = params(0,0);
    double c1y = params(1,0);
    double c2x = params(2,0);
    double c2y = params(3,0);
    double p2x = points[7].x;
    double p2y = points[7].y;

    double distances = 0;

    for (Point target : points)
    {
        double distance = _DMAX;

        for (double t = 0; t <= 1; t += 0.02)
        {
            Point pt = cubicBezier( p1x, p1y, c1x, c1y, c2x, c2y, p2x, p2y, t);

            double dCandidate = dist(pt.x, pt.y, target.x, target.y);

            distance = std::min(distance, dCandidate);
        }

        distances += distance;
    }

    // Thats the curve-fitting done.  Now incentivise slightly smoother curves.
    double p1c1 = dist(p1x, p1y, c1x, c1y);
    double p2c2 = dist(p2x, p2y, c2x, c2y);

    return distances + pow(p1c1, 0.6) + pow(p2c2, 0.6);
}

int main(int argc, char* argv[])
{
    column_vector params(4);
    params = points[7].x, points[7].y, points[0].x, points[0].y;

    dlib::find_min_using_approximate_derivatives(
            dlib::cg_search_strategy(),
            dlib::objective_delta_stop_strategy(1).be_verbose(),
            rateCurve,
            params,
            -1);

    printf("p1x = %f;\n", points[0].x);
    printf("p1y = %f;\n", points[0].y);
    printf("c1x = %f;\n", params(0,0));
    printf("c1y = %f;\n", params(1,0));
    printf("c2x = %f;\n", params(2,0));
    printf("c2y = %f;\n", params(3,0));
    printf("p2x = %f;\n", points[7].x);
    printf("p2y = %f;\n", points[7].y);

    return 0;
}

Aquí está la salida:

iteration: 0   objective: 9740.42
iteration: 1   objective: 4880.37
iteration: 2   objective: 2872.77
iteration: 3   objective: 2523.82
iteration: 4   objective: 2048.95
iteration: 5   objective: 1680.86
iteration: 6   objective: 1519.74
iteration: 7   objective: 1366.39
iteration: 8   objective: 1330.56
iteration: 9   objective: 1285.79
iteration: 10   objective: 1275.27
iteration: 11   objective: 1274.82
p1x = 254.000000;
p1y = 102.000000;
c1x = 127.524342;
c1y = -86.849427;
c2x = 146.034795;
c2y = 283.099363;
p2x = 40.000000;
p2y = 200.000000;

Y aquí hay una visualización que muestra mis puntos y el intento de ajustarlos:

Graphics Gems tiene un ejemplo de código C simple de ajuste de curva Bezier sin otras dependencias de biblioteca: https://github.com/erich666/GraphicsGems/blob/master/gems/FitCurves.c

(El código es de dominio público; consulte el archivo Léame en el repositorio).