# 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.
"""
Custom Ray Source with Per-Ray Properties
Provides a source where each ray's position, direction, intensity, and
wavelength can be individually specified. Ideal for:
- Chromatic dispersion simulations with different wavelengths per ray
- Custom ray distributions for specialized optical setups
- Importing ray data from external sources
- Testing and validation scenarios
Examples
--------
>>> from lsurf.sources import CustomRaySource
>>> import numpy as np
>>>
>>> # Create RGB rays from same position but different directions
>>> positions = np.array([
... [0, 0, 0],
... [0, 0, 0],
... [0, 0, 0],
... ])
>>> directions = np.array([
... [1, 0, 0.01], # Slight upward tilt
... [1, 0, 0], # Horizontal
... [1, 0, -0.01], # Slight downward tilt
... ])
>>> wavelengths = np.array([700e-9, 550e-9, 450e-9]) # RGB
>>> intensities = np.array([1.0, 1.0, 1.0])
>>>
>>> source = CustomRaySource(
... positions=positions,
... directions=directions,
... wavelengths=wavelengths,
... intensities=intensities,
... )
>>> rays = source.generate()
"""
import numpy as np
import numpy.typing as npt
from ..utilities.ray_data import RayBatch, create_ray_batch
[docs]
class CustomRaySource:
"""
Fully customizable ray source with per-ray properties.
Each ray can have its own position, direction, wavelength, and intensity.
This provides maximum flexibility for specialized simulations.
Parameters
----------
positions : array_like, shape (N, 3)
Starting position for each ray in meters.
directions : array_like, shape (N, 3)
Direction vector for each ray (will be normalized).
wavelengths : array_like, shape (N,)
Wavelength for each ray in meters.
intensities : array_like, shape (N,), optional
Intensity/power for each ray in watts. If None, defaults to
uniform distribution with total power of 1.0.
Attributes
----------
num_rays : int
Number of rays.
positions : ndarray, shape (N, 3)
Ray starting positions.
directions : ndarray, shape (N, 3)
Normalized ray directions.
wavelengths : ndarray, shape (N,)
Per-ray wavelengths.
intensities : ndarray, shape (N,)
Per-ray intensities.
power : float
Total power (sum of intensities).
Examples
--------
>>> # Chromatic dispersion study - same start, different wavelengths
>>> n_rays = 100
>>> positions = np.tile([0, 0, 1000], (n_rays, 1)) # All at same point
>>> directions = np.tile([1, 0, 0], (n_rays, 1)) # All same direction
>>> wavelengths = np.linspace(400e-9, 700e-9, n_rays) # Visible spectrum
>>> source = CustomRaySource(positions, directions, wavelengths)
>>> rays = source.generate()
>>> # Fan of rays for refraction analysis
>>> angles = np.linspace(-0.1, 0.1, 50) # +/- 5.7 degrees
>>> positions = np.zeros((50, 3))
>>> directions = np.column_stack([
... np.cos(angles),
... np.zeros(50),
... np.sin(angles)
... ])
>>> wavelengths = np.full(50, 550e-9)
>>> source = CustomRaySource(positions, directions, wavelengths)
"""
[docs]
def __init__(
self,
positions: npt.ArrayLike,
directions: npt.ArrayLike,
wavelengths: npt.ArrayLike,
intensities: npt.ArrayLike | None = None,
):
"""
Initialize custom ray source.
Parameters
----------
positions : array_like, shape (N, 3)
Starting position for each ray in meters.
directions : array_like, shape (N, 3)
Direction vector for each ray (will be normalized).
wavelengths : array_like, shape (N,)
Wavelength for each ray in meters.
intensities : array_like, shape (N,), optional
Intensity for each ray. If None, uniform distribution is used.
Raises
------
ValueError
If array shapes are inconsistent or invalid.
"""
# Convert and validate positions
self._positions = np.asarray(positions, dtype=np.float32)
if self._positions.ndim == 1:
self._positions = self._positions.reshape(1, 3)
if self._positions.ndim != 2 or self._positions.shape[1] != 3:
raise ValueError(
f"positions must have shape (N, 3), got {self._positions.shape}"
)
self._num_rays = self._positions.shape[0]
# Convert and validate directions
self._directions = np.asarray(directions, dtype=np.float32)
if self._directions.ndim == 1:
self._directions = self._directions.reshape(1, 3)
if self._directions.shape != (self._num_rays, 3):
raise ValueError(
f"directions must have shape ({self._num_rays}, 3), "
f"got {self._directions.shape}"
)
# Normalize directions
norms = np.linalg.norm(self._directions, axis=1, keepdims=True)
if np.any(norms < 1e-10):
raise ValueError("direction vectors cannot be zero")
self._directions = self._directions / norms
# Convert and validate wavelengths
self._wavelengths = np.asarray(wavelengths, dtype=np.float32).ravel()
if self._wavelengths.shape[0] != self._num_rays:
raise ValueError(
f"wavelengths must have shape ({self._num_rays},), "
f"got {self._wavelengths.shape}"
)
if np.any(self._wavelengths <= 0):
raise ValueError("wavelengths must be positive")
# Convert and validate intensities
if intensities is None:
# Default: uniform distribution with total power = 1.0
self._intensities = np.full(
self._num_rays, 1.0 / self._num_rays, dtype=np.float32
)
else:
self._intensities = np.asarray(intensities, dtype=np.float32).ravel()
if self._intensities.shape[0] != self._num_rays:
raise ValueError(
f"intensities must have shape ({self._num_rays},), "
f"got {self._intensities.shape}"
)
if np.any(self._intensities < 0):
raise ValueError("intensities cannot be negative")
@property
def num_rays(self) -> int:
"""Number of rays."""
return self._num_rays
@property
def positions(self) -> npt.NDArray[np.float32]:
"""Ray starting positions, shape (N, 3)."""
return self._positions
@property
def directions(self) -> npt.NDArray[np.float32]:
"""Normalized ray directions, shape (N, 3)."""
return self._directions
@property
def wavelengths(self) -> npt.NDArray[np.float32]:
"""Per-ray wavelengths in meters, shape (N,)."""
return self._wavelengths
@property
def intensities(self) -> npt.NDArray[np.float32]:
"""Per-ray intensities, shape (N,)."""
return self._intensities
@property
def power(self) -> float:
"""Total power (sum of intensities)."""
return float(np.sum(self._intensities))
[docs]
def generate(self) -> RayBatch:
"""
Generate ray batch with specified properties.
Creates a RayBatch with positions, directions, wavelengths, and
intensities as specified in the constructor.
Returns
-------
RayBatch
Ray batch ready for propagation.
"""
rays = create_ray_batch(num_rays=self._num_rays)
# Set all per-ray properties
rays.positions[:] = self._positions
rays.directions[:] = self._directions
rays.wavelengths[:] = self._wavelengths
rays.intensities[:] = self._intensities
rays.active[:] = True
return rays
[docs]
def __repr__(self) -> str:
"""Return string representation."""
wl_min = self._wavelengths.min()
wl_max = self._wavelengths.max()
if wl_min == wl_max:
wl_str = f"{wl_min:.2e}m"
else:
wl_str = f"[{wl_min:.2e}, {wl_max:.2e}]m"
return (
f"CustomRaySource(num_rays={self._num_rays}, "
f"wavelengths={wl_str}, power={self.power:.2e}W)"
)
[docs]
@classmethod
def from_spectral_fan(
cls,
origin: tuple[float, float, float],
direction: tuple[float, float, float],
wavelength_range: tuple[float, float],
num_rays: int,
total_power: float = 1.0,
) -> "CustomRaySource":
"""
Create rays with same position/direction but varying wavelengths.
Convenience factory for chromatic dispersion studies.
Parameters
----------
origin : tuple of float
Starting position for all rays (x, y, z) in meters.
direction : tuple of float
Direction for all rays (will be normalized).
wavelength_range : tuple of float
(min_wavelength, max_wavelength) in meters.
num_rays : int
Number of rays to create.
total_power : float, optional
Total power distributed uniformly. Default is 1.0 W.
Returns
-------
CustomRaySource
Source with spectral distribution.
Examples
--------
>>> # Visible spectrum from single point
>>> source = CustomRaySource.from_spectral_fan(
... origin=(0, 0, 0),
... direction=(1, 0, 0),
... wavelength_range=(400e-9, 700e-9),
... num_rays=100,
... )
"""
positions = np.tile(origin, (num_rays, 1)).astype(np.float32)
directions = np.tile(direction, (num_rays, 1)).astype(np.float32)
wavelengths = np.linspace(
wavelength_range[0], wavelength_range[1], num_rays
).astype(np.float32)
intensities = np.full(num_rays, total_power / num_rays, dtype=np.float32)
return cls(positions, directions, wavelengths, intensities)
[docs]
@classmethod
def from_angular_fan(
cls,
origin: tuple[float, float, float],
base_direction: tuple[float, float, float],
angle_range: tuple[float, float],
num_rays: int,
wavelength: float = 550e-9,
total_power: float = 1.0,
fan_axis: str = "vertical",
) -> "CustomRaySource":
"""
Create rays from same position with varying angles.
Convenience factory for refraction/reflection studies.
Parameters
----------
origin : tuple of float
Starting position for all rays (x, y, z) in meters.
base_direction : tuple of float
Central direction of the fan (will be normalized).
angle_range : tuple of float
(min_angle, max_angle) deviation from base direction in radians.
num_rays : int
Number of rays to create.
wavelength : float, optional
Wavelength for all rays. Default is 550 nm.
total_power : float, optional
Total power distributed uniformly. Default is 1.0 W.
fan_axis : str, optional
'vertical' for z-rotation, 'horizontal' for y-rotation.
Default is 'vertical'.
Returns
-------
CustomRaySource
Source with angular distribution.
Examples
--------
>>> # Fan of rays spanning +/- 10 degrees
>>> source = CustomRaySource.from_angular_fan(
... origin=(0, 0, 1000),
... base_direction=(1, 0, 0),
... angle_range=(-0.17, 0.17), # ~10 degrees
... num_rays=50,
... )
"""
positions = np.tile(origin, (num_rays, 1)).astype(np.float32)
# Normalize base direction
base = np.array(base_direction, dtype=np.float32)
base = base / np.linalg.norm(base)
# Create angles
angles = np.linspace(angle_range[0], angle_range[1], num_rays)
# Create directions by rotating base direction
directions = np.zeros((num_rays, 3), dtype=np.float32)
if fan_axis == "vertical":
# Rotate in xz plane (vertical fan)
cos_a = np.cos(angles)
sin_a = np.sin(angles)
directions[:, 0] = base[0] * cos_a - base[2] * sin_a
directions[:, 1] = base[1]
directions[:, 2] = base[0] * sin_a + base[2] * cos_a
else:
# Rotate in xy plane (horizontal fan)
cos_a = np.cos(angles)
sin_a = np.sin(angles)
directions[:, 0] = base[0] * cos_a - base[1] * sin_a
directions[:, 1] = base[0] * sin_a + base[1] * cos_a
directions[:, 2] = base[2]
wavelengths = np.full(num_rays, wavelength, dtype=np.float32)
intensities = np.full(num_rays, total_power / num_rays, dtype=np.float32)
return cls(positions, directions, wavelengths, intensities)