sis-commits mailing list archives

Site index · List index
Message view « Date » · « Thread »
Top « Date » · « Thread »
From desruisse...@apache.org
Subject [sis] branch geoapi-4.0 updated: GeodesicCalculator.toGeodesicPath2D() should build a chain of Bézier curves when a single one is not sufficient for the desired precision.
Date Mon, 20 May 2019 16:47:35 GMT
This is an automated email from the ASF dual-hosted git repository.

desruisseaux pushed a commit to branch geoapi-4.0
in repository https://gitbox.apache.org/repos/asf/sis.git


The following commit(s) were added to refs/heads/geoapi-4.0 by this push:
     new 1fc1ebc  GeodesicCalculator.toGeodesicPath2D() should build a chain of Bézier curves when a single one is not sufficient for the desired precision.
1fc1ebc is described below

commit 1fc1ebcd05d7beaf0d6e01e87847ddf5bee1cf23
Author: Martin Desruisseaux <martin.desruisseaux@geomatys.com>
AuthorDate: Mon May 20 18:45:02 2019 +0200

    GeodesicCalculator.toGeodesicPath2D() should build a chain of Bézier curves when a single one is not sufficient for the desired precision.
---
 .../sis/internal/referencing/j2d/Bezier.java       | 317 +++++++++++++++++++++
 .../internal/referencing/j2d/ShapeUtilities.java   | 174 -----------
 .../apache/sis/referencing/GeodeticCalculator.java | 151 +++++++---
 .../referencing/j2d/ShapeUtilitiesExt.java         | 149 ++++++++++
 .../referencing/j2d/ShapeUtilitiesTest.java        |   8 +-
 .../referencing/j2d/ShapeUtilitiesViewer.java      |   2 +-
 .../sis/referencing/GeodeticCalculatorTest.java    |  27 +-
 .../org/apache/sis/test/widget/ShapeViewer.java    |  42 ++-
 .../org/apache/sis/test/widget/VisualCheck.java    |  16 +-
 ide-project/NetBeans/nbproject/genfiles.properties |   2 +-
 ide-project/NetBeans/nbproject/project.xml         |   1 +
 11 files changed, 639 insertions(+), 250 deletions(-)

