I am trying to calculate a satellites Range Rate using Python and pyephem. Unfortunately pyephems result seems to be wrong.
After comparing the value with calculations made by other satellite tracking programs like GPredict or Ham Radio Deluxe the the difference goes up to 2km/sec.The calculated values for the Azemuth and Elevation ankle are almost the same thought. TLE's are new and the system clock is the same.
Do you see any mistake I made in my code or do you have an idea what else could cause the error?
Thank you very much!
Here is my Code:
import ephem
import time
#TLE Kepler elements
line1 = "ESTCUBE 1"
line2 = "1 39161U 13021C 13255.21187718 .00000558 00000-0 10331-3 0 3586"
line3 = "2 39161 98.1264 332.9982 0009258 190.0328 170.0700 14.69100578 18774"
satellite = ephem.readtle(line1, line2, line3) # create ephem object from tle information
while True:
city = ephem.Observer() # recreate Oberserver with current time
city.lon, city.lat, city.elevation = '52.5186' , '13.4080' , 100
satellite.compute(city)
RangeRate = satellite.range_velocity/1000 # get RangeRate in km/sec
print ("RangeRate: " + str(RangeRate))
time.sleep(1)
I recorded some Range Rate values from the script and from GPRedict to make the error reproducibly:
ESTCUBE 1
1 39161U 13021C 13255.96108453 .00000546 00000-0 10138-3 0 3602
2 39161 98.1264 333.7428 0009246 187.4393 172.6674 14.69101320 18883
date: 2013-09-13
time pyephem-Script Gpredict
14:07:02 -1.636 -3.204
14:12:59 -2.154 -4.355
14:15:15 -2.277 -4.747
14:18:48 -2.368 -5.291
And I added some lines to calculate the satellites elevation and coordinates:
elevation = satellite.elevation
sat_latitude = satellite.sublat
sat_longitude = satellite.sublong
The results with time stamp are:
2013-09-13 14:58:13
RangeRate: 2.15717797852 km/s
Range: 9199834.0
Sat Elevation: 660743.6875
Sat_Latitude: -2:22:27.3
Sat_Longitude: -33:15:15.4
2013-09-13 14:58:14
RangeRate: 2.15695092773 km/s
Range: 9202106.0
Sat Elevation: 660750.9375
Sat_Latitude: -2:26:05.8
Sat_Longitude: -33:16:01.7
Another important information might be that I am trying to calculate the Doppler Frequency for a satellite pass. And therefore I need the Range Rate:
f_Doppler_corrected = (c0/(c0 + RangeRate))*f0
Range Rate describes the velocity of a moving object on the visual axis to the observer. Maybe the range_velocity is something different?
It seems pyephem (libastro as a backend) and gpredict (predict) as a backend use different ways to calculate the satellite velocity. I am attaching detailed output of comparison for an actual reference observation. It can be seen that both output the correct position, while only gpredict outputs reasonable range_rate values. The error seems to occur in the satellite velocity vector. I would say that the reasons from gpredict are more reasonable (and the similar code is with question marks in libastro ..) therefore I will propose a fix in libastro to handle it as in gpredict, however maybe someone who understands the math behind it can add to this.
I added another tool, PyPredict (also predict based), to get some calculations here. However these values are off, so must be something else.
Pyephem
Pypredict
Debug output of Gpredict and PyEphem
PyPredict
Gpredict (sgp_math.h)
/------------------------------------------------------------------/
Ephem (Libastro)