Add cucim.skimage.morphology.medial_axis #342

grlee77 · 2022-07-26T14:43:53Z

closes #336

This PR adds a function for skeletonization of 2D images via the medial axis transform.

It should be reviewed after #318 is merged. The new commits only start from 19a6fed.

There is one sequential component to this algorithm that still must be run on the CPU, but the majority of the computations are on the GPU and acceleration is good. I will add benchmark results here soon.

fix automated setting of m1, m2, m3 parameters for 2D PBA+

When switching to 64-bit int can go up to size 2**20 on each axis.

grlee77 · 2022-07-27T21:09:35Z

Acceleration is around 10x relative to scikit-image v0.19.3. That is pretty good considering that there is still a sequential portion that needs to run on the GPU.

image.shape = (328, 400)
    Duration (cpu) = 0.0446 +/= 0.0006
    Duration (gpu) = 0.0041 +/= 0.0001
    Acceleration: 11.00
image.shape = (1312, 1600)
    Duration (cpu) = 0.3426 +/= 0.0013
    Duration (gpu) = 0.0284 +/= 0.0003
    Acceleration: 12.07
image.shape = (2624, 3200)
    Duration (cpu) = 1.7736 +/= 0.0006
    Duration (gpu) = 0.1583 +/= 0.0033
    Acceleration: 11.21
image.shape = (5248, 6400)
    Duration (cpu) = 8.8874 +/= 0.0000
    Duration (gpu) = 0.9180 +/= 0.0051
    Acceleration: 9.68

### Overview closes #319 This PR adds an implementation of `scipy.ndimage.distance_transform_edt`. For each foreground pixel in a binary image, this function computes the minimal Euclidean distance to reach a background pixel. This is a SciPy function rather than a scikit-image one so I have put it under `cucim.core` instead of `cucim.skimage`. Longer term we could move this upstream to CuPy, but I think we want to have a copy here so that we can use it in the near term to implement some of the missing scikit-image functionality. This function is used by the following functions in the scikit-image API, so this PR will help us implement these in future PRs: `skimage.morphology.medial_axis` (**update**: see #342) `skimage.segmentation.expand_labels` (**update**: see #341) `skimage.segmentation.chan_vese` (**update**: see #343) It is also used in examples related to `skimage.segmentation.watershed` and is needed for implementation of ITK's `SignedMaurerDistanceMapImageFilter` for [itk-cucim](InsightSoftwareConsortium/itk_cucim#13). The algorithm used here is adapted from the C++ code in the [PBA+ repository](https://github.com/orzzzjq/Parallel-Banding-Algorithm-plus) (MIT-licensed) ### extensions made to PBA+ kernel source - larger sizes can be supported due to runtime code generation of the integer datatypes and encode/decode functions. ### current known limitations - 2D and 3D only (this should cover the most common use cases) - The `sampling` argument is not fully supported. This can likely be done with a bit of effort, but will require modifications to the CUDA kernels - User-specified output arguments are not currently supported. We can potentially add this in the future - The indices returned by `return_indices=True` are equally valid, but not always identical to those returned by SciPy. This is because there can be multiple indices with an identical distance, so which one gets chosen in that case is implementation-dependent. ### initial benchmarks relative to SciPy Here `% true` is the percentage of the image that corresponds to foreground pixels. Even for fairly small images there is some benefit, with the benefit becoming two orders of magnitude at larger sizes. shape | % true | cuCIM | SciPy | acceleration -----------------|--------|---------|---------|-------------- (256, 256) | 5 | 0.00108 | 0.00353 | 3.25 (512, 512) | 5 | 0.00108 | 0.01552 | 14.42 (2048, 2048) | 5 | 0.00252 | 0.34434 | 136.86 (4096, 4096) | 5 | 0.00765 | 1.58948 | 207.88 (32, 32, 32) | 5 | 0.00103 | 0.00305 | 2.98 (128, 128, 128) | 5 | 0.00153 | 0.30103 | 196.26 (256, 256, 256) | 5 | 0.00763 | 3.17872 | 416.37 (384, 384, 384) | 5 | 0.02460 | 14.28779 | 580.89 (256, 256) | 50 | 0.00107 | 0.00430 | 4.01 (512, 512) | 50 | 0.00109 | 0.01878 | 17.30 (2048, 2048) | 50 | 0.00299 | 0.39304 | 131.60 (4096, 4096) | 50 | 0.00896 | 1.84686 | 206.19 (32, 32, 32) | 50 | 0.00102 | 0.00361 | 3.53 (128, 128, 128) | 50 | 0.00163 | 0.31657 | 194.66 (256, 256, 256) | 50 | 0.00914 | 3.35914 | 367.49 (384, 384, 384) | 50 | 0.03005 | 13.83219 | 460.30 (256, 256) | 95 | 0.00107 | 0.00344 | 3.22 (512, 512) | 95 | 0.00108 | 0.01638 | 15.17 (2048, 2048) | 95 | 0.00314 | 0.36996 | 117.81 (4096, 4096) | 95 | 0.00943 | 1.90475 | 202.05 (32, 32, 32) | 95 | 0.00102 | 0.00314 | 3.07 (128, 128, 128) | 95 | 0.00180 | 0.28843 | 159.80 (256, 256, 256) | 95 | 0.01073 | 3.23450 | 301.41 (384, 384, 384) | 95 | 0.03577 | 12.40526 | 346.82 ### other comments - can likely reduce memory overhead and improve performance a bit by refactoring some of the pre/post-processing code into elementwise kernels. (e.g. `encode3d`/`decode3d`, etc.) - (JK) This may be able to leverage CuPy in the future ( cupy/cupy#6919 ) Authors: - Gregory Lee (https://github.com/grlee77) Approvers: - Gigon Bae (https://github.com/gigony) - https://github.com/jakirkham URL: #318