diff --git a/core/sis-referencing/src/main/java/org/apache/sis/internal/referencing/j2d/Bezier.java b/core/sis-referencing/src/main/java/org/apache/sis/internal/referencing/j2d/Bezier.java
new file mode 100644
index 0000000..90f3999
--- /dev/null
+++ b/core/sis-referencing/src/main/java/org/apache/sis/internal/referencing/j2d/Bezier.java
@@ -0,0 +1,317 @@
+/*
+ * Licensed to the Apache Software Foundation (ASF) under one or more
+ * contributor license agreements.  See the NOTICE file distributed with
+ * this work for additional information regarding copyright ownership.
+ * The ASF licenses this file to You under the Apache License, Version 2.0
+ * (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.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+package org.apache.sis.internal.referencing.j2d;
+
+import java.awt.geom.Path2D;
+import org.opengis.referencing.operation.TransformException;
+
+import static java.lang.Math.abs;
+import static java.lang.Double.isInfinite;
+
+
+/**
+ * Helper class for appending Bézier curves to a path.
+ *
+ * @author  Martin Desruisseaux (Geomatys)
+ * @version 1.0
+ *
+ * @see <a href="https://pomax.github.io/bezierinfo/">A Primer on Bézier Curves</a>
+ *
+ * @since 1.0
+ * @module
+ */
+public abstract class Bezier {
+    /**
+     * The path where to append Bézier curves.
+     */
+    private final Path2D path;
+
+    /**
+     * Maximal distance on <var>x</var> and <var>y</var> axis between the cubic Bézier curve and quadratic
+     * or linear simplifications. Initial value is zero. Subclasses should set the values either in their
+     * constructor or at {@link #evaluateAt(double)} method call.
+     */
+    protected double εx, εy;
+
+    /**
+     * (<var>x</var> and <var>y</var>) coordinates and (∂y/∂x) derivative of starting point at <var>t</var>=0.
+     */
+    private double x0, y0, α0;
+
+    /**
+     * A buffer used by subclasses for storing results of {@link #evaluateAt(double)}.
+     * The two first elements are <var>x</var> and <var>y</var> coordinates respectively.
+     * Other elements (if any) are ignored.
+     */
+    protected final double[] point;
+
+    /**
+     * Creates a new builder.
+     *
+     * @param  dimension  length of the {@link #point} array. Must be at least 2.
+     */
+    protected Bezier(final int dimension) {
+        point = new double[dimension];
+        path  = new Path2D.Double();
+    }
+
+    /**
+     * Invoked for computing a new point on the Bézier curve. This method is invoked with a <var>t</var> value varying from
+     * 0 to 1 inclusive. Value 0 is for the starting point and value 1 is for the ending point. Other values are for points
+     * interpolated between the start and end points. In particular value ½ is for the point in the middle of the curve.
+     * This method will also be invoked at least for values ¼ and ¾, and potentially for other values too.
+     *
+     * <p>This method shall store the point coordinates in the {@link #point} array with <var>x</var> coordinate in the
+     * first element and <var>y</var> coordinate in the second element. The returned value is the derivative (∂y/∂x) at
+     * that location. If this method can not compute a coordinate, it can store {@link Double#NaN} values except for
+     * <var>t</var>=0, ½ and 1 where finite coordinates are mandatory.</p>
+     *
+     * <p>Subclasses can optionally update the {@link #εx} and {@link #εy} values if the tolerance thresholds change as a
+     * result of this method call, for example because we come closer to a pole. The tolerance values used for each Bézier
+     * curve are the ones computed at <var>t</var>=¼ and <var>t</var>=¾ of that curve.</p>
+     *
+     * @param  t  desired point on the curve, from 0 (start point) to 1 (end point) inclusive.
+     * @return derivative (∂y/∂x) at the point.
+     * @throws TransformException if the point coordinates can not be computed.
+     */
+    protected abstract double evaluateAt(double t) throws TransformException;
+
+    /**
+     * Creates a sequence of Bézier curves from the position given by {@code evaluateAt(0)} to the position given
+     * by {@code evaluateAt(0)}. This method determines the number of intermediate points required for achieving
+     * the precision requested by the {@code εx} and {@code εy} parameters given at construction time.
+     *
+     * @return the sequence of Bézier curves.
+     * @throws TransformException if the coordinates of a point can not be computed.
+     */
+    public final Path2D build()  throws TransformException {
+        α0 = evaluateAt(0);
+        x0 = point[0];
+        y0 = point[1];
+        path.moveTo(x0, y0);
+        final double α2 = evaluateAt(0.5);
+        final double x2 = point[0];
+        final double y2 = point[1];
+        final double α4 = evaluateAt(1);                                // `point` become (x₄,y₄) with this method call.
+        sequence(0, 1, x2, y2, α2, point[0], point[1], α4);             // Must be after the call to `evaluateAt(t₄)`.
+        return path;
+    }
+
+    /**
+     * Creates a sequence of Bézier curves from the position given by {@code evaluateAt(tmin)} to the position given by
+     * {@code evaluateAt(tmax)}. This method invokes itself recursively as long as more points are needed for achieving
+     * the precision requested by the {@code εx} and {@code εy} parameters given at construction time.
+     *
+     * <p>The {@link #x0}, {@link #y0} and {@link #α0} fields must be set to the start point before this method is invoked.
+     * On method completion those fields will be updated to the coordinates and derivative of end point, thus enabling
+     * chaining with the next {@code sequence(…)} call.</p>
+     *
+     * @param  tmin  <var>t</var> value of the start point (initially 0).
+     * @param  tmax  <var>t</var> value of the end point   (initially 1).
+     * @param  x2    <var>x</var> coordinate at <var>t</var>=½ (mid-point).
+     * @param  y2    <var>y</var> coordinate at <var>t</var>=½ (mid-point).
+     * @param  α2    the derivative (∂y/∂x) at mid-point.
+     * @param  x4    <var>x</var> coordinate at <var>t</var>=1 (end point).
+     * @param  y4    <var>y</var> coordinate at <var>t</var>=1 (end point).
+     * @param  α4    the derivative (∂y/∂x) at end point.
+     * @throws TransformException if the coordinates of a point can not be computed.
+     */
+    private void sequence(final double tmin, final double tmax,
+            final double x2, final double y2, final double α2,
+            final double x4, final double y4, final double α4) throws TransformException
+    {
+        final double α3 = evaluateAt(0.25*tmin + 0.75*tmax);
+        final double x3 = point[0];
+        final double y3 = point[1];
+        final double tx = εx;
+        final double ty = εy;
+        final double α1 = evaluateAt(0.75*tmin + 0.25*tmax);            // `point` become (x₁,y₁) with this method call.
+        if (tx < εx) εx = tx;                                           // Take smallest tolerance values.
+        if (ty < εy) εy = ty;
+        if (curve(point[0], point[1], x2, y2, x3, y3, x4, y4, α4)) {    // Must be after the call to `evaluateAt(t₁)`.
+            x0 = x4;
+            y0 = y4;
+            α0 = α4;
+        } else {
+            final double εxo = εx;
+            final double εyo = εy;
+            final double th = 0.5 * (tmin + tmax);
+            sequence(tmin, th, point[0], point[1], α1, x2, y2, α2);
+            sequence(th, tmax, x3,       y3,       α3, x4, y4, α4);
+            εx = εxo;
+            εy = εyo;
+        }
+    }
+
+    /**
+     * Computes a Bézier curve passing by the given points and with the given derivatives at end points.
+     * The cubic curve equation is:
+     *
+     * <blockquote>B(t) = (1-t)³⋅P₀ + 3(1-t)²t⋅A + 3(1-t)t²⋅B + t³⋅P₄</blockquote>
+     *
+     * where t ∈ [0…1], P₀ and P₄ are start and end points of the curve and A and B are control points generally not on the curve.
+     * For any point <var>P<sub>i</sub></var>, the index <var>i</var> gives the value of the parameter <var>t</var> where the point
+     * is located with <var>t</var> = <var>i</var>/4, and coordinates (<var>x</var>, <var>y</var>) follow the same indices.
+     * The (x₀,y₀) fields give the coordinates of point P₀ at <var>t</var>=0 (0/4).
+     * The (x₂,y₂) arguments give the coordinates of the point at <var>t</var>=½ (2/4).
+     * The (x₄,y₄) arguments give the coordinates of point P₄ at <var>t</var>=1 (4/4).
+     *
+     * <p>The P₁ and P₃ points are for quality control. They shall be points at <var>t</var>=¼ <var>t</var>=¾ respectively.
+     * Those points are not used for computing the curve, but are used for checking if the curve is an accurate approximation.
+     * If the curve is not accurate enough, then this method does nothing and return {@code false}. In that case, the caller
+     * should split the curve in two smaller parts and invoke this method again.</p>
+     *
+     * <p>If the full equation is required for representing the curve, then this method appends a {@link java.awt.geom.CubicCurve2D}.
+     * If the same curve can be represented by a quadratic curve, then this method appends a {@link java.awt.geom.QuadCurve2D}.
+     * If the curve is actually a straight line, then this method appends a {@link java.awt.geom.Line2D}.</p>
+     *
+     * @param  x1  <var>x</var> coordinate at <var>t</var>=¼ (quality control, may be NaN).
+     * @param  y1  <var>y</var> coordinate at <var>t</var>=¼ (quality control, may be NaN).
+     * @param  x2  <var>x</var> coordinate at <var>t</var>=½ (mid-point).
+     * @param  y2  <var>y</var> coordinate at <var>t</var>=½ (mid-point).
+     * @param  x3  <var>x</var> coordinate at <var>t</var>=¾ (quality control, may be NaN).
+     * @param  y3  <var>y</var> coordinate at <var>t</var>=¾ (quality control, may be NaN).
+     * @param  x4  <var>x</var> coordinate at <var>t</var>=1 (end point).
+     * @param  y4  <var>y</var> coordinate at <var>t</var>=1 (end point).
+     * @param  α4  the derivative (∂y/∂x) at ending point.
+     * @return {@code true} if the curve has been added to the path, or {@code false} if the approximation is not sufficient.
+     */
+    final boolean curve(double x1, double y1,
+                        double x2, double y2,
+                        double x3, double y3,
+                        final double x4, final double y4, final double α4)
+    {
+        /*
+         * Equations in this method are simplified as if (x₀,y₀) coordinates are (0,0).
+         * Adjust (x₁,y₁) to (x₄,y₄) consequently. If derivatives are equal, equation
+         * for cubic curve will not work (division by zero). But we can return a line
+         * instead if derivatives are equal to Δy/Δx slope and way-points are colinear.
+         */
+        x1 -= x0;  x2 -= x0;  x3 -= x0;
+        y1 -= y0;  y2 -= y0;  y3 -= y0;
+        final double Δx = x4 - x0;
+        final double Δy = y4 - y0;
+        final boolean isVertical = abs(Δx) <= εx;
+        if (isVertical || ((isInfinite(α0) || abs(Δx*α0 - Δy) <= εy)
+                       &&  (isInfinite(α4) || abs(Δx*α4 - Δy) <= εy)))
+        {
+            /*
+             * Following tests are partially redundant with above tests, but may detect a larger
+             * error than above tests did. They are also necessary if a derivative was infinite,
+             * or if `isVertical` was true.
+             */
+            final boolean isHorizontal = abs(Δy) <= εy;
+            if (isHorizontal || ((α0 == 0 || abs(Δy/α0 - Δx) <= εx)
+                             &&  (α4 == 0 || abs(Δy/α4 - Δx) <= εx)))
+            {
+                /*
+                 * Verify that all points are on the line joining P₀ to P₄. We test P₂ first since it is the most likely
+                 * to be away if the curve is not almost a line. If any coordinate tested below is NaN, ignore it since
+                 * they are not used in Line2D computation. This allows `evaluateAt(t)` method to return NaN if it can
+                 * not provide a point for a given `t` value.
+                 */
+                final double slope = Δy / Δx;
+                if ((!isVertical   && (abs(x2*slope - y2) > εy || abs(x1*slope - y1) > εy || abs(x3*slope - y3) > εy)) ||
+                    (!isHorizontal && (abs(y2/slope - x2) > εx || abs(y1/slope - x1) > εx || abs(y3/slope - x3) > εx)))
+                {
+                    return false;
+                }
+                path.lineTo(x4, y4);
+                return true;
+            }
+        }
+        /*
+         * Bezier curve equation for starting point P₀, ending point P₄ and control points A and B
+         *
+         * t ∈ [0…1]:  B(t)  =  (1-t)³⋅P₀ + 3(1-t)²t⋅A + 3(1-t)t²⋅B + t³⋅P₄
+         * Midpoint:   B(½)  =  ⅛⋅(P₀ + P₄) + ⅜⋅(A + B)
+         *
+         * Notation:   (x₀, y₀)   are coordinates of P₀ (same rule for P₄).
+         *             (x₂, y₂)   are coordinates of midpoint.
+         *             (Ax, Ay)   are coordinates of control point A (same rule for B).
+         *             α₀ and α₄  are derivative (∂y/∂x) at P₀ and P₄ respectively.
+         *
+         * Some relationships:
+         *
+         *     x₂ = ⅛⋅(x₀ + x₄) + ⅜⋅(Ax + Bx)
+         *     (Ay - y₀) / (Ax - x₀) = α₀
+         *     (y₄ - By) / (x₄ - Bx) = α₄
+         *
+         * Setting (x₀,y₀) = (0,0) for simplicity and rearranging above equations:
+         *
+         *     Ax = (8⋅x₂ - x₄)/3 - Bx              where    Bx = x₄ - (y₄ - By)/α₄
+         *        = (8⋅x₂ - 4⋅x₄)/3 + (y₄ - By)/α₄
+         *
+         * Doing similar rearrangement for y:
+         *
+         *     By = (8⋅y₂ - y₄)/3 - Ay    where    Ay = Ax⋅α₀
+         *        = (8⋅y₂ - y₄)/3 - Ax⋅α₀
+         *
+         * Putting together and isolating Ax:
+         *
+         *      Ax = (8⋅x₂ - 4⋅x₄)/3 + (Ax⋅α₀ - (8⋅y₂ - 4⋅y₄)/3)/α₄
+         *         = (8⋅x₂ - 4⋅x₄ - (8⋅y₂ - 4⋅y₄)/α₄) / 3(1 - α₀/α₄)
+         */
+        final double ax = ((8*x2 - 4*Δx)*α4 - (8*y2 - 4*Δy)) / (3*(α4 - α0));
+        final double ay = ax * α0;
+        final double bx = (8*x2 - Δx)/3 - ax;
+        final double by = Δy - (Δx - bx)*α4;
+        /*
+         * At this point we got the control points A and B. Verify if the curve is a good approximation.
+         * We compute the points on the curve at t=¼ and t=¾ as below (simplified with P₀ = (0,0)) and
+         * compare with the given values:
+         *
+         *     P(¼) = (27⋅P₀ + 27⋅A +  9⋅B +    P₄)/64
+         *     P(¾) = (   P₀ +  9⋅A + 27⋅B + 27⋅P₄)/64
+         *
+         * If any of (x₁,y₁) and (x₃,y₃) coordinates is NaN, then this method accepts the curve as valid.
+         * This allows `evaluateAt(t)` method to return NaN if it can not provide a point for a given `t` value.
+         */
+        if (abs(27./64*ax + 9./64*bx + 1./64*Δx - x1) > εx || abs(9./64*ax + 27./64*bx + 27./64*Δx - x3) > εx ||
+            abs(27./64*ay + 9./64*by + 1./64*Δy - y1) > εy || abs(9./64*ay + 27./64*by + 27./64*Δy - y3) > εy)
+        {
+            return false;
+        }
+        /*
+         * Verify if we can simplify cubic curve to a quadratic curve. If we were elevating the Bézier curve degree from
+         * quadratic to cubic, the control points A and B would be computed from control point C of the quadratic curve:
+         *
+         *     A  =  ⅓P₀ + ⅔C
+         *     B  =  ⅓P₄ + ⅔C
+         *
+         * We want C instead, which can be computed in two ways:
+         *
+         *     C   =  (3A - P₀)/2
+         *     C   =  (3B - P₄)/2
+         *
+         * We compute C both ways and check if they are close enough to each other:
+         *
+         *     ΔC  =  (3⋅(B - A) - (P₄ - P₀))/2
+         */
+        final double Δqx, Δqy;
+        if (abs(Δqx = (3*(bx - ax) - Δx)/2) <= εx &&        // P₀ is zero.
+            abs(Δqy = (3*(by - ay) - Δy)/2) <= εy)
+        {
+            final double qx = (3*ax + Δqx)/2;               // Take average of 2 control points.
+            final double qy = (3*ay + Δqy)/2;
+            path.quadTo(qx+x0, qy+y0, x4, y4);
+        } else {
+            path.curveTo(ax+x0, ay+y0, bx+x0, by+y0, x4, y4);
+        }
+        return true;
+    }
+}
diff --git a/core/sis-referencing/src/main/java/org/apache/sis/internal/referencing/j2d/ShapeUtilities.java b/core/sis-referencing/src/main/java/org/apache/sis/internal/referencing/j2d/ShapeUtilities.java
index 2226f4e..7915e2a 100644
--- a/core/sis-referencing/src/main/java/org/apache/sis/internal/referencing/j2d/ShapeUtilities.java
+++ b/core/sis-referencing/src/main/java/org/apache/sis/internal/referencing/j2d/ShapeUtilities.java
@@ -359,180 +359,6 @@ public final class ShapeUtilities extends Static {
     }
 
     /**
-     * Returns a Bézier curve passing by the given points and with the given derivatives at end points.
-     * The cubic curve equation is:
-     *
-     * <blockquote>B(t) = (1-t)³⋅P₁ + 3(1-t)²t⋅C₁ + 3(1-t)t²⋅C₂ + t³⋅P₂</blockquote>
-     *
-     * where t ∈ [0…1], P₁ and P₂ are end points of the curve and C₁ and C₂ are control points generally not on the curve.
-     * If the full equation is required for representing the curve, then this method builds a {@link CubicCurve2D}.
-     * If the same curve can be represented by a quadratic curve, then this method returns a {@link QuadCurve2D}.
-     * If the curve is actually a straight line, then this method returns a {@link Line2D}.
-     *
-     * <p>The (x₁,y₁) arguments give the coordinates of point P₁ at <var>t</var>=0.
-     * The (x<sub>m</sub>,y<sub>m</sub>) arguments give the coordinates of the point at <var>t</var>=½.
-     * The (x₂,y₂) arguments give the coordinates of point P₂ at <var>t</var>=1.</p>
-     *
-     * @param  x1  <var>x</var> value of the starting point.
-     * @param  y1  <var>y</var> value of the starting point.
-     * @param  xm  <var>x</var> value of the mid-point.
-     * @param  ym  <var>y</var> value of the mid-point.
-     * @param  x2  <var>x</var> value of the ending point.
-     * @param  y2  <var>y</var> value of the ending point.
-     * @param  α1  the derivative (∂y/∂x) at starting point.
-     * @param  α2  the derivative (∂y/∂x) at ending point.
-     * @param  εx  maximal distance on <var>x</var> axis between the cubic Bézier curve and quadratic or linear simplifications.
-     * @param  εy  maximal distance on <var>y</var> axis between the cubic Bézier curve and quadratic or linear simplifications.
-     * @return the Bézier curve passing by the 3 given points and having the given derivatives at end points.
-     *
-     * @see <a href="https://pomax.github.io/bezierinfo/">A Primer on Bézier Curves</a>
-     *
-     * @since 1.0
-     */
-    public static Shape bezier(final double x1, final double y1,
-                                     double xm,       double ym,
-                               final double x2, final double y2,
-                               final double α1, final double α2,
-                               final double εx, final double εy)
-    {
-        /*
-         * Equations in this method are simplified as if (x1,y1) coordinates are (0,0).
-         * Adjust (xm,ym) and (x2,y2) consequently. If derivatives are equal, equation
-         * for cubic curve will not work (division by zero). But we can return a line
-         * instead if derivatives are equal to Δy/Δx slope and mid-point is colinear.
-         */
-        xm -= x1;
-        ym -= y1;
-        final double Δx = x2 - x1;
-        final double Δy = y2 - y1;
-        final boolean isVertical = abs(Δx) <= εx;
-        if (isVertical || ((isInfinite(α1) || abs(Δx*α1 - Δy) <= εy)
-                       &&  (isInfinite(α2) || abs(Δx*α2 - Δy) <= εy)))
-        {
-            /*
-             * Following tests are partially redundant with above tests, but may detect a larger
-             * error than above tests did. They are also necessary if a derivative was infinite,
-             * or if `isVertical` was true.
-             */
-            final boolean isHorizontal = abs(Δy) <= εy;
-            if (isHorizontal || ((α1 == 0 || abs(Δy/α1 - Δx) <= εx)
-                             &&  (α2 == 0 || abs(Δy/α2 - Δx) <= εx)))
-            {
-                final double slope = Δy / Δx;
-                if ((isVertical   || abs(xm*slope - ym) <= εy) &&
-                    (isHorizontal || abs(ym/slope - xm) <= εx))
-                {
-                    return new Line2D.Double(x1, y1, x2, y2);
-                }
-            }
-        }
-        /*
-         * Bezier curve equation for starting point P₀, ending point P₃ and control points P₁ and P₂
-         * (note: not the same numbers than the ones in arguments and variables used in this method):
-         *
-         * t ∈ [0…1]:  B(t)  =  (1-t)³⋅P₀ + 3(1-t)²t⋅P₁ + 3(1-t)t²⋅P₂ + t³⋅P₃
-         * Midpoint:   B(½)  =  ⅛⋅(P₀ + P₃) + ⅜⋅(P₁ + P₂)
-         *
-         * Notation:   (x₀, y₀)   are coordinates of P₀ (same rule for P₁, P₂, P₃).
-         *             (xm, ym)   are coordinates of midpoint.
-         *             α₀ and α₃  are derivative (∂y/∂x) at P₀ and P₃ respectively.
-         *
-         * Some relationships:
-         *
-         *     xm = ⅛⋅(x₀ + x₃) + ⅜⋅(x₁ + x₂)
-         *     (y₁ - y₀) / (x₁ - x₀) = α₀
-         *     (y₃ - y₂) / (x₃ - x₂) = α₃
-         *
-         * Setting (x₀,y₀) = (0,0) for simplicity and rearranging above equations:
-         *
-         *     x₁ = (8⋅xm - x₃)/3 - x₂              where    x₂ = x₃ - (y₃ - y₂)/α₃
-         *     x₁ = (8⋅xm - 4⋅x₃)/3 + (y₃ - y₂)/α₃
-         *
-         * Doing similar rearrangement for y:
-         *
-         *     y₂ = (8⋅ym - y₃)/3 - y₁    where    y₁ = x₁⋅α₀
-         *     y₂ = (8⋅ym - y₃)/3 - x₁⋅α₀
-         *
-         * Putting together and isolating x₁:
-         *
-         *      x₁ = (8⋅xm - 4⋅x₃)/3 + (x₁⋅α₀ - (8⋅ym - 4⋅y₃)/3)/α₃
-         *      x₁ = (8⋅xm - 4⋅x₃ - (8⋅ym - 4⋅y₃)/α₃) / 3(1 - α₀/α₃)
-         *
-         * x₁ and x₂ are named cx1 and cx2 in the code below ("c" for "control").
-         * x₀ and x₃ are named x1 and x2 for consistency with Java2D usage.
-         * Same changes apply to y.
-         */
-        final double cx1 = ((8*xm - 4*Δx)*α2 - (8*ym - 4*Δy)) / (3*(α2 - α1));
-        final double cy1 = cx1 * α1;
-        final double cx2 = (8*xm - Δx)/3 - cx1;
-        final double cy2 = Δy - (Δx - cx2)*α2;
-        /*
-         * At this point we got the control points (cx1,cy1) and (cx2,cy2). Verify if we can simplify
-         * cubic curbe to a quadratic curve. If we were elevating the degree from quadratic to cubic,
-         * the control points C₁ and C₂ would be: (Q is the control point of the quadratic curve)
-         *
-         *     C₁  =  ⅓P₁ + ⅔Q
-         *     C₂  =  ⅓P₂ + ⅔Q
-         *
-         * We want Q instead, which can be computed in two ways:
-         *
-         *     Q   =  (3C₁ - P₁)/2
-         *     Q   =  (3C₂ - P₂)/2
-         *
-         * We compute Q both ways and check if they are close enough to each other:
-         *
-         *     ΔQ  =  (3⋅(C₂ - C₁) - (P₂ - P₁))/2
-         */
-        final double Δqx = (3*(cx2 - cx1) - Δx)/2;          // P₁ set to zero.
-        if (abs(Δqx) <= εx) {
-            final double Δqy = (3*(cy2 - cy1) - Δy)/2;
-            if (abs(Δqy) <= εy) {
-                final double qx = (3*cx1 + Δqx)/2;          // Take average of 2 control points.
-                final double qy = (3*cy1 + Δqy)/2;
-                return new QuadCurve2D.Double(x1, y1, qx+x1, qy+y1, x2, y2);
-            }
-        }
-        return new CubicCurve2D.Double(x1, y1, cx1+x1, cy1+y1, cx2+x1, cy2+y1, x2, y2);
-    }
-
-    /**
-     * Returns a point on the given linear, quadratic or cubic Bézier curve.
-     *
-     * @param  bezier  a {@link Line2D}, {@link QuadCurve2D} or {@link CubicCurve2D}.
-     * @param  t       a parameter from 0 to 1 inclusive.
-     * @return a point on the curve for the given <var>t</var> parameter.
-     *
-     * @see <a href="https://en.wikipedia.org/wiki/B%C3%A9zier_curve">Bézier curve on Wikipedia</a>
-     */
-    public static Point2D.Double pointOnBezier(final Shape bezier, final double t) {
-        final double x, y;
-        final double mt = 1 - t;
-        if (bezier instanceof Line2D) {
-            final Line2D z = (Line2D) bezier;
-            x = mt * z.getX1()  +  t * z.getX2();
-            y = mt * z.getY1()  +  t * z.getY2();
-        } else if (bezier instanceof QuadCurve2D) {
-            final QuadCurve2D z = (QuadCurve2D) bezier;
-            final double a = mt * mt;
-            final double b = mt * t * 2;
-            final double c =  t * t;
-            x = a * z.getX1()  +  b * z.getCtrlX()  +  c * z.getX2();
-            y = a * z.getY1()  +  b * z.getCtrlY()  +  c * z.getY2();
-        } else if (bezier instanceof CubicCurve2D) {
-            final CubicCurve2D z = (CubicCurve2D) bezier;
-            final double a = mt * mt * mt;
-            final double b = mt * mt * t  * 3;
-            final double c = mt * (t * t) * 3;
-            final double d =  t *  t * t;
-            x = a * z.getX1()  +  b * z.getCtrlX1()  +  c * z.getCtrlX2()  +  d * z.getX2();
-            y = a * z.getY1()  +  b * z.getCtrlY1()  +  c * z.getCtrlY2()  +  d * z.getY2();
-        } else {
-            throw new IllegalArgumentException();
-        }
-        return new Point2D.Double(x, y);
-    }
-
-    /**
      * Returns the center of a circle passing by the 3 given points. The distance between the returned
      * point and any of the given points will be constant; it is the circle radius.
      *
diff --git a/core/sis-referencing/src/main/java/org/apache/sis/referencing/GeodeticCalculator.java b/core/sis-referencing/src/main/java/org/apache/sis/referencing/GeodeticCalculator.java
index 2138605..1b88c27 100644
--- a/core/sis-referencing/src/main/java/org/apache/sis/referencing/GeodeticCalculator.java
+++ b/core/sis-referencing/src/main/java/org/apache/sis/referencing/GeodeticCalculator.java
@@ -40,6 +40,7 @@ import org.apache.sis.geometry.CoordinateFormat;
 import org.apache.sis.internal.referencing.PositionTransformer;
 import org.apache.sis.internal.referencing.ReferencingUtilities;
 import org.apache.sis.internal.referencing.j2d.ShapeUtilities;
+import org.apache.sis.internal.referencing.j2d.Bezier;
 import org.apache.sis.internal.referencing.Resources;
 import org.apache.sis.internal.referencing.Formulas;
 import org.apache.sis.internal.util.Numerics;
@@ -492,8 +493,8 @@ public class GeodeticCalculator {
         }
         final double Δλ = λ2 - λ1;
         final double Δφ = φ2 - φ1;
-        double factor;
-        if (abs(Δφ) < 1E-7) {
+        final double factor;
+        if (abs(Δφ) < Formulas.ANGULAR_TOLERANCE) {
             factor = Δλ * cos((φ1 + φ2)/2);
         } else {
             /*
@@ -628,19 +629,16 @@ public class GeodeticCalculator {
      * </ol>
      *
      * <b>Limitations:</b>
-     * current implementation builds a single linear, quadratic or cubic Bézier curve. It does not yet create a chain
-     * of Bézier curves. Consequently the errors may be larger than the given {@code tolerance} threshold. Another
-     * limitation is that this method depends on the presence of {@code java.desktop} module. Those limitations may be
-     * addressed in a future version (see <a href="https://issues.apache.org/jira/browse/SIS-453">SIS-453</a>).
+     * This method depends on the presence of {@code java.desktop} module. This limitations may be addressed
+     * in a future Apache SIS version (see <a href="https://issues.apache.org/jira/browse/SIS-453">SIS-453</a>).
      *
      * @param  tolerance  maximal error between the approximated curve and actual geodesic track
      *                    in the units of measurement given by {@link #getDistanceUnit()}.
-     *                    <em>See limitations above</em>.
      * @return an approximation of geodesic track as Bézier curves in a Java2D object.
      * @throws TransformException if the coordinates can not be transformed to {@linkplain #getPositionCRS() position CRS}.
      * @throws IllegalStateException if some required properties have not been specified.
      */
