Skip to main content
SZCZ

Best meeting point

Meeting points

I grew up in a small town (an ex-shtetl, actually) and lived very close to all of my friends. We had a few classic meet-up points and we'd sometimes discuss the "best" meeting point: the most fair place to meet given our individual starting points.

This exact scenario came up again recently, so I decided to attack it head on this time. I had an inkling that there wouldn't necessarily be a single answer. In fact each notion of "fair" lands on a different spot, and which one you'd root for turns out to depend on where you happen to live. And I always wanted to try out some OpenStreetMap[1] integrations; so I figured this would be a good occasion to do so. I'll go through a few different ways of looking at this problem.

Setup

For all the scenarios below I'll be assuming we have 3 points on a sphere: A, B, C. I'll assume we have us a proper sphere and disregard any equatorial bulges and other such nonsense.

Sub note: I focus on the 3-point case but two of the three centre points below generalise nicely to more than 3 points. One doesn't fully but there's a way to emulate it. More on that later.

Average

Let's start with the simplest meeting point: the average of the 3 points. There's one issue which is that if we blindly take A+B+C3 we're guaranteed that the result will put us below the surface of the sphere: the mean will be the centroid of the triangle we get from the 3 points. So unless we have A=B=C[2] we will get a point inside the Earth. Not good. To get around this we can normalise this point and then multiply it by the radius r to put us back on the surface. P=A+B+CA+B+Cr

Diagram showing the mean of three sphere points

This is a good start. E.g. if we have Alice, Bob and Charlie living in London, Rome and Warsaw, respectively, the mean point would land near Eichstätt in Germany. Very sensible. As long as we keep normalising we can guarantee this point lies on the sphere even as we add friend-points. It's very straightforward to calculate and... it feels like an okay spot at first.

Equidistance

One could argue that the mean isn't a very smart meeting point: just because it's the geometric centre doesn't make it a good meeting point. What if we require that everyone has to travel the same distance to meet? That's fair! Let's look for the circumcenter of the three points.[3]

Let's pretend for a minute we're in a 2D plane and not on a sphere. Finding a circumcenter between 3 points is then the same as finding the point where the three perpendicular bisectors meet. I.e. the perpendicular bisector of AB is the line of points equidistant between A and B; similarly for BC and AC. The point where those 3 meet is what we're looking for.

Diagram showing the meeting point of the perpendicular bisectors of 3 distinct plane points

This naturally extends to the sphere case if we replace perpendicular bisectors with a great circle that's equidistant from the two points. Then we find the point where the 3 great circles meet.

We can calculate this more easily by finding the normal to the plane going through our points.
n=(AB)×(AC)P=±n|n|r

Note that we will have two antipodal solutions on the sphere (hence the ±). Imagine three points laid out in a circle around the North Pole; one circumcenter is the North Pole and the other is the South Pole. If the three points lie on a common great circle, the circumcenters are the two poles of that great circle, each 90° from every point on it.

While it might seem like this is the most fair meeting point I think it's easy to show that it's not really. Let's say that Bob moves from Rome to just north of Prague. The circumcenter then will land somewhere near Porsgrunn in Norway. It would be a bit better if they moved to Zagreb; the circumcenter would then go into the centre of Germany.

The issue here is that the circumcenter will only stay inside the ABC triangle for points that don't form obtuse angles between each other. As soon as that happens (or for degenerate cases like cities on a common great circle) the result is not very practical. For example: if Bob moves yet again to Magdeburg in Germany the point equidistant to Alice, Bob and Charlie will be off the coast of Alaska :(

Notice that even our "silly" average is better for cases like this - it will always land inside the triangle.

There's one more catch worth flagging now: this exact equidistant point is really a three-point luxury. With four or more friends there's generally no single spot the same distance from all of them - the perpendicular bisectors stop meeting at one point. The usual fix is to stop demanding exactly equal and least-squares it instead: find the point that minimises the spread of the distances, i.e. the closest thing to equidistant you can manage. That's the method I warned didn't quite generalise back in the intro.

Median

Let's take a step back. So far we've considered ways of calculating our meeting point assuming that each point is equally important: either it forms a third of the answer OR we require each person to travel the same distance.

What if we instead looked at minimising the total distance[4] travelled? This is the eco-friendly way of thinking. Now, we could consider either the sum of the distances or the sum of the distances squared (we could use higher powers as well but I will ignore those).

