Calcul des phases lunaires pour la date du jour

Programmation d'applications complexes
Avatar de l’utilisateur
PhM
Messages : 137
Inscription : dim. 08/déc./2019 10:50

Calcul des phases lunaires pour la date du jour

Message par PhM »

Bonjour,

Petit programme fonctionnel de calcul des phases lunaire pour la date du jour pour l'hémisphère nord (Paris).

Code : Tout sélectionner

; ----------- Calcul des phases lunaire pour la date du jour ----------
;                      pour l'hémisphère nord (Paris)
;       Site de contrôle : https://starwalk.space/fr/moon-calendar
;          Ph. Mijon - septembre 2025 - PB 6.21 sur Windows 11
; ---------------------------------------------------------------------

; -------------------------- Déclaration des procédures de calcul et données fixes
Declare.f AgeLune(annee, mois, jour)
Declare.s PhaseLune(age.f)
Declare.f IlluminationLune(age.f)

Global RefJulien = 2451550.1     ; Nouvelle lune de référence du 6 janvier 2000 à 14h24
Global Lunaison  = 29.530588853  ; Longueur du mois lunaire 29j 12h 44m 2,9s

; -------------------------- Calcul de la phase lunaire du jour

annee = Year(Date())
mois = Month(Date())
jour = Day(Date())

age.f = AgeLune(annee, mois, jour)
illumination.f = IlluminationLune(age.f)
phase$ = PhaseLune(age.f)

MessageRequester("Phase Lunaire du " + Str(jour) + "/" + Str(mois) + "/" + Str(annee), "La Lune est en phase : " + phase$ + " (" + StrF(illumination.f, 1) + "% illuminée)")

End



; -------------------------- Calcul de l'âge lunaire avec correction

Procedure.f AgeLune(annee, mois, jour)
  
  If mois <= 2
    annee = annee - 1
    mois = mois + 12
  EndIf

  a = annee / 100
  b = 2 - a + a / 4
  JL = Int(365.25 * (annee + 4716)) + Int(30.6001 * (mois + 1)) + jour + b - 1524.5

  age.f = Mod(JL - RefJulien, Lunaison)
  
  If age.f < 0
    age.f = age.f + Lunaison
  EndIf

  ; Correction simple de l’anomalie moyenne
  M = 2 * #PI * age.f / Lunaison
  correction.f = 0.5 * Sin(M) - 0.25 * Sin(2 * M)
  age.f = age.f + correction.f

  ; Facteur de correction pour coller aux données observées
  ; aux coordonnées de Paris (Lat:48.86667 / Lon:2.33333)
  age.f = age.f + 0.56

  ProcedureReturn age.f
EndProcedure

; -------------------------- Détermination de la phase lunaire

Procedure.s PhaseLune(age.f)
  
  If age.f < 1.84566
    phase$ = "🌑 Nouvelle lune"
  ElseIf age.f < 5.53699
    phase$ = "🌒 Premier croissant"
  ElseIf age.f < 9.22831
    phase$ = "🌓 Premier quartier"
  ElseIf age.f < 12.91963
    phase$ = "🌔 Gibbeuse croissante"
  ElseIf age.f < 16.61096
    phase$ = "🌕 Pleine lune"
  ElseIf age.f < 20.30228
    phase$ = "🌖 Gibbeuse décroissante"
  ElseIf age.f < 23.99361
    phase$ = "🌗 Dernier quartier"
  ElseIf age.f < 27.68493
    phase$ = "🌘 Dernier croissant"
  Else
    phase$ = "🌑 Nouvelle lune"
  EndIf

  ProcedureReturn phase$
EndProcedure

; -------------------------- Calcul du pourcentage d'illumination