-    public Shape toGeodesicPath2D(double tolerance) throws TransformException {
+    public Shape toGeodesicPath2D(final double tolerance) throws TransformException {
         if (isInvalid(START_POINT | STARTING_AZIMUTH | END_POINT | ENDING_AZIMUTH | GEODESIC_DISTANCE)) {
             if (isInvalid(END_POINT)) {
                 computeEndPoint();
@@ -648,50 +646,113 @@ public class GeodeticCalculator {
                 computeDistance();
             }
         }
-        tolerance *= (180/PI) / radius;                                     // Angular tolerance in degrees.
-        double d1, x1, y1, d2, x2, y2;                                      // Parameters for the Bezier curve.
-        x2 = λ2;
-        final double sign = signum(α1);
-        if (sign == signum(λ1 - x2)) {
-            x2 += 2*PI * sign;                  // We need λ₁ < λ₂ if heading east, or λ₁ > λ₂ if heading west.
-        }
-        final double[] transformed = new double[ReferencingUtilities.getDimension(userToGeodetic.defaultCRS)];
-        d1 = slope(α1, geographic(φ1, λ1).inverseTransform(transformed)); x1 = transformed[0]; y1 = transformed[1];
-        d2 = slope(α2, geographic(φ2, x2).inverseTransform(transformed)); x2 = transformed[0]; y2 = transformed[1];
-        final double sφ2 = φ2;                                              // Save setting before modification.
-        final double sλ2 = λ2;
-        final double sα2 = α2;
-        final double sd  = geodesicDistance;
+        final PathBuilder bezier = new PathBuilder(toDegrees(tolerance / radius));
+        final Shape path;
         try {
-            geodesicDistance /= 2;
-            computeEndPoint();
-            final Matrix d = geographic(φ2, λ2).inverseTransform(transformed);      // Coordinates of midway point.
-            double εx;                                                              // Tolerance for φ (first coordinate).
-            double εy = tolerance / cos(φ2);                                        // Tolerance for λ (second coordinate).
-            εx = d.getElement(0,0)*tolerance + d.getElement(0,1)*εy;                // Tolerance for x in user CRS.
-            εy = d.getElement(1,0)*tolerance + d.getElement(1,1)*εy;                // Tolerance for y in user CRS.
-            return ShapeUtilities.bezier(x1, y1, transformed[0], transformed[1], x2, y2, d1, d2, εx, εy);
+            path = bezier.build();
         } finally {
-            φ2 = sφ2;                                                               // Restore the setting previously saved.
-            λ2 = sλ2;
-            α2 = sα2;
-            geodesicDistance = sd;
+            bezier.reset();
         }
+        return ShapeUtilities.toPrimitive(path);
     }
 
     /**
-     * Returns the tangent of the given angle converted to the user CRS space.
-     *
-     * @param  α  azimuth angle in radians, with 0 pointing toward north and values increasing clockwise.
-     * @param  d  Jacobian matrix from (φ,λ) to the user coordinate reference system.
-     * @return ∂y/∂x.
+     * Builds a geodesic path as a sequence of Bézier curves. The start point and end points are the points
+     * in enclosing {@link GeodeticCalculator} at the time this class is instantiated. The start coordinates
+     * given by {@link #φ1} and {@link #λ1} shall never change for this whole builder lifetime. However the
+     * end coordinates ({@link #φ2}, {@link #λ2}) will vary at each step.
      */
-    private double slope(final double α, final Matrix d) {
-        final double dx = cos(α);     // sin(π/2 - α) = -sin(α - π/2) = cos(α)
-        final double dy = sin(α);     // cos(π/2 - α) = +cos(α - π/2) = sin(α)
-        final double tx = d.getElement(0,0)*dx + d.getElement(0,1)*dy;
-        final double ty = d.getElement(1,0)*dx + d.getElement(1,1)*dy;
-        return ty / tx;
+    private final class PathBuilder extends Bezier {
+        /**
+         * End point coordinates and derivative, together with geodesic and loxodromic distances.
+         * Saved for later restoration by {@link #reset()}.
+         */
+        private final double φf, λf, αf, distance, length;
+
+        /**
+         * {@link #validity} flags at the time {@code PathBuilder} is instantiated.
+         * Saved for later restoration by {@link #reset()}.
+         */
+        private final int flags;
+
+        /**
+         * Angular tolerance at equator in degrees.
+         */
+        private final double tolerance;
+
+        /**
+         * Creates a builder for the given tolerance at equator in degrees.
+         */
+        PathBuilder(final double tolerance) {
+            super(ReferencingUtilities.getDimension(userToGeodetic.defaultCRS));
+            this.tolerance = tolerance;
+            φf       = φ2;
+            λf       = λ2;
+            αf       = α2;
+            distance = geodesicDistance;
+            length   = rhumblineLength;
+            flags    = validity;
+        }
+
+        /**
+         * Invoked for computing a new point on the Bézier curve. This method is invoked with a <var>t</var> value varying from
+         * 0 (start point) to 1 (end point) inclusive. This method stores the point coordinates in the {@link #point} array with
+         * <var>x</var> coordinate in the first element and <var>y</var> coordinate in the second element. The returned value is
+         * the derivative (∂y/∂x) at that location.
+         *
+         * @param  t  desired point on the curve, from 0 (start point) to 1 (end point) inclusive.
+         * @return derivative (∂y/∂x) at the point.
+         * @throws TransformException if the point coordinates can not be computed.
+         */
+        @Override
+        protected double evaluateAt(final double t) throws TransformException {
+            if (t == 0) {
+                φ2 = φ1;                        // Start point requested.
+                λ2 = λ1;
+                α2 = α1;
+            } else if (t == 1) {
+                φ2 = φf;                        // End point requested.
+                λ2 = λf;
+                α2 = αf;
+            } else {
+                geodesicDistance = distance * t;
+                computeEndPoint();
+            }
+            final double sign = signum(α1);
+            if (sign == signum(λ1 - λ2)) {
+                λ2 += 2*PI * sign;              // We need λ₁ < λ₂ if heading east, or λ₁ > λ₂ if heading west.
+            }
+            final Matrix d = geographic(φ2, λ2).inverseTransform(point);    // Coordinates and Jacobian of point.
+            final double m00 = d.getElement(0,0);
+            final double m01 = d.getElement(0,1);
+            final double m10 = d.getElement(1,0);
+            final double m11 = d.getElement(1,1);
+            εy = tolerance / cos(abs(φ2));                                  // Tolerance for λ (second coordinate).
+            εx = m00*tolerance + m01*εy;                                    // Tolerance for x in user CRS.
+            εy = m10*tolerance + m11*εy;                                    // Tolerance for y in user CRS.
+            /*
+             * Returns the tangent of the ending azimuth converted to the user CRS space.
+             * α2 is azimuth angle in radians, with 0 pointing toward north and values increasing clockwise.
+             * d  is the Jacobian matrix from (φ,λ) to the user coordinate reference system.
+             */
+            final double dx = cos(α2);          // sin(π/2 - α) = -sin(α - π/2) = cos(α)
+            final double dy = sin(α2);          // cos(π/2 - α) = +cos(α - π/2) = sin(α)
+            final double tx = m00*dx + m01*dy;
+            final double ty = m10*dx + m11*dy;
+            return ty / tx;
+        }
+
+        /**
+         * Restores the enclosing {@link GeodeticCalculator} to the state that it has at {@code PathBuilder} instantiation time.
+         */
+        void reset() {
+            φ2 = φf;
+            λ2 = λf;
+            α2 = αf;
+            geodesicDistance = distance;
+            rhumblineLength  = length;
+            validity         = flags;
+        }
     }
 
     /**
diff --git a/core/sis-referencing/src/test/java/org/apache/sis/internal/referencing/j2d/ShapeUtilitiesExt.java b/core/sis-referencing/src/test/java/org/apache/sis/internal/referencing/j2d/ShapeUtilitiesExt.java
new file mode 100644
index 0000000..b33a665
--- /dev/null
+++ b/core/sis-referencing/src/test/java/org/apache/sis/internal/referencing/j2d/ShapeUtilitiesExt.java
@@ -0,0 +1,149 @@
+/*
+ * Licensed to the Apache Software Foundation (ASF) under one or more
+ * contributor license agreements.  See the NOTICE file distributed with
+ * this work for additional information regarding copyright ownership.
+ * The ASF licenses this file to You under the Apache License, Version 2.0
+ * (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.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ */
+package org.apache.sis.internal.referencing.j2d;
+
+import java.awt.Shape;
+import java.awt.geom.Point2D;
+import java.awt.geom.Line2D;
+import java.awt.geom.QuadCurve2D;
+import java.awt.geom.CubicCurve2D;
+import org.opengis.referencing.operation.TransformException;
+
+
+/**
+ * Extensions to {@link ShapeUtilities} not yet needed in main API.
+ *
+ * @author  Martin Desruisseaux (Geomatys)
+ * @version 1.0
+ * @since   1.0
+ * @module
+ */
+public final class ShapeUtilitiesExt {
+    /**
+     * Do not allow instantiation of this class.
+     */
+    private ShapeUtilitiesExt() {
+    }
+
+    /**
+     * Returns a Bézier curve passing by the given points and with the given derivatives at end points.
+     * The cubic curve equation is:
+     *
+     * <blockquote>B(t) = (1-t)³⋅P₁ + 3(1-t)²t⋅C₁ + 3(1-t)t²⋅C₂ + t³⋅P₂</blockquote>
+     *
+     * where t ∈ [0…1], P₁ and P₂ are end points of the curve and C₁ and C₂ are control points generally not on the curve.
+     * If the full equation is required for representing the curve, then this method builds a {@link CubicCurve2D}.
+     * If the same curve can be represented by a quadratic curve, then this method returns a {@link QuadCurve2D}.
+     * If the curve is actually a straight line, then this method returns a {@link Line2D}.
+     *
+     * <p>The (x₁,y₁) arguments give the coordinates of point P₁ at <var>t</var>=0.
+     * The (x<sub>m</sub>,y<sub>m</sub>) arguments give the coordinates of the point at <var>t</var>=½.
+     * The (x₂,y₂) arguments give the coordinates of point P₂ at <var>t</var>=1.</p>
+     *
+     * @param  x1  <var>x</var> value of the starting point.
+     * @param  y1  <var>y</var> value of the starting point.
+     * @param  xm  <var>x</var> value of the mid-point.
+     * @param  ym  <var>y</var> value of the mid-point.
+     * @param  x2  <var>x</var> value of the ending point.
+     * @param  y2  <var>y</var> value of the ending point.
+     * @param  α1  the derivative (∂y/∂x) at starting point.
+     * @param  α2  the derivative (∂y/∂x) at ending point.
+     * @param  εx  maximal distance on <var>x</var> axis between the cubic Bézier curve and quadratic or linear simplifications.
+     * @param  εy  maximal distance on <var>y</var> axis between the cubic Bézier curve and quadratic or linear simplifications.
+     * @return the Bézier curve passing by the 3 given points and having the given derivatives at end points.
+     *
+     * @see <a href="https://pomax.github.io/bezierinfo/">A Primer on Bézier Curves</a>
+     */
+    public static Shape bezier(final double x1, final double y1,
+                               final double xm, final double ym,
+                               final double x2, final double y2,
+                               final double α1, final double α2,
+                               final double εx, final double εy)
+    {
+        final Bezier bezier = new Bezier(2) {
+            @Override protected double evaluateAt(final double t) {
+                final double x, y, α;
+                if (t == 0) {
+                    x = x1;
+                    y = y1;
+                    α = α1;
+                } else if (t == 1) {
+                    x = x2;
+                    y = y2;
+                    α = α2;
+                } else if (t == 0.5) {
+                    x = xm;
+                    y = ym;
+                    α = Double.NaN;
+                } else {
+                    x = Double.NaN;
+                    y = Double.NaN;
+                    α = Double.NaN;
+                }
+                point[0] = x;
+                point[1] = y;
+                return α;
+            }
+        };
+        bezier.εx = εx;
+        bezier.εy = εy;
+        final Shape path;
+        try {
+            path = bezier.build();
+        } catch (TransformException e) {
+            throw new IllegalStateException(e);         // Should never happen.
+        }
+        return ShapeUtilities.toPrimitive(path);
+    }
+
+    /**
+     * Returns a point on the given linear, quadratic or cubic Bézier curve.
+     *
+     * @param  bezier  a {@link Line2D}, {@link QuadCurve2D} or {@link CubicCurve2D}.
+     * @param  t       a parameter from 0 to 1 inclusive.
+     * @return a point on the curve for the given <var>t</var> parameter.
+     *
+     * @see <a href="https://en.wikipedia.org/wiki/B%C3%A9zier_curve">Bézier curve on Wikipedia</a>
+     */
+    public static Point2D.Double pointOnBezier(final Shape bezier, final double t) {
+        final double x, y;
+        final double mt = 1 - t;
+        if (bezier instanceof Line2D) {
+            final Line2D z = (Line2D) bezier;
+            x = mt * z.getX1()  +  t * z.getX2();
+            y = mt * z.getY1()  +  t * z.getY2();
+        } else if (bezier instanceof QuadCurve2D) {
+            final QuadCurve2D z = (QuadCurve2D) bezier;
+            final double a = mt * mt;
+            final double b = mt * t * 2;
+            final double c =  t * t;
+            x = a * z.getX1()  +  b * z.getCtrlX()  +  c * z.getX2();
+            y = a * z.getY1()  +  b * z.getCtrlY()  +  c * z.getY2();
+        } else if (bezier instanceof CubicCurve2D) {
+            final CubicCurve2D z = (CubicCurve2D) bezier;
+            final double a = mt * mt * mt;
+            final double b = mt * mt * t  * 3;
+            final double c = mt * (t * t) * 3;
+            final double d =  t *  t * t;
+            x = a * z.getX1()  +  b * z.getCtrlX1()  +  c * z.getCtrlX2()  +  d * z.getX2();
+            y = a * z.getY1()  +  b * z.getCtrlY1()  +  c * z.getCtrlY2()  +  d * z.getY2();
+        } else {
+            throw new IllegalArgumentException();
+        }
+        return new Point2D.Double(x, y);
+    }
+}
diff --git a/core/sis-referencing/src/test/java/org/apache/sis/internal/referencing/j2d/ShapeUtilitiesTest.java b/core/sis-referencing/src/test/java/org/apache/sis/internal/referencing/j2d/ShapeUtilitiesTest.java
index 21cf9ff..4d94f74 100644
--- a/core/sis-referencing/src/test/java/org/apache/sis/internal/referencing/j2d/ShapeUtilitiesTest.java
+++ b/core/sis-referencing/src/test/java/org/apache/sis/internal/referencing/j2d/ShapeUtilitiesTest.java
@@ -142,7 +142,7 @@ public final strictfp class ShapeUtilitiesTest extends TestCase {
                                                final double x2,  final double y2,
                                                final double α1,  final double α2)
     {
-        final CubicCurve2D p = (CubicCurve2D) ShapeUtilities.bezier(x1, y1, xm, ym, x2, y2, α1, α2, 1, 1);
+        final CubicCurve2D p = (CubicCurve2D) ShapeUtilitiesExt.bezier(x1, y1, xm, ym, x2, y2, α1, α2, 1, 1);
         assertPointEquals( x1,  y1, p.getP1());
         assertPointEquals( x2,  y2, p.getP2());
         assertPointEquals(cx1, cy1, p.getCtrlP1());
@@ -150,7 +150,7 @@ public final strictfp class ShapeUtilitiesTest extends TestCase {
     }
 
     /**
-     * Tests {@link ShapeUtilities#bezier(double, double, double, double, double, double, double, double, double, double)}.
+     * Tests {@link ShapeUtilitiesExt#bezier(double, double, double, double, double, double, double, double, double, double)}.
      * This is an anti-regression test with values computed by {@link ShapeUtilitiesViewer}.
      */
     @Test
@@ -159,7 +159,7 @@ public final strictfp class ShapeUtilitiesTest extends TestCase {
          * Case when the curve can be simplified to a straight line.
          * This test uses a line starting from (100,200) with a slope of 1.3333…
          */
-        final Line2D c1 = (Line2D) ShapeUtilities.bezier(
+        final Line2D c1 = (Line2D) ShapeUtilitiesExt.bezier(
                 100, 200,                           // Start point:  P1
                 175, 300,                           // Midle point:  Pm = P1 + (75,100)
                 250, 400,                           // End point:    P2 = Pm + (75,100)
@@ -183,7 +183,7 @@ public final strictfp class ShapeUtilitiesTest extends TestCase {
          *     C₁  =  ⅓P₁ + ⅔Q  = (300, 266.666…)
          *     C₂  =  ⅓P₂ + ⅔Q  = (450, 400.0)
          */
-        final QuadCurve2D c2 = (QuadCurve2D) ShapeUtilities.bezier(
+        final QuadCurve2D c2 = (QuadCurve2D) ShapeUtilitiesExt.bezier(
                 100,   200,                         // Start point
                 362.5, 350,                         // Midway point
                 550,   600,                         // End point
diff --git a/core/sis-referencing/src/test/java/org/apache/sis/internal/referencing/j2d/ShapeUtilitiesViewer.java b/core/sis-referencing/src/test/java/org/apache/sis/internal/referencing/j2d/ShapeUtilitiesViewer.java
index b2a3e4c..8d917a9 100644
--- a/core/sis-referencing/src/test/java/org/apache/sis/internal/referencing/j2d/ShapeUtilitiesViewer.java
+++ b/core/sis-referencing/src/test/java/org/apache/sis/internal/referencing/j2d/ShapeUtilitiesViewer.java
@@ -190,7 +190,7 @@ public final strictfp class ShapeUtilitiesViewer extends JPanel {
                 input.moveTo(x1, y1);
                 input.lineTo(x4, y4);
                 input.lineTo(x3, y3);
-                output.append(ShapeUtilities.bezier(x1, y1, x2, y2, x3, y3, α1, α2, 1, 1), false);
+                output.append(ShapeUtilitiesExt.bezier(x1, y1, x2, y2, x3, y3, α1, α2, 1, 1), false);
                 out.printf(Locale.ENGLISH, "fitCubicCurve(%d, %d, %d, %d, %d, %d, %g, %g)%n", x1, y1, x2, y2, x3, y3, α1, α2);
                 fillOutput = false;
                 fillInput = false;
diff --git a/core/sis-referencing/src/test/java/org/apache/sis/referencing/GeodeticCalculatorTest.java b/core/sis-referencing/src/test/java/org/apache/sis/referencing/GeodeticCalculatorTest.java
index 0b004fb..70751ec 100644
--- a/core/sis-referencing/src/test/java/org/apache/sis/referencing/GeodeticCalculatorTest.java
+++ b/core/sis-referencing/src/test/java/org/apache/sis/referencing/GeodeticCalculatorTest.java
@@ -17,6 +17,7 @@
 package org.apache.sis.referencing;
 
 import java.awt.Shape;
+import java.awt.geom.Path2D;
 import java.awt.geom.Point2D;
 import java.awt.geom.PathIterator;
 import java.util.Arrays;
@@ -26,11 +27,12 @@ import java.io.LineNumberReader;
 import org.opengis.geometry.DirectPosition;
 import org.opengis.referencing.cs.AxisDirection;
 import org.opengis.referencing.operation.TransformException;
-import org.apache.sis.internal.referencing.j2d.ShapeUtilities;
+import org.apache.sis.internal.referencing.j2d.ShapeUtilitiesExt;
 import org.apache.sis.internal.referencing.Formulas;
 import org.apache.sis.geometry.DirectPosition2D;
 import org.apache.sis.util.CharSequences;
 import org.apache.sis.measure.Units;
+import org.apache.sis.test.widget.VisualCheck;
 import org.apache.sis.test.OptionalTestData;
 import org.apache.sis.test.DependsOnMethod;
 import org.apache.sis.test.TestUtilities;
@@ -223,12 +225,23 @@ public final strictfp class GeodeticCalculatorTest extends TestCase {
     public void testGeodesicPath2D() throws TransformException {
         final GeodeticCalculator c = create(true);
         final double tolerance = 0.05;
-        c.setStartPoint(-33.0, -71.6);                  // Valparaíso
-        c.setEndPoint  ( 31.4, 121.8);                  // Shanghai
-        final Shape path = c.toGeodesicPath2D(1000);
-        assertPointEquals( -71.6, -33.0, ShapeUtilities.pointOnBezier(path, 0),   tolerance);
-        assertPointEquals(-238.2,  31.4, ShapeUtilities.pointOnBezier(path, 1),   tolerance);       // λ₂ = 121.8° - 360°
-        assertPointEquals(-159.2,  -6.8, ShapeUtilities.pointOnBezier(path, 0.5), tolerance);
+        c.setStartPoint(-33.0, -71.6);                                              // Valparaíso
+        c.setEndPoint  ( 31.4, 121.8);                                              // Shanghai
+        final Shape singleCurve = c.toGeodesicPath2D(Double.POSITIVE_INFINITY);
+        final Shape multiCurves = c.toGeodesicPath2D(10000);                        // 10 km tolerance.
+        /*
+         * The approximation done by a single curve is not very good, but is easier to test.
+         */
+        assertPointEquals( -71.6, -33.0, ShapeUtilitiesExt.pointOnBezier(singleCurve, 0),   tolerance);
+        assertPointEquals(-238.2,  31.4, ShapeUtilitiesExt.pointOnBezier(singleCurve, 1),   tolerance);       // λ₂ = 121.8° - 360°
+        assertPointEquals(-159.2,  -6.8, ShapeUtilitiesExt.pointOnBezier(singleCurve, 0.5), tolerance);
+        /*
+         * The more accurate curve can not be simplified to a Java2D primitive.
+         */
+        assertInstanceOf("Multicurves", Path2D.class, multiCurves);
+        if (VisualCheck.SHOW_WIDGET) {
+            VisualCheck.show(singleCurve, multiCurves);
+        }
     }
 
     /**
diff --git a/core/sis-referencing/src/test/java/org/apache/sis/test/widget/ShapeViewer.java b/core/sis-referencing/src/test/java/org/apache/sis/test/widget/ShapeViewer.java
index 7e2d8ef..dc471ea 100644
--- a/core/sis-referencing/src/test/java/org/apache/sis/test/widget/ShapeViewer.java
+++ b/core/sis-referencing/src/test/java/org/apache/sis/test/widget/ShapeViewer.java
@@ -26,7 +26,7 @@ import javax.swing.JPanel;
 
 
 /**
- * A simple viewer of {@link Shape} object. The shape is resized to fill most of the window,
+ * A simple viewer of {@link Shape} objects. The shapes are resized to fill most of the window,
  * with <var>y</var> axis oriented toward up. The bounding box is drawn in gray color behind.
  *
  * @author  Martin Desruisseaux (Geomatys)
@@ -44,14 +44,21 @@ final strictfp class ShapeViewer extends JPanel {
     /**
      * The shape to visualize.
      */
-    private final Shape shape;
+    private final Shape[] shapes;
+
+    /**
+     * Colors to use for drawing shapes.
+     */
+    private final Color[] colors = {
+        Color.RED, Color.GREEN, Color.BLUE
+    };
 
     /**
      * Creates a new panel for rendering the given shape.
      */
-    ShapeViewer(final Shape shape) {
+    ShapeViewer(final Shape[] shapes) {
         setBackground(Color.BLACK);
-        this.shape = shape;
+        this.shapes = shapes;
     }
 
     /**
@@ -61,14 +68,23 @@ final strictfp class ShapeViewer extends JPanel {
     protected void paintComponent(final Graphics graphics) {
         super.paintComponent(graphics);
         final Graphics2D g = (Graphics2D) graphics;
-        Rectangle2D bounds = shape.getBounds2D();
-        g.translate(MARGIN, MARGIN);
-        g.scale((getWidth() - 2*MARGIN) / bounds.getWidth(), (2*MARGIN - getHeight()) / bounds.getHeight());
-        g.translate(-bounds.getMinX(), -bounds.getMaxY());
-        g.setStroke(new BasicStroke(0));
-        g.setColor(Color.GRAY);
-        g.draw(bounds);
-        g.setColor(Color.RED);
-        g.draw(shape);
+        Rectangle2D bounds = null;
+        for (final Shape shape : shapes) {
+            final Rectangle2D b = shape.getBounds2D();
+            if (bounds == null) bounds = b;
+            else bounds.add(b);
+        }
+        if (bounds != null) {
+            g.translate(MARGIN, MARGIN);
+            g.scale((getWidth() - 2*MARGIN) / bounds.getWidth(), (2*MARGIN - getHeight()) / bounds.getHeight());
+            g.translate(-bounds.getMinX(), -bounds.getMaxY());
+            g.setStroke(new BasicStroke(0));
+            g.setColor(Color.GRAY);
+            g.draw(bounds);
+            for (int i=0; i<shapes.length; i++) {
+                g.setColor(colors[i % colors.length]);
+                g.draw(shapes[i]);
+            }
+        }
     }
 }
diff --git a/core/sis-referencing/src/test/java/org/apache/sis/test/widget/VisualCheck.java b/core/sis-referencing/src/test/java/org/apache/sis/test/widget/VisualCheck.java
index 32714c3..a47d5ad 100644
--- a/core/sis-referencing/src/test/java/org/apache/sis/test/widget/VisualCheck.java
+++ b/core/sis-referencing/src/test/java/org/apache/sis/test/widget/VisualCheck.java
@@ -21,6 +21,7 @@ import java.awt.event.WindowAdapter;
 import java.awt.event.WindowEvent;
 import javax.swing.JFrame;
 import javax.swing.JPanel;
+import org.apache.sis.test.TestConfiguration;
 
 
 /**
@@ -34,19 +35,24 @@ import javax.swing.JPanel;
  */
 public final strictfp class VisualCheck {
     /**
+     * Whether to show widgets.
+     */
+    public static final boolean SHOW_WIDGET = Boolean.getBoolean(TestConfiguration.SHOW_WIDGET_KEY);
+
+    /**
      * Do not allows instantiation of this class.
      */
     private VisualCheck() {
     }
 
     /**
-     * Visualizes the given shape in its own window. The shape is resized to fill most of the window,
-     * with <var>y</var> axis oriented toward up. The bounding box is drawn in gray color behind the shape.
+     * Visualizes the given shapes in a window. The shapes are resized to fill most of the window,
+     * with <var>y</var> axis oriented toward up. The bounding box is drawn in gray color behind the shapes.
      *
-     * @param  shape  the shape to visualize.
+     * @param  shapes  the shapes to visualize.
      */
-    public static void show(final Shape shape) {
-        show(new ShapeViewer(shape));
+    public static void show(final Shape... shapes) {
+        show(new ShapeViewer(shapes));
     }
 
     /**
diff --git a/ide-project/NetBeans/nbproject/genfiles.properties b/ide-project/NetBeans/nbproject/genfiles.properties
index 0e08e7c..b11d782 100644
--- a/ide-project/NetBeans/nbproject/genfiles.properties
+++ b/ide-project/NetBeans/nbproject/genfiles.properties
@@ -3,6 +3,6 @@
 build.xml.data.CRC32=58e6b21c
 build.xml.script.CRC32=462eaba0
 build.xml.stylesheet.CRC32=28e38971@1.53.1.46
-nbproject/build-impl.xml.data.CRC32=55754ab4
+nbproject/build-impl.xml.data.CRC32=790389f0
 nbproject/build-impl.xml.script.CRC32=a7689f96
 nbproject/build-impl.xml.stylesheet.CRC32=3a2fa800@1.89.1.48
diff --git a/ide-project/NetBeans/nbproject/project.xml b/ide-project/NetBeans/nbproject/project.xml
index 5b9530d..233ef9b 100644
--- a/ide-project/NetBeans/nbproject/project.xml
+++ b/ide-project/NetBeans/nbproject/project.xml
@@ -84,6 +84,7 @@
             <word>bilevel</word>
             <word>bitmask</word>
             <word>boolean</word>
+            <word>Bézier</word>
             <word>centimetre</word>
             <word>classname</word>
             <word>classnames</word>


Mime
View raw message