This document describes an algorithm for generating offset curves for X-splines, where the width of the X-spline is allowed to vary along its length. The implementation is provided by the grid.offsetXspline function in the 'vwline' package for R.

1. Introduction

X-splines (Blanc and Schlick, 1995) are a family of curves defined by control points (like a Bezier curve), with a shape parameter at each control point that allows the curve to vary between interpolation and approximation at each control point. In the diagram below, there are three X-splines, each drawn relative to five control points. The top X-spline has a shape parameter of 1 (except at the end points) so it approximates the control points. The bottom X-spline has a shape parameter of -1 (except at the end points) so it interpolates the control points. The middle X-spline has a shape parameter of 0, so it produces straight lines and sharp corners.

The R graphics system provides functions for drawing X-splines to produce smooth lines or closed shapes with smooth boundaries. The following code shows an example (used for the top X-spline in the diagram above).

x <- c(.1, .3, .5, .7, .9)
y <- c(.4, .6, .4, .6, .4)
grid.xspline(x, y, shape=1)

The R package 'vwline' provides several functions that draw curves based on X-splines, but with the width of the X-spline allowed to vary along the length of the curve. The following R code shows an example, where the width of the X-spline increases smoothly along the length of the line.

w <- c(0, .1)
grid.offsetXspline(x, y, w, shape=1)

In order to render a variable-width X-spline, we must be able to calculate locations along the boundary of the variable-width X-spline. This document focuses on the function grid.offsetXspline from the 'vwline' package, which uses a generalised offset curve approach to calculate the boundary of a variable-width X-spline.

In the Offset curves Section, we establish a general definition of offset curves and the X-spline curves Section describes the equations for X-spline curves. The X-spline tangent functions Section generates expressions for the tangent to an X-spline curve and the Rendering X-spline offset curves Section describes how the tangent is used to draw an offset curve.

2. Offset curves

Following Lin and Chen, 2014, if we have a planar (two-dimensional) parametric curve $\mathbf{r}=\mathbf{r}(t)$, $t \in [0, 1]$, then an offset curve with variable offset is given by $$\mathbf{r}_0(t)=\mathbf{r}(t) + d(t)\mathbf{n}(t)$$ where $\mathbf{n}(t)$ is the unit normal vector at each point on the original curve and $d(t)$ is a function that defines the offset at any point on the original curve.

In order to find the unit normal function, we must first obtain the unit tangent function. This requires differentiating the original curve and then dividing by the magnitude of the derivative: $$\mathbf{e}(t) = \frac{\mathbf{r}'(t)}{\left\Vert\mathbf{r}'(t)\right\Vert}$$

Because the unit tangent function has fixed length, its derivative is perpendicular to it, so the unit normal function is then found by (Apostol, 2007): $$\mathbf{n}(t) = \frac{\mathbf{e}'(t)}{\left\Vert\mathbf{e}'(t)\right\Vert}$$

We can also observe that the unit normal is just a 90 degree rotation of the unit tangent, so we can obtain the unit normal function simply as: $$\begin{array}{l} n_x(t) = e_y(t) \\ n_y(t) = -e_x(t) \\ \end{array}$$

3. X-spline curves

Following Blanc and Schlick, 1995, a (two-dimensional) spline is a planar parametric curve that is based on a discrete set of control points, $P_k \in \mathbb{R}^2$, and a discrete set of blending functions, $F_k : [0, 1] \rightarrow \mathbb{R}$, as follows: $$C(t) = \sum_{k=0}^n F_k(t)P_k, t \in [0, 1]$$