jakirkham · 2022-07-28T06:54:11Z

@grlee77, can you please take a look at the conflict here? Looks related to PR ( #318 )

…axis

grlee77 · 2022-07-28T14:12:31Z

fixed the conflict

jakirkham · 2022-08-01T23:07:04Z

@gigony, could you please review this as well? 🙂

gigony

Looks good to me!

gigony · 2022-08-02T22:27:57Z

python/cucim/src/cucim/skimage/morphology/tests/test_skeletonize.py

+
+class TestMedialAxis():
+    def test_00_00_zeros(self):
+        '''Test skeletonize on an array of all zeros'''


nit: just wondering if there is a reason to use ''' instead of """ in this file.

Yeah we should probably move to """. This is what black does

Happy to change this to """ and we can probably apply black at some point, but I didn't do it so far for two reasons
1.) It messes with readability of various small test array definitions unless # fmt: off, etc. are placed around them
2.) easier to compare diffs with upstream scikit-image without it.

jakirkham · 2022-08-03T10:28:21Z

Fixed a small conflict above. If that looks ok, please feel free to merge 🙂

grlee77 · 2022-08-03T11:08:43Z

@gpucibot merge

grlee77 · 2022-08-03T11:58:21Z

@jakirkham: I tried to merge via gpucibot in my previous comment, but it didn't seem to listen to me. Do I need to be granted some permissions here for it to work?

raydouglass · 2022-08-03T17:33:18Z

@gpucibot merge

jakirkham · 2022-08-03T17:41:20Z

Thanks Ray! 🙏

grlee77 added 13 commits June 27, 2022 11:17

add license info for the PBA+ algorithm

8e9b5a5

Add distance_transform_edt implementation based on PBA+

0462a2c

add tests for 2d distance transform

6a9997f

fix automated setting of m1, m2, m3 parameters for 2D PBA+

Update 3D kernels to allow use of size > 1024 per axis

3d44899

When switching to 64-bit int can go up to size 2**20 on each axis.

expand tests and extend to 3D cases

8509bff

test error on non-uniform sampling

0933d73

flake8/isort fixes

e98b611

minor update to kernel headers

e5a37a9

update outdated comment

e83d9a7

add fallback implementation of math.lcm for Python < 3.9

c67b7c4

make sure fallback lcm function is actually used if needed

cb06d78

add medial_axis

19a6fed

pep8

cce6185

grlee77 added feature request New feature or request non-breaking Introduces a non-breaking change labels Jul 26, 2022

grlee77 added this to the v22.08.00 milestone Jul 26, 2022

grlee77 requested a review from a team as a code owner July 26, 2022 14:43

pep8 fix

7adce8f

grlee77 mentioned this pull request Jul 27, 2022

Add Euclidean distance transform for images/volumes #318

Merged

Merge remote-tracking branch 'upstream/branch-22.08' into add-medial_…

3f4f981

…axis

jakirkham approved these changes Aug 1, 2022

View reviewed changes

jakirkham requested a review from gigony August 1, 2022 23:06

gigony approved these changes Aug 2, 2022

View reviewed changes

grlee77 and others added 2 commits August 2, 2022 22:07

consistently use double quotes

bcd24b4

Merge branch 'branch-22.08' into add-medial_axis

67d16fc

rapids-bot bot merged commit c834039 into rapidsai:branch-22.08 Aug 3, 2022

Provide feedback

Saved searches

Use saved searches to filter your results more quickly

Add cucim.skimage.morphology.medial_axis #342

Add cucim.skimage.morphology.medial_axis #342

grlee77 commented Jul 26, 2022 •

edited

Loading

grlee77 commented Jul 27, 2022

jakirkham commented Jul 28, 2022

grlee77 commented Jul 28, 2022

jakirkham commented Aug 1, 2022

gigony left a comment

gigony Aug 2, 2022

jakirkham Aug 2, 2022

grlee77 Aug 3, 2022

jakirkham commented Aug 3, 2022

grlee77 commented Aug 3, 2022

grlee77 commented Aug 3, 2022

raydouglass commented Aug 3, 2022

jakirkham commented Aug 3, 2022

Add cucim.skimage.morphology.medial_axis #342

Add cucim.skimage.morphology.medial_axis #342

Conversation

grlee77 commented Jul 26, 2022 • edited Loading

grlee77 commented Jul 27, 2022

jakirkham commented Jul 28, 2022

grlee77 commented Jul 28, 2022

jakirkham commented Aug 1, 2022

gigony left a comment

Choose a reason for hiding this comment

gigony Aug 2, 2022

Choose a reason for hiding this comment

jakirkham Aug 2, 2022

Choose a reason for hiding this comment

grlee77 Aug 3, 2022

Choose a reason for hiding this comment

jakirkham commented Aug 3, 2022

grlee77 commented Aug 3, 2022

grlee77 commented Aug 3, 2022

raydouglass commented Aug 3, 2022

jakirkham commented Aug 3, 2022

grlee77 commented Jul 26, 2022 •

edited

Loading