Code : Tout sélectionner
;+==========================================================+
;| |
;| Calcul de l'heure TU du lever et du coucher du Soleil |
;| |
;| Précision : La minute de temps |
;| Plateforme: Multiplateforme |
;| ToDo : Reste à optimiser le code (réduction du code)|
;| Auteur : Mesa 30 Janvier 2013 |
;| |
;+==========================================================+
;Grenoble +45° 11' 16" N (Nord = positif, Sud = négatif)
; -05° 43' 37" E (Est = négatif, Ouest = Positif)
;Altitude Min. 204 m Max. 475 m
; Bureau des longitudes : heure de lever du soleil
; http://www.imcce.fr/fr/ephemerides/phenomenes/rts/rts.php
; Soleil le 01/01/2013
; Lieu : 05°43’37" E / 45°11’16" N
; Date Lever Passage au méridien Coucher
; Temps Universel heure azimut heure hauteur heure azimut
; 2013-01-01 07:18 302.90 11:41 21.85 16:04 57.14
; Date Aube Crépuscule
; astronom. nautique civil civil nautique astronom.
; Temps Universel heure heure heure heure heure heure
; 2013-01-01 05:29 06:05 06:42 16:39 17:17 17:52
EnableExplicit
Procedure.d Frac(Xfrac.d)
ProcedureReturn Xfrac-Int(Xfrac)
EndProcedure
Procedure.d dms(Xdms.d)
;conversion decimales en degres minutes secondes
ProcedureReturn Int(Xdms)+Frac(Int(Frac(Xdms)*60)/100)+Frac(Frac(Frac(Xdms)*60)*0.006)
EndProcedure
Procedure.d dec(xdec.d)
;conversion degres minutes secondes en decimales
ProcedureReturn Int(Xdec)+Int(frac(Xdec)*100)/60+(frac(Xdec)*100-Int(frac(Xdec)*100))*100/3600
EndProcedure
Procedure.d jd2(Annee, Mois, jour, Heure, Minute, Seconde); Jour Julien à 0h TU (Temps Universel )
; --Ne fonctionne pas avant le 15 Octobre 15 1582
Protected j.d
Protected m
Protected a
Protected jj.d
a=Annee
m=Mois
j=jour+(Heure + Minute/60 + Seconde/3600)/24
jj=367*a -Int(1.75*(Int((m+9)/12)+a))
jj-Int(0.75*(1+Int((Int((m-9)/7)+a)*0.01)))
jj+Int((275*m)/9)+j+1721028.5
ProcedureReturn jj
EndProcedure
Procedure.d LeverCoucherSoleil(LongitudeT.d, LatitudeT.d, Jour, Mois, An, Lever=#True, Option.f=-0.61)
; Option
; -18<= Option <-12 : Crépuscule astronomique
; -12<= Option <-6 : Crépuscule nautique
; -6 <= Option <réfraction : Crépuscule civil
;
; réfraction = -36.6'= -0.61° : Ephémérides du Bureau des longitudes
; = -34' : Ephéméride nautique
; = centre du Soleil : Ephémérides astronomiques
; = bord inférieur du Soleil : Ephémérides pour navigateurs
;
; Lever
; Lever=#True : Calcul du lever du Soleil
; Lever=#False : Calcul du coucher du Soleil
Protected.d FA, H0, FI, LO, JD, T0, S0, TU, TS, T
Protected.d LL, M1, C, TL, EX, AL, X, DE, AH, TM
Protected IN0, i
FA=180/#PI
H0=0
; Radian
FI=dec(LatitudeT)/FA
LO=dec(LongitudeT)/FA
; Jour julien de la date à 00h:00'00" (ne pas changer l'heure)
JD=jd2(An,Mois,Jour,0,0,0)
;Debug JD
;Temps Sideral à Greenwitch
T0 = (JD - 2415020.0)/36525
S0= 6.6460656+2400.051262*T0 + 0.00002591*T0*T0
S0 = S0 - 24*Int(S0/24)
TU = 12
TS = S0 + TU *1.002737908
;Debug TS
;Position de la Terre sur son orbite (ou soleil)
T= T0+TU/24/36525
LL=279.69668+T*36000.76892+T*T*0.0003025
LL / FA
M1 = 358.47583+T*35999.04975-T*T*0.000151-T*T*T*0.0000033
M1 / FA
C = (1.91946-0.004789*T-T*T*0.000014)*Sin(M1)+(0.020094-T*0.0001)*Sin(2*M1)+0.000293*Sin(3*M1)
C / FA
;Longitue ecliptique
TL = LL +C
;Coordonnées écliptique
EX = 23.452294-0.0130125*T-0.00000164*T*T+T*T*T*0.000000503
EX / FA
AL = ATan(Cos(EX)*Tan(TL))
If Cos(TL)<0
AL + #PI
EndIf
X = Sin(EX)*Sin(TL)
DE= ATan(X/Sqr(-X*X+1))
;Angle Horaire
X=(Sin(H0)-Sin(FI)*Sin(DE))/Cos(FI)/Cos(DE)
If Abs(X)>1
MessageRequester("impossible",StrD(H0*FA))
End
EndIf
AH=-ATan(X/Sqr(-X*X+1))+ #PI/2
TM=TS-(LO+AL)*FA/15
TM=36-TM
TM=TM-24*Int(TM/24)
;Lever
If Lever=#True
IN0=-1
Else
IN0=1
EndIf
H0=Option/FA
For i= 1 To 3
;Position de la Terre sur son orbite (ou soleil)
T= T0+TU/24/36525
LL=279.69668+T*36000.76892+T*T*0.0003025
LL / FA
M1 = 358.47583+T*35999.04975-T*T*0.000151-T*T*T*0.0000033
M1 / FA
C = (1.91946-0.004789*T-T*T*0.000014)*Sin(M1)+(0.020094-T*0.0001)*Sin(2*M1)+0.000293*Sin(3*M1)
C / FA
;Longitue ecliptique
TL = LL +C
;Coordonnées écliptique
EX = 23.452294-0.0130125*T-0.00000164*T*T+T*T*T*0.000000503
EX / FA
AL = ATan(Cos(ex)*Tan(TL))
If Cos(TL)<0
AL + #PI
EndIf
X = Sin(EX)*Sin(TL)
DE= ATan(X/Sqr(-X*X+1))
;Angle Horaire
X=(Sin(H0)-Sin(FI)*Sin(DE))/Cos(FI)/Cos(DE)
If Abs(X)>1
MessageRequester("impossible",StrD(H0*FA))
End
EndIf
AH=-ATan(X/Sqr(-X*X+1))+ #PI/2
TS=(IN0*AH+AL+LO)*FA/15
While TS<S0
TS=TS+24
Wend
TU=(TS-S0)/1.002737908
TU=TU-24*Int(TU/24)
Next i
ProcedureReturn TU
EndProcedure
;Affichage
Debug "Le 01/01/2013 à Grenoble"
Debug ""
Debug "Lever/Coucher Civil: "
Debug dms(LeverCoucherSoleil(-5.4337,45.1116,1,1,2013)) ;07h17'37" ~ 07h 18'
Debug dms(LeverCoucherSoleil(-5.4337,45.1116,1,1,2013,#False))
Debug ""
Debug "Lever/CoucherNautique: "
Debug dms(LeverCoucherSoleil(-5.4337,45.1116,1,1,2013,#True,-12))
Debug dms(LeverCoucherSoleil(-5.4337,45.1116,1,1,2013,#False,-12))
Debug ""
Debug "Lever/CoucherAstronomique (Nuit noire): "
Debug dms(LeverCoucherSoleil(-5.4337,45.1116,1,1,2013,#True,-18))
Debug dms(LeverCoucherSoleil(-5.4337,45.1116,1,1,2013,#False,-18))
