3
$\begingroup$

I'm trying to get consistent normals along a 3D Bezier curve $B(t)$, where for any point I compute the normal as:

$$ \begin{align} \vec{a} &= B'(t) \\ \vec{b} &= B''(t) \\ \vec{c} &= \vec{a} + \vec{b} \\ \vec{r} &= \vec{c} × \vec{a} \\ \vec{n} &= \vec{r} × \vec{a} \\ \end{align} $$

So, get the derivative at a point for time value $t$, and implicitly get the plane of curvature at the point by computing the cross product of the derivative vector at the point, and the "next" derivative vector we get from moving the derivative by the amount dictated by the second derivative. The cross product yields the axis of rotation, so to then form the normal at the point for time value $t$ I take the cross product of the axis of rotation, and the original derivative vector, since these three vectors are by definition perpendicular.

The problem is that normals computed this way are not consistent: they will "flip" around inflections, and I'm not sure what the right way is to go about making sure that does not happen.

As visual illustration, consider the following 3D cubic Bezier curve:

$$ B(t) = \left[\begin{matrix}1&t&t^2&t^3\end{matrix}\right] \left[\begin{matrix}1&0&0&0\\-3&3&0&0\\3&-6&3&0\\-1&3&-3&1\end{matrix}\right] \left[\begin{matrix} 0 & 0 & 0\\ -0.38 & 2.68 & 0\\ -0.25 & 5.41 & 0\\ -0.15 & 8.21 & 0 \end{matrix}\right] $$

Now, this happens to be a 3D curve that lies entirely on the x/y plane, but it illustrates the problem rather well. The above procedure yields the following normals:

3d normals based on cross products

However, this is rather different from the 2D normals we get when taking advantage of the 2D plane, where a normal can be constructed by simply rotating the (normalised) derivative vector a quarter turn clockwise, setting $(x,y)$ as $(-y,x)$:

2D normals based on vector rotation in the x/y plane

I'd like to get something similar to the 2D case for the 3D case, but I don't know how to ensure that the cross products are unaffected by "which direction" the second derivative moves the derivative across its plane of curvature

(Effectively, how do I ensure that, when considering the triplet {normal,derivative,axis of rotation} that these always map to the local {x,y,z} axes, rather than sometimes mapping to {x,y,z} and somethings mapping to {y,x,z} axes)

Edit

While more "algorithmic" than I'd like, the only workable solution I've found so far is to compute the normals for two points $B(t)$ and $B(t+\varepsilon )$, then computing the angular difference in the plane for those two normals,

$$ \theta = \textit{acos} \left ( \frac{n_1 \cdot n_2 }{||n_1|| ||n_2||} \right ) $$

and then check whether that value is close to $\pi$ or not. Even in fast-changing curves, the angle between two "reasonable" normals is a relatively small value, so if the angle suddenly flips to "nearly $\pi$" then as of that time value the "desired normals" are negative actual normal.

While that works, it feels kind of hacky.

Without algorithmic flipping:

cabinet projection of the real curve normals

With algorithmic flipping:

cabinet projection of the algorithmically adjusted curve normals

Note this does not affect cuves with "reasonable twisting", e.g. when we set the $z$ values to $\{0,200,-200,600\}$ for the first, second, third and fourth control point respectively:

cabinet projection of a twisty curve, where the algorithmic approach does not change anything

$\endgroup$
11
  • $\begingroup$ It would be better if you showed the “2-D normals” for the same curve. $\endgroup$
    – amd
    Jul 6, 2018 at 23:40
  • $\begingroup$ Do you have some concrete way to characterize the “right” direction for the normal in 3-D? Which ones in your first illustration are the “correct” ones? Does the answer change if you look at it from the other side? $\endgroup$
    – amd
    Jul 6, 2018 at 23:42
  • $\begingroup$ In that case you should orient and scale them so that they look the same. The curve in the first illustration bends left as you move up the page, but the second one bends to the right. $\endgroup$
    – amd
    Jul 6, 2018 at 23:53
  • $\begingroup$ @Mike'Pomax'Kamermans: In my answer, I used the exact same matrices as you have, just combined them into vectors (using Maple), so they should match yours exactly. Even if they do not, don't mind the differences, as the numerical values are not the problem here. $\endgroup$ Jul 7, 2018 at 2:10
  • 2
    $\begingroup$ See this paper microsoft.com/en-us/research/wp-content/uploads/2016/12/… $\endgroup$
    – Oppenede
    Jul 8, 2018 at 9:11

5 Answers 5

4
$\begingroup$

Thanks to @Oppenede commenting on the question post, it turns out that what I was looking for in this case is called the "Rotation minimizing frame" of a point, also known as the "parallel transport frame", or "bishop frame".

