mirrors/gnss-sdr
1
0
Fork 0
mirror of https://github.com/gnss-sdr/gnss-sdr synced 2025-09-04 03:47:59 +00:00
Code Issues Releases Wiki Activity

skyplot improvements

Browse Source 
Add new --format optional argument. Allowed options: pdf, eps, png, svg. Default: pdf
Add new --version / -v optional arguments, show program version and exit
Improve --help / -h message
Improve error handling
Generate PDF files with embedded fonts, ready for journal submission
Improve argument handling
Improve documentation
Add explicit z-order control for cleaner display of labels in multiconstellation files
This commit is contained in:
Carles Fernandez
2025-08-22 12:42:36 +02:00
parent 15a7a9bae3
commit 1d31c8427d
2 changed files with 248 additions and 191 deletions
Download Patch File Download Diff File Expand all files Collapse all files

59
utils/skyplot/README.md
View File

@@ -16,18 +16,19 @@ showing satellite visibility over time.
## Features ## Features
- Processes RINEX navigation files. - Processes RINEX navigation files.
- Calculates satellite positions using broadcast ephemeris.
- Plots satellite tracks in azimuth-elevation coordinates.
- Customizable observer location.
- Color-codes satellites by constellation (GPS, Galileo, GLONASS, BeiDou).
- Elevation mask set to 5°, configurable via the `--elev-mask` option.
- Outputs high-quality image in PDF format. EPS, PNG, and SVG formats are also
available via the `--format` option.
- Non-interactive mode for CI jobs with the `--no-show` option.
- Constellations to plot can be configured via the `--system` option.
- Optionally uses an OBS file to limit plot to the receiver observation time - Optionally uses an OBS file to limit plot to the receiver observation time
(`--use-obs`). (`--use-obs`).
- When enabled, the tool looks for a matching file by replacing the last - When enabled, the tool looks for a matching file by replacing the last
character of the NAV filename with O/o and uses it if found. character of the NAV filename with `O`/`o` and uses it if found.
- Calculates satellite positions using broadcast ephemeris.
- Plots satellite tracks in azimuth-elevation coordinates.
- Elevation mask set to 5°, configurable via the `--elev-mask` flag.
- Color-codes satellites by constellation (GPS, Galileo, GLONASS, BeiDou).
- Constellations to plot can be configured via the `--system` flag.
- Customizable observer location.
- Outputs high-quality image in PDF format.
- Non-interactive mode for CI jobs (with `--no-show` flag).
## Requirements ## Requirements
@@ -41,21 +42,27 @@ showing satellite visibility over time.
### Basic Command ### Basic Command
``` ```
./skyplot.py <RINEX_FILE> [LATITUDE] [LONGITUDE] [ALTITUDE] [--use-obs] [--elev-mask] [--system ...] [--no-show] ./skyplot.py <RINEX_FILE> [LATITUDE] [LONGITUDE] [ALTITUDE]
[--elev-mask ELEV_MASK]
[--format {pdf,eps,png,svg}]
[--no-show]
[--system SYSTEM [SYSTEM ...]]
[--use-obs]
``` ```
### Arguments ### Arguments
| Argument | Type | Units | Description | Default | | Argument | Type | Units | Description | Default |
| ---------------- | -------- | ----------- | ---------------------- | -------- | | ---------------- | -------- | ----------- | ------------------------ | -------- |
| `RINEX_NAV_FILE` | Required | - | RINEX nav file path | - | | `RINEX_NAV_FILE` | Required | - | RINEX nav file path | - |
| `LATITUDE` | Optional | degrees (°) | North/South position | 41.275°N | | `LATITUDE` | Optional | degrees (°) | North/South position | 41.275°N |
| `LONGITUDE` | Optional | degrees (°) | East/West position | 1.9876°E | | `LONGITUDE` | Optional | degrees (°) | East/West position | 1.9876°E |
| `ALTITUDE` | Optional | meters (m) | Height above sea level | 80.0 m | | `ALTITUDE` | Optional | meters (m) | Height above sea level | 80.0 m |
| `--use-obs` | Optional | - | Use RINEX obs data | - | | `--elev-mask` | Optional | degrees (°) | Elevation mask | 5° |
| `--elev-mask` | Optional | degrees (°) | Elevation mask | 5° | | `--format` | Optional | - | Output {pdf,eps,png,svg} | pdf |
| `--system` | Optional | - | Systems to plot | All | | `--no-show` | Optional | - | Do not show plot | - |
| `--no-show` | Optional | - | Do not show plot | - | | `--system` | Optional | - | Systems to plot | All |
| `--use-obs` | Optional | - | Use RINEX obs data | - |
### Examples ### Examples
@@ -84,6 +91,10 @@ showing satellite visibility over time.
``` ```
./skyplot.py brdc0010.22n -33.4592 -70.6453 520.0 --system G E ./skyplot.py brdc0010.22n -33.4592 -70.6453 520.0 --system G E
``` ```
- Get a PNG file:
```
./skyplot.py brdc0010.22n -33.4592 -70.6453 520.0 --format png
```
## Output ## Output
@@ -91,8 +102,12 @@ The script generates a PDF file named `skyplot_<RINEX_FILE>.pdf` (with dots in
`<RINEX_FILE>` replaced by `_`) with: `<RINEX_FILE>` replaced by `_`) with:
- Satellite trajectories over all epochs in the file. - Satellite trajectories over all epochs in the file.
- NAV file - ephemeris time range (default) - NAV file - ephemeris time range (default).
- Receiver observation if `--use-obs` is specified and OBS file is found - Receiver observation if `--use-obs` is specified and OBS file is found.
- Color-coded by constellation. - Color-coded by constellation.
- Observer location in title. - Observer location in title.
- Time range in footer. - Time range in footer.
- Embedded fonts that display consistently across all systems and generate
publication-ready figures.
- EPS, PNG, and SVG output formats available via `--format eps`, `--format png`,
and `--format svg`.

