octave_packages/geometry-1.5.0/shape2d/shapearea.m

   1 %% Copyright (c) 2011 Juan Pablo Carbajal <carbajal@ifi.uzh.ch>
   2 %%
   3 %% This program is free software; you can redistribute it and/or modify
   4 %% it under the terms of the GNU General Public License as published by
   5 %% the Free Software Foundation; either version 3 of the License, or
   6 %% (at your option) any later version.
   7 %%
   8 %% This program is distributed in the hope that it will be useful,
   9 %% but WITHOUT ANY WARRANTY; without even the implied warranty of
  10 %% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
  11 %% GNU General Public License for more details.
  12 %%
  13 %% You should have received a copy of the GNU General Public License
  14 %% along with this program; if not, see <http://www.gnu.org/licenses/>.
  15
  16 %% -*- texinfo -*-
  17 %% @deftypefn {Function File} { @var{a} =} shapearea (@var{pp})
  18 %% Calculate the area of a 2D shape defined with piecewise smooth polynomials.
  19 %%
  20 %% Shape is defined with piecewise smooth polynomials. @var{pp} is a
  21 %% cell where each elements is a 2-by-(poly_degree+1) array containing a pair of
  22 %% polynomials.
  23 %%
  24 %% @code{px(i,:) = pp@{i@}(1,:)} and @code{py(i,:) = pp@{i@}(2,:)}.
  25 %%
  26 %% @seealso{shapecentroid, shape2polygon, shapeplot}
  27 %% @end deftypefn
  28
  29 function [A ccw] = shapearea (shape)
  30
  31   A = sum(cellfun (@Aint, shape));
  32   if A < 0
  33     warning ('geom2d:shapearea:InvalidResult', ...
  34     'Shape has negative area. Assuming this is due to a clockwise parametrization of the boundary');
  35     A = -A;
  36   end
  37
  38 endfunction
  39
  40 function dA = Aint (x)
  41
  42     px = x(1,:);
  43     py = x(2,:);
  44
  45     P = polyint (conv (px, polyder(py)));
  46
  47     dA = diff(polyval(P,[0 1]));
  48
  49 end
  50
  51 %!demo % non-convex piece-wise polynomial shape
  52 %! boomerang = {[ 0 -2 1; ...
  53 %!               -4  4 0]; ...
  54 %!              [0.25 -1; ...
  55 %!               0     0]; ...
  56 %!              [ 0 1.5 -0.75; ...
  57 %!               -3 3    0];
  58 %!              [0.25 0.75; ...
  59 %!               0 0]};
  60 %! A = shapearea (boomerang)
  61
  62 %!test
  63 %! triangle = {[1 0; 0 0]; [-0.5 1; 1 0]; [-0.5 0.5; -1 1]};
  64 %! A = shapearea (triangle);
  65 %! assert (0.5, A);
  66
  67 %!test
  68 %! circle = {[1.715729  -6.715729    0   5; ...
  69 %!            -1.715729  -1.568542   8.284271    0]; ...
  70 %!            [1.715729   1.568542  -8.284271    0; ...
  71 %!             1.715729  -6.715729    0   5]; ...
  72 %!            [-1.715729   6.715729    0  -5; ...
  73 %!             1.715729   1.568542  -8.284271    0]; ...
  74 %!            [-1.715729  -1.568542   8.284271    0; ...
  75 %!            -1.715729   6.715729    0  -5]};
  76 %! A = shapearea (circle);
  77 %! assert (pi*5^2, A, 5e-2);