[Octave-cvsupdate] octave-forge/main/image poly2mask.m,NONE,1.1 INDEX,1.8,1.9

[Josep M. i T. <jm...@us...>]

Update of /cvsroot/octave/octave-forge/main/image
In directory sc8-pr-cvs1.sourceforge.net:/tmp/cvs-serv9679

Modified Files:
	INDEX 
Added Files:
	poly2mask.m 
Log Message:
poly2mask added: creates filled polygon bw mask

--- NEW FILE: poly2mask.m ---
## Copyright (C) 2004 Josep Mones i Teixidor
##
## This program is free software; you can redistribute it and/or modify
## it under the terms of the GNU General Public License as published by
## the Free Software Foundation; either version 2 of the License, or
## (at your option) any later version.
##
## This program is distributed in the hope that it will be useful,
## but WITHOUT ANY WARRANTY; without even the implied warranty of
## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
## GNU General Public License for more details.
##
## You should have received a copy of the GNU General Public License
## along with this program; if not, write to the Free Software
## Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA

## -*- texinfo -*-
## @deftypefn {Function File} {@var{BW} = } poly2mask (@var{x},@var{y},@var{m},@var{n})
## Convert a polygon to a region mask
##
## BW=poly2mask(x,y,m,n) converts a polygon, specified by a list of
## vertices in @var{x} and @var{y} and returns in a @var{m}-by-@var{n}
## logical mask @var{BW} the filled polygon. Region inside the polygon
## is set to 1, values outside the shape are set to 0.
##
## @var{x} and @var{y} should always represent a closed polygon, first
## and last points should be coincident. If they are not poly2mask will
## close it for you. If @var{x} or @var{y} are fractional they are
## nearest integer.
##
## If all the polygon or part of it falls outside the masking area
## (1:m,1:n), it is discarded or clipped.
##
## This function uses scan-line polygon filling algorithm as described
## in http://www.cs.rit.edu/~icss571/ with some minor modifications:
## capability of clipping and scan order, which can affect the results
## of the algorithm (algorithm is described not to reach ymax, xmax
## border when filling to avoid enlarging shapes). In this function we
## scan the image backwards (we begin at ymax and end at ymin), and we
## don't reach ymin, xmin, which we believe should be compatibile with
## MATLAB.
## @end deftypefn


## TODO: convert BW to logical at creation time as soon as bool matrix
## TODO: arrays handle correctly assigning a bool (2.1.58?)

## Author:  Josep Mones i Teixidor <jm...@pu...>

