8000 [WIP] Adding @eldad-a's ridge directed ring detector by alexdesiqueira · Pull Request #4847 · scikit-image/scikit-image · GitHub
[go: up one dir, main page]
Skip to content

[WIP] Adding @eldad-a's ridge directed ring detector #4847

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Jump to bottom
Draft
wants to merge 35 commits into
base: main
Choose a base branch
from
Draft

[WIP] Adding @eldad-a's ridge directed ring detector #4847

wants to merge 35 commits into from

Conversation

alexdesiqueira
Copy link
Member

Description

Continues #1706, by @eldad-a:

an implementation of the algorithm in
Sci. Rep. 5, 13584; doi: 10.1038/srep13584 (2015).

Checklist

For reviewers

  • Check that the PR title is short, concise, and will make sense 1 year
    later.
  • Check that new functions are imported in corresponding __init__.py.
  • Check that new features, API changes, and deprecations are mentioned in
    doc/release/release_dev.rst.

@pep8speaks
Copy link
pep8speaks commented Jul 12, 2020 •
edited
Loading

Hello @alexdesiqueira! Thanks for updating this PR. We checked the lines you've touched for PEP 8 issues, and found:

Line 41:57: W605 invalid escape sequence '*'
Line 41:63: W605 invalid escape sequence '*'
Line 68:80: E501 line too long (143 > 79 characters)
Line 75:80: E501 line too long (153 > 79 characters)

Line 35:1: E302 expected 2 blank lines, found 1

Comment last updated at 2022-06-13 23:18:54 UTC

@alexdesiqueira alexdesiqueira mentioned this pull request Jul 12, 2020
Closed
@alexdesiqueira
Copy link
Member Author
alexdesiqueira commented Jul 12, 2020 •
edited
Loading

The following is based on the discussion I had with @clementkng on Slack during the SciPy 2020 sprint, and what I think we need to pursue for this PR.

What we have:

skimage/transform/ridge-directed-ring-detector/demo.ipynb
skimage/transform/ridge-directed-ring-detector/unprocessed_fig_46.png
skimage/transform/_ring_detector.pyx

I'd like to:
[OK] 1. create a function within transform — let's say ring_detector.py — that would have the API necessary to call ring_detector.pyx. Something on the line of skimage/transform/hough_transform.py or skimage/transform/radon_transform.py;
2. move demo.ipynb to doc/examples. The notebook would be converted to a .py file — as the other files in our documentation. It'd go to doc/examples/transform/. A possible name would be plot_ring_detector.py.
3. remove the ridge-directed-ring-detector folder.
[OK] 4. replace OpenCV functions on _ring_detector.pyx with skimage ones.

@eldad-a
Copy link
eldad-a commented Jul 13, 2020 •
edited
Loading

Hello @alexdesiqueira

Many thanks for moving forward with this PR, and sorry for my losing track of the process.

I should have a draft which is opencv independent - where all spatial derivatives are performed by skimage functions.

If that's helpful, please let me know where to push it to.

Thanking you,
Eldad

@alexdesiqueira alexdesiqueira mentioned this pull request Jul 14, 2020
Closed
@alexdesiqueira
Copy link
Member Author

Hey @eldad-a ,
don't worry at all — it happens to all of us 😄

I should have a draft which is opencv independent - where all spatial derivatives are performed by skimage functions.

If that's helpful, please let me know where to push it to.

That would be great, thank you! I will add you as a contributor, please feel free to push it here when you can.

mkcor
mkcor reviewed Jul 17, 2020
View reviewed changes
doc/examples/transform/plot_ring_detector.py Outdated

from skimage.transform import ring_detector
from skimage import data
from scipy import ndimage
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

import ndimage before skimage

doc/examples/transform/plot_ring_detector.py Outdated
using 2d microscope imaging.

Here's a quick description of the parameters. The explanations refer to the
terms as introduced in [Afik 2015](https://www.nature.com/articles/srep13584)
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Looks like you're using Markdown syntax. In my experience, it should be reStructuredText.

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

indeed :-)

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

you can take a look at existing gallery scripts for examples of rst syntax.

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Sorry about that @mkcor @emmanuelle, but this is too raw for a review yet. 😄 I asked @clementkng to move the code from a notebook to this file, starting the documentation process. We'll take care of it 🙂

skimage/transform/ring_detector.py Outdated
@@ -0,0 +1,30 @@
from ._ring_detector import RidgeHoughTransform

import matplotlib.pyplot as plt
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Please swap import order.
Reference: https://www.python.org/dev/peps/pep-0008/#imports

