Orbite elliptique d'un astre

[Fred_76]

Le bouquin est en anglais, 2nde édition de Willman-Bell, 477 pages. Il suffit de… 

['Bruno]

il y a une heure, LH44 a dit :

mais l'ouvrage reste très difficile à dénicher

On le trouve ici : https://www.astronome.fr/produit-saf30057-calculs-astronomiques-a-l-usage-des-amateurs-Prix-13-euro-id-2092.html .

C'est une version light ? 

[LH44]

Il y a 12 heures, 'Bruno a dit :

On le trouve ici : https://www.astronome.fr/produit-saf30057-calculs-astronomiques-a-l-usage-des-amateurs-Prix-13-euro-id-2092.html .

C'est une version light ? 

 

J'ai les deux ouvrages donc je peux en parler, l'original mais en 2eme édition de 1998 fait 477 pages et à l'époque on l'avais pour 40€ maintenant on le trouve assez facilement d'occasion pour 115€ environ de ce que j'ai vu, l'autre version plus récente de 2014 fait 160 pages et vaut 13€ soit directement chez l'éditeur donc à la boutique de la SAF sur leur site Web soit au magasin de Lorient puisque tu l'a aussi trouvé aussi la bas   Cette dernière édition même abrégée contient l'essentiel et surtout tous les exemples de calculs qui permettent de valider ses propres programmes.

 

Juste une remarque la théorie planétaire utilisée par Jean Meeus la VSOP87 est un peu datée, aujourd'hui il n'y a plus besoin de s'embêter à refaire ces calculs long et fastidieux car la NASA (JPL) nous met à disposition le site HORIZONS où l'on pourra télécharger toutes les positions (et d'autres données à jour) des planètes du systèmes solaire mais aussi celles de n'importe quel astéroïde ou "petit cailloux" répertorié dans le ciel ainsi que tous objets volants identifiés envoyés par l'homme depuis 1967 et qui rôdent encore au dessus de nos têtes et ça un progame Python ou autre saura très bien le faire en automatique.

 

PS: HORIZONS se base sur les toutes dernière théories planétaires DE-XXXX utilisées notamment pour envoyer des sondes dans le système solaire donc la précision fournie devrait largement suffire pour nous petits amateurs ...  

Modifié 28 août 2021 par LH44

['Bruno]

Il y a 3 heures, LH44 a dit :

car la NASA (JPL) nous met à disposition le site HORIZONS où l'on pourra télécharger toutes les positions (et d'autres données à jour)

 

Tout à fait, ce site est un must ! Pour les observateurs, il donne non seulement les positions des planètes mais aussi certains paramètres physiques comme la longitude du méridien central, la phase, etc.

 

L'intérêt de faire les calculs soi même, c'est d'essayer de comprendre comment ça marche. Cela dit, c'est compliqué. Disons que ça m'a appris à quel point c'est compliqué !

[LFiacre]

Bonjour, pour l'anomalie moyenne j'ai fait cette petite fonction.

def anomalie_moyenne(M, n, t, t_0):
    """ t: temps = 2459451.5 jours juliens
        n = moyen mouvement = 0.017204609895805228 radians/jours
        t_0 = instant passage périhélie = 2459217.0416666665 jours juliens
        M en radians, sortie en degrée"""
    M = n*(t-t_0)
    M = np.rad2deg(M) # radian_en_degree
    if M > 360:
        M-= 360
    elif M < 0:
        M+=360
    return M # 232.10341355880192 degrée

 

J'ai parcouru les codes sur github du livre de Jean Meeus. notamment pour l'anomalie excentrique qui m'aide énormément : https://github.com/pavolgaj/AstroAlgorithms4Python/blob/master/kepler.py

Il propose 4 méthodes différentes pour l'obtenir, pourquoi propose-t-il quatre méthodes ? qu'elles sont les enjeux de ces différentes formulations ?

 

J'ai repris ce que j'avais écrit pour l'anomalie excentrique:

def anomalie_excentrique(M, e):
    """ M: anomalie moyenne en degrée = 232.10341355880192
        e : excentricité = 0.01671123
        E : sortie en degrée
    """
    tolerance = 1e-09
    M = np.deg2rad(M) # degrée en radian
    E1 = M 
    answer = False
    while answer == False:
        E = M+e*np.sin(E1) # formule de la 1ère méthode proposé par Jean Meeus
        if abs(E - E1) < tolerance:
            answer = True
        E1 = E
    while E > (2*np.pi): # meme but que dans la fonction de anomalie moyenne mais pour les radians
        E = E - 2*np.pi
    return np.rad2deg(E) # 231.35558432523473 degrée

 

Vu la très faible excentricité de la terre, ca me semble être un résultat correct d'avoir moins de 1 degré de différence entre anomalie moyenne et excentrique ?

Je vous remercie du temps que vous m'accordez, bonne journée à vous !

['Bruno]

Ceci :

il y a une heure, LFiacre a dit :

    if M > 360:
        M-= 360
    elif M < 0:
        M+=360

ne marchera pas lorsque t sera suffisamment grand, ou petit. Tu risque d'avoir un M valant plus de 720, ou plus de 1080, ou moins de -720, etc. Une façon de faire serait d'itérer jusqu'à ce que M soit compris entre 0 et 360, mais ça peut être long. La bonne méthode, c'est d'utiliser math.fmod, la fonction dédiée.

 

Pour l'anomalie excentrique, il y a plusieurs méthodes mathématiquement équivalentes, mais pouvant donner plus ou moins de calculs selon les paramètres. Par exemple la méthode itérative (1ère méthode de Meeus) est d'autant plus rapide que l'excentricité est proche de 0 (du genre 3 ou 4 itérations suffisent). Mais si l'excentricité s'approche de 1, elle peut nécessiter des centaines voire des milliers d'itérations. Dans ce cas (fréquent avec les comètes), il faut utiliser la méthode de Newton (je ne sais pas si Meeus en parle).

 

Et là aussi tu pourrais utiliser math.fmod :

il y a une heure, LFiacre a dit :

while E > (2*np.pi): # meme but que dans la fonction de anomalie moyenne mais pour les radians
        E = E - 2*np.pi

 

Et pour la question

il y a une heure, LFiacre a dit :

Vu la très faible excentricité de la terre, ca me semble être un résultat correct d'avoir moins de 1 degré de différence entre anomalie moyenne et excentrique ?

 

Je pense que oui.

 

[LFiacre]

Merci Bruno, oui tu as raison, j'ai normalisé pour avoir une solution de modulo :

E = np.rad2deg(E)
    E = E%360

 

[LFiacre]

Enfin il me reste le calcul de l'anomalie vraie :

E = 231.35558432523473 # degrée, anomalie excentrique
e = 0.01671123 # excentricité
E = np.deg2rad(E) # radians

# http://www.heliodon.net/downloads/Beckers_2010_Helio_006_fr_2.pdf, début page 7
v_1 = np.arccos((np.cos(E)-e)/(1-e*np.cos(E)))
v_2 = np.arcsin((np.sqrt(1-e**2)*np.sin(E))/(1-e*np.cos(E)))

# http://clea-astro.eu/archives/cahiers-clairaut/CLEA_CahiersClairaut_028_10.pdf page 16
v_3 = np.arctan(np.sqrt(1+e/1-e)*np.tan(E/2))/2

v_1 = np.rad2deg(v_1) #  129.38840424533467 degrés
v_2 = np.rad2deg(v_2) # -50.61159575466535 degrés
v_3 = np.rad2deg(v_3) # -32.161103918691325 degrés

Il y a quelque chose qui m'échappe totalement dans ces calculs, j'en appelle encore à votre aide

['Bruno]

Ici il y a deux erreurs :

Il y a 14 heures, LFiacre a dit :

v_3 = np.arctan(np.sqrt(1+e/1-e)*np.tan(E/2))/2

 

1) 1+e/1-e signifie 1+(e/1)-e. Sur une ligne, il faut écrire (1+e)/(1-e). 

2) Ce n'est pas v_3 = arctan(truc) / 2 mais v_3 = 2 * arctan(truc).

 

L'arc tangente donne un angle entre -90° et +90°, donc v_3 est compris entre -180° et +180°, c'est qu'on voulait.

 

Concernant V_1 et v_2, c'est parce qu'il faut utiliser les deux formules, pas juste une seule. Le cosinus seul fournit (via arccos) un angle entre 0° et 180° ; le sinus seul fournit (via arcsin) un angle entre -90° et +90°. Ça, c'est une source d'erreur vraiment classique : il faut utiliser les deux valeurs, le cosinus et le sinus. (C'est d'ailleurs ce que dit l'article : les deux formules sont liées.)

 

Voici comment un procède :

cosv = (cos(E)-e)/(1-e*cos(E))            # on calcule le cosinus, qui est une abscisse x

sinv = (sqrt(1-e**2)*sin(E))/(1-e*cos(E)) # on calcule le sinus, qui est une ordonnée y

puis on utilise la fonction atan2 (arc tangente à 2 paramètres) :

v_1 = atan2(sinv, cosv)

Attention : y en premier paramètre, x en second. C'est parce que ça ressemble au calcul de atan(y/x).

 

J'essaie les deux méthodes :

>>> from math import cos, sin, tan, atan, atan2, pi, sqrt 
>>> E = 231.35558432523473 * pi/180
>>> e = 0.01671123
>>> cosv = (cos(E)-e)/(1-e*cos(E))
>>> sinv = (sqrt(1-e**2)*sin(E))/(1-e*cos(E))
>>> v_1 = atan2(sinv, cosv)
>>> print(v_1 * 180/pi)
-129.38840424533467  # il a le même cos que 129.38840424533467° et même sin que -50.61159575466535°
>>> v_3 = 2 * atan(sqrt((1+e)/(1-e))*tan(E/2))
>>> print(v_3 * 180/pi)                         
-129.38840424533467  # les maths, c'est une science exacte ! :-)

[LFiacre]

Merci beaucoup Bruno, je t'avoue qu'il me faut digérer ce que tu as écris, j'ai souvent l'impression de comprendre globalement, mais dans les faits y a plus personne. Merci encore !!