Ich bin kein SPICE-Experte, aber hier sind einige mögliche Lösungen (es sei denn, Sie möchten natürlich die Almanach-Methoden von Skyfield ausprobieren !)

Erste Möglichkeit, aber alle diese funktionieren möglicherweise nicht

diese Antwort verlinkt auf diesen Kommentar Links zu SPICE-Tutorials, aktualisiert am 11. Dezember 2019 Links zur 69-Folien-Präsentation SPICE Geometry Finder (GF) Subsystem; Suche nach Zeiten, in denen bestimmte geometrische Bedingungen vom Januar 2020 auftreten .

• Viel SPICE-Software berechnet einen Geometrieparameter zu einem gegebenen Zeitpunkt t, dh x = f(t). Beispiel: Wie hoch ist die Raumsonde am 30. März 2011 um 14:57:08 Uhr über dem Mars?
• Das Geometry Finder-Subsystem macht das Gegenteil: Es findet Zeiten, in denen bestimmte geometrische Ereignisse auftreten

GF bietet zwei Haupttypen von APIs zur Ereignissuche

• Boolesch: eine geometrische Bedingung (ein Ereignis) ist wahr oder falsch – Numerisch: eine geometrische Größe hat einen bestimmten Wert, liegt innerhalb eines bestimmten Bereichs oder hat ein lokales oder globales Maximum oder Minimum erreicht

In Ihrem Fall würden Sie eine Position auf einem Körper angeben und dann nach Zeiten suchen, in denen die Höhe 0 oder 0 +/- klein war.

Eine andere Möglichkeit, dies scheint genau das zu sein, wonach Sie suchen

• https://naif.jpl.nasa.gov/pub/naif/toolkit_docs/C/cspice/spkcpo_c.html

  Beispiele

  Die für diese Beispiele gezeigten numerischen Ergebnisse können je nach Plattform unterschiedlich sein. Die Ergebnisse hängen von den als Eingabe verwendeten SPICE-Kerneln, dem Compiler und den unterstützenden Bibliotheken sowie der maschinenspezifischen
  arithmetischen Implementierung ab.

  1. Berechnen Sie das scheinbare Azimut und die Elevation der Sonne, wie sie von einem bestimmten Oberflächenpunkt auf der Erde aus gesehen werden.

     Aufgabenbeschreibung

     In diesem Beispiel verwenden wir den Standort der DSN-Station
     DSS-14 als unseren Oberflächenpunkt.

     Wir führen die Berechnung des Sonnenazimuts und der Elevation auf zwei Arten durch:

     • Verwenden eines Stationsrahmenkerns, um die
       Spezifikation eines topozentrischen Referenzrahmens bereitzustellen, der bei DSS-14 zentriert ist.

     • Inline-Berechnung der Transformation vom erdfesten, erdzentrierten Rahmen ITRF93 zu einem topozentrischen Rahmen, der bei DSS-14 zentriert ist.

     Beachten Sie, dass sich die Ergebnisse der beiden Berechnungen geringfügig unterscheiden. Es gibt drei Quellen für die Unterschiede:

     1. Die Stationsposition ist aufgrund der tektonischen Plattenbewegung zeitabhängig, und die Epochen der Stationspositionen, die zur Angabe der Achsen des topozentrischen Rahmens verwendet werden, sind in den beiden Fällen unterschiedlich. Dadurch ergeben sich unterschiedliche Ausrichtungen der Rahmenachsen relativ zum Rahmen ITRF93.

     2. Die beiden Berechnungen verwenden unterschiedliche Erdradien; dies führt zur Berechnung verschiedener geodätischer Breiten der Station. Dieser Unterschied wirkt sich auch auf die topozentrische Rahmenorientierung relativ zu ITRF93 aus.

     3. Die Stationsbewegung zwischen ET und der Epoche, in der das DSS-14_TOPO-Frame spezifiziert ist, trägt einen sehr kleinen Versatz – in der Größenordnung von 10 cm – zum Stations-Sonnenpositionsvektor bei, der im ITRF93-Frame ausgedrückt wird.

und da ist noch viel mehr.

Es scheint, dass Sie dasselbe in Python mit SpicyPy tun können

Ermitteln Sie Azimut und Höhe der scheinbaren Position des Mondes, wie sie von der DSN-Station DSS-13 aus gesehen wird, indem Sie die folgenden Schritte ausführen:

11. Ermitteln Sie den scheinbaren Positionsvektor des Mondes relativ zur DSN-Station DSS-13 im topozentrischen Referenzrahmen DSS-13_TOPO zur Epoche ET. Verwenden Sie Lichtzeit- und Sternfehlerkorrekturen.

Für diesen Schritt müssen Sie eine Stations-SPK-Datei geladen haben, die geozentrische Stationspositionsvektoren bereitstellt, sowie einen Rahmenkern, der topozentrische Referenzrahmen angibt, die an den jeweiligen DSN-Stationen zentriert sind. (Andere Kernel werden ebenfalls benötigt; Sie müssen diese auswählen.)

12. Wandeln Sie den Positionsvektor in Breitenkoordinaten um. Verwenden Sie für diese Berechnung die Routine spiceypy.reclat.

13. Berechnen Sie Azimut und Höhe des Mondes wie folgt: Azimut ist das Negative der topozentrischen Länge und liegt im Bereich von 0 bis 360 Grad; Die Höhe ist gleich der topozentrischen Breite. Anzeige der Ergebnisse in Grad.