@alexdesiqueira
Copy link
Member Author

Thank you @mkcor! This is still very raw, not ready for a review, though.
Could we ping you when this is somewhat acceptable? 😅 Thank you again!

@mkcor
Copy link
Member
mkcor commented Jul 20, 2020

@alexdesiqueira sure, no worries! Your PR title says 'WIP' indeed :)

@eldad-a
Copy link
eldad-a commented Jul 24, 2020 •
edited
Loading

The following is based on the discussion I had with @clementkng on Slack during the SciPy 2020 sprint, and what I think we need to pursue for this PR.

What we have:
...
4. replace OpenCV functions on _ring_detector.pyx with skimage ones.

Thank you all for your patience with me!
Hope everyone is doing well during these different times.

  1. OpenCV dependency removed; currently replaced by skimage.feature.hessian_matrix.
  2. Updated curve_thresh parameter in the example plot_ring_detector.py accordingly.

Noticed that ring_detector.py returns the rings rather than the subpixel rings.
Where should the alternative be highlighted for users?
I would support having the subpixel rings as the default.

I've been exploring the outcome on the use-case for which this code was developed originally (an example was provided in ridge-directed-ring-detector/unprocessed_fig_46.png ; BTW, was this file removed?).
Included some observations and suggestions in the _ring_detector.pyx itself.
Happy to share the ipynb where I explored using scharr and farid; in case of interest, please let me know where to put them.
For convenience, I include the notes below as well; would be great to hear your take on these!

EA NOTES

  • The original opencv version shows better performance compared to the current skimage:
    • preprocessing takes x3 longer using hessian_matrix (from skimage), and
    • ring_detection takes about x1.2, showing higher sensitivity to choice of parameter (sigma and curv_thresh)
  • The results seem more precise and less sensitive to choice of parameters when derivatives are estimated using the filters.farid_v/h (or filters.scharr_v/h); however, these operators take x6 (x3) times longer to compute, w.r.t to hessian_matrix, which uses np.gradient; this could probably be improved by programming the more precise (but slower) operators to allow the 2nd order derivatives directly.

EA TO-Consider

  • store both principal curvatures and include verifying that |Lpp|>|Lqq| ; here Lpp denotes the least principal curvature (Lqq the largest principal curvature).
  • allow the user to choose between np.gradient (skimage hessian_matrix), scharr, and farid operators; the latter seem to offer more precise results, with less sensitivity to the choice of parameters. The former is faster.
  • add automated choice of curv_thresh based on percentiles.
  • verify results and efficiency are not significantly affected the transition from opencv to skimage (Hessian calculation)
  • Prepare PR for directed ridge detector ? There are several related algorithms in skimage; does any of them use non-maximum suppression in the direction of the eigen-vector associated with the least-principal curvature? (as performed in the ring_detection ). Explore the effect of separating the ridge detection function entirely; it may not affect the %%timeit results.

@alexdesiqueira
Copy link
Member Author

Hey @eldad-a,

  1. OpenCV dependency removed; currently replaced by skimage.feature.hessian_matrix.
  2. Updated curve_thresh parameter in the example plot_ring_detector.py accordingly.

Thank you very much, we appreciate your help!

Noticed that ring_detector.py returns the rings rather than the subpixel rings.
Where should the alternative be highlighted for users?
I would support having the subpixel rings as the default.

I'd ask your help to determine that technical features, Eldad. I still don't understand the algorithm that much, so please feel free to check what you consider doable on that subject — I would take a time to understand 🙂

I've been exploring the outcome on the use-case for which this code was developed originally (an example was provided in ridge-directed-ring-detector/unprocessed_fig_46.png ; BTW, was this file removed?).

I'd like to structure the example on the lines of the examples we have already — hence doc/examples/transform/plot_ring_detector.py. Would you like to base the example on that specific file? If yes, I'd suggest to check the license for that file. If it's permissible, we could add it to data and start the example from there; what do you think?

  • allow the user to choose between np.gradient (skimage hessian_matrix), scharr, and farid operators; the latter seem to offer more precise results, with less sensitivity to the choice of parameters. The former is faster.

I like this! I'd need to understand the other points better, though; what would be their implications, you know... I need to read your paper again 🙂
Thanks again!

@eldad-a
Copy link
eldad-a commented Aug 17, 2020 •
edited
Loading

Hi @alexdesiqueira :-)

  1. OpenCV dependency removed; currently replaced by skimage.feature.hessian_matrix.
  2. Updated curve_thresh parameter in the example plot_ring_detector.py accordingly.