This is an algorithmic procedure, where you compute the orthogonal vector triplet $\{\textit{tangent}, \textit{rotation axis}, \textit{normal}\}$ for the point at time value 0, and then compute subsequent frames based on "the previous frame", using an ever so slightly modified version of the procedure explained in section 4 or "Computation of Rotation Minimizing Frames" (Wenping Wang, Bert Jüttler, Dayue Zheng, and Yang Liu, 2008):

ArrayList<VectorFrame> getRMF(int steps) {
  ArrayList<VectorFrame> frames = new ArrayList<VectorFrame>();
  double c1, c2, step = 1.0/steps, t0, t1;
  PointVector v1, v2, riL, tiL, riN, siN;
  VectorFrame x0, x1;

  // Start off with the standard tangent/axis/normal frame
  // associated with the curve just prior the Bezier interval.
  t0 = -step;
  frames.add(getFrenetFrame(t0));

  // start constructing RM frames
  for (; t0 < 1.0; t0 += step) {
    // start with the previous, known frame
    x0 = frames.get(frames.size() - 1);

    // get the next frame: we're going to throw away its axis and normal
    t1 = t0 + step;
    x1 = getFrenetFrame(t1);

    // First we reflect x0's tangent and axis onto x1, through
    // the plane of reflection at the point midway x0--x1
    v1 = x1.o.minus(x0.o);
    c1 = v1.dot(v1);
    riL = x0.r.minus(v1.scale( 2/c1 * v1.dot(x0.r) ));
    tiL = x0.t.minus(v1.scale( 2/c1 * v1.dot(x0.t) ));

    // Then we reflection a second time, over a plane at x1
    // so that the frame tangent is aligned with the curve tangent:
    v2 = x1.t.minus(tiL);
    c2 = v2.dot(v2);
    riN = riL.minus(v2.scale( 2/c2 * v2.dot(riL) ));
    siN = x1.t.cross(riN);
    x1.n = siN;
    x1.r = riN;

    // we record that frame, and move on
    frames.add(x1);
  }

  // and before we return, we throw away the very first frame,
  // because it lies outside the Bezier interval.
  frames.remove(0);

  return frames;
}

(This uses Java syntax but should be easy enough to parse for porting to any other language)

The result of this procedure leads to rather aesthetically pleasing normals. For the original planar curve, we see the following, with the RMF normals in green and the original normals in blue:

RMF normals for the planar curve, in green

And for the non-planar curve, again with RMF normals in green and original normals in blue:

RMF normals for the non-planar curve, in green

So this works really well. The downside of course is that this means normals can no longer be computed "on demand", as each frame relies on the previous frame, necessitating a full RMF computation pass and then interpolating for missing normals. But, based on the literature available, there does not appear to be a way to get nice, consistent looking normals without an iterative approach like this.

So Rotation Minimizing Frames it is!

$\endgroup$
0
1
$\begingroup$

We have $$\begin{aligned} B(t) =& \left [ \begin{matrix} -0.54 t^3 + 1.53 t^2 - 1.14 t \\ 0.02 t^3 + 0.15 t^2 + 8.04 t \\ 0 \end{matrix} \right ] \\ B^\prime(t) =& \left [ \begin{matrix} -1.62 t^2 + 3.06 t - 1.14 \\ 0.06 t^2 + 0.30 t + 8.04 \\ 0 \end{matrix} \right ] \\ B^{\prime\prime}(t) =& \left [ \begin{matrix} -3.24 t + 3.06 \\ 0.12 t + 0.30 \\ 0 \end{matrix} \right ] \\ \; & 0 \le t \le 1\end{aligned}$$ with $B(t)$ the point on the curve, $B^{\prime}(t)$ the direction (velocity or tangent), and $B^{\prime\prime}(t)$ is curvature (or acceleration), at $t$.

If the curve curves (changes direction, even infinitesimally) at $t$, then the direction vector rotates around vector $n(t)$ at that point: $$n(t) = B^{\prime}(t) \times B^{\prime\prime}(t)$$

If the curve is straight at $t$, then $B^{\prime\prime}(t) = 0$.

A cubic curve can be straight, but "accelerate" or "decelerate". (If you consider a line parametrized as a cubic curve, putting the two control points along the line segment between the starting and ending points, will not change the shape of the curve at all, only how the curve is formed; i.e., the function describing the ratio $\lVert B(t) - B(0)\rVert / \lVert B(1) - B(0)\rVert$.) This "acceleration" and "deceleration" occurs when $B^{\prime} \parallel B^{\prime\prime}$. Because at such points the curve does not change its direction (the direction of $B^{\prime}(t)$ does not change, only its magnitude changes), $n(t) = 0$.

