Naive approach to average latitude and longitude positions

geobar.py

import webbrowser
from functools import cache

import numpy as np
import numpy.linalg as lin

from geopy.geocoders import Nominatim


@cache
def get_locator():
    return Nominatim(user_agent="cities_barycenter")


@cache
def get_position(city_name):
    loc = get_locator().geocode(city_name)
    return loc.latitude, loc.longitude


# Most of the code stolen from https://stackoverflow.com/a/68243962/1209937

# https://en.wikipedia.org/wiki/N-vector
E = np.array([
    [0, 0, 1],
    [0, 1, 0],
    [-1, 0, 0]
])


def lat_long_to_n_E(latitude, longitude):
    res = [
        np.sin(np.deg2rad(latitude)),
        np.sin(np.deg2rad(longitude)) * np.cos(np.deg2rad(latitude)),
        -np.cos(np.deg2rad(longitude)) * np.cos(np.deg2rad(latitude))
    ]
    return np.dot(E.T, np.array(res))


def n_E_to_lat_long(n_E):
    n_E = np.dot(E, n_E)
    longitude = np.arctan2(n_E[1], -n_E[2]);
    equatorial_component = np.sqrt(n_E[1] ** 2 + n_E[2] ** 2)
    latitude = np.arctan2(n_E[0], equatorial_component)
    return np.rad2deg(latitude), np.rad2deg(longitude)


def average(coords):
    res = []
    for lat, lon in coords:
        res.append(lat_long_to_n_E(lat, lon))
    res = np.array(res)
    m = np.mean(res, axis=0)
    m = m / lin.norm(m)
    return n_E_to_lat_long(m)


def show_location(lat, lon):
    url = f"https://www.google.com/maps/search/?api=1&query={lat},{lon}"
    webbrowser.open(url)


if __name__ == "__main__":
    from sys import argv
    cities = [x.strip() for x in argv[1].split(",")]
    positions = [get_position(city) for city in cities]
    centroid = average(positions)
    show_location(*centroid)