From c1eb37524da1ec9493f2941dd9aa44f5cd83856e Mon Sep 17 00:00:00 2001
From: "Michael J. Roberts" <michael@observerly.com>
Date: Thu, 18 May 2023 14:52:55 +0100
Subject: [PATCH] feat: Added get_angular_separation to celerity/astrometry
 module.

feat: Added get_angular_separation to celerity/astrometry module.
---
 src/celerity/astrometry.py | 34 +++++++++++++++++++++++++++++++++-
 tests/test_astrometry.py   | 10 ++++++++++
 2 files changed, 43 insertions(+), 1 deletion(-)

diff --git a/src/celerity/astrometry.py b/src/celerity/astrometry.py
index 62a35e2..61a66f3 100644
--- a/src/celerity/astrometry.py
+++ b/src/celerity/astrometry.py
@@ -7,7 +7,7 @@
 # *****************************************************************************************************************
 
 from datetime import datetime
-from math import atan2, cos, degrees, pow, radians, sin, tan
+from math import acos, atan2, cos, degrees, pow, radians, sin, tan
 
 from .common import EquatorialCoordinate, GeographicCoordinate
 from .temporal import get_julian_date, get_local_sidereal_time
@@ -15,6 +15,38 @@
 # *****************************************************************************************************************
 
 
+def get_angular_separation(A: EquatorialCoordinate, B: EquatorialCoordinate) -> float:
+    """
+    The angular separation between two objects in the sky is the angle between the two objects
+    as seen by an observer on Earth.
+
+    :param A: The equatorial coordinate of the observed object.
+    :param B: The equatorial coordinate of the observed object.
+    :return The angular separation in degrees between target A and target B.
+    """
+    # Calculate the angular separation between A and B (in degrees):
+    θ = (
+        degrees(
+            acos(
+                sin(radians(A["dec"])) * sin(radians(B["dec"]))
+                + cos(radians(A["dec"]))
+                * cos(radians(B["dec"]))
+                * cos(radians(A["ra"] - B["ra"]))
+            )
+        )
+        % 360
+    )
+
+    # Correct for negative angles:
+    if θ < 0:
+        θ += 360
+
+    return θ
+
+
+# *****************************************************************************************************************
+
+
 def get_hour_angle(date: datetime, ra: float, longitude: float) -> float:
     """
     Gets the hour angle for a particular object for a particular observer at a given datetime
diff --git a/tests/test_astrometry.py b/tests/test_astrometry.py
index 1b9c4f6..db1ecea 100644
--- a/tests/test_astrometry.py
+++ b/tests/test_astrometry.py
@@ -2,6 +2,7 @@
 from typing import TypedDict
 
 from src.celerity.astrometry import (
+    get_angular_separation,
     get_hour_angle,
     get_obliquity_of_the_ecliptic,
     get_parallactic_angle,
@@ -20,9 +21,18 @@
 
 betelgeuse: EquatorialCoordinate = {"ra": 88.7929583, "dec": 7.4070639}
 
+arcturus: EquatorialCoordinate = {"ra": 213.9153, "dec": 19.182409}
+
+spica: EquatorialCoordinate = {"ra": 201.2983, "dec": -11.1614}
+
 observer: GeographicCoordinate = {"lat": latitude, "lon": longitude}
 
 
+def test_get_angular_separation():
+    θ = get_angular_separation(arcturus, spica)
+    assert θ == 32.79290589269233
+
+
 def test_get_hour_angle():
     ha = get_hour_angle(date, betelgeuse["ra"], longitude)
     assert ha == 347.6988036852858