Sto provando a scrivere uno script Python di base che traccia un determinato satellite, definito con tle, da una determinata posizione. Non sono una persona asto/orbitale ma sto cercando di diventare più intelligente.Tracciamento satellitare Python con spg4, posizioni di pyephem non corrispondenti
Sto incontrando un problema in cui i diversi modelli che sto usando mi danno risposte di posizione molto diverse. Ho provato con: pyEphem SPG4 prevedere (chiamata di sistema exec dallo script)
I satelliti sto testando con sono l'ISS e directv10 (una fissa, una MOVING- con il monitoraggio internet a disposizione per la verifica):
0 Direct10
1 31862U 07032A 13099.15996183 -.00000126 00000-0 10000-3 0 1194
2 31862 000.0489 046.9646 0000388 001.7833 103.5813 01.00271667 21104
0 ISS
1 25544U 98067A 13112.50724749 .00016717 00000-0 10270-3 0 9148
2 25544 51.6465 24.5919 0009906 171.1474 188.9854 15.52429950 26067
Ho modificato la fonte di previsione per fornirmi la posizione eci in modo da poterlo utilizzare per conoscere la posizione reale. Ho anche dato l'intervallo az, el per l'utilizzo per verificare le osservazioni. Sto usando spg4 per ottenere la posizione reale. Per la posizione osservata, sto usando PyEphem.
sto ottenendo la posizione ECEF da SPG4 con:
def get_real(epoch, sv):
satellite = twoline2rv(sv.tle1, sv.tle2, wgs84)
#epoch = time.time()
obsTime = datetime.datetime.utcfromtimestamp(epoch)
position, velocity = satellite.propagate(obsTime.year,
obsTime.month,
obsTime.day,
obsTime.hour,
obsTime.minute,
obsTime.second)
x = position[0]
y = position[1]
z = position[2]
x *= 1000
y *= 1000
z *= 1000
Il mio codice per osservazioni basate pyephem è:
def get_ob(epoch, sv, obsLoc):
site = ephem.Observer()
site.lon = str(obsLoc.lat) # +E -104.77 here
site.lat = str(obsLoc.lon) # +N 38.95 here
site.elevation = obsLoc.alt # meters 0 here
#epoch = time.time()
site.date = datetime.datetime.utcfromtimestamp(epoch)
sat = ephem.readtle(sv.name,sv.tle1,sv.tle2)
sat.compute(site)
az = degrees(sat.az)
el = degrees(sat.alt)
#range in m
range = sat.range
sat_lat = degrees(sat.sublat)
sat_long = degrees(sat.sublong)
# elevation of sat in m
sat_elev = sat.elevation
#TODO: switch to using az,el,range for observed location calculation
#x, y, z = aer2ecef(az,el,range,38.95,-104.77,80/1000)
x,y,z = llh2ecef(sat_lat, sat_long, sat_elev)
La conversione llh2ecef:
def llh2ecef (flati,floni, altkmi):
# lat,lon,height to xyz vector
#
# input:
# flat geodetic latitude in deg
# flon longitude in deg
# altkm altitude in km
# output:
# returns vector x 3 long ECEF in km
dtr = pi/180.0;
flat = float(flati);
flon = float(floni);
altkm = float(altkmi);
clat = cos(dtr*flat);
slat = sin(dtr*flat);
clon = cos(dtr*flon);
slon = sin(dtr*flon);
rrnrm = radcur (flat);
rn = rrnrm[1];
re = rrnrm[0];
ecc1 = ecc;
esq1 = ecc1*ecc1
x = (rn + altkm) * clat * clon;
y = (rn + altkm) * clat * slon;
z = ((1-esq1)*rn + altkm) * slat;
return x,y,z
aer2ecef:
def aer2ecef(azimuthDeg, elevationDeg, slantRange, obs_lat, obs_long, obs_alt):
#site ecef in meters
sitex, sitey, sitez = llh2ecef(obs_lat,obs_long,obs_alt)
#some needed calculations
slat = sin(radians(obs_lat))
slon = sin(radians(obs_long))
clat = cos(radians(obs_lat))
clon = cos(radians(obs_long))
azRad = radians(azimuthDeg)
elRad = radians(elevationDeg)
# az,el,range to sez convertion
south = -slantRange * cos(elRad) * cos(azRad)
east = slantRange * cos(elRad) * sin(azRad)
zenith = slantRange * sin(elRad)
x = (slat * clon * south) + (-slon * east) + (clat * clon * zenith) + sitex
y = (slat * slon * south) + (clon * east) + (clat * slon * zenith) + sitey
z = (-clat * south) + (slat * zenith) + sitez
return x, y, z
Quando paragono e piazzano le posizioni su un globo 3D (usando le posizioni ecef), ricevo risposte dappertutto. La posizione eci di previsione (convertita in ecef) corrisponde a quello che vedo sui siti Web di tracciamento ISS (http://www.n2yo.com/?s=25544)
Il risultato di get_real() è molto diverso in scala e posizione. Il risultato da get_ob() sono proprio in scala, ma sbagliato in posizione sul globo
esempio i risultati:
prevedrebbe basano:
sv: ISS predict observed response @ epoch: 1365630559.000000 : [111.485527, -69.072949, 12351.471383]
sv: ISS predict aer2ecef position(m) @ epoch: 1365630559.000000 : [4731598.706291642, 1844098.7384999825, -4521102.9225004213]
sv: ISS predict ecef position(m) @ epoch: 1365630559.000000 : [-3207559.6840419229, -3937040.5048992992, -4521102.9110000003]
sv: ISS predict ecef2llh(m) @ epoch: 1365630559.000000 : [-41.67839724680753, -129.170165912171, 6792829.6884068651]
sv: Direct10 predict observed response @ epoch: 1365630559.000000 : [39.692138, -49.219935, 46791.914833]
sv: Direct10 predict aer2ecef position(m) @ epoch: 1365630559.000000 : [28401835.38849232, 31161334.784188181, 3419.5400331273049]
sv: Direct10 predict ecef position(m) @ epoch: 1365630559.000000 : [-9348629.6463202238, -41113211.570621684, 3419.8620000000005]
sv: Direct10 predict ecef2llh(m) @ epoch: 1365630559.000000 : [0.0046473273713214715, -102.81051792373036, 42156319.281573996]
basato su Python:
sv: ISS ephem observed response @ epoch: 1365630559.000000 : [344.067992722211, -72.38297754053431, 12587123.0][degrees(sat.az), degrees(sat.alt), sat.range]
sv: ISS ephem llh location(m) @ epoch: 1365630559.000000 : [-41.678271938092195, -129.16682754513502, 421062.90625][degrees(sat.sublat0, degrees(sat.sublong), sat.elevation]
sv: ISS ephem xyz location(m) @ epoch: 1365630559.000000 :[-201637.5647039332, -247524.53652043006, -284203.56557438202][llh2ecef(lat,long,elev)]
sv: ISS spg84 ecef position(m) @ epoch: 1365630559.000000 : [4031874.0758277094, 3087193.8810081254, -4521293.538866323]
sv: ISS spg84 ecef2llh(m) @ epoch: 1365630559.000000 : [-41.68067424524357, 37.4411722245808, 6792812.8704163525]
sv: Direct10 ephem observed response @ epoch: 1365630559.000000 : [320.8276456938389, -19.703680198781303, 43887572.0][degrees(sat.az), degrees(sat.alt), sat.range]
sv: Direct10 ephem llh location(m) @ epoch: 1365630559.000000 : [0.004647324660923812, -102.8070784813048, 35784688.0][degrees(sat.sublat0, degrees(sat.sublong), sat.elevation]
sv: Direct10 ephem xyz location(m) @ epoch: 1365630559.000000 :[-7933768.6901137345, -34900655.02490133, 2903.0498773286708][llh2ecef(lat,long,elev)]
sv: Direct10 spg84 ecef position(m) @ epoch: 1365630559.000000 : [18612307.532456037, 37832170.97306267, -14060.29781505302]
sv: Direct10 spg84 ecef2llh(m) @ epoch: 1365630559.000000 : [-0.019106864351793953, 63.80418030988552, 42156299.077687643]
L'az , el e range non corrispondono tra le due osservazioni. Le posizioni non corrispondono per le posizioni "vere". (Il lat e long non è all'altezza di una conversione ecef2llh
Rispetto al tracker basato sul web, noto che le "vere" posizioni corrispondono al sito Web. Per directv10, pyEphem corrisponde all'azimut e all'elevazione - ma non per ISS
Quando li ho tracciati su globo, la posizione di eci "vero" di previsione è nel punto giusto - corrisponde al sito Web di tracker). La posizione dell'ecef spg84 (che pensavo dovesse essere la stessa di quanto si possa prevedere, si trova dall'altra parte del globo: il predire posizione "osservata" è vicino alla posizione spg84.Il pyEphem è completamente fuori quota e non viene visualizzato (troppo basso, terra interna)
Quindi la mia domanda è dove sto usando i modelli Python errati?La mia comprensione era che spg84 propagate() call, dovrebbe restituire la posizione exec del satellite in metri. Avrei pensato che avrebbe dovuto corrispondere alla posizione di previsione dopo la conversione di eci2efec. Anch'io mi sarei aspettato che il match fosse poi llh2ecef() quando si usava sat.sublat, sat.sublong, sat.elevation.
Come ho già detto, sono nuovo a tutto ciò che orbita, quindi sono sicuro che sto facendo un semplice errore matematico o qualcosa di simile. Ho cercato di google e cercare risposte, esempi e tutorial il più possibile ma nulla ha aiutato finora (ho provato più metodi ecef2llh e llh2ecef per provare a risolvere quei bug
Qualsiasi suggerimento, consiglio, puntatori nella giusta direzione sarebbe molto apprezzato, posso postare/inviare il codice completo che sto usando se fosse utile per qualcuno, ho cercato di assicurarmi di aver postato le parti importanti qui e non volevo farlo (già molto) lungo post e più a lungo.
Grazie per l'aiuto.
Aaron
UPDATE:
Ho trovato almeno una parte del problema. spg84.propagate() restituisce la posizione in ECI, non ECEF. Esegui rapidamente eci2ecef e si allinea perfettamente con la risposta di previsione.
mi sembra sempre di trovare soluzioni dopo la pubblicazione aiuto;)
adesso ho bisogno di capire cosa sta succedendo con le posizioni di osservazione. Questo si riduce a: Come posso prendere il risultato da pyEphem.compute() e ottenere la posizione ecef per il satellite? Preferisci farlo con az, el, valori di intervallo, non latitudine, longitudine, elevazione.
Sto indovinando l'errore nella mia chiamata aer2ecef.
Grazie.
UPDATE 2:
Got l'osservazione per allinearsi con la posizione di "vero". Sembra che abbia avuto un problema con le unità. Il codice di lavoro:
az = degrees(sat.az)
el = degrees(sat.alt)
#range in km
range = sat.range
sat_lat = degrees(sat.sublat)
sat_long = degrees(sat.sublong)
# elevation of sat in km
sat_elev = sat.elevation
#x, y, z = aer2ecef(az,el,range,obsLoc.lat,obsLoc.long,obsLoc.alt)
x,y,z = llh2ecef(sat_lat, sat_long, sat_elev/1000)
x *= 1000
y *= 1000
z *= 1000
return x,y,z
Ora solo bisogno aer2ecef() per tornare posizione corretta ...
Sto per chiudere questa una nuova domanda aperta. Da quando ho risolto due dei problemi, immagino sia più facile fare una sola, più concisa domanda ora. Voglio lasciare questo però per essere utile a chiunque altro che potrebbe incontrare gli stessi problemi che ho fatto con le conversioni. Inoltre, se qualcuno ha più feedback, per favore fatemelo sapere. Sono sempre alla ricerca di miglioramenti. –
Potresti dirmi dove si ottiene la funzione radcur? – Brian