Herein lies the problem.

In two dimensions, the axis of rotation (around which you rotate $B^{\prime}(t)$ to get the perpendicular vector) stays constant; it is the "implicit" third axis. (It is implicit in the way the 2D analog of a cross product is defined.)

In three dimensions, the axis vanishes at $t \approx 0.9398$. Because all $z$ components for the curve are zero, when the rotation axis exists, it is parallel to the $z$ axis. Before $t \approx 0.9398$, the axis is towards $+z$, after it is towards $-z$.

A simple answer is to use $n(t) = N$ (so that the normal vector OP is looking for is $B^{\prime}(t) \times N$), with $N$ defined as the normal to the plane the curve lies in. (A cubic curve has four control points. If these lie in the same plane, they curve lies in that plane too. $N$ is the normal to this plane.)

That yields the exact same normal vector for the 3D-extended 2D curve.


The question of what to do with curves that are not planar, and actually curve in all three dimensions, remains.

If you switch to the $n(t)$ mentioned above, the normal vectors thus calculated will behave similarly to OP's image, even when the curve is otherwise planar. You could require it is in the same halfspace as $n(0)$, but it would still leave the point-like discontinuity at $t$ where $B^{\prime}(t) \parallel B^{\prime\prime}(t)$.

The proper answer to this remaining part is that it depends on what these normals are used for.

In a practical sense, when you extend a line or curve from 2D to 3D, you should get a surface to get analogous properties. If you keep the objects dimensionality unchanged, you get an additional degree of freedom. (If you extended the 2D cubic curve into 3D cubic patch, you'd have an analogous surface normal, too.)

This additional degree of freedom means that while a line in 2D has a normal vector, a line in 3D has a normal plane (as defined by its normal, $B^{\prime}(t)$).

Similarly, while a curve in 2D has a normal vector (the tangent, $B^{\prime}(t)$, rotated 90° clockwise or counterclockwise), a curve in 3D has a normal plane (as defined by $B^{\prime}(t)$).

$\endgroup$
10
  • $\begingroup$ It feels like there should be an approach here, somewhere, based on the fact that any non-planar curve can be described as a planar curve as long as we restrict the interval we're examining (in the same way that any curve can locally be considered equivalent to its tangent). While it is possible to catch the B''(t) transition from +z to -z, and flip the axis or rotation, that technically leaves a discontinuity at the point where B''(t)=0. I'm now trying to see if I can make things at least "work" by using the placement of the circle implied by the curve's radius of curvature. $\endgroup$ Jul 7, 2018 at 17:05
  • $\begingroup$ @Mike'Pomax'Kamermans: The issue is that a straight line, or a straight segment of a curve, has no normal vector; only a normal plane. If you extend from a 2D curve to a 3D curve, you need to extend the concept of a normal vector to a normal plane, to get analogous properties; that's what the extra degree of freedom does. To me, it feels like a desperate attempt of shoehorning a concept into a situation it does not fit in. Maybe you could expand in your question a bit about the purpose of these normal vectors? $\endgroup$ Jul 8, 2018 at 2:02
  • $\begingroup$ That is indeed the problem. Another is the nomenclature: the vector I'm looking for is not "the" normal vector as you point out, it's really more an aesthetically useful vector in the normal plane, for things like graphical design and visual artistry work - for instance, offsetting along these "normals" of a curve and getting something that yields the kind of contour designers and artists can actually make use of, rather than offsets that kink and flip wildly (as in this example). $\endgroup$ Jul 8, 2018 at 5:00
  • 1
    $\begingroup$ @Mike'Pomax'Kamermans: All I use are the Wikipedia articles on versors, axis-angle representation, Quaternions and spatial rotation, and Quaternion (Hamilton product). The key is that you can normalize a quaternion $(w, x, y, z)$ by dividing each component by $\sqrt{w^2 + x^2 + y^2 + z^2}$ after operations, [...] $\endgroup$ Jul 8, 2018 at 17:47
  • 1
    $\begingroup$ [...] and blend between two quaternions as you would scalars; just normalize the result. An unit quaternion or versor describes an orientation. You can think of it as describing an axis, and an angle around that axis. If you get your mind around those concepts, that's about it: quaternions are really much easier to use than for example Euler angles! $\endgroup$ Jul 8, 2018 at 17:49
0
$\begingroup$

The normal vector of a parametric curve $B(t)$ can be found in the following way:

1) compute the first derivative $B'(t)$ and the 2nd derivative $B''(t)$,
2) compute the unit tangent vector $t=B'(t)/|B'(t)|$,
3) compute the bi-normal vector $$b= \frac{B'(t) \times B''(t)}{\left|B'(t) \times B''(t)\right|}$$
4) compute the normal vector as $n= b \times t$.

