Source code for lsurf.detectors.constant_size_rings

# 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.

"""
Constant-Size Detector Rings with No Shadowing.

This module provides a detector ring geometry where each ring has constant
physical size (radial width) regardless of distance from origin. Adjacent rings
touch exactly (no shadowing) to provide complete angular coverage.

The detector centers lie on a sphere at fixed altitude above Earth's surface,
with normals pointing toward the origin (0,0,0).
"""

from dataclasses import dataclass

import numpy as np
from scipy.optimize import brentq

from ..surfaces import EARTH_RADIUS


[docs] @dataclass class ConstantSizeDetectorRings: """ Constant-size detector ring geometry with no shadowing. Creates a set of annular detector rings centered on a sphere at fixed altitude above Earth's surface. Each ring has constant physical radial width, with angular size varying based on distance from origin. Adjacent rings touch exactly (no overlap, no gaps) when viewed from origin, ensuring complete angular coverage with no shadowing. Parameters ---------- detector_radial_size : float Physical radial width of each detector in meters (default: 10 km) detector_altitude : float Altitude of detector sphere above Earth's surface in meters (default: 33 km) max_elevation_deg : float Maximum elevation angle from horizontal (90° = zenith) (default: 90°) min_elevation_deg : float Minimum elevation angle, where to stop generating rings (default: -2°) earth_radius : float Earth radius in meters (default: EARTH_RADIUS constant) Attributes ---------- ring_boundaries_deg : ndarray Elevation angles of ring boundaries (N+1 values for N rings) ring_centers_deg : ndarray Elevation angles of ring centers (N values) ring_distances : ndarray Distances from origin to ring centers in meters (N values) n_rings : int Number of detector rings detector_sphere_radius : float Radius of detector sphere from Earth center in meters Examples -------- >>> rings = ConstantSizeDetectorRings( ... detector_radial_size=10000.0, # 10 km ... detector_altitude=33000.0, # 33 km ... ) >>> print(f"Created {rings.n_rings} rings") Created 17 rings >>> print(f"Coverage: {rings.ring_boundaries_deg[-1]:.1f}° to {rings.ring_boundaries_deg[0]:.1f}°") Coverage: -2.3° to 90.0° """ detector_radial_size: float = 10000.0 # 10 km detector_altitude: float = 33000.0 # 33 km max_elevation_deg: float = 90.0 # Zenith min_elevation_deg: float = -2.0 # Stop 2 deg below horizontal earth_radius: float = EARTH_RADIUS # Computed attributes (set in __post_init__) ring_boundaries_deg: np.ndarray = None ring_centers_deg: np.ndarray = None ring_distances: np.ndarray = None n_rings: int = 0 detector_sphere_radius: float = 0.0
[docs] def __post_init__(self): """Compute ring geometry after initialization.""" self.detector_sphere_radius = self.earth_radius + self.detector_altitude self._compute_ring_geometry()
@property def detector_half_width(self) -> float: """Physical half-width of each detector in meters.""" return self.detector_radial_size / 2
[docs] def distance_at_elevation(self, elev_deg: float) -> float: """ Compute distance from origin (0,0,0) to detector sphere at given elevation. Uses the formula: d(θ) = -sin(θ)·R_E + √(R_d² - R_E²·cos²(θ)) where R_E = Earth radius, R_d = detector sphere radius, θ = elevation angle. Parameters ---------- elev_deg : float Elevation angle from horizontal in degrees (90° = zenith) Returns ------- float Distance from origin to intersection point (meters) Raises ------ ValueError If no intersection exists at the given elevation """ elev_rad = np.radians(elev_deg) cos_e, sin_e = np.cos(elev_rad), np.sin(elev_rad) discriminant = self.detector_sphere_radius**2 - self.earth_radius**2 * cos_e**2 if discriminant < 0: raise ValueError(f"No intersection at elevation {elev_deg}°") return -sin_e * self.earth_radius + np.sqrt(discriminant)
[docs] def point_at_elevation( self, elev_deg: float, azimuth_deg: float = 0.0 ) -> np.ndarray: """ Compute 3D point on detector sphere at given elevation and azimuth. Parameters ---------- elev_deg : float Elevation angle from horizontal in degrees azimuth_deg : float Azimuth angle in degrees (0 = +x direction) Returns ------- ndarray (x, y, z) position in meters """ dist = self.distance_at_elevation(elev_deg) elev_rad = np.radians(elev_deg) az_rad = np.radians(azimuth_deg) x = dist * np.cos(elev_rad) * np.cos(az_rad) y = dist * np.cos(elev_rad) * np.sin(az_rad) z = dist * np.sin(elev_rad) return np.array([x, y, z])
[docs] def find_detector_center(self, theta_top_deg: float) -> float: """ Find detector center elevation such that top edge is at theta_top. For a detector with physical half-width w and center at elevation θ_c, the angular half-width as seen from origin is: α = arctan(w / d(θ_c)) The top edge elevation is θ_top = θ_c + α. This function solves for θ_c given θ_top using brentq root finding. Parameters ---------- theta_top_deg : float Desired elevation angle of top edge (degrees) Returns ------- float Center elevation angle (degrees) """ hw = self.detector_half_width def residual(theta_c_deg): dist = self.distance_at_elevation(theta_c_deg) alpha_deg = np.degrees(np.arctan(hw / dist)) return theta_top_deg - (theta_c_deg + alpha_deg) # Search bounds: center must be below top edge lower_bound = max(theta_top_deg - 45.0, -10.0) upper_bound = theta_top_deg - 0.001 return brentq(residual, lower_bound, upper_bound)
[docs] def angular_width_at_ring(self, ring_index: int) -> float: """ Compute angular width of a ring as seen from origin. Parameters ---------- ring_index : int Index of the ring (0 = nearest to zenith) Returns ------- float Angular width in degrees """ if ring_index < 0 or ring_index >= self.n_rings: raise ValueError( f"Ring index {ring_index} out of range [0, {self.n_rings})" ) dist = self.ring_distances[ring_index] return 2 * np.degrees(np.arctan(self.detector_half_width / dist))
[docs] def horizontal_distance_at_elevation(self, elev_deg: float) -> float: """ Compute horizontal distance from origin to detector at given elevation. Parameters ---------- elev_deg : float Elevation angle in degrees Returns ------- float Horizontal distance in meters """ dist = self.distance_at_elevation(elev_deg) return dist * np.cos(np.radians(elev_deg))
[docs] def get_ring_horizontal_distances(self) -> np.ndarray: """ Get horizontal distances for all ring boundaries. Returns ------- ndarray Horizontal distances in meters for each ring boundary """ return np.array( [self.horizontal_distance_at_elevation(e) for e in self.ring_boundaries_deg] )
def _compute_ring_geometry(self): """ Build ring boundaries with constant physical size (no shadowing). Algorithm: 1. Start at θ_boundary[0] = max_elevation_deg (zenith) 2. Find center θ_c such that top edge is at θ_boundary[i] 3. Compute bottom edge: θ_bottom = θ_c - α = 2·θ_c - θ_top 4. Set θ_boundary[i+1] = θ_bottom 5. Repeat until θ_boundary < min_elevation_deg """ boundaries = [self.max_elevation_deg] centers = [] distances = [] current_top = self.max_elevation_deg while current_top > self.min_elevation_deg: # Find center such that top edge is at current_top theta_c = self.find_detector_center(current_top) dist_c = self.distance_at_elevation(theta_c) # Angular half-width at this distance alpha_deg = np.degrees(np.arctan(self.detector_half_width / dist_c)) # Bottom edge (no-shadowing: adjacent ring's top = this ring's bottom) theta_bottom = theta_c - alpha_deg centers.append(theta_c) distances.append(dist_c) boundaries.append(theta_bottom) current_top = theta_bottom self.ring_boundaries_deg = np.array(boundaries) self.ring_centers_deg = np.array(centers) self.ring_distances = np.array(distances) self.n_rings = len(centers)
[docs] def summary(self) -> str: """ Return a summary string of the detector ring configuration. Returns ------- str Multi-line summary of configuration """ lines = [ "Constant-Size Detector Ring Configuration:", f" Detector altitude: {self.detector_altitude/1000:.1f} km", f" Sphere radius: {self.detector_sphere_radius/1000:.1f} km", f" Detector radial size: {self.detector_radial_size/1000:.1f} km", f" Elevation range: {self.ring_boundaries_deg[-1]:.2f}° to {self.ring_boundaries_deg[0]:.2f}°", f" Number of rings: {self.n_rings}", ] if self.n_rings > 0: aw0 = self.angular_width_at_ring(0) aw_last = self.angular_width_at_ring(self.n_rings - 1) lines.extend( [ f" Ring 0 (zenith): {self.ring_distances[0]/1000:.1f} km, {aw0:.2f}° angular width", f" Ring {self.n_rings-1} (edge): {self.ring_distances[-1]/1000:.1f} km, {aw_last:.2f}° angular width", ] ) return "\n".join(lines)
[docs] def get_ring_info(self, ring_index: int) -> dict: """ Get detailed information about a specific ring. Parameters ---------- ring_index : int Index of the ring Returns ------- dict Dictionary with ring properties """ if ring_index < 0 or ring_index >= self.n_rings: raise ValueError( f"Ring index {ring_index} out of range [0, {self.n_rings})" ) return { "ring_index": ring_index, "center_elevation_deg": self.ring_centers_deg[ring_index], "top_elevation_deg": self.ring_boundaries_deg[ring_index], "bottom_elevation_deg": self.ring_boundaries_deg[ring_index + 1], "distance_m": self.ring_distances[ring_index], "distance_km": self.ring_distances[ring_index] / 1000, "angular_width_deg": self.angular_width_at_ring(ring_index), "horizontal_distance_m": self.horizontal_distance_at_elevation( self.ring_centers_deg[ring_index] ), }
[docs] def azimuth_bins_for_ring( self, ring_index: int, az_bin_size_m: float, az_range_deg: float = 10.0 ) -> tuple[int, np.ndarray, np.ndarray]: """ Compute azimuthal bins of constant physical size for a ring. At each ring distance, computes how many bins of the given physical width fit within the ±az_range azimuth range. Parameters ---------- ring_index : int Index of the ring az_bin_size_m : float Physical azimuthal bin size in meters az_range_deg : float Azimuth range in degrees (±this value from beam direction) Returns ------- n_bins : int Number of azimuthal bins for this ring az_edges_deg : ndarray Azimuth bin edges in degrees (n_bins + 1 values) az_centers_deg : ndarray Azimuth bin centers in degrees (n_bins values) """ if ring_index < 0 or ring_index >= self.n_rings: raise ValueError( f"Ring index {ring_index} out of range [0, {self.n_rings})" ) dist = self.ring_distances[ring_index] # Arc length for full azimuth range at this distance # Arc = distance * angle_radians (for small angles on a sphere from origin) total_arc_m = dist * np.radians(2 * az_range_deg) # Number of bins that fit n_bins = max(1, int(np.round(total_arc_m / az_bin_size_m))) # Compute bin edges in degrees az_edges_deg = np.linspace(-az_range_deg, az_range_deg, n_bins + 1) az_centers_deg = (az_edges_deg[:-1] + az_edges_deg[1:]) / 2 return n_bins, az_edges_deg, az_centers_deg
[docs] def get_constant_size_grid( self, az_bin_size_m: float, az_range_deg: float = 10.0 ) -> list[dict]: """ Get a grid of constant-size bins across all rings. Each bin has approximately constant physical size: - Radial size: detector_radial_size (same for all rings) - Azimuthal size: az_bin_size_m (variable number of bins per ring) Parameters ---------- az_bin_size_m : float Physical azimuthal bin size in meters az_range_deg : float Azimuth range in degrees (±this value) Returns ------- list of dict List of bin specifications with keys: - ring_idx: int - az_bin_idx: int - n_az_bins: int (total azimuth bins for this ring) - az_lo_deg, az_hi_deg: float - az_center_deg: float - distance_m: float - bin_area_m2: float (approximate) """ grid = [] for ring_idx in range(self.n_rings): n_az_bins, az_edges, az_centers = self.azimuth_bins_for_ring( ring_idx, az_bin_size_m, az_range_deg ) dist = self.ring_distances[ring_idx] # Approximate bin area (radial_size * azimuthal_arc) az_width_rad = np.radians(az_edges[1] - az_edges[0]) az_arc_m = dist * az_width_rad bin_area = self.detector_radial_size * az_arc_m for az_bin_idx in range(n_az_bins): grid.append( { "ring_idx": ring_idx, "az_bin_idx": az_bin_idx, "n_az_bins": n_az_bins, "az_lo_deg": az_edges[az_bin_idx], "az_hi_deg": az_edges[az_bin_idx + 1], "az_center_deg": az_centers[az_bin_idx], "distance_m": dist, "elev_center_deg": self.ring_centers_deg[ring_idx], "bin_area_m2": bin_area, } ) return grid
[docs] def create_default_detector_rings() -> ConstantSizeDetectorRings: """ Create detector rings with default parameters (10 km size, 33 km altitude). Returns ------- ConstantSizeDetectorRings Detector ring configuration with default parameters """ return ConstantSizeDetectorRings()