-
Notifications
You must be signed in to change notification settings - Fork 33
/
sun_position.jl
77 lines (60 loc) · 2.24 KB
/
sun_position.jl
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
# # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # #
#
# Description
#
# Compute the sun position.
#
# # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # #
#
# References
#
# [1] The Astronomical Almanac for the year 2000 (p. C24).
#
# [2] http://aa.usno.navy.mil/faq/docs/SunApprox.php
#
# [3] Vallado, D. A (2013). Fundamentals of Astrodynamics and Applications.
# Microcosm Press, Hawthorne, CA.
#
# [4] The Astronomical Almanac for the year 2006.
#
# # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # #
export sun_position_i
"""
sun_position_i(JD::Number)
Compute the Sun position represented in the Mean Equinox of Date (MOD) at the
Julian Day `JD`. The algorithm was adapted from [3, p. 277-279].
"""
function sun_position_i(JD::Number)
# Constants.
deg2rad = π/180
# Number of Julian centuries from J2000 epoch.
T_UT1 = (JD-JD_J2000)/36525.0
# Here, we will assume that T_TBD ≈ T_UT1. This can be done because we are
# using a low-precision algorithm [3].
T_TBD = T_UT1
# Mean anomaly of the Sun [deg].
Ms = 357.529_109_2 + 35_999.050_34T_TBD
# Convert Ms to [rad] and limit to the interval [0,2π].
Ms = mod2pi(Ms*deg2rad)
# Compute auxiliary variables.
sinMs, cosMs = sincos(Ms)
sin2Ms, cos2Ms = sincos(2Ms)
# Mean longitude of the Sun [deg].
λ_m = 280.460 + 36_000.771T_UT1
# Ecliptic latitude of the Sun [deg].
λ_e = λ_m + 1.914_666_471sinMs + 0.019_994_643sin2Ms
# Obliquity of the ecliptic [deg].
ϵ = 23.439_291 - 0.013_004_2T_TBD
# Convert λ_e and ϵ to [rad] and limit to the interval [0,2π].
λ_e = mod2pi(λ_e*deg2rad)
ϵ = mod2pi( ϵ*deg2rad)
# Auxiliary variables.
sinϵ , cosϵ = sincos(ϵ)
sinλ_e , cosλ_e = sincos(λ_e)
# Distance of the Sun from Earth [m].
r = ( 1.000_140_612 - 0.016_708_617cosMs - 0.000_139_589cos2Ms )*au2m
# Compute the Sun vector represented in the Mean Equinox of Date (MOD).
S_MOD = SVector{3,Float64}(r*cosλ_e, r*cosϵ*sinλ_e, r*sinϵ*sinλ_e)
# TODO: This vector must be transformed to the inertial reference frame that
# will be used in this toolbox (J2000, True of Date, etc.).
end