The 3 unit vectors $t$, $b$ and $n$ are mutually orthogonal and will form the so-called "Frenet-Serret frame".

$\endgroup$
10
  • $\begingroup$ Those are the normals I already have, and are super unreliable for graphics work (the first image shows them based on unnormalized vectors, but their alignment is identical). I need "normals" that lie in the normal plane, but oriented to fit a regular viewer's idea of "sticking out from the curve in a relatively uniform manner", so in the stated case, 2D and 3D showing the same result. If those can't be called normals, then that's fine: I don't care about nomenclature in this case, only that the actual vectors look right. $\endgroup$ Jul 7, 2018 at 21:17
  • $\begingroup$ For points lying on a plane, it may be possible to check the orientation of the binormal by examining the signs of elements of the curvature vector given by $\boldsymbol{\kappa} = (B'(t) \times B''(t))/|B'(t)|^3$. It may also to useful to examine the torsion for more general situations. $\endgroup$ Jul 7, 2018 at 22:22
  • $\begingroup$ A possibly better way is to determine whether the normal needs to be flipped is make sure that the Frenet-Serret frame is always right-handed, i.e., check whether $(t \times b) \cdot n > 0$. $\endgroup$ Jul 7, 2018 at 22:39
  • $\begingroup$ Upon further thought, the frame will always have the same orientation. So, the only robust way to check the flip is to take two adjacent points on the curve, check $t_1 \cdot t_2 > 0$ and $b_1 \cdot b_2 > 0$ and flip when the sign is negative. $\endgroup$ Jul 7, 2018 at 22:54
  • $\begingroup$ That's what I initially thought, too, but that's still not enough: both dot products are positive on both sides of the crossing point, only having a flipped sign where we actually cross from positive to negative second derivative. $\endgroup$ Jul 8, 2018 at 5:15
0
$\begingroup$

(I came across this question trying to debug my implementation of essentially the same problem. The stated answer did not work for my case. It could be that I implemented it wrong, or maybe the stated answer works for less curvy curves than I need [see below]. I did manage to fix my own initial implementation, and it is conceptually simpler than the procedure in the accepted answer, so I thought it would be a good additional answer, regardless of whether the accepted answer is sometimes inferior.)

Given a curve approximated by a sequence of points and tangents, find a sequence of normals that follows the curve and follows some intuition of non-drifting. (Actually, we need a complete "frame", which consists of a sequence of three mutually orthogonal vectors, made from the tangents, the "consistent" normals, and their cross product.)

Let $t_i$ be the sequence of tangents.

We compute the three axes $x_0$, $y_0$, and $z_0|$ of the first frame using Method 1 described in this answer. It is simply a way to find a good set of orthogonal vectors in a plane, and we use the plane at the first point of our curve with the same normal as our first tangent.

We then compute each frame iteratively as follows:

$z_i = t_i$

$y_i = t_i \times x_{i - 1}$

$x_i = y_i \times z_i$

The idea is basically we start with some normal. Then the next normal is a vector that 1) is perpendicular to the tangent, and 2) lies in the same plane as the tangent and previous normal. The cross products are just calculating such a vector.

I use this frame to calculate a "thin tube" around a curve (represented by a sequence of points and tangents).

Below on the left, you can see the result of the implementation in the accepted answer (and remember my implementation may be wrong rather than the answer itself), and on the right is the implementation described here. Notice how the seem (the black line) drifts from the top to the bottom in the left mesh, but stays on the bottom in the right mesh.

enter image description here

$\endgroup$
1
  • $\begingroup$ Missed this answers - let me implement this as well and see what it does for some of the weirder curves, like 3D curves with one dimension set to 0. $\endgroup$ Feb 19, 2022 at 18:56
0
$\begingroup$

I compared this solution from Herman Tulleken with RMF approximation from Wang et al 2008 and for my examples, it seems indeed to work better. But one should consider a special case where the angle between $t_i$ and $x_{i-1}$ is zero (sudden curve turn), so the cross-product will return a zero length vector. I guess in this case one can just switch vectors (keeping right-handed orientation): $$z_i = t_i$$ $$x_i=y_{i−1}×t_i$$ $$y_i=z_i×x_i$$

$\endgroup$
2
  • 1
    $\begingroup$ That should probably be a comment on that answer, not an actual answer post? (which I realize is hard to do when you don't have enough rep for comments yet, but using answers to post comments is considered "circumventing" the rules) $\endgroup$ Mar 22, 2023 at 15:02
  • $\begingroup$ I'm sorry, indeed, my reputation points are not enough. $\endgroup$ Mar 24, 2023 at 11:02

You must log in to answer this question.

Not the answer you're looking for? Browse other questions tagged .