过冷水带你设计闭区域多边面积计算算法

3年前浏览4576

过冷水最近遇到了这么一个问题，有一系列点组成了如上图所示的封闭图形，该如何求面积？

在过冷水印象中求面积=求积分，可是该问题不同于普通的积分问题。以往的是它不能用函数形式去表示啊！这可为难我胖虎了，在网上百度封闭MATLAB封闭图像的面积计算

有现成函数调用,于是就输入

S1=polyarea(x,y)
S1 =  
 4.6645e 03

轻松解决我的疑惑，之后有人问我这个求面积的方法靠谱吗？使用的是什么算法，我怎么知道使用什么算法，于是使用 open polyarea查看该函数的算法。

function area = polyarea(x,y,dim)
if nargin==1
  error(message('MATLAB:polyarea:NotEnoughInputs'));
end

if ~isequal(size(x),size(y))
  error(message('MATLAB:polyarea:XYSizeMismatch'));
end
if nargin==2
  [x,nshifts] = shiftdim(x);
  y = shiftdim(y);
elseif nargin==3
  if ~isscalar(dim) || ~isnumeric(dim) || (dim ~= floor(dim))
      error(message('MATLAB:getdimarg:dimensionMustBePositiveInteger'));
  end
  % Preserve existing errors for non-integer dim.
  dim = min(ndims(y) 1, dim);
  perm = [dim:max(length(size(x)),dim) 1:dim-1];
  x = permute(x,perm);
  y = permute(y,perm);
end
siz = size(x);
if ~isempty(x)
  area = reshape(abs(sum( (x([2:siz(1) 1],:) - x(:,:)).* ...
                 (y([2:siz(1) 1],:)   y(:,:)))/2),[1 siz(2:end)]);
else
  area = sum(x); % SUM produces the right value for all empty cases
end
if nargin==2
  area = shiftdim(area,-nshifts);
elseif nargin==3
  area = ipermute(area,perm);
end

经过过冷水的提炼出来有效语句

siz = size(x);
area = reshape(abs(sum( (x([2:siz(1) 1],:) - x(:,:)).* ...
                 (y([2:siz(1) 1],:)   y(:,:)))/2),[1 siz(2:end)]);

好了！求面积就是使用这个长的公式来完成计算的，我们得到了计算面积的底层公式，可是还是看不懂啊！所以的依据算法来设计程序帮我我们理解，根据小学知识知道，欲求多边形的面积可以将多边形转换成多个三角形

所以就转化成求三角形的面积，然而已知三点该如何求三角形的坐标？过冷水当然知道一定可以求，关键是要尽可能的简单，有这么一个公式可以用：

在平面直角坐标系中A(x1,y1)、B(x2,y2)、C(x2,y2)构成的三角形面积公式为：

观察公式S可以化简为

多边形面积就可以用该公式做计算

x=[252,251,250,249,248,247,246,245,245,244,244,241,241,240,240,239,239,238,238,239,239,240,240,241,241,242,242,243,243,244,244,246,246,247,248,248,251,251,252,252,253,256,257,259,260,261,263,264,266,268,269,270,277,278,280,281,288,289,294,303,303,304,305,307,307,308,308,309,309,308,308,307,307,306,306,305,305,304,304,302,302,301,301,300,300,299,299,298,298,296,296,290,289,287,286,284,283,282,280,277,276,275,273,270,269,268,267,266,264,263,252];
y=[580,581,581,582,582,583,583,584,586,587,588,591,595,596,602,603,621,622,629,630,632,633,635,636,637,638,640,641,642,643,644,646,648,648,649,650,653,654,655,656,656,659,659,661,661,662,662,663,663,665,665,666,666,667,667,666,666,665,665,656,655,654,654,652,649,648,645,644,635,634,629,628,626,625,623,622,620,619,618,616,614,613,611,610,609,608,607,606,604,602,601,595,595,593,593,591,591,590,590,587,587,586,586,583,583,582,582,581,581,580,580];
figure1 = figure;
axes1 = axes('Parent',figure1);
hold(axes1,'on');
plot(x,y,'LineWidth',3,'Color',[0 0 1]);
patch('Parent',axes1,'YData',y,'XData',x,'LineWidth',3,'FaceColor',[1 0 0]);
xlim(axes1,[200 360]);
ylim(axes1,[520 700]);
hold(axes1,'off');
n=length(x);
s=0
for i=1:n-1
    a=x(i)*y(i 1)-x(i 1)*y(i);
    s=s a;
end
S=0.5*s;

这就是一个完整的计算多边形面积的程序，逻辑性好的读者就会发现两个算面积的程序实际是一样的可以互换，算面积本来也不是很难，过冷水想和大家分享的是解决问题的思路，有时候编程问题其实的返回到数学问题，如果数学上你理解不好，那么好多程序就看不懂，让我们大家在学习数学的道路一起奔跑吧！

以上就是本期过冷水和大家分享的，有关于封闭多边形的面积算法知识，欢迎大家积极留言讨论。