diff --git a/esbo_etc/classes/SpectralQty.py b/esbo_etc/classes/SpectralQty.py index 2102244..b173976 100644 --- a/esbo_etc/classes/SpectralQty.py +++ b/esbo_etc/classes/SpectralQty.py @@ -14,7 +14,7 @@ class SpectralQty: A class to hold and work with spectral quantities """ - def __init__(self, wl: u.Quantity, qty: u.Quantity, extrapolate: bool = False): + def __init__(self, wl: u.Quantity, qty: u.Quantity, fill_value: Union[bool, int, float] = 0): """ Initialize a new spectral quantity @@ -25,8 +25,9 @@ class SpectralQty: qty : Quantity The quantity values corresponding to the binned wavelengths. If the values are supplied without a unit, they are assumed to be dimensionless. - extrapolate : bool - Whether extrapolation should be allowed. If disabled, the spectrum will be truncated and a warning given. + fill_value : Union[bool, int, float] + How to treat missing values. True enables extrapolation, False disables extrapolation and the spectrum will + be truncated. If a numeric value is given, the missing values will be filled with this value. Returns ------- sqty : SpectralQty @@ -45,11 +46,11 @@ class SpectralQty: self.qty = qty * u.dimensionless_unscaled else: error("Lengths not matching") - self._extrapolate = extrapolate + self._fill_value = fill_value @classmethod def fromFile(cls, file: str, wl_unit_default: u.Quantity = None, qty_unit_default: u.Quantity = None, - extrapolate: bool = False) -> "SpectralQty": + fill_value: Union[bool, int, float] = 0) -> "SpectralQty": """ Initialize a new spectral quantity and read the values from a file @@ -64,8 +65,9 @@ class SpectralQty: Default unit to be used for the wavelength column if no units are provided by the file. qty_unit_default : Quantity Default unit to be used for the quantity column if no units are provided by the file. - extrapolate : bool - Whether extrapolation should be allowed. If disabled, the spectrum will be truncated and a warning given. + fill_value : Union[bool, int, float] + How to treat missing values. True enables extrapolation, False disables extrapolation and the spectrum will + be truncated. If a numeric value is given, the missing values will be filled with this value. Returns ------- sqty : SpectralQty @@ -89,7 +91,7 @@ class SpectralQty: elif wl_unit_default is not None and qty_unit_default is not None: data[data.colnames[0]].unit = wl_unit_default data[data.colnames[1]].unit = qty_unit_default - return cls(data[data.colnames[0]].quantity, data[data.colnames[1]].quantity, extrapolate=extrapolate) + return cls(data[data.colnames[0]].quantity, data[data.colnames[1]].quantity, fill_value=fill_value) def __eq__(self, other) -> bool: """ @@ -141,12 +143,11 @@ class SpectralQty: # Summand is of type SpectralQty else: if other.wl.unit.is_equivalent(self.wl.unit) and other.qty.unit.is_equivalent(self.qty.unit): - # Wavelengths are matching, just add the quantities if len(self.wl) == len(other.wl) and (self.wl == other.wl).all(): + # Wavelengths are matching, just add the quantities return SpectralQty(self.wl, self.qty + other.qty) - # Wavelengths are not matching, rebinning needed else: - # Rebin addend + # Wavelengths are not matching, rebinning needed other_rebinned = other.rebin(self.wl) if len(self.wl) == len(other_rebinned.wl) and (self.wl == other_rebinned.wl).all(): return SpectralQty(self.wl, self.qty + other_rebinned.qty) @@ -264,10 +265,15 @@ class SpectralQty: if not wl.unit.is_equivalent(self.wl.unit): error("Mismatching units for rebinning: " + wl.unit + ", " + self.wl.unit) - if not self._extrapolate: - if min(wl) < min(self.wl) or max(wl) > max(self.wl): - logging.warning("Extrapolation disabled, bandwidth will be reduced.") - # Remove new wavelengths where extrapolation would have been necessary - wl = [x.value for x in wl if min(self.wl) <= x <= max(self.wl)] * wl.unit - f = interp1d(self.wl, self.qty.value, fill_value="extrapolate") + if min(wl) < min(self.wl) or max(wl) > max(self.wl): + if isinstance(self._fill_value, bool): + if not self._fill_value: + logging.warning("Extrapolation disabled, bandwidth will be reduced.") + # Remove new wavelengths where extrapolation would have been necessary + wl = [x.value for x in wl if min(self.wl) <= x <= max(self.wl)] * wl.unit + f = interp1d(self.wl, self.qty.value, fill_value="extrapolate") + else: + f = interp1d(self.wl, self.qty.value, fill_value=self._fill_value, bounds_error=False) + else: + f = interp1d(self.wl, self.qty.value) return SpectralQty(wl, f(wl.to(self.wl.unit)) * self.qty.unit) diff --git a/tests/optical_component/test_Atmosphere.py b/tests/optical_component/test_Atmosphere.py index 91bf156..59ceec3 100644 --- a/tests/optical_component/test_Atmosphere.py +++ b/tests/optical_component/test_Atmosphere.py @@ -12,11 +12,12 @@ class TestAtmosphere(TestCase): def test_calcSignal(self): self.assertEqual(self.atmosphere.calcSignal(), - SpectralQty(np.arange(200, 208) << u.nm, [1.10e-15, 1.20e-15, 1.30e-15, 1.26e-15, 1.20e-15, - 1.12e-15, 1.02e-15, 0.9e-15] << u.W / ( - u.m ** 2 * u.nm))) + SpectralQty(np.arange(200, 210) << u.nm, + np.array([1.10e-15, 1.20e-15, 1.30e-15, 1.26e-15, 1.20e-15, 1.12e-15, 1.02e-15, + 0.9e-15, 0, 0]) << u.W / (u.m ** 2 * u.nm))) def test_calcNoise(self): self.assertEqual(self.atmosphere.calcNoise(), - SpectralQty(np.arange(200, 208) << u.nm, [1.1e-16, 1.2e-16, 1.3e-16, 1.4e-16, 1.5e-16, 1.6e-16, - 1.7e-16, 1.8e-16] << u.W / (u.m ** 2 * u.nm * u.sr))) + SpectralQty(np.arange(200, 208) << u.nm, + np.array([1.1e-16, 1.2e-16, 1.3e-16, 1.4e-16, 1.5e-16, 1.6e-16, 1.7e-16, + 1.8e-16]) << u.W / (u.m ** 2 * u.nm * u.sr))) diff --git a/tests/optical_component/test_BeamSplitter.py b/tests/optical_component/test_BeamSplitter.py index 8b6061c..f6777a2 100644 --- a/tests/optical_component/test_BeamSplitter.py +++ b/tests/optical_component/test_BeamSplitter.py @@ -16,7 +16,9 @@ class TestBeamSplitter(TestCase): def test___init__(self): self.assertEqual(self.splitter.calcNoise(), - SpectralQty(self.wl, [4.31413931e-96, 1.37122214e-95, 4.30844544e-95, 1.33846280e-94] << u.W / - (u.m ** 2 * u.nm * u.sr))) + SpectralQty(self.wl, np.array([4.31413931e-96, 1.37122214e-95, 4.30844544e-95, + 1.33846280e-94]) << u.W / (u.m ** 2 * u.nm * u.sr))) self.assertEqual(self.splitter.calcSignal(), - SpectralQty(self.wl, [1.20e-15, 1.30e-15, 1.40e-15, 1.35e-15] << u.W / (u.m ** 2 * u.nm))) + SpectralQty(np.arange(200, 210, 1) << u.nm, + np.array([0, 1.20e-15, 1.30e-15, 1.40e-15, 1.35e-15, 0, 0, 0, 0, 0]) << u.W / + (u.m ** 2 * u.nm))) diff --git a/tests/optical_component/test_Lens.py b/tests/optical_component/test_Lens.py index edd5595..11e8553 100644 --- a/tests/optical_component/test_Lens.py +++ b/tests/optical_component/test_Lens.py @@ -15,7 +15,9 @@ class TestLens(TestCase): def test___init__(self): self.assertEqual(self.lens.calcNoise(), - SpectralQty(self.wl, [4.31413931e-96, 1.37122214e-95, 4.30844544e-95, 1.33846280e-94] << u.W / - (u.m ** 2 * u.nm * u.sr))) + SpectralQty(self.wl, np.array([4.31413931e-96, 1.37122214e-95, 4.30844544e-95, + 1.33846280e-94]) << u.W / (u.m ** 2 * u.nm * u.sr))) self.assertEqual(self.lens.calcSignal(), - SpectralQty(self.wl, [1.20e-15, 1.30e-15, 1.40e-15, 1.35e-15] << u.W / (u.m ** 2 * u.nm))) + SpectralQty(np.arange(200, 210, 1) << u.nm, + np.array([0, 1.20e-15, 1.30e-15, 1.40e-15, 1.35e-15, 0, 0, 0, 0, 0]) << u.W / + (u.m ** 2 * u.nm))) diff --git a/tests/optical_component/test_Mirror.py b/tests/optical_component/test_Mirror.py index 5669673..ce12a93 100644 --- a/tests/optical_component/test_Mirror.py +++ b/tests/optical_component/test_Mirror.py @@ -15,7 +15,9 @@ class TestMirror(TestCase): def test___init__(self): self.assertEqual(self.mirror.calcNoise(), - SpectralQty(self.wl, [4.31413931e-96, 1.37122214e-95, 4.30844544e-95, 1.33846280e-94] << u.W / - (u.m ** 2 * u.nm * u.sr))) + SpectralQty(self.wl, np.array([4.31413931e-96, 1.37122214e-95, 4.30844544e-95, + 1.33846280e-94]) << u.W / (u.m ** 2 * u.nm * u.sr))) self.assertEqual(self.mirror.calcSignal(), - SpectralQty(self.wl, [1.20e-15, 1.30e-15, 1.40e-15, 1.35e-15] << u.W / (u.m ** 2 * u.nm))) + SpectralQty(np.arange(200, 210, 1) << u.nm, np.array([0, 1.20e-15, 1.30e-15, 1.40e-15, + 1.35e-15, 0, 0, 0, 0, 0]) << u.W / + (u.m ** 2 * u.nm))) diff --git a/tests/test_SpectralQty.py b/tests/test_SpectralQty.py index ab3e4b5..189db32 100644 --- a/tests/test_SpectralQty.py +++ b/tests/test_SpectralQty.py @@ -28,20 +28,20 @@ class TestSpectralQty(TestCase): np.arange(2.2, 3.0, 2e-1) << u.W / (u.m ** 2 * u.nm))) # SpectralQty self.assertEqual(self.sqty * SpectralQty(self.wl, np.arange(1, 5, 1) << u.m), - SpectralQty(self.wl, [1.1, 2.4, 3.9, 5.6] << u.W / (u.m * u.nm))) + SpectralQty(self.wl, np.array([1.1, 2.4, 3.9, 5.6]) << u.W / (u.m * u.nm))) self.assertEqual(SpectralQty(self.wl, np.arange(1, 5, 1) << u.m) * self.sqty, - SpectralQty(self.wl, [1.1, 2.4, 3.9, 5.6] << u.W / (u.m * u.nm))) + SpectralQty(self.wl, np.array([1.1, 2.4, 3.9, 5.6]) << u.W / (u.m * u.nm))) # rebin without extrapolation and without reduction self.assertEqual( self.sqty * SpectralQty(np.arange(199.5, 204.5, 1) << u.nm, np.arange(1, 6, 1) << u.m), SpectralQty(self.wl, [1.65, 3.0, 4.55, 6.3] * u.W / (u.m * u.nm))) # rebin without extrapolation and with reduction self.assertEqual( - self.sqty * SpectralQty(np.arange(200.5, 204.5, 1) << u.nm, np.arange(1, 5, 1) << u.m), - SpectralQty(range(201, 204) << u.nm, [1.8, 3.25, 4.9] << u.W / (u.m * u.nm))) + self.sqty * SpectralQty(np.arange(200.5, 204.5, 1) << u.nm, np.arange(1, 5, 1) << u.m, fill_value=False), + SpectralQty(np.arange(201, 204) << u.nm, np.array([1.8, 3.25, 4.9]) << u.W / (u.m * u.nm))) # lambda - sqty_2 = lambda wl: 0.7 * u.dimensionless_unscaled - self.assertEqual(self.sqty * sqty_2, SpectralQty(self.wl, self.qty * 0.7)) + self.assertEqual(self.sqty * (lambda wl: 0.7 * u.dimensionless_unscaled), + SpectralQty(self.wl, self.qty * 0.7)) def test___sub__(self): # Quantity @@ -58,12 +58,12 @@ class TestSpectralQty(TestCase): SpectralQty(self.wl, [-0.4, -1.3, -2.2, -3.1] * u.W / (u.m ** 2 * u.nm))) # rebin without extrapolation and with reduction self.assertEqual( - self.sqty - SpectralQty(np.arange(200.5, 204.5, 1) << u.nm, np.arange(1, 5, 1) << u.W / (u.m ** 2 * u.nm)), - SpectralQty(range(201, 204) << u.nm, [-0.3, -1.2, -2.1] << u.W / (u.m ** 2 * u.nm))) + self.sqty - SpectralQty(np.arange(200.5, 204.5, 1) << u.nm, np.arange(1, 5, 1) << u.W / (u.m ** 2 * u.nm), + fill_value=False), + SpectralQty(np.arange(201, 204) << u.nm, np.array([-0.3, -1.2, -2.1]) << u.W / (u.m ** 2 * u.nm))) # lambda - sqty_2 = lambda wl: 1 * u.W / (u.m ** 2 * u.nm ** 2) * wl - self.assertEqual(self.sqty - sqty_2, - SpectralQty(self.wl, [-198.9, -199.8, -200.7, -201.6] << u.W / (u.m ** 2 * u.nm))) + self.assertEqual(self.sqty - (lambda wl: 1 * u.W / (u.m ** 2 * u.nm ** 2) * wl), + SpectralQty(self.wl, np.array([-198.9, -199.8, -200.7, -201.6]) << u.W / (u.m ** 2 * u.nm))) def test___add__(self): # Quantity @@ -80,33 +80,38 @@ class TestSpectralQty(TestCase): SpectralQty(self.wl, [2.6, 3.7, 4.8, 5.9] * u.W / (u.m ** 2 * u.nm))) # rebin without extrapolation and with reduction self.assertEqual( - self.sqty + SpectralQty(np.arange(200.5, 204.5, 1) << u.nm, np.arange(1, 5, 1) << u.W / (u.m ** 2 * u.nm)), - SpectralQty(range(201, 204) << u.nm, [2.7, 3.8, 4.9] << u.W / (u.m ** 2 * u.nm))) + self.sqty + SpectralQty(np.arange(200.5, 204.5, 1) << u.nm, np.arange(1, 5, 1) << u.W / (u.m ** 2 * u.nm), + fill_value=False), + SpectralQty(np.arange(201, 204) << u.nm, np.array([2.7, 3.8, 4.9]) << u.W / (u.m ** 2 * u.nm))) # lambda - sqty_2 = lambda wl: 1 * u.W / (u.m ** 2 * u.nm ** 2) * wl - self.assertEqual(self.sqty + sqty_2, - SpectralQty(self.wl, [201.1, 202.2, 203.3, 204.4] << u.W / (u.m ** 2 * u.nm))) + self.assertEqual(self.sqty + (lambda wl: 1 * u.W / (u.m ** 2 * u.nm ** 2) * wl), + SpectralQty(self.wl, np.array([201.1, 202.2, 203.3, 204.4]) << u.W / (u.m ** 2 * u.nm))) def test_rebinning(self): # Test interpolation wl_new = np.arange(200.5, 210.5, 1) << u.nm - sqty_res = SpectralQty(wl_new[:3], [1.15, 1.25, 1.35] << u.W / (u.m ** 2 * u.nm)) - sqty_rebin = self.sqty.rebin(wl_new) + sqty_res = SpectralQty(wl_new[:3], np.array([1.15, 1.25, 1.35]) << u.W / (u.m ** 2 * u.nm)) + sqty_rebin = SpectralQty(self.wl, self.qty, fill_value=False).rebin(wl_new) self.assertEqual(sqty_rebin, sqty_res) # Test extrapolation wl_new = np.arange(200.5, 210.5, 1) << u.nm - sqty_res = SpectralQty(wl_new, [1.15, 1.25, 1.35, 1.45, 1.55, 1.65, 1.75, 1.85, - 1.95, 2.05] << u.W / (u.m ** 2 * u.nm)) - sqty_extrapol = SpectralQty(self.wl, self.qty, extrapolate=True) - sqty_rebin = sqty_extrapol.rebin(wl_new) + sqty_res = SpectralQty(wl_new, np.array([1.15, 1.25, 1.35, 1.45, 1.55, 1.65, 1.75, 1.85, + 1.95, 2.05]) << u.W / (u.m ** 2 * u.nm)) + sqty_rebin = SpectralQty(self.wl, self.qty, fill_value=True).rebin(wl_new) + self.assertEqual(sqty_rebin, sqty_res) + + # Test fill value + wl_new = np.arange(200.5, 210.5, 1) << u.nm + sqty_res = SpectralQty(wl_new, np.array([1.15, 1.25, 1.35, 0, 0, 0, 0, 0, 0, 0]) << u.W / (u.m ** 2 * u.nm)) + sqty_rebin = SpectralQty(self.wl, self.qty, fill_value=0).rebin(wl_new) self.assertEqual(sqty_rebin, sqty_res) # Test binning self.setUp() wl_new = np.arange(200.5, 210, 2) << u.nm - sqty_res = SpectralQty(wl_new[:2], [1.15, 1.35] << u.W / (u.m ** 2 * u.nm)) - sqty_rebin = self.sqty.rebin(wl_new) + sqty_res = SpectralQty(wl_new[:2], np.array([1.15, 1.35]) << u.W / (u.m ** 2 * u.nm)) + sqty_rebin = SpectralQty(self.wl, self.qty, fill_value=False).rebin(wl_new) self.assertEqual(sqty_rebin, sqty_res) def test_fromFile(self):