Mercurial > hg > octave-nkf
changeset 15393:b99c52303d0b
voronoi.m: Fix bug when there are multiple identical input points (bug #37270)
Also clean up coding style to reflect core Octave style.
* voronoi.m: Eliminate empty facets from output of __voronoi__ call before
proceeding with the rest of processing.
author | Rik <rik@octave.org> |
---|---|
date | Sat, 15 Sep 2012 12:32:04 -0700 |
parents | c7fd43f5a89d |
children | 2b72343ed9c4 |
files | scripts/geometry/voronoi.m |
diffstat | 1 files changed, 14 insertions(+), 19 deletions(-) [+] |
line wrap: on
line diff
--- a/scripts/geometry/voronoi.m +++ b/scripts/geometry/voronoi.m @@ -105,10 +105,7 @@ linespec = varargin(narg); endif - lx = length (x); - ly = length (y); - - if (lx != ly) + if (length (x) != length (y)) error ("voronoi: X and Y must be vectors of the same length"); endif @@ -133,35 +130,33 @@ [[x(:) ; xbox(:)], [y(:); ybox(:)]], opts{:}); - idx = find (! infi); - ll = length (idx); - c = c(idx).'; - k = sum (cellfun ("length", c)); + c = c(! infi).'; + ## Delete null entries which cause problems in next cellfun function + c(cellfun ("isempty", c)) = []; edges = cell2mat (cellfun (@(x) [x ; [x(end), x(1:end-1)]], c, "uniformoutput", false)); ## Identify the unique edges of the Voronoi diagram edges = sortrows (sort (edges).').'; - edges = edges (:, [(edges(1, 1: end - 1) != edges(1, 2 : end) | ... - edges(2, 1 :end - 1) != edges(2, 2 : end)), true]); + edges = edges(:, [(edges(1, 1 :end - 1) != edges(1, 2 : end) | ... + edges(2, 1 :end - 1) != edges(2, 2 : end)), true]); ## Eliminate the edges of the diagram representing the box - poutside = (1 : rows (p)) ... - (p (:, 1) < xmin - xdelta | p (:, 1) > xmax + xdelta | ... - p (:, 2) < ymin - ydelta | p (:, 2) > ymax + ydelta); - edgeoutside = ismember (edges (1, :), poutside) & ... - ismember (edges (2, :), poutside); - edges (:, edgeoutside) = []; + poutside = (1:rows (p)) ... + (p(:, 1) < xmin - xdelta | p(:, 1) > xmax + xdelta | ... + p(:, 2) < ymin - ydelta | p(:, 2) > ymax + ydelta); + edgeoutside = ismember (edges(1, :), poutside) & ... + ismember (edges(2, :), poutside); + edges(:, edgeoutside) = []; ## Get points of the diagram Vvx = reshape (p(edges, 1), size (edges)); Vvy = reshape (p(edges, 2), size (edges)); if (nargout < 2) + h = plot (handl, Vvx, Vvy, linespec{:}, x, y, '+'); lim = [xmin, xmax, ymin, ymax]; - h = plot (handl, Vvx, Vvy, linespec{:}, x, y, '+'); - axis (lim + 0.1 * [[-1, 1] * (lim (2) - lim (1)), ... - [-1, 1] * (lim (4) - lim (3))]); + axis (lim + 0.1 * [[-1, 1] * xdelta, [-1, 1] * ydelta]); if (nargout == 1) vx = h; endif