Thank you very much, we appreciate your help!

It is my pleasure;
I thank you for your support and encouragement!

Noticed that ring_detector.py returns the rings rather than the subpixel rings.
Where should the alternative be highlighted for users?
I would support having the subpixel rings as the default.

I'd ask your help to determine that technical features, Eldad. I still don't understand the algorithm that much, so please feel free to check what you consider doable on that subject — I would take a time to understand slightly_smiling_face

Sure, I will.

I've been exploring the outcome on the use-case for which this code was developed originally (an example was provided in ridge-directed-ring-detector/unprocessed_fig_46.png ; BTW, was this file removed?).

I'd like to structure the example on the lines of the examples we have already — hence doc/examples/transform/plot_ring_detector.py. Would you like to base the example on that specific file? If yes, I'd suggest to check the license for that file. If it's permissible, we could add it to data and start the example from there; what do you think?

Which data image to use... Why not both?
#3384 calls for scientific examples, so this is an opportunity.

My understanding of licensing goes a very short distance (well, maybe none...).
Here's what the original publication reads in terms of license: Rights and permissions.
The figure is part of Supplementary Figure S2(a)
in Afik, E. Robust and highly performant ring detection algorithm for 3d particle tracking using 2d microscope imaging. Sci. Rep. 5, 13584; doi: 10.1038/srep13584 (2015).

  • allow the user to choose between np.gradient (skimage hessian_matrix), scharr, and farid operators; the latter seem to offer more precise results, with less sensitivity to the choice of parameters. The former is faster.

I like this! I'd need to understand the other points better, though; what would be their implications, you know... I need to read your paper again

If you have questions, please let me know; may be easier than going over the paper.

Below there's a TODO list I created for my future reference; some are added features, so I guess aren't a must.
Could you share what are the priorities to have this PR ready?

EA TODO:

  • ring_detector.py : return subpixel precision ring parameters
  • allow the user to choose between np.gradient (skimage hessian_matrix), scharr, and farid operators; the latter seem to offer more precise results, with less sensitivity to the choice of parameters. The former is faster.
  • test: store both principal curvatures and include verifying that |Lpp|>|Lqq| ; here Lpp denotes the least principal curvature (Lqq the largest principal curvature).
  • add automated choice of curv_thresh based on percentiles.
  • verify results and efficiency are not significantly affected the transition from opencv to skimage (Hessian calculation)

In case skimage community has interest:

  • Prepare PR for directed ridge detector ? There are several related algorithms in skimage; does any of them use non-maximum suppression in the direction of the eigen-vector associated with the least-principal curvature? (as performed in the ring_detection ).
    • Explore the effect of separating the ridge detection function entirely; it may not affect the %%timeit results.

Base automatically changed from master to main February 18, 2021 18:23
@mkcor mkcor mentioned this pull request Feb 22, 2021
Closed
@stefanv
Copy link
Member
stefanv commented Apr 26, 2021

I think this one is slipping again.

@alexdesiqueira
Copy link
Member Author

@stefanv I have this in my mind for a while, now. 🙂 I'll come back to that next week.

@alexdesiqueira
Copy link
Member Author
alexdesiqueira commented May 4, 2021 •
edited
Loading

Alright, (cracking knuckles 🤜🤛 )
let's get back!

@eldad-a:

Which data image to use... Why not both?
#3384 calls for scientific examples, so this is an opportunity.

Of course 🙂 I'd like to see a well-structured example, though. Maybe using these images to show different aspects of the algorithm? Would that be feasible?

My understanding of licensing goes a very short distance (well, maybe none...).
Here's what the original publication reads in terms of license: Rights and permissions.
The figure is part of Supplementary Figure S2(a)
in Afik, E. Robust and highly performant ring detection algorithm for 3d particle tracking using 2d microscope imaging. Sci. Rep. 5, 13584; doi: 10.1038/srep13584 (2015).

The link says that the license is CC-BY, Eldad; we'll be good if the repository has an attribution to the paper. we are not able to distribute these images with scikit-image.

Below there's a TODO list I created for my future reference; some are added features, so I guess aren't a must.
Could you share what are the priorities to have this PR ready?

I would be happy if we have, for this PR:

  • a running version of the algorithm — the first one you presented us is great already;
  • a documentation page on how to use the algorithm;
  • basic tests.

We can implement the TODO's you proposed in subsequent PRs, Eldad. How does that sound?

@stefanv
Copy link
Member
stefanv commented May 5, 2021

The link says that the license is CC-BY, Eldad; we'll be good if the repository has an attribution to the paper.

