-- AER 213 - Fondamentaux de conception spatiale --

TN n°2 Notions d'Astronomie et d'astrophysique

 - les bases partie 2 -

Le système solaire, le Soleil et l'heure légale

simon.marie@lecnam.net
In [ ]:
import numpy as np
import scipy as sc
import matplotlib
import matplotlib.pyplot as plt
import datetime
import pandas as pd
import plotly.express as px
import time 
#Option d'affichage et taille de police sur les figures:
fs=20
plt.style.use('ggplot')
plt.rc('xtick',labelsize=fs)
plt.rc('ytick',labelsize=fs)
plt.rc('text', usetex=True)
plt.rcParams["font.family"] = "Ubuntu"
In [ ]:
YEAR=2023
## Constantes utilisées:
c0=299792458 # Vitesse de la lumière dans le vide en m/s

# Diametre solaire en m
Dsol=1392680000
# Masse solaire en kg
Msol=1.9884e30

# Constante de Gravitation
G=6.6743015e-11

# Parametre solaire GM0:
mu0=G*Msol

I - Le système solaire

Cette première partie du TN présente les principales caractéristiques physiques et orbitales des 8 planétes du système solaire.

Source Wikipédia

1. Description des planètes

Le fichier Syst_solaire.csv contient les principales données orbitales et physiques des 8 planètes du système solaire que l'on peut mémoriser dans l'ordre d'éloignement par les première lettre de la phrase "Mercredi Viendra Tu Manger Avec Jean Sur Une Nappe ?" - Le A étant la ceinture d'Astéroïdes.

En regardant le système depuis le pôle nord terrestre, toutes les planètes tounent dans le sens antihoraire autour du soleil. Leur sens de rotation autour de leur axe est également antihoraire sauf pour Venus et Uranus pour lesquelles un observateur à leur surface verrait le soleil se lever à l'Ouest et se coucher à l'Est.

Comme dans le TN1 on utilise pandas pour explorer facilement les données du fichier:

In [ ]:
#! wget https://hpp.education/Lessons/AER213/Files/Syst_solaire.csv
In [ ]:
SystSol=pd.read_csv('Syst_solaire.csv',delimiter=" ")
In [ ]:
SystSol["Size"]=SystSol["Diamètre(km)"]/SystSol["Diamètre(km)"][2]
SystSol
In [ ]:
fig=px.scatter(SystSol, hover_name="Planete",x="Demi-GA(km)", y="Masse(kg)",color="Axial-Tilt(°)",size="Size")
fig.show()

Il est également important de connaitre les ordre de grandeur des données orbitale par rapport à celles de la terre:

In [ ]:
# Diamètres par rapport à celui de la terre:
SystSol["Diamètre(km)"]/SystSol["Diamètre(km)"][2]
In [ ]:
#Distances en UA:
SystSol["Demi-GA(km)"]/SystSol["Demi-GA(km)"][2]

Calculer les demi grand axes des planètes en minutes lumière:

In [ ]:
#Distance en minutes lumière:
dml=...
print(dml)

2. Orbites elliptiques

L'équation d'une conique (Cf Chap 3.) s'écrit en coordonnées polaires:

$$ \rho=\dfrac{p}{1+e\cos(\theta)} $$

Avec:

  • $e$ exentricité
  • $p=b^2/a$ paramètre de la conique
  • $a$ demi grand axe
  • $b=a\sqrt{1-e^2}$

soit donc:

$$ \rho=\dfrac{a(1-e^2)}{1+e\cos(\theta)} $$

En prenant un vecteur azimutale $\theta$ entre 0 et $2\pi$ on peut tracer l'orbite des planètes (on prendra ici les planètes telluriques):

In [ ]:
# Vecteur Azimutale:
theta=np.linspace(0,2*np.pi,100)

# Creation de la figure
fig=plt.figure(figsize=(16,12))
plt.axes(projection='polar')

# On représente les orbites des 4 premières planètes
for i,plan in SystSol[0:4].iterrows():
    a=SystSol["Demi-GA(km)"][i]/SystSol["Demi-GA(km)"][2]
    p=a*(1-SystSol["exentricity"][i]**2)
    plt.plot(theta,p/(1+SystSol["exentricity"][i]*np.cos(theta)),label=SystSol["Planete"][i])
plt.legend(fontsize=fs)

3 - Ephemerides

En astronomie, un éphéméride est un modèle mathématique permettant le calcul de la position des corps du système solaire (planètes, lune, astéroides) sur des intervales de temps régulier. La précision et la qualité d'un modèle d'éphéméride dépend à la fois de la précision des observations et de la complexité mathématique des équations utilisé pour coller aux observations. Ainsi les modèles sont constament corrigés et publiés par les organismes responsable comme l'IMCCE (Institut de mécanique céleste et de calcul des éphémérides) en France. Le modèle le plus utilisé est le modèle INPOP21a (datant de 2021) [1]