The definition of an X-spline reparameterises the curve and adds a set of shape parameters, $s_k \in [-1, 1]$, which are used to both select and influence the blending functions. The X-spline curve is calculated piecewise for each (overlapping) set of four control points. The section of the curve between control points $k+1$ and $k+2$ is calculated as follows: $$C(t) = \frac{A_0(t)P_k + A_1(t)P_{k+1} + A_2(t)P_{k+2} + A_3(t)P_{k+3}}{A_0(t) + A_1(t) + A_2(t) + A_3(t)}, t \in [0, 1]$$ where the value and choice of blending functons used for $A_0$ to $A_3$ depend on the shape parameters $s_{k+1}$ and $s_{k+2}$: $$\begin{array}{l} A_0 = & \begin{cases} h(-t, -s_{k+1}), \\ f(t - s_{k+1}, -1 - s_{k+1}), \\ 0, \end{cases} & \begin{array}{l} \text{if $s_{k+1} \lt 0$} \\ \text{if $s_{k+1} \ge 0$ and $t \lt s_{k+1}$} \\ \text{otherwise} \end{array} \\[1em] A_1 = & \begin{cases} g(1 - t, -s_{k+2}), \\ f(t - 1 - s_{k+2}, -1 - s_{k+2}), \end{cases} & \begin{array}{l} \text{if $s_{k+2} \lt 0$} \\ \text{otherwise} \end{array} \\[1em] A_2 = & \begin{cases} g(t, -s_{k+1}), \\ f(t + s_{k+1}, 1 + s_{k+1}), \end{cases} & \begin{array}{l} \text{if $s_{k+1} \lt 0$} \\ \text{otherwise} \end{array} \\[1em] A_3 = & \begin{cases} h(t - 1, -s_{k+2}), \\ f(t - 1 + s_{k+2}, 1 + s_{k+2}), \\ 0, \end{cases} & \begin{array}{l} \text{if $s_{k+2} \lt 0$} \\ \text{if $s_{k+2} \ge 0$ and $t \gt 1 - s_{k+2}$} \\ \text{otherwise} \end{array} \end{array}$$ and the blending functions $f$, $g$, and $h$ are defined as follows: $$\begin{array}{l} F(u, p) = u^3(10 - p + (2p - 15)u + (6 - p)u^2) \\ f(n, d) = F(n/d, 2d^2) \\ g(u, q) = u(q + u(2q + u(8 - 12q + u(14q - 11 + u(4 - 5q))))) \\ h(u, q) = u(q + u(2q + u^2(-2q - uq))) \\ \end{array}$$

The notation above, particularly the definition of the blending function $f$, is also based on the C code implementation of X-splines from the XFig drawing program (Smith, 2001).

4. X-spline tangent functions

In order to generate an offset curve for an X-spline, we need to differentiate the function $C(t)$ to find its tangent function, $C'(t)$.

In the absence of the requisite discipline and mathematical wit, one approach to calculating the X-spline tangent function is to brute force the result via computational methods. This approach also has the benefit of directly producing output in the form of R code that can be evaluated and, because the X-spline equation is piecewise, it efficiently replicates across multiple equations.

The first step is to generate an expression for $C(t)$ that is just in terms of $t$ (and $P_k$ and $s_k$, both of which are constants for a particular X-spline). The idea is to expand the X-spline equations so that, for example ... $$A0 = h(-t, -s_{k+1})$$ ... and ... $$h(u, q) = u(q + u(2q + u^2(-2q - uq)))$$ ... becomes ... $$A0 = (-t)((-s_{k+1}) + (-t)(2(-s_{k+1}) + (-t)^2(-2(-s_{k+1}) - (-t)(-s_{k+1}))))$$

The following R code performs the above transformation using text substitution and demonstrates that the result can be executed (for specific values of t and s1).

hblend <- "u*(q + u*(2*q + u*u*(-2*q - u*q)))"
A0 <- gsub("u", "(-t)",
       gsub("q", "(-s1)", hblend))
A0

  [1] "(-t)*((-s1) + (-t)*(2*(-s1) + (-t)*(-t)*(-2*(-s1) - (-t)*(-s1))))"

eval(parse(text=A0), list(t=0.5, s1=1))

  [1] 0.09375

