changeset 9820:e65110da131a octave-forge

geometry: Fixed bugs detected by Rafael (dev packager). Improved docstring of cbezier2poly. Verified all other shape2d functions. Moved test form shapetransfrom to shapecentriod
author jpicarbajal
date Thu, 22 Mar 2012 22:05:06 +0000
parents 537f5b239baa
children c869fff01da2
files main/geometry/NEWS main/geometry/inst/geom2d/cbezier2poly.m main/geometry/inst/shape2d/shape2polygon.m main/geometry/inst/shape2d/shapearea.m main/geometry/inst/shape2d/shapecentroid.m main/geometry/inst/shape2d/shapetransform.m
diffstat 6 files changed, 98 insertions(+), 71 deletions(-) [+]
line wrap: on
line diff
--- a/main/geometry/NEWS	Thu Mar 22 17:54:01 2012 +0000
+++ b/main/geometry/NEWS	Thu Mar 22 22:05:06 2012 +0000
@@ -14,11 +14,12 @@
  - @svg/path2polygon.m
  - Fix addpath/rmpath installation warnings
  - Fix octclip/src/Makefile
+ - Fix shapecentriod.m for piece-wise polynomial shapes.
 
 * Known issues
  - simplifypolygon.m returns empty polygons when points are repeated, i.e when
   the polygon is not correctly formed.
- - shapecentriod.m gives wierd results for Bezier curves.
+
 
 ===============================================================================
 geometry-1.4.0   Release Date: 2012-01-25   Release Manager: Juan Pablo Carbajal
--- a/main/geometry/inst/geom2d/cbezier2poly.m	Thu Mar 22 17:54:01 2012 +0000
+++ b/main/geometry/inst/geom2d/cbezier2poly.m	Thu Mar 22 22:05:06 2012 +0000
@@ -39,8 +39,8 @@
 %% With only one input argument, calculates the polynomial @var{pp} of the cubic
 %% Bezier curve defined by the 4 control points stored in @var{points}. The first
 %% point is the inital point of the curve. The segment joining the first point
-%% with the second point defines the tangent of the curve at the initial point.
-%% The segment that joints the third point with the fourth defines the tanget at
+%% with the second point (first center) defines the tangent of the curve at the initial point.
+%% The segment that joints the third point (second center) with the fourth defines the tanget at
 %% the end-point of the curve, which is defined in the fourth point.
 %% @var{points} is either a 4-by-2 array (vertical concatenation of point
 %% coordinates), or a 1-by-8 array (horizotnal concatenation of point
--- a/main/geometry/inst/shape2d/shape2polygon.m	Thu Mar 22 17:54:01 2012 +0000
+++ b/main/geometry/inst/shape2d/shape2polygon.m	Thu Mar 22 22:05:06 2012 +0000
@@ -29,21 +29,24 @@
 function polygon = shape2polygon (shape, N=16)
 
   polygon = cell2mat ( ...
-             cellfun(@(x) func (x,N), shape,'UniformOutput',false) );
+             cellfun (@(x) func (x,N), shape,'UniformOutput',false) );
 
-  if size(polygon, 1) == 1
-    polygon(2,1) = polyval(shape{1}(1,:),1);
-    polygon(2,2) = polyval(shape{1}(2,:),1);
+  %% TODO simply the polygon based on curvature
+%  polygon = unique (polygon, 'rows');
+
+  if size (polygon, 1) == 1
+    polygon(2,1) = polyval (shape{1}(1,:), 1);
+    polygon(2,2) = polyval (shape{1}(2,:), 1);
   end
 
 endfunction
 
-function y = func(x,N)
+function y = func (x,N)
 
-  if size(x,2) > 2
-    t = linspace(0,1-1/N,N).';
-    y(:,1) = polyval(x(1,:),t);
-    y(:,2) = polyval(x(2,:),t);
+  if size (x,2) > 2
+    t = linspace (0,1-1/N,N).';
+    y(:,1) = polyval (x(1,:), t);
+    y(:,2) = polyval (x(2,:), t);
   else
     y = x(:,2).';
   end
--- a/main/geometry/inst/shape2d/shapearea.m	Thu Mar 22 17:54:01 2012 +0000
+++ b/main/geometry/inst/shape2d/shapearea.m	Thu Mar 22 22:05:06 2012 +0000
@@ -26,10 +26,12 @@
 %% @seealso{shapecentroid, shape2polygon, shapeplot}
 %% @end deftypefn
 
-function A = shapearea (shape)
+function [A ccw] = shapearea (shape)
 
   A = sum(cellfun (@Aint, shape));
   if A < 0
+    warning ('geom2d:shapearea:InvalidResult', ...
+    'Shape has negative area. Assuming this is due to a clockwise parametrization of the boundary');
     A = -A;
   end
 
@@ -46,12 +48,16 @@
 
 end
 
