使用Matlab计算两条线的交点及三角形垂心

Matlab版本：R2016A

操作系统：Win-8

为了让整个流程比较完整，

1、我用了12个点，每两个点一条线，能组成6条线；

2、每两条线相交于一点，一共有3个点；

3、3个点构成了三角形的三个顶点，然后再根据三个顶点计算三角形的垂心。

1、现在先给出12个点的坐标（坐标可以随意设置，只要构成的线不是平行没有交点即可）

(294.94, 2959.12) (668.61, 1464.42)

(1678.70, 2848.19) (1731.25, 1914.00)
(540.16, 1978.22) (1731.25, 1919.83)

(335.80, 2778.12) (1696.22, 2678.86)
(1106.51, 892.23) (1404.28, 1832.25)

(1351.74, 880.56) (1591.12, 1820.58)

上述点是我从一幅图片里面点选出来的，将每个点的橫坐标放在一个数组里，纵坐标放在一个数组里，然后进行绘制：

X = [294.94, 668.61, 1678.70, 1731.25, 540.16, 1731.25, 335.80, 1696.22, 1106.51, 1404.28, 1351.74, 1591.12];  
Y = [2959.12, 1464.42, 2848.19, 1914.00, 1978.22, 1919.83, 2778.12, 2678.86, 892.23, 1832.25, 880.56, 1820.58];
plot(X, Y, '.');

绘制出来的结果如图所示：

将12个点组成6条线：

for i = 1: 6
    plot([X(2*i-1), X(2*i)], [Y(2*i-1), Y(2*i)]);
    hold on;
end

依次绘制出的线如图所示：

2、下面计算每两条线之间的交点（即上图所示的1,2一组，3,4一组，5,6一组计算其交点）
两条直线相交，必然需要求出两条直线的表达式，每条直线的表达式可以用 y = ax + b来表示，下面用matlab求解a和b两个变量，a即为直线的斜率，b为直线的截距。

a = zeros(6, 1);
b = zeros(6, 1);
for i = 1: 6
    a(i, 1) = (Y(2*i) - Y(2*i-1)) / (X(2*i) - X(2*i-1));
    b(i, 1) = Y(2*i) - a(i, 1) * X(2*i);
end

上述公式即求出了全部直线的方程式，下面计算它们两两之间的交点

points = zeros(3, 2);
for i = 1: 3
    % 两个直线方程之差的函数
    fun = @(x) (a(2*i-1, 1) * x + b(2*i-1, 1)) - (a(2*i, 1) * x + b(2*i, 1));
    
    % 两个直线方程之差函数值为0时的x值，由得到的x值计算y值
    points(i, 1) = fzero(fun, 2);
    points(i, 2) = a(2*i) * points(i, 1) + b(2*i, 1);
end

得到直线交点如图所示：

3、将三个点连线绘制，可以得到三角形

hold on;
plot([points(:, 1); points(1, 1)], [points(:, 2); points(1, 2)])

绘制得到结果如图所示：

计算三角形每条边的垂线，需要先计算每条线的斜率：

% 计算各个顶点连线的方程式
slope1_2 = (points(2, 2) - points(1, 2)) / (points(2, 1) - points(1, 1));
slope2_3 = (points(3, 2) - points(2, 2)) / (points(3, 1) - points(2, 1));
slope1_3 = (points(3, 2) - points(1, 2)) / (points(3, 1) - points(1, 1));
b1_2 = points(2, 2) - slope1_2 * points(2, 1);
b2_3 = points(3, 2) - slope2_3 * points(3, 1);
b1_3 = points(1, 2) - slope1_3 * points(1, 1);
slopevertical1_2 = -1 / slope1_2;
slopevertical2_3 = -1 / slope2_3;
slopevertical1_3 = -1 / slope1_3;

verticalb1_2 = points(3, 2) - slopevertical1_2 * points(3, 1);
verticalb2_3 = points(1, 2) - slopevertical2_3 * points(1, 1);
verticalb1_3 = points(2, 2) - slopevertical1_3 * points(2, 1);

这样就求出了三个边的垂线，下面对这三个垂线进行绘制，按照上面的方程计算出两个点，两个点的横坐标一个选择很小，一个选择很大，这样可以保证画出来的直线足够长，方便观察（这个方法似乎有些简单粗暴~~~~~~~~，不过还是能够绘制出来的）

% 设置预绘制直线的两个端点的X值
verticalLinePoints = zeros(6, 2)
for i = 1: 3
    verticalLinePoints(2*i-1, 1) = 0;
    verticalLinePoints(2*i, 1) = 10000;
    if i == 1
         slopevertical = slopevertical1_2;
         verticalb = verticalb1_2;
    elseif i == 2
         slopevertical = slopevertical2_3;
         verticalb = verticalb2_3;
    else
	 slopevertical = slopevertical1_3;
	 verticalb = verticalb1_3;
    end
    
   % 计算两个端点的y值
    verticalLinePoints(2*i-1, 2) = slopevertical * verticalLinePoints(2*i-1, 1) + verticalb;
    verticalLinePoints(2*i, 2) = slopevertical * verticalLinePoints(2*i, 1) + verticalb;
end

然后绘制所得结果（为了使图更加简洁，开始的点线都去掉了，仅保留了三角形）

hold on;	% 为了保留前面绘制的图像，使用hold on
plot(verticalLinePoints(1:2, 1), verticalLinePoints(1:2, 2));
plot(verticalLinePoints(3:4, 1), verticalLinePoints(3:4, 2));
plot(verticalLinePoints(5:6, 1), verticalLinePoints(5:6, 2));
axis equal;	% 使横纵轴单位统一

绘制得到结果如图所示：

由图或者基础几何可知，三条垂线交于一点，称之为垂心（中心），求垂心的方法与前面求两条线的交点的方法一样

fun = @(x) (slopevertical1_2 * x + verticalb1_2) - (slopevertical1_3 * x + verticalb1_3);
finalX = fzero(fun, 2);
finalY = slopevertical1_2 * finalX + verticalb1_2;
hold on;
plot(finalX, finalY, 'o');

至此，计算垂心的工作已经完成，对垂心进行绘制，绘制的结果如图所示：

代码下载地址：http://download.csdn.net/detail/zhangwellyear/9796085