3

I need to somehow compute the distance between a point and an Ellipse. I describe the Ellipse in my program as coordinates x = a cos phi and y = b sin phi (where a,b are constants and phi the changing angle).

I want to compute the shortest distance between a point P and my ellipse. My thought were to calculate the vector from the center of my ellipse and the point P and then find the vector that start from the center and reaches the end of the ellipse in the direction of the point P and at the end subtract both vectors to have the distance (thi may not give the shortest distance but it's still fine for what I need. The problem is I don't know how to compute the second vector. Does someone has a better Idea or can tell me how I can find the second vetor?

Thanks in advance!


EDIT1:

ISSUE:COMPUTED ANGLE DOESN'T SEEM TO GIVE RIGHT POINT ON ELLIPSE

Following the suggestion of MARTIN R, I get this result:

enter image description here

The white part is created by the program of how he calculates the distance. I compute the angle phi using the vector from the center P (of ellipse) to the center of the body. But as I use the angle in the equation of my ellipse to get the point that should stay on the ellipse BUT also having same direction of first calculated vector (if we consider that point as a vector) it actually gives the "delayed" vector shown above.

What could be the problem? I cannot really understand this behavior (could it have something to do with atan2??)

EDIT2: I show also that in the other half of the ellipse it gives this result:

enter image description here

So we can see that the only case where this works is when we have phi = -+pi/2 and phi = -+pi


IMPLEMENTATION FAILED

I tried using the implementation of MARTIN R but I still get the things wrong.

At first I thought it could be the center (that is not always the same) and I changed the implementation this way:

func pointOnEllipse(ellipse: Ellipse, p: CGPoint) -> CGPoint {

    let maxIterations = 10
    let eps = CGFloat(0.1/max(ellipse.a, ellipse.b))

    // Intersection of straight line from origin to p with ellipse
    // as the first approximation:
    var phi = atan2(ellipse.a*p.y, ellipse.b*p.x)

    // Newton iteration to find solution of
    // f(θ) := (a^2 − b^2) cos(phi) sin(phi) − x a sin(phi) + y b cos(phi) = 0:
    for _ in 0..<maxIterations {
        // function value and derivative at phi:
        let (c, s) = (cos(phi), sin(phi))
        let f = (ellipse.a*ellipse.a - ellipse.b*ellipse.b)*c*s - p.x*ellipse.a*s + p.y*ellipse.b*c - ellipse.center.x*ellipse.a*s + ellipse.center.y*ellipse.b*c
        //for the second derivative
        let f1 = (ellipse.a*ellipse.a - ellipse.b*ellipse.b)*(c*c - s*s) - p.x*ellipse.a*c - p.y*ellipse.b*s - ellipse.center.x*ellipse.a*c - ellipse.center.y*ellipse.b*s

        let delta = f/f1
        phi = phi - delta
        if abs(delta) < eps { break }
    }

  return CGPoint(x: (ellipse.a * cos(phi)) + ellipse.center.x, y: (ellipse.b * sin(phi)) + ellipse.center.y)

}

We can see what happens here:

enter image description here

This is pretty strange, all points stay in that "quadrant". But I also noticed when I move the green box far far away from the ellipse it seems to get the right vector for the distance.

What could it be?


END RESULT

Using updated version of MARTIN R (with 3 iterations)

enter image description here

DevX10
  • 483
  • 2
  • 6
  • 16
  • Did you see http://stackoverflow.com/questions/22959698/distance-from-given-point-to-given-ellipse ? – Sulthan Oct 23 '16 at 20:51

4 Answers4

3

x = a cos(phi), y = b sin (phi) is an ellipse with the center at the origin, and the approach described in your question can be realized like this:

// Point on ellipse in the direction of `p`:
let phi = atan2(a*p.y, b*p.x)
let p2 = CGPoint(x: a * cos(phi), y: b * sin(phi))

// Vector from `p2` to `p`:
let v = CGVector(dx: p.x - p2.x, dy: p.y - p2.y)

// Length of `v`:
let distance = hypot(v.dx, v.dy)

You are right that this does not give the shortest distance of the point to the ellipse. That would require to solve 4th degree polynomial equations, see for example distance from given point to given ellipse or Calculating Distance of a Point from an Ellipse Border.

Here is a possible implementation of the algorithm described in http://wwwf.imperial.ac.uk/~rn/distance2ellipse.pdf:

// From http://wwwf.imperial.ac.uk/~rn/distance2ellipse.pdf .

func pointOnEllipse(center: CGPoint, a: CGFloat, b: CGFloat, closestTo p: CGPoint) -> CGPoint {

    let maxIterations = 10
    let eps = CGFloat(0.1/max(a, b))

    let p1 = CGPoint(x: p.x - center.x, y: p.y - center.y)

    // Intersection of straight line from origin to p with ellipse
    // as the first approximation:
    var phi = atan2(a * p1.y, b * p1.x)

    // Newton iteration to find solution of
    // f(θ) := (a^2 − b^2) cos(phi) sin(phi) − x a sin(phi) + y b cos(phi) = 0:
    for i in 0..<maxIterations {
        // function value and derivative at phi:
        let (c, s) = (cos(phi), sin(phi))
        let f = (a*a - b*b)*c*s - p1.x*a*s + p1.y*b*c
        let f1 = (a*a - b*b)*(c*c - s*s) - p1.x*a*c - p1.y*b*s

        let delta = f/f1
        phi = phi - delta
        print(i)
        if abs(delta) < eps { break }
    }

    return CGPoint(x: center.x + a * cos(phi), y: center.y + b * sin(phi))
}

You may have to adjust the maximum iterations and epsilon according to your needs, but those values worked well for me. For points outside of the ellipse, at most 3 iterations were required to find a good approximation of the solution.

Using that you would calculate the distance as

let p2 = pointOnEllipse(a: a, b: b, closestTo: p)
let v = CGVector(dx: p.x - p2.x, dy: p.y - p2.y)
let distance = hypot(v.dx, v.dy)
Community
  • 1
  • 1
Martin R
  • 529,903
  • 94
  • 1,240
  • 1,382
  • I cannot believe that I tried everything and also a computation that was similar to yours but it was just so obvious.. thanks man, i will try it now and let you know if it worked :) – DevX10 Oct 22 '16 at 19:01
  • As I thought it isn't so easy. I'm pretty sure this would work perfectly with a circle, but in this case the angle `phi` doesn't actually give the point I need. is because there isn't a fixed radius but two different – DevX10 Oct 22 '16 at 20:41
  • I think it is because there isn't a fixed radius but the two constants that change the thing.. I show you what happens in my edited question – DevX10 Oct 22 '16 at 20:42
  • Thanks for this update, I tried and updated my question – DevX10 Oct 23 '16 at 19:20
  • @Ergo: The case that the ellipse center is not the origin is not handled correctly in your code. See updated answer. – Martin R Oct 23 '16 at 20:26
  • So you just use the same problem but just take point p relative to the center of the ellipse (point p1). It's awesome it finally starts working. Thanks a lot – DevX10 Oct 23 '16 at 21:23