Procedure.f IlluminationLune(age.f)

  illumination.f = (1 - Cos(2 * #PI * age.f / Lunaison)) / 2 * 100
  
  ProcedureReturn illumination.f
EndProcedure


Dernière modification par PhM le mar. 23/sept./2025 9:46, modifié 1 fois.
Avatar de l’utilisateur
Jacobus
Messages : 1577
Inscription : mar. 06/avr./2004 10:35
Contact :

Re: Calcul des phases lunaire pour la date du jour

Message par Jacobus »

Bonjour,
Sympa ton application :)
Merci pour le partage.
Quand tous les glands seront tombés, les feuilles dispersées, la vigueur retombée... Dans la morne solitude, ancré au coeur de ses racines, c'est de sa force maturité qu'il renaîtra en pleine magnificence...Jacobus.
Avatar de l’utilisateur
MLD
Messages : 1125
Inscription : jeu. 05/févr./2009 17:58
Localisation : Bretagne

Re: Calcul des phases lunaire pour la date du jour

Message par MLD »

Sympa
Merci du partage
Avatar de l’utilisateur
venom
Messages : 3153
Inscription : jeu. 29/juil./2004 16:33
Localisation : Klyntar
Contact :

Re: Calcul des phases lunaire pour la date du jour

Message par venom »

Bonjour.

Je n'en ai pas d'utilité, mais merci pour le partage c'est cool, ça peut toujours être utile. 8)






@++
Windows 10 x64, PureBasic 5.73 x86 & x64
GPU : radeon HD6370M, CPU : p6200 2.13Ghz
Avatar de l’utilisateur
PhM
Messages : 137
Inscription : dim. 08/déc./2019 10:50

Re: Calcul des phases lunaire pour la date du jour

Message par PhM »

Amélioration du premier programme qui gagne en précision dans cette version :

Code : Tout sélectionner


; ----------- Calcul des phases lunaires pour la date du jour ----------
;                      pour l'hémisphère nord (Paris)
;    Version améliorée PhM - septembre 2025 - PB 6.21 sur Windows 11
; ----------------------------------------------------------------------

Declare.f AgeLune(annee, mois, jour)
Declare.s PhaseLune(age.f)
Declare.f IlluminationLune(age.f)

Global RefJulien = 2451550.1     ; Nouvelle lune de référence du 6 janvier 2000 à 14h24 UTC
Global Lunaison  = 29.530588853  ; Longueur moyenne du mois lunaire

; -------------------------- Date du jour
annee = Year(Date())
mois = Month(Date())
jour = Day(Date())

age.f = AgeLune(annee, mois, jour)
illumination.f = IlluminationLune(age.f)
phase$ = PhaseLune(age.f)

MessageRequester("Phase Lunaire du " + Str(jour) + "/" + Str(mois) + "/" + Str(annee),
                 "Âge lunaire : " + StrF(age.f, 1) + " jours" + #LF$ +
                 "Illumination : " + StrF(illumination.f, 1) + "%" + #LF$ +
                 "Phase : " + phase$)

End

; -------------------------- Calcul de l'âge lunaire avec correction
Procedure.f AgeLune(annee, mois, jour)
  If mois <= 2
    annee - 1
    mois + 12
  EndIf

  a = annee / 100
  b = 2 - a + a / 4
  JL = Int(365.25 * (annee + 4716)) + Int(30.6001 * (mois + 1)) + jour + b - 1524.5

  age.f = Mod(JL - RefJulien, Lunaison)
  If age.f < 0
    age.f = age.f + Lunaison
  EndIf

  ; Correction de l’anomalie moyenne
  M = 2 * #PI * age.f / Lunaison
  correction.f = 0.5 * Sin(M) - 0.25 * Sin(2 * M)
  age.f = age.f + correction.f

  ; Ajustement empirique pour Paris
  age.f = age.f + 1.6 * Sin(2 * M)
  age.f = Mod(age.f, Lunaison)

  ProcedureReturn age.f
EndProcedure

; -------------------------- Calcul du pourcentage d'illumination
Procedure.f IlluminationLune(age.f)
  Protected correction.f = 0.953
  ProcedureReturn (1 - Cos(2 * #PI * age.f / Lunaison)) / 2 * 100 * correction
EndProcedure

; -------------------------- Détermination de la phase lunaire
Procedure.s PhaseLune(age.f)
  illumination.f = IlluminationLune(age.f)

  If illumination.f < 1
    phase$ = "🌑 Nouvelle lune"
  ElseIf illumination.f < 25 And age.f < 7.5
    phase$ = "🌒 Premier croissant"
  ElseIf illumination.f >= 25 And illumination.f < 50 And age.f < 9
    phase$ = "🌓 Premier quartier"
  ElseIf illumination.f >= 50 And illumination.f < 99 And age.f < 14.8
    phase$ = "🌔 Gibbeuse croissante"
  ElseIf illumination.f >= 99
    phase$ = "🌕 Pleine lune"
  ElseIf illumination.f >= 50 And age.f > 15
    phase$ = "🌖 Gibbeuse décroissante"
  ElseIf illumination.f >= 25 And age.f > 21
    phase$ = "🌗 Dernier quartier"
  ElseIf illumination.f < 25 And age.f > 22
    phase$ = "🌘 Dernier croissant"
  Else
    phase$ = "🌑 Nouvelle lune"
  EndIf

  ProcedureReturn phase$
EndProcedure


Avatar de l’utilisateur
SPH
Messages : 4963
Inscription : mer. 09/nov./2005 9:53

Re: Calcul des phases lunaire pour la date du jour

Message par SPH »

Merci pour le code. J'admire ton code bien structuré. Bravo :)

!i!i!i!i!i!i!i!i!i!
!i!i!i!i!i!i!
!i!i!i!
//// Informations ////
Intel Core i7 4770 64 bits - GTX 650 Ti
Version de PB : 6.12LTS- 64 bits
Avatar de l’utilisateur
PhM
Messages : 137
Inscription : dim. 08/déc./2019 10:50

Re: Calcul des phases lunaire pour la date du jour

Message par PhM »

Ultime modifications avec l'introduction du lieu d'observation exact (lat/lon) ainsi que la correction automatique pour l'heure d'été/hiver et l'indication de l'altitude et de l'azimut de la Lune.

Les résultats deviennent plus précis mais plus complexes à calculer.

Code : Tout sélectionner


; ----------- Calcul des phases lunaires pour la date du jour ----------
;                   Programme de contrôle : Stellarium
;    Version améliorée PhM - septembre 2025 - PB 6.21 pour Windows 11
;         🌕 par l'introduction du lieu d'observation (lat/lon)
;         🌕 ainsi que la correction automatique pour l'heure d'été/hiver
;         🌕 indication de l'altitude et de l'azimut de la Lune
; ----------------------------------------------------------------------

EnableExplicit

; ------------------ Déclaration des procédures

Declare.f JourJulien(annee, mois, jour, heure, minute, seconde)
Declare.f AgeLune(JD.f)
Declare.f IlluminationLune(age.f)
Declare.s PhaseLune(age.f)
Declare.s PositionLune(JD.f)
Declare.f CalculUTCOffset(annee, mois, jour)
Declare.s SaisonHoraire(offset.f)

; ------------------ Coordonnées de l'observation

Global Latitude.f  = 43.64     ; Latitude d'observation à la précision souhaitée
Global Longitude.f = 3.96     ; Longitude d'observation à la précision souhaitée

; ------------------ Constantes lunaires

Global RefJulien = 2451550.1        ; Nouvelle lune de référence du 6 janvier 2000 à 14h24
Global Lunaison  = 29.530588853     ; Longueur du mois lunaire 29j 12h 44m 2,9s

; ------------------ Date et heure UTC corrigée 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)
Define illumination.f = IlluminationLune(age)
Define phase$ = PhaseLune(age)
Define position$ = PositionLune(JD)

MessageRequester("Observation lunaire du " + Str(jo) + "/" + Str(mo) + "/" + Str(an) + " (" + SaisonHoraire(UTCOffset) + ")",
                 "Âge lunaire : " + StrF(age, 1) + " jours" + #LF$ +
                 "Illumination : " + StrF(illumination, 1) + "%" + #LF$ +
                 "Phase : " + phase$ + #LF$ + #LF$ +
                 "Position apparente de la Lune :" + #LF$ + position$)

End


; ------------------ Procédures de calculs

Procedure.f JourJulien(annee, mois, jour, heure, minute, seconde)

  If mois <= 2
    annee - 1
    mois + 12
  EndIf
  
  Protected A = annee / 100
  Protected B = 2 - A + A / 4
  Protected JD.f = Int(365.25 * (annee + 4716)) + Int(30.6001 * (mois + 1)) + jour + B - 1524.5
  JD + (heure + minute / 60.0 + seconde / 3600.0) / 24.0
  
  ProcedureReturn JD
  
EndProcedure

Procedure.f AgeLune(JD.f)
  
  Protected age.f = Mod(JD - RefJulien, Lunaison)
  
  If age.f < 0
    age.f + Lunaison
  EndIf
  
  Protected M = 2 * #PI * age.f / Lunaison
  Protected correction.f = 0.95 * Sin(M) - 0.25 * Sin(2 * M)
  age.f + correction.f
  age.f = Mod(age.f, Lunaison)
  
  ProcedureReturn age.f
  
EndProcedure

Procedure.f IlluminationLune(age.f)
  
  Protected correction.f = 0.855
  ProcedureReturn (1 - Cos(2 * #PI * age.f / Lunaison)) / 2 * 100 * correction
  
EndProcedure

Procedure.s PhaseLune(age.f)
  
  Protected illumination.f = IlluminationLune(age.f)
  Protected phase$
  
  If illumination.f < 1
    phase$ = "🌑 Nouvelle lune"
  ElseIf illumination.f < 25 And age.f < 7.5
    phase$ = "🌒 Premier croissant"
  ElseIf illumination.f >= 25 And illumination.f < 50 And age.f < 9
    phase$ = "🌓 Premier quartier"
  ElseIf illumination.f >= 50 And illumination.f < 99 And age.f < 14.8
    phase$ = "🌔 Gibbeuse croissante"
  ElseIf illumination.f >= 99
    phase$ = "🌕 Pleine lune"
  ElseIf illumination.f >= 50 And age.f > 15
    phase$ = "🌖 Gibbeuse décroissante"
  ElseIf illumination.f >= 25 And age.f > 21
    phase$ = "🌗 Dernier quartier"
  ElseIf illumination.f < 25 And age.f > 22
    phase$ = "🌘 Dernier croissant"
  Else
    phase$ = "🌑 Nouvelle lune"
  EndIf
  
  ProcedureReturn phase$
  
EndProcedure

Procedure.s PositionLune(JD.f)
  
  Protected T.f = (JD - 2451545.0) / 36525.0
  Protected L.f = Mod(218.316 + 13.176396 * (JD - 2451545.0), 360)
  Protected M.f = Mod(134.963 + 13.064993 * (JD - 2451545.0), 360)
  Protected F.f = Mod(93.272 + 13.229350 * (JD - 2451545.0), 360)
  Protected lambda.f = L + 6.289 * Sin(Radian(M))
  Protected beta.f = 5.128 * Sin(Radian(F))
  Protected epsilon.f = 23.439 - 0.0000004 * T
  Protected alpha.f = Degree(ATan2(Cos(Radian(epsilon)) * Sin(Radian(lambda)), Cos(Radian(lambda))))
  Protected delta.f = Degree(ASin(Sin(Radian(epsilon)) * Sin(Radian(lambda)) * Cos(Radian(beta)) + Sin(Radian(beta)) * Cos(Radian(epsilon))))
  Protected GMST.f = Mod(280.46061837 + 360.98564736629 * (JD - 2451545.0), 360)
  Protected LST.f = Mod(GMST + Longitude, 360)
  Protected HA.f = Mod(LST - alpha, 360)
  Protected altitude.f = Degree(ASin(Sin(Radian(Latitude)) * Sin(Radian(delta)) + Cos(Radian(Latitude)) * Cos(Radian(delta)) * Cos(Radian(HA))))
  Protected azimut.f = Degree(ATan2(-Sin(Radian(HA)), Tan(Radian(delta)) * Cos(Radian(Latitude)) - Sin(Radian(Latitude)) * Cos(Radian(HA))))
  
  ProcedureReturn "Altitude : " + StrF(altitude, 2) + "°" + #LF$ + "Azimut : " + StrF(azimut, 2) + "°"
  
EndProcedure

Procedure.f CalculUTCOffset(annee, mois, jour)
  
  ; ------------------ Détection automatique de l'heure d'été/hiver
  
  Protected jourETE, jourHIVER
  Protected UTCOffset.f

  ; Dernier dimanche de mars
  jourETE = 31 - ((Date(annee, 3, 31, 0, 0, 0) / 86400 + 3) % 7)

  ; Dernier dimanche d'octobre
  jourHIVER = 31 - ((Date(annee, 10, 31, 0, 0, 0) / 86400 + 3) % 7)

  If (mois > 3 And mois < 10) Or (mois = 3 And jour >= jourETE) Or (mois = 10 And jour < jourHIVER)
    UTCOffset = 2.0 ; Heure d'été en France
  Else
    UTCOffset = 1.0 ; Heure d'hiver en France
  EndIf

  ProcedureReturn UTCOffset

EndProcedure

Procedure.s SaisonHoraire(offset.f)
  
    ; Pour l'affichage
  If offset = 2.0
    ProcedureReturn "Heure d'été : UTC+2"
  Else
    ProcedureReturn "Heure d'hiver : UTC+1"
  EndIf

EndProcedure


Avatar de l’utilisateur
PhM
Messages : 137
Inscription : dim. 08/déc./2019 10:50

Re: Calcul des phases lunaire pour la date du jour

Message par PhM »

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 https://archive.org/details/astronomica ... 3/mode/2up).

Le contrôle et l'étalonnage ont été fait avec le programme Stellaium (https://stellarium.org/fr/).

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


Avatar de l’utilisateur
Ar-S
Messages : 9543
Inscription : dim. 09/oct./2005 16:51
Contact :

Re: Calcul des phases lunaire pour la date du jour

Message par Ar-S »

Merci pour ce partage. Je testerai à l'occasion.
~~~~Règles du forum ~~~~
⋅.˳˳.⋅ॱ˙˙ॱ⋅.˳Ar-S ˳.⋅ॱ˙˙ॱ⋅.˳˳.⋅
W11x64 PB 6.x
Section HORS SUJET : ICI
LDV MULTIMEDIA : Dépannage informatique & mes Logiciels PB
UPLOAD D'IMAGES : Uploader des images de vos logiciels
Avatar de l’utilisateur
Jacobus
Messages : 1577
Inscription : mar. 06/avr./2004 10:35
Contact :

Re: Calcul des phases lunaire pour la date du jour

Message par Jacobus »

Il serait intéressant de combiner ce programme avec celui de LSI :
viewtopic.php?t=12523

Code : Tout sélectionner

; https://www.purebasic.fr/french/viewtopic.php?t=12523
; par LSI
Enumeration 0
		#MoonPhase_None = -1
		#MoonPhase_NewMoon
		#MoonPhase_FirstQuarter
		#MoonPhase_FullMoon
		#MoonPhase_LastQuarter
EndEnumeration
Structure MoonPhase_Structure
		Phase.b
		Date.l
		An.l
		Mois.b
		Jour.b
		Heure.b
		Minute.b
EndStructure
Global NewList MoonPhase.MoonPhase_Structure()

Procedure MoonPhase_Calculation(Year, Month)
		Static MoonPhase_Calculation_Year.l, MoonPhase_Calculation_Month.b
		Protected.l Date, i, ii, Heure, Minute, Mois, An, Jour, JJ
		Protected.d K, T, T2, T3, J, M, MP, F
		
		Debug #PB_Compiler_Procedure
		
		If Year <> MoonPhase_Calculation_Year Or Month <> MoonPhase_Calculation_Month
				MoonPhase_Calculation_Year = Year
				MoonPhase_Calculation_Month = Month
				
				ClearList(MoonPhase())
				
				Repeat
						GetSystemTime_(DateUTC.SYSTEMTIME)
						GetLocalTime_(DateLocale.SYSTEMTIME)
				Until DateUTC\wSecond = DateLocale\wSecond ; Garantit que la lecture de la date s'est effectuée sur la même seconde
				DateUTC_Seconde.q = Date(DateUTC\wYear, DateUTC\wMonth, DateUTC\wDay, DateUTC\wHour, DateUTC\wMinute, DateUTC\wSecond)
				DateLocale_Seconde.q = Date(DateLocale\wYear, DateLocale\wMonth, DateLocale\wDay, DateLocale\wHour, DateLocale\wMinute, DateLocale\wSecond)
				
				DecalageHoraire_Seconde.q = DateLocale_Seconde - DateUTC_Seconde
				
				CompilerIf #PB_Compiler_Debugger
						DecalageHoraire = DecalageHoraire_Seconde / 3600
						If DecalageHoraire > 0
								Debug "Votre fuseau horaire : GMT+" + Str(DecalageHoraire)
						Else
								Debug "Votre fuseau horaire : GMT" + Str(DecalageHoraire)
						EndIf
				CompilerEndIf
				
				Date = Date(Year, Month, 1, 0, 0, 0)
				Date = AddDate(Date, #PB_Date_Month, -1)
				K.d = Year(Date)
				Select Month(Date)
						Case 1
								K + 0.041
						Case 2
								K + 0.126
						Case 3
								K + 0.203
						Case 4
								K + 0.288
						Case 5
								K + 0.370
						Case 6
								K + 0.455
						Case 7
								K + 0.537
						Case 8
								K + 0.622
						Case 9
								K + 0.707
						Case 10
								K + 0.789
						Case 11
								K + 0.874
						Case 12
								K + 0.956
				EndSelect
				K = (K - 1900) * 12.3685
				K = Int(K) - 0.25
				If K < 0
						K - 1
				EndIf
				
				For ii = 0 To 11
						
						K + 0.25
						T.d = K / 1236.85
						T2.d = T * T
						T3.d = T * T2
						J.d = 2415020.75933 + 29.5305888531 * K + 0.0001337 * T2 - 0.000000150 * T3 + 0.00033 * Sin((166.56 + 132.87 * T - 0.009 * T2) * #PI / 180)
						M.d = (359.2242 + 29.10535608 * K - 0.0000333 * T2 - 0.00000347 * T3) * #PI / 180
						M = M - Int(M / (2 * #PI)) * (2 * #PI)
						MP.d = (306.0253 + 385.81691806 * K + 0.0107306 * T2 + 0.00001236 * T3) * #PI / 180
						MP = MP - Int(MP / (2 * #PI)) * (2 * #PI)
						F.d = (21.2964 + 390.67050646 * K - 0.0016528 * T2 - 0.00000239 * T3) * #PI / 180
						F = F - Int(F / (2 * #PI)) * (2 * #PI)
						
						i = ii % 4
						If i = 0 Or i = 2
								J + (0.1734 - 0.000393 * T) * Sin(M)
								J + 0.0021 * Sin(2 * M) - 0.4068 * Sin(MP)
								J + 0.0161 * Sin(2 * MP) - 0.0004 * Sin(3 * MP)
								J + 0.0104 * Sin(2 * F) - 0.0051 * Sin(M + MP)
								J - 0.0074 * Sin(M - MP) + 0.0004 * Sin(2 * F + M)
								J - 0.0004 * Sin(2 * F - M) - 0.0006 * Sin(2 * F + MP)
								J + 0.001 * Sin(2 * F - MP) + 0.0005 * Sin(M + 2 * MP)
						Else
								J + (0.1721 - 0.0004 * T) * Sin(M)
								J + 0.0021 * Sin(2 * M) - 0.6280 * Sin(MP)
								J + 0.0089 * Sin(2 * MP) - 0.0004 * Sin(3 * MP)
								J + 0.0079 * Sin(2 * F) - 0.0119 * Sin(M + MP)
								J - 0.0047 * Sin(M - MP) + 0.0003 * Sin(2 * F + M)
								J - 0.0004 * Sin(2 * F - M) - 0.0006 * Sin(2 * F + MP)
								J + 0.0021 * Sin(2 * F - MP) + 0.0003 * Sin(M + 2 * MP)
								J + 0.0004 * Sin(M - 2 * MP) - 0.0003 * Sin(2 * M + MP)
								If i = 1
										J + 0.0028 - 0.0004 * Cos(M)
										J + 0.0003 * Cos(MP)
								Else
										J - 0.0028 + 0.0004 * Cos(M)
										J - 0.0003 * Cos(MP)
								EndIf
						EndIf
						
						J + 0.5
						JJ = Int(J)
						If JJ >= 2299160.5
								Alpha.d = Int((JJ - 1867216.25) / 36524.25)
								JJ = JJ + 1 + Alpha - Int(Alpha / 4)
						EndIf
						JJ + 1524
						Calcul_An = Int((JJ - 122.1) / 365.25)
						Calcul_Jour = Int(Calcul_An * 365.25)
						Calcul_Mois = Int((JJ - Calcul_Jour) / 30.6001)
						Jour = Int(JJ - Calcul_Jour - Int(Calcul_Mois * 30.6001))
						If Calcul_Mois < 13.5
								Mois = Int(Calcul_Mois - 1)
						Else
								Mois = Int(Calcul_Mois - 13)
						EndIf
						If Mois >= 3
								An = Int(Calcul_An - 4716)
						Else
								An = Int(Calcul_An - 4715)
						EndIf
						J - Int(J)
						Heure = Int(J * 24)
						Minute = Int((J - Heure / 24) * 1440)
						
						CompilerIf #PB_Compiler_Debugger
								Select i
										Case #MoonPhase_NewMoon
												Debug "Nouvelle lune"
										Case #MoonPhase_FirstQuarter
												Debug "Premier quartier"
										Case #MoonPhase_FullMoon
												Debug "Pleine lune"
										Case #MoonPhase_LastQuarter
												Debug "Dernier quartier"
								EndSelect
								Debug Str(Jour) + "/" + RSet(Str(Mois), 2, "0") + "/" + Str(An) + " à " + Str(Heure) + ":" + RSet(Str(Minute), 2, "0") + " (UTC)"
						CompilerEndIf
						
						AddElement(MoonPhase())
						MoonPhase()\Phase = i
						MoonPhase()\Date = Date(An, Mois, Jour, 0, 0, 0) + (Heure * 3600 + Minute * 60) + DecalageHoraire_Seconde
						MoonPhase()\An = Year(MoonPhase()\Date)
						MoonPhase()\Mois = Month(MoonPhase()\Date)
						MoonPhase()\Jour = Day(MoonPhase()\Date)
						MoonPhase()\Heure = Hour(MoonPhase()\Date)
						MoonPhase()\Minute = Minute(MoonPhase()\Date)
						
						Debug Str(MoonPhase()\Jour) + "/" + RSet(Str(MoonPhase()\Mois), 2, "0") + "/" + Str(MoonPhase()\An) + " à " + Str(MoonPhase()\Heure) + ":" + RSet(Str(MoonPhase()\Minute), 2, "0")
						
				Next
				
		EndIf
EndProcedure

Procedure NextMoonPhase(Year = 0, Month = 0, Day = 0) ; Get the next moon phase after the date in parameters or last result
		If Year And Month And Day
				MoonPhase_Calculation(Year, Month)
				Date = Date(Year, Month, Day, 0, 0, 0)
		ElseIf ListIndex(MoonPhase())
				Year = MoonPhase()\An
				Month = MoonPhase()\Mois
				Day = MoonPhase()\Jour
				MoonPhase_Calculation(Year, Month)
				Date = Date(Year, Month, Day, 0, 0, 0) + (24 * 60 * 60)
		EndIf
		If Date
				ForEach MoonPhase()
						If MoonPhase()\Date > Date
								Date = MoonPhase()\Date
								Break
						EndIf
				Next
		EndIf
		Debug #PB_Compiler_Procedure
		Debug Date
		ProcedureReturn Date
EndProcedure

Procedure GetMoonPhase(Year = 0, Month = 0, Day = 0) ; Get Moon phase of last result or of the date in parameters
		If Year And Month And Day
				MoonPhase_Calculation(Year, Month)
				Phase = #MoonPhase_None
				ForEach MoonPhase()
						If Year = MoonPhase()\An And Month = MoonPhase()\Mois And Day = MoonPhase()\Jour
								Phase = MoonPhase()\Phase
								Break
						EndIf
				Next
		ElseIf ListIndex(MoonPhase())
						Phase = MoonPhase()\Phase
				EndIf
		CompilerIf #PB_Compiler_Debugger
				Debug #PB_Compiler_Procedure
				Select Phase
						Case #MoonPhase_None
								Debug "Ce jour n'est pas un état spécifique de la Lune"
						Case #MoonPhase_NewMoon
								Debug "Nouvelle lune"
						Case #MoonPhase_FirstQuarter
								Debug "Premier quartier"
						Case #MoonPhase_FullMoon
								Debug "Pleine lune"
						Case #MoonPhase_LastQuarter
								Debug "Dernier quartier"
				EndSelect
		CompilerEndIf
		ProcedureReturn Phase
EndProcedure

Procedure GetMoonPhaseYear() ; Get year of moon phase of last result
		An = #MoonPhase_None
		If ListIndex(MoonPhase())
				An = MoonPhase()\An
		EndIf
		Debug #PB_Compiler_Procedure
		Debug An
		ProcedureReturn An
EndProcedure

Procedure GetMoonPhaseMonth() ; Get month of moon phase of last result
		Mois = #MoonPhase_None
		If ListIndex(MoonPhase())
				Mois = MoonPhase()\Mois
		EndIf
		Debug #PB_Compiler_Procedure
		Debug Mois
		ProcedureReturn Mois
EndProcedure

Procedure GetMoonPhaseDay() ; Get day of moon phase of last result
		Jour = #MoonPhase_None
		If ListIndex(MoonPhase())
				Jour = MoonPhase()\Jour
		EndIf
		Debug #PB_Compiler_Procedure
		Debug Jour
		ProcedureReturn Jour
EndProcedure

Procedure GetMoonPhaseHour(Year = 0, Month = 0, Day = 0) ; Get hour of Moon phase of last result or of the date in parameters
		If Year And Month And Day
				MoonPhase_Calculation(Year, Month)
				Heure = #MoonPhase_None
				ForEach MoonPhase()
						If Year = MoonPhase()\An And Month = MoonPhase()\Mois And Day = MoonPhase()\Jour
								Heure = MoonPhase()\Heure
								Break
						EndIf
				Next
		ElseIf ListIndex(MoonPhase())
						Heure = MoonPhase()\Heure
				EndIf
		Debug #PB_Compiler_Procedure
		Debug Heure
		ProcedureReturn Heure
EndProcedure

Procedure GetMoonPhaseMinute(Year = 0, Month = 0, Day = 0) ; Get minute of Moon phase of last result or of the date in parameters
		If Year And Month And Day
				MoonPhase_Calculation(Year, Month)
				Minute = #MoonPhase_None
				ForEach MoonPhase()
						If Year = MoonPhase()\An And Month = MoonPhase()\Mois And Day = MoonPhase()\Jour
								Minute = MoonPhase()\Minute
								Break
						EndIf
				Next
		ElseIf ListIndex(MoonPhase())
						Minute = MoonPhase()\Minute
				EndIf
		Debug #PB_Compiler_Procedure
		Debug Minute
		ProcedureReturn Minute
EndProcedure



; Test du programme

Enumeration
		#ListeMois
		#ListeAn
		#An
		#Mois
EndEnumeration

Procedure RemplirListe(An, Mois)
		
		; Pour tous les jours du mois
		ClearGadgetItems(#ListeMois)
		Date = Date(An, Mois, 1, 0, 0, 0) ; Premier jour du mois
		While Month(Date) = Mois
				; Jour
				Jour = Day(Date)
				Jour_Texte.s = RSet(Str(Jour), 2, "0") + "/" + RSet(Str(Mois), 2, "0") + "/" + Str(An)
				; Calcul de la phase de la Lune
				Phase = GetMoonPhase(An, Mois, Jour)
				Select Phase
						Case #MoonPhase_None
								Phase_Texte.s = ""
						Case #MoonPhase_NewMoon
								Phase_Texte.s = "Nouvelle lune"
						Case #MoonPhase_FirstQuarter
								Phase_Texte.s = "Premier quartier"
						Case #MoonPhase_FullMoon
								Phase_Texte.s = "Pleine lune"
						Case #MoonPhase_LastQuarter
								Phase_Texte.s = "Dernier quartier"
				EndSelect
				; Calcul de l'heure
				If Phase <> #MoonPhase_None
						Heure = GetMoonPhaseHour(An, Mois, Jour)
						Minute = GetMoonPhaseMinute(An, Mois, Jour)
						Heure_Texte.s = RSet(Str(Heure), 2, "0") + ":" + RSet(Str(Minute), 2, "0")
				Else
						Heure_Texte.s = ""
				EndIf
				; Affiche dans la liste
				AddGadgetItem(#ListeMois, -1, Jour_Texte + Chr(10) + Heure_Texte + Chr(10) + Phase_Texte)
				; Ajoute un jour
				Date = AddDate(Date, #PB_Date_Day, 1)
		Wend
		
		
		; Pour l'année complète
		ClearGadgetItems(#ListeAn)
		Date = Date(An, 1, 1, 0, 0, 0) ; Premier jour de l'année
		Date = NextMoonPhase(An, Mois, 1)
		While Year(Date) = An
				; Jour
				Mois = GetMoonPhaseMonth()
				Jour = GetMoonPhaseDay()
				Jour_Texte.s = RSet(Str(Jour), 2, "0") + "/" + RSet(Str(Mois), 2, "0") + "/" + Str(An)
				; Calcul de la phase de la Lune
				Phase = GetMoonPhase()
				Select Phase
						Case #MoonPhase_None
								Phase_Texte.s = ""
						Case #MoonPhase_NewMoon
								Phase_Texte.s = "Nouvelle lune"
						Case #MoonPhase_FirstQuarter
								Phase_Texte.s = "Premier quartier"
						Case #MoonPhase_FullMoon
								Phase_Texte.s = "Pleine lune"
						Case #MoonPhase_LastQuarter
								Phase_Texte.s = "Dernier quartier"
				EndSelect
				; Calcul de l'heure
				Heure = GetMoonPhaseHour()
				Minute = GetMoonPhaseMinute()
				Heure_Texte.s = RSet(Str(Heure), 2, "0") + ":" + RSet(Str(Minute), 2, "0")
				; Affiche dans la liste
				AddGadgetItem(#ListeAn, -1, Jour_Texte + Chr(10) + Heure_Texte + Chr(10) + Phase_Texte)
				; Phase suivante
				Date = NextMoonPhase()
		Wend
		
EndProcedure

; Lecture de la date actuelle
Date = Date()
An = Year(Date)
Mois = Month(Date)

If OpenWindow(0, 0, 0, 600, 600, "Phase de la Lune", #PB_Window_SystemMenu | #PB_Window_MinimizeGadget | #PB_Window_ScreenCentered) = 0
		End
EndIf
ListIconGadget(#ListeMois, 0, 0, 300, 560, "Jour du mois", 100, #PB_ListIcon_FullRowSelect)
AddGadgetColumn(#ListeMois, 1, "Heure", 50)
AddGadgetColumn(#ListeMois, 2, "Phase de la lune", 120)
ListIconGadget(#ListeAn, 300, 0, 300, 560, "Jour de l'année", 100, #PB_ListIcon_FullRowSelect)
AddGadgetColumn(#ListeAn, 1, "Heure", 50)
		AddGadgetColumn(#ListeAn, 2, "Phase de la lune", 120)
TextGadget(#PB_Any, 0, 560, 300, 16, "Année", #PB_Text_Center)
TextGadget(#PB_Any, 200, 560, 300, 16, "Mois", #PB_Text_Center)
StringGadget(#An, 0, 576, 300, 24, Str(An), #PB_String_Numeric | #ES_CENTER)
StringGadget(#Mois, 300, 576, 300, 24, Str(Mois), #PB_String_Numeric | #ES_CENTER)

RemplirListe(An, Mois)

Repeat
		Event = WaitWindowEvent()
		
		Select Event
				Case #PB_Event_Gadget
						Select EventGadget()
								Case #An, #Mois
										If EventType() = #PB_EventType_Change
												
												An = Val(GetGadgetText(#An))
												Mois = Val(GetGadgetText(#Mois))
												If An And Mois
														RemplirListe(An, Mois)
												EndIf
												
										EndIf
						EndSelect
		EndSelect
		
Until Event = #PB_Event_CloseWindow
Quand tous les glands seront tombés, les feuilles dispersées, la vigueur retombée... Dans la morne solitude, ancré au coeur de ses racines, c'est de sa force maturité qu'il renaîtra en pleine magnificence...Jacobus.
Avatar de l’utilisateur
PhM
Messages : 137
Inscription : dim. 08/déc./2019 10:50

Re: Calcul des phases lunaire pour la date du jour

Message par PhM »

Bonjour,

Oui, j'avais vu ce programme. Malheureusement, il n'est pas très précis, notamment car les résultats ne tiennent pas compte de l'heure.

Pour ma part, bien que mon dernier code était annoncé comme une " version finale", je continu à améliorer la précision. Le code continu d'évoluer toujours de manière autonome. C'est à dire, sans connexion web ou recours à des API (tables Nasa, JPL, etc.).

Pour la beauté du code, comme on dit...

Bonne journée.
Avatar de l’utilisateur
PhM
Messages : 137
Inscription : dim. 08/déc./2019 10:50

Re: Calcul des phases lunaires pour la date du jour

Message par PhM »

Jacobus,

Pour illustrer l'importance de l'heure dans les calculs de phases lunaires, j'ai modifié le programme pour qu'il calcule sur un mois entier (au choix) et sur 3 horaires à choisir par jour (au choix), un éphéméride mensuel.

Tu verras, par exemple avec les dates et heures de ce programme que le 22/10/205, la phase n'est pas la même à 6h00 du matin qu'à 22h00. Cela est normal en tenant compte de l'âge de la Lune et de son illumination qui évoluent sans cesse dans le temps.

La précision actuelle de ces deux programmes est de : ±0.3 jour / ±1.5% illum — transitions validées avec Stellarium

Je vous livre cette version du programme mensuelle maintenant et, à la suite, la toute dernière version du programme journalier, objet de mon premier post :

Code : Tout sélectionner


; ----------- Calcul des phases lunaires pour un mois entier (3x/jour) ---------------
;                Programme de contrôle : Stellarium (v25.2)
;           Version mensuelle 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, illumination.f)
Declare.f CalculUTCOffset(annee, mois, jour)
Declare.s SaisonHoraire(offset.f)
Declare.i JoursDansMois(annee, mois)

; ------------------ 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)
Global g_anomalie_vraie.f           ; Anomalie vraie de la Lune en degrés

; ------------------ Paramètres utilisateur : mois et année à analyser

Global an = 2025                                        ; ← MODIFIE ICI L'ANNEE
Global mo = 10                                          ; ← MODIFIE ICI LE MOIS (1-12)

; ------------------ Fenêtre de sortie textuelle

If OpenConsole("Phases Lunaires - " + Str(an) + "/" + Str(mo))
  ; Console ouverte, on peut écrire dedans
Else
  MessageRequester("Erreur", "Impossible d'ouvrir la console.", #PB_MessageRequester_Error)
  End
EndIf

; ------------------ Boucle sur tous les jours du mois

Define nbJours = JoursDansMois(an, mo)
Define jour, h, UTCOffset.f, JD.f, age.f, illumination.f, phase$, heure_locale$

PrintN("---------------------------------------------------------")
PrintN(" Calcul des phases lunaires pour un mois entier (3x/jour)")
PrintN("---------------------------------------------------------")
PrintN("")

For jour = 1 To nbJours
  
  ; Calcul du fuseau horaire AU DÉBUT DU JOUR
  UTCOffset = CalculUTCOffset(an, mo, jour)
  
  PrintN("---------------------------------------------------------")
  PrintN("Observation lunaire du " + Str(jour) + "/" + Str(mo) + "/" + Str(an) + " (" + SaisonHoraire(UTCOffset) + ")")
  PrintN("---------------------------------------------------------")
  
  For h = 0 To 2
    Select h
      Case 0 : Define heure_calc = 6                    ; ← MODIFIE ICI L'HEURE JOURNALIERE 1
      Case 1 : Define heure_calc = 12                   ; ← MODIFIE ICI L'HEURE JOURNALIERE 2
      Case 2 : Define heure_calc = 22                   ; ← MODIFIE ICI L'HEURE JOURNALIERE 3
    EndSelect
    
    ; Heure UTC = heure locale - décalage
    Define heure_utc = heure_calc - UTCOffset
    Define minute = 0
    Define seconde = 0
    
    ; Calcul du Jour Julien
    JD = JourJulien(an, mo, jour, heure_utc, minute, seconde)
    
    ; Calculs lunaires
    age = AgeLune(JD)
    illumination = IlluminationLune(JD)
    phase$ = PhaseLune(age, illumination)
    heure_locale$ = Str(heure_calc) + "h00"
    
    ; Affichage formaté
    PrintN("Heure locale : " + heure_locale$)
    PrintN("Âge lunaire : " + StrF(age, 2) + " jours")
    PrintN("Illumination : " + StrF(illumination, 2) + " %")
    PrintN("Phase : " + phase$)
    PrintN("")
    
  Next h
  
Next jour

PrintN("Fin du calcul pour le mois " + Str(mo) + "/" + Str(an))
PrintN("")
PrintN("---------------------------------------------------------")
PrintN("Appuyez sur Entrée pour quitter...")
Input()

End

; ------------------ Procédures inchangées (sauf ajout de JoursDansMois)

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 ---
  Protected age.f = diff / 360.0 * Lunaison
  
  ; --- CORRECTION DYNAMIQUE ---
  Protected v_corr.f = g_anomalie_vraie
  If v_corr > 180 : v_corr = 360 - v_corr : EndIf
  v_corr = (v_corr - 90) / 90.0  ; [-1, +1]
  
  Protected correction_dynamique.f = 0.05 * v_corr
  age = age + correction_dynamique
  
  ; ⚠️ LIGNE SUPPRIMÉE : heure_locale.f = heure + UTCOffset — INUTILE ET ILLÉGALE
  
  ; --- Normalisation modulo lunaison ---
  age = Mod(age, Lunaison)
  
  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    
  ; -------------------------------------------------------------
  
  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
  
  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)
  
  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)
  
  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
  
  g_anomalie_vraie = v * 180 / #PI
  
  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
  
  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)
  
  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
  
  i = i + 0.518 * Sin(2 * F * #PI / 180)
  If i < 0 : i = 0 : EndIf
  If i > 180 : i = 180 : EndIf
  
  elong = lambdaL - lambdaS
  If elong < 0 : elong + 360 : EndIf
  g_elongation = elong
  
  k = (1 + Cos(i * #PI / 180)) / 2 * 100
  k = k - 0.55 * (1 - Cos(i * #PI / 180)) * (1 - Cos(beta * #PI / 180)) * 100
  k = k - 0.04
  If k < 0 : k = 0 : EndIf
  
  ProcedureReturn k
  
EndProcedure

Procedure.s PhaseLune(age.f, illumination.f)
  
  ; ----------------------------------------------------------------------------
  ; Nomination de la phase en cours
  ; Algorithme hybride âge + illumination — calibré sur Stellarium v25.2
  ; Précision : ±0.3 jour / ±1.5% illum — transitions validées visuellement
  ; ----------------------------------------------------------------------------
  
  Protected phase$
  
  ; --- NOUVELLE LUNE ---
  If illumination < 0.5
    phase$ = "🌑 Nouvelle lune"
    
    ; --- PREMIER CROISSANT (très fin) ---
  ElseIf illumination < 5.0 And age < 7.4
    phase$ = "🌒 Premier croissant"
    
    ; --- PREMIER QUARTIER ---
  ElseIf illumination >= 48.0 And illumination <= 52.0 And age >= 7.4 And age <= 8.4
    phase$ = "🌓 Premier quartier"
    
    ; --- GIBBEUSE CROISSANTE ---
  ElseIf illumination > 52.0 And illumination < 99.5 And age < 15.8
    phase$ = "🌔 Gibbeuse croissante"
    
    ; --- PLEINE LUNE ---
  ElseIf illumination >= 99.5
    phase$ = "🌕 Pleine lune"
    
    ; --- GIBBEUSE DÉCROISSANTE ---
  ElseIf illumination > 52.0 And illumination < 99.5 And age >= 15.8 And age < 22.2
    phase$ = "🌖 Gibbeuse décroissante"
    
    ; --- DERNIER QUARTIER ---
  ElseIf illumination >= 48.0 And illumination <= 52.0 And age >= 22.2 And age <= 23.2
    phase$ = "🌗 Dernier quartier"
    
    ; --- DERNIER CROISSANT ---
  ElseIf illumination < 5.0 And age >= 23.2
    phase$ = "🌘 Dernier croissant"
    
    
    ; --- CAS PAR DÉFAUT : on utilise la logique plus simple basée sur l'âge
    ;   (au cas où un cas échappe aux conditions ci-dessus)
  Else
    If age < 7.4
      phase$ = "🌒 Premier croissant"
    ElseIf age < 8.4
      phase$ = "🌓 Premier quartier"
    ElseIf age < 14.8
      phase$ = "🌔 Gibbeuse croissante"
    ElseIf age < 15.8
      phase$ = "🌕 Pleine lune"
    ElseIf age < 22.2
      phase$ = "🌖 Gibbeuse décroissante"
    ElseIf age < 23.2
      phase$ = "🌗 Dernier quartier"
    Else
      phase$ = "🌘 Dernier croissant"
    EndIf
  EndIf
  
  ProcedureReturn phase$
  
EndProcedure

Procedure.f CalculUTCOffset(annee, mois, jour)
  
  Protected dimancheMars, dimancheOctobre, i
  
  For i = 31 To 25 Step -1
    If DayOfWeek(Date(annee, 3, i, 0, 0, 0)) = 0
      dimancheMars = i
      Break
    EndIf
  Next
  
  For i = 31 To 25 Step -1
    If DayOfWeek(Date(annee, 10, i, 0, 0, 0)) = 0
      dimancheOctobre = i
      Break
    EndIf
  Next
  
  If (mois > 3 And mois < 10) Or (mois = 3 And jour >= dimancheMars) Or (mois = 10 And jour < dimancheOctobre)
    ProcedureReturn 2.0
  Else
    ProcedureReturn 1.0
  EndIf
  
EndProcedure

Procedure.s SaisonHoraire(offset.f)
  
  If offset = 2.0
  ProcedureReturn "Heure d'été : UTC+2"
                  Else
    ProcedureReturn "Heure d'hiver : UTC+1"
  EndIf
  
EndProcedure

Procedure.i JoursDansMois(annee, mois)
  
  Select mois
    Case 1, 3, 5, 7, 8, 10, 12 : ProcedureReturn 31
    Case 4, 6, 9, 11          : ProcedureReturn 30
    Case 2
      ; Vérification bissextile
      If (annee % 4 = 0 And annee % 100 <> 0) Or (annee % 400 = 0)
        ProcedureReturn 29
      Else
        ProcedureReturn 28
      EndIf
  EndSelect
  
EndProcedure

Le programme d'origine dans sa dernière mouture :

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, illumination.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)
Global g_anomalie_vraie.f           ; Anomalie vraie de la Lune en degrés (pour corriger l'âge)

; ------------------ Date et heure UTC avec correction automatique heure d'été/hiver

Global an = Year(Date())                         ; Version date actuelle
Global mo = Month(Date())                        ;
Global jo = Day(Date())                          ;

;Global an = 2025                                  ; Version date fixe
;Global mo = 10                                    ;
;Global jo = 29                                     ;

Global UTCOffset.f = CalculUTCOffset(an, mo, jo)

Global heure = Hour(Date()) - UTCOffset           ; Version heure locale actuelle
Global minute = Minute(Date())                    ;
Global seconde = Second(Date())                   ;

;Global heure =  22 - UTCOffset                      ; Version heure locale fixe
;Global minute = 0                                  ;
;Global seconde = 0                                 ;

; ------------------ Boucle principale du programme ------------------

Define JD.f = JourJulien(an, mo, jo, heure, minute, seconde)
Define age.f = AgeLune(JD)

; Appel d'IlluminationLune
Define illumination.f = IlluminationLune(JD)
Define phase$ = PhaseLune(age, illumination)

; Affichage des résultats
MessageRequester("Observation lunaire du " + Str(jo) + "/" + Str(mo) + "/" + Str(an) + " (" + SaisonHoraire(UTCOffset) + ")",
                  "Heure locale : " + Str(heure + UTCOffset) + "h" + Str(minute) + "mn" + #LF$ +
                  "Â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 ---
  Protected age.f = diff / 360.0 * Lunaison
  
  ; --- CORRECTION DYNAMIQUE ---
  Protected v_corr.f = g_anomalie_vraie
  If v_corr > 180 : v_corr = 360 - v_corr : EndIf
  v_corr = (v_corr - 90) / 90.0  ; [-1, +1]

  Protected correction_dynamique.f = 0.05 * v_corr
  age = age + correction_dynamique
  
  ; Utilise l'heure locale = heure UTC + UTCOffset
  Protected heure_locale.f = heure + UTCOffset
  
  ; --- Normalisation modulo lunaison ---
  age = Mod(age, Lunaison)
  
  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)
  
  ; --- 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
  
  ; --- Exporte l'anomalie vraie en degrés ---
  g_anomalie_vraie = v * 180 / #PI   ; ← OUI, en degrés

  ; --- 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 Jean 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.04
  
  If k < 0 : k = 0 : EndIf    ; cas de % d'illumination négatif !
  
  ProcedureReturn k
  
EndProcedure

Procedure.s PhaseLune(age.f, illumination.f)
  
  ; ----------------------------------------------------------------------------
  ; Nomination de la phase en cours
  ; Algorithme hybride âge + illumination — calibré sur Stellarium v25.2
  ; Précision : ±0.3 jour / ±1.5% illum — transitions validées visuellement
  ; ----------------------------------------------------------------------------
  
  Protected phase$
  
  ; --- NOUVELLE LUNE ---
  If illumination < 0.5
    phase$ = "🌑 Nouvelle lune"
    
    ; --- PREMIER CROISSANT (très fin) ---
  ElseIf illumination < 5.0 And age < 7.4
    phase$ = "🌒 Premier croissant"
    
    ; --- PREMIER QUARTIER ---
  ElseIf illumination >= 48.0 And illumination <= 52.0 And age >= 7.4 And age <= 8.4
    phase$ = "🌓 Premier quartier"
    
    ; --- GIBBEUSE CROISSANTE ---
  ElseIf illumination > 52.0 And illumination < 99.5 And age < 15.8
    phase$ = "🌔 Gibbeuse croissante"
    
    ; --- PLEINE LUNE ---
  ElseIf illumination >= 99.5
    phase$ = "🌕 Pleine lune"
    
    ; --- GIBBEUSE DÉCROISSANTE ---
  ElseIf illumination > 52.0 And illumination < 99.5 And age >= 15.8 And age < 22.2
    phase$ = "🌖 Gibbeuse décroissante"
    
    ; --- DERNIER QUARTIER ---
  ElseIf illumination >= 48.0 And illumination <= 52.0 And age >= 22.2 And age <= 23.2
    phase$ = "🌗 Dernier quartier"
    
    ; --- DERNIER CROISSANT ---
  ElseIf illumination < 5.0 And age >= 23.2
    phase$ = "🌘 Dernier croissant"
    
    
    ; --- CAS PAR DÉFAUT : on utilise la logique plus simple basée sur l'âge
    ;   (au cas où un cas échappe aux conditions ci-dessus)
  Else
    If age < 7.4
      phase$ = "🌒 Premier croissant"
    ElseIf age < 8.4
      phase$ = "🌓 Premier quartier"
    ElseIf age < 14.8
      phase$ = "🌔 Gibbeuse croissante"
    ElseIf age < 15.8
      phase$ = "🌕 Pleine lune"
    ElseIf age < 22.2
      phase$ = "🌖 Gibbeuse décroissante"
    ElseIf age < 23.2
      phase$ = "🌗 Dernier quartier"
    Else
      phase$ = "🌘 Dernier croissant"
    EndIf
  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


Avatar de l’utilisateur
Philippe_GEORGES
Messages : 141
Inscription : mer. 28/janv./2009 13:28

Re: Calcul des phases lunaires pour la date du jour

Message par Philippe_GEORGES »

Bonjour à tous,

Je vais poser une question qui va vous paraitre complètement stupide.... Mais bon, je me lance !

Comment faites vous pour intégrer les petits symboles graphiques des lunaisons dans une chaine de caractères ????

Code : Tout sélectionner

phase$ = "🌓 Premier quartier"
Franchement là, c'est génial !!

J'ai dû rater un épisode !

Merci d'avance pour tout ce partage ! du très beau code ! Bien propre !

Amitiés,
Phil
Philippe GEORGES
"La simplicité est la sophistication suprême" (De Vinci)
assistance informatique, création de logiciels
Avatar de l’utilisateur
Jacobus
Messages : 1577
Inscription : mar. 06/avr./2004 10:35
Contact :

Re: Calcul des phases lunaires pour la date du jour

Message par Jacobus »

Bonjour, PhM,
Excellent, c'est bien d'avoir cette possibilité de connaître les phases de la lune par mois. Avec la sélection du jour cela fait un bon éphéméride en effet. A combiner avec d'autres options de calendrier, ça peut donner un programme sympa genre agenda :)
Philippe_GEORGES a écrit : mar. 23/sept./2025 17:47 Comment faites vous pour intégrer les petits symboles graphiques des lunaisons dans une chaine de caractères ????
Si tu as un traitement de texte comme Word, >Insertion>Symbole> cherche dans la table des caractères il y en a des tonnes, comme par exemple

Code : Tout sélectionner

 chr(25D0) et chr(25D1) " ◐ ◑ " 
Dans "Arial unicode ms" tu trouves ça par exemple :
Quand tous les glands seront tombés, les feuilles dispersées, la vigueur retombée... Dans la morne solitude, ancré au coeur de ses racines, c'est de sa force maturité qu'il renaîtra en pleine magnificence...Jacobus.
Avatar de l’utilisateur
venom
Messages : 3153
Inscription : jeu. 29/juil./2004 16:33
Localisation : Klyntar
Contact :

Re: Calcul des phases lunaires pour la date du jour

Message par venom »

Philippe_GEORGES a écrit : mar. 23/sept./2025 17:47Comment faites vous pour intégrer les petits symboles graphiques des lunaisons dans une chaine de caractères ????
C'est grâce à l'unicode comme le précise jacobus 8)

Un site d'exemple






@++
Windows 10 x64, PureBasic 5.73 x86 & x64
GPU : radeon HD6370M, CPU : p6200 2.13Ghz
Répondre