-%!demo % non-convex bezier shape
-%! weirdhearth ={[-17.6816  -34.3989    7.8580    3.7971; ...
-%!                15.4585  -28.3820  -18.7645    9.8519]; ...
-%!                 [-27.7359   18.1039  -34.5718    3.7878; ...
-%!                  -40.7440   49.7999  -25.5011    2.2304]};
-%! A = shapearea (weirdhearth)
+%!demo % non-convex piece-wise polynomial shape
+%! boomerang = {[ 0 -2 1; ...
+%!               -4  4 0]; ...
+%!              [0.25 -1; ...
+%!               0     0]; ...
+%!              [ 0 1.5 -0.75; ...
+%!               -3 3    0];
+%!              [0.25 0.75; ...
+%!               0 0]};
+%! A = shapearea (boomerang)
 
 %!test
 %! triangle = {[1 0; 0 0]; [-0.5 1; 1 0]; [-0.5 0.5; -1 1]};
--- a/main/geometry/inst/shape2d/shapecentroid.m	Thu Mar 22 17:54:01 2012 +0000
+++ b/main/geometry/inst/shape2d/shapecentroid.m	Thu Mar 22 22:05:06 2012 +0000
@@ -15,30 +15,39 @@
 
 %% -*- texinfo -*-
 %% @deftypefn {Function File} { @var{cm} =} shapecentroid (@var{pp})
-%%  Centroid of a plane shape defined with piecewise smooth polynomials.
+%%  Centroid of a simple plane shape defined with piecewise smooth polynomials.
 %%
 %% The shape is defined with piecewise smooth polynomials. @var{pp} is a
 %% cell where each elements is a 2-by-(poly_degree+1) matrix containing a pair
 %% of polynomials.
 %% @code{px(i,:) = pp@{i@}(1,:)} and @code{py(i,:) = pp@{i@}(2,:)}.
 %%
+%% The edges of the shape should not self-intersect. This function does not check for the
+%% sanity of the shape.
+%%
 %% @seealso{shapearea, shape2polygon}
 %% @end deftypefn
 
 function cm = shapecentroid (shape)
 
-  cm = sum( cell2mat ( cellfun (@CMint, shape, 'UniformOutput', false)));
+  cm = sum( cell2mat ( cellfun (@CMint, shape, 'UniformOutput', false)), 1);
   A = shapearea(shape);
   cm = cm / A;
 
+  [~,id] = lastwarn ('','');
+  if strcmp (id ,'geom2d:shapearea:InvalidResult')
+    lastwarn('Inverting centriod','geom2d:shapecentriod:InvalidResult');
+    cm = -cm;
+  end
+
 endfunction
 
 function dcm = CMint (x)
 
     px = x(1,:);
     py = x(2,:);
-    Px = polyint (conv(conv (px , px)/2 , polyder (py)));
-    Py = polyint (conv(-conv (py , py)/2 , polyder (px)));
+    Px = polyint (conv(conv (-px , py) , polyder (px)));
+    Py = polyint (conv(conv (px , py) , polyder (py)));
 
     dcm = zeros (1,2);
     dcm(1) = diff(polyval(Px,[0 1]));
@@ -47,17 +56,22 @@
 endfunction
 
 %!demo % non-convex bezier shape
-%! weirdhearth ={[34.81947,-63.60585 41.35964,1.61093; ...
-%!                73.22086,4.95439 7.1796,-34.7948]; ...
-%!                 [30.26599,-50.0316 77.6279,8.52058; ...
-%!                  -18.66371,58.02699 -168.20415,52.74819]};
-%! CoM = shapecentroid (weirdhearth)
-%! Gcentriod = centroid(shape2polygon(weirdhearth))
+%! boomerang = {[ 0 -2 1; ...
+%!               -4  4 0]; ...
+%!              [0.25 -1; ...
+%!               0     0]; ...
+%!              [ 0 1.5 -0.75; ...
+%!               -3 3    0];
+%!              [0.25 0.75; ...
+%!               0 0]};
+%! CoM = shapecentroid (boomerang)
+%! Gcentriod = centroid(shape2polygon(boomerang))
 %!
-%! shapeplot(weirdhearth);
+%! figure(1); clf;
+%! shapeplot(boomerang,10,'-o');
 %! hold on
-%! drawPoint(CoM,'ok');
-%! drawPoint(Gcentriod,'or');
+%! drawPoint(CoM,'xk;shape centriod;');
+%! drawPoint(Gcentriod,'xr;point centriod;');
 %! hold off
 %! axis equal
 
@@ -71,10 +85,10 @@
 %! CoM = shapecentroid (Lshape)
 %! Gcentriod = centroid (shape2polygon (Lshape))
 %!
-%! shapeplot(Lshape);
+%! shapeplot(Lshape,10,'-o');
 %! hold on
