Add unit tests for extraction edge cases

jwst/extract_1d/extract1d.py

@@ -683,7 +683,8 @@ def _extract_colpix(image_data, x, j, limits):
             continue
         maxval = min(ns, int(math.floor(interval[1] + 0.5)))
         minval = max(0, int(math.floor(interval[0] + 0.5)))
-        npts += maxval - minval + 1
+        if maxval - minval + 1 > 0:
+            npts += maxval - minval + 1
     if npts == 0:
         return [], [], []
 

jwst/extract_1d/tests/test_extract_src_flux.py

@@ -4,12 +4,13 @@
 import math
 
 import numpy as np
+import pytest
 
 from jwst.extract_1d import extract1d
 
 
-def test_extract_src_flux():
-
+@pytest.fixture
+def inputs_constant():
     shape = (9, 5)
     image = np.arange(shape[0] * shape[1], dtype=np.float32).reshape(shape)
     var_poisson = image.copy()
@@ -24,6 +25,14 @@ def test_extract_src_flux():
     upper = np.zeros(shape[1], dtype=np.float64) + 7.   # middle of pixel 7
     srclim = [[lower, upper]]
 
+    return (image, var_poisson, var_rnoise, var_rflat,
+            x, j, lam, srclim, weights, bkgmodels)
+
+
+def test_extract_src_flux(inputs_constant):
+    (image, var_poisson, var_rnoise, var_rflat,
+     x, j, lam, srclim, weights, bkgmodels) = inputs_constant
+
     (total_flux, f_var_poisson, f_var_rnoise, f_var_flat,
      bkg_flux, b_var_poisson, b_var_rnoise, b_var_flat,
      tarea, twht) = extract1d._extract_src_flux(
@@ -61,3 +70,48 @@ def test_extract_src_flux():
     assert np.isnan(total_flux)
 
     assert tarea == 0.
+
+
+@pytest.mark.parametrize('test_type', ['all_empty', 'all_equal'])
+def test_extract_src_flux_empty_interval(inputs_constant, test_type):
+    (image, var_poisson, var_rnoise, var_rflat,
+     x, j, lam, srclim, weights, bkgmodels) = inputs_constant
+
+    if test_type == 'all_empty':
+        # no limits provided
+        srclim = []
+    else:
+        # empty extraction range: upper equals lower
+        srclim[0][1] = srclim[0][0].copy()
+
+    (total_flux, f_var_poisson, f_var_rnoise, f_var_flat,
+     bkg_flux, b_var_poisson, b_var_rnoise, b_var_flat,
+     tarea, twht) = extract1d._extract_src_flux(
+        image, var_poisson, var_rnoise, var_rflat,
+        x, j, lam, srclim, weights, bkgmodels)
+
+    # empty interval, so no flux returned
+    assert np.isnan(total_flux)
+    assert bkg_flux == 0.
+    assert tarea == 0.
+
+
+@pytest.mark.parametrize('offset', [-100, 100])
+def test_extract_src_flux_interval_out_of_range(inputs_constant, offset):
+    (image, var_poisson, var_rnoise, var_rflat,
+     x, j, lam, srclim, weights, bkgmodels) = inputs_constant
+
+    # extraction limits out of range
+    srclim[0][0] += offset
+    srclim[0][1] += offset
+
+    (total_flux, f_var_poisson, f_var_rnoise, f_var_flat,
+     bkg_flux, b_var_poisson, b_var_rnoise, b_var_flat,
+     tarea, twht) = extract1d._extract_src_flux(
+        image, var_poisson, var_rnoise, var_rflat,
+        x, j, lam, srclim, weights, bkgmodels)
+
+    # empty interval, so no flux returned
+    assert np.isnan(total_flux)
+    assert bkg_flux == 0.
+    assert tarea == 0.

jwst/extract_1d/tests/test_fit_background_model.py

@@ -2,6 +2,8 @@
 Test for extract_1d._fit_background_model
 """
 import math
+from copy import deepcopy
+
 import numpy as np
 import pytest
 
@@ -138,3 +140,89 @@ def test_handles_one_value(inputs_constant):
     assert math.isclose(b_var_flat_model(8.), 0.0, rel_tol=1.e-8, abs_tol=1.e-8)
 
     assert npts == 0
+
+
+@pytest.mark.parametrize('test_type', ['all_empty', 'all_equal'])
+def test_handles_empty_interval(inputs_constant, test_type):
+    image, var_poisson, var_rnoise, var_rflat, x, j, bkglim, bkg_fit, bkg_order = inputs_constant
+
+    if test_type == 'all_empty':
+        # no limits provided
+        bkglim = []
+    else:
+        # empty extraction range: upper equals lower
+        bkglim[0][1] = bkglim[0][0].copy()
+
+    # No data available: background model is 0
+    (bkg_model, b_var_poisson_model, b_var_rnoise_model, b_var_flat_model, npts) = \
+        extract1d._fit_background_model(image, var_poisson, var_rnoise, var_rflat,
+                                        x, j, bkglim, bkg_fit, bkg_order)
+
+    assert math.isclose(bkg_model(0.), 0.0, rel_tol=1.e-8, abs_tol=1.e-8)
+    assert math.isclose(bkg_model(8.), 0.0, rel_tol=1.e-8, abs_tol=1.e-8)
+
+    assert math.isclose(b_var_poisson_model(0.), 0.0, rel_tol=1.e-8, abs_tol=1.e-8)
+    assert math.isclose(b_var_poisson_model(8.), 0.0, rel_tol=1.e-8, abs_tol=1.e-8)
+
+    assert math.isclose(b_var_rnoise_model(0.), 0.0, rel_tol=1.e-8, abs_tol=1.e-8)
+    assert math.isclose(b_var_rnoise_model(8.), 0.0, rel_tol=1.e-8, abs_tol=1.e-8)
+
+    assert math.isclose(b_var_flat_model(0.), 0.0, rel_tol=1.e-8, abs_tol=1.e-8)
+    assert math.isclose(b_var_flat_model(8.), 0.0, rel_tol=1.e-8, abs_tol=1.e-8)
+
+    assert npts == 0
+
+
+@pytest.mark.parametrize('offset', [-100, 100])
+def test_handles_interval_out_of_range(inputs_constant, offset):
+    image, var_poisson, var_rnoise, var_rflat, x, j, bkglim, bkg_fit, bkg_order = inputs_constant
+    bkglim[0][0] += offset
+    bkglim[0][1] += offset
+
+    # No data available: background model is 0
+    (bkg_model, b_var_poisson_model, b_var_rnoise_model, b_var_flat_model, npts) = \
+        extract1d._fit_background_model(image, var_poisson, var_rnoise, var_rflat,
+                                        x, j, bkglim, bkg_fit, bkg_order)
+
+    assert math.isclose(bkg_model(0.), 0.0, rel_tol=1.e-8, abs_tol=1.e-8)
+    assert math.isclose(bkg_model(8.), 0.0, rel_tol=1.e-8, abs_tol=1.e-8)
+
+    assert math.isclose(b_var_poisson_model(0.), 0.0, rel_tol=1.e-8, abs_tol=1.e-8)
+    assert math.isclose(b_var_poisson_model(8.), 0.0, rel_tol=1.e-8, abs_tol=1.e-8)
+
+    assert math.isclose(b_var_rnoise_model(0.), 0.0, rel_tol=1.e-8, abs_tol=1.e-8)
+    assert math.isclose(b_var_rnoise_model(8.), 0.0, rel_tol=1.e-8, abs_tol=1.e-8)
+
+    assert math.isclose(b_var_flat_model(0.), 0.0, rel_tol=1.e-8, abs_tol=1.e-8)
+    assert math.isclose(b_var_flat_model(8.), 0.0, rel_tol=1.e-8, abs_tol=1.e-8)
+
+    assert npts == 0
+
+
+def test_handles_one_empty_interval(inputs_constant):
+    image, var_poisson, var_rnoise, var_rflat, x, j, bkglim, bkg_fit, bkg_order = inputs_constant
+
+    # add an extra interval that is empty
+    bkglim.append(deepcopy(bkglim[0]))
+    bkglim[1][0] += 2.0
+    bkglim[1][1] = bkglim[1][0].copy()
+    print(bkglim)
+
+    # should ignore the second interval and return a valid answer for the first
+    (bkg_model, b_var_poisson_model, b_var_rnoise_model, b_var_flat_model, npts) = \
+        extract1d._fit_background_model(image, var_poisson, var_rnoise, var_rflat,
+                                        x, j, bkglim, bkg_fit, bkg_order)
+
+    assert math.isclose(bkg_model(0.), 22.0, rel_tol=1.e-8, abs_tol=1.e-8)
+    assert math.isclose(bkg_model(8.), 22.0, rel_tol=1.e-8, abs_tol=1.e-8)
+
+    assert math.isclose(b_var_poisson_model(0.), 22.0, rel_tol=1.e-8, abs_tol=1.e-8)
+    assert math.isclose(b_var_poisson_model(8.), 22.0, rel_tol=1.e-8, abs_tol=1.e-8)
+
+    assert math.isclose(b_var_rnoise_model(0.), 22.0, rel_tol=1.e-8, abs_tol=1.e-8)
+    assert math.isclose(b_var_rnoise_model(8.), 22.0, rel_tol=1.e-8, abs_tol=1.e-8)
+
+    assert math.isclose(b_var_flat_model(0.), 22.0, rel_tol=1.e-8, abs_tol=1.e-8)
+    assert math.isclose(b_var_flat_model(8.), 22.0, rel_tol=1.e-8, abs_tol=1.e-8)
+
+    assert npts == 2