Is this for an image that lives inside skimage.data? We only add CC0 and public domain images there, I think, since we cannot easily have our users display the correct copyright otherwise.

@alexdesiqueira
Copy link
Member Author

Is this for an image that lives inside skimage.data? We only add CC0 and public domain images there, I think, since we cannot easily have our users display the correct copyright otherwise.

Oh, sorry @stefanv; you're right! Do you have any images you didn't use in that paper — and you would consider donating — by any chance, @eldad-a?

@alexdesiqueira
Copy link
Member Author
alexdesiqueira commented May 5, 2021 •
edited
Loading

Also, @scikit-image/core may I ask some help on "beautifying" the cython code here? Thank you!
EDIT: I should say "refactoring," I know 🙂

@alexdesiqueira alexdesiqueira added this to the 0.19 milestone May 11, 2021
seberg
seberg reviewed May 11, 2021
View reviewed changes
Copy link
Contributor
@seberg seberg left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Sorry, this is super random. High level, I think mostly that the parameter handling should maybe be hidden away (either into the init or as function parameters to the second function). And some of the others, like the derivatives should probably just be cdef class attributes, so nobody can mess with them.

skimage/transform/_ring_detector.pyx Outdated
Comment on lines 162 to 170
#cdef:
# DTYPE_t trH
# DTYPE_t discriminant
# DTYPE_t curvature
#trH = Lrr+Lcc
# numerically stabilising for the cases of extrimly small Lrc:
#discriminant = (Lrr-Lcc)**2 + 4*Lrc*Lrc
#curvature = 0.5*(trH - sqrt( discriminant ))
#return curvature
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Seems like some unnecessary performance tries. I guess just delete (or use this instead of the actual one).

skimage/transform/_ring_detector.pyx Outdated
kernel = np.empty((ksize),DTYPE)
scale2x = -0.5/sigma**2

for i in xrange(ksize):
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
for i in xrange(ksize):
for i in range(ksize):

Use range throughout. Do you use the cython language level? Might set it to 3...

skimage/transform/_ring_detector.pyx Outdated
Comment on lines 330 to 331
assert Lrr is not None and Lcc is not None and Lrc is not None and\
curv is not None
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This seems wrong to me. If this is exposed to python, we probably need more than that. Even f it works, it would give a weird error i think?

I am not sure what the best option is here, maybe even create one cdef version, and then a very brief def ... version for Python that uses DType_t[:, :] Lrr not None.

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Honestly, there is probably no use for the cdef here, so just using def is maybe the best way.

skimage/transform/_ring_detector.pyx Outdated
Comment on lines 736 to 737
RIJ_t [:,:,:] sparse_3d_Hough
RIJ_t [:,:] voted4
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
RIJ_t [:,:,:] sparse_3d_Hough
RIJ_t [:,:] voted4
RIJ_t[:, :, ::1] sparse_3d_Hough
RIJ_t[:, ::1] voted4

It can be typed as contiguous.

(I prefer that spacing, but I never really coded reviewed cython, so I am not sure, whats typical.)

skimage/transform/_ring_detector.pyx Outdated
Comment on lines 798 to 799
DTYPE_t [:,:,:] smoothed_hough_array
DTYPE_t [:] kernel
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
DTYPE_t [:,:,:] smoothed_hough_array
DTYPE_t [:] kernel
DTYPE_t[:, :, ::1] smoothed_hough_array
DTYPE_t[::1] kernel

skimage/transform/_ring_detector.pyx Outdated
ring_mask[row_min-i+r+thickness:r+thickness+row_max-i,\
col_min-j+r+thickness:r+thickness+col_max-j]
))
if False: # eccentricity: ## excluding
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

never used branch? Should this be around?

skimage/transform/_ring_detector.pyx Outdated
RIJ_t [:,:,:] hough_slice # a.k.a. sparse_3d_Hough
DTYPE_t [:,:,:] smoothed_slice # a.k.a. smoothed_hough_array
RIJ_t [:,:,:] hough_hotspots, hough_modified
RIJ_t [:] hough_counter, hotspots_counter
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

There seem all allocated (and some below as well), so it probably makes sense to type them as RIJ_t[:, :, ::1] ... (i.e. end with ::1 to tell cython that they are contiguous).

skimage/transform/_ring_detector.pyx Outdated
#from collections import defaultdict

import numpy as np
cimport numpy as np
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Do you use this or cimport numpy as npc? My guess is, there is only a few things that would need changing (mainly the typedefs ending with _t)

