Closed tmamedzadeh closed 3 years ago
Hi @Tarlan0001 ,
No, I have not yet. This must be done because I want a numerical orbit propagator. I am working right now in porting another atmospheric density model into SatelliteToolbox (JB2008), which seems more accurate than NRLMSISE-00 that is already coded.
I'm also working on a propagator, comparing the results with GMAT. Currently, I'm taking into account the Newtonian gravity (Earth, Sun, Moon), harmonics (JGM-3 10th degree), SRP and relativistic effect. The difference for 70 days made 30 m. How do you think, may I neglect the atmospheric drag on the 700km, or should add it also? Which other effects do you think should I consider additionally?
The density at 700 km is pretty low. However, it seems more important than relativistic effects, for example, I think. Hence, if you really want high precision, then you must take this into account.
Have you already thought about the implementation? I think, the algorithm in Mathlab is complicated. Is there any other way?
For example, convert to ECEF, calculate the relative velocity and then acceleration.
I have already implemented a very simple one considering only the high-order gravitation. I am not sure how I will do this, but it is on the schedule. First, I need to get more atmospheric models, I need to figure out how to use the package JPLEphemeris.jl, and then how is the best integrator for the problem.
In the nrlmsise00 example you set the daily magnetic index to 12. How is it possible? Or, it's not the KP index?
Hi @Tarlan0001 ,
Actually the input to NRLMSISE00 is the AP and not the Kp.
Is there a way to use the function as MSISE90, not NRLMSISE00, or they are different? The point is, I want to compare results with GMAT, but there is only MSISE90 and Jacchia-Roberts.
MSISE90 is the old version of NRLMSISE00, and Jacchia-Roberts is the analytical version of Jacchia 1971 model. I found papers showing that Jacchia-Bowman (2008) is more accurate. Do you have any use case for those older models? I can add them to SatelliteToolbox.
@Tarlan0001
Is it possible, please, to provide some data from Jacchia-Roberts model? I need a table like this:
Altitude [km] | Density [kg/m^3]
The altitude should vary from 90 to 1000km. Can you do this for me please? I need to validate the implementation. You can select whatever latitude, longitude, F10.7, Kp you want, just let me know the inputs.
Actually, I don't have the implementation of Jacchia-Roberts. I used it in the General Mission Analysis Tool.
Is there a way to download the actual AP index and flux data?
Hi @Tarlan0001 ,
The AP I think you can get from GFZ website.
I have just pushed a commit with The Jacchia-Roberts atmosphere model! You can use it by:
julia> jr1971(DatetoJD(2017,1,1,0,0,0), 45*pi/180, 0*pi/180, 300e3, 76.5, 72.5, 7)
JR1971_Output{Float64}
nN2: Float64 4.0269432130021164e13
nO2: Float64 1.8040704522884658e12
nO: Float64 4.976154354226309e14
nAr: Float64 4.914191344437922e9
nHe: Float64 1.0731503643725797e13
nH: Float64 0.0
rho: Float64 1.526144428631946e-11
T_exo: Float64 871.4209785497403
Tz: Float64 858.5927825073542
The documentation is:
function jr1971(JD::Number, glat::Number, glon::Number, h::Number, F10::Number,
Compute the atmospheric density using the Jacchia-Roberts 1971 model.
Args
≡≡≡≡≡≡
• JD: Julian day.
• glat: Geodetic latitude [rad].
• glon: Geodetic longitude [rad].
• h: Altitude [m].
• F10: 10.7-cm solar flux [10⁻²² W/(M² Hz)].
• F10ₐ: 10.7-cm averaged solar flux, 81-day centered on input time.
• Kp: Kp geomagnetic index (with a delay of 3 hours).
Returns
≡≡≡≡≡≡≡≡≡
An instance of the structure JR1971_Output with the computed values.
Hi @Tarlan0001 ,
I have implemented a way to get the online indices in SatelliteToolbox. However, the services are down and you will not be able to download them. When everything is normalised, I will test and let you know.
In the satellite propagation algorithm, you have to use the actual daily data. This will be very useful!
Yes! It will be very handy. The code is already on master
. However, the SET server is going on a huge upgrade and the files are not available now.
If you want to try, you can do:
init_space_indices()
get_Ap(DatetoJD(2018,6,19,10,0,0); daily = true)
This will get the daily mean of Ap at 19 June 2018.
Thank you, I'll check. What about the daily and 81 day average solar flux? What if there is no data for the given date, would it acquire the previous existing data? You could also make it also a part of the nrlmsise00 algorithm.
If the data is not available, it will use from the nearest available date.
The 81 day average F10.7 can be obtained by get_F81a
.
There is already a function to call NRLMSISE00 in which it will obtain all the values automatically. However, the solar flux depends on the SET server that is still offline.
Could you please describe the usage of NRLMSISE00, which obtains the data automatically?
@Tarlan0001 sorry about the delay. The code is ready now! You need to do this in order to run NRLMSISE00 fetching the data automatically:
(remember to pull the latest code from master
)
# Initialize only the space indices.
init_space_indices(enabled_files = [:wdcfiles, :fluxtable])
julia> nrlmsise00(DatetoJD(2016,10,19,14,00,00), 750e3, -23*pi/180, -45*pi/180; output_si = false, dversion = false)
NRLMSISE00_Output{Float64}
den_N: Float64 2860.4183009695835
den_N2: Float64 12.141431187588672
den_O: Float64 67781.3862464794
den_aO: Float64 5044.235371497685
den_O2: Float64 0.04503651487434012
den_H: Float64 148814.26555193574
den_He: Float64 700462.0267827582
den_Ar: Float64 2.0371887656858425e-6
den_Total: Float64 6.7654160048712925e-18
T_exo: Float64 821.1537129687955
T_alt: Float64 821.1536995447245
flags: Dict{Symbol,Bool}
Compare with the online NRLMSISE-00 result:
Selected parameters are:
1 Height, km
2 O, cm-3
3 N2, cm-3
4 O2, cm-3
5 Mass_density, g/cm-3
6 Temperature_neutral, K
7 Temperature_exospheric, K
8 He, cm-3
9 Ar, cm-3
10 H, cm-3
11 N, cm-3
12 Anomalous_Oxygen, cm-3
1 2 3 4 5 6 7 8 9 10 11 12
750.0 6.820E+04 1.227E+01 4.549E-02 6.785E-18 821.6 822 7.017E+05 2.067E-06 1.486E+05 2.875E+03 5.046E+03
Notice that if you want to compute something in a data more than 3 years ago, then you must pass the option wdcfiles_oldest_year
to the function init_space_indices
. You can look at the documentation of the function to get more information.
Hi @Tarlan0001
Could you compute the drag using those functions? Can we close this bug?
Closing due to inactivity.
Have you developed the calculation of the acceleration due to the atmospheric drag?
I've found a function in Mathlab https://www.mathworks.com/matlabcentral/fileexchange/54877-satellite-orbit-modeling?focused=7116317&tab=function.