-%! drawPoint(CoM,'ok');
-%! drawPoint(Gcentriod,'or');
+%! drawPoint(CoM,'xk;shape centriod;');
+%! drawPoint(Gcentriod,'xr;point centriod;');
 %! hold off
 %! axis equal
 
@@ -82,10 +96,12 @@
 %! square = {[1 -0.5; 0 -0.5]; [0 0.5; 1 -0.5]; [-1 0.5; 0 0.5]; [0 -0.5; -1 0.5]};
 %! CoM = shapecentroid (square);
 %! assert (CoM, [0 0], sqrt(eps));
+
+%!test
+%! square = {[1 -0.5; 0 -0.5]; [0 0.5; 1 -0.5]; [-1 0.5; 0 0.5]; [0 -0.5; -1 0.5]};
 %! square_t = shapetransform (square,[1;1]);
-%! CoM_t = shapecentroid (square_t);
-%! assert (CoM, [0 0], sqrt(eps));
-%! assert (CoM_t, [1 1], sqrt(eps));
+%! CoM = shapecentroid (square_t);
+%! assert (CoM, [1 1], sqrt(eps));
 
 %!test
 %! circle = {[1.715729  -6.715729    0   5; ...
@@ -98,3 +114,37 @@
 %!            -1.715729   6.715729    0  -5]};
 %! CoM = shapecentroid (circle);
 %! assert (CoM , [0 0], 5e-3);
+
+%!shared shape
+%! shape = {[-93.172   606.368  -476.054   291.429; ...
+%!          -431.196   637.253    11.085   163.791]; ...
+%!         [-75.3626  -253.2337   457.1678   328.5714; ...
+%!           438.7659  -653.6278    -7.9953   380.9336]; ...
+%!         [-89.5841   344.9716  -275.3876   457.1429; ...
+%!          -170.3613   237.8858     1.0469   158.0765];...
+%!         [32.900  -298.704   145.804   437.143; ...
+%!         -243.903   369.597   -34.265   226.648]; ...
+%!         [-99.081   409.127  -352.903   317.143; ...
+%!           55.289  -114.223   -26.781   318.076]; ...
+%!         [-342.231   191.266   168.108   274.286; ...
+%!           58.870   -38.083   -89.358   232.362]};
+
+%!test % x-Reflection
+%! v = shapecentroid (shape)(:);
+%! T = createLineReflection([0 0 1 0]);
+%! nshape = shapetransform (shape, T);
+%! vn = shapecentroid (nshape)(:);
+%! assert(vn,T(1:2,1:2)*v);
+
+%!test % Rotation
+%! v = shapecentroid (shape)(:);
+%! T = createRotation(v.',pi/2);
+%! nshape = shapetransform (shape, T);
+%! vn = shapecentroid (nshape)(:);
+%! assert(vn,v,1e-2);
+
+%!test % Translation
+%! v = shapecentroid (shape)(:);
+%! nshape = shapetransform (shape, -v);
+%! vn = shapecentroid (nshape)(:);
+%! assert(vn,[0; 0],1e-2);
--- a/main/geometry/inst/shape2d/shapetransform.m	Thu Mar 22 17:54:01 2012 +0000
+++ b/main/geometry/inst/shape2d/shapetransform.m	Thu Mar 22 22:05:06 2012 +0000
@@ -109,36 +109,3 @@
 %! axis square
 
 
-%!shared shape
-%! shape = {[-93.172   606.368  -476.054   291.429; ...
-%!          -431.196   637.253    11.085   163.791]; ...
-%!         [-75.3626  -253.2337   457.1678   328.5714; ...
-%!           438.7659  -653.6278    -7.9953   380.9336]; ...
-%!         [-89.5841   344.9716  -275.3876   457.1429; ...
-%!          -170.3613   237.8858     1.0469   158.0765];...
-%!         [32.900  -298.704   145.804   437.143; ...
-%!         -243.903   369.597   -34.265   226.648]; ...
-%!         [-99.081   409.127  -352.903   317.143; ...
-%!           55.289  -114.223   -26.781   318.076]; ...
-%!         [-342.231   191.266   168.108   274.286; ...
-%!           58.870   -38.083   -89.358   232.362]};
-
-%!test
-%! v = shapecentroid (shape)(:);
-%! nshape = shapetransform (shape, -v);
-%! vn = shapecentroid (nshape)(:);
-%! assert(vn,[0; 0],1e-2);
-
-%!test
-%! v = shapecentroid (shape)(:);
-%! T = createLineReflection([0 0 1 0]);
-%! nshape = shapetransform (shape, T);
-%! vn = shapecentroid (nshape)(:);
-%! assert(vn,T(1:2,1:2)*v);
-
-%!test
-%! v = shapecentroid (shape)(:);
-%! T = createRotation(v.',pi/2);
-%! nshape = shapetransform (shape, T);
-%! vn = shapecentroid (nshape)(:);
-%! assert(vn,v,1e-2);