Numerically computing normal, binormal, and tangent directions of non-parametric curve in $\mathbb{R}^{3}$

Question

Let's say that I have numerical data for a curve in $\mathbb{R}^{3}$, but I do not have the parametric equations of the curve; all that I have is a sampling of $N$-many points that lie on this curve, $R_{i}=[x_{i} \, y_{i} \, z_{i}]$ for rows $1\leq i \leq N$. How would I numerically determine the normal, binormal, and tangent vectors for each point along this curve?

I do know the direction that points on the curve travel, and so I can numerically compute/approximate the tangent vector using, say, a central difference scheme:

$$\vec{T}v_{i} \approx \frac{v_{i+1}-v_{i-1}}{2\Delta t},$$

where $\Delta t$ is the time delay in the sampling of this parametric curve (and $i$ runs $2 \leq N-1$).

Conventionally, the binormal and directions require the parametric equation(s) of motion in order to be computed (e.g., the Frenet-Serret Equations). Are there methods for numerically approximating the normal $\vec{N}$ and binormal $\vec{B}$, just as I have numerically approximated the tangent vector $\vec{T}$ at $v_{i}$?

The key, I think, is the normal vector because once $\vec{N}v_{i}$ and $\vec{T}v_{i}$ have been found, $\vec{B}v_{i}$ can then be computed by its definition as $\vec{B} = \vec{T} \times \vec{N}$. In particular, $\vec{N}$ should be orthogonal to $\vec{T}$ and $\vec{B}$, and it should point in the direction that the curve bends. This refers to the direction in which the curve deviates from straight line motion. The binormal vector points in the direction around which the tangent vector turns.

Answer 1

This question is at the heart of the general field of "Discrete Differential Geometry". It's also the kind of thing my colleague Tom Banchoff has studied for more than 50 years.

The first question you have to ask is "Does the underlying curve have a Frenet frame?" If, for instance, with finer and finer sampling you get a sequence of tangent vectors that don't converge, you've got a problem. (A good example of this: look at the graph of $y = |x|$ near the origin.) Even if there is a limiting tangent vector, it's possible that the curvature ends up zero, in which case there's no Frenet frame.

The second question (which maybe should have been the first) is "is there an underlying curve at all?" If I draw, say, 6 random points from a gaussian distribution around the origin, you can "connect the dots" to form a curve, but the next point I draw from that distribution is unlikely to lie on it. Indeed, the underlying set isn't a curve at all -- it's all of $\Bbb R^2$. That's an extreme case, but even things as simple as the solutions of polynomial equations like $x^2 + y^2 = 1$ might have problems, as $y^3 = x^2$ will show you.

Having asked those questions, which I beg you not to ignore, you can pretend everything is fine and then there's a nice solution: to compute the normal at point $i$, look at the points $P_{i-1}, P_i, P_{i+1}$; these three lie on a plane, and the curve normal is the vector in that plane that's closest to being orthogonal to the tangent. Letting $v = P_i - P_{i-1}$ and $u = P_{i+1} - P_i$m the normal to the plane (not the curve normal!) is $n = u \times v$, and you can then approximate $N$ (the curve normal) as $N = \pm n \times T$ (I leave you to work out the sign). By the way, the vector $\pm n$ (again, you get to work out the sign) is a good approximation of the binormal (once you make it a unit vector).

To be honest, the approach in the previous paragraph will work OK, but near a point of inflection, the three points will be collinear (or nearly so), and then you're screwed. Then again, at such a point, you don't actually have a Frenet frame either, so maybe you're lucky that the computation will tell you that. More likely, you'll just get a really short vector for $n$, and thus a short vector for $N$, you'll normalize it without thinking, get a vector in some random direction due to numerical glitches and then claim that my reasoning was wrong. (Sigh.)

A slightly improved approach is to estimate the best-fit-plane through a sequence of $5$ or $7$ or more points, i.e., $P_{i-k},\ldots, P_{i-1}, P_i, P_{i+1}, \ldots, P_{i+k}$; if your sampling is fine enough, this'll give a great approximation of $n$, and hence $N$, at least when $N$ is well-defined. If your sampling is too coarse, then a large $k$ will lead to garbage. The problem is that "how big should $k$ be" depends on your sampling rate, and on the unknown curvature of the underlying curve. I have no advice on how to address this particular tradeoff.

Why does this work? Well, it helps to go back to classical differential geometry, where the "osculating plane" (the "TN plane" in Frenet terms!) is defined by some weird phrase like "it's the plane passing through three successive points of the curve". Once you understand that kind of phrase, translating to numerical computation isn't too tough.

As another (and better in my opinion) approach, you could regard the sequence of points you've got as a polygon, and ask what the frenet frame of a polygon actually "should be". Once you think about it, you find things like "the tangent vector along an edge should be that edge", and "the binormal vector at a vertex should be the cross product of the two adjacent edge vectors." What about the binormal vector along the edge? Well, there's a pretty good case for extending from each end towards the edge-midpoint as a constant, so that the binormal becomes a locally-constant function (rather than a continuous one), and then filling in, at the edge-midpoint, with an arc of vectors that transitions from one binormal to another...but that's too tough to describe here. This idea belongs (as far as I know) to Tom Banchoff, who told it to me over a cup of coffee one day. I'm pretty sure that Eitan Grinspun has used this idea (and related stuff for surfaces) in some of his computer graphics work as well, probably having discovered it independently. You might want to take a look at that.