Formula or Algorithm to Draw curved lines between points

Question

I'm developing a script to connect point with curved lines.

The points are in ascendent order in x-axis like this:

I'm studying Bezier Curves, but I don't think it's the best solution (see https://www.geogebra.org/m/qcnExXbn):

But, instead of straight lines, I just would like to connect these points in a smooth way.

Could anyone help me with a formula or algorithm?

Answer 1

Given an arbitrary set of $n$ points, it's possible to find the equation of a unique polynomial of degree at most $n-1$ that passes through all $n$ points.

---

Let's consider a very simple example; say we we wanted to find the parabola that passes through the points $(-4,2)$, $(-2,-1)$, and $(1,1)$.

The general equation of a parabola is $$y=ax^2+bx+c$$ so we can substitute the $x$- and $y$-coordinates of the three points we have to create a system of three equations with three unknowns, like so: $$\begin{align} 2&=16a-4b+c \\ -1&=4a-2b+c \\ 1&=a+b+c \end{align}$$ Solving this system, we find that $a=\frac{13}{30}$, $b=\frac{11}{10}$, and $c=\frac{-8}{15}$. The parabola $y=\frac{13}{30}x^2+\frac{11}{10}x-\frac{8}{15}$ does indeed pass cleanly through all three points, so we're done.

---

And of course, you could do this with a system of six equations with six unknowns as well, but that would just take a little longer than this example.

Edit:

Mathematica reveals that the interpolating polynomial for those six points is approximately $$y=1.51771x^5-10.5198x^4+23.6478x^3-17.2688x^2-0.696339x+2.81946$$ The graph of that function looks like this:

Answer 2

I've actually done a lot of curve drawing "by hand" - all the curves in the figures in Complex Made Simple were drawn using PostScript's curveto function, which is Bezier curves. If you try polynomial interpolation I predict you won't like the results - the curve passing though $p_1,\dots p_n$ has funny wiggles near $p_1$ due to the exact location of the other points.

I don't know what you're talking about when you say Bezier curves don't touch the points. A bezier curve is specified by four points; it passes through two of those points exactly, and the tangent vectors at those two endpoints are determined by the two "control points". The interface in terms of "control points" makes no sense to me - I wrote code to convert endpoints and tangent vectors at endpoints to endpoints and control points.

You don't say in what sense you want to "draw" this curve. If you want to write PostScript, here's what to do.

Say you want a curve $c:[0,1]\to\Bbb R^2$ with $c(0)=(x,y)$, $c(1)=(xx,yy)$, $c'(0)=(dx,dy)$ and $c'(1)=(dxx,dyy)$. The following two lines of PostScript give you exactly that:

x y moveto

x+dx/3 y+dy/3 xx-dxx/3 y-dyy/3 xx yy curveto

Not that text literally - you want strings giving the numeric value of expressions. For exammple if $(x,y)=(0,0)$, $(xx,yy)=(1,1)$, $(dx,dy)=(1,1)$ and $(dxx,dyy)=(2,3)$ you'd say

0 0 moveto
0.3333 0.3333 0.3333 0 1 1 curveto

