Matlab研究小问题：如何计算一条线段所经过的网格区域和各区域内的长度

问题描述

在一个10*10的网格里，判断一条线段经过的网格位置（判断序号）。并计算经过的每个网格内线段的长度。

出于问题简单化的需要，当线段正好在网格线上时，不计算经过的网格和长度。

分析

这个问题的复杂度在于，需要计算线段与网格线的交点，进而判断线段经过网格的序号。由于线段是以两个端点的形式给出，除了垂直于横轴以外，均可以通过斜截式来表示直线。因此需要对这两种情况分别进行处理，加上不计算网格和长度的情形，一共需要分三种情况进行处理。

设计一个calLength函数来进行相关计算，设计segs结构体来存储结果。

segs.index —— 经过网格的序号；

segs.index_x —— 经过网格的横轴序号；

segs.index_y —— 经过网格的纵轴序号；

segs.length —— 经过网格内的线段长度。

分段结果以结构数组形式存储。

代码

function segs = calLength(P1,P2)

segs.length = [];
segs.index_x = [];
segs.index_y = [];
segs.index = [];

P1_x = P1(1);
P1_y = P1(2);
P2_x = P2(1);
P2_y = P2(2);

xmin = min(P1_x,P2_x);
ymin = min(P1_y,P2_y);
xmax = max(P1_x,P2_x);
ymax = max(P1_y,P2_y);

if (P1_x == P2_x) && (round(P1_x) == P1_x)...
|| (P1_y == P2_y) && (round(P1_y) == P1_y)
return;
end

pos = reshape(1:100,10,10);

if (P1_x == P2_x) && (round(P1_x) ~= P1_x)
SP = unique([ymin,ymax,ceil(ymin):floor(ymax)]);
for t = 1 : size(SP,2)-1
segs(t).length = SP(t+1)-SP(t);
segs(t).index_x = ceil(P1_x);
segs(t).index_y = max(ceil(SP(t+1)),ceil(SP(t)));
segs(t).index = pos(segs(t).index_x, segs(t).index_y);
end
end

if (P1_x ~= P2_x)
K = polyfit([P1_x,P2_x],[P1_y,P2_y],1);
xpx = []; xpy = [];
for i = ceil(xmin):floor(xmax)
xpx(i-ceil(xmin)+1) = i;
xpy(i-ceil(xmin)+1) = K(1)*i+K(2);
end

ypx = []; ypy = [];
for j = ceil(ymin):floor(ymax)
ypy(j-ceil(ymin)+1) = j;
syms x;
ypx(j-ceil(ymin)+1) = double(solve(K(1)*x + K(2) - j, x));
end

SP = unique([P1_x,P2_x,xpx,ypx;P1_y,P2_y,xpy,ypy]','rows');

L = @(x) sqrt((SP(x+1,1) - SP(x,1)).^2 + (SP(x+1,2) - SP(x,2)).^2);

for t = 1 : size(SP,1)-1
segs(t).length = L(t);
segs(t).index_x = max(ceil(SP(t+1,1)),ceil(SP(t,1)));
segs(t).index_y = max(ceil(SP(t+1,2)),ceil(SP(t,2)));
segs(t).index = pos(segs(t).index_x, segs(t).index_y);
end
end

n = [];
for i = 1 : size(segs,2)
if segs(i).length < eps
n = [n i];
end
end
segs(n) = [];

clear;clc;clf;
x = linspace(0,10,11);
y = linspace(0,10,11);
[X,Y] = meshgrid(x,y);
line(X,Y,'color','b');
line(X',Y','color','b');
axis equal;
axis([0 10 0 10]);
set(gca,'xtick',0:10);

gridindex = reshape(1:100,10,10)';
numposx = 0.5*(X(1:end-1,2:end)+X(1:end-1,1:end-1))-0.1;
numposy = 0.5*(Y(2:end,1:end-1)+Y(1:end-1,1:end-1));
for i = 1 : 10
for j = 1 : 10
text(numposx(i,j),numposy(i,j),num2str(gridindex(i,j)));
end
end

P1 = input('P1=');
P2 = input('P2=');
segs = calLength(P1,P2);
line([P1(1) P2(1)],[P1(2) P2(2)],'color','r');

display('所经过的网格序号\长度分别为:');
for i = 1 : size(segs,2)
display(['序号: ' num2str(segs(i).index)]);
display(['长度: ' num2st
af5a
r(segs(i).length)]);
end

示例结果

P1=[8.7,9]

P2=[3.5,2]

所经过的网格序号\长度分别为:

序号: 24

长度: 0.83847

序号: 25

长度: 0.40726

序号: 35

长度: 1.2457

序号: 45

长度: 0.023956

序号: 46

长度: 1.2218

序号: 56

长度: 0.45517

序号: 57

长度: 0.79056

序号: 67

长度: 0.88638

序号: 68

长度: 0.35934

序号: 78

长度: 1.2457

序号: 88

长度: 0.071869

序号: 89

长度: 1.1739

经核对结果是正确的。

注：本例中对特殊情况也进行了测试，详细结果限于篇幅不在此做验证。