bezier.cpp

Go to the documentation of this file.

/*
 * bezier.cpp
 *
 * Copyright (C) 2003 Thomas A. Vaughan
 * All rights reserved.
 *
 *
 * Redistribution and use in source and binary forms, with or without
 * modification, are permitted provided that the following conditions are met:
 *     * Redistributions of source code must retain the above copyright
 *       notice, this list of conditions and the following disclaimer.
 *     * Redistributions in binary form must reproduce the above copyright
 *       notice, this list of conditions and the following disclaimer in the
 *       documentation and/or other materials provided with the distribution.
 *     * Neither the name of the <organization> nor the
 *       names of its contributors may be used to endorse or promote products
 *       derived from this software without specific prior written permission.
 *
 * THIS SOFTWARE IS PROVIDED BY THOMAS A. VAUGHAN ''AS IS'' AND ANY
 * EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED
 * WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
 * DISCLAIMED. IN NO EVENT SHALL THOMAS A. VAUGHAN BE LIABLE FOR ANY
 * DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES
 * (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
 * LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND
 * ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
 * (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
 * SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
 *
 *
 * simple bezier APIs
 */

// includes --------------------------------------------------------------------
#include "bezier.h"             // include our own header first

#include <sstream>

#include <math.h>
#include <stdio.h>

#include "common/wave_ex.h"
#include "perf/perf.h"
#include "util/token_stream.h"