To get a complete expression for $C(t)$, we must repeat this sort of expansion for $A1$, $A2$, and $A3$ as well, and there are 16 different scenarios to create equations for: $s_{k+1}$ can be negative or non-negative, $s_{k+2}$ can be negative or non-negative, $A0$ can be zero or non-zero, and $A3$ can be zero or non-zero. Furthermore, we need to generate both $x$ and $y$ versions to provide the two components of the vector function $C(t)$. Not all of those scenarios can actually occur, but because we are programmatically generating the equations, it is easier, and costs no more, to generate all 16 scenarios.

The xsplineFunGenerator function was written to generate an expression for a specific scenario. The following code demonstrates its use for the scenario where $s_{k+1}$ is negative, $s_{k+2}$ is negative, $A0$ is zero, and $A3$ is zero. The expression that is generated is the equation for $C_x(t)$ in that scenario (the px terms in the expression are the x-values from the four control points that are controlling the curve, e.g., px0 corresponds to the x-component of $P_k$; s1 and s2 correspond to $s_{k+1}$ and $s_{k+2}$).

xsplineFunGenerator("s1neg", "s2neg", "noA0", "noA3")("x")

  expression(((-t) * ((-s1) + (-t) * (2 * (-s1) + (-t) * (-t) * 
      (-2 * (-s1) - (-t) * (-s1)))) * px0 + (1 - t) * ((-s2) + 
      (1 - t) * (2 * (-s2) + (1 - t) * (8 - 12 * (-s2) + (1 - t) * 
          (14 * (-s2) - 11 + (1 - t) * (4 - 5 * (-s2)))))) * px1 + 
      t * ((-s1) + t * (2 * (-s1) + t * (8 - 12 * (-s1) + t * (14 * 
          (-s1) - 11 + t * (4 - 5 * (-s1)))))) * px2 + (t - 1) * 
      ((-s2) + (t - 1) * (2 * (-s2) + (t - 1) * (t - 1) * (-2 * 
          (-s2) - (t - 1) * (-s2)))) * px3)/((-t) * ((-s1) + (-t) * 
      (2 * (-s1) + (-t) * (-t) * (-2 * (-s1) - (-t) * (-s1)))) + 
      (1 - t) * ((-s2) + (1 - t) * (2 * (-s2) + (1 - t) * (8 - 
          12 * (-s2) + (1 - t) * (14 * (-s2) - 11 + (1 - t) * (4 - 
          5 * (-s2)))))) + t * ((-s1) + t * (2 * (-s1) + t * (8 - 
      12 * (-s1) + t * (14 * (-s1) - 11 + t * (4 - 5 * (-s1)))))) + 
      (t - 1) * ((-s2) + (t - 1) * (2 * (-s2) + (t - 1) * (t - 
          1) * (-2 * (-s2) - (t - 1) * (-s2))))))

A full set of equations for $C(t)$ for all scenarios can be generated with a simple loop. A link to the complete definition of xsplineFunGenerator is given in the Resources.

Now that we have an equation for $C(t)$ in terms of $t$, we need to find its derivative. Again, we can take a computational approach using the D (deriv) function from the 'Ryacas' package (Goedman et al., 2016), which provides an interface to the Yacas computer algebra system (Pinkus and Winitzki, 2002). The xsplineTangentExpr function was written to take the output of xsplineFunGenerator and differentiate it. The following code demonstrates (the x-component of) the result from xsplineTangentExpr for the scenario above.

