diff --git a/src/sorcha/ephemeris/simulation_geometry.py b/src/sorcha/ephemeris/simulation_geometry.py
index 9e3f07e4..dc226525 100644
--- a/src/sorcha/ephemeris/simulation_geometry.py
+++ b/src/sorcha/ephemeris/simulation_geometry.py
@@ -195,9 +195,16 @@ def barycentricObservatoryRates(et, obsCode, observatories, Rearth=RADIUS_EARTH_
     pos = posvel[0:3]
     vel = posvel[3:6]
     # Get the matrix that rotates from the Earth's equatorial body fixed frame to the J2000 equatorial frame.
-    m = spice.pxform("ITRF93", "J2000", et)
-    mp = spice.pxform("ITRF93", "J2000", et + delta_et)
-    mm = spice.pxform("ITRF93", "J2000", et - delta_et)
+    # We don't need to compute this repeatedly
+    et_1962 = spice.str2et("1962-Jan-20")
+    if et >= et_1962:
+        m = spice.pxform("ITRF93", "J2000", et)
+        mp = spice.pxform("ITRF93", "J2000", et + delta_et)
+        mm = spice.pxform("ITRF93", "J2000", et - delta_et)
+    else:
+        m = spice.pxform("IAU_EARTH", "J2000", et)
+        mp = spice.pxform("IAU_EARTH", "J2000", et + delta_et)
+        mm = spice.pxform("IAU_EARTH", "J2000", et - delta_et)
     # Get the MPC's unit vector from the geocenter to
     # the observatory
     obsVec = observatories.ObservatoryXYZ[obsCode]
diff --git a/tests/ephemeris/test_simulation_parsing.py b/tests/ephemeris/test_simulation_parsing.py
index 46717123..3ba77fed 100644
--- a/tests/ephemeris/test_simulation_parsing.py
+++ b/tests/ephemeris/test_simulation_parsing.py
@@ -2,6 +2,8 @@
 import sorcha.ephemeris.simulation_parsing as sp
 from sorcha.utilities.dataUtilitiesForTests import get_test_filepath
 from sorcha.utilities.sorchaConfigs import auxiliaryConfigs
+import sorcha.ephemeris.simulation_geometry as sg
+import sorcha.ephemeris.simulation_setup as ss
 
 
 def test_observatory_compared_to_original():
@@ -27,3 +29,22 @@ def test_observatory_for_moving_observatories():
     obs = observatory.ObservatoryXYZ
 
     assert obs["250"] == (None, None, None)
+
+
+def test_observatory_before_1962():
+    auxconfigs = auxiliaryConfigs()
+    observatory = sp.Observatory(
+        auxconfigs=auxconfigs, args=None, oc_file=get_test_filepath("ObsCodes_test.json")
+    )
+
+    et = -1229083166.3680687  # thanks to Rahil Makadia
+
+    pos, _ = sg.barycentricObservatoryRates(et, "Z20", observatory)
+
+    x, y, z = -7.426461821563405e07, 1.177545544454091e08, 5.105719899396534e07  # taken from JPL
+
+    # note these should not match precisely because pre-1962 we're using an approximation
+    # but this is good enough
+    assert np.isclose(pos[0], x, 5)
+    assert np.isclose(pos[1], y, 5)
+    assert np.isclose(pos[2], z, 5)