2
$\begingroup$

I'm trying to calculate the position of each satellite of the GPS constalation, using TLE files and library pyephem, I want to know the position Elevation-Azimuth. But when I compare with the results with the page calsky with the results of my code are diferents, why?

The initial information are: Bogotá, Colombia
Longitude = -74.0817 Latitude = 4.60971 Alt: 2578.0
Date: 2017-04-18 Hour: 10:08:51.61

Results in CALSKY >45º

PRN 27/GPS BIIF-(39166 2013-023-A) az: 152.7° SSE h: 61.1°
PRN 16/GPS BIIR-(27663 2003-005-A) az: 341.8° NNW h: 52.5°
PRN 14/GPS BIIR-(26605 2000-071-A) az: 109.5° ESE h: 49.5°

Results in my code > 45º

GPS BIIR-7 (PRN 18) - azim: 27.68 elev: 58.34
GPS BIIR-9 (PRN 21) - azim: 200.54 elev: 55.05
GPS BIIF-1 (PRN 25) - azim: 111.35 elev: 46.92

def GetPositionELAZ(Time,Longitude,Latitude,constelation,AllTheLines):
'''
- time: Time for look the satellite
- longitude: observer longitude
- latitude: observer latitude
- constelation: 
            1 - just GPS 
            2 - GPS - GLONASS
            3 - GPS - GLONASS - GALILEO
            
RETURN: (Tuble of lists...)
Satellite Name, Satellite Elevation, Satellite Azimuth, Type Satellite (1:GPS, 2: GLO, 3:GAL)
'''

contain_GPS = 'GPS'
contain_GLONASS = 'COSMOS'
contain_GALILEO = 'GSAT'

sat_alt, sat_az, sat_name, sat_type = [], [], [], []
original_name = []

longlat = str(Longitude) + "," + str(Latitude)
g = geocoder.elevation(longlat) #lon,lat
#print g.meters
import ephem
observer = ephem.Observer()
observer.long = Longitude
observer.lat = Latitude
observer.date = Time
observer.elevation = 2600. #float(str(g.meters))

GPS_list = 'https://celestrak.org/NORAD/elements/gps-ops.txt'
GLONASS_list = 'https://celestrak.org/NORAD/elements/glo-ops.txt'
GALILEO_list = 'https://celestrak.org/NORAD/elements/galileo.txt'

#'http://www.amsat.org/amsat/ftp/keps/current/nasabare.txt').readlines()
GPS2_list = 'http://www.tle.info/data/gps-ops.txt'
GLONASS2_list = 'http://www.tle.info/data/glo-ops.txt'

if constelation == 1:
    list_TLEs = [GPS_list]
if constelation == 2:
    list_TLEs = [GPS_list,GLONASS_list]
if constelation == 3:
    list_TLEs = [GPS_list,GLONASS_list,GALILEO_list]

for satlist in list_TLEs:
    tles = AllTheLines #urllib2.urlopen(satlist).readlines() 
    tles = [item.strip() for item in tles]
    tles = [(tles[i],tles[i+1],tles[i+2]) for i in xrange(0,len(tles)-2,3)]
    
    for tle in tles:

        try:
            sat = ephem.readtle(tle[0], tle[1], tle[2])
            rt, ra, tt, ta, st, sa = observer.next_pass(sat)

            if rt is not None and st is not None:
                sat.compute(observer)

                if Time >= ephem.localtime(st) and Time <= ephem.localtime(rt):
                    Title = tle[0]
                    sat_alt.append(np.rad2deg(sat.alt))
                    sat_az.append(np.rad2deg(sat.az))
                    original_name.append(tle[0])

                    text2 = Title.rsplit(')', 1)[0]
                    NamePRN = text2.rsplit('(', 1)[1]
                    sat_name.append(NamePRN)

                    if contain_GPS in Title:
                        sat_type.append(1)
                    elif contain_GLONASS in Title:
                        sat_type.append(2)
                    elif contain_GALILEO in Title:
                        sat_type.append(3)

        except ValueError as e:
            print e
return sat_name, sat_alt, sat_az, sat_type, original_name

Calling the function

GPS_list = 'https://celestrak.org/NORAD/elements/gps-ops.txt'
Downloaded = urllib2.urlopen(GPS_list).readlines()
TimeNow = datetime.datetime.now()
print Time2SecondsOfDay(TimeNow), 'time: ',TimeNow

names,elevs,azs,types,originalName = 
GetPositionELAZ(TimeNow,Longitude,Latitude,1,Downloaded)
$\endgroup$
8
  • $\begingroup$ I hope see the same vehicles in the same time, but I only see others vehicles, but I don't know why is not the same. for compare a few vehicles only I show the satellites over 45º. $\endgroup$ Apr 18, 2017 at 15:54
  • $\begingroup$ This probably isn't helpful, but you have your latitude and longitude reversed. Let me know if that does/doesn't help. $\endgroup$
    – user7073
    Apr 18, 2017 at 16:24
  • $\begingroup$ thank you @barrycarter I had the wrong latitude/longitude, I updated the results (and the question), although I find some differences in the results yet. $\endgroup$ Apr 18, 2017 at 17:18
  • $\begingroup$ If I'm not mistaken, the TLE information will give you the orbit. From that, you get compute the position and velocity of the spacecraft in the ECI frame. However, you need to convert that to the ECEF frame in order to compute the altitude and elevation information, and I don't see that happening in your code. That said, if there is a library that works well, why do you need to recode such a feature, and not just use the library that works? $\endgroup$
    – ChrisR
    Apr 18, 2017 at 20:03
  • $\begingroup$ Also, what's the code from sat.compute()? Could the reference frame conversion happen there? $\endgroup$
    – ChrisR
    Apr 18, 2017 at 20:04

1 Answer 1

3
$\begingroup$

I figured out,It was something wrong in the code: The function observer.next_pass(sat), return a vector of 5 values

0 Rise time
1 Rise azimuth
2 Maximum altitude time
3 Maximum altitude
4 Set time
5 Set azimuth

So the instruction if Time >= ephem.localtime(st) and Time <= ephem.localtime(rt): I was using a innecesary comparation because you already compute the time with the instruction observer.date = Time also I was comparing with the localtime, the best option is have all in UTC time.
and a advice with satellites don't use datetime.now(), use datetime.utcnow(), and later you can add the GTM of your region.

$\endgroup$
0

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.