Exercice: En utilisant les éphémérides fournis, extraire l'évolution INPOP des coordonnées équatoriales de Jupiter entre 2023 et 2033 et représenter sa trajectoire sur la carte stellaire du TN1.

II - Le soleil et sa trajectoire apparente à l'horizon

Le soleil est un des foyers des orbites planétaires. On s'interresse maintenant à l'orbite terrestre et à la trajectoire apparente du soleil vu depuis la terre:

In [ ]:
terre=SystSol.iloc[2]

1. Jour solaire et jour sidéral

Le jour solaire est défini comme la durée que met le soleil pour revenir à la même position dans le ciel d'un jour à l'autre (1->3 sur la figure). La durée du jour solaire est fixé à 86400 secondes (soit 24 heures).

Source Wikipédia

Le jour sidéral est défini comme la durée que met la terre pour faire un tour sur elle même (1->2 sur la figure). En pratique ce n'est pas tout à fait la même chose que le jour solaire car d'un jour à l'autre elle a parcouru un petit arc de cercle dans sa révolution autour du soleil:

On peut donc calculer la durée du jour sidéral en remarquant qu'au bout d'un an la terre aura bien fait un tour de plus sur elle même soit $366.25$ jours. Ainsi on a:

In [ ]:
Tsid=86400*(365.25/366.25)
print("Tsid =",Tsid)

La durée d'un jour sidéral est donc de 86164 seconde soit 23h 56 min et 4 secondes.

2. Durée de l'année sidérale

L'orbite terrestre permet de décrire sa position autour du soleil. Cette orbite ne fixe pas de calendrier (on ne peut pas a priori positioner le 1er Janvier ou le 13 Mai sur la figure précédente). Le calendrier est principalement fixé par deux paramètres:

  1. La durée de l'année en nombre de jour solaire
  2. La position apparente du soleil dans le ciel à midi

A partir de la troisième loi de Kepler, la periode de révolution d'une planète en fonction de son demi-grand axe $a$ est donnée par:

$$ T^2=\dfrac{4\pi^2 a^3}{GM_{sol}} $$

soit:

In [ ]:
T=np.sqrt(4*np.pi**2*(1000*terre["Demi-GA(km)"])**3/mu0)
print(T)

Soit en terme de jour solaire et de jour sidérale :

In [ ]:
T0=T/(86400);print(T0)
In [ ]:
T0=T/Tsid;print(T0)

En pratique la durée réelle de l'année tropique (interval de temps pour que la lattitude ecliptique du soleil augmente de 360°) est de $365,2421905166$ jours solaires. On distingue également l'année sidérale (interval de temps pour que la soleil revienne au même point par rapport aux étoiles lointaine) qui est de $365,256363004$ jours. En pratique l'année sidérale est légèrement plus longue car en 1 an, les astres éloignés ont légèrement bougé à cause de la precession des équinoxes.

Un calendrier consiste donc à fixer un nombre de jour entier dans une année tropique en s'approchant au plus de la valeur réelle:

  • Le calendrier Julien déffinissait une année de $365.25$ jours, soit 365 jours avec un jour supplémentaire tous les 4 ans. Cela introduisait donc un décalage de $~11$ min par an (soit 11.4 jours en 1500 ans).
  • Le calendrier Grégorien défini une année de $365.2425$ jours, soit 365 jours avec un jour supplémentaire tous les 4 ans sauf pour les nombres de centaines divisibles par 4 (1600,2000,2400...). Cela introduit un décalage de $~26$ secondes par an (soit 11.4 jours en 38076 ans).

Les accords entre l'orbite et le calendrier sont alors fixés par la postion de l'axe de rotation de la terre par rapport à l'éclitptique:

  1. Les équinoxes sont les positions correspondant aux instants d'égale durée du jour et de la nuit (la déclinaison du soleil est nulle).
  2. Les solstices sont les positions correspondant aux durées du jour maximales/minimales en fonction de l'hémisphére (la déclinaison du soleil est maximal).

Il est important de noter que les positions des solstices ne correspondent pas au périhélie et à l'aphélie qui sont respectivement atteint le 03 janvier et le 03 juillet:

Source Wikipédia

Si l'on considère le système de coordonnée écliptique $(\beta,\lambda)$ défini dans le TN1 dans sa version héliocentrique, les coordonnées de la Terre sont donc $(\beta=0,\lambda)$. Par définition $\lambda=0$ lors de l'équinoxe de printemps. On peut alors définir un calendrier à l'aide d'un entier $N$ variant de $1$ (premier jour du calendrier) à $365$ ou $366$. On considère dans toute la suite une année à 365 jours. La longitude écliptique de la terre est donc donnée en degré par:

