Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Compared to the position of the sun in STK, there is a significant difference. Why does this situation occur? #3

Open
blackbutton opened this issue Aug 12, 2024 · 1 comment

Comments

@blackbutton
Copy link

 jd = date_to_jd(2024,8,6,4,0,0)
s_mod = sun_position_mod(jd)
r_mod_j2000 = r_eci_to_eci(MOD(),J2000(),jd)
r_mod_j2000 * s_mod
3-element StaticArraysCore.SVector{3, Float64} with indices SOneTo(3):
 -1.0512053177080336e11
  1.0039274359027356e11
  4.351950513418951e10

STK:
-105129703.817743 100395407.099930 43520964.349647

@ronisbr
Copy link
Member

ronisbr commented Aug 12, 2024

Hi @blackbutton !

The algorithm we currently have is a "low-precision" implemented from the Astronomical Almanac (and described on Vallado's book). We still lack a high precision Sun position algorithm. Nevertheless, there is a minor problem in your code because sun_position_mod expects the time in TBD (Barycentric Dynamical Time) instead of UTC. The difference between the Terrestrial time and TBD is period by Earth's orbit. Hence, we can assume that both are the same here:

julia> jd = date_to_jd(2024, 8, 6, 4, 0, 0)
2.4605286666666665e6

julia> jd_utc = date_to_jd(2024, 8, 6, 4, 0, 0)
2.4605286666666665e6

julia> jd_tt = jd_utc_to_tt(jd_utc)
2.4605286674674074e6

julia> s_mod = sun_position_mod(jd_tt)
3-element StaticArraysCore.SVector{3, Float64} with indices SOneTo(3):
 -1.05776297878576e11
  9.981141741230861e10
  4.326692615766075e10

julia> s_j2000 = r_eci_to_eci(MOD(), J2000(), jd) * s_mod
3-element StaticArraysCore.SVector{3, Float64} with indices SOneTo(3):
 -1.0512198446655511e11
  1.003914406336404e11
  4.351894032716577e10

Now, if we take the angular difference between the STK answer and the one we just obtained, we get:

julia> s_j2000_stk = [-105129703.817743, 100395407.099930, 43520964.349647]
3-element Vector{Float64}:
 -1.05129703817743e8
  1.0039540709993e8
  4.3520964349647e7

julia> using LinearAlgebra

julia> acosd(s_j2000 / norm(s_j2000)  s_j2000_stk / norm(s_j2000_stk)) * 3600
3.4026739070005374

The result indicates that we differ from STK by less than 3.5 arcs, which is more that enough for most satellite applications. If you want a really accurate algorithm for satellite purposes, we can implement it here.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

2 participants