380
utils/skyplot/skyplot.py
View File

@@ -2,10 +2,17 @@
""" """
skyplot.py skyplot.py
Reads a RINEX navigation file and generates a skyplot. Optionally, a RINEX observation file can Reads a RINEX navigation file and generates a skyplot. Optionally, a RINEX
also be read to match the skyplot to the receiver processing time. observation file can also be read to match the skyplot to the receiver
processing time.
Usage: python skyplot.py <RINEX_NAV_FILE> [observer_lat] [observer_lon] [observer_alt] [--use-obs] Usage:
skyplot.py <RINEX_NAV_FILE> [observer_lat] [observer_lon] [observer_alt]
[--elev-mask ELEV_MASK]
[--format {pdf,eps,png,svg}]
[--no-show]
[--system SYSTEM [SYSTEM ...]]
[--use-obs]
----------------------------------------------------------------------------- -----------------------------------------------------------------------------
@@ -25,9 +32,15 @@ from datetime import datetime, timedelta
from math import atan2, cos, sin, sqrt from math import atan2, cos, sin, sqrt
from pathlib import Path from pathlib import Path
import matplotlib.pyplot as plt try:
import numpy as np import matplotlib.pyplot as plt
import numpy as np
except ImportError:
print("Error: This script requires matplotlib and numpy.")
print("Install them with: pip install matplotlib numpy")
sys.exit(1)
__version__ = "1.0.0"
def _read_obs_time_bounds(obs_path): def _read_obs_time_bounds(obs_path):
"""Return (start_time, end_time) from a RINEX 2/3 OBS file by scanning epoch lines. """Return (start_time, end_time) from a RINEX 2/3 OBS file by scanning epoch lines.
@@ -360,9 +373,15 @@ def ecef_to_az_el(x, y, z, obs_lat, obs_lon, obs_alt):
def plot_satellite_tracks(satellites, obs_lat, obs_lon, obs_alt, def plot_satellite_tracks(satellites, obs_lat, obs_lon, obs_alt,
footer_text=None, filename=None, footer_text=None, filename=None,
show_plot=True, start_time=None, show_plot=True, start_time=None,
end_time=None, elev_mask=5.0): end_time=None, elev_mask=5.0,
output_format="pdf"):
"""Plot trajectories for all visible satellites""" """Plot trajectories for all visible satellites"""
plt.rcParams["font.family"] = "Times New Roman" plt.rcParams['pdf.fonttype'] = 42 # TrueType fonts
plt.rcParams['ps.fonttype'] = 42 # TrueType fonts
plt.rcParams['font.family'] = 'serif'
plt.rcParams['font.serif'] = ['Times New Roman', 'Times', 'DejaVu Serif']
plt.rcParams['mathtext.fontset'] = 'dejavuserif' # For math text
plt.rcParams['svg.fonttype'] = 'none' # Make SVG text editable
fig = plt.figure(figsize=(8, 8)) fig = plt.figure(figsize=(8, 8))
ax = fig.add_subplot(111, projection='polar') ax = fig.add_subplot(111, projection='polar')
ax.tick_params(labelsize=16, pad=7) ax.tick_params(labelsize=16, pad=7)
@@ -442,13 +461,13 @@ def plot_satellite_tracks(satellites, obs_lat, obs_lon, obs_alt,
for az_seg, el_seg in segments: for az_seg, el_seg in segments:
theta = np.radians(az_seg) theta = np.radians(az_seg)
r = 90 - np.array(el_seg) r = 90 - np.array(el_seg)
ax.plot(theta, r, '-', color=color, alpha=0.7, linewidth=2.5) ax.plot(theta, r, '-', color=color, alpha=0.7, linewidth=2.5, zorder=1)
# Arrow at end # Arrow at end
if len(theta) >= 2: if len(theta) >= 2:
dx = theta[-1] - theta[-2] dx = theta[-1] - theta[-2]
dy = r[-1] - r[-2] dy = r[-1] - r[-2]
arrow_length_factor = 1.3 arrow_length_factor = 1.8
extended_theta = theta[-2] + dx * arrow_length_factor extended_theta = theta[-2] + dx * arrow_length_factor
extended_r = r[-2] + dy * arrow_length_factor extended_r = r[-2] + dy * arrow_length_factor
ax.annotate('', ax.annotate('',
@@ -461,13 +480,15 @@ def plot_satellite_tracks(satellites, obs_lat, obs_lon, obs_alt,
'linewidth': 1.5, 'linewidth': 1.5,
'shrinkA': 0, 'shrinkA': 0,
'shrinkB': 0 'shrinkB': 0
}) },
zorder=2)
# Label at midpoint of the segment # Label at midpoint of the segment
mid_idx = len(theta)//2 mid_idx = len(theta)//2
ax.text(theta[mid_idx], r[mid_idx], prn, ax.text(theta[mid_idx], r[mid_idx], prn,
fontsize=12, ha='center', va='center', fontsize=12, ha='center', va='center',
bbox={"facecolor": "white", "alpha": 0.8, "pad": 2}) bbox={"facecolor": "white", "alpha": 0.8, "pad": 2},
zorder=3)
# Legend for present systems # Legend for present systems
legend_elements = [ legend_elements = [
@@ -503,11 +524,11 @@ def plot_satellite_tracks(satellites, obs_lat, obs_lon, obs_alt,
if filename: if filename:
filename_no_path = Path(filename).name filename_no_path = Path(filename).name
filename_no_dots = filename_no_path.replace('.', '_') filename_no_dots = filename_no_path.replace('.', '_')
output_name = f"skyplot_{filename_no_dots}.pdf" output_name = f"skyplot_{filename_no_dots}.{output_format}"
else: else:
output_name = "skyplot.pdf" output_name = f"skyplot.{output_format}"
plt.savefig(output_name, format='pdf', bbox_inches='tight') plt.savefig(output_name, format=output_format, bbox_inches='tight')
print(f"Image saved as {output_name}") print(f"Image saved as {output_name}")
if show_plot: if show_plot:
plt.show() plt.show()
@@ -517,177 +538,198 @@ def plot_satellite_tracks(satellites, obs_lat, obs_lon, obs_alt,
def main(): def main():
"""Generate the skyplot""" """Generate the skyplot"""
# Set system names and codes
system_name_to_code = {
'GPS': 'G',
'GLONASS': 'R',
'GALILEO': 'E',
'BEIDOU': 'C',
'QZSS': 'J',
'IRNSS': 'I',
'SBAS': 'S'
}
# Set up argument parser
parser = argparse.ArgumentParser(
description='Generate GNSS skyplot from RINEX navigation file',
add_help=False
)
# Add the no-show flag
parser.add_argument(
'--no-show',
action='store_true',
help='Run without displaying plot window'
)
# Add the use-obs flag
parser.add_argument(
'--use-obs',
action='store_true',
help='Use corresponding RINEX observation file to bound the skyplot to the receiver time window'
)
# Add the elev-mask flag.
parser.add_argument(
'--elev-mask', type=float, default=5.0,
help='Elevation mask in degrees for plotting satellite tracks (default: 5°)'
)
# Add the system flag.
parser.add_argument(
'--system',
nargs='+',
help='Only plot satellites from these systems (e.g., G R E or GPS Galileo GLONASS)'
)
# Parse known args (this ignores other positional args)
args, remaining_args = parser.parse_known_args()
# Handle help manually
if '-h' in remaining_args or '--help' in remaining_args:
print("""
Usage: python skyplot.py <RINEX_FILE> [LATITUDE] [LONGITUDE] [ALTITUDE] [--use-obs] [--elev-mask] [--system ...] [--no-show]
Example:
python skyplot.py brdc0010.22n 41.275 1.9876 80.0 --use-obs --no-show --elev-mask=10 --system GPS Galileo
""")
sys.exit(0)
if len(remaining_args) < 1:
print("Error: RINEX file required")
sys.exit(1)
filename = remaining_args[0]
# Default observer location (Castelldefels, Barcelona, Spain)
obs_lat = np.radians(41.2750)
obs_lon = np.radians(1.9876)
obs_alt = 80.0
# Override with command line arguments if provided
if len(remaining_args) >= 4:
try:
obs_lat = np.radians(float(remaining_args[1]))
obs_lon = np.radians(float(remaining_args[2]))
if len(remaining_args) >= 5:
obs_alt = float(remaining_args[3])
except ValueError:
print("Invalid observer coordinates. Using defaults.")
# Read RINEX file
print(f"Reading {filename}...")
try: try:
satellites = read_rinex_nav(filename) # Set system names and codes
except FileNotFoundError: system_name_to_code = {
print(f"Error: File '{filename}' not found.") 'GPS': 'G',
return 'GLONASS': 'R',
'GALILEO': 'E',
'BEIDOU': 'C',
'QZSS': 'J',
'IRNSS': 'I',
'SBAS': 'S'
}
if not satellites: # Set up argument parser
print("No satellite data found in the file.") parser = argparse.ArgumentParser(
return description='Generate a GNSS skyplot from a RINEX navigation file',
epilog="Example: skyplot.py brdc0010.22n -33.4592 -70.6453 520.0 --format png --system G E --elev-mask 10 --no-show"
)
if args.system: # Positional arguments
systems_upper = set() parser.add_argument(
for s in args.system: 'filename',
s_upper = s.upper() help='RINEX navigation file path'
if s_upper in system_name_to_code: )
systems_upper.add(system_name_to_code[s_upper])
else:
systems_upper.add(s_upper) # Assume user passed the code
satellites = {prn: eph_list for prn, eph_list in satellites.items() if prn[0].upper() in systems_upper} parser.add_argument(
'lat', nargs='?', type=float, default=41.2750,
help='Observer latitude in degrees (default: 41.275° N)'
)
parser.add_argument(
'lon', nargs='?', type=float, default=1.9876,
help='Observer longitude in degrees (default: 1.9876° E)'
)
parser.add_argument(
'alt', nargs='?', type=float, default=80.0,
help='Observer altitude in meters (default: 80.0 m)'
)
# Optional arguments
parser.add_argument(
'--elev-mask',
type=float,
default=5.0,
help='Elevation mask in degrees for plotting satellite tracks (default: 5°)'
)
parser.add_argument(
'--format',
type=str,
default="pdf",
choices=["pdf", "eps", "png", "svg"],
help='Output file format for plot (default: pdf)'
)
parser.add_argument(
'--no-show',
action='store_true',
help='Run without displaying plot window'
)
parser.add_argument(
'--system',
nargs='+',
help='Only plot satellites from these systems (e.g., G R E or GPS GLONASS Galileo)'
)
parser.add_argument(
'--use-obs',
action='store_true',
help='Use corresponding RINEX observation file to bound the skyplot to the receiver time window'
)
parser.add_argument(
"-v", "--version",
action="version",
version=f"%(prog)s {__version__}",
help="Show program version and exit"
)
# Parse all arguments with full validation
args = parser.parse_args()
# Convert coordinates to radians
obs_lat = np.radians(args.lat)
obs_lon = np.radians(args.lon)
obs_alt = args.alt
filename = args.filename
# Read RINEX file
print(f"Reading {filename}...")
try:
satellites = read_rinex_nav(filename)
except FileNotFoundError:
print(f"Error: File '{filename}' not found.")
return 1
if not satellites: if not satellites:
print(f"No satellites found for systems: {', '.join(sorted(systems_upper))}") print("No satellite data found in the file.")
return return 1
# Print summary information if args.system:
all_epochs = sorted(list(set( systems_upper = set()
e['epoch'] for prn, ephemerides in satellites.items() for e in ephemerides for s in args.system:
))) s_upper = s.upper()
print("\nFile contains:") if s_upper in system_name_to_code:
print(f"- {len(satellites)} unique satellites") systems_upper.add(system_name_to_code[s_upper])
print(f"- {len(all_epochs)} unique epochs") else:
print(f"- From {all_epochs[0]} to {all_epochs[-1]}") systems_upper.add(s_upper) # Assume user passed the code
# Calculate and print satellite counts by system satellites = {prn: eph_list for prn, eph_list in satellites.items() if prn[0].upper() in systems_upper}
system_counts = {}
for prn in satellites:
system = prn[0]
system_counts[system] = system_counts.get(system, 0) + 1
print("\nSatellite systems found:") if not satellites:
for system, count in sorted(system_counts.items()): print(f"No satellites found for systems: {', '.join(sorted(systems_upper))}")
system_name = { return 1
'G': 'GPS',
'R': 'GLONASS',
'E': 'Galileo',
'C': 'BeiDou'
}.get(system, 'Unknown')
print(f"- {system_name} ({system}): {count} satellites")
# Generate the combined skyplot # Print summary information
# Time window: OBS bounds if provided; else NAV span all_epochs = sorted(list(set(
use_start, use_end = all_epochs[0], all_epochs[-1] e['epoch'] for prn, ephemerides in satellites.items() for e in ephemerides
if args.use_obs: )))
tried = [] print("\nFile contains:")
obs_path = None print(f"- {len(satellites)} unique satellites")
stem = filename[:-1] print(f"- {len(all_epochs)} unique epochs")
print(f"- From {all_epochs[0]} to {all_epochs[-1]}")
for s in ('O', 'o'): # Try uppercase then lowercase # Calculate and print satellite counts by system
candidate = stem + s system_counts = {}
tried.append(candidate) for prn in satellites:
if Path(candidate).exists(): system = prn[0]
obs_path = candidate system_counts[system] = system_counts.get(system, 0) + 1
break
if obs_path: print("\nSatellite systems found:")
obs_start, obs_end = _read_obs_time_bounds(obs_path) for system, count in sorted(system_counts.items()):
if obs_start and obs_end: system_name = {
use_start, use_end = obs_start, obs_end 'G': 'GPS',
print(f"\nObservation window detected from {obs_path}: {use_start} → {use_end}") 'R': 'GLONASS',
'E': 'Galileo',
'C': 'BeiDou'
}.get(system, 'Unknown')
print(f"- {system_name} ({system}): {count} satellites")
# Generate the combined skyplot
# Time window: OBS bounds if provided; else NAV span
use_start, use_end = all_epochs[0], all_epochs[-1]
if args.use_obs:
tried = []
obs_path = None
stem = filename[:-1]
for s in ('O', 'o'): # Try uppercase then lowercase
candidate = stem + s
tried.append(candidate)
if Path(candidate).exists():
obs_path = candidate
break
if obs_path:
obs_start, obs_end = _read_obs_time_bounds(obs_path)
if obs_start and obs_end:
use_start, use_end = obs_start, obs_end
print(f"\nObservation window detected from {obs_path}: {use_start} → {use_end}")
else:
print(f"\nWarning: Could not read valid times from {obs_path}. Using NAV span instead.")
else: else:
print(f"\nWarning: Could not read valid times from {obs_path}. Using NAV span instead.") print(f"\nOBS file not found. Tried: {', '.join(tried)}. Using NAV span instead.")
else:
print(f"\nOBS file not found. Tried: {', '.join(tried)}. Using NAV span instead.")
# Ensure at least two samples with the default 5-minute step # Ensure at least two samples with the default 5-minute step
if (use_end - use_start) < timedelta(minutes=5): if (use_end - use_start) < timedelta(minutes=5):
use_end = use_start + timedelta(minutes=5) use_end = use_start + timedelta(minutes=5)
# Generate the combined skyplot # Generate the combined skyplot
print("\nGenerating skyplot...") print("\nGenerating skyplot...")
footer = f"From {use_start} to {use_end} UTC" footer = f"From {use_start} to {use_end} UTC"
plot_satellite_tracks( plot_satellite_tracks(
satellites, satellites,
obs_lat, obs_lat,
obs_lon, obs_lon,
obs_alt, obs_alt,
footer_text=footer, footer_text=footer,
filename=filename, filename=filename,
show_plot=not args.no_show, show_plot=not args.no_show,
start_time=use_start, start_time=use_start,
end_time=use_end, end_time=use_end,
elev_mask=args.elev_mask elev_mask=args.elev_mask,
) output_format=args.format
)
except Exception as e:
print(f"Error: {str(e)}")
return 1
return 0
if __name__ == "__main__": if __name__ == "__main__":