skimage/transform/_ring_detector.pyx Outdated
return 1./denominator, tangent / denominator


#@cython.profile(False)
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Delete these, at least when they are commented out?

skimage/transform/_ring_detector.pyx Outdated

dr = Rmax - Rmin
sr=1
coords = np.empty((dr+1,2), dtype=INDX)
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

If we go so far as to inline this function... this function should just not return anything (or if it can error int, but I don't think it can). And coord should be an input and modified in-place.

The empty call should be early on (before the loop) in the function where it is used. (Also type coord as [:, ::1] so cython can optimize it as being contiguous.

@alexdesiqueira alexdesiqueira self-assigned this May 14, 2021
@eldad-a
Copy link
eldad-a commented Jun 16, 2021

@alexdesiqueira

Could you please guide me through the following:

  1. Where should example images be added to (folder)?
  2. Attempting to compile commit cbab232 using the command (found at #installation-from-source ):
python setup.py build_ext -i

results in :

ValueError: '.../scikit-image/skimage/transform/_ring_detector.pyx' doesn't match any files

Any suggestions? Should I test using my original PR and commit the test files here?

@alexdesiqueira
Copy link
Member Author

Hey @eldad-a,
sorry about the silence.

  1. Where should example images be added to (folder)?

It would be best to add these in our GitLab repo.

  1. Attempting to compile commit cbab232 using the command...
    Any suggestions? Should I test using my original PR and commit the test files here?

Sorry about that, Eldad 😅 I broke the PR while converting from procedural to OO, and I didn't have time to fix it properly yet. We need to think on an API that could agree with the other transforms, so we can write the tests on that scheme. Minimally, what do you think that the function would need as input and output?
Thank you!

@eldad-a
Copy link
eldad-a commented Jun 16, 2021 via email

@alexdesiqueira
Copy link
Member Author

Hey @eldad-a,

Would skimage circle_hough transform be a good reference?

Sure! I think all our hough_ transforms would be great for comparison.
Thank you!

clementkng and others added 24 commits May 24, 2022 16:23
@clementkng
DOC: Refactor ridge directed ring detector example
2822a36
@mkcor
Update doc/examples/transform/plot_ring_detector.py
2827dbe 
Co-authored-by: Marianne Corvellec <marianne.corvellec@ens-lyon.org>
@mkcor
Update doc/examples/transform/plot_ring_detector.py
d8b15a5 
Co-authored-by: Marianne Corvellec <marianne.corvellec@ens-lyon.org>
@eldad-a
removed opencv dependencies from ring_detector
06b9432 
1. Replace OpenCV Hessian estimation by skimage ones

2. Comment-out 'eccentricity', which relies on OpenCV for EllipseFit

3. Rename, following skimage conventions:
Lxx -> Lrr , Lyy -> Lcc , Lxy -> Lrc

3. Add performance observations and TODO to consider
@eldad-a
adjust 'curv_thresh' following the transition from OpenCV to skimage
6afb9e4
@clementkng
DOC: Fix import order and convert Markdown to RST
46591e1
Starting to un-OO _ring_detector
4af6f55
Rewriting img_preprocess -> preprocessing
c901cb5
RidgeHoughTransform -> hough_transform_ridge
c6244b5
Deleting unused: least_principal_curvature
51c939f
Trimming comments, xrange -> range
1d6e4b3
Adding spaces for readability, removing comments
6d0a6d0
Removing debugging function
027a9f1
DTYPE_t -> np_floats
5962f9e
pi -> M_PI
b56a426
RIJ, DTYPE, INDX -> np_uints, np_ints, np_floats
fb818a6
Fixing 'Fused types not allowed here'
13dc170
Fixing cnp_float, int, uint -> double, cnp.uint32
4e04596
Move ring_detector -> hough_..._ridge; remove .py file
d7ec00e
Removing a small todo to check if my copy still works
52740fd
@tibuch
Add test for central moment computation with spacing and fix implemen…
53e29ff 
…tation.
@tibuch
Add moments test with spacing.
91da300
@rfezzani
Fix first pass declaration type
7723fac
colour -> color on examples/transform/plot_ring_detector
ad5d7b4
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Reviewers

@seberg seberg seberg left review comments

@emmanuelle emmanuelle emmanuelle left review comments

@mkcor mkcor mkcor left review comments

At least 1 approving review is required to merge this pull request.

Assignees

@alexdesiqueira alexdesiqueira

Labels
💪 Work in progress 🙏 Feature request
Projects
None yet
Milestone
No milestone
Development

Successfully merging this pull request may close these issues.
0