Skip to content

Commit

Permalink
BUGFIX for irregular NAV3 files
Browse files Browse the repository at this point in the history
add galileo nav plots
  • Loading branch information
scivision committed Oct 3, 2018
1 parent 9c1df83 commit 5508d12
Show file tree
Hide file tree
Showing 4 changed files with 61 additions and 18 deletions.
12 changes: 10 additions & 2 deletions georinex/keplerian.py
Original file line number Diff line number Diff line change
Expand Up @@ -31,15 +31,23 @@ def keplerian2ecef(sv: xarray.DataArray) -> Tuple[np.ndarray, np.ndarray, np.nda
n = n0 + sv['DeltaN'] # corrected mean motion

# from GPS Week 0
t0 = datetime(1980, 1, 6) + timedelta(weeks=sv['GPSWeek'][0].astype(int).item())
if sv.svtype[0] == 'E':
weeks = sv['GALWeek'].values - 1024
# TODO: need to add N*4096 for year 2058 and beyond
elif sv.svtype[0] == 'G':
weeks = sv['GPSWeek'].values
else:
raise ValueError(f'Unknown system type {sv.svtype[0]}')

T0 = [datetime(1980, 1, 6) + timedelta(weeks=week) for week in weeks]

tk = np.empty(sv['time'].size, dtype=float)

# FIXME: so ugly...
# time elapsed since reference epoch
# seems to be a bug in MyPy, this line computes "correctly"

for i, (t1, t2) in enumerate(zip(sv['time'], sv['Toe'])):
for i, (t0, t1, t2) in enumerate(zip(T0, sv['time'], sv['Toe'])):
tsv = datetime.utcfromtimestamp(t1.item()/1e9)
toe = timedelta(seconds=t2.values.astype(int).item()) + t0 # type: ignore # noqa
tk[i] = (tsv - toe).total_seconds() # type: ignore # noqa
Expand Down
56 changes: 41 additions & 15 deletions georinex/nav3.py
Original file line number Diff line number Diff line change
Expand Up @@ -3,6 +3,7 @@
import xarray
import logging
import numpy as np
import math
from io import BytesIO
from datetime import datetime
from .io import opener, rinexinfo
Expand All @@ -11,6 +12,7 @@
# constants
STARTCOL3 = 4 # column where numerical data starts for RINEX 3
Nl = {'C': 7, 'E': 7, 'G': 7, 'J': 7, 'R': 3, 'S': 3} # number of additional SV lines
Lf = 19 # string length per field


def rinexnav3(fn: Path,
Expand All @@ -24,7 +26,6 @@ def rinexnav3(fn: Path,
The "eof" stuff is over detection of files that may or may not have a trailing newline at EOF.
"""
Lf = 19 # string length per field

fn = Path(fn).expanduser()

Expand Down Expand Up @@ -96,19 +97,8 @@ def rinexnav3(fn: Path,
if tu.size != t[svi].size:
logging.warning(f'duplicate times detected on SV {sv}, using first of duplicate times')
""" I have seen that the data rows match identically when times are duplicated"""
# %% check for optional GPS "fit interval" presence
cf = fields[sv[0]]
testread = np.genfromtxt(BytesIO(raws[svi[0]].encode('ascii')), delimiter=Lf)
# %% patching for Spare entries, some receivers include, and some don't include...
if sv[0] == 'G' and len(cf) == testread.size + 1:
cf = cf[:-1]
elif sv[0] == 'C' and len(cf) == testread.size - 1:
cf.insert(20, 'spare')
elif sv[0] == 'E' and len(cf) == testread.size - 1:
cf.insert(22, 'spare')

if testread.size != len(cf):
raise ValueError(f'{sv[0]} NAV data @ {tu} is not the same length as the number of fields.')
cf = _sparefields(fields[sv[0]], sv[0], raws[svi[0]])

darr = np.empty((svi.size, len(cf)))

Expand Down Expand Up @@ -162,11 +152,42 @@ def _time(ln: str) -> Optional[datetime]:
return None


def _sparefields(cf: List[str], sys: str, raw: str) -> List[str]:
"""
check for optional spare fields, or GPS "fit interval" field
"""
numval = math.ceil(len(raw) / Lf) # need this for irregularly defined files
# %% patching for Spare entries, some receivers include, and some don't include...
if sys == 'G' and len(cf) == numval + 1:
cf = cf[:-1]
elif sys == 'C' and len(cf) == numval - 1:
cf.insert(20, 'spare')
elif sys == 'E':
if numval == 28:
cf = cf[:-3]
elif numval == 27:
cf = cf[:22] + cf[23:-3]

if numval != len(cf):
raise ValueError(f'System {sys} NAV data is not the same length as the number of fields.')

return cf


def _newnav(ln: str, sv: str) -> List[str]:

if sv.startswith('G'):
"""
ftp://igs.org/pub/data/format/rinex303.pdf page A-23 - A-24
ftp://igs.org/pub/data/format/rinex303.pdf
pages:
G: A23-A24
E: A25-A28
R: A29-A30
J: A31-A32
C: A33-A34
S: A35-A36
I: A37-A39
"""
fields = ['SVclockBias', 'SVclockDrift', 'SVclockDriftRate',
'IODE', 'Crs', 'DeltaN', 'M0',
Expand Down Expand Up @@ -211,8 +232,13 @@ def _newnav(ln: str, sv: str) -> List[str]:
'Toe', 'Cic', 'Omega0', 'Cis',
'Io', 'Crc', 'omega', 'OmegaDot',
'IDOT', 'DataSrc', 'GALWeek',
'spare0',
'SISA', 'health', 'BGDe5a', 'BGDe5b',
'TransTime']
'TransTime',
'spare1', 'spare2', 'spare3']
assert len(fields) == 31
elif sv.startswith('I'):
raise NotImplementedError('please raise GitHub issue to request IRNSS NAV3')
else:
raise ValueError(f'Unknown SV type {sv[0]}')

Expand Down
9 changes: 9 additions & 0 deletions georinex/plots.py
Original file line number Diff line number Diff line change
Expand Up @@ -112,6 +112,15 @@ def navtimeseries(nav: xarray.Dataset):

if ((lat < -57) | (lat > 57)).any():
logging.warning('GPS inclination ~ 55 degrees')
elif sv[0] == 'E':
ecef = keplerian2ecef(nav.sel(sv=sv))
lat, lon, alt = pm.ecef2geodetic(*ecef)

if ((alt < 23e6) | (alt > 24e6)).any():
logging.warning('unrealistic Galileo satellite altitudes')

if ((lat < -57) | (lat > 57)).any():
logging.warning('Galileo inclination ~ 56 degrees')

else:
continue
Expand Down
2 changes: 1 addition & 1 deletion setup.cfg
Original file line number Diff line number Diff line change
@@ -1,6 +1,6 @@
[metadata]
name = georinex
version = 1.6.7.2
version = 1.6.7.3
author = Michael Hirsch, Ph.D.
author_email = [email protected]
description = Python RINEX 2/3 NAV/OBS reader with speed and simplicity.
Expand Down

0 comments on commit 5508d12

Please sign in to comment.