namespace bezier {


static const float s_eps                        = 1.0e-8;

/*
 * static helper methods
 */

/*
 * finds roots of 3at^2 + 2bt + c = 0
 *
 * roots t(i) must be in the range 0 < t < 1
 *
 * TODO: should generalize this into a general quadratic root finder
 */
static int
findRoots
(
IN float a,
IN float b,
IN float c,
OUT float * roots
)
throw()
{
        perf::Timer timer("bezier::findRoots");

        int nroots = 0;

        ASSERT(roots, "null roots array");

        float t;

        if (a < s_eps && a > -s_eps) {
                if (b < s_eps && b > -s_eps) {
                        // DPRINTF("no roots: a=b=0");
                        return 0;
                }
                t = (-1.0 * c) / (2.0 * b);
                if (t > 0.0 && t < 1.0) {
                        // DPRINTF("one root: a=0, t = %f", t);
                        roots[0] = t;
                        ++nroots;
                }
                return nroots;
        }

        float factor = b * b - 3.0 * a * c;
        if (factor < 0.0) {
                // DPRINTF("no roots--imaginary");
                return 0;       // imaginary roots
        }

        factor = sqrt(factor);
        t = (1.0 / (3.0 * a)) * (factor - b);
        if (t > 0.0 && t < 1.0) {
                // DPRINTF("Found a root: %f", t);
                roots[nroots] = t;
                ++nroots;
        }

        // skip degenerate roots
        if (0.0 == factor)
                return nroots;

        t = (1.0 / (3.0 * a)) * (-factor - b);
        if (t > 0.0 && t < 1.0) {
                // DPRINTF("Found a root: %f", t);
                roots[nroots] = t;
                ++nroots;
        }
        return nroots;
}



static void
checkBoundaries
(
IO float& min,
IO float& max,
IN float value
)
throw()
{
        if (value < min)
                min = value;
        if (value > max)
                max = value;
}



/*
 * getRect() -- gets the bounding rectangle for the bezier.
 *
 * solutions for maxima are
 *
 *  t = (1/3a) * (-b +/- sqrt(b^2 - 3ac))
 *
 *   1) see if roots exist
 *   2) see if any roots are between 0 and 1
 *   3) maxima are then endpoints + valid roots
 */
void
getRect
(
IN const curve_t& curve,
OUT rect_t& rect
)
throw()
{
        perf::Timer timer("bezier::getRect");

        // see if we have any roots for x,y
        float tx[3], ty[3];
        int xroots = findRoots(curve.ax, curve.bx, curve.cx, tx);
        int yroots = findRoots(curve.ay, curve.by, curve.cy, ty);

        // DPRINTF("Found %d x-roots...", xroots);
        //DPRINTF("Found %d y-roots...", yroots);

        // add endpoint as another root
        tx[xroots++] = 1.0;
        ty[yroots++] = 1.0;
        ASSERT(xroots <= 3, "Bad xroot count: %d", xroots);
        ASSERT(yroots <= 3, "Bad yroot count: %d", yroots);

        // check roots for boundary
        float xmin = curve.x0;
        float xmax = curve.x0;
        float ymin = curve.y0;
        float ymax = curve.y0;

        // check roots
        for (int i = 0; i < xroots; ++ i) {
                float x = curve.x0 + curve.cx * tx[i] +
                    curve.bx * tx[i] * tx[i] +
                    curve.ax * tx[i] * tx[i] * tx[i];
                //DPRINTF("x-root %d: %f", i, x);
                checkBoundaries(xmin, xmax, x);
        }
        for (int i = 0; i < yroots; ++i) {
                float y = curve.y0 + curve.cy * ty[i] +
                    curve.by * ty[i] * ty[i] +
                    curve.ay * ty[i] * ty[i] * ty[i];
                //DPRINTF("y-root: %d: %f", i, y);
                checkBoundaries(ymin, ymax, y);
        }

        // return bounding box
        rect.left = xmin;
        rect.right = xmax;
        rect.top = ymin;
        rect.bottom = ymax;
}



void
scale
(
IO curve_t& curve,
IN float factor
)
throw()
{
        curve.ax *= factor;
        curve.bx *= factor;
        curve.cx *= factor;

        curve.ay *= factor;
        curve.by *= factor;
        curve.cy *= factor;
}



static float
writeVal
(
IN float x,
IN float eps
)
throw()
{
        if (x < eps && x > -eps) {
                // within +/- epsilon -- effectively zero
                return 0.0;
        }
        return x;
}



void
write
(
IN const curve_t& curve,
OUT std::string& value,
IN float eps
)
throw()
{
        ASSERT(eps >= 0.0, "Bad epsilon value: %f", eps);
        value = "";

        std::ostringstream out;
        out << writeVal(curve.x0, eps) << " ";
        out << writeVal(curve.ax, eps) << " ";
        out << writeVal(curve.bx, eps) << " ";
        out << writeVal(curve.cx, eps) << " ";

        out << writeVal(curve.y0, eps) << " ";
        out << writeVal(curve.ay, eps) << " ";
        out << writeVal(curve.by, eps) << " ";
        out << writeVal(curve.cy, eps);

        value = out.str();
}


#define READ_COMPONENT(x)                                       \
        if (!*p) {                                              \
                WAVE_EX(wex);                                   \
                wex << "Could not read bezier component" << #x ; \
        }                                                       \
        curve.x = atof(p);                                      \
        while (*p && !isspace(*p)) { ++p; }                     \
        while (*p && isspace(*p)) { ++p; }


void
read
(
IN const char * value,
OUT curve_t& curve
)
{
        // read all components in order x0, ax, bx, cx, y0, ay, by, cy
        // (see write() routine above)

        const char * p = value;
        READ_COMPONENT(x0)
        READ_COMPONENT(ax)
        READ_COMPONENT(bx)
        READ_COMPONENT(cx)
        READ_COMPONENT(y0)
        READ_COMPONENT(ay)
        READ_COMPONENT(by)
        READ_COMPONENT(cy)
}



void getCurveFromControlPoints
(
IN const control_points_t& cp,
OUT curve_t& b
)
throw()
{
        // initial point
        b.x0 = cp.p0.x;
        b.y0 = cp.p0.y;

        // linear with t
        b.cx = 3.0 * (cp.p1.x - cp.p0.x);
        b.cy = 3.0 * (cp.p1.y - cp.p0.y);

        // quadradic with t
        b.bx = 3.0 * (cp.p0.x - 2.0 * cp.p1.x + cp.p2.x);
        b.by = 3.0 * (cp.p0.y - 2.0 * cp.p1.y + cp.p2.y);

        // cubic with t
        b.ax = -cp.p0.x + 3.0 * cp.p1.x - 3.0 * cp.p2.x + cp.p3.x;
        b.ay = -cp.p0.y + 3.0 * cp.p1.y - 3.0 * cp.p2.y + cp.p3.y;
}



void
getControlPointsFromCurve
(
IN const curve_t& b,
OUT control_points_t& cp
)
throw()
{
        // see 8 march 07 journal entry
        cp.p0.set(b.x0, b.y0);

        float inv3 = 1.0 / 3.0;

        cp.p1.set(b.x0 + inv3 * b.cx,
                  b.y0 + inv3 * b.cy);

        cp.p2.set(inv3 * (b.bx + 2.0 * b.cx + 3.0 * b.x0),
                  inv3 * (b.by + 2.0 * b.cy + 3.0 * b.y0));

        cp.p3.set(b.ax + b.bx + b.cx + b.x0,
                  b.ay + b.by + b.cy + b.y0);
}



};      // bezier namespace