Je reviens vers vous une dernière fois pour vous donnez la version finale de ce programme de calculs des phases de la Lune.
Cette ultime version est encore plus précise que les précédentes car elle se base sur de nouveaux calculs et, notamment, ceux de Jean Meeus (Algorithmes astronomiques de Jean Meeus et théorie ELP2000/82
J'ai commenté le plus possible pour essayer de rendre lisible et compréhensible les nombreux calculs autonomes car sans liaisons avec des ressources extérieures (API ou web).
Code : Tout sélectionner
; ----------- Calcul des phases lunaires pour la date du jour ---------------
; Programme de contrôle : Stellarium (v25.2)
; Version finale PhM - septembre 2025 - PB 6.21 (x64)
; ---------------------------------------------------------------------------
EnableExplicit
; ------------------ Déclaration des procédures
Declare.f JourJulien(annee, mois, jour, heure, minute, seconde)
Declare.f AgeLune(JD.f)
Declare.f IlluminationLune(JD.f)
Declare.s PhaseLune(age.f)
Declare.f CalculUTCOffset(annee, mois, jour)
Declare.s SaisonHoraire(offset.f)
; ------------------ Constantes lunaires
Global Lunaison.f = 29.530588853 ; Longueur moyenne du mois synodique
; ------------------ Variables globales temporaires pour IlluminationLune
Global g_elongation.f ; Élongation Lune-Soleil (degrés)
Global g_distance.f ; Distance Terre-Lune (km)
; ------------------ Date et heure UTC avec correction automatique heure d'été/hiver
Global an = Year(Date())
Global mo = Month(Date())
Global jo = Day(Date())
Global UTCOffset.f = CalculUTCOffset(an, mo, jo)
Global heure = Hour(Date()) - UTCOffset
Global minute = Minute(Date())
Global seconde = Second(Date())
; ------------------ Boucle principale du programme ------------------
Define JD.f = JourJulien(an, mo, jo, heure, minute, seconde)
Define age.f = AgeLune(JD)
; Appel d'IlluminationLune → remplit automatiquement g_elongation et g_distance
Define illumination.f = IlluminationLune(JD)
Define phase$ = PhaseLune(age)
; Affichage des résultats
MessageRequester("Observation lunaire du " + Str(jo) + "/" + Str(mo) + "/" + Str(an) + " (" + SaisonHoraire(UTCOffset) + ")",
"Âge lunaire : " + StrF(age, 2) + " jours" + #LF$ +
"Illumination : " + StrF(illumination, 2) + " %" + #LF$ +
"Phase : " + phase$)
; --------------------------------------------------------------------
End
; ------------------ Procédures de calculs
Procedure.f JourJulien(annee, mois, jour, heure, minute, seconde)
Protected a, b, y, m
y = annee
m = mois
If m <= 2
y = annee - 1
m = mois + 12
EndIf
a = y / 100
b = 2 - a + a / 4
Protected JD.f = Int(365.25 * (y + 4716)) + Int(30.6001 * (m + 1)) + jour + b - 1524.5
JD = JD + (heure + minute / 60.0 + seconde / 3600.0) / 24.0
ProcedureReturn JD
EndProcedure
Procedure.f AgeLune(JD.f)
Protected T.f = (JD - 2451545.0) / 36525.0
; --- Calcul de la longitude vraie du Soleil ---
Protected M.f = 357.5291092 + 35999.0502909 * T - 0.0001536 * T * T + Pow(T, 3) / 24490000
M = Mod(M, 360)
Protected Ls.f = 280.46646 + 36000.76983 * T + 0.0003032 * T * T
Ls = Mod(Ls, 360)
Protected lambdaS.f = Ls + 1.914602 * Sin(M * #PI / 180) - 0.004817 * Sin(2 * M * #PI / 180) - 0.019993 * Sin(3 * M * #PI / 180)
lambdaS = Mod(lambdaS, 360)
; --- Calcul de la longitude vraie de la Lune (avec terme de Variation) ---
Protected L.f = 218.3164591 + 481267.88134236 * T - 0.0013268 * T * T + Pow(T, 3) / 538841 - Pow(T, 4) / 65194000
Protected Mp.f = 134.9633964 + 477198.8675055 * T + 0.0087414 * T * T + Pow(T, 3) / 69699 - Pow(T, 4) / 147120000
Mp = Mod(Mp, 360)
Protected D.f = L - Ls : D = Mod(D, 360)
Protected correctionLune.f = 6.289 * Sin(Mp * #PI / 180) - 1.274 * Sin((2*L - lambdaS) * #PI / 180) + 0.658 * Sin(2*(L - lambdaS) * #PI / 180) + 0.214 * Sin(2 * D * #PI / 180)
Protected lambdaL.f = L + correctionLune
lambdaL = Mod(lambdaL, 360)
; --- Élongation apparente ---
Protected diff.f = lambdaL - lambdaS
If diff < 0 : diff + 360 : EndIf
; --- Calcul de l'âge lunaire (jours depuis dernière nouvelle lune) ---
Protected age.f = diff / 360.0 * Lunaison
; --- Correction de +0.15 jour : optimale pour la stabilité et la précision globale ---
; très proche de Stellarium (version 25.2)
age = age + 0.05
ProcedureReturn age.f
EndProcedure
Procedure.f IlluminationLune(JD.f)
; -------------------------------------------------------------
; Algorithmes astronomiques de Jean Meeus et théorie ELP2000/82
; https://archive.org/details/astronomicalalgorithmsjeanmeeus1991/page/n323/mode/2up
; -------------------------------------------------------------
; 357.5291092 Anomalie moyenne du Soleil (M)
; 35999.0502909 Vitesse angulaire de M
; 134.9633964 Anomalie moyenne de la Lune (Mp)
; 477198.8675055 Vitesse angulaire de Mp
; --- Déclaration de TOUTES les variables locales ---
Protected T.f, ecc.f, M_rad.f, E.f, tanHalfE.f, v.f
Protected M.f, Ls.f, lambdaS.f, rS.f, L.f, Mp.f, F.f, beta.f, D.f, rL.f
Protected cosElong.f, i.f, elong.f, k.f, lambdaL.f
; Temps en siècles depuis J2000.0
T = (JD - 2451545.0) / 36525.0
; --- Soleil ---
M = 357.5291092 + 35999.0502909 * T - 0.0001536 * T * T + Pow(T, 3) / 24490000
M = Mod(M, 360)
Ls = 280.46646 + 36000.76983 * T + 0.0003032 * T * T
Ls = Mod(Ls, 360)
lambdaS = Ls + 1.914602 * Sin(M * #PI / 180) - 0.004817 * Sin(2 * M * #PI / 180) - 0.019993 * Sin(3 * M * #PI / 180)
lambdaS = Mod(lambdaS, 360)
rS = 149597870 * (1.00014 - 0.01671 * Cos(M * #PI / 180) - 0.00014 * Cos(2 * M * #PI / 180))
; --- Lune ---
L = 218.3164591 + 481267.88134236 * T - 0.0013268 * T * T + Pow(T, 3) / 538841 - Pow(T, 4) / 65194000
Mp = 134.9633964 + 477198.8675055 * T + 0.0087414 * T * T + Pow(T, 3) / 69699 - Pow(T, 4) / 147120000
Mp = Mod(Mp, 360)
; --- Argument de latitude F ---
F = 93.2720950 + 483202.0175233 * T - 0.0036539 * T * T - Pow(T, 3) / 3526000 + Pow(T, 4) / 863310000
F = Mod(F, 360)
beta = 5.128 * Sin(F * #PI / 180)
; --- [CORRECTION] Calcul de l'anomalie vraie v ---
ecc = 0.0549006
M_rad = Mp * #PI / 180.0
E = M_rad + ecc * Sin(M_rad) + (ecc * ecc / 2.0) * Sin(2.0 * M_rad)
tanHalfE = Tan(E / 2.0)
v = 2.0 * ATan(Sqr((1.0 + ecc) / (1.0 - ecc)) * tanHalfE)
If v < 0 : v + 2.0 * #PI : EndIf
; --- Distance Terre-Lune ---
D = L - Ls
D = Mod(D, 360)
rL = 382700.0 * (1.0 - ecc * Cos(v))
rL = rL * (1.0 - 0.00029 * Cos((2*D - Mp) * #PI / 180) - 0.00019 * Cos(2*D * #PI / 180) + 0.00011 * Cos(Mp * #PI / 180) + 0.00009 * Cos((2*D + Mp) * #PI / 180) - 0.00005 * Cos(2*Mp * #PI / 180))
g_distance = rL
; --- Calcul de lambdaL ---
Protected correctionLune.f = 6.289 * Sin(Mp * #PI / 180) - 1.274 * Sin((2*L - lambdaS) * #PI / 180) + 0.658 * Sin(2*(L - lambdaS) * #PI / 180) + 0.214 * Sin(2 * D * #PI / 180)
lambdaL = L + correctionLune
lambdaL = Mod(lambdaL, 360)
; --- Angle de phase i ---
cosElong = Cos(beta * #PI / 180) * Cos((lambdaL - lambdaS) * #PI / 180)
If cosElong > 1 : cosElong = 1 : ElseIf cosElong < -1 : cosElong = -1 : EndIf
i = 180 - ACos(cosElong) * 180 / #PI
; --- Correction empirique de Meeus ---
i = i + 0.518 * Sin(2 * F * #PI / 180)
If i < 0 : i = 0 : EndIf
If i > 180 : i = 180 : EndIf
; --- Élongation ---
elong = lambdaL - lambdaS
If elong < 0 : elong + 360 : EndIf
g_elongation = elong
; --- Illumination ---
k = (1 + Cos(i * #PI / 180)) / 2 * 100
k = k - 0.55 * (1 - Cos(i * #PI / 180)) * (1 - Cos(beta * #PI / 180)) * 100
; --- Calibration pour alignement précis avec Stellarium ---
k = k + 0.1
ProcedureReturn k
EndProcedure
Procedure.s PhaseLune(age.f)
; ------------------------------------------------------------------------------
; Détermination de la phase lunaire en texte
; Seuils ajustés pour coller aux observations réelles
; ------------------------------------------------------------------------------
Protected phase$
If age.f < 1.8
phase$ = "🌑 Nouvelle lune"
ElseIf age.f < 7.4
phase$ = "🌒 Premier croissant"
ElseIf age.f < 8.4
phase$ = "🌓 Premier quartier"
ElseIf age.f < 14.8
phase$ = "🌔 Gibbeuse croissante"
ElseIf age.f < 15.8
phase$ = "🌕 Pleine lune"
ElseIf age.f < 22.2
phase$ = "🌖 Gibbeuse décroissante"
ElseIf age.f < 23.2
phase$ = "🌗 Dernier quartier"
Else
phase$ = "🌘 Dernier croissant"
EndIf
ProcedureReturn phase$
EndProcedure
Procedure.f CalculUTCOffset(annee, mois, jour)
; ------------------------------------------------------------------------------
; Calcul du décalage UTC (France métropolitaine)
; Correction : recherche du dernier dimanche par DayOfWeek()
; ------------------------------------------------------------------------------
Protected dimancheMars, dimancheOctobre, i ; ← i déclarée ici !
; Dernier dimanche de mars
For i = 31 To 25 Step -1
If DayOfWeek(Date(annee, 3, i, 0, 0, 0)) = 0 ; 0 = dimanche
dimancheMars = i
Break
EndIf
Next
; Dernier dimanche d'octobre
For i = 31 To 25 Step -1
If DayOfWeek(Date(annee, 10, i, 0, 0, 0)) = 0
dimancheOctobre = i
Break
EndIf
Next
; Heure d'été ?
If (mois > 3 And mois < 10) Or (mois = 3 And jour >= dimancheMars) Or (mois = 10 And jour < dimancheOctobre)
ProcedureReturn 2.0 ; UTC+2
Else
ProcedureReturn 1.0 ; UTC+1
EndIf
EndProcedure
Procedure.s SaisonHoraire(offset.f)
; ------------------------------------------------------------------------------
; Libellé de la saison horaire
; ------------------------------------------------------------------------------
If offset = 2.0
ProcedureReturn "Heure d'été : UTC+2"
Else
ProcedureReturn "Heure d'hiver : UTC+1"
EndIf
EndProcedure