xsplineTangentExpr(xsplineFunGenerator("s1neg", "s2neg", "noA0", "noA3"))$x

  (((-t) * ((-t) * ((-t) * (-t) * (-s1) - ((-t) + (-t)) * (-2 * 
      (-s1) - (-t) * (-s1))) - (2 * (-s1) + (-t) * (-t) * (-2 * 
      (-s1) - (-t) * (-s1)))) - ((-s1) + (-t) * (2 * (-s1) + (-t) * 
      (-t) * (-2 * (-s1) - (-t) * (-s1))))) * px0 - ((1 - t) * 
      ((1 - t) * ((1 - t) * ((1 - t) * (4 - 5 * (-s2)) + (14 * 
          (-s2) - 11 + (1 - t) * (4 - 5 * (-s2)))) + (8 - 12 * 
          (-s2) + (1 - t) * (14 * (-s2) - 11 + (1 - t) * (4 - 5 * 
          (-s2))))) + (2 * (-s2) + (1 - t) * (8 - 12 * (-s2) + 
          (1 - t) * (14 * (-s2) - 11 + (1 - t) * (4 - 5 * (-s2)))))) + 
      ((-s2) + (1 - t) * (2 * (-s2) + (1 - t) * (8 - 12 * (-s2) + 
          (1 - t) * (14 * (-s2) - 11 + (1 - t) * (4 - 5 * (-s2))))))) * 
      px1 + (((-s1) + t * (2 * (-s1) + t * (8 - 12 * (-s1) + t * 
      (14 * (-s1) - 11 + t * (4 - 5 * (-s1)))))) + t * ((2 * (-s1) + 
      t * (8 - 12 * (-s1) + t * (14 * (-s1) - 11 + t * (4 - 5 * 
          (-s1))))) + t * ((8 - 12 * (-s1) + t * (14 * (-s1) - 
      11 + t * (4 - 5 * (-s1)))) + t * ((14 * (-s1) - 11 + t * 
      (4 - 5 * (-s1))) + t * (4 - 5 * (-s1)))))) * px2 + (((-s2) + 
      (t - 1) * (2 * (-s2) + (t - 1) * (t - 1) * (-2 * (-s2) - 
          (t - 1) * (-s2)))) + (t - 1) * ((2 * (-s2) + (t - 1) * 
      (t - 1) * (-2 * (-s2) - (t - 1) * (-s2))) + (t - 1) * (((t - 
      1) + (t - 1)) * (-2 * (-s2) - (t - 1) * (-s2)) - (t - 1) * 
      (t - 1) * (-s2)))) * px3)/((-t) * ((-s1) + (-t) * (2 * (-s1) + 
      (-t) * (-t) * (-2 * (-s1) - (-t) * (-s1)))) + (1 - t) * ((-s2) + 
      (1 - t) * (2 * (-s2) + (1 - t) * (8 - 12 * (-s2) + (1 - t) * 
          (14 * (-s2) - 11 + (1 - t) * (4 - 5 * (-s2)))))) + t * 
      ((-s1) + t * (2 * (-s1) + t * (8 - 12 * (-s1) + t * (14 * 
          (-s1) - 11 + t * (4 - 5 * (-s1)))))) + (t - 1) * ((-s2) + 
      (t - 1) * (2 * (-s2) + (t - 1) * (t - 1) * (-2 * (-s2) - 
          (t - 1) * (-s2))))) - ((-t) * ((-s1) + (-t) * (2 * (-s1) + 
      (-t) * (-t) * (-2 * (-s1) - (-t) * (-s1)))) * px0 + (1 - 
      t) * ((-s2) + (1 - t) * (2 * (-s2) + (1 - t) * (8 - 12 * 
      (-s2) + (1 - t) * (14 * (-s2) - 11 + (1 - t) * (4 - 5 * (-s2)))))) * 
      px1 + t * ((-s1) + t * (2 * (-s1) + t * (8 - 12 * (-s1) + 
      t * (14 * (-s1) - 11 + t * (4 - 5 * (-s1)))))) * px2 + (t - 
      1) * ((-s2) + (t - 1) * (2 * (-s2) + (t - 1) * (t - 1) * 
      (-2 * (-s2) - (t - 1) * (-s2)))) * px3) * ((-t) * ((-t) * 
      ((-t) * (-t) * (-s1) - ((-t) + (-t)) * (-2 * (-s1) - (-t) * 
          (-s1))) - (2 * (-s1) + (-t) * (-t) * (-2 * (-s1) - (-t) * 
      (-s1)))) - ((-s1) + (-t) * (2 * (-s1) + (-t) * (-t) * (-2 * 
      (-s1) - (-t) * (-s1)))) - ((1 - t) * ((1 - t) * ((1 - t) * 
      ((1 - t) * (4 - 5 * (-s2)) + (14 * (-s2) - 11 + (1 - t) * 
          (4 - 5 * (-s2)))) + (8 - 12 * (-s2) + (1 - t) * (14 * 
      (-s2) - 11 + (1 - t) * (4 - 5 * (-s2))))) + (2 * (-s2) + 
      (1 - t) * (8 - 12 * (-s2) + (1 - t) * (14 * (-s2) - 11 + 
          (1 - t) * (4 - 5 * (-s2)))))) + ((-s2) + (1 - t) * (2 * 
      (-s2) + (1 - t) * (8 - 12 * (-s2) + (1 - t) * (14 * (-s2) - 
      11 + (1 - t) * (4 - 5 * (-s2))))))) + (((-s1) + t * (2 * 
      (-s1) + t * (8 - 12 * (-s1) + t * (14 * (-s1) - 11 + t * 
      (4 - 5 * (-s1)))))) + t * ((2 * (-s1) + t * (8 - 12 * (-s1) + 
      t * (14 * (-s1) - 11 + t * (4 - 5 * (-s1))))) + t * ((8 - 
      12 * (-s1) + t * (14 * (-s1) - 11 + t * (4 - 5 * (-s1)))) + 
      t * ((14 * (-s1) - 11 + t * (4 - 5 * (-s1))) + t * (4 - 5 * 
          (-s1)))))) + (((-s2) + (t - 1) * (2 * (-s2) + (t - 1) * 
      (t - 1) * (-2 * (-s2) - (t - 1) * (-s2)))) + (t - 1) * ((2 * 
      (-s2) + (t - 1) * (t - 1) * (-2 * (-s2) - (t - 1) * (-s2))) + 
      (t - 1) * (((t - 1) + (t - 1)) * (-2 * (-s2) - (t - 1) * 
          (-s2)) - (t - 1) * (t - 1) * (-s2)))))/((-t) * ((-s1) + 
      (-t) * (2 * (-s1) + (-t) * (-t) * (-2 * (-s1) - (-t) * (-s1)))) + 
      (1 - t) * ((-s2) + (1 - t) * (2 * (-s2) + (1 - t) * (8 - 
          12 * (-s2) + (1 - t) * (14 * (-s2) - 11 + (1 - t) * (4 - 
          5 * (-s2)))))) + t * ((-s1) + t * (2 * (-s1) + t * (8 - 
      12 * (-s1) + t * (14 * (-s1) - 11 + t * (4 - 5 * (-s1)))))) + 
      (t - 1) * ((-s2) + (t - 1) * (2 * (-s2) + (t - 1) * (t - 
          1) * (-2 * (-s2) - (t - 1) * (-s2)))))^2

