Adding upstream version 1.3.upstream/1.3upstream

Signed-off-by: Daniel Baumann <daniel.baumann@progress-linux.org>

1 files changed, 680 insertions, 0 deletions

diff --git a/tests/bezier-test.cpp b/tests/bezier-test.cpp
new file mode 100644
index 0000000..0799393
--- /dev/null
+++ b/tests/bezier-test.cpp
@@ -0,0 +1,680 @@
+/** @file
+ * @brief Unit tests for Affine.
+ * Uses the Google Testing Framework
+ *//*
+ * Authors:
+ *   Nathan Hurst <njh@njhurst.com>
+ *   Krzysztof Kosiński <tweenk.pl@gmail.com>
+ *   Johan Engelen <j.b.c.engelen@alumnus.utwente.nl>
+ *
+ * Copyright 2010 Authors
+ *
+ * This library is free software; you can redistribute it and/or
+ * modify it either under the terms of the GNU Lesser General Public
+ * License version 2.1 as published by the Free Software Foundation
+ * (the "LGPL") or, at your option, under the terms of the Mozilla
+ * Public License Version 1.1 (the "MPL"). If you do not alter this
+ * notice, a recipient may use your version of this file under either
+ * the MPL or the LGPL.
+ *
+ * You should have received a copy of the LGPL along with this library
+ * in the file COPYING-LGPL-2.1; if not, write to the Free Software
+ * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
+ * You should have received a copy of the MPL along with this library
+ * in the file COPYING-MPL-1.1
+ *
+ * The contents of this file are subject to the Mozilla Public License
+ * Version 1.1 (the "License"); you may not use this file except in
+ * compliance with the License. You may obtain a copy of the License at
+ * http://www.mozilla.org/MPL/
+ *
+ * This software is distributed on an "AS IS" basis, WITHOUT WARRANTY
+ * OF ANY KIND, either express or implied. See the LGPL or the MPL for
+ * the specific language governing rights and limitations.
+ */
+
+#include "testing.h"
+#include <iostream>
+
+#include <2geom/bezier.h>
+#include <2geom/polynomial.h>
+#include <2geom/basic-intersection.h>
+#include <2geom/bezier-curve.h>
+#include <vector>
+#include <iterator>
+#include <glib.h>
+
+using std::vector, std::min, std::max;
+using namespace Geom;
+
+Poly lin_poly(double a, double b) { // ax + b
+    Poly p;
+    p.push_back(b);
+    p.push_back(a);
+    return p;
+}
+
+bool are_equal(Bezier A, Bezier B) {
+    int maxSize = max(A.size(), B.size());
+    double t = 0., dt = 1./maxSize;
+
+    for(int i = 0; i <= maxSize; i++) {
+        EXPECT_FLOAT_EQ(A.valueAt(t), B.valueAt(t));// return false;
+        t += dt;
+    }
+    return true;
+}
+
+class BezierTest : public ::testing::Test {
+protected:
+
+    BezierTest()
+        : zero(fragments[0])
+        , unit(fragments[1])
+        , hump(fragments[2])
+        , wiggle(fragments[3])
+    {
+        zero = Bezier(0.0,0.0);
+        unit = Bezier(0.0,1.0);
+        hump = Bezier(0,1,0);
+        wiggle = Bezier(0,1,-2,3);
+    }
+
+    Bezier fragments[4];
+    Bezier &zero, &unit, &hump, &wiggle;
+};
+
+TEST_F(BezierTest, Basics) {
+
+    //std::cout << unit <<std::endl;
+    //std::cout << hump <<std::endl;
+
+    EXPECT_TRUE(Bezier(0,0,0,0).isZero());
+    EXPECT_TRUE(Bezier(0,1,2,3).isFinite());
+
+    EXPECT_EQ(3u, Bezier(0,2,4,5).order());
+
+    ///cout << " Bezier::Bezier(const Bezier& b);\n";
+    //cout << Bezier(wiggle) << " == " << wiggle << endl;
+
+    //cout << "explicit Bezier(unsigned ord);\n";
+    //cout << Bezier(10) << endl;
+
+    //cout << "Bezier(Coord c0, Coord c1);\n";
+    //cout << Bezier(0.0,1.0) << endl;
+
+    //cout << "Bezier(Coord c0, Coord c1, Coord c2);\n";
+    //cout << Bezier(0,1, 2) << endl;
+
+    //cout << "Bezier(Coord c0, Coord c1, Coord c2, Coord c3);\n";
+    //cout << Bezier(0,1,2,3) << endl;
+
+    //cout << "unsigned degree();\n";
+    EXPECT_EQ(2u, hump.degree());
+
+    //cout << "unsigned size();\n";
+    EXPECT_EQ(3u, hump.size());
+}
+
+TEST_F(BezierTest, ValueAt) {
+    EXPECT_EQ(0.0, wiggle.at0());
+    EXPECT_EQ(3.0, wiggle.at1());
+
+    EXPECT_EQ(0.0, wiggle.valueAt(0.5));
+
+    EXPECT_EQ(0.0, wiggle(0.5));
+
+    //cout << "SBasis toSBasis();\n";
+    //cout << unit.toSBasis() << endl;
+    //cout << hump.toSBasis() << endl;
+    //cout << wiggle.toSBasis() << endl;
+}
+
+TEST_F(BezierTest, Casteljau) {
+    unsigned N = wiggle.order() + 1;
+    std::vector<Coord> left(N), right(N);
+    std::vector<Coord> left2(N), right2(N);
+    double eps = 1e-15;
+
+    for (unsigned i = 0; i < 10000; ++i) {
+        double t = g_random_double_range(0, 1);
+        double vok = bernstein_value_at(t, &wiggle[0], wiggle.order());
+        double v = casteljau_subdivision<double>(t, &wiggle[0], &left[0], &right[0], wiggle.order());
+        EXPECT_near(v, vok, eps);
+        EXPECT_EQ(left[0], wiggle.at0());
+        EXPECT_EQ(left[wiggle.order()], right[0]);
+        EXPECT_EQ(right[wiggle.order()], wiggle.at1());
+
+        double vl = casteljau_subdivision<double>(t, &wiggle[0], &left2[0], NULL, wiggle.order());
+        double vr = casteljau_subdivision<double>(t, &wiggle[0], NULL, &right2[0], wiggle.order());
+        EXPECT_EQ(vl, vok);
+        EXPECT_near(vr, vok, eps);
+        EXPECT_vector_near(left2, left, eps);
+        EXPECT_vector_equal(right2, right);
+
+        double vnone = casteljau_subdivision<double>(t, &wiggle[0], NULL, NULL, wiggle.order());
+        EXPECT_near(vnone, vok, 1e-12);
+    }
+}
+
+TEST_F(BezierTest, Portion) {
+    constexpr Coord eps{1e-12};
+
+    for (unsigned i = 0; i < 10000; ++i) {
+        double from = g_random_double_range(0, 1);
+        double to = g_random_double_range(0, 1);
+        for (auto & input : fragments) {
+            Bezier result = portion(input, from, to);
+
+            // the endpoints must correspond exactly
+            EXPECT_near(result.at0(), input.valueAt(from), eps);
+            EXPECT_near(result.at1(), input.valueAt(to), eps);
+        }
+    }
+}
+
+TEST_F(BezierTest, Subdivide) {
+    std::vector<std::pair<Bezier, double> > errors;
+    for (unsigned i = 0; i < 10000; ++i) {
+        double t = g_random_double_range(0, 1e-6);
+        for (auto & input : fragments) {
+            std::pair<Bezier, Bezier> result = input.subdivide(t);
+
+            // the endpoints must correspond exactly
+            // moreover, the subdivision point must be exactly equal to valueAt(t)
+            EXPECT_DOUBLE_EQ(result.first.at0(), input.at0());
+            EXPECT_DOUBLE_EQ(result.first.at1(), result.second.at0());
+            EXPECT_DOUBLE_EQ(result.second.at0(), input.valueAt(t));
+            EXPECT_DOUBLE_EQ(result.second.at1(), input.at1());
+
+            // ditto for valueAt
+            EXPECT_DOUBLE_EQ(result.first.valueAt(0), input.valueAt(0));
+            EXPECT_DOUBLE_EQ(result.first.valueAt(1), result.second.valueAt(0));
+            EXPECT_DOUBLE_EQ(result.second.valueAt(0), input.valueAt(t));
+            EXPECT_DOUBLE_EQ(result.second.valueAt(1), input.valueAt(1));
+
+            if (result.first.at1() != result.second.at0()) {
+                errors.emplace_back(input, t);
+            }
+        }
+    }
+    if (!errors.empty()) {
+        std::cout << "Found " << errors.size() << " subdivision errors" << std::endl;
+        for (unsigned i = 0; i < errors.size(); ++i) {
+            std::cout << "Error #" << i << ":\n"
+                          << errors[i].first << "\n"
+                          << "t: " << format_coord_nice(errors[i].second) << std::endl;
+        }
+    }
+}
+
+TEST_F(BezierTest, Mutation) {
+//Coord &operator[](unsigned ix);
+//Coord const &operator[](unsigned ix);
+//void setCoeff(unsigned ix double val);
+    //cout << "bigun\n";
+    Bezier bigun(Bezier::Order(30));
+    bigun.setCoeff(5,10.0);
+    for(unsigned i = 0; i < bigun.size(); i++) {
+        EXPECT_EQ((i == 5) ? 10 : 0, bigun[i]);
+    }
+
+    bigun[5] = -3;
+    for(unsigned i = 0; i < bigun.size(); i++) {
+        EXPECT_EQ((i == 5) ? -3 : 0, bigun[i]);
+    }
+}
+
+TEST_F(BezierTest, MultiDerivative) {
+    vector<double> vnd = wiggle.valueAndDerivatives(0.5, 5);
+    expect_array((const double[]){0,0,12,72,0,0}, vnd);
+}
+
+TEST_F(BezierTest, DegreeElevation) {
+    EXPECT_TRUE(are_equal(wiggle, wiggle));
+    Bezier Q = wiggle;
+    Bezier P = Q.elevate_degree();
+    EXPECT_EQ(P.size(), Q.size()+1);
+    //EXPECT_EQ(0, P.forward_difference(1)[0]);
+    EXPECT_TRUE(are_equal(Q, P));
+    Q = wiggle;
+    P = Q.elevate_to_degree(10);
+    EXPECT_EQ(10u, P.order());
+    EXPECT_TRUE(are_equal(Q, P));
+    //EXPECT_EQ(0, P.forward_difference(10)[0]);
+    /*Q = wiggle.elevate_degree();
+    P = Q.reduce_degree();
+    EXPECT_EQ(P.size()+1, Q.size());
+    EXPECT_TRUE(are_equal(Q, P));*/
+}
+//std::pair<Bezier, Bezier > subdivide(Coord t);
+
+// Constructs a linear Bezier with root at t
+Bezier linear_root(double t) {
+    return Bezier(0-t, 1-t);
+}
+
+// Constructs a Bezier with roots at the locations in x
+Bezier array_roots(vector<double> x) {
+    Bezier b(1);
+    for(double i : x) {
+        b = multiply(b, linear_root(i));
+    }
+    return b;
+}
+
+TEST_F(BezierTest, Deflate) {
+    Bezier b = array_roots(vector_from_array((const double[]){0,0.25,0.5}));
+    EXPECT_FLOAT_EQ(0, b.at0());
+    b = b.deflate();
+    EXPECT_FLOAT_EQ(0, b.valueAt(0.25));
+    b = b.subdivide(0.25).second;
+    EXPECT_FLOAT_EQ(0, b.at0());
+    b = b.deflate();
+    const double rootposition = (0.5-0.25) / (1-0.25);
+    constexpr Coord eps{1e-12};
+    EXPECT_near(0.0, b.valueAt(rootposition), eps);
+    b = b.subdivide(rootposition).second;
+    EXPECT_near(0.0, b.at0(), eps);
+}
+
+TEST_F(BezierTest, Roots) {
+    expect_array((const double[]){0, 0.5, 0.5}, wiggle.roots());
+
+    /*Bezier bigun(Bezier::Order(30));
+    for(unsigned i = 0; i < bigun.size(); i++) {
+        bigun.setCoeff(i,rand()-0.5);
+    }
+    cout << bigun.roots() << endl;*/
+
+    // The results of our rootfinding are at the moment fairly inaccurate.
+    double eps = 5e-4;
+
+    vector<vector<double> > tests;
+    tests.push_back(vector_from_array((const double[]){0}));
+    tests.push_back(vector_from_array((const double[]){1}));
+    tests.push_back(vector_from_array((const double[]){0, 0}));
+    tests.push_back(vector_from_array((const double[]){0.5}));
+    tests.push_back(vector_from_array((const double[]){0.5, 0.5}));
+    tests.push_back(vector_from_array((const double[]){0.1, 0.1}));
+    tests.push_back(vector_from_array((const double[]){0.1, 0.1, 0.1}));
+    tests.push_back(vector_from_array((const double[]){0.25,0.75}));
+    tests.push_back(vector_from_array((const double[]){0.5,0.5}));
+    tests.push_back(vector_from_array((const double[]){0, 0.2, 0.6, 0.6, 1}));
+    tests.push_back(vector_from_array((const double[]){.1,.2,.3,.4,.5,.6}));
+    tests.push_back(vector_from_array((const double[]){0.25,0.25,0.25,0.75,0.75,0.75}));
+
+    for(auto & test : tests) {
+        Bezier b = array_roots(test);
+        //std::cout << tests[test_i] << ": " << b << std::endl;
+        //std::cout << b.roots() << std::endl;
+        EXPECT_vector_near(test, b.roots(), eps);
+    }
+}
+
+TEST_F(BezierTest, BoundsExact) {
+    OptInterval unit_bounds = bounds_exact(unit);
+    EXPECT_EQ(unit_bounds->min(), 0);
+    EXPECT_EQ(unit_bounds->max(), 1);
+
+    OptInterval hump_bounds = bounds_exact(hump);
+    EXPECT_EQ(hump_bounds->min(), 0);
+    EXPECT_FLOAT_EQ(hump_bounds->max(), hump.valueAt(0.5));
+
+    OptInterval wiggle_bounds = bounds_exact(wiggle);
+    EXPECT_EQ(wiggle_bounds->min(), 0);
+    EXPECT_EQ(wiggle_bounds->max(), 3);
+}
+
+TEST_F(BezierTest, Operators) {
+    // Test equality operators
+    EXPECT_EQ(zero, zero);
+    EXPECT_EQ(hump, hump);
+    EXPECT_EQ(wiggle, wiggle);
+    EXPECT_EQ(unit, unit);
+
+    EXPECT_NE(zero, hump);
+    EXPECT_NE(hump, zero);
+    EXPECT_NE(wiggle, hump);
+    EXPECT_NE(zero, wiggle);
+    EXPECT_NE(wiggle, unit);
+
+    // Recall that hump == Bezier(0,1,0);
+    EXPECT_EQ(hump + 3, Bezier(3, 4, 3));
+    EXPECT_EQ(hump - 3, Bezier(-3, -2, -3));
+    EXPECT_EQ(hump * 3, Bezier(0, 3, 0));
+    EXPECT_EQ(hump / 3, Bezier(0, 1.0/3.0, 0));
+    EXPECT_EQ(-hump, Bezier(0, -1, 0));
+
+    Bezier reverse_wiggle = reverse(wiggle);
+    EXPECT_EQ(reverse_wiggle.at0(), wiggle.at1());
+    EXPECT_EQ(reverse_wiggle.at1(), wiggle.at0());
+    EXPECT_TRUE(are_equal(reverse(reverse_wiggle), wiggle));
+
+    //cout << "Bezier portion(const Bezier & a, double from, double to);\n";
+    //cout << portion(Bezier(0.0,2.0), 0.5, 1) << endl;
+
+// std::vector<Point> bezier_points(const D2<Bezier > & a) {
+
+    /*cout << "Bezier derivative(const Bezier & a);\n";
+    std::cout << derivative(hump) <<std::endl;
+    std::cout << integral(hump) <<std::endl;*/
+
+    EXPECT_TRUE(are_equal(derivative(integral(wiggle)), wiggle));
+    //std::cout << derivative(integral(hump)) <<std::endl;
+    expect_array((const double []){0.5}, derivative(hump).roots());
+
+    EXPECT_TRUE(bounds_fast(hump)->contains(Interval(0,hump.valueAt(0.5))));
+
+    EXPECT_EQ(Interval(0,hump.valueAt(0.5)), *bounds_exact(hump));
+
+    Interval tight_local_bounds(min(hump.valueAt(0.3),hump.valueAt(0.6)),
+             hump.valueAt(0.5));
+    EXPECT_TRUE(bounds_local(hump, Interval(0.3, 0.6))->contains(tight_local_bounds));
+
+    Bezier Bs[] = {unit, hump, wiggle};
+    for(auto B : Bs) {
+        Bezier product = multiply(B, B);
+        for(int i = 0; i <= 16; i++) {
+            double t = i/16.0;
+            double b = B.valueAt(t);
+            EXPECT_near(b*b, product.valueAt(t), 1e-12);
+        }
+    }
+}
+
+struct XPt {
+    XPt(Coord x, Coord y, Coord ta, Coord tb)
+        : p(x, y), ta(ta), tb(tb)
+    {}
+    XPt() {}
+    Point p;
+    Coord ta, tb;
+};
+
+struct XTest {
+    D2<Bezier> a;
+    D2<Bezier> b;
+    std::vector<XPt> s;
+};
+
+struct CILess {
+    bool operator()(CurveIntersection const &a, CurveIntersection const &b) const {
+        if (a.first < b.first) return true;
+        if (a.first == b.first && a.second < b.second) return true;
+        return false;
+    }
+};
+
+TEST_F(BezierTest, Intersection) {
+    /* Intersection test cases taken from:
+     * Dieter Lasser (1988), Calculating the Self-Intersections of Bezier Curves
+     * https://archive.org/stream/calculatingselfi00lass
+     *
+     * The intersection points are not actually calculated to a high precision
+     * in the paper. The most relevant tests are whether the curves actually
+     * intersect at the returned time values (i.e. whether a(ta) = b(tb))
+     * and whether the number of intersections is correct.
+     */
+    typedef D2<Bezier> D2Bez;
+    std::vector<XTest> tests;
+
+    // Example 1
+    tests.emplace_back();
+    tests.back().a = D2Bez(Bezier(-3.3, -3.3, 0, 3.3, 3.3), Bezier(1.3, -0.7, 2.3, -0.7, 1.3));
+    tests.back().b = D2Bez(Bezier(-4.0, -4.0, 0, 4.0, 4.0), Bezier(-0.35, 3.0, -2.6, 3.0, -0.35));
+    tests.back().s.resize(4);
+    tests.back().s[0] = XPt(-3.12109, 0.76362, 0.09834, 0.20604);
+    tests.back().s[1] = XPt(-1.67341, 0.60298, 0.32366, 0.35662);
+    tests.back().s[2] = XPt(1.67341, 0.60298, 0.67634, 0.64338);
+    tests.back().s[3] = XPt(3.12109, 0.76362, 0.90166, 0.79396);
+
+    // Example 2
+    tests.emplace_back();
+    tests.back().a = D2Bez(Bezier(0, 0, 3, 3), Bezier(0, 14, -9, 5));
+    tests.back().b = D2Bez(Bezier(-1, 13, -10, 4), Bezier(4, 4, 1, 1));
+    tests.back().s.resize(9);
+    tests.back().s[0] = XPt(0.00809, 1.17249, 0.03029, 0.85430);
+    tests.back().s[1] = XPt(0.02596, 1.97778, 0.05471, 0.61825);
+    tests.back().s[2] = XPt(0.17250, 3.99191, 0.14570, 0.03029);
+    tests.back().s[3] = XPt(0.97778, 3.97404, 0.38175, 0.05471);
+    tests.back().s[4] = XPt(1.5, 2.5, 0.5, 0.5);
+    tests.back().s[5] = XPt(2.02221, 1.02596, 0.61825, 0.94529);
+    tests.back().s[6] = XPt(2.82750, 1.00809, 0.85430, 0.96971);
+    tests.back().s[7] = XPt(2.97404, 3.02221, 0.94529, 0.38175);
+    tests.back().s[8] = XPt(2.99191, 3.82750, 0.96971, 0.14570);
+
+    // Example 3
+    tests.emplace_back();
+    tests.back().a = D2Bez(Bezier(-5, -5, -3, 0, 3, 5, 5), Bezier(0, 3.555, -1, 4.17, -1, 3.555, 0));
+    tests.back().b = D2Bez(Bezier(-6, -6, -3, 0, 3, 6, 6), Bezier(3, -0.555, 4, -1.17, 4, -0.555, 3));
+    tests.back().s.resize(6);
+    tests.back().s[0] = XPt(-3.64353, 1.49822, 0.23120, 0.27305);
+    tests.back().s[1] = XPt(-2.92393, 1.50086, 0.29330, 0.32148);
+    tests.back().s[2] = XPt(-0.77325, 1.49989, 0.44827, 0.45409);
+    tests.back().s[3] = XPt(0.77325, 1.49989, 0.55173, 0.54591);
+    tests.back().s[4] = XPt(2.92393, 1.50086, 0.70670, 0.67852);
+    tests.back().s[5] = XPt(3.64353, 1.49822, 0.76880, 0.72695);
+
+    // Example 4
+    tests.emplace_back();
+    tests.back().a = D2Bez(Bezier(-4, -10, -2, -2, 2, 2, 10, 4), Bezier(0, 6, 6, 0, 0, 6, 6, 0));
+    tests.back().b = D2Bez(Bezier(-8, 0, 8), Bezier(1, 6, 1));
+    tests.back().s.resize(4);
+    tests.back().s[0] = XPt(-5.69310, 2.23393, 0.06613, 0.14418);
+    tests.back().s[1] = XPt(-2.68113, 3.21920, 0.35152, 0.33243);
+    tests.back().s[2] = XPt(2.68113, 3.21920, 0.64848, 0.66757);
+    tests.back().s[3] = XPt(5.69310, 2.23393, 0.93387, 0.85582);
+
+    //std::cout << std::setprecision(5);
+
+    for (unsigned i = 0; i < tests.size(); ++i) {
+        BezierCurve a(tests[i].a), b(tests[i].b);
+        std::vector<CurveIntersection> xs;
+        xs = a.intersect(b, 1e-8);
+        std::sort(xs.begin(), xs.end(), CILess());
+        //xs.erase(std::unique(xs.begin(), xs.end(), XEqual()), xs.end());
+
+        std::cout << "\n\n"
+                  << "===============================\n"
+                  << "=== Intersection Testcase " << i+1 << " ===\n"
+                  << "===============================\n" << std::endl;
+
+        EXPECT_EQ(xs.size(), tests[i].s.size());
+        //if (xs.size() != tests[i].s.size()) continue;
+
+        for (unsigned j = 0; j < std::min(xs.size(), tests[i].s.size()); ++j) {
+            std::cout << xs[j].first << " = " << a.pointAt(xs[j].first) << "   "
+                      << xs[j].second << " = " << b.pointAt(xs[j].second) << "\n"
+                      << tests[i].s[j].ta << " = " << tests[i].a.valueAt(tests[i].s[j].ta) << "   "
+                      << tests[i].s[j].tb << " = " << tests[i].b.valueAt(tests[i].s[j].tb) << std::endl;
+        }
+
+        EXPECT_intersections_valid(a, b, xs, 1e-6);
+    }
+
+    #if 0
+    // these contain second-order intersections
+    Coord a5x[] = {-1.5, -1.5, -10, -10, 0, 10, 10, 1.5, 1.5};
+    Coord a5y[] = {0, -8, -8, 9, 9, 9, -8, -8, 0};
+    Coord b5x[] = {-3, -12, 0, 12, 3};
+    Coord b5y[] = {-5, 8, 2.062507, 8, -5};
+    Coord p5x[] = {-3.60359, -5.44653, 0, 5.44653, 3.60359};
+    Coord p5y[] = {-4.10631, -0.76332, 4.14844, -0.76332, -4.10631};
+    Coord p5ta[] = {0.01787, 0.10171, 0.5, 0.89829, 0.98213};
+    Coord p5tb[] = {0.12443, 0.28110, 0.5, 0.71890, 0.87557};
+
+    Coord a6x[] = {5, 14, 10, -12, -12, -2};
+    Coord a6y[] = {1, 6, -6, -6, 2, 2};
+    Coord b6x[] = {0, 2, -10.5, -10.5, 3.5, 3, 8, 6};
+    Coord b6y[] = {0, -8, -8, 9, 9, -4.129807, -4.129807, 3};
+    Coord p6x[] = {6.29966, 5.87601, 0.04246, -4.67397, -3.57214};
+    Coord p6y[] = {1.63288, -0.86192, -2.38219, -2.17973, 1.91463};
+    Coord p6ta[] = {0.03184, 0.33990, 0.49353, 0.62148, 0.96618};
+    Coord p6tb[] = {0.96977, 0.85797, 0.05087, 0.28232, 0.46102};
+    #endif
+}
+
+/** Basic test for intersecting a quadratic Bézier with a line segment. */
+TEST_F(BezierTest, QuadraticIntersectLineSeg)
+{
+    double const EPS = 1e-12;
+    auto const bow = QuadraticBezier({0, 0}, {1, 1}, {2, 0});
+    auto const highhoriz  = LineSegment(Point(0, 0), Point(2, 0));
+    auto const midhoriz   = LineSegment(Point(0, 0.25), Point(2, 0.25));
+    auto const lowhoriz   = LineSegment(Point(0, 0.5), Point(2, 0.5));
+    auto const noninters  = LineSegment(Point(0, 0.5 + EPS), Point(2, 0.5 + EPS));
+    auto const noninters2 = LineSegment(Point(1, 0), Point(1, 0.5 - EPS));
+
+    auto const endpoint_intersections = bow.intersect(highhoriz, EPS);
+    EXPECT_EQ(endpoint_intersections.size(), 2);
+    EXPECT_intersections_valid(bow, highhoriz, endpoint_intersections, EPS);
+    for (auto const &ex : endpoint_intersections) {
+        EXPECT_DOUBLE_EQ(ex.point()[Y], 0.0);
+    }
+
+    auto const mid_intersections = bow.intersect(midhoriz, EPS);
+    EXPECT_EQ(mid_intersections.size(), 2);
+    EXPECT_intersections_valid(bow, midhoriz, mid_intersections, EPS);
+    for (auto const &mx : mid_intersections) {
+        EXPECT_DOUBLE_EQ(mx.point()[Y], 0.25);
+    }
+
+    auto const tangent_intersection = bow.intersect(lowhoriz, EPS);
+    EXPECT_EQ(tangent_intersection.size(), 1);
+    EXPECT_intersections_valid(bow, lowhoriz, tangent_intersection, EPS);
+    for (auto const &tx : tangent_intersection) {
+        EXPECT_DOUBLE_EQ(tx.point()[Y], 0.5);
+    }
+
+    auto no_intersections = bow.intersect(noninters, EPS);
+    EXPECT_TRUE(no_intersections.empty());
+
+    no_intersections = bow.intersect(noninters2, EPS);
+    EXPECT_TRUE(no_intersections.empty());
+}
+
+TEST_F(BezierTest, QuadraticIntersectLineRandom)
+{
+    g_random_set_seed(0xB747A380);
+    auto const diagonal = LineSegment(Point(0, 0), Point(1, 1));
+    double const EPS = 1e-12;
+
+    for (unsigned i = 0; i < 10'000; i++) {
+        auto q = QuadraticBezier({0, 1}, {g_random_double_range(0.0, 1.0), g_random_double_range(0.0, 1.0)}, {1, 0});
+        auto xings = q.intersect(diagonal, EPS);
+        ASSERT_EQ(xings.size(), 1);
+        auto pt = xings[0].point();
+        EXPECT_TRUE(are_near(pt[X], pt[Y], EPS));
+        EXPECT_intersections_valid(q, diagonal, xings, EPS);
+    }
+}
+
+/** Basic test for intersecting a cubic Bézier with a line segment. */
+TEST_F(BezierTest, CubicIntersectLine)
+{
+    double const EPS = 1e-12;
+    auto const wavelet = CubicBezier({0, 0}, {1, 2}, {0, -2}, {1, 0});
+
+    auto const unit_seg = LineSegment(Point(0, 0), Point(1, 0));
+    auto const expect3 = wavelet.intersect(unit_seg, EPS);
+    EXPECT_EQ(expect3.size(), 3);
+    EXPECT_intersections_valid(wavelet, unit_seg, expect3, EPS);
+
+    auto const half_seg = LineSegment(Point(0, 0), Point(0.5, 0));
+    auto const expect2 = wavelet.intersect(half_seg, EPS);
+    EXPECT_EQ(expect2.size(), 2);
+    EXPECT_intersections_valid(wavelet, half_seg, expect2, EPS);
+
+    auto const less_than_half = LineSegment(Point(0, 0), Point(0.5 - EPS, 0));
+    auto const expect1 = wavelet.intersect(less_than_half, EPS);
+    EXPECT_EQ(expect1.size(), 1);
+    EXPECT_intersections_valid(wavelet, less_than_half, expect1, EPS);
+
+    auto const dollar_stroke = LineSegment(Point(0, 0.5), Point(1, -0.5));
+    auto const dollar_xings = wavelet.intersect(dollar_stroke, EPS);
+    EXPECT_EQ(dollar_xings.size(), 3);
+    EXPECT_intersections_valid(wavelet, dollar_stroke, dollar_xings, EPS);
+}
+
+TEST_F(BezierTest, CubicIntersectLineRandom)
+{
+    g_random_set_seed(0xCAFECAFE);
+    auto const diagonal = LineSegment(Point(0, 0), Point(1, 1));
+    double const EPS = 1e-8;
+
+    for (unsigned i = 0; i < 10'000; i++) {
+        double a1 = g_random_double_range(0.0, 1.0);
+        double a2 = g_random_double_range(a1, 1.0);
+        double b1 = g_random_double_range(0.0, 1.0);
+        double b2 = g_random_double_range(0.0, b1);
+
+        auto c = CubicBezier({0, 1}, {a1, a2}, {b1, b2}, {1, 0});
+        auto xings = c.intersect(diagonal, EPS);
+        ASSERT_EQ(xings.size(), 1);
+        auto pt = xings[0].point();
+        EXPECT_TRUE(are_near(pt[X], pt[Y], EPS));
+        EXPECT_intersections_valid(c, diagonal, xings, EPS);
+    }
+}
+
+/** Regression test for issue https://gitlab.com/inkscape/lib2geom/-/issues/47 . */
+TEST_F(BezierTest, Balloon)
+{
+    auto const loop = CubicBezier({0, 0}, {4, -2}, {4, 2}, {0, 0});
+    auto const seghoriz = LineSegment(Point(-1, 0), Point(0, 0));
+
+    for (double EPS : {1e-6, 1e-9, 1e-12}) {
+        // We expect that 2 intersections are found: one at each end of the loop,
+        // both at the coordinates (0, 0).
+        auto xings_horiz = loop.intersect(seghoriz, EPS);
+        EXPECT_EQ(xings_horiz.size(), 2);
+        EXPECT_intersections_valid(loop, seghoriz, xings_horiz, EPS);
+    }
+}
+
+TEST_F(BezierTest, ExpandToTransformedTest)
+{
+    auto test_curve = [] (Curve const &c) {
+        constexpr int N = 50;
+        for (int i = 0; i < N; i++) {
+            auto angle = 2 * M_PI * i / N;
+            auto transform = Affine(Rotate(angle));
+            
+            auto copy = std::unique_ptr<Curve>(c.duplicate());
+            *copy *= transform;
+            auto box1 = copy->boundsExact();
+            
+            auto pt = c.initialPoint() * transform;
+            auto box2 = Rect(pt, pt);
+            c.expandToTransformed(box2, transform);
+            
+            for (auto i : { X, Y }) {
+                EXPECT_DOUBLE_EQ(box1[i].min(), box2[i].min());
+                EXPECT_DOUBLE_EQ(box1[i].max(), box2[i].max());
+            }
+        }
+    };
+    
+    test_curve(LineSegment(Point(-1, 0), Point(1, 2)));
+    test_curve(QuadraticBezier(Point(-1, 0), Point(1, 1), Point(3, 0)));
+    test_curve(CubicBezier(Point(-1, 0), Point(1, 1), Point(2, -2), Point(3, 0)));
+}
+
+TEST_F(BezierTest, ForwardDifferenceTest)
+{
+    auto b = Bezier(3, 4, 2, -5, 7);
+    EXPECT_EQ(b.forward_difference(1), Bezier(19, 34, 22, 5));
+    EXPECT_EQ(b.forward_difference(2), Bezier(-3, 2, 2));
+}
+
+/*
+  Local Variables:
+  mode:c++
+  c-file-style:"stroustrup"
+  c-file-offsets:((innamespace . 0)(inline-open . 0)(case-label . +))
+  indent-tabs-mode:nil
+  fill-column:99
+  End:
+*/
+// vim: filetype=cpp:expandtab:shiftwidth=4:tabstop=8:softtabstop=4:fileencoding=utf-8:textwidth=99 :