1

Create new coordinate system, which transforms ellipse into circle https://math.stackexchange.com/questions/79842/is-an-ellipse-a-circle-transformed-by-a-simple-formula, then find distance of point to circle, and convert distance

Community
  • 1
  • 1
Alex Aparin
  • 4,393
  • 5
  • 25
  • 51
  • it's interesting could you explain better how you would do it? – DevX10 Oct 22 '16 at 21:54
  • @DevX10 This is _by far_ the simplest solution; confirmed to work (my use case) on an HTML5 Canvas. The linked question/answer is important: it shows that an ellipse is indeed "just" a circle plus affine transformation of the coordinate system. Instead of computing the distance to the ellipse, you reduce the problem to computing the distance to the circle. You need the transformation matrix which transforms a circle into your ellipse (in case of HTML5 Canvas, it will give it to you after you translate and skew to get ellipse from a circle). Then get inverse matrix transform... (TBC) – rawpower Jan 03 '22 at 11:08
  • ... and then transform the point you're computing the distance from into the coordinate system of the circle, then compute the distance. Finally, convert the distance into your normal coordinate system using the same transform which converts circle to the ellipse. You can liberally pick the circle, it doesn't matter what is radius is, as long as it can be transformed to the ellipse. I chose smaller of the two "radiuses" of the ellipse as the radius of the circle. Everything works well and fast, just a few matrix multiplications (performed by the Canvas for me). – rawpower Jan 03 '22 at 11:12
  • Finally, inverting a matrix is not complex, but most graphics systems - including the HTML5 Canvas - will likely do it for you, so you don't have to fret about that either. Pedantic: `DOMMatrix` which is returned by HTML5 Canvas predates the Canvas itself and can be used independently. – rawpower Jan 03 '22 at 11:15
1

I wrote up an explanation using Latex so it could be more readable and just took some screen shots. The approach I am sharing is one using a Newton step based optimization approach to the problem.

Note that for situations where you have an ellipse with a smaller ratio between the major and minor axis lengths, you only need a couple iterations, at most, to get pretty good accuracy. For smaller ratios, you could even probably get away with just the initial guess's result, which is essentially what Martin R shows. But if your ellipses can be any shape, you may want to add in some code to improve the approximation.

enter image description here

enter image description here

enter image description here

enter image description here

spektr
  • 777
  • 5
  • 14
  • This is really interesting, I think I will try this. It seems really nice but what would you suggest me to do in case the ratio of the ellipses I use is almost never bigger than 4 (considering a/b). Would you still suggest me to use newtons iteration or go with something simplier and less computation effort – DevX10 Oct 23 '16 at 09:40
  • @Ergo Based on the algorithm above, I would recommend 1 or 2 Newton iterations, depending on the level of accuracy you want. If you are okay with some error, not using any Newton steps could be fine. I also notice Martin R updated his answer with some Swift code to do the Newton iteration refinement, so looks like you don't need to worry much about the implementation. – spektr Oct 23 '16 at 14:30
0

You have the Ellipsis center of (a, b) and an arbitrary point of P(Px, Py). The equation of the line defined by these two points looks like this:

(Y - Py) / (b - Py) = (X - Px) / (a - Px)

The other form you have is an ellipse. You need to find out which are the (X, Y) points which are both on the ellipse and on the line between the center and the point. There will be two such points and you need to calculate both their distance from P and choose the smaller distance.

Lajos Arpad
  • 64,414
  • 37
  • 100
  • 175