求空间直线的交点.docx

1、求空间直线的交点因为一个小任务需要用到求空间直线交点的 MATLAB函数和求空间中两个平面的相交线的函数，但是在网上找了一下没有，只好自己写了几个函数，自 己觉得还能用，在这里给大家分享一下。1.MATLAB求空间中的两个平面的相交线functionflag,n,p = Plane2Line(n1,p1,n2,p2)% % calulate the line shared by two intersecting plane% input:% n1 normal vector of plane one% p1 any point on plane one% n2 normal vector of

2、 plane two% p2 any point on plane two% output:% flag whether the two planes are intersecting ( 1 or 0)% n the direction vector of the expected line% p any point in the expected line% author:Lai Zhenzhou from Harbin Institute of Technology% email:% date:2014.1.16% if(isvector(n1) & isvector(p1) & isv

3、ector(n2) &isvector(p2)error(Plane2Line:the parameter is not vector);endif(length(n1)=3)|(length(p1)=3)|(length(n2)=3)|(length(p2)=3) error(Plane2Line:the parameter is not 3d vector);endA = n1(1)n1(2)n1(3);n2(1)n2(2)n2(3);if(rank(A)2)flag = 0;elseflag = 1;endif(flag=1)% calculate the normal vectorn

4、= cross(n1,n2);c1 = n1(1)n1(2)n1(3)-dot(n1,p1);c2 = n2(1)n2(2)n2(3)-dot(n2,p2);% calculate the simplest Row echelon matrix temp1 = rref(A);temp2 = 1 2 3;index(1)= find(temp1(1,:),1,first);% the first nonzero elementindex(2)= find(temp1(2,:),1,first);% the first nonzero elementD = A(: ,index(1) A(:,i

5、ndex(2);Y = dot(n1,p1);dot(n2,p2);X = inv(D)*Y; fori=1:3 if(i=index(1)& i=index(2) index(3)= i;endend p(index(1) = X(1); p(index(2) = X(2); p(index(3) = 0; else n=;p=; endfor test% flag n p = Plane2Line(1 2 3,1 0 1,2 3 4,0 -1 0)% flag n p =Pla ne2Li ne(O 0 1,0.5 0.5 0.5,1 0 0,0 0 1)2.MATLAB求空间中两条直线的

6、交点functionflag,p = Line2Point(n1,p1,n2,p2)% % determine the relation between two straight lines and% calulate the intersection point if they are intersecting% input:% n1 direction vector of line one% p1 any point in line one% n2 direction vector of line two% p2 any point in line two% output:% flag t

7、he relation of the two line% flag = 0 the two line are on different plane% flag=1 thetwo line are on the same plane and they are parallel%flag = 2 the two line are on the same plane and they areintersecting% p the point shared by the two intersecting line% author:Lai Zhenzhou from Harbin Institute o

8、f Technology% email:% date:2014.1.17% reference:http:% if(isvector(n1) & isvector(p1) & isvector(n2) & isvector(p2)error(Line2Point:the parameter is not vector);endif(length(n1)=3)|(length(p1)=3)|(length(n2)=3)|(length(p2)=3)error(Line2Point:the parameter is not 3d vector);endA = p2(1)-p1(1)p2(2)-p1

9、(2)p2(3)-p1(3);n1(1)n1(2)n1(3);n2(1)n2n2(2)(3); if(det(A)=0) flag = 0; else if(rank(A(2:3,: )2) flag = 1; else flag = 2; end end if(flag = 2) B = rref(A(2:3,: );index(1)= find(B(1,: ),1,first); index(2)= find(B(2,: ),1,first);fori=1:3if(i=index(1)& i=index(2)index(3)= i;endendY(1,1) = -p1(index(1) +

10、 p2(index(1);Y(2,1) = -p1(index(2) + p2(index(2);D = n1(index(1) -n2(index(1);n1(index(2) -n2(index(2);t = inv(D)*Y;p(1) = p1(1) + n1(1)*t(1);p(2) = p1(2) + n1(2)*t(1);p(3) = p1(3)+ n1(3)*t(1);elsep = ;endend% for test % flag p = Line2Point(1 0 0,0 0 0,0 1 0,1 1 0)% flag p = Line2Point(1 1 0,0 0 0,0

11、 1 0,1 1 0)% flag p = Line2Point(1 1 -1,0 0 1,0 1 0,1 1 0)% flag p = Line2Point(1 1 - 1,0 0 1,1 1 -1,1 1 0)% flag p = Line2Point(1 1 1,0 0 0,2 0 0,2 1 1)% flagp=Line2Point(0-42050 -21025,-100,00145,00 0)3.MATLAB求空间中直线和某一线段的交点functionflag,p = LineSegment2Point(n1,p1,p21,p22)% % determine the relation

12、 between a straight line and% a line segment and calulate the intersecting point% if they are intersecting% input:% n1 direction vector of the straight line% p1 any point in the straight line% p21 any point in the line segment% p22 any another different point in the line segment% output:% flag the r

13、elation of the two line% flag=0 thestraight line and line segment are not intersecting% flag = 1 the straight line and line segment are intersecting% p the point shared by the two intersecting line% author:Lai Zhenzhou from Harbin Institute of Technology% email:% date:2014.1.17% if(isvector(n1) & is

14、vector(p1) & isvector(p21) & isvector(p22)error(LineSegment2Point:the parameter is not vector);endif(length(n1)=3)|(length(p1)=3)|(length(p21)=3)|(length(p22)=3)error(LineSegment2Point:the parameter is not 3d vector);endif(p21(1)=p22(1)&(p21(2)=p22(2)&(p21(3)=p22(3)error(LineSegment2Point:the two ap

15、pointed points of line segment are the same);end%calculate the direction vector of line segmentn2 = p22(1)-p21(1),p22(2)-p21(2),p22(3)-p21(3);flag0,p0 = Line2Point(n1,p1,n2,p21);if(flag0 = 2)if(norm(p0-p21)=norm(p22-p21) & norm(p0-p22)=norm(p22-p21)flag = 1; p = p0;elseflag = 0;p=;endelseflag = 0;p=;endend% for test % flag p = LineSegment2Point(0 0 1,0 0 0,0 1 1,0 2 1)% flag p = LineSegment2Point(1 1 1,0 0 0,0 1 1,0 2 1)% flag p = LineSegment2Point(1 1 1,0 0 0,0 1 1,1 1 1)% flag p = LineSegment2Point(1 1 1,0 0 0,0 1 1,2 1 1)