wassname/scikit-image
1
0
Fork 0
mirror of https://github.com/wassname/scikit-image.git synced 2026-07-08 00:34:48 +08:00
Code Issues Packages Projects Releases Wiki Activity

Merge pull request #597 from sciunto/hough_ellipse

Browse Source 
Hough transform for ellipses
This commit is contained in:
Johannes Schönberger
2013-06-21 12:47:50 -07:00
parent d8d5e30920 33bbda17c1
commit 360f5b55e3
4 changed files with 168 additions and 9 deletions
Download Patch File Download Diff File Expand all files Collapse all files
CONTRIBUTORS.txt
+1 -1
View File
@@ -132,7 +132,7 @@
- François Boulogne
Drawing: Andres Method for circle perimeter, ellipse perimeter drawing, Bezier curve.
Circular Hough Transform
Circular and elliptical Hough Transforms
Various fixes
- Thouis Jones
skimage/transform/__init__.py
+2 -1
View File
@@ -1,4 +1,4 @@
from ._hough_transform import (hough_circle, hough_line,
from ._hough_transform import (hough_circle, hough_ellipse, hough_line,
probabilistic_hough_line)
from .hough_transform import (hough, probabilistic_hough, hough_peaks,
hough_line_peaks)
@@ -15,6 +15,7 @@ from .pyramids import (pyramid_reduce, pyramid_expand,
__all__ = ['hough_circle',
'hough_ellipse',
'hough_line',
'probabilistic_hough_line',
'hough',
skimage/transform/_hough_transform.pyx
+127 -5
View File
@@ -101,7 +101,129 @@ def hough_circle(cnp.ndarray img,
return acc
def hough_line(cnp.ndarray img, cnp.ndarray[ndim=1, dtype=cnp.double_t] theta=None):
def hough_ellipse(cnp.ndarray img, int threshold=4, double accuracy=1,
int min_size=4, max_size=None):
"""Perform an elliptical Hough transform.
Parameters
----------
img : (M, N) ndarray
Input image with nonzero values representing edges.
threshold: int, optional (default 4)
Accumulator threshold value.
accuracy : double, optional (default 1)
Bin size on the minor axis used in the accumulator.
min_size : int, optional (default 4)
Minimal major axis length.
max_size : int, optional
Maximal minor axis length. (default None)
If None, the value is set to the half of the smaller
image dimension.
Returns
-------
res : list of tuples [(x0, y0, a, b, angle, accumulator)]
Where (x0, y0) is the center, (a, b) major and minor axis.
The angle value follows `draw.ellipse_perimeter()` convention.
Examples
--------
>>> img = np.zeros((25, 25), dtype=int)
>>> rr, cc = draw.ellipse_perimeter(10, 10, 6, 8)
>>> img[rr, cc] = 1
>>> result = hough_ellipse(img, threshold=6)
[(10.0, 10.0, 8.0, 6.0474292058692187, 0.0, 8)]
Notes
-----
The accuracy must be chosen to produce a peak in the accumulator
distribution. In other words, a flat accumulator distribution with low
values may be caused by a too low bin size.
References
----------
.. [1] Xie, Yonghong, and Qiang Ji. "A new efficient ellipse detection
method." Pattern Recognition, 2002. Proceedings. 16th International
Conference on. Vol. 2. IEEE, 2002
"""
if img.ndim != 2:
raise ValueError('The input image must be 2D.')
cdef long[:, :] pixels = np.transpose(np.nonzero(img))
cdef Py_ssize_t num_pixels = pixels.shape[0]
cdef list acc = list()
cdef list results = list()
cdef bin_size = accuracy**2
cdef int max_b_squared
if max_size is None:
if img.shape[0] < img.shape[1]:
max_b_squared = np.round(0.5 * img.shape[0])**2
else:
max_b_squared = np.round(0.5 * img.shape[1])**2
else:
max_b_squared = max_size**2
cdef Py_ssize_t p1, p2, p3, p1x, p1y, p2x, p2y, p3x, p3y
cdef double x0, y0, a, b, d, k
cdef double cos_tau_squared, b_squared, f_squared, angle
for p1 in range(num_pixels):
p1x = pixels[p1, 1]
p1y = pixels[p1, 0]
for p2 in range(p1):
p2x = pixels[p2, 1]
p2y = pixels[p2, 0]
# Candidate: center (x0, y0) and main axis a
a = 0.5 * sqrt((p1x - p2x)**2 + (p1y - p2y)**2)
if a > 0.5 * min_size:
x0 = 0.5 * (p1x + p2x)
y0 = 0.5 * (p1y + p2y)
for p3 in range(num_pixels):
p3x = pixels[p3, 1]
p3y = pixels[p3, 0]
d = sqrt((p3x - x0)**2 + (p3y - y0)**2)
if d > min_size:
f_squared = (p3x - p1x)**2 + (p3y - p1y)**2
cos_tau_squared = ((a**2 + d**2 - f_squared) \
/ (2 * a * d))**2
# Consider b2 > 0 and avoid division by zero
k = a**2 - d**2 * cos_tau_squared
if k > 0 and cos_tau_squared < 1:
b_squared = a**2 * d**2 * (1 - cos_tau_squared) / k
# b2 range is limited to avoid histogram memory
# overflow
if b_squared <= max_b_squared:
acc.append(b_squared)
if len(acc) > 0:
bins = np.arange(0, np.max(acc) + bin_size, bin_size)
hist, bin_edges = np.histogram(acc, bins=bins)
hist_max = np.max(hist)
if hist_max > threshold:
angle = np.arctan2(p1x - p2x, p1y - p2y)
# pi - angle to keep ellipse_perimeter() convention
if angle != 0:
angle = np.pi - angle
b = sqrt(bin_edges[hist.argmax()])
results.append((x0,
y0,
a,
b,
angle,
hist_max, # Accumulator
))
acc = []
return results
def hough_line(cnp.ndarray img,
cnp.ndarray[ndim=1, dtype=cnp.double_t] theta=None):
"""Perform a straight line Hough transform.
Parameters
@@ -211,8 +333,8 @@ def probabilistic_hough_line(cnp.ndarray img, int threshold=10,
Returns
-------
lines : list
List of lines identified, lines in format ((x0, y0), (x1, y0)), indicating
line start and end.
List of lines identified, lines in format ((x0, y0), (x1, y0)),
indicating line start and end.
References
----------
@@ -334,14 +456,14 @@ def probabilistic_hough_line(cnp.ndarray img, int threshold=10,
y1 = py >> shift
else:
x1 = px >> shift
y1 = py;
y1 = py
# check when line exits image boundary
if x1 < 0 or x1 >= width or y1 < 0 or y1 >= height:
break
gap += 1
# if non-zero point found, continue the line
if mask[y1, x1]:
gap = 0;
gap = 0
line_end[k, 1] = y1
line_end[k, 0] = x1
# if gap to this point was too large, end the line
skimage/transform/tests/test_hough_transform.py
+38 -2
View File
@@ -2,7 +2,7 @@ import numpy as np
from numpy.testing import *
import skimage.transform as tf
from skimage.draw import circle_perimeter, line
from skimage.draw import line, circle_perimeter, ellipse_perimeter
def append_desc(func, description):
@@ -126,6 +126,7 @@ def test_hough_circle():
assert_equal(x[0], x_0)
assert_equal(y[0], y_0)
def test_hough_circle_extended():
# Prepare picture
# The circle center is outside the image
@@ -133,7 +134,7 @@ def test_hough_circle_extended():
radius = 20
x_0, y_0 = (-5, 50)
y, x = circle_perimeter(y_0, x_0, radius)
img[x[np.where(x>0)], y[np.where(x>0)]] = 1
img[x[np.where(x > 0)], y[np.where(x > 0)]] = 1
out = tf.hough_circle(img, np.array([radius]), full_output=True)
@@ -142,5 +143,40 @@ def test_hough_circle_extended():
assert_equal(x[0], x_0 + radius)
assert_equal(y[0], y_0 + radius)
def test_hough_ellipse_zero_angle():
img = np.zeros((25, 25), dtype=int)
a = 6
b = 8
x0 = 12
y0 = 12
angle = 0
rr, cc = ellipse_perimeter(x0, x0, b, a)
img[rr, cc] = 1
result = tf.hough_ellipse(img, threshold=9)
assert_equal(result[0][0], x0)
assert_equal(result[0][1], y0)
assert_almost_equal(result[0][2], b, decimal=1)
assert_almost_equal(result[0][3], a, decimal=1)
assert_equal(result[0][4], angle)
def test_hough_ellipse_non_zero_angle():
img = np.zeros((20, 20), dtype=int)
a = 6
b = 9
x0 = 10
y0 = 10
angle = np.pi/1.35
rr, cc = ellipse_perimeter(x0, x0, b, a, orientation=angle)
img[rr, cc] = 1
result = tf.hough_ellipse(img, threshold=15, accuracy=3)
print(result)
assert_almost_equal(result[0][0]/100., x0/100., decimal=1)
assert_almost_equal(result[0][1]/100., y0/100., decimal=1)
assert_almost_equal(result[0][2]/100., b/100., decimal=1)
assert_almost_equal(result[0][3]/100., a/100., decimal=1)
assert_almost_equal(result[0][4], angle, decimal=1)
if __name__ == "__main__":
run_module_suite()