pengmun / destopy

Thermospheric Mass Density Estimation Toolbox in python

DESTOPy - Atmospheric Density Estimation Toolbox in Python

The Atmospheric Density Estimation Toolbox in python (DESTOPy) is a python toolbox for estimation of the global thermospheric mass density based on the Matlab Density Estimation Toolbox (DESTO) by David Gondelach and Richard Linares.

Currently, only a single reduced-order density model (JB2008 reduced-order model) is available and the python toolbox only supports assimilation of Two-Line Element (TLE) sets.

Copyright © 2021 by Peng Mun Siew and Richard Linares

License

The code is licensed under the MIT License - see the LICENSE file for details.

Acknowledgments

The code is developed based upon work supported by the National Science Foundation under award NSF-PHY-2028125.

The python scripts are coded based on work taken from the Matlab Density Estimation Toolbox (DESTO) by David Gondelach and Richard Linares.

The Jacchia-Bowman 2008 model was downloaded from the ATMOS python package by Chunxiao Li.

The solar radiation pressure and third-body perturbations models are adapted from the Matlab code by Meysam Mahooti (copyright 2018) and were downloaded from Matlab File Exchange.

The SGP4 model used was based on the python-sgp4 module by Brandon Rhodes.

Several time and reference frame routines are adapted from the Matlab code developed by David Vallado (and others) and were downloaded from https://celestrak.com/software/vallado-sw.php.

The Earth Gravitational Model 2008 (EGM2008) coefficients were obtained from the NGA's Office of Geomatics: https://earth-info.nga.mil.

The toolbox also makes use of NASA's SPICE toolkit (N66) via the SpiceyPy module by Annex et al. available here.

Some of the scripts also take advantage of Numba for faster computations.

References

The thermospheric density modeling and estimation techniques using Two-Line Element (TLE) sets are described in:

@article{gondelach2020tle,
  author = {Gondelach, David J. and Linares, Richard},
  title = {Real-Time Thermospheric Density Estimation Via Two-Line-Element Data Assimilation},
  journal = {Space Weather},
  volume = {18},
  number = {2},
  pages = {e2019SW002356},
  doi = {10.1029/2019SW002356},
  url = {https://doi.org/10.1029/2019SW002356}
}

see https://doi.org/10.1029/2019SW002356.

Dependencies (partial list)

• python
• numpy
• scipy
• matplotlib
• spiceypy
• datetime
• sgp4
• numba
• h5py

Installation instructions

1. Download the DESTOPy python code.
2. Install the necessary dependencies by running 'pip install -r requirements.txt'.
3. Download SPICE kernels (i.e. ephemeris files) from https://naif.jpl.nasa.gov/pub/naif/generic_kernels/ and put them in the folder Data. See Ephemeris files section below.
4. Download space weather file from Celestrak and put in folder Data: https://www.celestrak.com/SpaceData/SW-All.txt
5. Download Earth orientation data file from Celestrak and put in folder Data: https://www.celestrak.com/SpaceData/EOP-All.txt
6. Download 2 space weather files needed for the JB2008 model and put in folder Data: http://sol.spacenvironment.net/jb2008/indices/SOLFSMY.TXT and http://sol.spacenvironment.net/jb2008/indices/DTCFILE.TXT
7. (Optional) Download new TLE data and put them in the folder TLEdata.
8. For each space object used for estimation, ensure that their ballistic coefficient (BC) are specified in the text file: Data/BCdata.txt

Run instructions

1. Open the DensityEstimationTLE.ipynb or SO_propagation.ipynb using jupyter notebook.
2. (Optional) Under section 1.0.1 Settings, select the estimation window (start year, start month, start day, number of days, and number of hours).
3. (Optional) Under section 1.0.1 Settings, select the reduced-order model dimensions.
4. (Optional) Under section 1.0.1 Settings, select the space objects to use for density estimation.
5. Run all cells.

Ephemeris files

Download the following ephemeris files and put them in the Data folder (Note: Update the ephemeris filenames in Data/kernel.txt to match the downloaded files.):

• latest_leapseconds.tls: https://naif.jpl.nasa.gov/pub/naif/generic_kernels/lsk/
• de430.bsp: https://naif.jpl.nasa.gov/pub/naif/generic_kernels/spk/planets/
• earthstns_itrf93_201023.bsp: https://naif.jpl.nasa.gov/pub/naif/generic_kernels/spk/stations/
• pck00010.tpc, earth_fixed.tf, earth_200101_990628_predict.bpc, earth_000101_210530_210308.bpc, earth_720101_070426.bpc, earth_latest_high_prec.bpc: https://naif.jpl.nasa.gov/pub/naif/generic_kernels/pck/

Python 3.7 was used to develop the code and was tested on Python 3.11.

Peng Mun Siew, Oct 2023 email: [email protected]

destopy's Issues

object of type 'NoneType' has no len()

Unable to perform rom density output

Dear sir or madam
I couldn't input my own latitude and longitude information to get the correct rom density output . I'm not sure how to use getDensityROM.I need your help.pls!!

I can't output the same curve

wanna run the program

Dear sir/madam,
I am an engineering student and now calculated the TLE in many aspects. That I have downloaded all the needed Dependencies and files and put them into the Data file.
When I start the program everything is fine, but when running the cell of 'Run unscented Kalman filter', it prompts an error!

0 24

RemoteTraceback Traceback (most recent call last)
RemoteTraceback:
"""
Traceback (most recent call last):
File "C:\Users\WNT\anaconda3\lib\multiprocessing\pool.py", line 125, in worker
result = (True, func(*args, **kwds))
File "C:\Users\WNT\anaconda3\lib\multiprocessing\pool.py", line 48, in mapstar
return list(map(*args))
File "C:\Users\WNT\DESTOPy-ae771d5140d20cb35dd9ef156055ff8bf49e3091\python_utils\astro_function_jit.py", line 2001, in propagateState_MeeBcRom
xf_out = ode_ope.integrate(tf)
File "C:\Users\WNT\anaconda3\lib\site-packages\scipy\integrate_ode.py", line 433, in integrate
self._y, self.t = mth(self.f, self.jac or (lambda: None),
File "C:\Users\WNT\anaconda3\lib\site-packages\scipy\integrate_ode.py", line 1009, in run
y1, t, istate = self.runner(*args)
File "C:\Users\WNT\DESTOPy-ae771d5140d20cb35dd9ef156055ff8bf49e3091\python_utils\astro_function_jit.py", line 1756, in computeDerivative_PosVelBcRom
x2 = x[b].transpose(1, 0, 2).reshape(3,-1)
IndexError: index 1146574368 is out of bounds for axis 0 with size 150
"""

The above exception was the direct cause of the following exception:

IndexError Traceback (most recent call last)
~\AppData\Local\Temp\ipykernel_29436\3028581531.py in
1 # Run Unscented Kalman filter estimation
----> 2 X_est,Pv = UKF(X_est,Meas,time,P,RM,Q,nop,svs,r, AC, BC, SWinputs, F_U, M_U, maxAtmAlt, et0, jd0, highFidelity, GM, Re, C_gravmodel,S_gravmodel, gravdegree, sF_gravmod)
3 Pv[:,0] = np.diag(P); # Add initial state variance to state variance history

~\DESTOPy-ae771d5140d20cb35dd9ef156055ff8bf49e3091\python_utils\astro_function_jit.py in UKF(X_est, Meas, time, P, RM, Q, nop, svs, r, AC, BC, SWinputs, F_U, M_U, maxAtmAlt, et0, jd0, highFidelity, GM, Re, C_gravmodel, S_gravmodel, gravdegree, sF_gravmod)
2110 # Time Update
2111 with poolcontext(processes=multiprocessing.cpu_count()-1 or 1) as pool:
-> 2112 Xp = pool.map(partial(propagateState_MeeBcRom, t0=time[i], tf=time[i+1], AC=AC, BC=BC, SWinputs=SWinputs, r=r, nop=nop, svs=svs, F_U=F_U, M_U=M_U,maxAtmAlt=maxAtmAlt,et0=et0,jdate0=jd0,highFidelity=highFidelity,GM=GM, Re=Re, C_gravmodel=C_gravmodel, S_gravmodel=S_gravmodel, gravdegree=gravdegree, sF_gravmod=sF_gravmod), xx.T)
2113
2114 Xp = np.squeeze(np.asarray(Xp).T)

~\anaconda3\lib\multiprocessing\pool.py in map(self, func, iterable, chunksize)
362 in a list that is returned.
363 '''
--> 364 return self._map_async(func, iterable, mapstar, chunksize).get()
365
366 def starmap(self, func, iterable, chunksize=None):

~\anaconda3\lib\multiprocessing\pool.py in get(self, timeout)
769 return self._value
770 else:
--> 771 raise self._value
772
773 def _set(self, i, obj):

IndexError: index 1146574368 is out of bounds for axis 0 with size 150

I have checked all of the parameters, and could not find the reason. I am wondering you could help me to solve the problem!
Thank you a lot!!!

IndexError: index 0 is out of bounds for axis 0 with size 0

Dear sir/madam,
I am an engineering student and I'm doing some work on the atmosphere. That I have downloaded all the needed Dependencies and files and put them into the Data file. However, there was an error in “Generate reduced-order density model ”section.

display IndexError: index 0 is out of bounds for axis 0 with size 0

I checked all the parameters, but I haven't found the problem yet，I'd like some guidance. I'd appreciate it。

error：
IndexError Traceback (most recent call last)
~\AppData\Local\Temp\ipykernel_5920\3705748063.py in
1 ## Load reduced-order density models
----> 2 AC,BC,Uh,F_U,Dens_Mean,M_U,SLTm,LATm,ALTm,maxAtmAlt,SWinputs,Qrom = generateROMdensityModel(ROMmodel,r,jd0,jdf,spice)

D:\Code2024\py1\DESTOPy-main\DESTOPy-main\python_utils\astro_function_jit.py in generateROMdensityModel(ROMmodel, r, jd0, jdf, spice)
1271 # Read geomagnetic storm DTC values
1272 DTCdata = readDTCFILE('Data/DTCFILE.txt');
-> 1273 SWinputs = computeSWinputs_JB2008(jd0,jdf+1,eopdata,SOLdata,DTCdata,spice);
1274
1275 # Maximum altitude of ROM-JB2008 density model (for higher altitudes

D:\Code2024\py1\DESTOPy-main\DESTOPy-main\python_utils\astro_function_jit.py in computeSWinputs_JB2008(jd0, jdf, eopdata, SOLdata, DTCdata, spice)
1169 UThrs = hhUTC + mnmnUTC/60 + ssUTC/3600;
1170 # Get JB2008 space weather data
-> 1171 _,GWRAS,SUN,F10,F10B,S10,S10B,XM10,XM10B,Y10,Y10B,DSTDTC = computeJB2000SWinputs(yyUTC,doyUTC,hhUTC,mnmnUTC,ssUTC,SOLdata,DTCdata,eopdata,spice);
1172
1173 # [jdate; doy; UThrs; F10; F10B; S10; S10B; XM10; XM10B; Y10; Y10B; DSTDTC; GWRAS; SUN(1); SUN(2)]

D:\Code2024\py1\DESTOPy-main\DESTOPy-main\python_utils\astro_function_jit.py in computeJB2000SWinputs(year, doy, hour, minute, sec, SOLdata, DTCdata, eopdata, spice)
1087 # USE 1 DAY LAG FOR F10 AND S10 FOR JB2008
1088 JD = np.rint(JD_0)-1; # 1 day lag (JD=MJD+2400000.5)
-> 1089 i = np.where(JD==SOLdata[2,:])[0][0];
1090 SOL = SOLdata[:,i];
1091 F10 = SOL[3];

IndexError: index 0 is out of bounds for axis 0 with size 0