I will argue that minimising the squared distances is not the point we want. By squaring the distances we are weighting the far off point more heavily... and I don't think that's very fair (unless you're that friend who lives very far from all your other friends...). I think the most sensible way to minimise this total travel distance is to look at the non-squared sum instead. I will refer to this point as the median; the one minimising the sum of distances squared is referred to as the Fréchet mean.

We had nice, closed-form formulas for everything so far. Not quite the case with the median and the Fréchet mean. The way you calculate these is that you start somewhere on the sphere (e.g. at the mean point) and iteratively optimise towards the solution.

Also since we'll be using gradient descent it means we are actually doing AI!! We're using AI to compute the Fréchet mean!!!!
import numpy as np

def log_map(mu, x, tol=1e-12):
    dot = np.clip(np.dot(mu, x), -1, 1)
    theta = np.arccos(dot)
    if theta < tol:
        return np.zeros(3)
    return theta / np.sin(theta) * (x - dot*mu)

def exp_map(mu, v, tol=1e-12):
    theta = np.linalg.norm(v)
    if theta < tol:
        return mu
    return np.cos(theta)*mu + np.sin(theta)*v/theta

def frechet_mean(points, tol=1e-10, max_iter=100):
    mu = points.mean(axis=0)
    mu /= np.linalg.norm(mu)
    for _ in range(max_iter):
        v = np.mean([log_map(mu, x) for x in points], axis=0)
        if np.linalg.norm(v) < tol:
            break
        mu = exp_map(mu, v)
    return mu
    
# naive first attempt - we'll see below why this one misbehaves
def median(points, tol=1e-10, max_iter=100, step=0.25):
    mu = points.mean(axis=0)
    mu /= np.linalg.norm(mu)
    for _ in range(max_iter):
        v = np.zeros(3)
        for x in points:
            w = log_map(mu, x)
            norm = np.linalg.norm(w)
            if norm > tol:
                v += w / norm
        v /= len(points)
        if np.linalg.norm(v) < tol:
            break
        mu = exp_map(mu, step * v)
    return mu
Here's how you calculate this Assume we are dealing with the general case of n points on a sphere x1,x2,,xnS2. We'll use the unit sphere - remember that we can always scale up by the real radius. The Fréchet mean minimises id(p,xi)2 while our median minimises id(p,xi)

Our "great circle" distance is just d(p,x)=arccos(px) since the dot product of our two sphere vectors will be cos(θ) (where θ is the angle between them) and that angle (in radians) is exactly the great circle distance.

Diagram showing the great circle distance is just the angle between sphere points

We seed our search with the mean μ=ixiixi.

Now we need to work out how to map our sphere points to tangent vectors. This is the log map. Let's take our mean μ and x (one of the xis). If we look at xμ we'll get a vector pointing from μ to x BUT also pointing slightly inward. We can fix it by removing the projection of x onto μ.

So, we want to calculate x(μx)μ as you can see on this expertly made drawing.

Diagram showing how we arrive at the tangent vector to a point in the direction of a different point

It's easiest to compute the squared norm x(μx)μ2 where we use
x=1μ=11cos(θ)2=sin(θ)2xμ=cos(θ)
to get x(μx)μ2=sin2(θ) and hence x(μx)μ=sin(θ)

So now we need to normalise by sin(θ) to get the right length (equal to the great circle distance - θ)
So in the end we have:
θsin(θ)(x(μx)μ)
or in code

import numpy as np

def log_map(mu, x, tol=1e-12):
    dot = np.clip(np.dot(mu, x), -1, 1)
    theta = np.arccos(dot)

    if theta < tol:
        return np.zeros(3)

    return theta / np.sin(theta) * (x - dot*mu)

We will also want the opposite: something to turn a tangent vector BACK to a point on a sphere. Suppose we have a tangent vector v with length (as calculated previously) θ.

We normalise u=vθ. Now u is a unit tangent vector. We also have our mean point μ. We can use these two points to make a great circle on the sphere, similar to a way we generate a circle on the plane.

On the plane we can start with the point (1,0) and move an angle θ to arrive at (cos(θ),sin(θ)). And in our case, instead of (1,0) we have (μ,u) and so we get cos(θ),μ+sin(θ)vθ.

In code

def exp_map(mu, v, tol=1e-12):
    theta = np.linalg.norm(v)

    if theta < tol:
        return mu

    return np.cos(theta)*mu + np.sin(theta)*v/theta

Then the glue is just iterating for up to max_iter times. In the Fréchet case we use the mean of the tangent vectors to move our mean point μ. In the median case we instead normalise each tangent vector to unit length first (that's the naive version - we'll fix its convergence below).

What this means in practice is that if we have a cluster of points with some outlier point further out, the median point will be positioned very close to the cluster; it's on the outliers to travel to the cluster.

In a situation where there is a "central" point X and two points on either side the median tends to lie very close to X - that's what minimises total travel after all. For example, in the crazy case from earlier with Bob moving to Magdeburg we'd get a median point situated at Bob's house: Alice and Charlie live on opposite sides in London and Warsaw and it makes most sense for them to fly in.

In general, any sort of clustering will weigh very heavily for this one.

Median or Fréchet? A worked example

To feel the difference between the two, keep Alice in London and put Bob and Charlie up in Liverpool and Newcastle upon Tyne. Now Bob and Charlie are the northern pair and Alice is the lone southerner.

The median lands just north of Warrington - it tucks in near the northern pair and leaves Alice to make the trek down from London. The Fréchet mean gets dragged noticeably further south-east, out toward the Peak District (it lands almost exactly on the plain mean[5] - over distances this small the sphere is near enough flat that the two barely differ). That pull is the squaring at work: penalising Alice's large distance squared drags the answer toward London so the burden is shared, exactly the behaviour I grumbled about earlier.

It gets dramatic if Charlie swaps Newcastle for Carlisle. Now Liverpool sits right between Carlisle and London, and the median snaps onto Liverpool itself - Bob doesn't move an inch, and it's on Alice and Charlie to come to him. The Fréchet mean refuses to let Bob off so lightly: it only edges west and still strands everyone out in open country. Which of those two you'd call "fair" is precisely the argument you'd be having with your friends.

Wait, does the code actually agree?

If you take the median function from above and run it on this exact London / Magdeburg / Warsaw setup, you will not get Magdeburg. With step=0.25 it swings wildly back and forth between about 47° N and 58° N; drop the step down to 0.1 and it will settle into a tighter cycle around the 52nd parallel. Changing max_iter won't help since we fall into a cycle in both cases.

The direction of the update is okay but since we have a constant step size it's easy for us to move too far and never actually settle properly in some edge cases.

E.g. London and Warsaw are on opposite sides, so their tangent vectors very nearly cancel, but not quite, and whatever is left over keeps v away from zero. Our fixed step then keeps kicking us off Magdeburg back and forth for ever and ever.

The standard fix is Weiszfeld's algorithm. Instead of a hand-tuned step, we move to the inverse-distance-weighted average of the points. Each point gets a weight of 1/di (so nearby points weigh harder) and we divide by the total weight:

μnext=expμ!(i1di,logμ(xi)i1di)

I.e. we scale the step dynamically instead of using a constant one. We also add a short-circuit for when we're already sitting on one of the input points (that point is the median in that case).

def median(points, tol=1e-9, max_iter=300):
    mu = points.mean(axis=0)
    mu /= np.linalg.norm(mu)
    for _ in range(max_iter):
        num = np.zeros(3)
        den = 0.0
        for x in points:
            d = np.arccos(np.clip(np.dot(mu, x), -1, 1))
            if d < tol:           # already on a data point => that's the median
                return x
            num += (x - np.dot(mu, x) * mu) / np.sin(d)  # unit tangent toward x
            den += 1.0 / d
        v = num / den
        if np.linalg.norm(v) < tol:
            break
        mu = exp_map(mu, v)
    return mu

With that we don't get the cycling behaviour anymore.

Demo

Here's a little demo where you can move the 3 points around and see the mean, circumcenter, median and Fréchet mean adjust accordingly (toggle any of them on or off in the legend). Some things to try:

  • with 2 points held in place move the 3rd around and notice when the circumcenter jumps outside the triangle
  • try lining up the points on the same great circle; what happens to the circumcenter? or the median?
  • drag one point far away from the other two and watch the median cling to the pair while the Fréchet mean drifts out toward the outlier
  • check the sum row in the table: the median always has the smallest total distance (highlighted), since that's exactly what it minimises
  • try the Liverpool / London / Newcastle setup from above, then shift Newcastle to Carlisle and watch the median snap onto Liverpool while the Fréchet mean stays out in the open
  • try the best meeting point between London, Berlin and Warsaw

Wrap up

I guess the end result is that none of the four solutions we considered is ideal. If you're a person that has your friends on either side you should be advocating for the median: you won't have to move at all. If you're the outlier friend who lives far from everyone else you should advocate for the Fréchet mean OR for the circumcenter (if one exists); you'll have a bit to go anyway but so will the others. And if you don't like the idea of the circumcenter and are disgusted by gradient descent based optimisation you should go for the mean. It's the simplest - and still a good starting point.


  1. Love it. ↩︎

  2. OR the Flat Earth Theory is correct; in that case the mean works out the box! ↩︎

  3. Also, I've calculated the circumcircle analytically before; see Monte Carlo ↩︎

  4. This is distance travelled via "great circle" ↩︎

  5. They're not quite the same thing - this is the difference between the extrinsic and the intrinsic mean. Our plain mean minimises the squared straight-line (chord) distances that cut through the planet; the Fréchet mean minimises the squared great-circle distances along the surface. But for points close together the chord and the arc are nearly the same length (the chord only really starts cutting the corner once the angle gets big), so the two land almost on top of each other. You need points scattered across half the globe before they visibly come apart. ↩︎