I actually did that to get a curve passing through those points. Then I found this thing doesn't recognize eps as an image format. Here's a screen shot. (This was tedious enough - your points $A,B,\dots$ are just marked with X's.):

Answer 3

Cubic Bezier spline is perfectly suitable to smoothly connect the points.

Given an ordered sequence of $n$ points $p_0,\dots p_{n-1}$, we need to define $n$ cubic Bezier segments. The points on the $k$-th segment are defined parametrically as

\begin{align} s_k(t)=&a_k(1-t)^3+3b_k(1-t)^2t+3c_k(1-t)t^2+d_k t^3 ,\quad t\in[0,1] \tag{1}\label{1} . \end{align}

The first and second derivatives of \eqref{1} are given by \begin{align} s_k'(t)=&3((b_k-a_k)(1-t)^2+2(c_k-b_k)(1-t)t+(d_k-c_k)t^2) \tag{2}\label{2} ,\\ s_k''(t)=&6((a_k-2b_k+c_k)(1-t)+(b_k-2c_k+d_k)t) \tag{3}\label{3} . \end{align}

So, we need to define $a_k,b_k,c_k,d_k$ that satisfy (all indices here are taken $\mod\,n$): \begin{align} s_{k}'(0)=s_{k-1}'(1)&\Rightarrow &\quad b_k-a_k =&d_{k-1}-c_{k-1} \\ s_{k+1}'(0)=s_{k}'(1)&\Rightarrow &\quad b_{k+1}-a_{k+1}=& d_k-c_k \\ s_{k}''(0)=s_{k-1}''(1)&\Rightarrow &\quad c_k-2b_k+a_k=& d_{k-1}-2c_{k-1}+b_{k-1} \end{align} Excluding $c_{k-1}$ and using $d_k=a_{k+1}$, we arrive at $n\times n$ linear system for $b_k$: \begin{align} b_{k-1}+4b_k+b_{k+1}=4a_k+2a_{k+1} \tag{4}\label{4} \\ \text{for }k=0,\dots,n-1 . \end{align}

Then $c_k$ can be easily found:

\begin{align} c_k&=2a_{k+1}-b_{k+1} . \end{align}

Example:

In case if the curve is not closed, you can try this answer.

Answer 4

I created my own code to generate what you're looking for using the math in this video: https://www.youtube.com/watch?v=WCGKqJrf4N4

I required this code for a flutter project using Dart.

class LinearInterpolationPainter extends CustomPainter {
  final List<double> data;
  final double spacing;

  LinearInterpolationPainter({
    this.data = const [],
    this.spacing,
  });

  @override
  void paint(Canvas canvas, Size size) {
    
    //Used to determine the pixels between each data point
    double _spacing = spacing ?? size.width / (data.length - 1);

    var paint = Paint();
    // TODO: Set properties to paint

    paint.color = Colors.green.shade800;
    paint.style = PaintingStyle.stroke;
    paint.strokeWidth = 2.0;

    var path = Path();

    /// Set the path to initial data point
    path.moveTo(0, size.height - data[0]);

    // Draw path
    for (int x = 0; x < data.length - 1; x++) {
      for (int space = 0; space < _spacing; space++) {
        path.lineTo(
          (x * _spacing) + space,
          size.height - interpolate(x + (space / _spacing)),
        );
      }
    }

    /// Draw path to final data point
    path.lineTo(
      (data.length - 1) * _spacing,
      size.height - data.last,
    );

    canvas.drawPath(path, paint);
  }

  @override
  bool shouldRepaint(CustomPainter oldDelegate) {
    return true;
  }

  /// Calculates the y-value for an x-value given a set of points
  double interpolate(double a) {
    double result = 0;

    for (int j = 0; j < data.length; j++) {
      double product = data[j];

      for (int i = 0; i < data.length; i++) {
        if (i != j) product *= (a - i) / (j - i);
      }

      result += product;
    }

    return result;
  }
}

Answer 5

A path built from a series of cubic Bézier curves is fine for this application, you just need to add control points between your main points.

A single cubic Bézier curve is defined by its 2 endpoints and 2 control points. The control points lie on the tangents which pass through the endpoints; the distance between an endpoint and its associated control point affects the curvature.

To make a path passing through multiple points, we use a series of Bézier curves, connected at the endpoints, with control points chosen so that the tangents match where the curves connect (otherwise, we get cusps at the connections).

For the example points in the question, we need 5 curves, one from A to B, one from B to C, etc. Since we just have a set of points, with no tangent (or curvature) data, we construct a "default" Bézier path, using each main point's immediate neighbours to determine the tangents. Like this:

The tangent at B is parallel to AC, the tangent at C is parallel to BD, etc. The point A only has 1 neighbour, so we use AB for its tangent. Similarly, the tangent at F is EF.

The distance from a main point to its associated control point is one third of the (straight line) distance between that curve's endpoints. Eg, the distance from B to its left control point is AB/3, and the distance from B to its right control point is BC/3.

See How to find a bezier curve between two points and two tangents without anchor points for an explanation of the "1/3 rule".

---

Here's the Sage / Python script I used to create that diagram. SageMath is a vast mathematics system built on top of Python. It has nice plotting capabilities (courtesy of matplotlib), and its built-in vector type is handy for this kind of task. If you can read plain Python, (hopefully) you shouldn't have problems reading my script.

from itertools import chain 

def bezier_fd(pos):
    "Bezier path plot, using finite differences"
    # Make the normalized tangent vectors
    tgt = [v - u for u, v in zip(pos, pos[2:])]
    if pos[0] == pos[-1]:
        end = [pos[1] - pos[-2]]
        tgt = end + tgt + end
    else:
        tgt = [pos[1] - pos[0]] + tgt + [pos[-1] - pos[-2]]
    tgt = [u.normalized() for u in tgt] 

    # Build Bezier curves
    pos_tgt = zip(pos, tgt)
    p0, t0 = next(pos_tgt)
    row = [p0]
    curves = []
    for p, t in pos_tgt:
        s = abs(p - p0) / 3
        row.extend([p0 + t0*s, p - t*s, p])
        curves.append(row)
        row = []
        p0, t0 = p, t
    return curves 

pos = [
  (0.5, 0.5),
  (1, -0.5),
  (1.5, 1),
  (2.25, 1.1),
  (2.6, -0.5),
  (3, 0.5),
] 

# Plot the main points
pos = [vector(p) for p in pos]
P = points(pos, color="blue") 

# Label the points
ne = vector((0.1, 0.1))
nw = vector((-0.1, 0.1))
offsets = (nw, nw, nw, ne, ne, ne)
P += sum(text(c, p + delta)
  for c, p, delta in zip("ABCDEF", pos, offsets)) 

# Compute & plot the bezier curves
curves = bezier_fd(pos)
P += bezier_path(curves, color="red") 

# The main points & the control points in a flat list
pts = list(chain.from_iterable(curves)) 

# The tangents
P += sum(line(t, color="#333", linestyle="--")
  for t in chain([pts[:2]], zip(pts[2::3], pts[4::3]), [pts[-2:]])) 

# The rest of the hull
P += sum(line(t, color="#333", linestyle=":")
  for t in zip(pts[1::3], pts[2::3])) 

# The control points
P += points([p for i, p in enumerate(pts) if i%3], color="magenta") 

P.show(xmin=0, aspect_ratio=1, ticks=(0.5, 0.5), gridlines="minor")

Here's a live version of the script, running on the SageMathCell server.

And here's an SVG version of the basic path, without visible control points, etc.

Source code

Just for fun, here's a single interactive Bézier curve. You can drag points around and watch what happens. (It should work with a mouse or touchscreen). You can run it on the server, or download it to run it locally.