$$ \lambda=\dfrac{360}{365} (N-80) $$

ou en radian par:

$$ \lambda=\dfrac{2\pi}{365} (N-80) $$

2. Déclinaison du soleil en fonction du temps

En considérant une orbite circulaire de vitesse constante, la hauteur du soleil $h$ dans le système de coordonnée horizontal définis dans le TN1 dépend de la déclinaison $\delta$ du soleil (qui correspond aussi à l'angle entre les rayons du soleil et le plan de l'équateur). La déclinaison du soleil est nulle aux équinoxes et maximum aux solstices. Elle s'exprime simplement:

$$ \delta=\epsilon \cos\left(\dfrac{2\pi}{365}(N+10)\right) $$

La formule exacte tenant compte de l'execentricité de l'orbite est donnée par:

$$ \delta=\arcsin\left(\sin(-\epsilon)\cos\left(\dfrac{360}{365.2425}(N+10)+\dfrac{360}{\pi} e \sin\left(\dfrac{360}{365.2425}(N-3) \right) \right) \right) $$

Avec $\epsilon = 23.44$, l'obliquité de la terre et $e$ l'excentricité de son orbite.

La hauteur $h$ en degrés du soleil dans le ciel au midi solaire s'exprime donc pour une latitude donnée $L$:

$$ h=90-L-\delta $$

Application: Quelle est la hauteur du soleil à Paris lors du midi solaire le jour de l'équinoxe de printemps ?

In [ ]:
h=90-48.8
print(h)q
In [ ]:
def Hsol(N,L):
    long_ec=2*np.pi*(N+10)/365
    delta=23.44*np.cos(long_ec)
    return 90-L-delta

N=np.arange(1,366)
fig=plt.figure(figsize=(16,8))
plt.plot(N,Hsol(N,23.44),label="Lat=23.44°")
plt.plot(N,Hsol(N,48.8),label="Lat=48.8° (Paris)")
plt.plot(N,Hsol(N,60),label="Lat=60°")
plt.xlabel("Jour de l'année",fontsize=fs)
plt.ylabel("Hauteur du soleil",fontsize=fs)
plt.legend(fontsize=fs)

3. Equation du temps

L'équation du temps désigne l'ensemble des corrections qu'il faut apporter à l'azimute du soleil dans le système de coordonée azimutale pour passer de la position moyenne (basée sur une orbite circulaire uniforme) à la position réelle (basée sur une orbite elliptique).

Ces corrections sont souvent données en minutes et permettent alors de passer de l'heure solaire moyenne à l'heure solaire réelle. Ces corrections sont dues d'une part à l'orbite elliptique de la terre mais également à l'inclinaison de l'axe de rotation de la terre.

  1. La correction due à l'orbite elliptique de la terre s'écrit dans sa version simplifiée:
$$C_E=7.678\sin(B+1.374)$$
  1. La correction due à l'oblicité de la terre s'écrit dans sa version simplifiée:
$$ C_O=-9.87\sin(2*C) $$

L'équation du temps s'écrit alors: $$ \Delta z=C_E+C_O $$

Que l'on exprime souvent en minute en considérant que 24H (24*60 minutes) correspondent à $360°$ on a donc 4 minutes par degré:

