Skip to content

Commit

Permalink
Fixed bug in `named_arrays._scalars.uncertainties.uncertainties_named…
Browse files Browse the repository at this point in the history
…_array_functions.optimize_root_secant()` where some arguments weren't being broadcasted along the distribution axis.
  • Loading branch information
byrdie committed Aug 31, 2023
1 parent 1ff3d8a commit 9b0513f
Show file tree
Hide file tree
Showing 3 changed files with 222 additions and 0 deletions.
1 change: 1 addition & 0 deletions optika/__init__.py
Original file line number Diff line number Diff line change
Expand Up @@ -2,3 +2,4 @@
from . import rays
from . import sags
from . import apertures
from . import rulings
82 changes: 82 additions & 0 deletions optika/_tests/test_rulings.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,82 @@
import pytest
import numpy as np
import astropy.units as u
import named_arrays as na
import optika
import optika._tests.test_transforms


class AbstractTestAbstractRulings(
optika._tests.test_transforms.AbstractTestTransformable,
):
def test_diffraction_order(self, a: optika.rulings.AbstractRulings):
assert na.get_dtype(a.diffraction_order) == int

def test_ruling_normal(self, a: optika.rulings.AbstractRulings):
position = na.Cartesian3dVectorArray() * u.mm
result = a.ruling_normal(position)
assert isinstance(result, na.AbstractCartesian3dVectorArray)
assert np.all(result.length == 1)

@pytest.mark.parametrize(
argnames="rays",
argvalues=[
optika.rays.RayVectorArray(
wavelength=500 * u.nm,
position=na.Cartesian3dVectorArray(
x=na.linspace(-5, 5, axis="x", num=3) * u.mm,
y=na.linspace(-5, 5, axis="y", num=4) * u.mm,
z=0 * u.mm,
),
direction=na.Cartesian3dVectorArray(0, 0, 1),
)
],
)
def test_rays_apparent(
self,
a: optika.rulings.AbstractRulings,
rays: optika.rays.RayVectorArray,
):
result = a.rays_apparent(rays, index_refraction=1)
assert isinstance(result, optika.rays.AbstractRayVectorArray)
assert np.all(result.wavelength == rays.wavelength)
assert np.all(result.position == rays.position)
assert np.any(result.direction.z != rays.direction.z)
assert np.allclose(result.direction.length, 1)
assert np.all(result.intensity == rays.intensity)
assert np.any(result.index_refraction != rays.index_refraction)


class AbstractTestAbstractConstantDensityRulings(
AbstractTestAbstractRulings,
):
def test_ruling_density(
self,
a: optika.rulings.AbstractConstantDensityRulings,
):
assert a.ruling_density.unit.is_equivalent(1 / u.mm)

def test_ruling_spacing(
self,
a: optika.rulings.AbstractConstantDensityRulings,
):
assert a.ruling_spacing.unit.is_equivalent(u.mm)


@pytest.mark.parametrize(
argnames="a",
argvalues=[
optika.rulings.ConstantDensityRulings(
ruling_density=5000 / u.mm,
diffraction_order=1,
),
optika.rulings.ConstantDensityRulings(
ruling_density=na.linspace(1, 5, axis="rulings", num=4) / u.mm,
diffraction_order=na.ScalarArray(np.array([-1, 0, 1]), axes="m"),
),
],
)
class TestConstantDensityRulings(
AbstractTestAbstractConstantDensityRulings,
):
pass
139 changes: 139 additions & 0 deletions optika/rulings.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,139 @@
import abc
import dataclasses
import astropy.units as u
import named_arrays as na
import optika
import optika.mixins

__all__ = [
"AbstractRulings",
"AbstractConstantDensityRulings",
"ConstantDensityRulings",
]


@dataclasses.dataclass(eq=False, repr=False)
class AbstractRulings(
optika.transforms.Transformable,
):
"""
Interface for the interaction of a ruled surface with incident light
"""

@property
@abc.abstractmethod
def diffraction_order(self) -> na.ScalarLike:
"""
the diffraction order to simulate
"""

@abc.abstractmethod
def ruling_normal(
self,
position: na.AbstractCartesian3dVectorArray,
) -> na.AbstractCartesian3dVectorArray:
"""
vector normal to the plane of the rulings
Parameters
----------
position
location to evaluate the normal vector
"""

@abc.abstractmethod
def rays_apparent(
self,
rays: optika.rays.RayVectorArray,
index_refraction: float,
) -> optika.rays.RayVectorArray:
"""
the apparent input rays from the given actual input rays
Parameters
----------
rays
actual input rays
index_refraction
index of refraction of the output rays
"""


@dataclasses.dataclass(eq=False, repr=False)
class AbstractConstantDensityRulings(
AbstractRulings,
):
@property
@abc.abstractmethod
def ruling_density(self) -> na.ScalarLike:
"""
the frequency of the ruling pattern
"""

@property
def ruling_spacing(self) -> na.ScalarLike:
"""
the distance between successive rulings
"""
return 1 / self.ruling_density

def ruling_normal(
self,
position: na.AbstractCartesian3dVectorArray,
) -> na.AbstractCartesian3dVectorArray:
"""
unit vector normal to the planes of the rulings
Parameters
----------
position
the location to evalulate the ruling normal
Returns
-------
"""
return na.Cartesian3dVectorArray(1, 0, 0)
# ruling_density = self.ruling_density
# return na.Cartesian3dVectorArray(
# x=ruling_density,
# y=0 * ruling_density.unit,
# z=0 * ruling_density.unit,
# )

def rays_apparent(
self,
rays: optika.rays.RayVectorArray,
index_refraction: float,
) -> optika.rays.RayVectorArray:
"""
apparent direction of incoming rays if there was no diffraction
Parameters
----------
rays
incoming rays that will be diffracted by the rulings
index_refraction
the index of refraction after the rulings
"""
a = rays.index_refraction * rays.direction
diffraction_order = self.diffraction_order
ruling_density = self.ruling_density
ruling_normal = ruling_density * self.ruling_normal(rays.position)
a = a + index_refraction * diffraction_order * rays.wavelength * ruling_normal
length_a = a.length
return optika.rays.RayVectorArray(
wavelength=rays.wavelength,
position=rays.position,
direction=a / length_a,
intensity=rays.intensity,
index_refraction=length_a,
)


@dataclasses.dataclass(eq=False, repr=False)
class ConstantDensityRulings(
AbstractConstantDensityRulings,
):
ruling_density: na.ScalarLike = 0 / u.mm
diffraction_order: na.ScalarLike = 1
transform: None | optika.transforms.AbstractTransform = None

0 comments on commit 9b0513f

Please sign in to comment.