From 1ff2c410a15a85ffc1f77f5b84e72a6c752064d9 Mon Sep 17 00:00:00 2001 From: Matt Einhorn Date: Fri, 7 Jun 2024 00:46:14 -0400 Subject: [PATCH 1/3] Prevent zero underflow in volume. --- .../filters/volume/structure_splitting.py | 16 ++++++--- .../test_detection_structure_splitting.py | 35 +++++++++++++++++++ 2 files changed, 46 insertions(+), 5 deletions(-) diff --git a/cellfinder/core/detect/filters/volume/structure_splitting.py b/cellfinder/core/detect/filters/volume/structure_splitting.py index f0018df8..7480e8ec 100644 --- a/cellfinder/core/detect/filters/volume/structure_splitting.py +++ b/cellfinder/core/detect/filters/volume/structure_splitting.py @@ -74,6 +74,7 @@ def ball_filter_imgs( good_tiles_mask = np.ones((1, 1, volume.shape[2]), dtype=np.bool_) plane_width, plane_height = volume.shape[:2] + current_z = ball_z_size // 2 bf = BallFilter( plane_width, @@ -86,9 +87,7 @@ def ball_filter_imgs( threshold_value=threshold_value, soma_centre_value=soma_centre_value, ) - cell_detector = CellDetector( - plane_width, plane_height, start_z=ball_z_size // 2 - ) + cell_detector = CellDetector(plane_width, plane_height, start_z=current_z) # FIXME: hard coded type ball_filtered_volume = np.zeros(volume.shape, dtype=np.uint32) @@ -98,7 +97,11 @@ def ball_filter_imgs( if bf.ready: bf.walk() middle_plane = bf.get_middle_plane() - ball_filtered_volume[:, :, z] = middle_plane[:] + + # first valid middle plane is the current_z, not z + ball_filtered_volume[:, :, current_z] = middle_plane[:] + current_z += 1 + # DEBUG: TEST: transpose previous_plane = cell_detector.process( middle_plane.copy(), previous_plane @@ -134,7 +137,10 @@ def iterative_ball_filter( vol, cell_centres = ball_filter_imgs( vol, threshold_value, soma_centre_value ) - vol -= 1 + + # vol is unsigned, so can't let zeros underflow to max value + vol[:, :, :] = np.where(vol != 0, vol - 1, 0) + n_structures = len(cell_centres) ns.append(n_structures) centres.append(cell_centres) diff --git a/tests/core/test_integration/test_detection_structure_splitting.py b/tests/core/test_integration/test_detection_structure_splitting.py index 6fa83d07..d45027f6 100644 --- a/tests/core/test_integration/test_detection_structure_splitting.py +++ b/tests/core/test_integration/test_detection_structure_splitting.py @@ -8,9 +8,13 @@ import os +import numpy as np import pytest from brainglobe_utils.IO.image.load import read_with_dask +from cellfinder.core.detect.filters.volume.structure_splitting import ( + split_cells, +) from cellfinder.core.main import main data_dir = os.path.join( @@ -46,3 +50,34 @@ def test_structure_splitting(signal_array, background_array): voxel_sizes, n_free_cpus=0, ) + + +def test_underflow_issue_435(): + # two cells centered at (9, 10, 10), (19, 10, 10) with radius 5 + p1 = np.array([9, 10, 10]) + p2 = np.array([19, 10, 10]) + radius = 5 + + bright_voxels = np.zeros((20, 20, 30), dtype=np.bool_) + + pos = np.empty((20, 20, 30, 3)) + pos[:, :, :, 0] = np.arange(20).reshape((-1, 1, 1)) + pos[:, :, :, 1] = np.arange(20).reshape((1, -1, 1)) + pos[:, :, :, 2] = np.arange(30).reshape((1, 1, -1)) + + dist1 = pos - p1.reshape((1, 1, 1, 3)) + dist1 = np.sqrt(np.sum(np.square(dist1), axis=3)) + inside1 = dist1 <= radius + dist2 = pos - p2.reshape((1, 1, 1, 3)) + dist2 = np.sqrt(np.sum(np.square(dist2), axis=3)) + inside2 = dist2 <= radius + + bright_voxels[np.logical_or(inside1, inside2)] = True + bright_indices = np.argwhere(bright_voxels) + + centers = split_cells(bright_indices) + + # for some reason, same with pytorch, it's shifted by 1. Probably rounding + expected = {(10, 11, 11), (18, 11, 11)} + got = set(map(tuple, centers.tolist())) + assert expected == got From c21a9caab9d8db453851528f09cf2ca4926f6922 Mon Sep 17 00:00:00 2001 From: Matt Einhorn Date: Fri, 7 Jun 2024 12:14:44 -0400 Subject: [PATCH 2/3] Update tests/core/test_integration/test_detection_structure_splitting.py Co-authored-by: Alessandro Felder --- .../test_detection_structure_splitting.py | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/tests/core/test_integration/test_detection_structure_splitting.py b/tests/core/test_integration/test_detection_structure_splitting.py index d45027f6..be6f3b53 100644 --- a/tests/core/test_integration/test_detection_structure_splitting.py +++ b/tests/core/test_integration/test_detection_structure_splitting.py @@ -58,12 +58,12 @@ def test_underflow_issue_435(): p2 = np.array([19, 10, 10]) radius = 5 - bright_voxels = np.zeros((20, 20, 30), dtype=np.bool_) + bright_voxels = np.zeros((30, 20, 20), dtype=np.bool_) - pos = np.empty((20, 20, 30, 3)) - pos[:, :, :, 0] = np.arange(20).reshape((-1, 1, 1)) + pos = np.empty((30, 20, 20, 3)) + pos[:, :, :, 0] = np.arange(30).reshape((-1, 1, 1)) pos[:, :, :, 1] = np.arange(20).reshape((1, -1, 1)) - pos[:, :, :, 2] = np.arange(30).reshape((1, 1, -1)) + pos[:, :, :, 2] = np.arange(20).reshape((1, 1, -1)) dist1 = pos - p1.reshape((1, 1, 1, 3)) dist1 = np.sqrt(np.sum(np.square(dist1), axis=3)) From aa192d89fe536eefea1daddcf71928698b0e14fc Mon Sep 17 00:00:00 2001 From: Matt Einhorn Date: Fri, 7 Jun 2024 12:16:43 -0400 Subject: [PATCH 3/3] Now that point is fully inside, it'll all shifted upward by one. --- .../core/test_integration/test_detection_structure_splitting.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/tests/core/test_integration/test_detection_structure_splitting.py b/tests/core/test_integration/test_detection_structure_splitting.py index be6f3b53..31f04672 100644 --- a/tests/core/test_integration/test_detection_structure_splitting.py +++ b/tests/core/test_integration/test_detection_structure_splitting.py @@ -78,6 +78,6 @@ def test_underflow_issue_435(): centers = split_cells(bright_indices) # for some reason, same with pytorch, it's shifted by 1. Probably rounding - expected = {(10, 11, 11), (18, 11, 11)} + expected = {(10, 11, 11), (20, 11, 11)} got = set(map(tuple, centers.tolist())) assert expected == got