Merge pull request #104 from spacesys-finch/feature/snr

Created SNR vs Wavelength notebook

payload_designer/components/foreoptics.py

@@ -23,6 +23,7 @@ class Foreoptic(Component):
     def __init__(
         self,
         diameter=None,
+        image_diameter=None,
         focal_length=None,
         mass=None,
         length=None,
@@ -31,10 +32,19 @@ def __init__(
         super().__init__(mass=mass, dimensions=(diameter, diameter, length))
 
         self.diameter = diameter
+        self.image_diameter = image_diameter
         self.focal_length = focal_length
         self.length = length
         self.transmittance = transmittance
 
+    def get_image_area(self):
+        """Calculate the image area from the image diameter."""
+        assert self.image_diameter is not None, "image_diameter is not set."
+
+        area = np.pi * (self.image_diameter / 2) ** 2
+
+        return area
+
     def get_magnification(self):
         """Calculate the magnification of the foreoptics.
 
@@ -56,12 +66,7 @@ def get_f_number(self):
             float: f/# (unitless).
 
         """
-        if self.na is not None:
-            n = np.divide(1, 2 * self.na)
-        elif self.ds_i is not None and self.dm_a is not None:
-            n = np.divide(self.ds_i, self.dm_a)
-        else:
-            raise ValueError("ds_i and dm_a or na must be set.")
+        n = self.focal_length / self.diameter
 
         return n
 
@@ -132,19 +137,6 @@ def get_radiant_flux(self):
 
         return f
 
-    def get_image_area(self):
-        """Calculate the image area.
-
-        Returns:
-            float: image area [mm^2].
-
-        """
-        assert self.d_i is not None, "d_i is not set."
-
-        a_i = math.pi * (self.d_i / 2) ** 2
-
-        return a_i
-
 
 class Chromar(Foreoptic):
     def __init__(self):

payload_designer/components/masks.py

@@ -27,7 +27,7 @@ def __init__(self, mass=None, diameter=None, thickness=None, size: tuple = None)
         self.thickness = thickness
         self.diameter = diameter
 
-    def get_area(self):
+    def get_clear_area(self):
         """Get the clear aperture slit area."""
         assert self.size is not None, "Size must be specified."
 

payload_designer/components/sensors.py

@@ -65,7 +65,7 @@ def get_size(self) -> tuple:
         cubesat frame."""
         assert self.n_px is not None, "n_px must be specified."
         assert self.pitch is not None, "Pitch must be specified."
-        
+
         size = (self.n_px * self.pitch, self.n_px * self.pitch)
 
         return size
@@ -95,7 +95,7 @@ def get_dark_noise(self):
         """Get the dark noise of the sensor."""
         assert self.i_dark is not None, "i_dark must be specified."
         assert self.dt is not None, "dt amplitude must be specified."
-        
+
         dark_noise = self.i_dark * self.dt
 
         return dark_noise
@@ -104,16 +104,16 @@ def get_quantization_noise(self):
         """Get the quantization noise of the sensor."""
         assert self.n_well is not None, "n_well must be specified."
         assert self.n_bit is not None, "n_bit must be specified."
-        
-        quant_noise = (1 / math.sqrt(12)) * self.n_well / 2**self.n_bit
+
+        quant_noise = (1 / math.sqrt(12)) * self.n_well / 2**self.n_bit.value
 
         return quant_noise
 
     def get_noise(self, signal):
         """Get the net noise of the sensor."""
         assert self.n_bin is not None, "n_bin must be specified."
         assert self.noise_read is not None, "noise_read must be specified."
-        
+
         noise = np.sqrt(
             signal
             + self.n_bin * self.get_dark_noise() ** 2

payload_designer/systems/payloads.py

@@ -1,3 +1,6 @@
+# stdlib
+import math
+
 # external
 import astropy.constants as const
 import astropy.units as unit
@@ -69,6 +72,18 @@ def get_transmittance(self):
 
         return transmittance
 
+    def get_ratio_cropped_light_through_slit(self):
+        """Get the ratio of the light area passing through the slit to the area of the image of the foreoptic."""
+        assert (
+            self.foreoptic.image_diameter is not None
+        ), "Foreoptic image diameter must be set."
+
+        effective_width = min(self.slit.size[0], self.foreoptic.image_diameter)
+        effective_slit_area = effective_width * self.slit.size[1]
+        ratio = effective_slit_area / self.foreoptic.get_image_area()
+
+        return ratio
+
     def get_signal_to_noise(self, radiance: LUT, wavelength):
         """Get the signal to noise ratio of the system.
 
@@ -77,37 +92,63 @@ def get_signal_to_noise(self, radiance: LUT, wavelength):
             wavelength: Wavelength(s) at which to evaluate SNR.
 
         """
-        assert self.sensor is not None, "A sensor component must be specified." 
+        assert self.sensor is not None, "A sensor component must be specified."
         assert self.foreoptic is not None, "A foreoptic component must be specified."
         assert self.slit is not None, "A slit component must be specified."
 
-        signal = (
-            (const.pi / 4)
-            * (wavelength / (const.h * const.c))
-            * (self.sensor.get_pixel_area() / self.foreoptic.f_number**2)
-            * self.sensor.efficiency(wavelength)
-            * self.get_transmittance()
-            * self.slit.get_aperture_area()
+        signal1 = (
+            (wavelength / (const.h * const.c))
             * radiance(wavelength)
+            * (math.pi / 4)
             * self.sensor.dt
-        )
+        )  # [sr-1 m-3]
+        signal2 = self.sensor.get_pixel_area() / (
+            ((self.foreoptic.get_f_number()).decompose()) ** 2 * 1 / unit.sr
+        )  # [m2 sr]
+        signal3 = (
+            (self.sensor.efficiency(wavelength)).decompose()
+            * unit.electron
+            * self.get_transmittance()
+        )  # [electron]
+        signal4 = (
+            self.get_ratio_cropped_light_through_slit()
+            * (800 * 10 ** (-9))
+            * unit.meter
+        )  # [dimensionless * m]
 
-        noise = self.sensor.get_noise(signal)
+        signal = signal1 * signal2 * signal3 * signal4
+
+        print("signal", signal.decompose())
+
+        print("shot noise sqr", signal.decompose() * unit.electron)
+        print(
+            "dark noise sqr",
+            self.sensor.n_bin * (self.sensor.get_dark_noise() * unit.pix) ** 2,
+        )
+        print("quantization noise sqr", (self.sensor.get_quantization_noise()) ** 2)
+        print("read noise sqr", self.sensor.n_bin * (self.sensor.noise_read) ** 2)
+
+        noise = np.sqrt(
+            (signal * unit.electron)
+            + self.sensor.n_bin * (self.sensor.get_dark_noise() * unit.pix) ** 2
+            + self.sensor.get_quantization_noise() ** 2
+            + self.sensor.n_bin * self.sensor.noise_read**2
+        )
 
         snr = signal / noise
 
-        return snr
+        return snr.decompose()
 
     def get_FOV(self) -> np.ndarray[float, float]:
         """Get the field of view vector.
-        
+
         A vector that defines the angular extent that can be imaged by the payload in
         the along-track and the across-track directions.
 
         """
         assert self.slit is not None, "A slit component must be specified."
         assert self.foreoptic is not None, "A foreoptic component must be specified."
-        
+
         fov = 2 * np.arctan(self.slit.size / (2 * self.foreoptic.focal_length))
 
         return fov
@@ -124,7 +165,9 @@ def get_iFOV(self) -> np.ndarray[float, float]:
     def get_sensor_spatial_resolution(self, target_distance):
         """Get the sensor-limited spatial resolution."""
 
-        spatial_resolution = target_distance * self.sensor.pitch / self.foreoptic.focal_length
+        spatial_resolution = (
+            target_distance * self.sensor.pitch / self.foreoptic.focal_length
+        )
 
         return spatial_resolution.decompose()
 
@@ -148,7 +191,11 @@ def get_swath(
         return swath
 
     def get_optical_spatial_resolution(self, wavelength, target_distance, skew_angle=0):
-        """Get the optically-limited spatial resolution. Aka GRD (ground-resolved distance)"""
+        """Get the optically-limited spatial resolution.
+
+        Aka GRD (ground-resolved distance)
+
+        """
         assert self.foreoptic is not None, "A foreoptic component must be specified."
 
         optical_spatial_resolution = (
@@ -165,7 +212,8 @@ def get_spatial_resolution(self, wavelength, target_distance, skew_angle=0):
             return self.spatial_resolution
 
         sensor_spatial_resolution = self.get_sensor_spatial_resolution(
-            target_distance=target_distance)
+            target_distance=target_distance
+        )
 
         optical_spatial_resolution = self.get_optical_spatial_resolution(
             wavelength=wavelength,
@@ -240,7 +288,6 @@ def get_pointing_accuracy_constraint(
 
         return constraint_angle
 
-
 
 class FINCHEye(HyperspectralImager):
     def __init__(

tests/test_components/test_slits/test_Slit.py

@@ -1,81 +0,0 @@
-"""Tests for Slit component."""
-# stdlib
-import logging
-
-# external
-import numpy as np
-import pytest
-
-# project
-from payload_designer.components import slits
-
-LOG = logging.getLogger(__name__)
-
-
-def test_get_horizontal_field_of_view():
-    """Test ThinFocuser.get_horizontal_field_of_view()."""
-
-    # parameters
-    l_s = 1
-    f = 70
-
-    # component instantiation
-    slit = slits.Slit(l_s=l_s, f=f)
-
-    # evaluation
-    fov = slit.get_horizontal_field_of_view()
-    LOG.info(f"Horizontal field of view: {fov}")
-
-    assert fov == pytest.approx(360 * np.arctan(1 / 140) / np.pi)
-
-
-def test_get_vertical_field_of_view():
-    """Test ThinFocuser.get_vertical_field_of_view()."""
-
-    # parameters
-    w_s = 0.01
-    f = 70
-
-    # component instantiation
-    slit = slits.Slit(w_s=w_s, f=f)
-
-    # evaluation
-    fov = slit.get_vertical_field_of_view()
-    LOG.info(f"Vertical field of view: {fov}")
-
-    assert fov == pytest.approx(360 * np.arctan(1 / 14000) / np.pi)
-
-
-def test_get_image_width():
-    """Test ThinFocuser.get_image_width()."""
-    # w_i = np.power(np.multiply(np.power(self.m, 2), np.power(self.w_s, 2)) + np.power(self.w_o, 2), 0.2)
-    # parameters
-    m = 1
-    w_s = 0.01
-    w_o = 0.05
-
-    # component instantiation
-    slit = slits.Slit(m=m, w_s=w_s, w_o=w_o)
-
-    # evaluation
-    w = slit.get_image_width()
-    LOG.info(f"Image width: {w}")
-
-    assert w == pytest.approx(np.sqrt(13 / 5000))
-
-
-def test_get_slit_width():
-    """Test ThinFocuser.get_slit_width()."""
-    # w_s = np.divide(self.w_d, self.m)
-    # parameters
-    w_d = 5
-    m = 1
-
-    # component instantiation
-    slit = slits.Slit(w_d=w_d, m=m)
-
-    # evaluation
-    w = slit.get_slit_width()
-    LOG.info(f"Slit width: {w}")
-
-    assert w == 5