Unify computation of distances in spherical and Cartesian coordinates

[santisoler]

Description of the desired feature

To prevent repeated code, would be nice to have single functions for computing distances between two points in both spherical and Cartesian coordinates. They must be Numba jitted functions to speed things up.

For the spherical coordinates case, the points are given in longitude, latitude and radius coordinates, but several trigonometric functions must be applied in order to ultimately get the distance between these two points. Some forward models precompute these quantities in order to save time. So, would be nice to have a private function that takes all these precomputed quantities as arguments and returns the distance, and a public function that just takes the points' coordinates, compute this quantities and use this private function to get the distance between them.

Moreover, some forward models just needs the distance square, so we could make another private function that just returns the distance square, i.e. doesn't apply the square root.

For example:

@jit(nopython=True)
def distance_spherical(point_a, point_b):
    # Get coordinates of the two points
    longitude, latitude, radius = point_a[:]
    longitude_p, latitude_p, radius_p = point_b[:]
    # Convert angles to radians
    longitude, latitude = np.radians(longitude), np.radians(latitude)
    longitude_p, latitude_p = np.radians(longitude_p), np.radians(latitude_p)
    # Compute trigonometric quantities
    cosphi_p = np.cos(latitude_p)
    sinphi_p = np.sin(latitude_p)
    cosphi = np.cos(latitude)
    sinphi = np.sin(latitude)
    return _distance_spherical(
        longitude, cosphi, sinphi, radius, longitude_p, cosphi_p, sinphi_p, radius_p
    )


@jit(nopython=True)
def _distance_spherical(
    longitude, cosphi, sinphi, radius, longitude_p, cosphi_p, sinphi_p, radius_p
):
    return np.sqrt(
        _distance_sq_spherical(
            longitude, cosphi, sinphi, radius, longitude_p, cosphi_p, sinphi_p, radius_p
        )
    )


@jit(nopython=True)
def _distance_sq_spherical(
    longitude, cosphi, sinphi, radius, longitude_p, cosphi_p, sinphi_p, radius_p
):
    coslambda = np.cos(longitude_p - longitude)
    cospsi = sinphi_p * sinphi + cosphi_p * cosphi * coslambda
    return (radius - radius_p) ** 2 + 2 * radius * radius_p * (1 - cospsi)

Are you willing to help implement and maintain this feature? Yes