Convert hyperbola to Bézier curve for graphing orbital paths

Question

I am writing a 2D simulator and game using the HTML canvas which involves orbital mechanics. One feature of the program is to take the position and velocity vector of a satellite at one point and return the semi-major axis, eccentricity, argument of periapsis, etc of a 2D orbit around one planet. When the eccentricity is less than one, I can easily graph the orbit as an ellipse using ctx.ellipse(). However, for eccentricities greater than one, the correct shape of the orbit is a hyperbola. At the moment, my program just draws nothing if the eccentricity is greater than one, but I would like it to graph the correct hyperbolic orbit. Since there is no built in "hyperbola" function, I need to convert my orbit into a Bézier curve. I am at a bit of a loss as to how to do this. The inputs would be the location of one focus, semi-major axis, eccentricity, and argument of periapsis (basically how far the orbit is rotated) and it should return the correct control points to graph a Bézier curve approximation of a hyperbola. It does not have to be exactly perfect, as long as it is a close enough fit. How can I approach this problem?

Answer 1

In terms of conic sections, hyperbola are unfortunately the one class of curves that the Canvas cannot natively render, so you're stuck with approximating the curve you need. There are some options here:

1. Flatten your curve, by sampling the hyperbola at one or two points in the distance and lots of points near the extrema so that you can draw a simple polygon that looks like a curve.
2. Model the hyperbola with a single "best approximation" quadratic or cubic curve.
3. As @fang mentions: sample the curve at a few points and convert the Catmull-Rom spline through those points to Bezier form.
4. Combine approaches 1 and 2. using a single Bezier curve to approximate the part of the hyperbola that actually looks curved, and using straight lines for the parts that don't.
5. Combine approaches 1 and 3, using a Catmull-Rom spline for the curvy bit, and straight lines for the straight bits.

1: Curve flattening

Curve flattening is basically trivial. Rotate your curve until it's axis-aligned, and then just compute y given x using the standard hyperbolic function, where a is half the distance between the extrema, and b is the semi-minor axis:

x²/a² - y²/b² = 1
x²/a² = 1 + y²/b² 
x²/a² - 1 = y²/b² 
b²(x²/a² - 1) = y²
b²(x²/a² - 1) = y²
± sqrt(b²(x²/a² - 1)) = y

plug in your values, iterate over x to get a sequence of (x,y) coordinates (remembering to generate more coordinates near the extrema), then turn those into a moveTo() for the first coordinate, followed by however many lineTo() calls you need for the rest. As long as your point density is high enough for the scale you're presenting, this should look fine:

function flattenHyperbola(a, b, inf=1000) {
  const points = [],
        a2 = a**2,
        b2 = b**2;

  let x, y, x2;

  for (x=inf; x>0.1; x/=2) {
    x2 = (a+x)**2;
    y = -Math.sqrt(b2*x2/a2 - b2);
    points.push({x: a+x, y});
  }

  points.push({x:a, y:0});

  for (x=0.1; x<inf; x*=2) {
    x2 = (a+x)*(a+x);
    y = Math.sqrt(b2*x2/a2 - b2);
    points.push({x:  a+x, y});
  }

  return points;
}

Let's plot the hyperbola in red and the approximation in blue:

Of course the downside with this approach is that you will need to create a separate flattened curve for every scale the user may view your graphics at. Or, you need to generate a flattened curve with lots and lots of points, and then draw it by skipping over coordinates depending on how zoomed in/out things are.

2: Bezier approximation