$$ \Delta T (')= 4\Delta z $$

La version plus précise de l'équation du temps est développé ici.

In [ ]:
def correction_excentrique(N):
    B=2*np.pi*(N-81)/365
    E=7.678*np.sin(B+1.374)
    return E

def correction_oblique(N):
    B=2*np.pi*(N-81)/365
    O=-9.87*np.sin(2*B)
    return O

def equation_du_tps(N):
    return correction_excentrique(N)+correction_oblique(N)

def Real_time_eq(L,N):
    rad=np.pi/180
    # Ecart au 1er Janvier de l'époque J2000:
    d0=(datetime.date(YEAR,1,1)-datetime.date(2000,1,1)).days
    N=N+d0
    M=357.5291+0.98560028*N
    #Excentricité de l'orbite  Terrestre
    e=0.01671
    #Obliquité
    eps=23.43929
    C=(2*e-0.25*e**3)*np.sin(M*rad)+1.25*e**2*np.sin(2*M*rad)+(13/12)*e**3*np.sin(3*M*rad)
    l=280.4665+0.98564736*N+C
    y=np.tan(eps*rad/2)
    R=-y**2*np.sin(2*l*rad)+0.5*y**4*np.sin(4*l*rad)-(y**6/3)*np.sin(6*l*rad)
    T=60*(C+R)
    #Hauteur solaire
    P=Hsol(N,L)
    return 4*T,P
In [ ]:
#Lattitude:
Lat=48.8667

fig=plt.figure(figsize=(16,8))
plt.plot(N,equation_du_tps(N),'--',label='Version Simplifiée')
plt.plot(N,Real_time_eq(Lat,N)[0],'-',label='Version Précise')
plt.xlabel('Jour de l\'année',fontsize=fs)
plt.ylabel('Décallage en min',fontsize=fs)
plt.grid(True)
plt.legend(fontsize=fs)
plt.show()

On peut également représenter l'équation du temps dans les coordonnées polaires en fonction de la longitude ecliptique (nulle le 21 Mars):

In [ ]:
fig=plt.figure(figsize=(14,14))
long_eq=2*np.pi*(N-80)/365
plt.polar(long_eq,equation_du_tps(N),'--',long_eq,Real_time_eq(Lat,N)[0],'-')

4. Analemme

En représentant sur la même figure la déclinaison du soleil en fonction de l'équation du temps, on obtient une Analemme. Cette figure peut également être obtenue en marquant quotidiennement sur le sol la position de l'ombre d'un gnomon à la même heure ou en faisant un cliché de la position d'un astre tous les jours à la même heure.

En pratique on utilise l'analemme solaire pour passer de l'heure solaire vrai à l'heure légale qui s'obtient en ajoutant à l'équation du temps une correction de longitude pour tenir compte de l'heure au méridien local. En effet, l'heure solaire réelle change avec la longitude. On considérant toujours une correction de 4 minutes par degrés on a donc:

$$ H_{legale}=H_{solaire} + \delta_f + 4l + \Delta T $$

Ou $\delta_f$ est la correction due au fuseau horaire (+1 pour l'Europe en hivers, +2 en été), $l$ est la longitude en degrés et $\Delta T$ est l'équation du temps.

In [ ]:
# Longitude (CNAM Paris)
long=2.3559

T,P=Real_time_eq(Lat,N)

# On repère le premier jour du calendrier
T0,P0=Real_time_eq(Lat,1)

# On repère le premier jour de chaque mois
N1=np.cumsum([1,31,28,31,30,31,30,31,31,30,31,30]);
T1,P1=Real_time_eq(Lat,N1)

# Liste des mois
month=['Jan','Fev','Mar','Avr','Mai','Juin','Juil','Aout','Sep','Oct','Nov','Dec']
In [ ]:
def dtr(d):
    return 1/np.tan(d*np.pi/180)

def rtd(r):
    return np.arctan2(1,r)
In [ ]:
from matplotlib.ticker import (MultipleLocator, AutoMinorLocator)
fig=plt.figure(figsize=(16,22))
ax=fig.add_subplot(111)
plt.plot(T,P,'.k')
plt.plot(T0,P0,'om',markersize=12)
plt.plot(T1,P1,'or')
plt.title('Analemme solaire '+str(YEAR)+' (Lattitude: $'+str(Lat)+'$°) \n Correction Longitude cl: '+str(-4*long)+' min \n  \n Heure légale = Heure solaire + 1h(hiver ou 2h été) + cl(min) + $\Delta T$(min)',fontsize=fs)
plt.xlabel('$\Delta T$ (min)',fontsize=fs)
plt.ylabel('Elevation au midi solaire (°)',fontsize=fs)
ax.xaxis.set_minor_locator(AutoMinorLocator())
ax.yaxis.set_minor_locator(MultipleLocator(1))
secaxy = ax.secondary_yaxis('right', functions=(dtr, rtd))
secaxy.set_ylabel(r'$l_o=l/tan(\delta)$',fontsize=fs)
plt.grid(True,which='major',color='k',linestyle='-')
plt.grid(True,which='minor',color='gray',linestyle='--')
for i,n in enumerate(month):
    ax.annotate('%s' % n, xy=(T1[i],P1[i]),xytext=(T1[i],P1[i]), xycoords='data',textcoords='data',fontsize=fs)
    
fig.savefig('Analeme_CNAM_'+str(YEAR)+'.png',bbox_inches='tight',dpi=120)

Applications:

  1. Quelle est l'heure légale le 14 Décembre au midi solaire ?
  2. Le soleil a une hauteur nulle à 07h32 du matin le 20 Novembre, quelle est l'heure légale du lever de soleil ?
  3. Le soleil est à 34° de hauteur dans le ciel au midi solaire, il est 12h37 à la montre. Quel jour sommes-nous ?
  4. Quel heure est-il au midi solaire le 12 juillet Mai ?

Références

[1] INPOP21a planetary ephemerides, IMCCE
[2] INPOP19a planetary ephemerides, IMCCE

In [26]:
from IPython.core.display import HTML
style=open('notebooks.css', "r").read()
HTML(style)
Out[26]:

Retour en haut de la page