假设AB1, AB2, CD1, CD2
.AB 1&AB 2和CD 1&CD 23D点构成一条线段。并且所述线段 * 不在同一平面内 。AP
是点线段AB 1&AB 2,BP
是点线段CD 1&CD 2。Point1
和Point2
彼此最接近( 两条线段之间的最短距离 *)
现在,我怎样才能找到这两个点Point1
和Point2
?我应该用什么方法?
这只是部分解决...因为这个功能不工作时,两条线是在同一平面...感谢@MBo我遇到了Geometry GoldMine of Code and Explanations!他们有许多源代码的贡献者!我从那里挑选了一个在这里,它是干净和伟大的!
bool CalculateLineLineIntersection(Vector3D p1, Vector3D p2, Vector3D p3, Vector3D p4, Vector3D& resultSegmentPoint1, Vector3D& resultSegmentPoint2)
{
// Algorithm is ported from the C algorithm of
// Paul Bourke at http://local.wasp.uwa.edu.au/~pbourke/geometry/lineline3d/
resultSegmentPoint1 = { 0,0,0 };
resultSegmentPoint2 = { 0,0,0 };
Vector3D p13 = VectorMinus(p1, p3);
Vector3D p43 = VectorMinus(p4, p3);
/*if (p43.LengthSq() < Math.Epsilon) {
return false;
}*/
Vector3D p21 = VectorMinus(p2, p1);
/*if (p21.LengthSq() < Math.Epsilon) {
return false;
}*/
double d1343 = p13.x * (double)p43.x + (double)p13.y * p43.y + (double)p13.z * p43.z;
double d4321 = p43.x * (double)p21.x + (double)p43.y * p21.y + (double)p43.z * p21.z;
double d1321 = p13.x * (double)p21.x + (double)p13.y * p21.y + (double)p13.z * p21.z;
double d4343 = p43.x * (double)p43.x + (double)p43.y * p43.y + (double)p43.z * p43.z;
double d2121 = p21.x * (double)p21.x + (double)p21.y * p21.y + (double)p21.z * p21.z;
double denom = d2121 * d4343 - d4321 * d4321;
/*if (Math.Abs(denom) < Math.Epsilon) {
return false;
}*/
double numer = d1343 * d4321 - d1321 * d4343;
double mua = numer / denom;
double mub = (d1343 + d4321 * (mua)) / d4343;
resultSegmentPoint1.x = (float)(p1.x + mua * p21.x);
resultSegmentPoint1.y = (float)(p1.y + mua * p21.y);
resultSegmentPoint1.z = (float)(p1.z + mua * p21.z);
resultSegmentPoint2.x = (float)(p3.x + mub * p43.x);
resultSegmentPoint2.y = (float)(p3.y + mub * p43.y);
resultSegmentPoint2.z = (float)(p3.z + mub * p43.z);
return true;
}
到目前为止,我已经尝试了所有这些下面的工作,只有当两个线段有相同的幅度...
Link 1Link 2
我尝试计算两条线段的质心,并计算线段上距离中点最近的点。(* 我知道如何计算距离另一个点最近的点线段 *)
但这只适用于当两个线段是相等的长度和每个两个线段的中点是相互垂直和质心...
注:Visual Geometry Geogbra3D用于这些点的视觉表示注:AB1CD
表示从点AB 1到线CD(非线段)
AB1 = (6.550000, -7.540000, 0.000000 )
AB2 = (4.540000, -3.870000, 6.000000 )
CD1 = (0.000000, 8.000000, 3.530000 )
CD2 = (0.030000, -7.240000, -1.340000 )
PointCD1AB = (3.117523, -1.272742, 10.246199 )
PointCD2AB = (6.318374, -7.117081, 0.691420 )
PointAB1CD = (0.029794, -7.135321, -1.306549 )
PointAB2CD = (0.019807, -2.062110, 0.314614 )
Magntidue of PointCD1AB - P1LineSegmentCD = 11.866340
Magntidue of PointCD2AB - P2LineSegmentCD = 6.609495
Magntidue of PointAB1CD - P1LineSegmentAB = 6.662127
Magntidue of PointAB2CD - P2LineSegmentAB = 9.186399
Magntidue of PointCD1AB - PointAB1CD = 13.318028
Magntidue of PointCD2AB - PointAB2CD = 8.084965
Magntidue of PointCD1AB - PointAB2CD = 10.433375
Magntidue of PointCD2AB - PointAB1CD = 6.598368
Actual Shortest Point are
Point1 = (0.01, 1.59, 1.48 )
Point2 = (-1.23, 1.11, 3.13 )
Magnitude of Point1 And Point2 = 2.1190799890518526
对于上面的数据,我使用下面的函数
void NearestPointBetweenTwoLineSegmentOfVariedLength(Vector3D P1LineSegmentAB, Vector3D P2LineSegmentAB, Vector3D P1LineSegmentCD, Vector3D P2LineSegmentCD, Vector3D Testing)
{
/* float Line1Mag = Magnitude(VectorMinus(P1LineSegmentAB, P2LineSegmentAB));
float Line2Mag = Magnitude(VectorMinus(P1LineSegmentCD, P2LineSegmentCD));
P2LineSegmentAB = VectorMinus(P2LineSegmentAB, P1LineSegmentAB);
P1LineSegmentCD = VectorMinus(P1LineSegmentCD, P1LineSegmentAB);
P2LineSegmentCD = VectorMinus(P2LineSegmentCD, P1LineSegmentAB);
P1LineSegmentAB = VectorMinus(P1LineSegmentAB, P1LineSegmentAB);
Vector3D P1P2UnitDirection = GetUnitVector(P2LineSegmentAB, { 0,0,0 });
AngleBetweenTwoVectorsWithCommonUnitVectorAngleOfSecondArgument(P1LineSegmentAB, P2LineSegmentAB, P1P2UnitDirection);*/
Vector3D ReturnVal;
Vector3D PointCD1AB;
Vector3D PointCD2AB;
Vector3D PointAB1CD;
Vector3D PointAB2CD;
NearestPointOnLineFromPoint(P1LineSegmentCD, P1LineSegmentAB, P2LineSegmentAB, PointCD1AB, false);
PrintVector3Dfor(VectorMinus(PointCD1AB, Testing), "PointCD1AB", true);
NearestPointOnLineFromPoint(P2LineSegmentCD, P1LineSegmentAB, P2LineSegmentAB, PointCD2AB, false);
PrintVector3Dfor(VectorMinus(PointCD2AB, Testing), "PointCD2AB", true);
NearestPointOnLineFromPoint(P1LineSegmentAB, P1LineSegmentCD, P2LineSegmentCD, PointAB1CD, false);
PrintVector3Dfor(VectorMinus(PointAB1CD, Testing), "PointAB1CD", true);
NearestPointOnLineFromPoint(P2LineSegmentAB, P1LineSegmentCD, P2LineSegmentCD, PointAB2CD, false);
PrintVector3Dfor(VectorMinus(PointAB2CD, Testing), "PointAB2CD", true);
float m1 = Magnitude(VectorMinus(PointCD1AB, P1LineSegmentCD));
float m2 = Magnitude(VectorMinus(PointCD2AB, P2LineSegmentCD));
float m3 = Magnitude(VectorMinus(PointAB1CD, P1LineSegmentAB));
float m4 = Magnitude(VectorMinus(PointAB1CD, P2LineSegmentAB));
float m5 = Magnitude(VectorMinus(PointCD1AB, PointAB1CD));
float m6 = Magnitude(VectorMinus(PointCD2AB, PointAB2CD));
float m7 = Magnitude(VectorMinus(PointCD1AB, PointAB2CD));
float m8 = Magnitude(VectorMinus(PointCD2AB, PointAB1CD));
Printfloatfor(m1, "Magntidue of PointCD1AB - P1LineSegmentCD");
Printfloatfor(m2, "Magntidue of PointCD2AB - P2LineSegmentCD");
Printfloatfor(m3, "Magntidue of PointAB1CD - P1LineSegmentAB");
Printfloatfor(m4, "Magntidue of PointAB2CD - P2LineSegmentAB");
Printfloatfor(m5, "Magntidue of PointCD1AB - PointAB1CD");
Printfloatfor(m6, "Magntidue of PointCD2AB - PointAB2CD");
Printfloatfor(m7, "Magntidue of PointCD1AB - PointAB2CD");
Printfloatfor(m8, "Magntidue of PointCD2AB - PointAB1CD");
//NearestPointBetweenTwoLineSegmentOfSameLength1(P1LineSegmentAB, P2LineSegmentAB, P1LineSegmentCD, P2LineSegmentCD);
//NearestPointBetweenTwoLineSegmentOfSameLength2(P1LineSegmentAB, P2LineSegmentAB, P1LineSegmentCD, P2LineSegmentCD);
//NearestPointBetweenTwoLineSegmentOfSameLength3(P1LineSegmentAB, P2LineSegmentAB, P1LineSegmentCD, P2LineSegmentCD);
}
void NearestPointOnLineFromPoint(Vector3D Point, Vector3D LineSegmentStart, Vector3D LineSegmentEnd, Vector3D& ReturnVector, bool ClampTheValue)
{
//Get Heading Direction of Capsule from Origin To End
Vector3D CapsuleHeading = VectorMinus(LineSegmentEnd, LineSegmentStart);
float MagnitudeOfLineSegment = Magnitude(CapsuleHeading);
CapsuleHeading = VectorDivide(CapsuleHeading, MagnitudeOfLineSegment);
// Project From Point to Origin
Vector3D Projection = VectorMinus(Point, LineSegmentStart);
float DotProd = DotProduct(Projection, CapsuleHeading);
if (ClampTheValue)
{
DotProd = Clamp(DotProd, 0.0f, MagnitudeOfLineSegment);
}
ReturnVector = VectorAdd(LineSegmentStart, VectorMultiply(CapsuleHeading, DotProd));
}
我已经将此代码从C#转换为C++,但它没有按预期工作...我不知道是代码转换有问题还是代码本身有问题?
Vector3D ClampPointToLine(Vector3D pointToClamp, Vector3D LineStart, Vector3D LineEnd)
{
Vector3D clampedPoint = {0,0,0};
double minX, minY, minZ, maxX, maxY, maxZ;
if (LineStart.x <= LineEnd.x)
{
minX = LineStart.x;
maxX = LineEnd.x;
}
else
{
minX = LineEnd.x;
maxX = LineStart.x;
}
if (LineStart.y <= LineEnd.y)
{
minY = LineStart.y;
maxY = LineEnd.y;
}
else
{
minY = LineEnd.y;
maxY = LineStart.y;
}
if (LineStart.z <= LineEnd.z)
{
minZ = LineStart.z;
maxZ = LineEnd.z;
}
else
{
minZ = LineEnd.z;
maxZ = LineStart.z;
}
clampedPoint.x = (pointToClamp.x < minX) ? minX : (pointToClamp.x > maxX) ? maxX : pointToClamp.x;
clampedPoint.y = (pointToClamp.y < minY) ? minY : (pointToClamp.y > maxY) ? maxY : pointToClamp.y;
clampedPoint.z = (pointToClamp.z < minZ) ? minZ : (pointToClamp.z > maxZ) ? maxZ : pointToClamp.z;
return clampedPoint;
}
void distBetweenLines(Vector3D p1, Vector3D p2, Vector3D p3, Vector3D p4, Vector3D& ClosestPointOnLineP1P2, Vector3D& ClosestPointOnLineP3P4)
{
Vector3D d1;
Vector3D d2;
d1 = VectorMinus(p2, p1);
d2 = VectorMinus(p4, p3);
double eq1nCoeff = (d1.x * d2.x) + (d1.y * d2.y) + (d1.z * d2.z);
double eq1mCoeff = (-(powf(d1.x, 2)) - (powf(d1.y, 2)) - (powf(d1.z, 2)));
double eq1Const = ((d1.x * p3.x) - (d1.x * p1.x) + (d1.y * p3.y) - (d1.y * p1.y) + (d1.z * p3.z) - (d1.z * p1.z));
double eq2nCoeff = ((powf(d2.x, 2)) + (powf(d2.y, 2)) + (powf(d2.z, 2)));
double eq2mCoeff = -(d1.x * d2.x) - (d1.y * d2.y) - (d1.z * d2.z);
double eq2Const = ((d2.x * p3.x) - (d2.x * p1.x) + (d2.y * p3.y) - (d2.y * p2.y) + (d2.z * p3.z) - (d2.z * p1.z));
double M[2][3] = { { eq1nCoeff, eq1mCoeff, -eq1Const }, { eq2nCoeff, eq2mCoeff, -eq2Const } };
int rowCount = 2;
// pivoting
for (int col = 0; col + 1 < rowCount; col++) if (M[col, col] == 0)
// check for zero coefficients
{
// find non-zero coefficient
int swapRow = col + 1;
for (; swapRow < rowCount; swapRow++) if (M[swapRow, col] != 0) break;
if (M[swapRow, col] != 0) // found a non-zero coefficient?
{
// yes, then swap it with the above
double tmp[2];
for (int i = 0; i < rowCount + 1; i++)
{
tmp[i] = M[swapRow][i];
M[swapRow][i] = M[col][i];
M[col][i] = tmp[i];
}
}
else
{
std::cout << "\n the matrix has no unique solution";
return; // no, then the matrix has no unique solution
}
}
// elimination
for (int sourceRow = 0; sourceRow + 1 < rowCount; sourceRow++)
{
for (int destRow = sourceRow + 1; destRow < rowCount; destRow++)
{
double df = M[sourceRow][sourceRow];
double sf = M[destRow][sourceRow];
for (int i = 0; i < rowCount + 1; i++)
M[destRow][i] = M[destRow][i] * df - M[sourceRow][i] * sf;
}
}
// back-insertion
for (int row = rowCount - 1; row >= 0; row--)
{
double f = M[row][row];
if (f == 0) return;
for (int i = 0; i < rowCount + 1; i++) M[row][i] /= f;
for (int destRow = 0; destRow < row; destRow++)
{
M[destRow][rowCount] -= M[destRow][row] * M[row][rowCount]; M[destRow][row] = 0;
}
}
double n = M[0][2];
double m = M[1][2];
Vector3D i1 = { p1.x + (m * d1.x), p1.y + (m * d1.y), p1.z + (m * d1.z) };
Vector3D i2 = { p3.x + (n * d2.x), p3.y + (n * d2.y), p3.z + (n * d2.z) };
Vector3D i1Clamped = ClampPointToLine(i1, p1, p2);
Vector3D i2Clamped = ClampPointToLine(i2, p3, p4);
ClosestPointOnLineP1P2 = i1Clamped;
ClosestPointOnLineP3P4 = i2Clamped;
return;
}
2条答案
按热度按时间n3ipq98p1#
你的问题是要找出两条线段
AB
和CD
之间的最短连线P1P2
,我们定义l1
为通过点A
和B
的直线,l2
为通过点C
和D
的直线。您可以将此问题分解为多个子任务:
l1
和l2
之间的最短连接。A
、B
中的任一个到线段CD
的最短连接(对于C
、D
到线段AB
同样)。让我们从第一个子任务开始,通过
A
和B
的l1
线可以用一个标量来参数化,比如sc
,与方向矢量
u=(B-A)
。因此,我们还有l1(0) = A
和l(1) = B
。现在,我们要找出这条线与另一条经过C
和D
的线之间的最小距离,即与
v = D-C
类似,我们有l2(0) = C
和l(1) = D
。现在,我们定义也就是两条直线之间距离平方的一半。如果我们现在想最小化这个函数,我们需要满足
你会发现
引入
w=A-C
并按向量和矩阵排列,得到:线性系统的解为
result = -m^(⁻1)* rhs
,其中m^(⁻1)
是m
的倒数。如果a
和c
小于0或大于1,则直线的最近点位于线段AB
和CD
之外。您也可以返回这些值。第二个子任务与此问题密切相关,但我们将
直接产生
如果
sc < 0
,我们设置sc = 0
,如果sc > 1
,我们设置sc = 1
(对于tc
也是如此),以获得线段上的点。这是我从这里得到并修改的实现。首先,我们定义一些辅助对象,即向量和一些基本的数学关系。
此函数执行实际工作:
用法:
这个打印
live demo .
请注意,此代码仍需要更多测试。
eyh26e7m2#
下面是@StefanKssmr的代码的“压缩”版本,它是Here,这个“压缩”版本可以很容易地移植到OpenCL
非常感谢@StefanKssmr发布正确答案,