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

MAINT: Refactor ellipsoid generator into skimage.draw

Browse Source
This commit is contained in:
Josh Warner (Mac)
2013-09-01 20:51:20 -05:00
parent d0d9fee36e
commit 14a0685838
3 changed files with 128 additions and 109 deletions
Download Patch File Download Diff File Expand all files Collapse all files
skimage/draw/__init__.py
+3
View File
@@ -1,11 +1,14 @@
from .draw import circle, ellipse, set_color
from ._draw import line, polygon, ellipse_perimeter, circle_perimeter, \
bezier_segment
from .draw3d import ellipsoid, ellipsoid_stats
__all__ = ['line',
'polygon',
'ellipse',
'ellipse_perimeter',
'ellipsoid',
'ellipsoid_stats',
'circle',
'circle_perimeter',
'set_color']
skimage/draw/draw3d.py
+117
View File
@@ -0,0 +1,117 @@
# coding: utf-8
import numpy as np
from scipy.special import (ellipkinc as ellip_F, ellipeinc as ellip_E)
def ellipsoid(a, b, c, sampling=(1., 1., 1.), levelset=False):
"""
Generates ellipsoid with semimajor axes aligned with grid dimensions
on grid with specified `sampling`.
Parameters
----------
a : float
Length of semimajor axis aligned with x-axis.
b : float
Length of semimajor axis aligned with y-axis.
c : float
Length of semimajor axis aligned with z-axis.
sampling : tuple of floats, length 3
Sampling in (x, y, z) spatial dimensions.
levelset : bool
If True, returns the level set for this ellipsoid (signed level
set about zero, with positive denoting interior) as np.float64.
False returns a binarized version of said level set.
Returns
-------
ellip : (N, M, P) array
Ellipsoid centered in a correctly sized array for given `sampling`.
Boolean dtype unless `levelset=True`, in which case a float array is
returned with the level set above 0.0 representing the ellipsoid.
"""
if (a <= 0) or (b <= 0) or (c <= 0):
raise ValueError('Parameters a, b, and c must all be > 0')
offset = np.r_[1, 1, 1] * np.r_[sampling]
# Calculate limits, and ensure output volume is odd & symmetric
low = np.ceil((- np.r_[a, b, c] - offset))
high = np.floor((np.r_[a, b, c] + offset + 1))
for dim in range(3):
if (high[dim] - low[dim]) % 2 == 0:
low[dim] -= 1
num = np.arange(low[dim], high[dim], sampling[dim])
if 0 not in num:
low[dim] -= np.max(num[num < 0])
# Generate (anisotropic) spatial grid
x, y, z = np.mgrid[low[0]:high[0]:sampling[0],
low[1]:high[1]:sampling[1],
low[2]:high[2]:sampling[2]]
if not levelset:
arr = ((x / float(a)) ** 2 +
(y / float(b)) ** 2 +
(z / float(c)) ** 2) <= 1
else:
arr = ((x / float(a)) ** 2 +
(y / float(b)) ** 2 +
(z / float(c)) ** 2) - 1
return arr
def ellipsoid_stats(a, b, c, sampling=(1., 1., 1.)):
"""
Calculates analytical surface area and volume for ellipsoid with
semimajor axes aligned with grid dimensions of specified `sampling`.
Parameters
----------
a : float
Length of semimajor axis aligned with x-axis.
b : float
Length of semimajor axis aligned with y-axis.
c : float
Length of semimajor axis aligned with z-axis.
sampling : tuple of floats, length 3
Sampling in (x, y, z) spatial dimensions.
Returns
-------
vol : float
Analytically calculated volume of ellipsoid.
surf : float
Analytically calculated surface area of ellipsoid.
"""
if (a <= 0) or (b <= 0) or (c <= 0):
raise ValueError('Parameters a, b, and c must all be > 0')
# Calculate volume & surface area
# Surface calculation requires a >= b >= c and a != c.
abc = [a, b, c]
abc.sort(reverse=True)
a = abc[0]
b = abc[1]
c = abc[2]
# Volume
vol = 4 / 3. * np.pi * a * b * c
# Analytical ellipsoid surface area
phi = np.arcsin((1. - (c ** 2 / (a ** 2.))) ** 0.5)
d = float((a ** 2 - c ** 2) ** 0.5)
m = (a ** 2 * (b ** 2 - c ** 2) /
float(b ** 2 * (a ** 2 - c ** 2)))
F = ellip_F(phi, m)
E = ellip_E(phi, m)
surf = 2 * np.pi * (c ** 2 +
b * c ** 2 / d * F +
b * d * E)
return vol, surf
skimage/measure/tests/test_marching_cubes.py
+8 -109
View File
@@ -1,114 +1,13 @@
import numpy as np
from numpy.testing import assert_raises
from scipy.special import (ellipkinc as ellip_F, ellipeinc as ellip_E)
from skimage.draw import ellipsoid, ellipsoid_stats
from skimage.measure import marching_cubes, mesh_surface_area
def _ellipsoid(a, b, c, sampling=(1., 1., 1.), info=False, tight=False,
levelset=False):
"""
Generates ellipsoid with semimajor axes aligned with grid dimensions,
on grid with specified `sampling`.
Parameters
----------
a : float
Length of semimajor axis aligned with x-axis
b : float
Length of semimajor axis aligned with y-axis
c : float
Length of semimajor axis aligned with z-axis
sampling : tuple of floats, length 3
Sampling in each spatial dimension
info : bool
If False, only `bool_arr` returned.
If True, (`bool_arr`, `vol`, `surf`) returned; the additional
values are analytical volume and surface area calculated for
this ellipsoid.
tight : bool
Controls if the ellipsoid will precisely be contained within
the returned volume (tight=True) or if each dimension will be
2 longer than necessary (tight=False). For algorithms which
need both sides of a contour, use False.
levelset : bool
If True, returns the level set for this ellipsoid (signed level
set about zero, with positive denoting interior) as np.float64.
False returns a binarized version of said level set.
Returns
-------
bool_arr : (N, M, P) array
Sphere in an appropriately sized boolean array.
vol : float
Analytically calculated volume of ellipsoid. Only returned if
`info` is True.
surf : float
Analytically calculated surface area of ellipsoid. Only returned
if `info` is True.
"""
if not tight:
offset = np.r_[1, 1, 1] * np.r_[sampling]
else:
offset = np.r_[0, 0, 0]
# Calculate limits, and ensure output volume is odd & symmetric
low = np.ceil((-np.r_[a, b, c] - offset))
high = np.floor((np.r_[a, b, c] + offset + 1))
for dim in range(3):
if (high[dim] - low[dim]) % 2 == 0:
low[dim] -= 1
num = np.arange(low[dim], high[dim], sampling[dim])
if 0 not in num:
low[dim] -= np.max(num[num < 0])
# Generate (anisotropic) spatial grid
x, y, z = np.mgrid[low[0]:high[0]:sampling[0],
low[1]:high[1]:sampling[1],
low[2]:high[2]:sampling[2]]
if not levelset:
arr = ((x / float(a)) ** 2 +
(y / float(b)) ** 2 +
(z / float(c)) ** 2) <= 1
else:
arr = ((x / float(a)) ** 2 +
(y / float(b)) ** 2 +
(z / float(c)) ** 2) - 1
if not info:
return arr
else:
# Surface calculation requires a >= b >= c and a != c.
abc = [a, b, c]
abc.sort(reverse=True)
a = abc[0]
b = abc[1]
c = abc[2]
# Volume
vol = 4 / 3. * np.pi * a * b * c
# Analytical ellipsoid surface area
phi = np.arcsin((1. - (c ** 2 / (a ** 2.))) ** 0.5)
d = float((a ** 2 - c ** 2) ** 0.5)
m = (a ** 2 * (b ** 2 - c ** 2) /
float(b ** 2 * (a ** 2 - c ** 2)))
F = ellip_F(phi, m)
E = ellip_E(phi, m)
surf = 2 * np.pi * (c ** 2 +
b * c ** 2 / d * F +
b * d * E)
return arr, vol, surf
def test_marching_cubes_isotropic():
ellipsoid_isotropic, _, surf = _ellipsoid(6, 10, 16,
levelset=True,
info=True)
ellipsoid_isotropic = ellipsoid(6, 10, 16, levelset=True)
_, surf = ellipsoid_stats(6, 10, 16, levelset=True)
verts, faces = marching_cubes(ellipsoid_isotropic, 0.)
surf_calc = mesh_surface_area(verts, faces)
@@ -118,13 +17,13 @@ def test_marching_cubes_isotropic():
def test_marching_cubes_anisotropic():
sampling = (1., 10 / 6., 16 / 6.)
ellipsoid_isotropic, _, surf = _ellipsoid(6, 10, 16,
sampling=sampling,
levelset=True,
info=True)
verts, faces = marching_cubes(ellipsoid_isotropic, 0.,
ellipsoid_anisotropic, _, surf = ellipsoid(6, 10, 16, sampling=sampling,
levelset=True)
_, surf = ellipsoid_stats(6, 10, 16, sampling=sampling, levelset=True)
verts, faces = marching_cubes(ellipsoid_anisotropic, 0.,
sampling=sampling)
surf_calc = mesh_surface_area(verts, faces)
# Test within 1.5% tolerance for anisotropic. Will always underestimate.
assert surf > surf_calc and surf_calc > surf * 0.985