function BW = poly2mask(x, y, m, n)
  if (nargin!=4)
    usage("BW = poly2mask(x, y, m, n)");
  endif

  ## check x and y
  x=round(x(:)');
  y=round(y(:)');
  if(length(x)<3)
    error("poly2mask: polygon must have at least 3 vertices.");
  endif
  if(length(x)!=length(y))
    error("poly2mask: length of x doesn't match length of y.");
  endif

  ## create output matrix
  BW=zeros(m,n);

  ## close polygon if needed
  if((x(1)!=x(length(x)))||(y(1)!=y(length(y))))
    x=horzcat(x,x(1));
    y=horzcat(y,y(1));
  endif

  ## build global edge table
  ex=[x(1:length(x)-1);x(1,2:length(x))]; ## x values for each edge
  ey=[y(1:length(y)-1);y(1,2:length(y))]; ## y values for each edge
  idx=find(ey(1,:)!=ey(2,:));             ## eliminate horizontal edges
  ex=ex(:,idx);
  ey=ey(:,idx);
  eminy=min(ey);                          ## minimum y for each edge
  emaxy=max(ey);                          ## maximum y for each edge
  t=ey==[eminy;eminy];                    ## values associated to miny
  exminy=ex(:)(find(t));                  ## x values associated to min y
  exmaxy=ex(:)(find(!t));                 ## x values associated to max y
  emaxy=emaxy';                           ## we want them vertical now...
  eminy=eminy';
  m_inv=(exmaxy-exminy)./(emaxy-eminy);   ## calculate inverse slope
  ge=[emaxy, eminy, exmaxy, m_inv];       ## build global edge table
  ge=sortrows(ge,[1,3]);                  ## sort on eminy and exminy
  
  ## we add an extra dummy edge at the end just to avoid checking
  ## while indexing it
  ge=[-Inf,-Inf,-Inf,-Inf;ge];
  
  ## initial parity is even (0)
  parity=0;
  
  ## init scan line set to bottom line
  sl=ge(size(ge,1),1);

  ## init active edge table
  ## we use a loop because the table is sorted and edge list could be
  ## huge
  ae=[];
  gei=size(ge,1);
  while sl==ge(gei,1)
    ae=[ge(gei,2:4);ae];
    gei-=1;
  endwhile

  ## calc minimum y to draw
  miny=min(y);
  if (miny<1)
    miny=1;
  endif

  while sl>=miny
    ## check vert clipping
    if(sl<=m)
      ## draw current scan line
      ## we have to round because 1/m is fractional
      ie=round(reshape(ae(:,2),2,size(ae)/2));

      ## this discards left border of image (this differs from version at
      ## http://www.cs.rit.edu/~icss571/filling which discards right
      ## border) but keeps an exception when the point is a vertex.
      ie(1,:)+=ie(1,:)!=ie(2,:);

      ## we'll clip too, just in case m,n is not big enough
      ie(1,find(ie(1,:)<1))=1;
      ie(2,find(ie(2,:)>n))=n;

      ## we eliminate segments outside window
      ie=ie(:,find(ie(1,:)<=n));
      ie=ie(:,find(ie(2,:)>=1));

      for i=1:columns(ie)
	BW(sl,ie(1,i):ie(2,i))=1;
      endfor
    endif
    
    ## decrement scan line
    sl-=1;

    ## eliminate edges that eymax==sl
    ## this discards ymin border of image (this differs from version at
    ## http://www.cs.rit.edu/~icss571/filling which discards ymax).
    ae=ae(find(ae(:,1)!=sl),:);
    
    ## update x (x1=x0-1/m)
    ae(:,2)-=ae(:,3);

    ## update ae with new values
    while sl==ge(gei,1)
      ae=vertcat(ae,ge(gei,2:4));
      gei-=1;
    endwhile
    
    ## order the edges in ae by x value
    if(rows(ae)>0)
      ae=sortrows(ae,2);
    endif
  endwhile
  
  ## In the future, we'll work directly on a logical array.
  BW=logical(BW);
endfunction


%!demo
%! s=[0:pi/4:2*pi];
%! x=cos(s)*90+101;
%! y=sin(s)*90+101;
%! bw=poly2mask(x,y,200,200);
%! imshow(bw);
%! %Creates a filled octagon

%!demo
%! s=[0:2*pi/5:pi*4];
%! s=s([1,3,5,2,4,6]);
%! x=cos(s)*90+101;
%! y=sin(s)*90+101;
%! bw=poly2mask(x,y,200,200);
%! imshow(bw);
%! %Creates a 5-vertex star

%!# Convex polygons

%!shared xs, ys, Rs, xt, yt, Rt
%! xs=[3,3,10,10];
%! ys=[4,12,12,4];
%! Rs=zeros(16,14);
%! Rs(5:12,4:10)=1;
%! Rs=logical(Rs);
%! xt=[1,4,7];
%! yt=[1,4,1];
%! Rt=[0,0,0,0,0,0,0;
%!     0,0,1,1,1,1,0;
%!     0,0,0,1,1,0,0;
%!     0,0,0,1,0,0,0;
%!     0,0,0,0,0,0,0];
%! Rt=logical(Rt);

%!assert(poly2mask(xs,ys,16,14),Rs);          # rectangle
%!assert(poly2mask(xs,ys,8,7),Rs(1:8,1:7));   # clipped
%!assert(poly2mask(xs-7,ys-8,8,7),Rs(9:16,8:14)); # more clipping

%!assert(poly2mask(xt,yt,5,7),Rt);            # triangle
%!assert(poly2mask(xt,yt,3,3),Rt(1:3,1:3));   # clipped


%!# Concave polygons

%!test
%! x=[3,3,5,5,8,8,10,10];
%! y=[4,12,12,8,8,11,11,4];
%! R=zeros(16,14);
%! R(5:12,4:5)=1;
%! R(5:8,6:8)=1;
%! R(5:11,9:10)=1;
%! R=logical(R);
%! assert(poly2mask(x,y,16,14), R);

%!# Complex polygons
%!test
%! x=[1,5,1,5];
%! y=[1,1,4,4];
%! R=[0,0,0,0,0,0;
%!    0,0,1,1,0,0;
%!    0,0,1,1,0,0;
%!    0,1,1,1,1,0;
%!    0,0,0,0,0,0];
%! R=logical(R);
%! assert(poly2mask(x,y,5,6), R);





%
% $Log: poly2mask.m,v $
% Revision 1.1  2004/08/11 17:34:11  jmones
% poly2mask added: creates filled polygon bw mask
%
%
Index: INDEX
===================================================================
RCS file: /cvsroot/octave/octave-forge/main/image/INDEX,v
retrieving revision 1.8
retrieving revision 1.9
diff -u -d -r1.8 -r1.9
--- INDEX	11 Aug 2004 00:57:02 -0000	1.8
+++ INDEX	11 Aug 2004 17:34:11 -0000	1.9
@@ -95,4 +95,5 @@
  roicolor
  roifill
  roifilt2
- roipoly
\ No newline at end of file
+ roipoly
+ poly2mask