The script I'm wanting to develop uses the cartesian coordinates (XYZ) from a satellite, and in conjunction with the range, elevation and azimuth from a location, I then take a satellite’s orbital information and get the ground longitude/latitude under that satellite at a given time.
One step further from this: imagine the signal from a satellite piercing the atmosphere at exactly 300km above sea level. At this particular point when altitude is 300km, I need to calculate the ground longitude/latitude.
In the pyephem module there appears to be already a method (ephem.readtle) that can achieve this, but for TLE (two line element) data only. I'd like to use a satellite's cartesian coordinates to develop this. Is there such a method already out there? Or perhaps somebody with experience in this domain can point me in the right direction.
A similar question already exists referring to ECEF from Azimuth, Elevation, Range and Observer Lat,Lon,Alt, but it's not the same problem.
Here's what I have developed already:
- satellite cartesian coordinates, XYZ
- azimuth, elevation and range of satellite from ground station
- ground station coordinates in lat, long, height above sea level
Here's what I need: ground longitude/latitude under a satellite at a specific epoch, and in particular where the piercing point in the atmosphere (the point which the signal from the satellite pierces the atmosphere) is 300km altitude.
Ok, so I didn't resolve yet the issue with how to resolve for ground tracks for altitudes of 300km, but I believe the method I've wrote converting XYZ to ellipsoidal is complete:
def cartesian_to_ellipsoidal(self): x = 4433469.9438 y = 362672.7267 z = 4556211.6409 r = np.sqrt(x**2 + y**2 + z**2) # WGS-84 PARAMETERS, semimajor and semiminor axis a = 6378137.0 b = 6356752.314 # Eccentricity e_squared = (a**2 - b**2) / a**2 # Auxiliary quantities p = np.sqrt(x**2 + y**2) # Latitude (phi) & Longitude (lam) phi = np.rad2deg(np.arctan(z / ((1- e_squared) * p))) lam = np.rad2deg(np.arctan(y/x)) # Radius of curvature in prime vertical N = a / np.sqrt(1 - e_squared * (np.sin(np.deg2rad(phi)))**2) # Altitude h = (p / np.cos(np.deg2rad(phi))) - N return lam, phi, h
The XYZ coordinates are taken as a sample from this Matlab worked example Matlab worked example. The results are to an accuracy I'm satisfied with. What I don't understand and perhaps a dumb question, but when a ground track is calculated along with the altitude, why is the altitude not zero? One would expect since a ground track is mapped to the surface of the earth, the altitude should be zero.
One last point: the desired ground track was for where the line betweent the satellite and ground station pierces the amosphere 300km above the earth. Considering distances of 20000km (26000 km radius), then would an adjustment to my code - just to compensate for the desired 300km altitude scenario - make much a difference to outcome? If not, ground track data may just suffice.