perf: num_transitions rather than num_foreground + return_N + automatic out_dtype + remove sparse (#55)

* feat(voxel_connectivity_graph): extracts directed voxel connectivity graph

Each voxel of the output is a bitfield that represents the allowed
transit directions.

* docs: added description of bitfields to python help

* refactor: move region_graph tests to bottom

* test: add simple tests for 4 and 8 graph connectivity

* test: skip big memory test for 32-bit windows

* fix: voxel size + fix compiler warnings

* test: adds test for 3d voxel connectivity graph

* fix: the masks dont make sense in row-major format so just drop it for now

* docs: update README with voxel_connectivity_graph

* docs: describe how to use the voxel graph in C++

* docs: correction to data width of black cube

* perf(sparse): use num transitions rather than num foreground voxels

* fix: if first voxel is foreground, add one to transitions

* BREAKING: radical revision that eliminates max_labels, sparse, and out_dtype

Uses num_transitions to select appropriate out_dtype and obviates sparse.

* fix: compiler warnings around printf w/ templates

* test: fix and add tests, drop sparse, out_dtype, add zeroth_pass

* fix: needed to account for boundaries in transition estimation

* docs: describe automatic out dtype assignment

* chore: add fastremap to setup.py

* fix: actually cc3d doesn't have a limitation on the input labels...

* feat: add return_N argument

* docs: show how to use the new library features

* test: remove py3.9 from appveyor

* test: add py3.8 to travis

* fix: wrong invocation for C++

Author: william-silversmith

.travis.yml

@@ -3,7 +3,8 @@ python:
 - '2.7'
 - '3.5'
 - '3.6'
-- '3.7-dev'
+- '3.7'
+- '3.8'
 before_install:
 - PYTHON_MAJOR_VERSION=`echo $TRAVIS_PYTHON_VERSION | head -c 1`
 - if [[ $PYTHON_MAJOR_VERSION == 3 ]]; then sudo apt-get install python3-pip; fi

README.md

@@ -41,8 +41,6 @@ python setup.py develop
 
 ## Python Use
 
-Important limitation: Only label values less than or equal to the size of the image in voxels (pixels) are supported currently. If you want to use larger values, consider using [fastremap.renumber](https://github.com/seung-lab/fastremap/).
-
 ```python
 import cc3d
 import numpy as np
@@ -53,20 +51,20 @@ labels_out = cc3d.connected_components(labels_in) # 26-connected
 connectivity = 6 # only 4,8 (2D) and 26, 18, and 6 (3D) are allowed
 labels_out = cc3d.connected_components(labels_in, connectivity=connectivity)
 
-# You can adjust the bit width of the output to accomodate
-# different expected image statistics with memory usage tradeoffs.
-# uint16, uint32 (default), and uint64 are supported. Note that
-# uint16 will limit you to 2^16 provisional labels.
-labels_out = cc3d.connected_components(labels_in, out_dtype=np.uint16)
+# We use a third pass to estimate the number of provisional labels prior 
+# to running CCL. This is usually faster and lower memory for large datasets. 
+# If you have a special situation, you can switch this feature off. Note that
+# np.bool datatypes can use fewer labels without this estimation scan.
+labels_out = cc3d.connected_components(labels_in, zeroth_pass=False)
 
-# If you know that the number of foreground voxels is relatively
-# low, you can save memory and sometimes time by enabling sparse mode
-# which will shrink the memory allocation of the Union-Find datastructure
-# to exactly match the number of foreground voxels if not lower.
-labels_out = cc3d.connected_components(labels_in, sparse=True)
+# You can extract the number of labels (which is also the maximum 
+# label value) like so:
+labels_out, N = cc3d.connected_components(labels_in, return_N=True) # free
+# -- OR -- 
+labels_out = cc3d.connected_components(labels_in) 
+N = np.max(labels_out) # costs a full read
 
 # You can extract individual components like so:
-N = np.max(labels_out)
 for segid in range(1, N+1):
   extracted_image = labels_out * (labels_out == segid)
   process(extracted_image)
@@ -84,7 +82,7 @@ graph = cc3d.voxel_connectivity_graph(labels, connectivity=connectivity)
 
 ```
 
-If you know approximately how many labels you are going to generate, you can save some memory by specifying a number a safety factor above that range. The max label ID in your input labels must be less than `max_labels`.
+If you know approximately how many labels you are going to generate, you can save some memory by specifying a number a safety factor above that range. The max label ID in your input labels must be less than `max_labels`.  *This option is only recommended for expert users.*
 
 ```python
 labels_out = connected_components(labels_in, max_labels=20000)

appveyor.yml

@@ -4,11 +4,11 @@ environment:
     - PYTHON: "C:\\Python36"
     - PYTHON: "C:\\Python37"
     - PYTHON: "C:\\Python38"
-    - PYTHON: "C:\\Python39"
+    # - PYTHON: "C:\\Python39"
     - PYTHON: "C:\\Python36-x64"
     - PYTHON: "C:\\Python37-x64"
     - PYTHON: "C:\\Python38-x64"
-    - PYTHON: "C:\\Python39-x64"
+    # - PYTHON: "C:\\Python39-x64"
 
 install:
   # We need wheel installed to build wheels

automated_test.py

@@ -21,8 +21,7 @@ def gt_c2f(gt):
 
 @pytest.mark.parametrize("connectivity", (4, 6, 8, 18, 26))
 @pytest.mark.parametrize("dtype", TEST_TYPES)
-@pytest.mark.parametrize("out_dtype", OUT_TYPES)
-def test_2d_square(out_dtype, dtype, connectivity):
+def test_2d_square(dtype, connectivity):
   def test(order, ground_truth):
     input_labels = np.zeros( (16,16), dtype=dtype, order=order )
     input_labels[:8,:8] = 8
@@ -31,7 +30,7 @@ def test(order, ground_truth):
     input_labels[8:,8:] = 11
 
     output_labels = cc3d.connected_components(
-      input_labels, out_dtype=out_dtype, connectivity=connectivity
+      input_labels, connectivity=connectivity
     ).astype(dtype)
 
     print(output_labels)
@@ -62,8 +61,7 @@ def test(order, ground_truth):
 
 @pytest.mark.parametrize("connectivity", (4, 6, 8, 18, 26))
 @pytest.mark.parametrize("dtype", TEST_TYPES)
-@pytest.mark.parametrize("out_dtype", OUT_TYPES)
-def test_2d_rectangle(out_dtype, dtype, connectivity):
+def test_2d_rectangle(dtype, connectivity):
   def test(order, ground_truth):
     input_labels = np.zeros( (16,13), dtype=dtype, order=order )
     input_labels[:8,:8] = 8
@@ -72,7 +70,7 @@ def test(order, ground_truth):
     input_labels[8:,8:] = 11
 
     output_labels = cc3d.connected_components(
-      input_labels, out_dtype=out_dtype, connectivity=connectivity
+      input_labels, connectivity=connectivity
     ).astype(dtype)
     print(output_labels.shape)
     output_labels = output_labels[:,:]
@@ -412,9 +410,7 @@ def test_max_labels_nonsensical():
   assert np.all(real_labels == zero_labels)
   assert np.all(real_labels == negative_labels)
 
-@pytest.mark.parametrize("out_dtype", OUT_TYPES)
-@pytest.mark.parametrize("sparse", (True, False))
-def test_compare_scipy_26(out_dtype, sparse):
+def test_compare_scipy_26():
   import scipy.ndimage.measurements
 
   sx, sy, sz = 128, 128, 128
@@ -426,16 +422,16 @@ def test_compare_scipy_26(out_dtype, sparse):
     [[1,1,1], [1,1,1], [1,1,1]]
   ]
 
-  cc3d_labels = cc3d.connected_components(labels, connectivity=26, out_dtype=out_dtype, sparse=sparse)
-  scipy_labels, wow = scipy.ndimage.measurements.label(labels, structure=structure)
+  cc3d_labels, Ncc3d = cc3d.connected_components(labels, connectivity=26, return_N=True)
+  scipy_labels, Nscipy = scipy.ndimage.measurements.label(labels, structure=structure)
 
   print(cc3d_labels)
   print(scipy_labels)
 
+  assert Ncc3d == Nscipy
   assert np.all(cc3d_labels == scipy_labels)
 
-@pytest.mark.parametrize("sparse", (True, False))
-def test_compare_scipy_18(sparse):
+def test_compare_scipy_18():
   import scipy.ndimage.measurements
 
   sx, sy, sz = 256, 256, 256
@@ -447,27 +443,28 @@ def test_compare_scipy_18(sparse):
     [[0,1,0], [1,1,1], [0,1,0]]
   ]
 
-  cc3d_labels = cc3d.connected_components(labels, connectivity=18, sparse=sparse)
-  scipy_labels, wow = scipy.ndimage.measurements.label(labels, structure=structure)
+  cc3d_labels, Ncc3d = cc3d.connected_components(labels, connectivity=18, return_N=True)
+  scipy_labels, Nscipy = scipy.ndimage.measurements.label(labels, structure=structure)
 
   print(cc3d_labels)
   print(scipy_labels)
 
+  assert Ncc3d == Nscipy
   assert np.all(cc3d_labels == scipy_labels)
 
-@pytest.mark.parametrize("sparse", (True, False))
-def test_compare_scipy_6(sparse):
+def test_compare_scipy_6():
   import scipy.ndimage.measurements
 
   sx, sy, sz = 256, 256, 256
   labels = np.random.randint(0,2, (sx,sy,sz), dtype=np.bool)
 
-  cc3d_labels = cc3d.connected_components(labels, connectivity=6, sparse=sparse)
-  scipy_labels, wow = scipy.ndimage.measurements.label(labels)
+  cc3d_labels, Ncc3d = cc3d.connected_components(labels, connectivity=6, return_N=True)
+  scipy_labels, Nscipy = scipy.ndimage.measurements.label(labels)
 
   print(cc3d_labels)
   print(scipy_labels)
 
+  assert Ncc3d == Nscipy
   assert np.all(cc3d_labels == scipy_labels)
 
 @pytest.mark.skipif("sys.maxsize <= 2**33")
@@ -477,22 +474,24 @@ def test_sixty_four_bit():
   cc3d.connected_components(input_labels, max_labels=3)  
 
 @pytest.mark.parametrize("size", (255,256))
-def test_stress_upper_bound_for_binary_6(size):
+@pytest.mark.parametrize("zeroth_pass", (True,False))
+def test_stress_upper_bound_for_binary_6(size, zeroth_pass):
   labels = np.zeros((size,size,size), dtype=np.bool)
   for z in range(labels.shape[2]):
     for y in range(labels.shape[1]):
       off = (y + (z % 2)) % 2
       labels[off::2,y,z] = True
 
-  out = cc3d.connected_components(labels, connectivity=6)
+  out = cc3d.connected_components(labels, connectivity=6, zeroth_pass=zeroth_pass)
   assert np.max(out) + 1 <= (256**3) // 2 + 1
 
 @pytest.mark.parametrize("size", (255,256))
-def test_stress_upper_bound_for_binary_8(size):
+@pytest.mark.parametrize("zeroth_pass", (True,False))
+def test_stress_upper_bound_for_binary_8(size, zeroth_pass):
   labels = np.zeros((size,size), dtype=np.bool)
   labels[0::2,0::2] = True
 
-  out = cc3d.connected_components(labels, connectivity=8)
+  out = cc3d.connected_components(labels, connectivity=8, zeroth_pass=zeroth_pass)
   assert np.max(out) + 1 <= (256**2) // 4 + 1
 
   for _ in range(10):
@@ -501,49 +500,51 @@ def test_stress_upper_bound_for_binary_8(size):
     assert np.max(out) + 1 <= (256**2) // 4 + 1    
 
 @pytest.mark.parametrize("size", (255,256))
-def test_stress_upper_bound_for_binary_18(size):
+@pytest.mark.parametrize("zeroth_pass", (True,False))
+def test_stress_upper_bound_for_binary_18(size, zeroth_pass):
   labels = np.zeros((size,size,size), dtype=np.bool)
   labels[::2,::2,::2] = True
   labels[1::2,1::2,::2] = True
 
-  out = cc3d.connected_components(labels, connectivity=26)
+  out = cc3d.connected_components(labels, connectivity=18, zeroth_pass=zeroth_pass)
   assert np.max(out) + 1 <= (256**3) // 4 + 1
 
   for _ in range(10):
     labels = np.random.randint(0,2, (256,256,256), dtype=np.bool)
-    out = cc3d.connected_components(labels, connectivity=26)
-    assert np.max(out) + 1 <= (256**2) // 4 + 1
+    out = cc3d.connected_components(labels, connectivity=18, zeroth_pass=zeroth_pass)
+    assert np.max(out) + 1 <= (256**3) // 4 + 1
 
 @pytest.mark.parametrize("size", (255,256))
-def test_stress_upper_bound_for_binary_26(size):
+@pytest.mark.parametrize("zeroth_pass", (True,False))
+def test_stress_upper_bound_for_binary_26(size, zeroth_pass):
   labels = np.zeros((size,size,size), dtype=np.bool)
   labels[::2,::2,::2] = True
 
+  out = cc3d.connected_components(labels, connectivity=26, zeroth_pass=zeroth_pass)
+  assert np.max(out) + 1 <= (256**3) // 8 + 1
+
+  for _ in range(10):
+    labels = np.random.randint(0,2, (256,256,256), dtype=np.bool)
+    out = cc3d.connected_components(labels, connectivity=26, zeroth_pass=zeroth_pass)
+    assert np.max(out) + 1 <= (256**3) // 8 + 1
+
 @pytest.mark.parametrize("connectivity", (8, 18, 26))
 @pytest.mark.parametrize("dtype", TEST_TYPES)
-@pytest.mark.parametrize("out_dtype", OUT_TYPES)
 @pytest.mark.parametrize("order", ("C", "F"))
-@pytest.mark.parametrize("sparse", (True, False))
-def test_all_zeros_3d(connectivity, dtype, out_dtype, order, sparse):
+def test_all_zeros_3d(connectivity, dtype, order):
   labels = np.zeros((128,128,128), dtype=dtype, order=order)
-  out = cc3d.connected_components(labels, out_dtype=out_dtype, sparse=sparse)
+  out = cc3d.connected_components(labels)
   assert np.all(out == 0)
 
 @pytest.mark.parametrize("connectivity", (8, 18, 26))
 @pytest.mark.parametrize("dtype", TEST_TYPES)
-@pytest.mark.parametrize("out_dtype", OUT_TYPES)
 @pytest.mark.parametrize("order", ("C", "F"))
 @pytest.mark.parametrize("lbl", (1, 100, 7))
-@pytest.mark.parametrize("sparse", (True, False))
-def test_all_single_foreground(connectivity, dtype, out_dtype, order, lbl, sparse):
+def test_all_single_foreground(connectivity, dtype, order, lbl):
   labels = np.zeros((64,64,64), dtype=dtype, order=order) + lbl
-  out = cc3d.connected_components(labels, out_dtype=out_dtype, sparse=sparse)
+  out = cc3d.connected_components(labels)
   assert np.all(out == 1)
 
-# def test_sixty_four_bit():
-#   input_labels = np.ones((1626,1626,1626), dtype=np.uint8)
-#   cc3d.connected_components(input_labels, max_labels=0)  
-
 def test_region_graph_26():
   labels = np.zeros( (10, 10, 10), dtype=np.uint32 )