Again, a simple loop can be used to generate tangent expressions for all X-spline scenarios.

It is easy to generate an expression for the unit tangent from the X-spline tangent expression by dividing by the square-root of the sum of squares of the x-component and y-component; the xsplineUnitTangentExpr function can do this. We can then repeat the process, differentiating the unit tangent to get an expression for the normal function and then generating an expression for the unit normal. This is provided by the xsplineUnitNormal function. However, these expressions become uncomfortably large. For the scenario demonstrated above, where the tangent expression consists of 65 lines of R code, the unit tangent expression is 193 lines of R code and the unit normal expression is 4,458 lines of R code. Consequently, the implementation in the 'vwline' package only uses the R expressions that were generated for X-spline tangents and then divides by the tangent lengths and performs a 90 degree rotation to obtain unit normals.

5. Rendering X-spline offset curves

An X-spline offset curve is a continuous mathematical function. Rendering such a curve is most easily achieved by drawing a series of many short straight line segments that give the appearance of a smooth curve. More segments are used for longer curves and for curves with sharper bends. For example, the image below shows an X-spline drawn in R with alternating straight line segments coloured black and white; the segments are shorter on the sharper bends.

The offset curve is drawn by evaluating the tangent at the end of each of these straight line segments, dividing by the length of the tangent, rotating 90 degrees to get the unit normal, then multiplying by the appropriate offset. The diagram below shows the calculation of the left offset for a fixed offset (the unit normals are green lines and the left offset curve is red).