The parametric representation of a hyperbola is f(t)=(a*sec(t), b*tan(t)) (or rather, that's the representation for the y-axis aligned hyperbola - we can get any other variant by applying the standard rotation transform). We can have a quick look at the Taylor series for these functions to see which order of Bezier curve we can use:

sec(t) = 1 + t²/2 + 5t⁴/15 + ...
tan(t) = x + t³/3 + 2t⁵/15 + ...

So we might be able to get away with just the first two terms for each dimension in which case we can use a cubic Bezier (as the highest order is t³):

Turns out, that won't do: It's just way too inaccurate, so we're going to have to better approximate: we create a Bezier curve with start and end points "well off into the distance", with the control points set such that the Bezier midpoint coincides with the hyperbola's extrema. If we try this, we might be fooled into thinking that'll work:

But if we pick x far enough away, we see this approximation quickly stops working:

function touchingParabolicHyperbola(a, b, inf=1000) {
  const beziers = [],
        a2 = a**2,
        b2 = b**2;

  let x, x2, y, A, CA;

  for(x=50; x<inf; x+=50) {
    x2 = x**2;
    y = sqrt(b2*x2/a2 - b2);

    // Hit up https://pomax.github.io/bezierinfo/#abc
    // and model the hyperbola in the cubic graphic to
    // understand why the next, very simple-looking,
    // line actually works:
    A = a - (x-a)/3;

    // We want the control points for this A to lie on
    // the asymptote, but for small x we want it to be 0,
    // otherwise the curve won't run parallel to the
    // hyperbola at the start and end points.
    CA = lerp(0, A*b/a, x/inf);

    beziers.push([
      {x,    y: -y}, 
      {x: A, y:-CA}, 
      {x: A, y: CA}, 
      {x,    y}, 
    ]);
  }

  return beziers;
}

This shows us a sequence of curves that starts off looking decent but becomes completely useless pretty fast:

One obvious problem is that the curves end up going past the asymptotes. We can fix that by forcing the control points to (0,0) so that the Bezier hull is a triangle and the curve will always lie inside of that.

function tangentialParabolicHyperbola(a, b, inf=1000) {
  const beziers = [],
        a2 = a**2,
        b2 = b**2;

  let x, x2, y;

  for(x=50; x<inf; x+=50) {
    x2 = x**2;
    y = sqrt(b2*x2/a2 - b2);  
    beziers.push([
      {x, y:-y}, 
      {x: 0, y:0}, 
      {x: 0, y:0}, 
      {x, y}, 
    ]);
  }

  return beziers;
}

This leads to a series of curves that go from useless on one side, to useless on the other side:

So single curve approximations are not all that great. What if we use more curves?

3: Poly-Bezier using a Catmull-Rom spline

We can overcome the above problem by using multiple Bezier curves along the hyperbola, which we can (almost trivially) compute by picking a few coordinates on the hyperbola, and then constructing a Catmull-Rom spline through those points. Since a Catmull-Rom spline through N points is equivalent to a poly-Bezier made up of N-3 segments, this could be the winning strategy.

function hyperbolaToPolyBezier(a, b, inf=1000) {
  const points = [],
        a2 = a**2,
        b2 = b**2,
        step = inf/10;

  let x, y, x2,
   
  for (x=a+inf; x>a; x-=step) {
    x2 = x**2;
    y = -Math.sqrt(b2*x2/a2 - b2);
    points.push({x, y});
  }

  for (x=a; x<a+inf; x+=step) {
    x2 = x**2;
    y = Math.sqrt(b2*x2/a2 - b2);
    points.push({x, y});
  }

  return crToBezier(points);
}

With the conversion function being:

function crToBezier(points) {
  const beziers = [];

  for(let i=0; i<points.length-3; i++) {
    //  NOTE THE i++ HERE! We're performing a sliding window conversion.
    let [p1, p2, p3, p4] = points.slice(i);
    beziers.push({
      start: p2,
      end: p3,
      c1: { x: p2.x + (p3.x-p1.x)/6, y: p2.y + (p3.y-p1.y)/6 },
      c2: { x: p3.x - (p4.x-p2.x)/6, y: p3.y - (p4.y-p2.y)/6 }
    })
  }

  return beziers;
}

Let's plot that:

We have to do a bit more work up front compared to flattening, but the upside is that we now have a curve that actually "looks like a curve" at any scale.

4: Combining (1) and (2)

Now, most of the hyperbola actually "looks straight", so using lots of Bezier curves for those parts does feel a bit silly: why not only model the curvy bit with a curve, and model the straight bits with straight lines?

We already saw that if we fix the control point to (0,0), there might be a curve that's at least decent enough, so let's combine approaches 1 and 2, where we can create a single Bezier curve with start and end points "close enough" to the curve, and tacking two line segments onto the ends that join up the bezier curves to two distant points on the asymptotes (which are at y=±b/a * x, so any large value for x will yield a usable-enough y)

Of course the trick is to find the distance at which the single curve still captures the curvature, while also making our lines to infinity look like they smoothly join up to our single curve. The Bezier projection identity comes in handy again: we want A to be at (0,0) and we want the Bezier midpoint to be at (a,0), which means our start and end points should have an x coordinate of 4a:

function hyperbolicallyFitParabolica(a, b, inf=1000) {
  const a2 = a**2,
        b2 = b**2,
        x = 4*a,
        x2 = x**2,
        y = sqrt(b2*x2/a2 - b2)
        bezier = [
          {x: x, y:-y}, 
          {x: 0, y: 0}, 
          {x: 0, y: 0}, 
          {x: x, y: y}, 
        ],
        start = { x1:x, y1:-y, x2:inf, y2: -inf * b/a},
        end   = { x1:x, y1: y, x2:inf, y2:  inf * b/a};

  return [start, bezier, end];
}

Which gives us the following result (Bezier in blue, line segments in black):

So that's not great, but it's not terrible either. It's certainly good enough if the audience doesn't scrutinize the render, and it's definitely cheap, but we can do quite a bit better with just a little more work, so: let's also look at the best approximation we can probably come up with here:

5: Combining (1) and (3)

If a single Bezier didn't work, and we already saw that using a Catmull-Rom spline instead of a single curve works way better, then we can of course also just combine approaches 1 and 3. We can form a much better fit around the extrema by constructing two Bezier curves rather than one, by generating five points centered on the extrema and converting the resulting Catmull-Rom spline through those points to Bezier form:

function probablyTheBestHyperbola(a, b, inf=1000) {
  let curve = [],
      a2 = a**2,
      b2 = b**2,
      x, y, x2,
      cover = 100;

  // generate two points approaching the midpoint
  for (x=a+cover; x>a; x-=cover/2) {
    x2 = x**2;
    y = -Math.sqrt(b2*x2/a2 - b2);
    curve.add(new Vec2(x, y));
  }

  // generate three points departing at the midpoint
  for (x=a; x<=a+cover; x+=cover/2) {
    x2 = x*x;
    y = sqrt(b2*x2/a2 - b2);
    curve.add(new Vec2(x, y));
  }

  const beziers = crToBezier(curve),
        start = {
          x1: points.get(1).x, y1: points.get(1).y,
          x2: inf, y2: -inf * b/a
        },
        end = {
          x1: points.get(3).x, y1: points.get(3).y,
          x2: inf, y2: inf * b/a
        };

  return { start, beziers, end };
}

Which gives us the following result (CR in blue, line segments in black):

And that's probably the best we're going to get in a trade-off between "cheap to compute", "easy to scale", and "looks rights".