Space
The great thing about the #30DayMapChallenge is that you can choose how many challenges to get involved in. Some people go all out and manage every day. I tend to be more selective. Day 09 - Space seemed like a great opportunity to dig out an old globe model and look outwards from planet earth. Researching possibilities, I stumbled across the amazing Skyfield that is able to predict the positions of Earth satellites.
Behind the Scenes
Skyfield from https://rhodesmill.org/skyfield/earth-satellites.html
Satellite positions from http://celestrak.org/NORAD/elements/
Working thorugh the code snippets on the Skyfield pages, I ended up with this code to get the ISS positions:
from skyfield.api import load, wgs84
import math
stations_url = 'http://celestrak.com/NORAD/elements/stations.txt'
#stations_url = 'https://celestrak.org/NORAD/elements/gp.php?GROUP=starlink&FORMAT=tle'
#stations_url = 'https://celestrak.org/NORAD/elements/gp.php?GROUP=active&FORMAT=tle'
satellites = load.tle_file(stations_url)
print('Loaded', len(satellites), 'satellites')
#for sat in satellites:
# print(sat.epoch.utc_jpl())
#sat = satellites[0]
ts = load.timescale()
for hour in range(0, 25, 1):
for min in range(0, 61, 1):
t = ts.utc(2022, 11, 9, hour, min, 0)
# print("lat, lon, height")
# for sat in satellites:
sat = satellites[0]
geocentric = sat.at(t)
# print(geocentric.position.km)
lat, lon = wgs84.latlon_of(geocentric)
height = wgs84.height_of(geocentric)
height_km = height.km
if (math.isnan(lat.degrees) == False and math.isnan(lon.degrees) == False and math.isnan(height_km) == False):
print(f"({lat.degrees},{lon.degrees},{height_km}),")
This gives a set of points that we can load into a python list
points = [
(30.197778074499244,-41.69238013378573,414.3551803386221),
...
many entries excluded
...
(51.44188187480004,9.04164498351871,419.61570652925167)
]
In Blender, we can use this points list to create a point cloud
def point_cloud(ob_name, coords, edges=[], faces=[]):
"""Create point cloud object based on given coordinates and name.
Keyword arguments:
ob_name -- new object name
coords -- float triplets eg: [(-1.0, 1.0, 0.0), (-1.0, -1.0, 0.0)]
"""
# Create new mesh and a new object
me = bpy.data.meshes.new(ob_name + "Mesh")
ob = bpy.data.objects.new(ob_name, me)
# Make a mesh from a list of vertices/edges/faces
me.from_pydata(coords, edges, faces)
# Display name and update the mesh
ob.show_name = True
me.update()
return ob
Stringing everything together we end up with a short function to create the ISS trails and positions.
getXYZfromLatLon() was covered in Globe Flights and addCurve() will be covered in my next post...
def showPoints():
R = 6371 / 1000
coords = []
for p in points:
newPoint = getXYZfromLatLon(R+p[2]/1000, p[0], p[1])
coords.append(newPoint)
pc = point_cloud('satellites', coords)
bpy.context.collection.objects.link(pc)
addCurve(coords)
showPoints()
