diff --git a/optika/rays/_ray_vectors.py b/optika/rays/_ray_vectors.py index cc3176b2..49fb0bfc 100644 --- a/optika/rays/_ray_vectors.py +++ b/optika/rays/_ray_vectors.py @@ -70,6 +70,22 @@ def unvignetted(self) -> na.AbstractScalar: and :obj:`False` indicates the ray has been blocked by an aperture. """ + @property + def n(self) -> complex | na.AbstractScalar: + """ + The complex index of refraction of the medium that the ray is + traveling in. + + The real part is :attr:`index_refraction` and the imaginary part is the + extinction coefficient, :math:`k = \\alpha \\lambda / 4 \\pi`, computed + from the :attr:`attenuation` coefficient :math:`\\alpha` and the + :attr:`wavelength` :math:`\\lambda`. + """ + return ( + self.index_refraction + + self.attenuation * self.wavelength / (4 * np.pi) * 1j + ) + @property def type_abstract(self) -> type[na.AbstractArray]: return AbstractRayVectorArray diff --git a/optika/rays/_tests/test_ray_vectors.py b/optika/rays/_tests/test_ray_vectors.py index c5aa626a..432ef85a 100644 --- a/optika/rays/_tests/test_ray_vectors.py +++ b/optika/rays/_tests/test_ray_vectors.py @@ -80,6 +80,21 @@ def test_attenuation(self, array: optika.rays.AbstractRayVectorArray): def test_index_refraction(self, array: optika.rays.AbstractRayVectorArray): assert isinstance(na.as_named_array(array.index_refraction), na.AbstractScalar) + def test_n(self, array: optika.rays.AbstractRayVectorArray): + result = array.n + assert isinstance(na.as_named_array(result), na.AbstractScalar) + + # The real part is the real index of refraction. + assert np.allclose(np.real(result), array.index_refraction) + + # The imaginary part is the (non-negative) extinction coefficient, + # k = alpha * lambda / (4 pi). + extinction = (array.attenuation * array.wavelength / (4 * np.pi)).to( + u.dimensionless_unscaled + ) + assert np.allclose(np.imag(result), extinction) + assert np.all(np.imag(result) >= 0) + def test_unvignetted(self, array: optika.rays.AbstractRayVectorArray): assert isinstance(na.as_named_array(array.unvignetted), na.AbstractScalar) diff --git a/optika/sensors/__init__.py b/optika/sensors/__init__.py index 42f8d154..d5d89bba 100644 --- a/optika/sensors/__init__.py +++ b/optika/sensors/__init__.py @@ -17,6 +17,7 @@ fano_factor_inf, transmittance, absorbance, + absorption_effective, charge_collection_efficiency, quantum_efficiency_effective, probability_measurement, @@ -43,6 +44,7 @@ "fano_factor_inf", "transmittance", "absorbance", + "absorption_effective", "charge_collection_efficiency", "quantum_efficiency_effective", "probability_measurement", diff --git a/optika/sensors/_sensors.py b/optika/sensors/_sensors.py index 58adf079..43586b3a 100644 --- a/optika/sensors/_sensors.py +++ b/optika/sensors/_sensors.py @@ -81,70 +81,54 @@ def aperture(self): half_width=self.width_pixel * self.num_pixel / 2, ) - def readout( + def collect( self, rays: optika.rays.RayVectorArray, wavelength: na.AbstractScalar, - timedelta: None | u.Quantity | na.AbstractScalar = None, axis: None | str | Sequence[str] = None, where: bool | na.AbstractScalar = True, - noise: bool = True, - ) -> na.FunctionArray[ - na.SpectralPositionalVectorArray, + ) -> tuple[ + na.FunctionArray[na.SpectralPositionalVectorArray, na.AbstractScalar], na.AbstractScalar, ]: """ - Given a set of rays incident on the sensor surface, - where each ray represents an expected number of photons per unit time, - simulate the number of electrons that would be measured by the sensor. + Bin a cloud of rays onto the pixel grid. + + Returns the per-pixel photon image (the binned ray intensity) and the + flux-weighted mean cosine of the *refracted* angle inside the + light-sensitive region in each pixel: the two quantities :meth:`expose` + needs. Refracting each ray here (with its own ambient index of + refraction) and binning the result is what lets :meth:`expose` and the + material's :meth:`~optika.sensors.materials.AbstractSensorMaterial.signal` + model be shared with systems that have no rays, + without threading a separate ambient-index argument through them. Parameters ---------- rays - A set of incident rays in local coordinates to measure. + A set of incident rays in local coordinates to bin. wavelength The edges of the wavelength bins to sample. - timedelta - The exposure time of the measurement. - If :obj:`None` (the default), the value in :attr:`timedelta_exposure` - will be used. axis The logical axes along which to collect photons. where - A boolean mask used to indicate which photons should be considered - when calculating the signal measured by the sensor. - noise - Whether to add noise to the result + A boolean mask used to indicate which rays should be considered. """ - if timedelta is None: - timedelta = self.timedelta_exposure - where = where & rays.unvignetted - sgn = np.sign(rays.intensity) - - rays = dataclasses.replace( - rays, - intensity=np.abs(rays.intensity) * timedelta, - ) - normal = self.sag.normal(rays.position) - rays = self.material.signal( - rays=rays, + # Cosine of the refracted angle inside the light-sensitive region, + # folding in each ray's ambient index of refraction. This is generally + # complex, so its real and imaginary parts are binned separately below. + direction = self.material.direction_refracted( + wavelength=rays.wavelength, + direction=rays.direction, + n=rays.n, normal=normal, - noise=noise, - ) - - rays = dataclasses.replace( - rays, - intensity=sgn * rays.intensity, ) - rays = self.material.charge_diffusion( - rays=rays, - normal=normal, - ) + flux = rays.intensity * where bins = na.SpectralPositionalVectorArray( wavelength=wavelength, @@ -155,15 +139,149 @@ def readout( num=self.num_pixel + 1, ), ) + a = na.SpectralPositionalVectorArray( + wavelength=rays.wavelength, + position=rays.position.xy, + ) - return na.histogram( - a=na.SpectralPositionalVectorArray( - wavelength=rays.wavelength, - position=rays.position.xy, - ), - bins=bins, + image = na.histogram(a, bins=bins, axis=axis, weights=flux) + moment_real = na.histogram( + a, bins=bins, axis=axis, weights=flux * np.real(direction) + ) + moment_imag = na.histogram( + a, bins=bins, axis=axis, weights=flux * np.imag(direction) + ) + + # flux-weighted mean refracted cosine; empty pixels carry no flux, so + # assume normal incidence (cosine of 1). + nonempty = image.outputs > 0 + with np.errstate(invalid="ignore", divide="ignore"): + direction_real = np.where(nonempty, moment_real.outputs / image.outputs, 1) + direction_imag = np.where(nonempty, moment_imag.outputs / image.outputs, 0) + direction = direction_real + direction_imag * 1j + + return image, direction + + def expose( + self, + image: na.FunctionArray[ + na.SpectralPositionalVectorArray, + na.AbstractScalar, + ], + direction: float | na.AbstractScalar = 1, + axis_wavelength: None | str = None, + timedelta: None | u.Quantity | na.AbstractScalar = None, + noise: bool = True, + ) -> na.FunctionArray[ + na.SpectralPositionalVectorArray, + na.AbstractScalar, + ]: + """ + Convert a per-pixel photon image into the electrons measured by the sensor. + + This is the detector-physics step shared by every optical system: it + operates on a pixel grid (a photon image plus a refracted-cosine map), + so it can be driven by :meth:`collect` for a ray-traced system, or by + any model that produces a per-pixel photon image directly. + + The photon flux is multiplied by the exposure time and converted to + electrons using + :meth:`~optika.sensors.materials.AbstractSensorMaterial.signal`, which + applies the quantum efficiency, noise, and charge-diffusion models. + + Parameters + ---------- + image + The expected photon flux incident on each pixel, as a function of + wavelength and pixel position. + The wavelength inputs (``image.inputs.wavelength``) must be the + bin *edges*, not the centers. + direction + The cosine of the refracted angle inside the light-sensitive region + in each pixel, as produced by :meth:`collect`. + axis_wavelength + The logical axis of `image` corresponding to changing wavelength. + If :obj:`None` (the default), ``image.inputs.wavelength`` must have + only one logical axis. + timedelta + The exposure time of the measurement. + If :obj:`None` (the default), the value in :attr:`timedelta_exposure` + will be used. + noise + Whether to add shot noise and intrinsic sensor noise to the result. + """ + if axis_wavelength is None: + shape_wavelength = na.shape(image.inputs.wavelength) + if len(shape_wavelength) != 1: # pragma: nocover + raise ValueError( + f"if `axis_wavelength` is `None`, `image.inputs.wavelength` " + f"must have exactly one logical axis, got {shape_wavelength}." + ) + (axis_wavelength,) = shape_wavelength + + if timedelta is None: + timedelta = self.timedelta_exposure + + photons = image.outputs * timedelta + + electrons = self.material.signal( + photons=photons, + wavelength=image.inputs.wavelength.cell_centers(axis_wavelength), + direction=direction, + width_pixel=self.width_pixel, + axis_xy=(self.axis_pixel.x, self.axis_pixel.y), + noise=noise, + ) + + return dataclasses.replace(image, outputs=electrons) + + def measure( + self, + rays: optika.rays.RayVectorArray, + wavelength: na.AbstractScalar, + axis: None | str | Sequence[str] = None, + where: bool | na.AbstractScalar = True, + timedelta: None | u.Quantity | na.AbstractScalar = None, + noise: bool = True, + ) -> na.FunctionArray[ + na.SpectralPositionalVectorArray, + na.AbstractScalar, + ]: + """ + Bin a set of rays onto the pixel grid and convert them to the electrons + measured by the sensor. + + This composes :meth:`collect` (gather rays into the pixel grid) with + :meth:`expose` (apply the detector physics). + + Parameters + ---------- + rays + A set of incident rays in local coordinates to measure. + wavelength + The edges of the wavelength bins to sample. + axis + The logical axes along which to collect photons. + where + A boolean mask used to indicate which rays should be considered. + timedelta + The exposure time of the measurement. + If :obj:`None` (the default), the value in :attr:`timedelta_exposure` + will be used. + noise + Whether to add shot noise and intrinsic sensor noise to the result. + """ + image, direction = self.collect( + rays=rays, + wavelength=wavelength, axis=axis, - weights=rays.intensity * where, + where=where, + ) + return self.expose( + image, + direction=direction, + timedelta=timedelta, + noise=noise, ) diff --git a/optika/sensors/_sensors_test.py b/optika/sensors/_sensors_test.py index 4256e0be..a48e8f88 100644 --- a/optika/sensors/_sensors_test.py +++ b/optika/sensors/_sensors_test.py @@ -26,6 +26,7 @@ def test_timedelta_exposure(self, a: optika.sensors.AbstractImagingSensor): intensity=100 * u.photon / u.s, wavelength=500 * u.nm, position=na.Cartesian3dVectorArray() * u.mm, + direction=na.Cartesian3dVectorArray(0, 0, 1), ), optika.rays.RayVectorArray( intensity=na.random.poisson(100, shape_random=dict(t=11)) * u.erg / u.s, @@ -35,6 +36,7 @@ def test_timedelta_exposure(self, a: optika.sensors.AbstractImagingSensor): y=na.random.uniform(-1, 1, shape_random=dict(t=11)) * u.mm, z=0 * u.mm, ), + direction=na.Cartesian3dVectorArray(0, 0, 1), ), ], ) @@ -44,13 +46,13 @@ def test_timedelta_exposure(self, a: optika.sensors.AbstractImagingSensor): na.linspace(500, 600, axis="wavelength", num=11) * u.nm, ], ) - def test_readout( + def test_measure( self, a: optika.sensors.AbstractImagingSensor, rays: optika.rays.RayVectorArray, wavelength: u.Quantity | na.AbstractScalar, ): - result = a.readout(rays, wavelength) + result = a.measure(rays, wavelength) assert isinstance(result, na.FunctionArray) assert isinstance(result.inputs, na.SpectralPositionalVectorArray) assert isinstance(result.outputs, na.AbstractScalar) diff --git a/optika/sensors/materials/_e2v_ccd203/_e2v_ccd203.py b/optika/sensors/materials/_e2v_ccd203/_e2v_ccd203.py index 8d5511bc..d33193ec 100644 --- a/optika/sensors/materials/_e2v_ccd203/_e2v_ccd203.py +++ b/optika/sensors/materials/_e2v_ccd203/_e2v_ccd203.py @@ -88,10 +88,8 @@ def e2v_ccd203( # Evaluate the fitted QE using the given wavelengths eqe_fit = material.quantum_efficiency_effective( - rays=optika.rays.RayVectorArray( - wavelength=wavelength_fit, - direction=na.Cartesian3dVectorArray(0, 0, 1), - ), + wavelength=wavelength_fit, + direction=na.Cartesian3dVectorArray(0, 0, 1), normal=na.Cartesian3dVectorArray(0, 0, -1), ) @@ -154,11 +152,7 @@ def e2v_ccd203( # Compute the width of the charge diffusion kernel # for each wavelength. width = material.width_charge_diffusion( - rays=optika.rays.RayVectorArray( - wavelength=wavelength_fit, - direction=na.Cartesian3dVectorArray(0, 0, 1), - ), - normal=na.Cartesian3dVectorArray(0, 0, -1), + wavelength=wavelength_fit, ) # Plot the results diff --git a/optika/sensors/materials/_e2v_ccd97/_e2v_ccd97.py b/optika/sensors/materials/_e2v_ccd97/_e2v_ccd97.py index c445dbd8..b6e0cb21 100644 --- a/optika/sensors/materials/_e2v_ccd97/_e2v_ccd97.py +++ b/optika/sensors/materials/_e2v_ccd97/_e2v_ccd97.py @@ -95,10 +95,8 @@ def e2v_ccd97( # Evaluate the fitted QE using the given wavelengths eqe_fit = material.quantum_efficiency_effective( - rays=optika.rays.RayVectorArray( - wavelength=wavelength_fit, - direction=na.Cartesian3dVectorArray(0, 0, 1), - ), + wavelength=wavelength_fit, + direction=na.Cartesian3dVectorArray(0, 0, 1), normal=na.Cartesian3dVectorArray(0, 0, -1), ) @@ -161,10 +159,8 @@ def e2v_ccd97( material_ccd_aia = optika.sensors.materials.e2v_ccd203() eqe_fit_aia = material_ccd_aia.quantum_efficiency_effective( - rays=optika.rays.RayVectorArray( - wavelength=wavelength_fit, - direction=na.Cartesian3dVectorArray(0, 0, 1), - ), + wavelength=wavelength_fit, + direction=na.Cartesian3dVectorArray(0, 0, 1), normal=na.Cartesian3dVectorArray(0, 0, -1), ) @@ -205,11 +201,7 @@ def e2v_ccd97( # Compute the width of the charge diffusion kernel # for each wavelength. width = material.width_charge_diffusion( - rays=optika.rays.RayVectorArray( - wavelength=wavelength_fit, - direction=na.Cartesian3dVectorArray(0, 0, 1), - ), - normal=na.Cartesian3dVectorArray(0, 0, -1), + wavelength=wavelength_fit, ) # Plot the results diff --git a/optika/sensors/materials/_materials.py b/optika/sensors/materials/_materials.py index 38606333..770905b8 100644 --- a/optika/sensors/materials/_materials.py +++ b/optika/sensors/materials/_materials.py @@ -13,6 +13,7 @@ _thickness_oxide, _thickness_implant, _thickness_substrate, + _width_pixel, _cce_backsurface, ) from ._ramanathan_2020 import ( @@ -33,6 +34,8 @@ "quantum_yield_ideal", "fano_factor", "fano_factor_inf", + "absorption_effective", + "transmittance", "absorbance", "charge_collection_efficiency", "quantum_efficiency_effective", @@ -49,10 +52,51 @@ ] +def absorption_effective( + wavelength: u.Quantity | na.AbstractScalar, + n_substrate: complex | na.AbstractScalar, + direction_substrate: complex | na.AbstractScalar, +) -> na.AbstractScalar: + r""" + The absorption coefficient per unit perpendicular depth for a wave + refracted into the (absorbing) light-sensitive region at oblique incidence. + + This is the rigorous generalization of :math:`4 \pi \, \text{Im}(n) / \lambda` + to oblique incidence: it uses the complex refracted cosine + ``direction_substrate`` so that the geometric :math:`\sec\theta` path-length + factor is no longer needed downstream. Reduces to the normal-incidence + coefficient when ``direction_substrate`` is 1. + + Parameters + ---------- + wavelength + The wavelength of the incident light in vacuum. + n_substrate + The complex index of refraction of the light-sensitive region. + direction_substrate + The (generally complex) cosine of the refracted angle inside the + light-sensitive region, e.g. from :func:`optika.materials.snells_law_scalar`. + + Notes + ----- + For a wave refracted from a transparent ambient medium into an absorbing + medium, the planes of constant amplitude remain parallel to the interface, + so the intensity attenuates along the surface normal with coefficient + + .. math:: + + \alpha_\perp = \frac{4 \pi}{\lambda} \, \text{Im}(n_2 \cos\theta_2), + + where :math:`n_2` is the complex index of the light-sensitive region and + :math:`\cos\theta_2` is its (complex) refracted cosine. + """ + return 4 * np.pi * np.imag(n_substrate * direction_substrate) / wavelength + + def transmittance( wavelength: u.Quantity | na.AbstractScalar, direction: float | na.AbstractScalar = 1, - n: float | na.AbstractScalar = 1, + n: complex | na.AbstractScalar = 1, thickness_oxide: u.Quantity | na.AbstractScalar = _thickness_oxide, thickness_substrate: u.Quantity | na.AbstractScalar = _thickness_substrate, chemical_oxide: str | optika.chemicals.AbstractChemical = "SiO2", @@ -69,11 +113,10 @@ def transmittance( wavelength The wavelength of the incident light in vacuum. direction - The component of the incident light's propagation direction antiparallel - to the surface normal of the sensor. + The cosine of the incidence angle. Default is normal incidence. n - The index of refraction in the ambient medium. + The complex index of refraction in the ambient medium. thickness_oxide The thickness of the oxide layer on the illuminated surface of the sensor. Default is the value given in :cite:t:`Stern1994`. @@ -157,13 +200,14 @@ def transmittance( def absorbance( wavelength: u.Quantity | na.AbstractScalar, direction: float | na.AbstractScalar = 1, - n: float | na.AbstractScalar = 1, + n: complex | na.AbstractScalar = 1, thickness_oxide: u.Quantity | na.AbstractScalar = _thickness_oxide, thickness_substrate: u.Quantity | na.AbstractScalar = _thickness_substrate, chemical_oxide: str | optika.chemicals.AbstractChemical = "SiO2", chemical_substrate: str | optika.chemicals.AbstractChemical = "Si", roughness_oxide: u.Quantity | na.AbstractScalar = 0 * u.nm, roughness_substrate: u.Quantity | na.AbstractScalar = 0 * u.nm, + method: Literal["exact", "Beer-Lambert"] = "Beer-Lambert", ) -> optika.vectors.PolarizationVectorArray: """ The fraction of incident energy absorbed by the light-sensitive @@ -174,11 +218,10 @@ def absorbance( wavelength The wavelength of the incident light in vacuum. direction - The component of the incident light's propagation direction antiparallel - to the surface normal of the sensor. + The cosine of the incidence angle. Default is normal incidence. n - The index of refraction in the ambient medium. + The complex index of refraction in the ambient medium. thickness_oxide The thickness of the oxide layer on the illuminated surface of the sensor. Default is the value given in :cite:t:`Stern1994`. @@ -195,6 +238,13 @@ def absorbance( The RMS roughness the oxide layer surface. roughness_substrate The RMS roughness of the substrate surface. + method + The method to use to compute the absorbance. + If ``exact``, this method allows thin-film interference effects + inside the light-sensitive region. + If ``Beer-Lambert``, this method assumes no interference effects. + These methods only differ in the infrared, where the wavelength is + commensurate with the thickness of the light-sensitive region. Examples -------- @@ -218,8 +268,14 @@ def absorbance( ) # Compute the absorbance vs wavelength - absorbance = optika.sensors.absorbance( + absorbance_exact = optika.sensors.absorbance( + wavelength=wavelength, + method="exact", + ) + + absorbance_beer = optika.sensors.absorbance( wavelength=wavelength, + method="Beer-Lambert", ) # Plot the average absorbance vs. wavelength @@ -232,9 +288,15 @@ def absorbance( ); na.plt.plot( wavelength, - absorbance.average, + absorbance_exact.average, + ax=ax, + label="exact absorbance", + ); + na.plt.plot( + wavelength, + absorbance_beer.average, ax=ax, - label="absorbance", + label="Beer-Lambert absorbance", ); ax.set_xscale("log"); ax.set_xlabel(f"wavelength ({wavelength.unit:latex_inline})"); @@ -247,28 +309,67 @@ def absorbance( if not isinstance(chemical_substrate, optika.chemicals.AbstractChemical): chemical_substrate = optika.chemicals.Chemical(chemical_substrate) - result = optika.materials.layer_absorbance( - index=1, - wavelength=wavelength, - direction=direction, - n=n, - layers=[ - optika.materials.Layer( - chemical=chemical_oxide, - thickness=thickness_oxide, - interface=optika.materials.profiles.ErfInterfaceProfile( - width=roughness_oxide, + if method == "exact": + + result = optika.materials.layer_absorbance( + index=1, + wavelength=wavelength, + direction=direction, + n=n, + layers=[ + optika.materials.Layer( + chemical=chemical_oxide, + thickness=thickness_oxide, + interface=optika.materials.profiles.ErfInterfaceProfile( + width=roughness_oxide, + ), ), - ), - optika.materials.Layer( - chemical=chemical_substrate, - thickness=thickness_substrate, - interface=optika.materials.profiles.ErfInterfaceProfile( - width=roughness_substrate, + optika.materials.Layer( + chemical=chemical_substrate, + thickness=thickness_substrate, + interface=optika.materials.profiles.ErfInterfaceProfile( + width=roughness_substrate, + ), ), - ), - ], - ) + ], + ) + + elif method == "Beer-Lambert": + + _transmittance = transmittance( + wavelength=wavelength, + direction=direction, + n=n, + thickness_oxide=thickness_oxide, + thickness_substrate=thickness_substrate, + chemical_oxide=chemical_oxide, + chemical_substrate=chemical_substrate, + roughness_oxide=roughness_oxide, + roughness_substrate=roughness_substrate, + ) + + n_substrate = chemical_substrate.n(wavelength) + + direction_substrate = optika.materials.snells_law_scalar( + cos_incidence=direction, + index_refraction=n, + index_refraction_new=n_substrate, + ) + + absorption = absorption_effective( + wavelength=wavelength, + n_substrate=n_substrate, + direction_substrate=direction_substrate, + ) + + _transmittance_total = _transmittance * np.exp( + -absorption * thickness_substrate + ) + + result = _transmittance - _transmittance_total + + else: # pragma: nocover + raise ValueError(f"Method {method} not implemented") return np.real(result) @@ -277,7 +378,6 @@ def charge_collection_efficiency( absorption: u.Quantity | na.AbstractScalar, thickness_implant: u.Quantity | na.AbstractScalar = _thickness_implant, cce_backsurface: u.Quantity | na.AbstractScalar = _cce_backsurface, - cos_incidence: float | na.AbstractScalar = 1, ) -> na.AbstractScalar: r""" Compute the average charge collection efficiency using the differential @@ -286,8 +386,11 @@ def charge_collection_efficiency( Parameters ---------- absorption - The absorption coefficient of the light-sensitive material for the - wavelength of interest. + The absorption coefficient of the light-sensitive material per unit + perpendicular depth. + For oblique incidence, supply the effective coefficient from + :func:`absorption_effective`, which folds in the refracted angle, so no + separate angle argument is needed. thickness_implant The thickness of the implant layer, the layer where recombination can occur. @@ -296,9 +399,6 @@ def charge_collection_efficiency( The differential charge collection efficiency on the back surface of the sensor. Default is the value given in :cite:t:`Stern1994`. - cos_incidence - The cosine of the angle of the incident light's propagation direction - inside the substrate with the surface normal Examples -------- @@ -381,27 +481,14 @@ def charge_collection_efficiency( with the addition of an :math:`e^{-\alpha W}` term which represents photons absorbed inside the epitaxial layer but outside the implant layer. - Equation :eq:`cce` is only valid for normally-incident light. - We can generalize it to obliquely-incident light by making the substitution - - .. math:: - :label: x-oblique - - x \rightarrow \frac{x}{\cos \theta} - - where :math:`\theta` is the angle between the propagation direction - inside the silicon substrate and the normal vector. - - Substituting :eq:`x-oblique` into Equation :eq:`cce` and solving yields - - .. math:: - :label: eqe-oblique - - \text{CCE}(\lambda, \theta) = - \eta_0 - + \left( \frac{1 - \eta_0}{\alpha W \sec \theta} \right) (1 - e^{-\alpha W \sec \theta}) + In Equation :eq:`cce`, :math:`\alpha` and :math:`z` are measured along the + surface normal. Oblique incidence is therefore handled entirely by the + absorption coefficient: passing the perpendicular-depth coefficient + :math:`\alpha_\perp` from :func:`absorption_effective` (which folds in the + refracted angle, and reduces to :math:`\alpha` at normal incidence) makes + Equation :eq:`cce` valid at any incidence angle without further modification. """ - z0 = absorption * thickness_implant / np.real(cos_incidence) + z0 = absorption * thickness_implant exp_z0 = np.exp(-z0) term_1 = cce_backsurface @@ -412,8 +499,8 @@ def charge_collection_efficiency( def quantum_efficiency_effective( wavelength: u.Quantity | na.AbstractScalar, - direction: None | na.AbstractCartesian3dVectorArray = None, - n: float | na.AbstractScalar = 1, + direction: float | na.AbstractScalar = 1, + n: complex | na.AbstractScalar = 1, thickness_oxide: u.Quantity | na.AbstractScalar = _thickness_oxide, thickness_implant: u.Quantity | na.AbstractScalar = _thickness_implant, thickness_substrate: u.Quantity | na.AbstractScalar = _thickness_substrate, @@ -422,7 +509,6 @@ def quantum_efficiency_effective( chemical_substrate: str | optika.chemicals.AbstractChemical = "Si", roughness_oxide: u.Quantity | na.AbstractScalar = 0 * u.nm, roughness_substrate: u.Quantity | na.AbstractScalar = 0 * u.nm, - normal: None | na.AbstractCartesian3dVectorArray = None, ) -> na.AbstractScalar: r""" Calculate the effective quantum efficiency of a backilluminated detector. @@ -432,8 +518,8 @@ def quantum_efficiency_effective( wavelength The wavelength of the incident light in vacuum. direction - The propagation direction of the incident light in the ambient medium. - If :obj:`None` (default), normal incidence (:math:`\hat{z}`) is assumed. + The cosine of the incidence angle. + Default is normal incidence. n The complex index of refraction of the ambient medium. thickness_oxide @@ -459,9 +545,6 @@ def quantum_efficiency_effective( The RMS roughness the oxide layer surface. roughness_substrate The RMS roughness of the substrate surface. - normal - The vector perpendicular to the surface of the sensor. - If :obj:`None`, then the normal is assumed to be :math:`-\hat{z}` Examples -------- @@ -522,13 +605,9 @@ def quantum_efficiency_effective( # Define an array of wavelengths with which to sample the EQE wavelength = na.geomspace(10, 10000, axis="wavelength", num=1001) * u.AA - # Define the incidence directions + # Define the cosines of the incidence angles angle = na.linspace(0, 30, axis="angle", num=2) * u.deg - direction = na.Cartesian3dVectorArray( - x=np.sin(angle), - y=0, - z=np.cos(angle), - ) + direction = np.cos(angle) eqe = optika.sensors.quantum_efficiency_effective( wavelength=wavelength, @@ -564,8 +643,6 @@ def quantum_efficiency_effective( and :math:`\text{CCE}(\lambda)` is the charge collection efficiency (computed by :func:`charge_collection_efficiency`). """ - if direction is None: - direction = na.Cartesian3dVectorArray(0, 0, 1) if not isinstance(chemical_oxide, optika.chemicals.AbstractChemical): chemical_oxide = optika.chemicals.Chemical(chemical_oxide) @@ -573,11 +650,6 @@ def quantum_efficiency_effective( if not isinstance(chemical_substrate, optika.chemicals.AbstractChemical): chemical_substrate = optika.chemicals.Chemical(chemical_substrate) - if normal is None: - normal = na.Cartesian3dVectorArray(0, 0, -1) - - direction = -direction @ normal - absorbance_substrate = absorbance( wavelength=wavelength, direction=direction, @@ -591,8 +663,6 @@ def quantum_efficiency_effective( ) n_substrate = chemical_substrate.n(wavelength) - wavenumber_substrate = np.imag(n_substrate) - absorption_substrate = 4 * np.pi * wavenumber_substrate / wavelength direction_substrate = optika.materials.snells_law_scalar( cos_incidence=direction, @@ -600,11 +670,16 @@ def quantum_efficiency_effective( index_refraction_new=n_substrate, ) + absorption_substrate = absorption_effective( + wavelength=wavelength, + n_substrate=n_substrate, + direction_substrate=direction_substrate, + ) + cce = charge_collection_efficiency( absorption=absorption_substrate, thickness_implant=thickness_implant, cce_backsurface=cce_backsurface, - cos_incidence=direction_substrate, ) result = absorbance_substrate.average * cce @@ -750,6 +825,7 @@ def electrons_measured_approx( wavelength: u.Quantity | na.ScalarArray, absorption: None | u.Quantity | na.AbstractScalar = None, thickness_implant: u.Quantity | na.AbstractScalar = _thickness_implant, + thickness_substrate: u.Quantity | na.AbstractScalar = _thickness_substrate, cce_backsurface: u.Quantity | na.AbstractScalar = _cce_backsurface, temperature: u.Quantity | na.ScalarArray = 300 * u.K, iqy: None | u.Quantity | na.AbstractScalar = None, @@ -758,8 +834,8 @@ def electrons_measured_approx( ) -> na.AbstractScalar: r""" A random sample from an approximate distribution of measured electrons - given the number of photons absorbed by the light-sensitive layer of the - sensor. + given the number of photons absorbed by the light-sensitive layer + of the sensor. This function accounts for both Fano noise and recombination noise due to partial-charge collection. @@ -767,21 +843,27 @@ def electrons_measured_approx( Parameters ---------- photons_absorbed - The number of photons absorbed by the light-sensitive layer of the sensor. + The number of photons absorbed by the light-sensitive layer + of the sensor. wavelength The vacuum wavelength of the absorbed photons. absorption - The absorption coefficient of the light-sensitive material for the - wavelength of interest. + The absorption coefficient of silicon per unit perpendicular depth. + For oblique incidence, supply the effective coefficient from + :func:`optika.sensors.absorption_effective`, which folds in the + refracted angle, so no separate angle argument is needed. thickness_implant - The thickness of the implant layer. - Default is the value given in :cite:t:`Stern1994`. + The thickness of the implant layer, where partial-charge collection occurs. + thickness_substrate + The thickness of the entire light-sensitive region of the device. + The absorbed photons are distributed within this region, which sets + the fraction that land in the implant layer. cce_backsurface The differential charge collection efficiency on the back surface of the sensor. - Default is the value given in :cite:t:`Stern1994`. temperature - The temperature of the light-sensitive silicon layer. + The temperature of the silicon detector. + Default is room temperature. iqy The ideal quantum yield of the sensor in electrons per photon. If :obj:`None` (the default), the result of :func:`quantum_yield_ideal` @@ -798,7 +880,7 @@ def electrons_measured_approx( Examples -------- - Plot the energy spectrum of 100 6 keV photons emitted from an Fe-55 + Plot the energy spectrum of twenty 6 keV photons emitted from an Fe-55 radioactive source and compare it to the exact spectrum .. jupyter-execute:: @@ -814,7 +896,7 @@ def electrons_measured_approx( # Define the expected number of photons # for each experiment - photons_absorbed = (100 * u.photon).astype(int) + photons_absorbed = (20 * u.photon).astype(int) # Define the wavelength at which to sample the distribution wavelength = 5.9 * u.keV @@ -890,9 +972,9 @@ def electrons_measured_approx( shape = na.shape_broadcasted( photons_absorbed, absorption, - absorbance, iqy, thickness_implant, + thickness_substrate, cce_backsurface, fano_factor, ) @@ -902,11 +984,18 @@ def electrons_measured_approx( a = absorption W = thickness_implant + D = thickness_substrate n0 = cce_backsurface aW = (a * W).to(u.dimensionless_unscaled).value - - fraction_complete = np.exp(-aW) - fraction_partial = 1 - fraction_complete + aD = (a * D).to(u.dimensionless_unscaled).value + + # The absorbed photons are distributed within the substrate `[0, D)`, so the + # depth follows an exponential truncated to that region (matching the exact + # `electrons_measured` kernel). Conditioned on absorption within the + # substrate, a photon lands in the implant region (partial charge collection) + # with probability (1 - exp(-alpha * W)) / (1 - exp(-alpha * D)); otherwise + # it is absorbed deeper in the substrate (complete charge collection). + fraction_partial = (1 - np.exp(-aW)) / (1 - np.exp(-aD)) photons_absorbed_partial = na.random.binomial( n=photons_absorbed, p=fraction_partial, @@ -945,22 +1034,32 @@ def electrons_measured_approx( return result +_transmittance = transmittance _absorbance = absorbance def signal( photons_expected: u.Quantity | na.AbstractScalar, wavelength: u.Quantity | na.ScalarArray, + direction: float | na.AbstractScalar = 1, + n: complex | na.AbstractScalar = 1, + n_substrate: None | complex | na.AbstractScalar = None, absorbance: None | float | na.AbstractScalar = None, - absorption: None | u.Quantity | na.AbstractScalar = None, thickness_implant: u.Quantity | na.AbstractScalar = _thickness_implant, + thickness_depletion: u.Quantity | na.AbstractScalar = _thickness_substrate, + thickness_substrate: None | na.AbstractScalar = _thickness_substrate, + width_pixel: ( + u.Quantity | na.AbstractScalar | na.AbstractCartesian2dVectorArray + ) = _width_pixel, cce_backsurface: u.Quantity | na.AbstractScalar = _cce_backsurface, temperature: u.Quantity | na.ScalarArray = 300 * u.K, - method: Literal["exact", "approx", "expected"] = "exact", + method: Literal["monte-carlo", "expected"] = "monte-carlo", + axis_xy: None | tuple[str, str] = None, + wrap: bool = False, shape_random: None | dict[str, int] = None, ) -> na.AbstractScalar: r""" - A random sample from the approximate distribution of measured electrons + A random sample from the distribution of measured electrons given the expected number of photons incident on the front surface of the sensor. @@ -974,19 +1073,32 @@ def signal( The `expected` number of photons incident on the detector surface. wavelength The vacuum wavelength of the absorbed photons. + direction + The cosine of the incidence angle. + n + The complex index of refraction of the ambient medium. + n_substrate + The complex index of refraction of the light-sensitive material + If :obj:`None` (the default), the result of + :meth:`optika.chemicals.Chemical.n` for silicon will be used. absorbance - The fraction of incident energy absorbed by the light-sensitive layer - of the detector computed using the average of :func:`absorbance`. + The fraction of incident energy absorbed by the light-sensitive region + of the detector. If :obj:`None` (the default), the result of :func:`absorbance` called with default values will be used. - absorption - The absorption coefficient of the light-sensitive material for the - wavelength of interest. - If :obj:`None` (the default), the result of - :meth:`optika.chemicals.Chemical.absorption` for silicon will be used. thickness_implant The thickness of the implant layer. Default is the value given in :cite:t:`Stern1994`. + thickness_depletion + The thickness of the depletion region, the region with significant electric + field. + If :obj:`None` (the default), this is set to the same value as + `thickness_substrate`. + thickness_substrate + The thickness of the entire light-sensitive region of the device. + The default + width_pixel + The size of a single pixel on the sensor. cce_backsurface The differential charge collection efficiency on the back surface of the sensor. @@ -994,13 +1106,23 @@ def signal( temperature The temperature of the light-sensitive silicon layer. method - The method used to generate random samples of measured electrons. - The `exact` method simulates every photon so it is accurate for low - photon counts but slow for high photon counts. - The `approx` method is much faster, but is only accurate if the number - of photons is high. - The `expected` method does not add any noise to the signal and just - returns the expected number of electrons. + The method used to generate samples of measured electrons. + The `monte-carlo` method draws a random sample by simulating every + photon using :func:`electrons_measured`, including shot, Fano, and + recombination noise as well as charge diffusion. + The `expected` method adds no noise and just returns the expected + number of electrons in each pixel; since it is a per-pixel expectation, + it does not apply charge diffusion. + axis_xy + The two logical axes corresponding to the pixel grid of the sensor + along which electrons will diffuse. + If :obj:`None` (the default), there is no charge diffusion. + wrap + Controls how diffused charge is treated at the edges of the pixel grid. + If :obj:`False` (the default), charge that diffuses past the edge of the + grid is lost, as it would be at the physical edge of a sensor. + If :obj:`True`, the grid is treated as periodic and the charge re-enters + the opposite edge (a toroidal boundary). shape_random Additional shape used to specify the number of samples to draw. @@ -1049,16 +1171,34 @@ def signal( """ if absorbance is None: - absorbance = _absorbance(wavelength).average + absorbance = _absorbance( + wavelength=wavelength, + direction=direction, + n=n, + thickness_substrate=thickness_substrate, + ).average - if absorption is None: - absorption = optika.chemicals.Chemical("Si").absorption(wavelength) + if n_substrate is None: + n_substrate = optika.chemicals.Chemical("Si").n(wavelength) + + direction_substrate = optika.materials.snells_law_scalar( + cos_incidence=direction, + index_refraction=n, + index_refraction_new=n_substrate, + ) + + absorption = absorption_effective( + wavelength=wavelength, + n_substrate=n_substrate, + direction_substrate=direction_substrate, + ) if method == "expected": iqy = quantum_yield_ideal( wavelength=wavelength, temperature=temperature, ) + cce = charge_collection_efficiency( absorption=absorption, thickness_implant=thickness_implant, @@ -1066,33 +1206,37 @@ def signal( ) return iqy * absorbance * cce * photons_expected.to(u.ph) - photons_absorbed_expected = absorbance * photons_expected.to(u.ph) - photons_absorbed = na.random.poisson( - lam=photons_absorbed_expected, - shape_random=shape_random, - ).astype(int) + elif method == "monte-carlo": - kwargs = dict( - photons_absorbed=photons_absorbed, - wavelength=wavelength, - absorption=absorption, - thickness_implant=thickness_implant, - cce_backsurface=cce_backsurface, - temperature=temperature, - shape_random=shape_random, - ) + photons = na.random.poisson( + lam=absorbance * photons_expected.to(u.ph), + shape_random=shape_random, + ).astype(int) + + return electrons_measured( + photons_absorbed=photons, + wavelength=wavelength, + absorption=absorption, + thickness_implant=thickness_implant, + thickness_depletion=thickness_depletion, + thickness_substrate=thickness_substrate, + width_pixel=width_pixel, + cce_backsurface=cce_backsurface, + temperature=temperature, + axis_xy=axis_xy, + wrap=wrap, + shape_random=shape_random, + ) - if method == "exact": - return electrons_measured(**kwargs) - elif method == "approx": - return electrons_measured_approx(**kwargs) else: # pragma: nocover raise ValueError(f"Unrecognized method: {method}") def vmr_signal( wavelength: u.Quantity | na.ScalarArray, - absorption: None | u.Quantity | na.AbstractScalar = None, + direction: float | na.AbstractScalar = 1, + n: complex | na.AbstractScalar = 1, + n_substrate: None | complex | na.AbstractScalar = None, thickness_implant: u.Quantity | na.AbstractScalar = _thickness_implant, cce_backsurface: u.Quantity | na.AbstractScalar = _cce_backsurface, temperature: u.Quantity | na.ScalarArray = 300 * u.K, @@ -1108,11 +1252,14 @@ def vmr_signal( ---------- wavelength The vacuum wavelength of the absorbed photons. - absorption - The absorption coefficient of the light-sensitive material for the - wavelength of interest. + direction + The cosine of the incidence angle. + n + The complex index of refraction of the ambient medium. + n_substrate + The complex index of refraction of the light-sensitive material. If :obj:`None` (the default), the result of - :meth:`optika.chemicals.Chemical.absorption` for silicon will be used. + :meth:`optika.chemicals.Chemical.n` for silicon will be used. thickness_implant The thickness of the implant layer. Default is the value given in :cite:t:`Stern1994`. @@ -1205,10 +1352,22 @@ def vmr_signal( and :math:`\eta_0` is the CCE at the back surface. """ - if absorption is None: - absorption = optika.chemicals.Chemical("Si").absorption(wavelength) + if n_substrate is None: + n_substrate = optika.chemicals.Chemical("Si").n(wavelength) - n = quantum_yield_ideal( + direction_substrate = optika.materials.snells_law_scalar( + cos_incidence=direction, + index_refraction=n, + index_refraction_new=n_substrate, + ) + + absorption = absorption_effective( + wavelength=wavelength, + n_substrate=n_substrate, + direction_substrate=direction_substrate, + ) + + iqy = quantum_yield_ideal( wavelength=wavelength, temperature=temperature, ) @@ -1224,7 +1383,7 @@ def vmr_signal( result = 0 if shot: - F_shot = n * cce + F_shot = iqy * cce result = result + F_shot @@ -1235,12 +1394,13 @@ def vmr_signal( if pcc: n0 = cce_backsurface - aW = (absorption * thickness_implant).to(u.dimensionless_unscaled).value + aW = absorption * thickness_implant + aW = aW.to(u.dimensionless_unscaled).value F_cce = 2 * np.exp(-aW) * np.square((n0 - 1) / aW) * (np.sinh(aW) - aW) / cce unit = u.electron / u.photon - F_pcc = 1 * unit - cce * unit - F_cce * unit + n * F_cce + F * F_cce + F_pcc = 1 * unit - cce * unit - F_cce * unit + iqy * F_cce + F * F_cce result = result + F_pcc @@ -1258,65 +1418,97 @@ class AbstractSensorMaterial( @abc.abstractmethod def signal( self, - rays: optika.rays.RayVectorArray, - normal: na.AbstractCartesian3dVectorArray, + photons: u.Quantity | na.AbstractScalar, + wavelength: u.Quantity | na.AbstractScalar, + direction: float | na.AbstractScalar = 1, + width_pixel: ( + u.Quantity | na.AbstractScalar | na.AbstractCartesian2dVectorArray + ) = 0 + * u.um, + axis_xy: None | tuple[str, str] = None, noise: bool = True, - ) -> optika.rays.RayVectorArray: + wrap: bool = False, + ) -> na.AbstractScalar: """ - Given a set of absorbed rays, compute the number of electrons - measured by the sensor using :func:`signal`. + Given the photons incident on each pixel, compute the number of + electrons measured by the sensor using :func:`signal`. Parameters ---------- - rays - The rays absorbed by the light-sensitive silicon layer. - The :attr:`optika.rays.RayVectorArray.intensity` field should - either be in units of photons or energy. - normal - The vector perpendicular to the surface of the sensor. + photons + The number of photons incident on each pixel. + wavelength + An assumed grid of wavelengths for the incident photons. + direction + The cosine of the refracted angle inside the light-sensitive region, + as produced by :meth:`direction_refracted`. + width_pixel + The physical size of each pixel, used by the charge-diffusion model. + axis_xy + The two logical axes corresponding to the pixel grid of the sensor. + If provided, charge diffusion will occur along these two axes. + If :obj:`None` (the default), no diffusion is performed. noise Whether to add noise to the result. + wrap + Controls how diffused charge is treated at the edges of the pixel grid. + If :obj:`False` (the default), charge that diffuses past the edge of + the grid is lost, as it would be at the physical edge of a sensor. + If :obj:`True`, the grid is treated as periodic and the charge + re-enters the opposite edge (a toroidal boundary). """ @abc.abstractmethod - def photons_incident( + def direction_refracted( self, - electrons: u.Quantity | na.AbstractScalar, wavelength: u.Quantity | na.AbstractScalar, - direction: na.AbstractCartesian3dVectorArray, - normal: na.AbstractCartesian3dVectorArray, - ) -> na.AbstractScalar: + direction: None | na.AbstractCartesian3dVectorArray = None, + n: complex | na.AbstractScalar = 1, + normal: None | na.AbstractCartesian3dVectorArray = None, + ) -> complex | na.AbstractScalar: """ - Given the number of electrons measured by the sensor, - and a grid of wavelengths, compute the expected number of - photons incident on the sensor. + The cosine of the refracted propagation angle inside the light-sensitive + region of the sensor. + + This is the quantity :meth:`signal` expects as its ``direction`` + argument. Performing the refraction here lets the sensor's ``collect`` + method fold the ambient index of refraction into the per-pixel cosine, + so :meth:`signal` does not need a separate ambient-index argument. Parameters ---------- - electrons - The number of electrons measured by the sensor. wavelength - An assumed grid of wavelengths for the incident photons. + The wavelength of the incident light in vacuum. direction - An assumed propagation direction for the incident photons. + The propagation direction of the incident light in the ambient medium. + If :obj:`None` (default), normal incidence (:math:`\\hat{z}`) is assumed. + n + The complex index of refraction of the ambient medium. normal The vector perpendicular to the surface of the sensor. """ @abc.abstractmethod - def charge_diffusion( + def photons_incident( self, - rays: optika.rays.RayVectorArray, + electrons: u.Quantity | na.AbstractScalar, + wavelength: u.Quantity | na.AbstractScalar, + direction: na.AbstractCartesian3dVectorArray, normal: na.AbstractCartesian3dVectorArray, - ) -> optika.rays.RayVectorArray: + ) -> na.AbstractScalar: """ - Given a set of incident rays, compute how much the position of each ray - is perturbed due to charge diffusion within the sensor. + Given the number of electrons measured by the sensor, + and a grid of wavelengths, compute the expected number of + photons incident on the sensor. Parameters ---------- - rays - The rays incident on the sensor surface. + electrons + The number of electrons measured by each pixel. + wavelength + An assumed grid of wavelengths for the incident photons. + direction + An assumed propagation direction for the incident photons. normal The vector perpendicular to the surface of the sensor. """ @@ -1334,25 +1526,45 @@ class IdealSensorMaterial( def signal( self, - rays: optika.rays.RayVectorArray, - normal: na.AbstractCartesian3dVectorArray, - noise: bool = False, - ) -> optika.rays.RayVectorArray: - intensity = rays.intensity - if not intensity.unit.is_equivalent(u.photon): + photons: u.Quantity | na.AbstractScalar, + wavelength: u.Quantity | na.AbstractScalar, + direction: float | na.AbstractScalar = 1, + width_pixel: ( + u.Quantity | na.AbstractScalar | na.AbstractCartesian2dVectorArray + ) = 0 + * u.um, + axis_xy: None | tuple[str, str] = None, + noise: bool = True, + wrap: bool = False, + ) -> na.AbstractScalar: + + if not photons.unit.is_equivalent(u.photon): h = astropy.constants.h c = astropy.constants.c - intensity = intensity / (h * c / rays.wavelength) * u.photon + photons = photons / (h * c / wavelength) * u.photon if noise: - intensity = na.random.poisson(intensity.to(u.ph)).astype(int) + photons = na.random.poisson(photons.to(u.ph)).astype(int) - electrons = intensity * u.electron / u.photon + electrons = photons * u.electron / u.photon electrons = electrons.to(u.electron) - result = dataclasses.replace(rays, intensity=electrons) + return electrons - return result + def direction_refracted( + self, + wavelength: u.Quantity | na.AbstractScalar, + direction: None | na.AbstractCartesian3dVectorArray = None, + n: complex | na.AbstractScalar = 1, + normal: None | na.AbstractCartesian3dVectorArray = None, + ) -> complex | na.AbstractScalar: + # An ideal material does not refract, so the refracted cosine is just + # the cosine of the incidence angle. + if direction is None: + direction = na.Cartesian3dVectorArray(0, 0, 1) + if normal is None: + normal = na.Cartesian3dVectorArray(0, 0, -1) + return -direction @ normal def photons_incident( self, @@ -1363,13 +1575,6 @@ def photons_incident( ) -> na.AbstractScalar: return electrons * u.photon / u.electron - def charge_diffusion( - self, - rays: optika.rays.RayVectorArray, - normal: na.AbstractCartesian3dVectorArray, - ) -> optika.rays.RayVectorArray: - return rays - @dataclasses.dataclass(eq=False, repr=False) class AbstractSiliconSensorMaterial( @@ -1499,8 +1704,7 @@ def depletion(self) -> AbstractDepletionModel: def width_charge_diffusion( self, - rays: optika.rays.RayVectorArray, - normal: na.AbstractCartesian3dVectorArray, + wavelength: u.Quantity | na.AbstractScalar, ) -> na.AbstractScalar: """ The standard deviation of the charge diffusion kernel for this sensor. @@ -1508,37 +1712,89 @@ def width_charge_diffusion( Parameters ---------- - rays - The rays incident on the sensor surface. - normal - The vector perpendicular to the surface of the sensor. + wavelength + The wavelength of the incident light in vacuum. """ return optika.sensors.charge_diffusion( - self._chemical.absorption(rays.wavelength), + self._chemical.absorption(wavelength), thickness_substrate=self.thickness_substrate, thickness_depletion=self.depletion.thickness, ) - def absorbance( + def transmittance( self, - rays: optika.rays.AbstractRayVectorArray, - normal: na.AbstractCartesian3dVectorArray, + wavelength: u.Quantity | na.AbstractScalar, + direction: None | na.AbstractCartesian3dVectorArray = None, + n: complex | na.AbstractScalar = 1, + normal: None | na.AbstractCartesian3dVectorArray = None, ) -> optika.vectors.PolarizationVectorArray: + r""" + Compute the fraction of energy transmitted to the light-sensitive region + of the sensor. + + Parameters + ---------- + wavelength + The wavelength of the incident light in vacuum. + direction + The propagation direction of the incident light in the ambient medium. + If :obj:`None` (default), normal incidence (:math:`\hat{z}`) is assumed. + n + The complex index of refraction of the ambient medium. + normal + The vector perpendicular to the surface of the CCD sensor. """ + if direction is None: + direction = na.Cartesian3dVectorArray(0, 0, 1) + + if normal is None: + normal = na.Cartesian3dVectorArray(0, 0, -1) + + return transmittance( + wavelength=wavelength, + direction=-direction @ normal, + n=n, + thickness_oxide=self.thickness_oxide, + thickness_substrate=self.thickness_substrate, + chemical_oxide=self._chemical_oxide, + chemical_substrate=self._chemical, + roughness_oxide=self.roughness_oxide, + roughness_substrate=self.roughness_substrate, + ) + + def absorbance( + self, + wavelength: u.Quantity | na.AbstractScalar, + direction: None | na.AbstractCartesian3dVectorArray = None, + n: complex | na.AbstractScalar = 1, + normal: None | na.AbstractCartesian3dVectorArray = None, + ) -> optika.vectors.PolarizationVectorArray: + r""" Compute the fraction of energy absorbed by the light-sensitive region of the sensor. Parameters ---------- - rays - The light rays incident on the CCD surface. + wavelength + The wavelength of the incident light in vacuum. + direction + The propagation direction of the incident light in the ambient medium. + If :obj:`None` (default), normal incidence (:math:`\hat{z}`) is assumed. + n + The complex index of refraction of the ambient medium. normal The vector perpendicular to the surface of the CCD sensor. """ + if direction is None: + direction = na.Cartesian3dVectorArray(0, 0, 1) + + if normal is None: + normal = na.Cartesian3dVectorArray(0, 0, -1) + return absorbance( - wavelength=rays.wavelength, - direction=-rays.direction @ normal, - n=rays.index_refraction, + wavelength=wavelength, + direction=-direction @ normal, + n=n, thickness_oxide=self.thickness_oxide, thickness_substrate=self.thickness_substrate, chemical_oxide=self._chemical_oxide, @@ -1549,49 +1805,88 @@ def absorbance( def charge_collection_efficiency( self, - rays: optika.rays.AbstractRayVectorArray, - normal: na.AbstractCartesian3dVectorArray, + wavelength: u.Quantity | na.AbstractScalar, + direction: None | na.AbstractCartesian3dVectorArray = None, + n: complex | na.AbstractScalar = 1, + normal: None | na.AbstractCartesian3dVectorArray = None, ) -> na.AbstractScalar: - """ + r""" Compute the charge collection efficiency of this CCD sensor material using :func:`charge_collection_efficiency`. Parameters ---------- - rays - The light rays incident on the CCD surface. + wavelength + The wavelength of the incident light in vacuum. + direction + The propagation direction of the incident light in the ambient medium. + If :obj:`None` (default), normal incidence (:math:`\hat{z}`) is assumed. + n + The complex index of refraction of the ambient medium. normal The vector perpendicular to the surface of the CCD sensor. """ + + if direction is None: + direction = na.Cartesian3dVectorArray(0, 0, 1) + + if normal is None: + normal = na.Cartesian3dVectorArray(0, 0, -1) + + direction = -direction @ normal + + n_substrate = self._chemical.n(wavelength) + + direction_substrate = optika.materials.snells_law_scalar( + cos_incidence=direction, + index_refraction=n, + index_refraction_new=n_substrate, + ) + return charge_collection_efficiency( - absorption=self._chemical.absorption(rays.wavelength), + absorption=absorption_effective( + wavelength=wavelength, + n_substrate=n_substrate, + direction_substrate=direction_substrate, + ), thickness_implant=self.thickness_implant, cce_backsurface=self.cce_backsurface, - cos_incidence=-rays.direction @ normal, ) def quantum_efficiency_effective( self, - rays: optika.rays.AbstractRayVectorArray, - normal: na.AbstractCartesian3dVectorArray, + wavelength: u.Quantity | na.AbstractScalar, + direction: None | na.AbstractCartesian3dVectorArray = None, + n: complex | na.AbstractScalar = 1, + normal: None | na.AbstractCartesian3dVectorArray = None, ) -> na.AbstractScalar: - """ + r""" Compute the effective quantum efficiency of this CCD material using :func:`quantum_efficiency_effective`. Parameters ---------- - rays - The light rays incident on the CCD surface + wavelength + The wavelength of the incident light in vacuum. + direction + The propagation direction of the incident light in the ambient medium. + If :obj:`None` (default), normal incidence (:math:`\hat{z}`) is assumed. + n + The complex index of refraction of the ambient medium. normal The vector perpendicular to the surface of the CCD. """ - k = rays.attenuation * rays.wavelength / (4 * np.pi) - n = rays.index_refraction + k * 1j + if direction is None: + direction = na.Cartesian3dVectorArray(0, 0, 1) + + if normal is None: + normal = na.Cartesian3dVectorArray(0, 0, -1) + + direction = -direction @ normal return quantum_efficiency_effective( - wavelength=rays.wavelength, - direction=rays.direction, + wavelength=wavelength, + direction=direction, n=n, thickness_oxide=self.thickness_oxide, thickness_implant=self.thickness_implant, @@ -1601,117 +1896,188 @@ def quantum_efficiency_effective( chemical_substrate=self._chemical, roughness_oxide=self.roughness_oxide, roughness_substrate=self.roughness_substrate, - normal=normal, ) def quantum_efficiency( self, - rays: optika.rays.AbstractRayVectorArray, - normal: na.AbstractCartesian3dVectorArray, + wavelength: u.Quantity | na.AbstractScalar, + direction: None | na.AbstractCartesian3dVectorArray = None, + n: complex | na.AbstractScalar = 1, + normal: None | na.AbstractCartesian3dVectorArray = None, ) -> na.AbstractScalar: - """ + r""" Compute the quantum efficiency of this CCD material using :meth:`quantum_efficiency_effective` and :meth:`quantum_yield_ideal`. Parameters ---------- - rays - The light rays incident on the CCD surface + wavelength + The wavelength of the incident light in vacuum. + direction + The propagation direction of the incident light in the ambient medium. + If :obj:`None` (default), normal incidence (:math:`\hat{z}`) is assumed. + n + The complex index of refraction of the ambient medium. normal The vector perpendicular to the surface of the CCD. """ - iqy = self.quantum_yield_ideal(rays.wavelength) - eqe = self.quantum_efficiency_effective(rays, normal) + iqy = self.quantum_yield_ideal(wavelength) + eqe = self.quantum_efficiency_effective( + wavelength=wavelength, + direction=direction, + n=n, + normal=normal, + ) return iqy * eqe def probability_measurement( self, - rays: optika.rays.AbstractRayVectorArray, - normal: na.AbstractCartesian3dVectorArray, + wavelength: u.Quantity | na.AbstractScalar, + direction: None | na.AbstractCartesian3dVectorArray = None, + normal: None | na.AbstractCartesian3dVectorArray = None, ) -> na.AbstractScalar: - """ + r""" Compute the probability of measuring an absorbed photon for this sensor using :func:`probability_measurement`. Parameters ---------- - rays - The light rays incident on the CCD surface + wavelength + The wavelength of the incident light in vacuum. + direction + The propagation direction of the incident light in the ambient medium. + If :obj:`None` (default), normal incidence (:math:`\hat{z}`) is assumed. normal The vector perpendicular to the surface of the CCD. """ return probability_measurement( - iqy=self.quantum_yield_ideal(rays.wavelength), - cce=self.charge_collection_efficiency(rays, normal), + iqy=self.quantum_yield_ideal(wavelength), + cce=self.charge_collection_efficiency( + wavelength=wavelength, + direction=direction, + normal=normal, + ), + ) + + def direction_refracted( + self, + wavelength: u.Quantity | na.AbstractScalar, + direction: None | na.AbstractCartesian3dVectorArray = None, + n: complex | na.AbstractScalar = 1, + normal: None | na.AbstractCartesian3dVectorArray = None, + ) -> complex | na.AbstractScalar: + if direction is None: + direction = na.Cartesian3dVectorArray(0, 0, 1) + + if normal is None: + normal = na.Cartesian3dVectorArray(0, 0, -1) + + return optika.materials.snells_law_scalar( + cos_incidence=-direction @ normal, + index_refraction=n, + index_refraction_new=self._chemical.n(wavelength), ) def electrons_measured( self, - rays: optika.rays.RayVectorArray, - normal: na.AbstractCartesian3dVectorArray, - ) -> optika.rays.RayVectorArray: + photons_absorbed: na.AbstractScalar, + wavelength: u.Quantity | na.AbstractScalar, + direction: None | na.AbstractCartesian3dVectorArray = None, + n: complex | na.AbstractScalar = 1, + normal: None | na.AbstractCartesian3dVectorArray = None, + width_pixel: ( + u.Quantity | na.AbstractScalar | na.AbstractCartesian2dVectorArray + ) = 0 + * u.um, + axis_xy: None | tuple[str, str] = None, + wrap: bool = False, + ) -> na.AbstractScalar: """ Randomly sample the number of measured electrons given the number of absorbed photons using :func:`electrons_measured`. """ - intensity = rays.intensity - wavelength = rays.wavelength + direction_substrate = self.direction_refracted( + wavelength=wavelength, + direction=direction, + n=n, + normal=normal, + ) - electrons = electrons_measured( - photons_absorbed=intensity, + return electrons_measured( + photons_absorbed=photons_absorbed, wavelength=wavelength, - absorption=self._chemical.absorption(wavelength), + absorption=absorption_effective( + wavelength=wavelength, + n_substrate=self._chemical.n(wavelength), + direction_substrate=direction_substrate, + ), thickness_implant=self.thickness_implant, + thickness_depletion=self.depletion.thickness, + thickness_substrate=self.thickness_substrate, + width_pixel=width_pixel, cce_backsurface=self.cce_backsurface, + axis_xy=axis_xy, + wrap=wrap, temperature=self.temperature, ) - result = dataclasses.replace(rays, intensity=electrons) - - return result - def signal( self, - rays: optika.rays.RayVectorArray, - normal: na.AbstractCartesian3dVectorArray, + photons: u.Quantity | na.AbstractScalar, + wavelength: u.Quantity | na.AbstractScalar, + direction: float | na.AbstractScalar = 1, + width_pixel: ( + u.Quantity | na.AbstractScalar | na.AbstractCartesian2dVectorArray + ) = 0 + * u.um, + axis_xy: None | tuple[str, str] = None, noise: bool = True, - ) -> optika.rays.RayVectorArray: - intensity = rays.intensity - wavelength = rays.wavelength + wrap: bool = False, + ) -> na.AbstractScalar: - if not intensity.unit.is_equivalent(u.photon): + if not photons.unit.is_equivalent(u.photon): h = astropy.constants.h c = astropy.constants.c - intensity = intensity / (h * c / wavelength) * u.photon + photons = photons / (h * c / wavelength) * u.photon if noise: - method = "exact" + method = "monte-carlo" else: method = "expected" - electrons = signal( - photons_expected=intensity, + # `direction` is the cosine of the refracted angle *inside* the substrate + # (computed by the sensor's `collect` method), so use equal ambient and + # substrate indices to make the Snell refraction inside `signal` a no-op + # while still computing the effective absorption with the correct index. + n_substrate = self._chemical.n(wavelength) + + return signal( + photons_expected=photons, wavelength=wavelength, + direction=direction, + n=n_substrate, + n_substrate=n_substrate, absorbance=1, - absorption=self._chemical.absorption(wavelength), thickness_implant=self.thickness_implant, + thickness_depletion=self.depletion.thickness, + thickness_substrate=self.thickness_substrate, + width_pixel=width_pixel, cce_backsurface=self.cce_backsurface, temperature=self.temperature, method=method, + axis_xy=axis_xy, + wrap=wrap, ) - result = dataclasses.replace(rays, intensity=electrons) - - return result - def photons_incident( self, electrons: u.Quantity | na.AbstractScalar, wavelength: u.Quantity | na.AbstractScalar, - direction: na.AbstractCartesian3dVectorArray, - normal: na.AbstractCartesian3dVectorArray, + direction: None | na.AbstractCartesian3dVectorArray = None, + n: complex | na.AbstractScalar = 1, + normal: None | na.AbstractCartesian3dVectorArray = None, ) -> na.AbstractScalar: """ Compute the expected number of incident photons for a given number @@ -1725,47 +2091,30 @@ def photons_incident( The assumed wavelength of the incident photons. direction The assumed direction of the incident photons. + n + The complex index of refraction of the ambient medium. normal The vector perpendicular to the surface of the sensor. """ qe = self.quantum_efficiency( - rays=optika.rays.RayVectorArray( - wavelength=wavelength, - direction=direction, - ), + wavelength=wavelength, + direction=direction, + n=n, normal=normal, ) return electrons / qe - def charge_diffusion( - self, - rays: optika.rays.RayVectorArray, - normal: na.AbstractCartesian3dVectorArray, - ) -> optika.rays.RayVectorArray: - width = self.width_charge_diffusion(rays, normal) - - position = dataclasses.replace( - rays.position, - x=na.random.normal(rays.position.x, width), - y=na.random.normal(rays.position.y, width), - ) - - rays = dataclasses.replace( - rays, - position=position, - ) - - return rays - def efficiency( self, rays: optika.rays.AbstractRayVectorArray, normal: na.AbstractCartesian3dVectorArray, ) -> na.ScalarLike: result = self.absorbance( - rays=rays, + wavelength=rays.wavelength, + direction=rays.direction, + n=rays.n, normal=normal, ) @@ -1848,7 +2197,7 @@ def eqe_rms_difference(x: tuple[float, float, float, float]): ) = x qe_fit = quantum_efficiency_effective( wavelength=eqe_measured.inputs, - direction=na.Cartesian3dVectorArray(0, 0, 1), + direction=1, thickness_oxide=thickness_oxide << unit_thickness_oxide, thickness_implant=thickness_implant << unit_thickness_implant, thickness_substrate=thickness_substrate, diff --git a/optika/sensors/materials/_materials_test.py b/optika/sensors/materials/_materials_test.py index 1f662c2c..4367359f 100644 --- a/optika/sensors/materials/_materials_test.py +++ b/optika/sensors/materials/_materials_test.py @@ -89,12 +89,20 @@ def test_transmittance( 10 * u.um, ], ) +@pytest.mark.parametrize( + argnames="method", + argvalues=[ + "Beer-Lambert", + "exact", + ], +) def test_absorbance( wavelength: u.Quantity | na.AbstractScalar, direction: float | na.AbstractScalar, n: float | na.AbstractScalar, thickness_oxide: u.Quantity | na.AbstractScalar, thickness_substrate: u.Quantity | na.AbstractScalar, + method: str, ): result = optika.sensors.absorbance( wavelength=wavelength, @@ -102,6 +110,7 @@ def test_absorbance( n=n, thickness_oxide=thickness_oxide, thickness_substrate=thickness_substrate, + method=method, ) assert np.all(result >= 0) @@ -127,23 +136,15 @@ def test_absorbance( 1, ], ) -@pytest.mark.parametrize( - argnames="cos_incidence", - argvalues=[ - 1, - ], -) def test_charge_collection_efficiency( absorption: u.Quantity | na.AbstractScalar, thickness_implant: u.Quantity | na.AbstractScalar, cce_backsurface: u.Quantity | na.AbstractScalar, - cos_incidence: float | na.AbstractScalar, ): result = optika.sensors.charge_collection_efficiency( absorption=absorption, thickness_implant=thickness_implant, cce_backsurface=cce_backsurface, - cos_incidence=cos_incidence, ) assert np.all(result >= 0) @@ -160,8 +161,8 @@ def test_charge_collection_efficiency( @pytest.mark.parametrize( argnames="direction", argvalues=[ - None, - na.Cartesian3dVectorArray(0, 0, 1), + 1, + 0.5, ], ) @pytest.mark.parametrize( @@ -191,7 +192,7 @@ def test_charge_collection_efficiency( ) def test_quantum_efficiency_effective( wavelength: u.Quantity | na.AbstractScalar, - direction: na.AbstractCartesian3dVectorArray, + direction: float | na.AbstractScalar, thickness_oxide: u.Quantity | na.AbstractScalar, thickness_implant: u.Quantity | na.AbstractScalar, thickness_substrate: u.Quantity | na.AbstractScalar, @@ -299,8 +300,8 @@ def test_electrons_measured_approx( @pytest.mark.parametrize( argnames="method", argvalues=[ - "exact", - "approx", + "monte-carlo", + "expected", ], ) def test_signal( @@ -336,19 +337,22 @@ class AbstractTestAbstractSensorMaterial( AbstractTestAbstractMaterial, ): @pytest.mark.parametrize( - argnames="rays", + argnames="photons", argvalues=[ - optika.rays.RayVectorArray( - intensity=1e-6 * u.erg, - wavelength=100 * u.AA, - direction=na.Cartesian3dVectorArray(0, 0, 1), - ), + 1e-6 * u.erg, ], ) @pytest.mark.parametrize( - argnames="normal", + argnames="wavelength", argvalues=[ - na.Cartesian3dVectorArray(0, 0, -1), + 100 * u.AA, + ], + ) + @pytest.mark.parametrize( + argnames="direction", + argvalues=[ + 1, + 0.5, ], ) @pytest.mark.parametrize( @@ -358,13 +362,19 @@ class AbstractTestAbstractSensorMaterial( def test_signal( self, a: optika.sensors.materials.AbstractSensorMaterial, - rays: optika.rays.RayVectorArray, - normal: na.AbstractCartesian3dVectorArray, + photons: u.Quantity | na.AbstractScalar, + wavelength: u.Quantity | na.AbstractScalar, + direction: float | na.AbstractScalar, noise: bool, ): - result = a.signal(rays, normal, noise=noise) - assert isinstance(result, optika.rays.RayVectorArray) - assert np.all(result.intensity >= 0 * u.electron) + result = a.signal( + photons=photons, + wavelength=wavelength, + direction=direction, + noise=noise, + ) + assert isinstance(na.as_named_array(result), na.AbstractScalar) + assert np.all(result >= 0 * u.electron) @pytest.mark.parametrize( argnames="electrons", @@ -408,30 +418,21 @@ def test_photons_incident( assert result.unit.is_equivalent(u.photon) @pytest.mark.parametrize( - argnames="rays", - argvalues=[ - optika.rays.RayVectorArray( - intensity=10 * u.electron, - wavelength=100 * u.AA, - position=na.Cartesian3dVectorArray() * u.mm, - direction=na.Cartesian3dVectorArray(0, 0, 1), - ), - ], - ) - @pytest.mark.parametrize( - argnames="normal", + argnames="wavelength", argvalues=[ - na.Cartesian3dVectorArray(0, 0, -1), + 100 * u.AA, ], ) - def test_charge_diffusion( + def test_direction_refracted( self, a: optika.sensors.materials.AbstractSensorMaterial, - rays: optika.rays.RayVectorArray, - normal: na.AbstractCartesian3dVectorArray, + wavelength: u.Quantity | na.AbstractScalar, ): - result = a.charge_diffusion(rays, normal) - assert isinstance(result, optika.rays.RayVectorArray) + # with the default direction and normal (normal incidence) the + # refracted cosine is unity + result = a.direction_refracted(wavelength=wavelength) + assert isinstance(na.as_named_array(result), na.AbstractScalar) + assert np.all(np.real(result) > 0) @pytest.mark.parametrize( @@ -522,154 +523,223 @@ def test_depletion( ) @pytest.mark.parametrize( - argnames="rays", + argnames="wavelength", argvalues=[ - optika.rays.RayVectorArray( - wavelength=100 * u.AA, - position=na.Cartesian3dVectorArray() * u.mm, - direction=na.Cartesian3dVectorArray(0, 0, 1), - ), + 100 * u.AA, + ], + ) + def test_width_charge_diffusion( + self, + a: optika.sensors.materials.AbstractBackIlluminatedSiliconSensorMaterial, + wavelength: u.Quantity | na.AbstractScalar, + ): + result = a.width_charge_diffusion(wavelength=wavelength) + assert np.all(result >= 0 * u.um) + + @pytest.mark.parametrize( + argnames="wavelength", + argvalues=[ + 100 * u.AA, + ], + ) + @pytest.mark.parametrize( + argnames="direction", + argvalues=[ + None, ], ) @pytest.mark.parametrize( argnames="normal", argvalues=[ - na.Cartesian3dVectorArray(0, 0, -1), + None, ], ) - def test_width_charge_diffusion( + def test_transmittance( self, a: optika.sensors.materials.AbstractBackIlluminatedSiliconSensorMaterial, - rays: optika.rays.AbstractRayVectorArray, - normal: na.AbstractCartesian3dVectorArray, + wavelength: u.Quantity | na.AbstractScalar, + direction: None | na.AbstractCartesian3dVectorArray, + normal: None | na.AbstractCartesian3dVectorArray, ): - result = a.width_charge_diffusion(rays, normal) - assert np.all(result >= 0 * u.um) + result = a.transmittance( + wavelength=wavelength, + direction=direction, + normal=normal, + ) + assert np.all(result >= 0) + assert np.all(result <= 1) @pytest.mark.parametrize( - argnames="rays", + argnames="wavelength", argvalues=[ - optika.rays.RayVectorArray( - wavelength=100 * u.AA, - direction=na.Cartesian3dVectorArray(0, 0, 1), - ), + 100 * u.AA, + ], + ) + @pytest.mark.parametrize( + argnames="direction", + argvalues=[ + None, ], ) @pytest.mark.parametrize( argnames="normal", argvalues=[ - na.Cartesian3dVectorArray(0, 0, -1), + None, ], ) def test_absorbance( self, a: optika.sensors.materials.AbstractBackIlluminatedSiliconSensorMaterial, - rays: optika.rays.AbstractRayVectorArray, - normal: na.AbstractCartesian3dVectorArray, + wavelength: u.Quantity | na.AbstractScalar, + direction: None | na.AbstractCartesian3dVectorArray, + normal: None | na.AbstractCartesian3dVectorArray, ): - result = a.absorbance(rays, normal) + result = a.absorbance( + wavelength=wavelength, + direction=direction, + normal=normal, + ) assert np.all(result >= 0) assert np.all(result <= 1) @pytest.mark.parametrize( - argnames="rays", + argnames="wavelength", argvalues=[ - optika.rays.RayVectorArray( - wavelength=100 * u.AA, - direction=na.Cartesian3dVectorArray(0, 0, 1), - ), + 100 * u.AA, + ], + ) + @pytest.mark.parametrize( + argnames="direction", + argvalues=[ + None, ], ) @pytest.mark.parametrize( argnames="normal", argvalues=[ - na.Cartesian3dVectorArray(0, 0, -1), + None, ], ) def test_charge_collection_efficiency( self, a: optika.sensors.materials.AbstractBackIlluminatedSiliconSensorMaterial, - rays: optika.rays.AbstractRayVectorArray, - normal: na.AbstractCartesian3dVectorArray, + wavelength: u.Quantity | na.AbstractScalar, + direction: None | na.AbstractCartesian3dVectorArray, + normal: None | na.AbstractCartesian3dVectorArray, ): - result = a.charge_collection_efficiency(rays, normal) + result = a.charge_collection_efficiency( + wavelength=wavelength, + direction=direction, + normal=normal, + ) assert np.all(result >= 0) assert np.all(result <= 1) @pytest.mark.parametrize( - argnames="rays", + argnames="wavelength", + argvalues=[ + 100 * u.AA, + ], + ) + @pytest.mark.parametrize( + argnames="direction", argvalues=[ - optika.rays.RayVectorArray( - wavelength=100 * u.AA, - direction=na.Cartesian3dVectorArray(0, 0, 1), - ), + None, ], ) @pytest.mark.parametrize( argnames="normal", argvalues=[ - na.Cartesian3dVectorArray(0, 0, -1), + None, ], ) def test_quantum_efficiency( self, a: optika.sensors.materials.AbstractBackIlluminatedSiliconSensorMaterial, - rays: optika.rays.AbstractRayVectorArray, - normal: na.AbstractCartesian3dVectorArray, + wavelength: u.Quantity | na.AbstractScalar, + direction: None | na.AbstractCartesian3dVectorArray, + normal: None | na.AbstractCartesian3dVectorArray, ): - result = a.quantum_efficiency(rays, normal) + result = a.quantum_efficiency( + wavelength=wavelength, + direction=direction, + normal=normal, + ) assert result > 0 * u.electron / u.photon @pytest.mark.parametrize( - argnames="rays", + argnames="wavelength", argvalues=[ - optika.rays.RayVectorArray( - wavelength=100 * u.AA, - direction=na.Cartesian3dVectorArray(0, 0, 1), - ), + 100 * u.AA, + ], + ) + @pytest.mark.parametrize( + argnames="direction", + argvalues=[ + None, ], ) @pytest.mark.parametrize( argnames="normal", argvalues=[ - na.Cartesian3dVectorArray(0, 0, -1), + None, ], ) def test_probability_measurement( self, a: optika.sensors.materials.AbstractBackIlluminatedSiliconSensorMaterial, - rays: optika.rays.AbstractRayVectorArray, - normal: na.AbstractCartesian3dVectorArray, + wavelength: u.Quantity | na.AbstractScalar, + direction: None | na.AbstractCartesian3dVectorArray, + normal: None | na.AbstractCartesian3dVectorArray, ): - result = a.probability_measurement(rays, normal) + result = a.probability_measurement( + wavelength=wavelength, + direction=direction, + normal=normal, + ) assert np.all(result >= 0) assert np.all(result <= 1) @pytest.mark.parametrize( - argnames="rays", + argnames="photons_absorbed", argvalues=[ - optika.rays.RayVectorArray( - intensity=100 * u.photon, - wavelength=100 * u.AA, - direction=na.Cartesian3dVectorArray(0, 0, 1), - ), + (100 * u.photon).astype(int), + ], + ) + @pytest.mark.parametrize( + argnames="wavelength", + argvalues=[ + 100 * u.AA, + ], + ) + @pytest.mark.parametrize( + argnames="direction", + argvalues=[ + None, ], ) @pytest.mark.parametrize( argnames="normal", argvalues=[ - na.Cartesian3dVectorArray(0, 0, -1), + None, ], ) def test_electrons_measured( self, a: optika.sensors.materials.AbstractBackIlluminatedSiliconSensorMaterial, - rays: optika.rays.RayVectorArray, - normal: na.AbstractCartesian3dVectorArray, + photons_absorbed: u.Quantity | na.AbstractScalar, + wavelength: u.Quantity | na.AbstractScalar, + direction: None | na.AbstractCartesian3dVectorArray, + normal: None | na.AbstractCartesian3dVectorArray, ): - result = a.electrons_measured(rays, normal) - assert isinstance(result, optika.rays.RayVectorArray) - assert np.all(result.intensity >= 0 * u.electron) + result = a.electrons_measured( + photons_absorbed=photons_absorbed, + wavelength=wavelength, + direction=direction, + normal=normal, + ) + assert isinstance(na.as_named_array(result), na.AbstractScalar) + assert np.all(result >= 0 * u.electron) @pytest.mark.parametrize( diff --git a/optika/sensors/materials/_ramanathan_2020/_ramanathan_2020.py b/optika/sensors/materials/_ramanathan_2020/_ramanathan_2020.py index af58fa3c..b7432f33 100644 --- a/optika/sensors/materials/_ramanathan_2020/_ramanathan_2020.py +++ b/optika/sensors/materials/_ramanathan_2020/_ramanathan_2020.py @@ -8,6 +8,8 @@ import optika from .._stern_1994 import ( _thickness_implant, + _thickness_substrate, + _width_pixel, _cce_backsurface, ) @@ -468,8 +470,15 @@ def electrons_measured( wavelength: u.Quantity | na.ScalarArray, absorption: None | u.Quantity | na.AbstractScalar = None, thickness_implant: u.Quantity | na.AbstractScalar = _thickness_implant, + thickness_depletion: None | u.Quantity | na.AbstractScalar = None, + thickness_substrate: u.Quantity | na.AbstractScalar = _thickness_substrate, + width_pixel: ( + u.Quantity | na.AbstractScalar | na.AbstractCartesian2dVectorArray + ) = _width_pixel, cce_backsurface: u.Quantity | na.AbstractScalar = _cce_backsurface, temperature: u.Quantity | na.ScalarArray = 300 * u.K, + axis_xy: None | tuple[str, str] = None, + wrap: bool = False, shape_random: None | dict[str, int] = None, ) -> na.AbstractScalar: r""" @@ -483,26 +492,53 @@ def electrons_measured( Parameters ---------- photons_absorbed - The number of photons absorbed by the light-sensitive layer of the sensor. + The number of photons absorbed by the light-sensitive layer of the + sensor. wavelength The vacuum wavelength of the absorbed photons. absorption - The absorption coefficient of silicon at the given wavelength. + The absorption coefficient of silicon per unit perpendicular depth. + For oblique incidence, supply the effective coefficient from + :func:`optika.sensors.absorption_effective`, which folds in the + refracted angle, so no separate angle argument is needed. thickness_implant The thickness of the implant layer, where partial-charge collection occurs. + thickness_depletion + The thickness of the depletion region, the region with significant electric + field. + If :obj:`None` (the default), this is set to the same value as + `thickness_substrate`. + thickness_substrate + The thickness of the entire light-sensitive region of the device. + width_pixel + The size of a single pixel on the sensor. + A scalar gives square pixels; a + :class:`named_arrays.AbstractCartesian2dVectorArray` whose ``x``/``y`` + components are the pixel widths along ``axis_xy[0]``/``axis_xy[1]`` + gives rectangular pixels. cce_backsurface The differential charge collection efficiency on the back surface of the sensor. temperature The temperature of the silicon detector. Default is room temperature. + axis_xy + The two logical axes corresponding to the pixel grid of the sensor + along which electrons will diffuse. + If :obj:`None` (the default), there is no charge diffusion. + wrap + Controls how diffused charge is treated at the edges of the pixel grid. + If :obj:`False` (the default), charge that diffuses past the edge of the + grid is lost, as it would be at the physical edge of a sensor. + If :obj:`True`, the grid is treated as periodic and the charge re-enters + the opposite edge (a toroidal boundary). shape_random Additional shape used to specify the number of samples to draw. Examples -------- - Plot the energy spectrum of 100 6 keV photons emitted from an Fe-55 + Plot the energy spectrum of twenty 6 keV photons emitted from an Fe-55 radioactive source. .. jupyter-execute:: @@ -518,7 +554,7 @@ def electrons_measured( # Define the expected number of photons # for each experiment - photons_absorbed = (100 * u.photon).astype(int) + photons_absorbed = (20 * u.photon).astype(int) # Define the wavelength at which to sample the distribution wavelength = 5.9 * u.keV @@ -561,22 +597,49 @@ def electrons_measured( if absorption is None: absorption = optika.chemicals.Chemical("Si").absorption(wavelength) + if thickness_depletion is None: + thickness_depletion = thickness_substrate + if shape_random is None: shape_random = dict() + if not isinstance(width_pixel, na.AbstractCartesian2dVectorArray): + width_pixel = na.Cartesian2dVectorArray(width_pixel, width_pixel) + + width_pixel_x = width_pixel.x + width_pixel_y = width_pixel.y + shape = na.broadcast_shapes( na.shape(photons_absorbed), na.shape(wavelength), na.shape(absorption), na.shape(thickness_implant), + na.shape(thickness_depletion), + na.shape(thickness_substrate), + na.shape(width_pixel_x), + na.shape(width_pixel_y), na.shape(cce_backsurface), na.shape(temperature), shape_random, ) + if axis_xy is not None: + axis_x, axis_y = axis_xy + shape[axis_x] = shape.pop(axis_x) + shape[axis_y] = shape.pop(axis_y) + else: + axis_x = "__dummy_x__" + axis_y = "__dummy_y__" + shape[axis_x] = 1 + shape[axis_y] = 1 + photons_absorbed = na.broadcast_to(photons_absorbed, shape) absorption = na.broadcast_to(absorption, shape) thickness_implant = na.broadcast_to(thickness_implant, shape) + thickness_depletion = na.broadcast_to(thickness_depletion, shape) + thickness_substrate = na.broadcast_to(thickness_substrate, shape) + width_pixel_x = na.broadcast_to(width_pixel_x, shape) + width_pixel_y = na.broadcast_to(width_pixel_y, shape) cce_backsurface = na.broadcast_to(cce_backsurface, shape) if not isinstance(cce_backsurface, u.Quantity): @@ -604,58 +667,89 @@ def electrons_measured( wavelength=wavelength.ndarray, absorption=absorption.ndarray, thickness_implant=thickness_implant.ndarray, + thickness_depletion=thickness_depletion.ndarray, + thickness_substrate=thickness_substrate.ndarray, + width_pixel_x=width_pixel_x.ndarray, + width_pixel_y=width_pixel_y.ndarray, cce_backsurface=cce_backsurface.ndarray, p_n=p_n.ndarray, n=n.ndarray, energy_pair_inf=energy_pair_inf.ndarray, fano_inf=fano_inf.ndarray, + wrap=wrap, ) - return na.ScalarArray( + result = na.ScalarArray( ndarray=result, axes=tuple(shape), ) + if axis_xy is None: + result = result[{axis_x: 0, axis_y: 0}] + + return result + def _electrons_measured_quantity( photons_absorbed: u.Quantity, wavelength: u.Quantity, absorption: u.Quantity, thickness_implant: u.Quantity, + thickness_depletion: u.Quantity, + thickness_substrate: u.Quantity, + width_pixel_x: u.Quantity, + width_pixel_y: u.Quantity, cce_backsurface: u.Quantity, p_n: np.ndarray, n: np.ndarray, energy_pair_inf: u.Quantity, fano_inf: u.Quantity, + wrap: bool, ) -> u.Quantity: shape = np.broadcast_shapes( photons_absorbed.shape, absorption.shape, thickness_implant.shape, + thickness_depletion.shape, + thickness_substrate.shape, cce_backsurface.shape, + width_pixel_x.shape, + width_pixel_y.shape, ) + num_x = shape[~1] + num_y = shape[~0] + unit_length = u.mm photons_absorbed = photons_absorbed.to_value(u.photon) energy = wavelength.to(u.eV, equivalencies=u.spectral()) absorption = absorption.to_value(1 / unit_length) thickness_implant = thickness_implant.to_value(unit_length) + thickness_depletion = thickness_depletion.to_value(unit_length) + thickness_substrate = thickness_substrate.to_value(unit_length) + width_pixel_x = width_pixel_x.to_value(unit_length) + width_pixel_y = width_pixel_y.to_value(unit_length) cce_backsurface = cce_backsurface.to_value(u.dimensionless_unscaled) energy_pair_inf = energy_pair_inf.to_value(u.eV) fano_inf = fano_inf.to_value(u.electron / u.photon) result = _electrons_measured_numba( - photons_absorbed=photons_absorbed.reshape(-1), - energy=energy.reshape(-1), - absorption=absorption.reshape(-1), - thickness_implant=thickness_implant.reshape(-1), - cce_backsurface=cce_backsurface.reshape(-1), - p_n=p_n.reshape(-1, p_n.shape[~0]), - n=n.reshape(-1, n.shape[~0]), - energy_pair_inf=energy_pair_inf.reshape(-1), - fano_inf=fano_inf.reshape(-1), + photons_absorbed=photons_absorbed.reshape(-1, num_x, num_y), + energy=energy.reshape(-1, num_x, num_y), + absorption=absorption.reshape(-1, num_x, num_y), + thickness_implant=thickness_implant.reshape(-1, num_x, num_y), + thickness_depletion=thickness_depletion.reshape(-1, num_x, num_y), + thickness_substrate=thickness_substrate.reshape(-1, num_x, num_y), + width_pixel_x=width_pixel_x.reshape(-1, num_x, num_y), + width_pixel_y=width_pixel_y.reshape(-1, num_x, num_y), + cce_backsurface=cce_backsurface.reshape(-1, num_x, num_y), + p_n=p_n.reshape(-1, num_x, num_y, p_n.shape[~0]), + n=n.reshape(-1, num_x, num_y, n.shape[~0]), + energy_pair_inf=energy_pair_inf.reshape(-1, num_x, num_y), + fano_inf=fano_inf.reshape(-1, num_x, num_y), + wrap=wrap, ) result = result.reshape(shape) @@ -665,75 +759,118 @@ def _electrons_measured_quantity( return result -@numba.njit(fastmath=True, parallel=True) -def _electrons_measured_numba( +@numba.njit( + cache=True, + fastmath=True, + parallel=True, +) +def _electrons_measured_numba( # pragma: nocover photons_absorbed: np.ndarray, energy: np.ndarray, absorption: np.ndarray, thickness_implant: np.ndarray, + thickness_depletion: np.ndarray, + thickness_substrate: np.ndarray, + width_pixel_x: np.ndarray, + width_pixel_y: np.ndarray, cce_backsurface: np.ndarray, p_n: np.ndarray, n: np.ndarray, energy_pair_inf: np.ndarray, fano_inf: np.ndarray, -) -> np.ndarray: # pragma: nocover + wrap: bool, +) -> np.ndarray: - num_i, num_n = p_n.shape + num_i, num_x, num_y, num_n = p_n.shape - result = np.empty(num_i) + result = np.zeros((num_i, num_x, num_y)) for i in numba.prange(num_i): - - num_photon = int(photons_absorbed[i]) - energy_i = energy[i] - a = absorption[i] - W = thickness_implant[i] - h_0 = cce_backsurface[i] - cmf_i = np.cumsum(p_n[i]) - n_i = n[i] - energy_pair_inf_i = energy_pair_inf[i] - fano_inf_i = fano_inf[i] - - d = 1 / a - - num_electron_i = 0 - - mean_inf = energy_i / energy_pair_inf_i - std_inf = math.sqrt(fano_inf_i * mean_inf) - - low_energy = energy_i <= 50 - - for j in range(num_photon): - - if low_energy: - - x_ij = random.uniform(0, 1) - - for k_ij in range(num_n): - if cmf_i[k_ij] > x_ij: - break - - n_ij = n_i[k_ij] - - else: - n_ij = random.normalvariate( - mu=mean_inf, - sigma=std_inf, - ) - n_ij = round(n_ij) - - y_ij = random.uniform(0, 1) - z_ij = -d * math.log(1 - y_ij) - - if z_ij < W: - h_ij = h_0 + (1 - h_0) * z_ij / W - else: - h_ij = 1 - - m_ij = np.random.binomial(n=n_ij, p=h_ij) - - num_electron_i = num_electron_i + m_ij - - result[i] = num_electron_i + for x in range(num_x): + for y in range(num_y): + num_photon = int(photons_absorbed[i, x, y]) + energy_i = energy[i, x, y] + a = absorption[i, x, y] + W = thickness_implant[i, x, y] + h_0 = cce_backsurface[i, x, y] + cmf_i = np.cumsum(p_n[i, x, y]) + n_i = n[i, x, y] + energy_pair_inf_i = energy_pair_inf[i, x, y] + fano_inf_i = fano_inf[i, x, y] + z_substrate = thickness_substrate[i, x, y] + z_ff = z_substrate - thickness_depletion[i, x, y] + wp_x = width_pixel_x[i, x, y] + wp_y = width_pixel_y[i, x, y] + + # guard against a non-absorbing medium (a == 0), where the + # reciprocal absorption length is undefined + d = 1 / a if a > 0 else 0.0 + + # every photon is absorbed within the substrate, so sample the + # depth from the exponential truncated to [0, z_substrate) + fraction_absorbed = 1 - math.exp(-a * z_substrate) + + mean_inf = energy_i / energy_pair_inf_i + std_inf = math.sqrt(fano_inf_i * mean_inf) + + low_energy = energy_i <= 50 + + for j in range(num_photon): + if low_energy: + x_ij = random.uniform(0, 1) + + for k_ij in range(num_n): + if cmf_i[k_ij] > x_ij: + break + + n_ij = n_i[k_ij] + + else: + n_ij = random.normalvariate( + mu=mean_inf, + sigma=std_inf, + ) + n_ij = round(n_ij) + + y_ij = random.uniform(0, 1) + if a > 0: + z_ij = -d * math.log(1 - y_ij * fraction_absorbed) + else: + # with no absorption the truncated-exponential depth + # distribution reduces to a uniform one over [0, z_substrate) + z_ij = y_ij * z_substrate + + if z_ij < W: + h_ij = h_0 + (1 - h_0) * z_ij / W + else: + h_ij = 1 + + m_ij = np.random.binomial(n=n_ij, p=h_ij) + + u = random.uniform(-0.5, 0.5) + v = random.uniform(-0.5, 0.5) + + for e in range(m_ij): + if z_ij < z_ff and wp_x > 0 and wp_y > 0: + w = z_ff * math.sqrt(1 - z_ij / z_ff) + + p = random.gauss(u, w / wp_x) + q = random.gauss(v, w / wp_y) + + p = round(p) + q = round(q) + + else: + p = q = 0 + + x_e = x + p + y_e = y + q + + if wrap: + # toroidal boundary: charge re-enters the opposite edge + result[i, x_e % num_x, y_e % num_y] += 1 + elif (0 <= x_e < num_x) and (0 <= y_e < num_y): + result[i, x_e, y_e] += 1 + # otherwise the electron diffused off the sensor and is lost return result diff --git a/optika/sensors/materials/_ramanathan_2020/_ramanathan_2020_test.py b/optika/sensors/materials/_ramanathan_2020/_ramanathan_2020_test.py index 6ef10eb3..0dda2fc0 100644 --- a/optika/sensors/materials/_ramanathan_2020/_ramanathan_2020_test.py +++ b/optika/sensors/materials/_ramanathan_2020/_ramanathan_2020_test.py @@ -2,6 +2,7 @@ import numpy as np import astropy.units as u import named_arrays as na +import optika from . import _ramanathan_2020 _wavelength = [ @@ -100,7 +101,14 @@ def test_fano_factor_inf( @pytest.mark.parametrize( argnames="photons_absorbed", argvalues=[ - 100 * u.photon, + na.broadcast_to((100 * u.photon).astype(int), dict(pixel_x=2, pixel_y=2)), + ], +) +@pytest.mark.parametrize( + argnames="axis_xy", + argvalues=[ + None, + ("pixel_x", "pixel_y"), ], ) @pytest.mark.parametrize( @@ -125,6 +133,12 @@ def test_fano_factor_inf( 0.5, ], ) +@pytest.mark.parametrize( + argnames="width_pixel", + argvalues=[ + 15 * u.um, + ], +) @pytest.mark.parametrize("temperature", _temperture) def test_electrons_measured( photons_absorbed: u.Quantity | na.AbstractScalar, @@ -133,14 +147,18 @@ def test_electrons_measured( thickness_implant: u.Quantity | na.AbstractScalar, cce_backsurface: u.Quantity | na.AbstractScalar, temperature: u.Quantity | na.ScalarArray, + width_pixel: u.Quantity | na.AbstractCartesian2dVectorArray, + axis_xy: None | tuple[str, str], ): result = _ramanathan_2020.electrons_measured( photons_absorbed=photons_absorbed, wavelength=wavelength, absorption=absorption, thickness_implant=thickness_implant, + width_pixel=width_pixel, cce_backsurface=cce_backsurface, temperature=temperature, + axis_xy=axis_xy, ) assert np.all(result >= 0 * u.electron) @@ -155,3 +173,90 @@ def test_electrons_measured( ) assert result.shape == shape + + +def test_electrons_measured_diffusion(): + """ + Inject many photons into a single pixel and check that the spatial spread + of the diffused charge matches the analytic charge-diffusion width given by + :func:`optika.sensors.charge_diffusion`. + """ + num = 41 + axis_xy = ("pixel_x", "pixel_y") + + absorption = 1 / u.um + thickness_substrate = 14 * u.um + thickness_depletion = 0 * u.um + width_pixel = 3 * u.um + + # Place all of the photons in the central pixel of an otherwise empty grid. + photons = np.zeros((num, num)) + photons[num // 2, num // 2] = 20000 + photons = na.ScalarArray(photons << u.photon, axes=axis_xy).astype(int) + + electrons = _ramanathan_2020.electrons_measured( + photons_absorbed=photons, + wavelength=500 * u.nm, + absorption=absorption, + thickness_implant=0 * u.um, + thickness_depletion=thickness_depletion, + thickness_substrate=thickness_substrate, + width_pixel=width_pixel, + cce_backsurface=1, + axis_xy=axis_xy, + ) + + # The diffused charge should spread beyond the central pixel. + assert electrons[{axis_xy[0]: num // 2, axis_xy[1]: num // 2}] < electrons.sum( + axis_xy + ) + + # Physical offset of each pixel from the center of the grid. + offset_x = (na.arange(0, num, axis=axis_xy[0]) - num // 2) * width_pixel + offset_y = (na.arange(0, num, axis=axis_xy[1]) - num // 2) * width_pixel + + total = electrons.sum(axis_xy) + mean_x = (electrons * offset_x).sum(axis_xy) / total + mean_y = (electrons * offset_y).sum(axis_xy) / total + var_x = (electrons * np.square(offset_x - mean_x)).sum(axis_xy) / total + var_y = (electrons * np.square(offset_y - mean_y)).sum(axis_xy) / total + std_measured = np.sqrt((var_x + var_y) / 2) + + std_expected = optika.sensors.charge_diffusion( + absorption=absorption, + thickness_substrate=thickness_substrate, + thickness_depletion=thickness_depletion, + ) + + assert np.allclose(std_measured, std_expected, rtol=0.05) + + +def test_electrons_measured_wrap(): + """ + On a grid small compared to the diffusion width, charge that diffuses off + the grid is lost when ``wrap=False`` but retained (toroidally) when + ``wrap=True``, so the wrapped grid holds strictly more charge. + """ + num = 3 + axis_xy = ("pixel_x", "pixel_y") + + photons = np.zeros((num, num)) + photons[num // 2, num // 2] = 5000 + photons = na.ScalarArray(photons << u.photon, axes=axis_xy).astype(int) + + kwargs = dict( + photons_absorbed=photons, + wavelength=500 * u.nm, + absorption=1 / u.um, + thickness_implant=0 * u.um, + thickness_depletion=0 * u.um, + thickness_substrate=14 * u.um, + width_pixel=2 * u.um, + cce_backsurface=1, + axis_xy=axis_xy, + ) + + total_drop = _ramanathan_2020.electrons_measured(**kwargs, wrap=False).sum(axis_xy) + total_wrap = _ramanathan_2020.electrons_measured(**kwargs, wrap=True).sum(axis_xy) + + assert total_wrap > total_drop diff --git a/optika/sensors/materials/_stern_1994.py b/optika/sensors/materials/_stern_1994.py index 99f66151..ee336e19 100644 --- a/optika/sensors/materials/_stern_1994.py +++ b/optika/sensors/materials/_stern_1994.py @@ -3,4 +3,5 @@ _thickness_oxide = 50 * u.AA _thickness_substrate = 7 * u.um _thickness_implant = 2317 * u.AA +_width_pixel = 27 * u.um _cce_backsurface = 0.21 diff --git a/optika/sensors/materials/_tektronix_tk512cb/_tektronix_tk512cb.py b/optika/sensors/materials/_tektronix_tk512cb/_tektronix_tk512cb.py index 7cd263cf..6f167983 100644 --- a/optika/sensors/materials/_tektronix_tk512cb/_tektronix_tk512cb.py +++ b/optika/sensors/materials/_tektronix_tk512cb/_tektronix_tk512cb.py @@ -101,10 +101,8 @@ def tektronix_tk512cb( # Evaluate the fitted QE using the given wavelengths eqe_fit = material.quantum_efficiency_effective( - rays=optika.rays.RayVectorArray( - wavelength=wavelength_fit, - direction=na.Cartesian3dVectorArray(0, 0, 1), - ), + wavelength=wavelength_fit, + direction=na.Cartesian3dVectorArray(0, 0, 1), normal=na.Cartesian3dVectorArray(0, 0, -1), ) @@ -167,11 +165,7 @@ def tektronix_tk512cb( # Compute the width of the charge diffusion kernel # for each wavelength. width = material.width_charge_diffusion( - rays=optika.rays.RayVectorArray( - wavelength=wavelength_fit, - direction=na.Cartesian3dVectorArray(0, 0, 1), - ), - normal=na.Cartesian3dVectorArray(0, 0, -1), + wavelength=wavelength_fit, ) # Plot the results diff --git a/optika/surfaces.py b/optika/surfaces.py index 07ec6a86..ae25f1d5 100644 --- a/optika/surfaces.py +++ b/optika/surfaces.py @@ -123,7 +123,6 @@ def is_stop(self) -> bool: def propagate_rays( self, rays: optika.rays.RayVectorArray, - # material: None | optika.materials.AbstractMaterial = None, ) -> optika.rays.RayVectorArray: r""" Refract, reflect, and/or diffract the given rays off of this surface @@ -131,7 +130,7 @@ def propagate_rays( Parameters ---------- rays - a set of input rays that will interact with this surface + A set of input rays that will interact with this surface. """ sag = self.sag material = self.material diff --git a/optika/systems/_sequential.py b/optika/systems/_sequential.py index a304f3b2..971e17a3 100644 --- a/optika/systems/_sequential.py +++ b/optika/systems/_sequential.py @@ -915,7 +915,7 @@ def image( axis=axis_wavelength, ) - return self.sensor.readout( + return self.sensor.measure( rays=rayfunction.outputs, wavelength=wavelength, axis=(axis_wavelength,) + axis_field + axis_pupil,