判断线条是否相交
判断线条是否相交，3位立体点是否相交，平行等
Delphi
Author: Arash Partow
Orientation in 2D of a point relative to two other points?
function Orientation(x1, y1, x2, y2, Px, Py: Double): Integer;
var
Orin: Double;
begin
(* Linear determinant of the 3 points *)
Orin := (x2 - x1) * (py - y1) - (px - x1) * (y2 - y1);
if Orin > 0.0 then Result := +1 (* Orientaion is to the right-hand side *)
else if Orin < 0.0 then Result := -1 (* Orientaion is to the left-hand side *)
else
Result := 0; (* Orientaion is neutral if result is 0 *)
end;
(* End Of Orientation *)
Determine if two 2d segments intersect?
function Intersect(x1, y1, x2, y2, x3, y3, x4, y4: Double): Boolean;
begin
Result := (Orientation(x1, y1, x2, y2, x3, y3) <> Orientation(x1, y1, x2, y2, x4, y4))
and
(Orientation(x3, y3, x4, y4, x1, y1) <> Orientation(x3, y3, x4, y4, x2, y2));
end;
(* End Of SegmentIntersect *)
Orientation in 3D of a point relative to a plane?
function Orientation(x1, y1, z1, x2, y2, z2, x3, y3, z3, Px, Py, Pz: Double): Integer;
var
Px1, Px2, Px3: Double;
Py1, Py2, Py3: Double;
Pz1, Pz2, Pz3: Double;
Orin: Double;
begin
Px1 := x1 - px;
Px2 := x2 - px;
Px3 := x3 - px;
Py1 := y1 - py;
Py2 := y2 - py;
Py3 := y3 - py;
Pz1 := z1 - pz;
Pz2 := z2 - pz;
Pz3 := z3 - pz;
Orin := Px1 * (Py2 * Pz3 - Pz2 * Py3) +
Px2 * (Py3 * Pz1 - Pz3 * Py1) +
Px3 * (Py1 * Pz2 - Pz1 * Py2);
if Orin < 0.0 then Result := -1 (* Orientaion is below plane *)
else if Orin > 0.0 then Result :=
+1 (* Orientaion is above plane *)
else
Result := 0; (* Orientaion is coplanar to plane if result is 0 *)
end;
(* End Of Orientation *)
Signed Volume in 2D?
function Signed(x1, y1, x2, y2, Px, Py: Double): Double;
begin
Result := (x2 - x1) * (py - y1) - (px - x1) * (y2 - y1);
end;
(* End Of Signed *)
Signed Volume in 3D?
function Signed(x1, y1, z1, x2, y2, z2, x3, y3, z3, Px, Py, Pz: Double): Double;
var
Px1, Px2, Px3: Double;
Py1, Py2, Py3: Double;
Pz1, Pz2, Pz3: Double;
begin
Px1 := x1 - px;
Px2 := x2 - px;
Px3 := x3 - px;
Py1 := y1 - py;
Py2 := y2 - py;
Py3 := y3 - py;
Pz1 := z1 - pz;
Pz2 := z2 - pz;
Pz3 := z3 - pz;
Result := Px1 * (Py2 * Pz3 - Pz2 * Py3) +
Px2 * (Py3 * Pz1 - Pz3 * Py1) +
Px3 * (Py1 * Pz2 - Pz1 * Py2);
end;
(* End Of Signed *)
Determine if 3 points are collinear in 2D?
function Collinear(x1, y1, x2, y2, x3, y3: Double): Boolean;
begin
Result := (((x2 - x1) * (y3 - y1) - (x3 - x1) * (y2 - y1)) = 0);
end;
(* End Of Collinear *)
Determine if 3 points are collinear in 3D?
function Collinear(x1, y1, z1, x2, y2, z2, x3, y3, z3: Double): Boolean;
var
Dx1, Dx2: Double;
Dy1, Dy2: Double;
Dz1, Dz2: Double;
Cx, Cy, Cz: Double;
//Var AB,AC,BC:Double;
begin
{find the difference between the 2 points P2 and P3 to P1 }
Dx1 := x2 - x1;
Dy1 := y2 - y1;
Dz1 := z2 - z1;
Dx2 := x3 - x1;
Dy2 := y3 - y1;
Dz2 := z3 - z1;
{perform a 3d cross product}
Cx := Dy1 * Dz2 - Dy2 * Dz1;
Cy := Dx2 * Dz1 - Dx1 * Dz2;
Cz := Dx1 * Dy2 - Dx2 * Dy1;
Result := IsEqual(Cx * Cx + Cy * Cy + Cz * Cz, 0.0);
{
Note:
The method below is very stable and logical, however at the same time
it is "VERY" inefficient, it requires 3 SQRTs which is not acceptable...
Result:=False;
AB:=Distance(x1,y1,z1,x2,y2,z2);
AC:=Distance(x1,y1,z1,x3,y3,z3);
BC:=Distance(x2,y2,z2,x3,y3,z3);
If (AB+AC) = BC Then Result:=True
Else
If (AB+BC) = AC Then Result:=True
Else
If (AC+BC) = AB Then Result:=True;
}
end;
(* End Of Collinear *)
Determine the point at which two 2D segments intersect?
procedure IntersectPoint(x1, y1, x2, y2, x3, y3, x4, y4: Double; var Nx, Ny: Double);
var
R: Double;
dx1, dx2, dx3: Double;
dy1, dy2, dy3: Double;
begin
dx1 := x2 - x1;
dx2 := x4 - x3;
dx3 := x1 - x3;
dy1 := y2 - y1;
dy2 := y1 - y3;
dy3 := y4 - y3;
R := dx1 * dy3 - dy1 * dx2;
if R <> 0 then
begin
R := (dy2 * (x4 - x3) - dx3 * dy3) / R;
Nx := x1 + R * dx1;
Ny := y1 + R * dy1;
end
else
begin
if Collinear(x1, y1, x2, y2, x3, y3) then
begin
Nx := x3;
Ny := y3;
end
else
begin
Nx := x4;
Ny := y4;
end;
end;
end;
function Collinear(x1, y1, x2, y2, x3, y3: Double): Boolean;
begin
Result := (((x2 - x1) * (y3 - y1) - (x3 - x1) * (y2 - y1)) = 0);
end;
(* End Of Collinear *)
Determine the vertex angle made by two 3D segments?
function VertexAngle(x1, y1, z1, x2, y2, z2, x3, y3, z3: Double): Double;
var
Dist: Double;
begin
(* Quantify coordinates *)
x1 := x1 - x2;
x3 := x3 - x2;
y1 := y1 - y2;
y3 := y3 - y2;
z1 := z1 - z2;
z3 := z3 - z2;
(* Calculate Lay Distance *)
Dist := (x1 * x1 + y1 * y1 + z1 * z1) * (x3 * x3 + y3 * y3 + z3 * z3);
if IsEqual(Dist, 0) then Result := 0.0
else
Result := ArcCos((x1 * x3 + y1 * y3 + z1 * z3) / sqrt(Dist)) * _180DivPI;
end;
(* End Of VertexAngle *)
Determine the point at which two 2D segments intersect?
procedure IntersectPoint(x1, y1, x2, y2, x3, y3, x4, y4: Double; var Nx, Ny: Double);
var
R: Double;
dx1, dx2, dx3: Double;
dy1, dy2, dy3: Double;
begin
dx1 := x2 - x1;
dx2 := x4 - x3;
dx3 := x1 - x3;
dy1 := y2 - y1;
dy2 := y1 - y3;
dy3 := y4 - y3;
R := dx1 * dy3 - dy1 * dx2;
if R <> 0 then
begin
R := (dy2 * (x4 - x3) - dx3 * dy3) / R;
Nx := x1 + R * dx1;
Ny := y1 + R * dy1;
end
else
begin
if Collinear(x1, y1, x2, y2, x3, y3) then
begin
Nx := x3;
Ny := y3;
end
else
begin
Nx := x4;
Ny := y4;
end;
end;
end;
function Collinear(x1, y1, x2, y2, x3, y3: Double): Boolean;
begin
Result := (((x2 - x1) * (y3 - y1) - (x3 - x1) * (y2 - y1)) = 0);
end;
(* End Of Collinear *)
...Perpendicular from Point to Segment in 2D?
procedure PerpendicularPntToSegment(x1, y1, x2, y2, Px, Py: Double; var Nx, Ny: Double);
var
R: Double;
Dx: Double;
Dy: Double;
begin
Dx := x2 - x1;
Dy := y2 - y1;
R := ((Px - x1) * Dx + (Py - y1) * Dy) / Sqr(Dx * Dx + Dy * Dy);
Nx := x1 + R * Dx;
Ny := y1 + R * Dy;
end;
(* End PerpendicularPntSegment *)
Perpendicular distance from Point to Segment in 2D?
function PntToSegmentDistance(Px, Py, x1, y1, x2, y2: Double): Double;
var
Ratio: Double;
Dx: Double;
Dy: Double;
begin
if IsEqual(x1, x2) and IsEqual(y1, y2) then
begin
Result := Distance(Px, Py, x1, y1);
end
else
begin
Dx := x2 - x1;
Dy := y2 - y1;
Ratio := ((Px - x1) * Dx + (Py - y1) * Dy) / (Dx * Dx + Dy * Dy);
if Ratio < 0 then Result := Distance(Px, Py, x1, y1)
else if Ratio > 1 then Result := Distance(Px, Py, x2, y2)
else
Result := Distance(Px, Py, (1 - Ratio) * x1 + Ratio * x2,
(1 - Ratio) * y1 + Ratio * y2);
end;
end;
(* End PntToSegmentDistance *)
Determine if two 2D segments are parallel?
function SegmentsParallel(x1, y1, x2, y2, x3, y3, x4, y4: Double): Boolean;
begin
Result := (((y1 - y2) * (x1 - x2)) = ((y3 - y4) * (x3 - x4)));
end;
(* End Of SegmentsParallel *)
Determine if two 3D segments are parallel?
function SegmentsParallel(x1, y1, z1, x2, y2, z2, x3, y3, z3, x4, y4, z4: Double): Boolean;
var
Dx1, Dx2: Double;
Dy1, Dy2: Double;
Dz1, Dz2: Double;
begin
{
Theory:
If the gradients in the following planes x-y, y-z, z-x are equal then it can be
said that the segments are parallel in 3D, However as of yet I haven't been able
to prove this "mathematically".
Worst case scenario: 6 floating point divisions and 9 floating point subtractions
}
Result := False;
{
There is a division-by-zero problem that needs attention.
My initial solution to the problem is to check divisor of the divisions.
}
Dx1 := x1 - x2;
Dx2 := x3 - x4;
//If (IsEqual(dx1,0.0) Or IsEqual(dx2,0.0)) And NotEqual(dx1,dx2) Then Exit;
Dy1 := y1 - y2;
Dy2 := y3 - y4;
//If (IsEqual(dy1,0.0) Or IsEqual(dy2,0.0)) And NotEqual(dy1,dy2) Then Exit;
Dz1 := z1 - z2;
Dz2 := z3 - z4;
//If (IsEqual(dy1,0.0) Or IsEqual(dy2,0.0)) And NotEqual(dy1,dy2) Then Exit;
if NotEqual(Dy1 / Dx1, Dy2 / Dx2) then Exit;
if NotEqual(Dz1 / Dy1, Dz2 / Dy2) then Exit;
if NotEqual(Dx1 / Dz1, Dx2 / Dz2) then Exit;
Result := True;
end;
(* End Of SegmentsParallel*)
const
Epsilon = 1.0E-12;
function IsEqual(Val1, Val2: Double): Boolean;
var
Delta: Double;
begin
Delta := Abs(Val1 - Val2);
Result := (Delta <= Epsilon);
end;
(* End Of Is Equal *)
function NotEqual(Val1, Val2: Double): Boolean;
var
Delta: Double;
begin
Delta := Abs(Val1 - Val2);
Result := (Delta > Epsilon);
end;
(* End Of Not Equal *)
determine if two 2D segments are perpendicular?
function SegmentsPerpendicular(x1, y1, x2, y2, x3, y3, x4, y4: Double): Boolean;
begin
Result := (((y2 - y1) * (x3 - x4)) = ((y4 - y3) * (x2 - x1) * -1));
end;
(* End Of SegmentsPerpendicular *)
Determine if two 3D segments are perpendicular?
function SegmentsPerpendicular(x1, y1, z1, x2, y2, z2, x3, y3, z3, x4, y4, z4: Double): Boolean;
var
Dx1, Dx2: Double;
Dy1, Dy2: Double;
Dz1, Dz2: Double;
begin
{
The dot product of the vector forms of the segments will be
0 if the segments are perpendicular
}
Dx1 := x1 - x2;
Dx2 := x3 - x4;
Dy1 := y1 - y2;
Dy2 := y3 - y4;
Dz1 := z1 - z2;
Dz2 := z3 - z4;
Result := (((Dx1 * Dx2) + (Dy1 * Dy2) + (Dz1 * Dz2)) = 0)
end;
(* End Of *)
Determine if a 2D point exists within a 2D triangle?
function PntInTriangle(Px, Py, x1, y1, x2, y2, x3, y3: Double): Boolean;
var
Or1, Or2, Or3: Double;
begin
Or1 := Orientation(x1, y1, x2, y2, Px, Py);
Or2 := Orientation(x2, y2, x3, y3, Px, Py);
Or3 := Orientation(x3, y3, x1, y1, Px, Py);
Result := (Or1 = Or2) and (Or2 = Or3);
end;
(* End Of PntInTriangle *)
function Orientation(x1, y1, x2, y2, Px, Py: Double): Integer;
var
Orin: Double;
begin
(* Linear determinant of the 3 points *)
Orin := (x2 - x1) * (py - y1) - (px - x1) * (y2 - y1);
if Orin > 0.0 then Result := +1 (* Orientaion is to the right-hand side *)
else if Orin < 0.0 then Result := -1 (* Orientaion is to the left-hand side *)
else
Result := 0; (* Orientaion is neutral if result is 0 *)
end;
(* End Of Orientation *)
Determine if a 2D point exists within a 2D triangle?
判断一个点是否在三角形内
function PntInTriangle(Px, Py, x1, y1, x2, y2, x3, y3: Double): Boolean;
var
Or1, Or2, Or3: Double;
begin
Or1 := Orientation(x1, y1, x2, y2, Px, Py);
Or2 := Orientation(x2, y2, x3, y3, Px, Py);
Or3 := Orientation(x3, y3, x1, y1, Px, Py);
Result := (Or1 = Or2) and (Or2 = Or3);
end;
(* End Of PntInTriangle *)
function Orientation(x1, y1, x2, y2, Px, Py: Double): Integer;
var
Orin: Double;
begin
(* Linear determinant of the 3 points *)
Orin := (x2 - x1) * (py - y1) - (px - x1) * (y2 - y1);
if Orin > 0.0 then Result := +1 (* Orientaion is to the right-hand side *)
else if Orin < 0.0 then Result := -1 (* Orientaion is to the left-hand side *)
else
Result := 0; (* Orientaion is neutral if result is 0 *)
end;
(* End Of Orientation *)
determine the 3rd point of an Equilateral Triangle from the other 2 points?
procedure CreateEquilateralTriangle(x1, y1, x2, y2: Double; var x3, y3: Double);
const
Sin60 = 0.86602540378443864676372317075294;
const
Cos60 = 0.50000000000000000000000000000000;
begin
{ Translate for x1,y1 to be origin }
x2 := x2 - x1;
y2 := y2 - y1;
{ Rotate 60 degrees and translate back }
x3 := ((x2 * Cos60) - (y2 * Sin60)) + x1;
y3 := ((y2 * Cos60) + (x2 * Sin60)) + y1;
end;
(* End Of Create Equilateral Triangle *)