我的交集检查算法有什么问题?
What is wrong with my intersection checking algorithm?
我知道有很多网站都解释了如何检查两条线的交点,但我发现为这样一个简单的数学任务复制和粘贴代码是非常无聊的。我越让我沮丧的是我无法让我的代码工作。我知道 "What is wrong in my code?" 的问题很愚蠢,但我不知道我的数学/代码到底出了什么问题,而且我的代码记录得很好(除了公认的错误变量命名),所以我想应该有对其背后的数学感兴趣的人:
bool segment::checkforIntersection(QPointF a, QPointF b) { //line 1: a+bx, line 2: c+dx, note that a and c are called offset and bx and dx are called gradients in this code
QPointF bx = b-a;
double firstGradient = bx.y() / bx.x(); //gradient of line 1
//now we have to calculate the offset of line 1: we have b from a+bx. Since QPointF a is on that line, it is:
//a + b * a.x = a.y with a as free variable, which yields a = a.y - b*a.x.
//One could also use the second point b for this calculation.
double firstOffset = a.y() - firstGradient * a.x();
double secondGradient, secondOffset;
for (int i = 0; i < poscount-3; i++) { //we dont check with the last line, because that could be the same line, as the one that emited intersection checking
QPointF c = pos[i];
QPointF d = pos[i+1];
QPointF dx = d-c;
secondGradient = dx.y() / dx.x(); //same formula as above
secondOffset = c.y() - secondGradient * c.x();
//a+bx=c+dx <=> a-c = (d-b)x <=> (a-c)/(d-b) = x
double x = (firstOffset - secondOffset) / (secondGradient - firstGradient);
//we have to check, if those lines intersect with a x \in [a.x,b.x] and x \in [c.x,d.x]. If this is the case, we have a collision
if (x >= a.x() && x <= b.x() && x >= c.x() && x <= d.x()) {
return true;
}
}
return false;
}
那么这是做什么的,它有 4 个点 a、b、c、d(第 1 行:a--b,第 2 行:c--d)(忽略 for 循环),它们具有绝对 x 和值。首先,它通过计算 deltay/deltax 来计算线条的梯度。然后它通过使用点 a(或分别为 c)在线上这一事实来计算偏移量。通过这种方式,我们将这 4 个点转换为这些直线的数学表示形式,如等式 a+bx,其中 x 为 0 表示我们在第一个点 (a / c),x 为 1 表示我们在第二个点 (b/d).接下来我们计算这两条线的交点(基础代数)。之后我们检查交集的 x 值是否有效。据我了解,这一切都是正确的。有人看到错误了吗?
根据经验检查这是不正确的。该代码没有给出任何误报(说有交集,但实际上没有),但它给出了假否定(说没有交集,但实际上有)。所以说有交集是正确的,但是如果说没有交集,你就不能一直相信我的算法。
再次,我在网上查了一下,但是算法不同(有一些方向技巧和其他东西),我只是想提出我自己的算法,如果有人能提供帮助,我会很高兴。 :)
编辑:这是一个最小的可重现的不工作示例,这次没有 Qt,但只有 C++:
#include <iostream>
#include <math.h>
using namespace std;
class Point {
private:
double xval, yval;
public:
// Constructor uses default arguments to allow calling with zero, one,
// or two values.
Point(double x = 0.0, double y = 0.0) {
xval = x;
yval = y;
}
// Extractors.
double x() { return xval; }
double y() { return yval; }
Point sub(Point b)
{
return Point(xval - b.xval, yval - b.yval);
}
};
bool checkforIntersection(Point a, Point b, Point c, Point d) { //line 1: a+bx, line 2: c+dx, note that a and c are called offset and bx and dx are called gradients in this code
Point bx = b.sub(a);
double firstGradient = bx.y() / bx.x(); //gradient of line 1
//now we have to calculate the offset of line 1: we have b from a+bx. Since Point a is on that line, it is:
//a + b * a.x = a.y with a as free variable, which yields a = a.y - b*a.x.
//One could also use the second point b for this calculation.
double firstOffset = a.y() - firstGradient * a.x();
double secondGradient, secondOffset;
Point dx = d.sub(c);
secondGradient = dx.y() / dx.x(); //same formula as above
secondOffset = c.y() - secondGradient * c.x();
//a+bx=c+dx <=> a-c = (d-b)x <=> (a-c)/(d-b) = x
double x = (firstOffset - secondOffset) / (secondGradient - firstGradient);
//we have to check, if those lines intersect with a x \in [a.x,b.x] and x \in [c.x,d.x]. If this is the case, we have a collision
if (x >= a.x() && x <= b.x() && x >= c.x() && x <= d.x()) {
return true;
}
return false;
}
int main(int argc, char const *argv[]) {
if (checkforIntersection(Point(310.374,835.171),Point(290.434,802.354), Point(333.847,807.232), Point(301.03,827.172)) == true) {
cout << "These lines do intersect so I should be printed out\n";
} else {
cout << "The algorithm does not work, so instead I do get printed out\n";
}
return 0;
}
因此,作为示例,我采用了点 ~ (310,835) -- (290,802) 和 (333,807) -- (301,827)。这些线相交:
\documentclass[crop,tikz]{standalone}
\begin{document}
\begin{tikzpicture}[x=.1cm,y=.1cm]
\draw (310,835) -- (290,802);
\draw (333,807) -- (301,827);
\end{tikzpicture}
\end{document}
Proof of intersection
然而,当 运行 上面的 C++ 代码时,它说它们不相交
First it calculates the gradient of the lines by calculating deltay/deltax.
当 deltax 非常接近于零时会发生什么?
看,你正在做的是让自己暴露在病态的情况下 - 在计算几何方面总是害怕分裂和与 0.0
的直接比较。
备选方案:
- 两条线不平行就会相交
- 如果两条 distinct 条线的定义向量的叉积为零,则它们将平行。
您的 (a,b) x (c,d) = (ax-bx)*(cy-dy)-(ay-by)*(cx-dx)
的叉积 - 如果这足够接近于零,则出于所有实际目的,您的线之间没有交点(交点太远并不重要)。
现在,剩下的就是:
- 在计算叉积之前需要进行 "are those line distinct?" 测试。更重要的是,您将需要处理退化的情况(一条或两条线通过重合端减少到一个点 - 如
a==b
and/or c==d
)
- 如果您不对定义向量进行归一化,"close enough to zero" 测试将不明确 - 想象一个 1 光秒长度的向量定义第一行,另一个 1 秒差距向量(什么测试 'proximity to zero' 你应该在这种情况下使用吗?)
要对两个向量进行归一化,只需应用除法... (你说除法?我已经害怕得发抖了) ...嗯..我说的是将得到的叉积除以hypot(ax-bx, ay-by)*hypot(cx-dx,cy-dy)
(你明白为什么退化案例需要提前处理了吗?)
- 归一化后,再一次,什么是对结果叉积的良好 'proximity to zero' 测试?好吧,我想我可以再进行一个小时左右的分析(例如,与 {a,b,c,d} 多边形的范围相比,交叉点有多远),但是...因为交叉-两个酉向量的乘积(归一化后)是
sin(angle-between-versors)
,你可以用你的常识说'如果角度小于 1 度,这是否足以认为两条线平行?不? 1 角秒怎么样?
(你可以说我是书呆子,但术语很重要)
如果你想看看直线段是否相交,那么依靠你的两个段的参数表示,在两个参数中求解系统,看看是否两个两个参数的解都在 [0,1]
范围内。
分段 [a
、b
]、组件方面的参数表示
{x, y}(t) = {(1-t)*ax+t*bx, (1-t)*ay+t*by}
with t
in the [0,1]
range
快速检查 - 在 t=0
你得到 a,在 t=1
你得到 b,表达式在 t
中是线性的,所以你知道了。
因此,您的 (a,b) (c,d) 交集问题变为:
// for some t1 and t2, the x coordinate is the same
(1-t1)*ax+t*bx=(1-t2)*cx+t2*dx;
(1-t1)*ay+t*by=(1-t2)*cy+t2*dy; // and so is the y coordinate
求解 t1
和 t2
中的系统。如果t1
在[0,1]
范围内,那么交集在a
和b
之间,t2
对于c
也是如此和 d
.
它留作 reader 的练习,研究在以下条件下会对系统产生什么影响,以及应该为稳健的算法实施哪些检查:
- 片段退化 - 一个或两个片段的重合末端
- 具有非空白重叠的共线段。有一个重叠点的特殊情况(必要的,那个点将是末端之一)
- 没有重叠的共线段
- 平行线段
我知道有很多网站都解释了如何检查两条线的交点,但我发现为这样一个简单的数学任务复制和粘贴代码是非常无聊的。我越让我沮丧的是我无法让我的代码工作。我知道 "What is wrong in my code?" 的问题很愚蠢,但我不知道我的数学/代码到底出了什么问题,而且我的代码记录得很好(除了公认的错误变量命名),所以我想应该有对其背后的数学感兴趣的人:
bool segment::checkforIntersection(QPointF a, QPointF b) { //line 1: a+bx, line 2: c+dx, note that a and c are called offset and bx and dx are called gradients in this code
QPointF bx = b-a;
double firstGradient = bx.y() / bx.x(); //gradient of line 1
//now we have to calculate the offset of line 1: we have b from a+bx. Since QPointF a is on that line, it is:
//a + b * a.x = a.y with a as free variable, which yields a = a.y - b*a.x.
//One could also use the second point b for this calculation.
double firstOffset = a.y() - firstGradient * a.x();
double secondGradient, secondOffset;
for (int i = 0; i < poscount-3; i++) { //we dont check with the last line, because that could be the same line, as the one that emited intersection checking
QPointF c = pos[i];
QPointF d = pos[i+1];
QPointF dx = d-c;
secondGradient = dx.y() / dx.x(); //same formula as above
secondOffset = c.y() - secondGradient * c.x();
//a+bx=c+dx <=> a-c = (d-b)x <=> (a-c)/(d-b) = x
double x = (firstOffset - secondOffset) / (secondGradient - firstGradient);
//we have to check, if those lines intersect with a x \in [a.x,b.x] and x \in [c.x,d.x]. If this is the case, we have a collision
if (x >= a.x() && x <= b.x() && x >= c.x() && x <= d.x()) {
return true;
}
}
return false;
}
那么这是做什么的,它有 4 个点 a、b、c、d(第 1 行:a--b,第 2 行:c--d)(忽略 for 循环),它们具有绝对 x 和值。首先,它通过计算 deltay/deltax 来计算线条的梯度。然后它通过使用点 a(或分别为 c)在线上这一事实来计算偏移量。通过这种方式,我们将这 4 个点转换为这些直线的数学表示形式,如等式 a+bx,其中 x 为 0 表示我们在第一个点 (a / c),x 为 1 表示我们在第二个点 (b/d).接下来我们计算这两条线的交点(基础代数)。之后我们检查交集的 x 值是否有效。据我了解,这一切都是正确的。有人看到错误了吗?
根据经验检查这是不正确的。该代码没有给出任何误报(说有交集,但实际上没有),但它给出了假否定(说没有交集,但实际上有)。所以说有交集是正确的,但是如果说没有交集,你就不能一直相信我的算法。
再次,我在网上查了一下,但是算法不同(有一些方向技巧和其他东西),我只是想提出我自己的算法,如果有人能提供帮助,我会很高兴。 :)
编辑:这是一个最小的可重现的不工作示例,这次没有 Qt,但只有 C++:
#include <iostream>
#include <math.h>
using namespace std;
class Point {
private:
double xval, yval;
public:
// Constructor uses default arguments to allow calling with zero, one,
// or two values.
Point(double x = 0.0, double y = 0.0) {
xval = x;
yval = y;
}
// Extractors.
double x() { return xval; }
double y() { return yval; }
Point sub(Point b)
{
return Point(xval - b.xval, yval - b.yval);
}
};
bool checkforIntersection(Point a, Point b, Point c, Point d) { //line 1: a+bx, line 2: c+dx, note that a and c are called offset and bx and dx are called gradients in this code
Point bx = b.sub(a);
double firstGradient = bx.y() / bx.x(); //gradient of line 1
//now we have to calculate the offset of line 1: we have b from a+bx. Since Point a is on that line, it is:
//a + b * a.x = a.y with a as free variable, which yields a = a.y - b*a.x.
//One could also use the second point b for this calculation.
double firstOffset = a.y() - firstGradient * a.x();
double secondGradient, secondOffset;
Point dx = d.sub(c);
secondGradient = dx.y() / dx.x(); //same formula as above
secondOffset = c.y() - secondGradient * c.x();
//a+bx=c+dx <=> a-c = (d-b)x <=> (a-c)/(d-b) = x
double x = (firstOffset - secondOffset) / (secondGradient - firstGradient);
//we have to check, if those lines intersect with a x \in [a.x,b.x] and x \in [c.x,d.x]. If this is the case, we have a collision
if (x >= a.x() && x <= b.x() && x >= c.x() && x <= d.x()) {
return true;
}
return false;
}
int main(int argc, char const *argv[]) {
if (checkforIntersection(Point(310.374,835.171),Point(290.434,802.354), Point(333.847,807.232), Point(301.03,827.172)) == true) {
cout << "These lines do intersect so I should be printed out\n";
} else {
cout << "The algorithm does not work, so instead I do get printed out\n";
}
return 0;
}
因此,作为示例,我采用了点 ~ (310,835) -- (290,802) 和 (333,807) -- (301,827)。这些线相交:
\documentclass[crop,tikz]{standalone}
\begin{document}
\begin{tikzpicture}[x=.1cm,y=.1cm]
\draw (310,835) -- (290,802);
\draw (333,807) -- (301,827);
\end{tikzpicture}
\end{document}
Proof of intersection 然而,当 运行 上面的 C++ 代码时,它说它们不相交
First it calculates the gradient of the lines by calculating deltay/deltax.
当 deltax 非常接近于零时会发生什么?
看,你正在做的是让自己暴露在病态的情况下 - 在计算几何方面总是害怕分裂和与 0.0
的直接比较。
备选方案:
- 两条线不平行就会相交
- 如果两条 distinct 条线的定义向量的叉积为零,则它们将平行。
您的 (a,b) x (c,d) = (ax-bx)*(cy-dy)-(ay-by)*(cx-dx)
的叉积 - 如果这足够接近于零,则出于所有实际目的,您的线之间没有交点(交点太远并不重要)。
现在,剩下的就是:
- 在计算叉积之前需要进行 "are those line distinct?" 测试。更重要的是,您将需要处理退化的情况(一条或两条线通过重合端减少到一个点 - 如
a==b
and/orc==d
) - 如果您不对定义向量进行归一化,"close enough to zero" 测试将不明确 - 想象一个 1 光秒长度的向量定义第一行,另一个 1 秒差距向量(什么测试 'proximity to zero' 你应该在这种情况下使用吗?)
要对两个向量进行归一化,只需应用除法... (你说除法?我已经害怕得发抖了) ...嗯..我说的是将得到的叉积除以hypot(ax-bx, ay-by)*hypot(cx-dx,cy-dy)
(你明白为什么退化案例需要提前处理了吗?) - 归一化后,再一次,什么是对结果叉积的良好 'proximity to zero' 测试?好吧,我想我可以再进行一个小时左右的分析(例如,与 {a,b,c,d} 多边形的范围相比,交叉点有多远),但是...因为交叉-两个酉向量的乘积(归一化后)是
sin(angle-between-versors)
,你可以用你的常识说'如果角度小于 1 度,这是否足以认为两条线平行?不? 1 角秒怎么样?
(你可以说我是书呆子,但术语很重要)
如果你想看看直线段是否相交,那么依靠你的两个段的参数表示,在两个参数中求解系统,看看是否两个两个参数的解都在 [0,1]
范围内。
分段 [a
、b
]、组件方面的参数表示
{x, y}(t) = {(1-t)*ax+t*bx, (1-t)*ay+t*by}
witht
in the[0,1]
range
快速检查 - 在 t=0
你得到 a,在 t=1
你得到 b,表达式在 t
中是线性的,所以你知道了。
因此,您的 (a,b) (c,d) 交集问题变为:
// for some t1 and t2, the x coordinate is the same
(1-t1)*ax+t*bx=(1-t2)*cx+t2*dx;
(1-t1)*ay+t*by=(1-t2)*cy+t2*dy; // and so is the y coordinate
求解 t1
和 t2
中的系统。如果t1
在[0,1]
范围内,那么交集在a
和b
之间,t2
对于c
也是如此和 d
.
它留作 reader 的练习,研究在以下条件下会对系统产生什么影响,以及应该为稳健的算法实施哪些检查:
- 片段退化 - 一个或两个片段的重合末端
- 具有非空白重叠的共线段。有一个重叠点的特殊情况(必要的,那个点将是末端之一)
- 没有重叠的共线段
- 平行线段