# The Clear BSD License
#
# Copyright (c) 2026 Tobias Heibges
# All rights reserved.
#
# Redistribution and use in source and binary forms, with or without
# modification, are permitted (subject to the limitations in the disclaimer
# below) provided that the following conditions are met:
#
# * Redistributions of source code must retain the above copyright notice,
# this list of conditions and the following disclaimer.
#
# * Redistributions in binary form must reproduce the above copyright
# notice, this list of conditions and the following disclaimer in the
# documentation and/or other materials provided with the distribution.
#
# * Neither the name of the copyright holder nor the names of its
# contributors may be used to endorse or promote products derived from this
# software without specific prior written permission.
#
# NO EXPRESS OR IMPLIED LICENSES TO ANY PARTY'S PATENT RIGHTS ARE GRANTED BY
# THIS LICENSE. THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND
# CONTRIBUTORS "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
# LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A
# PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR
# CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
# EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
# PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR
# BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER
# IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
# ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
# POSSIBILITY OF SUCH DAMAGE.
"""
Time Spread Visualization
Plotting functions for geometric time spread estimation and analysis.
"""
import matplotlib.pyplot as plt
import numpy as np
from ..utilities.time_spread import TimeSpreadResult
from .common import GRID_ALPHA, save_figure
# =============================================================================
# Single-Axis Functions
# =============================================================================
[docs]
def plot_time_spread_geometry_side(
ax: plt.Axes,
result: TimeSpreadResult,
source_position: tuple[float, float, float],
detector_position: tuple[float, float, float],
scale: float = 1000.0,
) -> None:
"""
Plot side view (X-Z plane) of time spread geometry.
Parameters
----------
ax : Axes
Matplotlib axes to plot on.
result : TimeSpreadResult
Time spread estimation result.
source_position : tuple
(x, y, z) source position in meters.
detector_position : tuple
(x, y, z) detector position in meters.
scale : float
Scale factor for display (default 1000 for km).
"""
edge_points = result.edge_points
center_point = result.center_point
min_pt = result.min_path_point
max_pt = result.max_path_point
det = np.array(detector_position)
src = np.array(source_position)
# Ocean surface
x_range = max(
abs(edge_points[:, 0].min()), abs(edge_points[:, 0].max()), abs(det[0])
)
x_ocean = np.linspace(min(src[0] / scale - 2, -5), x_range / scale + 5, 200)
ax.fill_between(x_ocean, 0, -2, color="lightblue", alpha=0.5, label="Ocean")
ax.axhline(y=0, color="blue", linewidth=2)
# Source
ax.plot(
src[0] / scale, src[2] / scale, "ro", markersize=12, label="Source", zorder=10
)
# Detector point
ax.plot(
det[0] / scale, det[2] / scale, "g*", markersize=15, label="Detector", zorder=10
)
# Detector sphere arc (faint)
det_r = np.linalg.norm(det) / scale
theta_arc = np.linspace(0, np.pi / 2, 100)
ax.plot(
det_r * np.cos(theta_arc),
det_r * np.sin(theta_arc),
"g-",
alpha=0.3,
linewidth=1,
)
# Center ray path
ax.plot(
[src[0] / scale, center_point[0] / scale],
[src[2] / scale, center_point[2] / scale],
"k-",
linewidth=2,
alpha=0.5,
)
ax.plot(
[center_point[0] / scale, det[0] / scale],
[center_point[2] / scale, det[2] / scale],
"k-",
linewidth=2,
alpha=0.5,
label="Center ray",
)
# Shortest path (min time)
ax.plot(
[src[0] / scale, min_pt[0] / scale],
[src[2] / scale, min_pt[2] / scale],
"b-",
linewidth=2,
)
ax.plot(
[min_pt[0] / scale, det[0] / scale],
[min_pt[2] / scale, det[2] / scale],
"b--",
linewidth=2,
label=f"Shortest path ({result.min_path/1000:.2f} km)",
)
ax.plot(min_pt[0] / scale, min_pt[2] / scale, "bo", markersize=8)
# Longest path (max time)
ax.plot(
[src[0] / scale, max_pt[0] / scale],
[src[2] / scale, max_pt[2] / scale],
"r-",
linewidth=2,
)
ax.plot(
[max_pt[0] / scale, det[0] / scale],
[max_pt[2] / scale, det[2] / scale],
"r--",
linewidth=2,
label=f"Longest path ({result.max_path/1000:.2f} km)",
)
ax.plot(max_pt[0] / scale, max_pt[2] / scale, "rs", markersize=8)
ax.set_xlabel("X (km)", fontsize=12)
ax.set_ylabel("Z (km)", fontsize=12)
ax.set_title("Time Spread Geometry - Side View", fontsize=14, fontweight="bold")
ax.legend(loc="upper left", fontsize=9)
ax.grid(True, alpha=GRID_ALPHA)
ax.set_xlim(src[0] / scale - 2, det[0] / scale + 2)
ax.set_ylim(-1, det[2] / scale + 2)
[docs]
def plot_arrival_time_histogram(
ax: plt.Axes,
result: TimeSpreadResult,
threshold_ns: float = 10.0,
) -> None:
"""
Plot histogram of arrival times for edge rays.
Parameters
----------
ax : Axes
Matplotlib axes to plot on.
result : TimeSpreadResult
Time spread estimation result.
threshold_ns : float
Time threshold to mark on plot.
"""
speed_of_light = 299792458.0
arrival_times_ns = (result.path_lengths / speed_of_light) * 1e9
relative_times_ns = arrival_times_ns - arrival_times_ns.min()
ax.hist(relative_times_ns, bins=30, color="steelblue", edgecolor="black", alpha=0.7)
ax.axvline(x=0, color="blue", linestyle="-", linewidth=2, label="Earliest arrival")
ax.axvline(
x=result.time_spread_ns,
color="red",
linestyle="-",
linewidth=2,
label="Latest arrival",
)
ax.axvline(
x=threshold_ns,
color="orange",
linestyle="--",
linewidth=1.5,
label=f"{threshold_ns:.0f} ns threshold",
)
ax.set_xlabel("Relative Arrival Time (ns)", fontsize=12)
ax.set_ylabel("Number of Edge Rays", fontsize=12)
ax.set_title("Arrival Time Distribution", fontsize=14, fontweight="bold")
ax.legend(fontsize=9)
ax.grid(True, alpha=GRID_ALPHA)
# Stats box
textstr = "\n".join(
[
f"Mean: {np.mean(relative_times_ns):.2f} ns",
f"Std: {np.std(relative_times_ns):.2f} ns",
f"Range: {result.time_spread_ns:.2f} ns",
]
)
props = {"boxstyle": "round", "facecolor": "wheat", "alpha": 0.8}
ax.text(
0.97,
0.97,
textstr,
transform=ax.transAxes,
fontsize=10,
verticalalignment="top",
horizontalalignment="right",
bbox=props,
)
[docs]
def plot_time_spread_vs_divergence(
ax: plt.Axes,
divergences_deg: np.ndarray,
time_spreads_ns: np.ndarray,
threshold_ns: float = 10.0,
source_altitude: float | None = None,
grazing_angle_deg: float | None = None,
detector_altitude: float | None = None,
) -> None:
"""
Plot time spread as function of beam divergence.
Parameters
----------
ax : Axes
Matplotlib axes to plot on.
divergences_deg : ndarray
Array of divergence angles in degrees.
time_spreads_ns : ndarray
Corresponding time spreads in nanoseconds.
threshold_ns : float
Time threshold to mark on plot.
source_altitude, grazing_angle_deg, detector_altitude : float, optional
Parameters for title display.
"""
ax.plot(divergences_deg, time_spreads_ns, "b-", linewidth=2)
ax.axhline(
y=threshold_ns,
color="r",
linestyle="--",
linewidth=1.5,
label=f"{threshold_ns:.0f} ns threshold",
)
# Find and mark crossing point
if time_spreads_ns[0] < threshold_ns < time_spreads_ns[-1]:
cross_idx = np.where(time_spreads_ns > threshold_ns)[0][0]
cross_div = divergences_deg[cross_idx]
ax.axvline(x=cross_div, color="orange", linestyle=":", linewidth=1.5)
ax.annotate(
f"{cross_div:.2f}°",
(cross_div, threshold_ns),
xytext=(cross_div + 0.3, threshold_ns * 1.5),
fontsize=10,
arrowprops={"arrowstyle": "->", "color": "orange"},
)
ax.set_xlabel("Beam Divergence Half-Angle (degrees)", fontsize=12)
ax.set_ylabel("Maximum Time Spread (ns)", fontsize=12)
# Build title
title_parts = ["Time Spread vs Beam Divergence"]
if source_altitude is not None or grazing_angle_deg is not None:
subtitle_parts = []
if source_altitude is not None:
subtitle_parts.append(f"Altitude: {source_altitude/1000:.1f} km")
if grazing_angle_deg is not None:
subtitle_parts.append(f"Grazing: {grazing_angle_deg}°")
if detector_altitude is not None:
subtitle_parts.append(f"Detector: {detector_altitude/1000:.0f} km")
title_parts.append(f"({', '.join(subtitle_parts)})")
ax.set_title("\n".join(title_parts), fontsize=12, fontweight="bold")
ax.legend(fontsize=10)
ax.grid(True, alpha=GRID_ALPHA)
ax.set_xlim(divergences_deg.min(), divergences_deg.max())
ax.set_yscale("log")
# =============================================================================
# Composite Figure Builders
# =============================================================================
[docs]
def create_time_spread_schematic(
result: TimeSpreadResult,
source_position: tuple[float, float, float],
detector_position: tuple[float, float, float],
grazing_angle_deg: float | None = None,
save_path: str | None = None,
) -> plt.Figure:
"""
Create composite figure showing time spread geometry.
Parameters
----------
result : TimeSpreadResult
Time spread estimation result.
source_position : tuple
(x, y, z) source position in meters.
detector_position : tuple
(x, y, z) detector position in meters.
grazing_angle_deg : float, optional
Grazing angle for display.
save_path : str, optional
Path to save figure.
Returns
-------
fig : Figure
Matplotlib figure.
"""
fig, axes = plt.subplots(1, 2, figsize=(16, 7))
plot_time_spread_geometry_side(axes[0], result, source_position, detector_position)
plot_beam_footprint_top(
axes[1], result, source_position, detector_position, grazing_angle_deg
)
plt.tight_layout()
if save_path:
save_figure(fig, save_path)
print(f" Saved: {save_path}")
return fig