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

Merge pull request #865 from jni/profile-line

Browse Source 
Move `profile_line()` out of viewer and refactor
This commit is contained in:
Stefan van der Walt
2014-01-30 11:39:42 -08:00
parent 305ef0dde0 01967e5060
commit 73bf701b5c
4 changed files with 228 additions and 79 deletions
Download Patch File Download Diff File Expand all files Collapse all files
skimage/measure/__init__.py
+3 -1
View File
@@ -4,6 +4,7 @@ from ._regionprops import regionprops, perimeter
from ._structural_similarity import structural_similarity
from ._polygon import approximate_polygon, subdivide_polygon
from ._moments import moments, moments_central, moments_normalized, moments_hu
from .profile import profile_line
from .fit import LineModel, CircleModel, EllipseModel, ransac
from .block import block_reduce
@@ -24,4 +25,5 @@ __all__ = ['find_contours',
'moments_normalized',
'moments_hu',
'marching_cubes',
'mesh_surface_area']
'mesh_surface_area',
'profile_line']
skimage/measure/profile.py
+87
View File
@@ -0,0 +1,87 @@
import numpy as np
import scipy.ndimage as nd
def profile_line(img, src, dst, linewidth=1,
order=1, mode='constant', cval=0.0):
"""Return the intensity profile of an image measured along a scan line.
Parameters
----------
img : numeric array, shape (M, N[, C])
The image, either grayscale (2D array) or multichannel
(3D array, where the final axis contains the channel
information).
src : 2-tuple of numeric scalar (float or int)
The start point of the scan line.
dst : 2-tuple of numeric scalar (float or int)
The end point of the scan line.
linewidth : int, optional
Width of the scan, perpendicular to the line
order : int in {0, 1, 2, 3, 4, 5}, optional
The order of the spline interpolation to compute image values at
non-integer coordinates. 0 means nearest-neighbor interpolation.
mode : string, one of {'constant', 'nearest', 'reflect', 'wrap'}, optional
How to compute any values falling outside of the image.
cval : float, optional
If `mode` is 'constant', what constant value to use outside the image.
Returns
-------
return_value : array
The intensity profile along the scan line. The length of the profile
is the ceil of the computed length of the scan line.
Examples
--------
>>> x = np.array([[1, 1, 1, 2, 2, 2]])
>>> img = np.vstack([np.zeros_like(x), x, x, x, np.zeros_like(x)])
>>> img
array([[0, 0, 0, 0, 0, 0],
[1, 1, 1, 2, 2, 2],
[1, 1, 1, 2, 2, 2],
[1, 1, 1, 2, 2, 2],
[0, 0, 0, 0, 0, 0]])
>>> profile_line(img, (2, 1), (2, 4))
array([ 1., 1., 2., 2.])
Notes
-----
The destination point is included in the profile, in contrast to
standard numpy indexing.
"""
src_row, src_col = src = np.asarray(src, dtype=float)
dst_row, dst_col = dst = np.asarray(dst, dtype=float)
d_row, d_col = dst - src
theta = np.arctan2(d_row, d_col)
length = np.ceil(np.hypot(d_row, d_col) + 1)
# we add one above because we include the last point in the profile
# (in contrast to standard numpy indexing)
line_col = np.linspace(src_col, dst_col, length)
line_row = np.linspace(src_row, dst_row, length)
# we subtract 1 from linewidth to change from pixel-counting
# (make this line 3 pixels wide) to point distances (the
# distance between pixel centers)
col_width = (linewidth - 1) * np.sin(-theta) / 2
row_width = (linewidth - 1) * np.cos(theta) / 2
perp_rows = np.array([np.linspace(row_i - row_width, row_i + row_width,
linewidth) for row_i in line_row])
perp_cols = np.array([np.linspace(col_i - col_width, col_i + col_width,
linewidth) for col_i in line_col])
perp_lines = np.array([perp_rows, perp_cols])
if img.ndim == 3:
pixels = [nd.map_coordinates(img[..., i], perp_lines,
order=order, mode=mode, cval=cval)
for i in range(img.shape[2])]
pixels = np.transpose(np.asarray(pixels), (1, 2, 0))
else:
pixels = nd.map_coordinates(img, perp_lines,
order=order, mode=mode, cval=cval)
intensities = pixels.mean(axis=1)
return intensities
skimage/measure/tests/test_profile.py
+110
View File
@@ -0,0 +1,110 @@
from numpy.testing import assert_equal, assert_almost_equal
import numpy as np
from skimage.measure import profile_line
image = np.arange(100).reshape((10, 10)).astype(np.float)
def test_horizontal_rightward():
prof = profile_line(image, (0, 2), (0, 8), order=0)
expected_prof = np.arange(2, 9)
assert_equal(prof, expected_prof)
def test_horizontal_leftward():
prof = profile_line(image, (0, 8), (0, 2), order=0)
expected_prof = np.arange(8, 1, -1)
assert_equal(prof, expected_prof)
def test_vertical_downward():
prof = profile_line(image, (2, 5), (8, 5), order=0)
expected_prof = np.arange(25, 95, 10)
assert_equal(prof, expected_prof)
def test_vertical_upward():
prof = profile_line(image, (8, 5), (2, 5), order=0)
expected_prof = np.arange(85, 15, -10)
assert_equal(prof, expected_prof)
def test_45deg_right_downward():
prof = profile_line(image, (2, 2), (8, 8), order=0)
expected_prof = np.array([22, 33, 33, 44, 55, 55, 66, 77, 77, 88])
# repeats are due to aliasing using nearest neighbor interpolation.
# to see this, imagine a diagonal line with markers every unit of
# length traversing a checkerboard pattern of squares also of unit
# length. Because the line is diagonal, sometimes more than one
# marker will fall on the same checkerboard box.
assert_almost_equal(prof, expected_prof)
def test_45deg_right_downward_interpolated():
prof = profile_line(image, (2, 2), (8, 8), order=1)
expected_prof = np.linspace(22, 88, 10)
assert_almost_equal(prof, expected_prof)
def test_45deg_right_upward():
prof = profile_line(image, (8, 2), (2, 8), order=1)
expected_prof = np.arange(82, 27, -6)
assert_almost_equal(prof, expected_prof)
def test_45deg_left_upward():
prof = profile_line(image, (8, 8), (2, 2), order=1)
expected_prof = np.arange(88, 21, -22. / 3)
assert_almost_equal(prof, expected_prof)
def test_45deg_left_downward():
prof = profile_line(image, (2, 8), (8, 2), order=1)
expected_prof = np.arange(28, 83, 6)
assert_almost_equal(prof, expected_prof)
def test_pythagorean_triangle_right_downward():
prof = profile_line(image, (1, 1), (7, 9), order=0)
expected_prof = np.array([11, 22, 23, 33, 34, 45, 56, 57, 67, 68, 79])
assert_equal(prof, expected_prof)
def test_pythagorean_triangle_right_downward_interpolated():
prof = profile_line(image, (1, 1), (7, 9), order=1)
expected_prof = np.linspace(11, 79, 11)
assert_almost_equal(prof, expected_prof)
pyth_image = np.zeros((6, 7), np.float)
line = ((1, 2, 2, 3, 3, 4), (1, 2, 3, 3, 4, 5))
below = ((2, 2, 3, 4, 4, 5), (0, 1, 2, 3, 4, 4))
above = ((0, 1, 1, 2, 3, 3), (2, 2, 3, 4, 5, 6))
pyth_image[line] = 1.8
pyth_image[below] = 0.6
pyth_image[above] = 0.6
def test_pythagorean_triangle_right_downward_linewidth():
prof = profile_line(pyth_image, (1, 1), (4, 5), linewidth=3, order=0)
expected_prof = np.ones(6)
assert_almost_equal(prof, expected_prof)
def test_pythagorean_triangle_right_upward_linewidth():
prof = profile_line(pyth_image[::-1, :], (4, 1), (1, 5),
linewidth=3, order=0)
expected_prof = np.ones(6)
assert_almost_equal(prof, expected_prof)
def test_pythagorean_triangle_transpose_left_down_linewidth():
prof = profile_line(pyth_image.T[:, ::-1], (1, 4), (5, 1),
linewidth=3, order=0)
expected_prof = np.ones(6)
assert_almost_equal(prof, expected_prof)
if __name__ == "__main__":
from numpy.testing import run_module_suite
run_module_suite()
skimage/viewer/plugins/lineprofile.py
+28 -78
View File
@@ -1,8 +1,9 @@
import warnings
import numpy as np
import scipy.ndimage as ndi
from skimage.util.dtype import dtype_range
from skimage import draw
from skimage import measure
from .plotplugin import PlotPlugin
from ..canvastools import ThickLineTool
@@ -70,7 +71,11 @@ class LineProfile(PlotPlugin):
on_change=self.line_changed)
self.line_tool.end_points = np.transpose([x, y])
scan_data = profile_line(image, self.line_tool.end_points)
scan_data = measure.profile_line(image,
*self.line_tool.end_points[:, ::-1])
self.scan_data = scan_data
if scan_data.ndim == 1:
scan_data = scan_data[:, np.newaxis]
self.reset_axes(scan_data)
@@ -104,8 +109,12 @@ class LineProfile(PlotPlugin):
def line_changed(self, end_points):
x, y = np.transpose(end_points)
self.line_tool.end_points = end_points
scan = profile_line(self.image_viewer.original_image, end_points,
linewidth=self.line_tool.linewidth)
scan = measure.profile_line(self.image_viewer.original_image,
*end_points[:, ::-1],
linewidth=self.line_tool.linewidth)
self.scan_data = scan
if scan.ndim == 1:
scan = scan[:, np.newaxis]
if scan.shape[1] != len(self.profile):
self.reset_axes(scan)
@@ -131,79 +140,20 @@ class LineProfile(PlotPlugin):
scan_data[:, 1], 'g-',
scan_data[:, 2], 'b-')
def output(self):
"""Return the drawn line and the resulting scan.
def _calc_vert(img, x1, x2, y1, y2, linewidth):
# Quick calculation if perfectly horizontal
pixels = img[min(y1, y2): max(y1, y2) + 1,
x1 - linewidth / 2: x1 + linewidth / 2 + 1]
Returns
-------
line_image : (M, N) uint8 array, same shape as image
An array of 0s with the scanned line set to 255.
scan : (P,) or (P, 3) array of int or float
The line scan values across the image.
"""
(x1, y1), (x2, y2) = self.line_tool.end_points
line_image = np.zeros(self.image_viewer.original_image.shape[:2],
np.uint8)
rr, cc = draw.line(y1, x1, y2, x2)
line_image[rr, cc] = 255
return line_image, self.scan_data
# Reverse index if necessary
if y2 > y1:
pixels = pixels[::-1, :]
return pixels.mean(axis=1)[:, np.newaxis]
def profile_line(img, end_points, linewidth=1):
"""Return the intensity profile of an image measured along a scan line.
Parameters
----------
img : 2d or 3d array
The image, in grayscale (2d) or RGB (3d) format.
end_points: (2, 2) list
End points ((x1, y1), (x2, y2)) of scan line.
linewidth: int
Width of the scan, perpendicular to the line
Returns
-------
return_value : array
The intensity profile along the scan line. The length of the profile
is the ceil of the computed length of the scan line.
"""
point1, point2 = end_points
x1, y1 = point1 = np.asarray(point1, dtype=float)
x2, y2 = point2 = np.asarray(point2, dtype=float)
dx, dy = point2 - point1
channels = 1
if img.ndim == 3:
channels = 3
# Quick calculation if perfectly vertical; shortcuts div0 error
if x1 == x2:
if channels == 1:
img = img[:, :, np.newaxis]
img = np.rollaxis(img, -1)
intensities = np.hstack([_calc_vert(im, x1, x2, y1, y2, linewidth)
for im in img])
return intensities
theta = np.arctan2(dy, dx)
a = dy / dx
b = y1 - a * x1
length = np.hypot(dx, dy)
line_x = np.linspace(x2, x1, np.ceil(length))
line_y = line_x * a + b
y_width = abs(linewidth * np.cos(theta) / 2)
perp_ys = np.array([np.linspace(yi - y_width,
yi + y_width, linewidth) for yi in line_y])
perp_xs = - a * perp_ys + (line_x + a * line_y)[:, np.newaxis]
perp_lines = np.array([perp_ys, perp_xs])
if img.ndim == 3:
pixels = [ndi.map_coordinates(img[..., i], perp_lines)
for i in range(3)]
pixels = np.transpose(np.asarray(pixels), (1, 2, 0))
else:
pixels = ndi.map_coordinates(img, perp_lines)
pixels = pixels[..., np.newaxis]
intensities = pixels.mean(axis=1)
if intensities.ndim == 1:
return intensities[..., np.newaxis]
else:
return intensities