feat: init

test/model/analysis/test_boundary_detect.py

@@ -0,0 +1,112 @@
+# Copyright (c) 2021-2025  The University of Texas Southwestern Medical Center.
+# All rights reserved.
+
+# Redistribution and use in source and binary forms, with or without
+# modification, are permitted for academic and research use only (subject to the
+# limitations in the disclaimer below) provided that the following conditions are met:
+
+#      * Redistributions of source code must retain the above copyright notice,
+#      this list of conditions and the following disclaimer.
+
+#      * Redistributions in binary form must reproduce the above copyright
+#      notice, this list of conditions and the following disclaimer in the
+#      documentation and/or other materials provided with the distribution.
+
+#      * Neither the name of the copyright holders nor the names of its
+#      contributors may be used to endorse or promote products derived from this
+#      software without specific prior written permission.
+
+# NO EXPRESS OR IMPLIED LICENSES TO ANY PARTY'S PATENT RIGHTS ARE GRANTED BY
+# THIS LICENSE. THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND
+# CONTRIBUTORS "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
+# LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A
+# PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR
+# CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
+# EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
+# PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR
+# BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER
+# IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
+# ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
+# POSSIBILITY OF SUCH DAMAGE.
+
+import math
+
+import numpy as np
+
+
+def im_circ(r=1, N=128):
+    X, Y = np.meshgrid(range(-N // 2, N // 2), range(-N // 2, N // 2))
+    return (X * X + Y * Y) < r * r
+
+
+def test_has_tissue():
+    from navigate.model.analysis.boundary_detect import has_tissue
+
+    for _ in range(100):
+        N = 2 ** np.random.randint(5, 9)
+        r = np.random.randint(math.ceil(0.2 * N), int(0.4 * N))
+        ds = np.random.randint(1, 6)
+
+        print(N, r, ds)
+
+        im = im_circ(r, N) * 1001
+
+        mu, sig = 100 * np.random.rand() + 1, 10 * np.random.rand() + 1
+        print(mu, sig)
+        offsets = [None, np.ones((N, N)) * mu]
+        variances = [None, np.ones((N, N)) * sig]
+
+        for off, var in zip(offsets, variances):
+            assert has_tissue(im, 0, 0, N, off, var) and not has_tissue(
+                im, 0, 0, N // 2 - r, off, var
+            )
+
+
+def test_find_tissue_boundary_2d():
+    from skimage.transform import downscale_local_mean
+
+    from navigate.model.analysis.boundary_detect import find_tissue_boundary_2d
+
+    for _ in range(100):
+        N = 2 ** np.random.randint(5, 9)
+        r = np.random.randint(1, int(0.4 * N))
+        ds = np.random.randint(1, 6)
+
+        print(N, r, ds)
+
+        im = im_circ(r, N)
+        b = find_tissue_boundary_2d(im, ds)
+        b = np.vstack([x for x in b if x is not None])
+
+        idx_x, idx_y = np.where(downscale_local_mean(im, (ds, ds)))
+        iixy = (np.unique(idx_x)[:, None] == idx_x[None, :]) * idx_y
+        low, high = idx_y[np.argmax(iixy != 0, 1)], np.max(iixy, 1)
+
+        np.testing.assert_equal(b, np.vstack([low, high]).T)
+
+
+def test_binary_detect():
+    from navigate.model.analysis.boundary_detect import (
+        find_tissue_boundary_2d,
+        binary_detect,
+    )
+
+    for _ in range(100):
+        N = 2 ** np.random.randint(5, 9)
+        r = np.random.randint(1, int(0.4 * N))
+        ds = np.random.randint(1, 6)
+
+        print(N, r, ds)
+
+        im = im_circ(r, N)
+        b = find_tissue_boundary_2d(im, ds)
+
+        assert binary_detect(im * 1001, b, ds) == b
+
+
+def test_map_boundary():
+    from navigate.model.analysis.boundary_detect import map_boundary
+
+    assert map_boundary([[1, 2]]) == [(0, 1), (0, 2)]
+    assert map_boundary([None, [1, 2]]) == [(1, 1), (1, 2)]
+    assert map_boundary([None, [1, 2], None]) == [(1, 1), (1, 2)]

test/model/analysis/test_camera.py

@@ -0,0 +1,56 @@
+import numpy as np
+import pytest
+
+
+@pytest.mark.skip("volatile")
+def test_compute_scmos_offset_and_variance_map():
+    from navigate.model.analysis.camera import compute_scmos_offset_and_variance_map
+
+    mu, sig = 100 * np.random.rand() + 1, 100 * np.random.rand() + 1
+    im = sig * np.random.randn(256, 256, 256) + mu
+
+    offset, variance = compute_scmos_offset_and_variance_map(im)
+
+    print(mu, sig)
+    # TODO: 1 is a bit high?
+    np.testing.assert_allclose(offset, mu, rtol=1)
+    np.testing.assert_allclose(variance, sig * sig, rtol=1)
+
+
+@pytest.mark.parametrize("local", [True, False])
+def test_compute_flatfield_map(local):
+    from navigate.model.analysis.camera import compute_flatfield_map
+
+    image = np.ones((256, 256))
+    offset = np.zeros((256, 256))
+    ffmap = compute_flatfield_map(image, offset, local)
+
+    np.testing.assert_allclose(ffmap, 0.5)
+
+
+def test_compute_noise_sigma():
+    from navigate.model.analysis.camera import compute_noise_sigma
+
+    Fn = np.random.rand()
+    qe = np.random.rand()
+    S = np.random.rand(256, 256)
+    Ib = np.random.rand()
+    Nr = np.random.rand()
+    M = np.random.rand()
+    sigma = compute_noise_sigma(Fn=Fn, qe=qe, S=S, Ib=Ib, Nr=Nr, M=M)
+    sigma_true = np.sqrt(Fn * Fn * qe * (S + Ib) + (Nr / M) ** 2)
+
+    np.testing.assert_allclose(sigma, sigma_true)
+
+
+def test_compute_signal_to_noise():
+    from navigate.model.analysis.camera import compute_signal_to_noise
+
+    A = np.random.rand() * 100 + 10
+    image = A * np.ones((256, 256))
+    offset = np.zeros((256, 256))
+    variance = 3 * A * A * np.ones((256, 256))
+
+    snr = compute_signal_to_noise(image, offset, variance)
+
+    np.testing.assert_allclose(snr, 0.5, rtol=0.2)