imsegm.ellipse_fitting module¶
Framework for ellipse fitting
Copyright (C) 2014-2018 Jiri Borovec <jiri.borovec@fel.cvut.cz>
- class imsegm.ellipse_fitting.EllipseModelSegm(*args: Any, **kwargs: Any)[source]¶
Bases:
skimage.measure.fit.EllipseModel
Total least squares estimator for 2D ellipses.
The functional model of the ellipse is:
xt = xc + a*cos(theta)*cos(t) - b*sin(theta)*sin(t) yt = yc + a*sin(theta)*cos(t) + b*cos(theta)*sin(t) d = sqrt((x - xt)**2 + (y - yt)**2)
where
(xt, yt)
is the closest point on the ellipse to(x, y)
. Thus d is the shortest distance from the point to the ellipse.The estimator is based on a least squares minimization. The optimal solution is computed directly, no iterations are required. This leads to a simple, stable and robust fitting method.
The
params
attribute contains the parameters in the following order:xc, yc, a, b, theta
Example
>>> from imsegm.utilities.drawing import ellipse_perimeter >>> params = 20, 30, 12, 16, np.deg2rad(30) >>> rr, cc = ellipse_perimeter(*params) >>> xy = np.array([rr, cc]).T >>> ellipse = EllipseModelSegm() >>> ellipse.estimate(xy) True >>> np.round(ellipse.params, 2) array([ 19.5 , 29.5 , 12.45, 16.52, 0.53]) >>> xy = EllipseModelSegm().predict_xy(np.linspace(0, 2 * np.pi, 25), params) >>> ellipse = EllipseModelSegm() >>> ellipse.estimate(xy) True >>> np.round(ellipse.params, 2) array([ 20. , 30. , 12. , 16. , 0.52]) >>> np.round(abs(ellipse.residuals(xy)), 5) array([ 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.]) >>> ellipse.params[2] += 2 >>> ellipse.params[3] += 2 >>> np.round(abs(ellipse.residuals(xy))) array([ 2., 2., 2., 2., 2., 2., 2., 2., 2., 2., 2., 2., 2., 2., 2., 2., 2., 2., 2., 2., 2., 2., 2., 2., 2.])
- criterion(points, weights, labels, table_prob=(0.1, 0.9))[source]¶
Determine residuals of data to model.
- Parameters
points – points coordinates
weights – weight for each point represent the region size
labels – vector of labels for each point
table_prob – is a vector or foreground probabilities for each class and being background is supplement to 1. Another option is define a matrix with number of columns related to number of classes and the first row denote probability being foreground and second being background
- Returns
Example
>>> seg = np.zeros((10, 15), dtype=int) >>> r, c = np.meshgrid(range(seg.shape[1]), range(seg.shape[0])) >>> el = EllipseModelSegm() >>> el.params = [4, 7, 3, 6, np.deg2rad(10)] >>> weights = np.ones(seg.ravel().shape) >>> seg[4:5, 6:8] = 1 >>> table_prob = [[0.1, 0.9]] >>> el.criterion(np.array([r.ravel(), c.ravel()]).T, weights, seg.ravel(), table_prob) 87.888... >>> seg[2:7, 4:11] = 1 >>> el.criterion(np.array([r.ravel(), c.ravel()]).T, weights, seg.ravel(), table_prob) 17.577... >>> seg[1:9, 1:14] = 1 >>> el.criterion(np.array([r.ravel(), c.ravel()]).T, weights, seg.ravel(), table_prob) -70.311...
- imsegm.ellipse_fitting.add_overlap_ellipse(segm, ellipse_params, label, thr_overlap=1.0)[source]¶
add to existing image ellipse with specific label if the new ellipse does not ouvelap with already existing object / ellipse
- Parameters
- Return ndarray
>>> seg = np.zeros((15, 20), dtype=int) >>> ell_params = 7, 10, 5, 8, np.deg2rad(30) >>> ell = add_overlap_ellipse(seg, ell_params, 1) >>> ell array([[0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0], [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0], [0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0], [0, 0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0], [0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0], [0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0], [0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0], [0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0], [0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0], [0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0, 0], [0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0], [0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 0], [0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0], [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0], [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]]) >>> ell2_params = 4, 5, 2, 3, np.deg2rad(-30) >>> add_overlap_ellipse(ell, ell2_params, 2) array([[0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0], [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0], [0, 0, 0, 2, 2, 2, 0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0], [0, 0, 0, 2, 2, 2, 2, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0], [0, 0, 0, 2, 2, 2, 2, 2, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0], [0, 0, 0, 0, 2, 2, 2, 2, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0], [0, 0, 0, 0, 1, 2, 2, 2, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0], [0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0], [0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0], [0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0, 0], [0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0], [0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 0], [0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0], [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0], [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0]])
- imsegm.ellipse_fitting.get_slic_points_labels(segm, img=None, slic_size=20, slic_regul=0.1)[source]¶
run SLIC on image or supepixels and return superpixels, their centers and also lebels (label from segmentation in position of superpixel centre)
- imsegm.ellipse_fitting.prepare_boundary_points_close(seg, centers, sp_size=25, relative_compact=0.3)[source]¶
extract some point around foreground boundaries
- Parameters
- Return list(ndarray)
>>> seg = np.zeros((100, 200), dtype=int) >>> ell_params = 50, 100, 40, 60, np.deg2rad(30) >>> seg = add_overlap_ellipse(seg, ell_params, 1) >>> pts = prepare_boundary_points_close(seg, [(40, 90)]) >>> pts [array([[ 6, 85], [ 8, 150], ... [ 92, 118]])]
- imsegm.ellipse_fitting.prepare_boundary_points_ray_dist(seg, centers, close_points=1, sel_bg=15, sel_fg=5)[source]¶
extract some point around foreground boundaries
- Parameters
- Return list(ndarray)
>>> seg = np.zeros((10, 20), dtype=int) >>> ell_params = 5, 10, 4, 6, np.deg2rad(30) >>> seg = add_overlap_ellipse(seg, ell_params, 1) >>> pts = prepare_boundary_points_ray_dist(seg, [(4, 9)], 2, sel_bg=0, sel_fg=0) >>> np.round(pts, 2).tolist() [[[4.0, 16.0], [6.8, 15.0], [9.0, 5.5], [4.35, 5.0], [1.0, 6.9], [1.0, 9.26], [0.0, 11.31], [0.5, 14.0], [1.45, 16.0]]]
- imsegm.ellipse_fitting.prepare_boundary_points_ray_edge(seg, centers, close_points=5, min_diam=25.0, sel_bg=15, sel_fg=5)[source]¶
extract some point around foreground boundaries
- Parameters
- Return list(ndarray)
>>> seg = np.zeros((10, 20), dtype=int) >>> ell_params = 5, 10, 4, 6, np.deg2rad(30) >>> seg = add_overlap_ellipse(seg, ell_params, 1) >>> pts = prepare_boundary_points_ray_edge(seg, [(4, 9)], 2.5, 3, sel_bg=1, sel_fg=0) >>> np.round(pts).tolist() [[[4.0, 16.0], [7.0, 15.0], [9.0, 5.0], [4.0, 5.0], [1.0, 7.0], [0.0, 14.0]]]
- imsegm.ellipse_fitting.prepare_boundary_points_ray_join(seg, centers, close_points=5, min_diam=25.0, sel_bg=15, sel_fg=5)[source]¶
extract some point around foreground boundaries
- Parameters
- Return list(ndarray)
>>> seg = np.zeros((10, 20), dtype=int) >>> ell_params = 5, 10, 4, 6, np.deg2rad(30) >>> seg = add_overlap_ellipse(seg, ell_params, 1) >>> pts = prepare_boundary_points_ray_join(seg, [(4, 9)], 5., 3, sel_bg=1, sel_fg=0) >>> np.round(pts).tolist() [[[4.0, 16.0], [7.0, 10.0], [9.0, 5.0], [4.0, 16.0], [7.0, 10.0]]]
- imsegm.ellipse_fitting.prepare_boundary_points_ray_mean(seg, centers, close_points=5, min_diam=25.0, sel_bg=15, sel_fg=5)[source]¶
extract some point around foreground boundaries
- Parameters
- Return list(ndarray)
>>> seg = np.zeros((10, 20), dtype=int) >>> ell_params = 5, 10, 4, 6, np.deg2rad(30) >>> seg = add_overlap_ellipse(seg, ell_params, 1) >>> pts = prepare_boundary_points_ray_mean(seg, [(4, 9)], 2.5, 3, sel_bg=1, sel_fg=0) >>> np.round(pts).tolist() [[[4.0, 16.0], [7.0, 15.0], [9.0, 5.0], [4.0, 5.0], [1.0, 7.0], [0.0, 14.0]]]
- imsegm.ellipse_fitting.ransac_segm(points, model_class, points_all, weights, labels, table_prob, min_samples, residual_threshold=1, max_trials=100)[source]¶
Fit a model to points with the RANSAC (random sample consensus).
- Parameters
points ([list, tuple of] (N, D) array) – Data set to which the model is fitted, where N is the number of points points and D the dimensionality of the points. If the model class requires multiple input points arrays (e.g. source and destination coordinates of
skimage.transform.AffineTransform
), they can be optionally passed as tuple or list. Note, that in this case the functionsestimate(*points)
,residuals(*points)
,is_model_valid(model, *random_data)
andis_data_valid(*random_data)
must all take each points array as separate arguments.model_class (class) –
Object with the following object methods:
success = estimate(*points)
residuals(*points)
where success indicates whether the model estimation succeeded (True or None for success, False for failure).
points_all (list) –
weights (list) –
labels (list) –
table_prob (list) –
min_samples (int float) – The minimum number of points points to fit a model to.
residual_threshold (float) – Maximum distance for a points point to be classified as an inlier.
max_trials (int, optional) – Maximum number of iterations for random sample selection.
- Returns
model (object) – Best model with largest consensus set.
inliers ((N, ) array) – Boolean mask of inliers classified as
True
.
Examples
>>> seg = np.zeros((120, 150), dtype=int) >>> ell_params = 60, 75, 40, 65, np.deg2rad(30) >>> seg = add_overlap_ellipse(seg, ell_params, 1) >>> slic, points_all, labels = get_slic_points_labels(seg, slic_size=10, slic_regul=0.3) >>> points = prepare_boundary_points_ray_dist(seg, [(40, 90)], 2, sel_bg=1, sel_fg=0)[0] >>> table_prob = [[0.01, 0.75, 0.95, 0.9], [0.99, 0.25, 0.05, 0.1]] >>> weights = np.bincount(slic.ravel()) >>> ransac_model, _ = ransac_segm( ... points, EllipseModelSegm, points_all, weights, labels, table_prob, 0.6, 3, max_trials=15) >>> np.round(ransac_model.params[:4]).astype(int) array([60, 75, 41, 65]) >>> np.round(ransac_model.params[4], 1) 0.5
- imsegm.ellipse_fitting.split_segm_background_foreground(seg, sel_bg=15, sel_fg=5)[source]¶
smoothing segmentation with morphological operation
- Parameters
seg (ndarray) – input segmentation
sel_bg (int|float) – smoothing background with morphological operation
sel_fg (int) – smoothing foreground with morphological operation
- Return tuple(ndarray,ndarray)
>>> seg = np.zeros((10, 20), dtype=int) >>> ell_params = 5, 10, 4, 6, np.deg2rad(30) >>> seg = add_overlap_ellipse(seg, ell_params, 1) >>> seg_bg, seg_fc = split_segm_background_foreground(seg, 1.5, 0) >>> seg_bg.astype(int) array([[1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1], [1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1], [1, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 1, 1], [1, 1, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 1, 1], [1, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 1, 1], [1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 1, 1], [1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1], [1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 1], [1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 1], [1, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1]]) >>> seg_fc.astype(int) array([[0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0], [0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0], [0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0, 0], [0, 0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0, 0], [0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0, 0], [0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0, 0], [0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0], [0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 0], [0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 0, 0], [0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0]])