Skip to content

Finding The Intersection Of Two Arcs That Lie On A Sphere

Published On:
Feb 15, 2018
Last Updated:
Nov 10, 2022

Finding The Intersection Of Two Arcs That Lie On A Sphere

This is a common problem which usually arises because you want to find if two great circles (or segments of two great circles) on the earth’s surface intersect (obviously with are approximating the earth as a sphere).

Let’s Start With Our Inputs To The Problem

We will choose to define each of the two arcs with two points on the surface of the sphere/earth. We could use either the geocentric (also called ECEF, coordinates are provided in terms of x, y, z) or geodetic (latitude/longitude) coordinate system to define these points. Let’s use geodetic for now:

Points for arc 1 (a1a1):

Pa11=[λ11ϕ11]Pa12=[λ12ϕ12]\begin{align} \b{P_{a11}} &= \left[ {\begin{array}{c} \lambda_{11} \\ \phi_{11} \end{array} } \right] \\ \b{P_{a12}} &= \left[ {\begin{array}{c} \lambda_{12} \\ \phi_{12} \end{array} } \right] \\ \end{align}

Points for arc 2 (a2a2):

Pa21=[λ21ϕ21]Pa22=[λ22ϕ22]\begin{align} \b{P_{a21}} &= \left[ {\begin{array}{c} \lambda_{21} \\ \phi_{21} \end{array} } \right] \\ \b{P_{a22}} &= \left[ {\begin{array}{c} \lambda_{22} \\ \phi_{22} \end{array} } \right] \\ \end{align}

where:
λ\lambda = latitude, in degrees
ϕ\phi = longitude, in degrees

We will then convert these into spherical coordinates with the following equations:

x=cos(λ)sin(ϕ)y=cos(λ)cos(ϕ)z=sin(λ)\begin{align} x = \cos(\lambda) \, \sin(\phi) \\ y = \cos(\lambda) \, \cos(\phi) \\ z = \sin(\lambda) \\ \end{align} Pa11=[x11y11z11]Pa12=[x12y12z12]Pa21=[x21y21z21]Pa22=[x22y22z22]\begin{align} \b{P_{a11}} &= \left[ {\begin{array}{c} x_{11} \\ y_{11} \\ z_{11} \end{array} } \right] \\ \b{P_{a12}} &= \left[ {\begin{array}{c} x_{12} \\ y_{12} \\ z_{12} \end{array} } \right] \\ \b{P_{a21}} &= \left[ {\begin{array}{c} x_{21} \\ y_{21} \\ z_{21} \end{array} } \right] \\ \b{P_{a22}} &= \left[ {\begin{array}{c} x_{22} \\ y_{22} \\ z_{22} \end{array} } \right] \\ \end{align}

Calculate Cross-Products

We then need to calculate a cross-product for each arc, using the two vectors which define the arc as the inputs to the cross-product:

N1=Pa11×Pa12N2=Pa21×Pa22\begin{align} \b{N_1} = \b{P_{a11}} \times \b{P_{a12}} \\ \b{N_2} = \b{P_{a21}} \times \b{P_{a22}} \\ \end{align}

Note how we switched from using “point” to “vector” terminology here. Before we were talking about two “points” that lie on the sphere and define the arc. This is the easiest way to intuitively understand how two points can fully define an arc on the surface of a sphere. But now we can describe them as vectors. When we talk about them as vectors, we are referring to the vector which starts at [0,0,0][0, 0, 0] and passes through this point.

These cross-products will be perpendicular to the plane that the arc lies in. Hence these cross-products are called plane normals.

Calculate The Cross-Product of the Cross-Product

Yup, you read that correctly. We now calculate the cross-product of the two cross-products we calculated above. Because the two arc cross-products define two plane normals, the cross-product of two plane normals gives us a vector that is in the same direction as the intersecting line.

L=N1×N2\b{L} = \b{N_1} \times \b{N_2}

Because the line must pass through [0,0,0]T[0, 0, 0]^T, this vector is also the line itself!

Find The Intersecting Points

Now that we have the line of intersection, we can now find the two points of intersection between the arc planes by normalizing L (remember, we are using a unit sphere, so normalizing the line reduces it to a vector to one of the intersecting points):

I1=LLI2=I1\b{I_1} = \frac{\b{L}}{||\b{L}||} \\ \b{I_2} = -\b{I_1}

To get the other point of intersection, just take the negative of the first one.

Check If These Intersecting Points Are Within The Original Arc Segments

So we now have the two points of intersection between the arcs. If these arcs were complete great circles, then job done, you are guaranteed that the two above points lie within your great circles. But what if you had great-circle arc segments? These intersection points may not be within both arc segments.

We can determine if these intersection points lies within arc segments by doing a angle test, using the dot product:

ab=abcos(θ)\b{a} \cdot \b{b} = ||\b{a}|| ||\b{b}|| cos (\theta)

We will re-arrange for θ\theta:

θ=arccos(abab)\theta = \arccos \left( \frac{\b{a} \cdot \b{b}}{||\b{a}|| ||\b{b}|| } \right)

The angle from the start of arc 1 to intersecting point 1, and the angle from the end of arc 1 to intersecting point 1, and compare it against the angle between the start and end of arc 1.

The intersecting point is within arc 1 if:

θa11,i1+θa12,i1=θa11,a12\theta_{a11,i1} + \theta_{a12,i1} = \theta_{a11,a12}

Take note that if calculating this result using any data type that can lose precision (e.g. floats, doubles), you will have to check it is close to equal rather than exactly equal. This can be done by adding some epsilon. Usually 1101^{-10} is sufficient.

This test has to be repeated between potential intersecting point 1 and arc 2. If both arcs intersect this potential intersecting point, then with have just confirmed that both arcs do intersect at this point!

The same process has to be applied to potential intersecting point 2 (remember, normalizing the line of intersection between the two planes gives us TWO potential points of intersection).

Further Reading

https://stackoverflow.com/questions/2954337/great-circle-rhumb-line-intersection is a great page showing how to calculate the intersection between a arc and a Rhumb line (similar to above).