Specifying the amount of offset

The offset curve equation consists of the original curve plus the unit normal function mulitplied by a function, $d(t)$, that specifies that amount of offset at any point on the original curve: $$\mathbf{r}_0(t)=\mathbf{r}(t) + d(t)\mathbf{n}(t)$$

The function $d(t)$ can in theory be any function that takes a single numeric argument in the range 0 to 1. In the 'vwline' package, the offset can be specified using the widthSpline function. This describes the offset as an X-spline, where the x-values for the control points are distances along the main curve and the y-values for the control points are amounts of offset. The following code gives an example where the offset rises smoothly from 0 to 1cm at half-way along the main curve, then back to 0 (the main curve is shown as a white line).

x <- 1:4/5
y <- c(.2, .7, .2, .7)
w <- widthSpline(unit(c(0, 1, 0), "cm"))
grid.offsetXspline(x, y, w, shape=-1)

Line ends

The offset curve only produces left and right offsets for the main curve, as shown below (the red lines).

If we want to produce a closed shape, we need to connect the ends of the offset curves. In the 'vwline' package, there are options for "butt", "square", "round", and "mitre" endings. The example below shows "round" ends.

Loops in the offset curve

The offset curve on the inside of a sharp bend in the main curve can form loops, as shown in the example below, where the left offset curve (the red line) forms a loop when the main curve takes a sharp left turn.

The 'vwline' package eliminates those loops by calling polysimplify from the 'polyclip' package (Johnson and Baddeley, 2017) on the final closed shape (after adding line ends). In the example below, the left image shows offset curves that contain loops and the right image shows the final result with round ends added and the loops eliminated.

Sharp corners and wide offsets

The offset curve is obtained by breaking the main curve into straight line segments, calculating offset points at the ends of those segments, and then joining up those offset points to create segments on the offset curve. This means that, if there is a sharp bend in the main curve and the offset is large, the segments in the offset curve may become long. This may result in a visibly non-smooth offset curve, as shown below.

This problem could be avoided if the offset curve function was evaluated directly, but this approach has not even been attempted because of the cost of evaluating the enormous offset curve expressions.

Byte compiling

The expressions for the X-spline tangent functions, even though they are "only" 65 lines of R code long, can be slow to evaluate. For this reason, the 'vwline' package has the 'ByteCompile' option set so that the package code is byte-compiled at installation to improve performance.

6. Summary

This document has described an algorithm for generating offset curves for X-splines, where the amount of offset can vary along the length of the X-spline. The algorithm relies on having expressions for the X-spline tangent function and this was generated computationally. The algorithm involves flattening the X-spline to a series of straight line segments, evaluating the tangent function for an X-spline at the end of each segment, calculating the unit normal at each point, and finally multiplying the unit normals by an offset. The algorithm is implemented as the grid.offsetXspline function in the 'vwline' package for R.

7. Technical requirements

8. Resources

How to cite this document

Murrell, P. (2017). Offset Curves for Variable-Width X-splines. Technical Report 2017-03, University of Auckland. [ bib ]

Offset Curves for Variable-Width X-splines

Table of Contents: