View Javadoc
1   package org.opentrafficsim.core.geometry;
2   
3   import java.awt.geom.Line2D;
4   
5   import org.djutils.draw.point.OrientedPoint2d;
6   import org.djutils.draw.point.Point2d;
7   import org.djutils.exceptions.Throw;
8   import org.opentrafficsim.base.geometry.OtsGeometryException;
9   import org.opentrafficsim.base.geometry.OtsLine2d;
10  
11  /**
12   * Generation of B&eacute;zier curves. <br>
13   * The class implements the cubic(...) method to generate a cubic B&eacute;zier curve using the following formula: B(t) = (1 -
14   * t)<sup>3</sup>P<sub>0</sub> + 3t(1 - t)<sup>2</sup> P<sub>1</sub> + 3t<sup>2</sup> (1 - t) P<sub>2</sub> + t<sup>3</sup>
15   * P<sub>3</sub> where P<sub>0</sub> and P<sub>3</sub> are the end points, and P<sub>1</sub> and P<sub>2</sub> the control
16   * points. <br>
17   * For a smooth movement, one of the standard implementations of the cubic(...) function offered is the case where P<sub>1</sub>
18   * is positioned halfway between P<sub>0</sub> and P<sub>3</sub> starting from P<sub>0</sub> in the direction of P<sub>3</sub>,
19   * and P<sub>2</sub> is positioned halfway between P<sub>3</sub> and P<sub>0</sub> starting from P<sub>3</sub> in the direction
20   * of P<sub>0</sub>.<br>
21   * Finally, an n-point generalization of the B&eacute;zier curve is implemented with the bezier(...) function.
22   * <p>
23   * Copyright (c) 2013-2024 Delft University of Technology, PO Box 5, 2600 AA, Delft, the Netherlands. All rights reserved. <br>
24   * BSD-style license. See <a href="https://opentrafficsim.org/docs/license.html">OpenTrafficSim License</a>.
25   * </p>
26   * @author <a href="https://github.com/averbraeck">Alexander Verbraeck</a>
27   * @author <a href="https://github.com/peter-knoppers">Peter Knoppers</a>
28   */
29  public final class Bezier
30  {
31      /** The default number of points to use to construct a B&eacute;zier curve. */
32      private static final int DEFAULT_NUM_POINTS = 64;
33  
34      /** Cached factorial values. */
35      private static long[] fact = new long[] {1L, 1L, 2L, 6L, 24L, 120L, 720L, 5040L, 40320L, 362880L, 3628800L, 39916800L,
36              479001600L, 6227020800L, 87178291200L, 1307674368000L, 20922789888000L, 355687428096000L, 6402373705728000L,
37              121645100408832000L, 2432902008176640000L};
38  
39      /** Utility class. */
40      private Bezier()
41      {
42          // do not instantiate
43      }
44  
45      /**
46       * Construct a cubic B&eacute;zier curve from start to end with two control points.
47       * @param numPoints the number of points for the B&eacute;zier curve
48       * @param start the start point of the B&eacute;zier curve
49       * @param control1 the first control point
50       * @param control2 the second control point
51       * @param end the end point of the B&eacute;zier curve
52       * @return a cubic B&eacute;zier curve between start and end, with the two provided control points
53       * @throws OtsGeometryException in case the number of points is less than 2 or the B&eacute;zier curve could not be
54       *             constructed
55       */
56      public static OtsLine2d cubic(final int numPoints, final Point2d start, final Point2d control1, final Point2d control2,
57              final Point2d end) throws OtsGeometryException
58      {
59          Throw.when(numPoints < 2, OtsGeometryException.class, "Number of points too small (got %d; minimum value is 2)",
60                  numPoints);
61          Point2d[] points = new Point2d[numPoints];
62          for (int n = 0; n < numPoints; n++)
63          {
64              double t = n / (numPoints - 1.0);
65              double x = B3(t, start.x, control1.x, control2.x, end.x);
66              double y = B3(t, start.y, control1.y, control2.y, end.y);
67              points[n] = new Point2d(x, y);
68          }
69          return new OtsLine2d(points);
70      }
71  
72      /**
73       * Construct a cubic B&eacute;zier curve from start to end with two generated control points at half the distance between
74       * start and end. The z-value is interpolated in a linear way.
75       * @param numPoints the number of points for the B&eacute;zier curve
76       * @param start the directed start point of the B&eacute;zier curve
77       * @param end the directed end point of the B&eacute;zier curve
78       * @return a cubic B&eacute;zier curve between start and end, with the two provided control points
79       * @throws OtsGeometryException in case the number of points is less than 2 or the B&eacute;zier curve could not be
80       *             constructed
81       */
82      public static OtsLine2d cubic(final int numPoints, final OrientedPoint2d start, final OrientedPoint2d end)
83              throws OtsGeometryException
84      {
85          return cubic(numPoints, start, end, 1.0);
86      }
87  
88      /**
89       * Construct a cubic B&eacute;zier curve from start to end with two generated control points at half the distance between
90       * start and end. The z-value is interpolated in a linear way.
91       * @param numPoints the number of points for the B&eacute;zier curve
92       * @param start the directed start point of the B&eacute;zier curve
93       * @param end the directed end point of the B&eacute;zier curve
94       * @param shape 1 = control points at half the distance between start and end, &gt; 1 results in a pointier shape, &lt; 1
95       *            results in a flatter shape, value should be above 0
96       * @return a cubic B&eacute;zier curve between start and end, with the two determined control points
97       * @throws OtsGeometryException in case the number of points is less than 2 or the B&eacute;zier curve could not be
98       *             constructed
99       */
100     public static OtsLine2d cubic(final int numPoints, final OrientedPoint2d start, final OrientedPoint2d end,
101             final double shape) throws OtsGeometryException
102     {
103         return cubic(numPoints, start, end, shape, false);
104     }
105 
106     /**
107      * Construct a cubic B&eacute;zier curve from start to end with two generated control points at half the distance between
108      * start and end. The z-value is interpolated in a linear way.
109      * @param numPoints the number of points for the B&eacute;zier curve
110      * @param start the directed start point of the B&eacute;zier curve
111      * @param end the directed end point of the B&eacute;zier curve
112      * @param shape shape factor; 1 = control points at half the distance between start and end, &gt; 1 results in a pointier
113      *            shape, &lt; 1 results in a flatter shape, value should be above 0
114      * @param weighted control point distance relates to distance to projected point on extended line from other end
115      * @return a cubic B&eacute;zier curve between start and end, with the two determined control points
116      * @throws OtsGeometryException in case the number of points is less than 2 or the B&eacute;zier curve could not be
117      *             constructed
118      */
119     public static OtsLine2d cubic(final int numPoints, final OrientedPoint2d start, final OrientedPoint2d end,
120             final double shape, final boolean weighted) throws OtsGeometryException
121     {
122         return bezier(cubicControlPoints(start, end, shape, weighted));
123     }
124 
125     /**
126      * Construct control points for a cubic B&eacute;zier curve from start to end with two generated control points at half the
127      * distance between start and end.
128      * @param start the directed start point of the B&eacute;zier curve
129      * @param end the directed end point of the B&eacute;zier curve
130      * @param shape shape factor; 1 = control points at half the distance between start and end, &gt; 1 results in a pointier
131      *            shape, &lt; 1 results in a flatter shape, value should be above 0
132      * @param weighted control point distance relates to distance to projected point on extended line from other end
133      * @return a cubic B&eacute;zier curve between start and end, with the two determined control points
134      */
135     public static Point2d[] cubicControlPoints(final OrientedPoint2d start, final OrientedPoint2d end, final double shape,
136             final boolean weighted)
137     {
138         Throw.when(shape <= 0.0, IllegalArgumentException.class, "Shape factor must be above 0.0.");
139         Point2d control1;
140         Point2d control2;
141 
142         if (weighted)
143         {
144             // each control point is 'w' * the distance between the end-points away from the respective end point
145             // 'w' is a weight given by the distance from the end point to the extended line of the other end point
146             double dx = end.x - start.x;
147             double dy = end.y - start.y;
148             double distance = shape * Math.hypot(dx, dy);
149             double cosEnd = Math.cos(end.getDirZ());
150             double sinEnd = Math.sin(end.getDirZ());
151             double dStart = Line2D.ptLineDist(end.x, end.y, end.x + cosEnd, end.y + sinEnd, start.x, start.y);
152             double cosStart = Math.cos(start.getDirZ());
153             double sinStart = Math.sin(start.getDirZ());
154             double dEnd = Line2D.ptLineDist(start.x, start.y, start.x + cosStart, start.y + sinStart, end.x, end.y);
155             double wStart = dStart / (dStart + dEnd);
156             double wEnd = dEnd / (dStart + dEnd);
157             double wStartDistance = wStart * distance;
158             double wEndDistance = wEnd * distance;
159             control1 = new Point2d(start.x + wStartDistance * cosStart, start.y + wStartDistance * sinStart);
160             // - (minus) as the angle is where the line leaves, i.e. from shape point to end
161             control2 = new Point2d(end.x - wEndDistance * cosEnd, end.y - wEndDistance * sinEnd);
162         }
163         else
164         {
165             // each control point is half the distance between the end-points away from the respective end point
166             double dx = end.x - start.x;
167             double dy = end.y - start.y;
168             double distance2 = shape * .5 * Math.hypot(dx, dy);
169             control1 = new Point2d(start.x + distance2 * Math.cos(start.getDirZ()),
170                     start.y + distance2 * Math.sin(start.getDirZ()));
171             control2 = new Point2d(end.x - distance2 * Math.cos(end.getDirZ()), end.y - distance2 * Math.sin(end.getDirZ()));
172         }
173 
174         return new Point2d[] {start, control1, control2, end};
175     }
176 
177     /**
178      * Construct a cubic B&eacute;zier curve from start to end with two generated control points at half the distance between
179      * start and end. The z-value is interpolated in a linear way.
180      * @param start the directed start point of the B&eacute;zier curve
181      * @param end the directed end point of the B&eacute;zier curve
182      * @return a cubic B&eacute;zier curve between start and end, with the two provided control points
183      * @throws OtsGeometryException in case the number of points is less than 2 or the B&eacute;zier curve could not be
184      *             constructed
185      */
186     public static OtsLine2d cubic(final OrientedPoint2d start, final OrientedPoint2d end) throws OtsGeometryException
187     {
188         return cubic(DEFAULT_NUM_POINTS, start, end);
189     }
190 
191     /**
192      * Calculate the cubic B&eacute;zier point with B(t) = (1 - t)<sup>3</sup>P<sub>0</sub> + 3t(1 - t)<sup>2</sup>
193      * P<sub>1</sub> + 3t<sup>2</sup> (1 - t) P<sub>2</sub> + t<sup>3</sup> P<sub>3</sub>.
194      * @param t the fraction
195      * @param p0 the first point of the curve
196      * @param p1 the first control point
197      * @param p2 the second control point
198      * @param p3 the end point of the curve
199      * @return the cubic bezier value B(t)
200      */
201     @SuppressWarnings("checkstyle:methodname")
202     private static double B3(final double t, final double p0, final double p1, final double p2, final double p3)
203     {
204         double t2 = t * t;
205         double t3 = t2 * t;
206         double m = (1.0 - t);
207         double m2 = m * m;
208         double m3 = m2 * m;
209         return m3 * p0 + 3.0 * t * m2 * p1 + 3.0 * t2 * m * p2 + t3 * p3;
210     }
211 
212     /**
213      * Construct a B&eacute;zier curve of degree n.
214      * @param numPoints the number of points for the B&eacute;zier curve to be constructed
215      * @param points the points of the curve, where the first and last are begin and end point, and the intermediate ones are
216      *            control points. There should be at least two points.
217      * @return the B&eacute;zier value B(t) of degree n, where n is the number of points in the array
218      * @throws OtsGeometryException in case the number of points is less than 2 or the B&eacute;zier curve could not be
219      *             constructed
220      */
221     public static OtsLine2d bezier(final int numPoints, final Point2d... points) throws OtsGeometryException
222     {
223         Point2d[] result = new Point2d[numPoints];
224         double[] px = new double[points.length];
225         double[] py = new double[points.length];
226         for (int i = 0; i < points.length; i++)
227         {
228             px[i] = points[i].x;
229             py[i] = points[i].y;
230         }
231         for (int n = 0; n < numPoints; n++)
232         {
233             double t = n / (numPoints - 1.0);
234             double x = Bn(t, px);
235             double y = Bn(t, py);
236             result[n] = new Point2d(x, y);
237         }
238         return new OtsLine2d(result);
239     }
240 
241     /**
242      * Construct a B&eacute;zier curve of degree n.
243      * @param points the points of the curve, where the first and last are begin and end point, and the intermediate ones are
244      *            control points. There should be at least two points.
245      * @return the B&eacute;zier value B(t) of degree n, where n is the number of points in the array
246      * @throws OtsGeometryException in case the number of points is less than 2 or the B&eacute;zier curve could not be
247      *             constructed
248      */
249     public static OtsLine2d bezier(final Point2d... points) throws OtsGeometryException
250     {
251         return bezier(DEFAULT_NUM_POINTS, points);
252     }
253 
254     /**
255      * Calculate the B&eacute;zier point of degree n, with B(t) = Sum(i = 0..n) [C(n, i) * (1 - t)<sup>n-i</sup> t<sup>i</sup>
256      * P<sub>i</sub>], where C(n, k) is the binomial coefficient defined by n! / ( k! (n-k)! ), ! being the factorial operator.
257      * @param t the fraction
258      * @param p the points of the curve, where the first and last are begin and end point, and the intermediate ones are control
259      *            points
260      * @return the B&eacute;zier value B(t) of degree n, where n is the number of points in the array
261      */
262     @SuppressWarnings("checkstyle:methodname")
263     static double Bn(final double t, final double... p)
264     {
265         double b = 0.0;
266         double m = (1.0 - t);
267         int n = p.length - 1;
268         double fn = factorial(n);
269         for (int i = 0; i <= n; i++)
270         {
271             double c = fn / (factorial(i) * (factorial(n - i)));
272             b += c * Math.pow(m, n - i) * Math.pow(t, i) * p[i];
273         }
274         return b;
275     }
276 
277     /**
278      * Calculate factorial(k), which is k * (k-1) * (k-2) * ... * 1. For factorials up to 20, a lookup table is used.
279      * @param k the parameter
280      * @return factorial(k)
281      */
282     private static double factorial(final int k)
283     {
284         if (k < fact.length)
285         {
286             return fact[k];
287         }
288         double f = 1;
289         for (int i = 2; i <= k; i++)
290         {
291             f = f * i;
292         }
293         return f;
294     }
295 
296 }