我想检查三次贝塞尔曲线是否是另一个贝塞尔曲线的子曲线。
我想我基本上了解如何做到这一点,将贝塞尔曲线表示为 x 和 y 中的两个三次方,然后测试三次方是否相互缩放或平移。如果缩放和平移匹配,则表明曲线是同一曲线的子段,并在曲线 As 空间中为我们提供曲线 B 的 t0 素数和 t1 素数。
但我不太清楚如何检查三次方是否相等。
根据以下评论回答:
假设我们采用贝塞尔曲线,并使用 de Casteljau 算法将其拆分。显然结果是原始曲线的很多子曲线。问题是如何返回并恢复 t 值,以及曲线是同一曲线的一部分这一事实,仅给定它们的 4 个控制点
简短的回答:除非你有一台无限精密的机器,否则你不能。
所以我们坚持使用“错误阈值”测试。给定主曲线 A 和“希望是子曲线”曲线 B,如果 B是A 的子曲线,则遍历需要为真的事物:
如果所有这些都成立,我们可以合理地确定 B 是 A 的子曲线。
此外,由于我们需要得出t
值来检查一个点是否位于 A 上,以及 A 在该点的导数是多少,我们已经知道t
定义 A 上映射到完整曲线 B 的区间的值。
这是工作代码。(你可以很容易地找到三次根查找器)
/*
A = p3 + 3.0 * p1 - 3.0 * p2 - p0;
B = 3.0 * p0 - 6.0 * p1 + 3.0 * p2;
C = 3.0 * p1 - 3.0 * p0;
D = p0;
*/
bool CurveIsSubCurve(BezierCurve bez, BezierCurve sub, double epsilon, double *t)
{
int Nr;
double tcand[6];
int i, ii;
double ts[6], te[6];
int Ns = 0;
int Ne = 0;
Vector2 p;
/*
Take two bites at the cherry. The points may have slight errors, and a small error in x or y could represent a big error in
t. However with any luck either x or y will be close
*/
Nr = cubic_roots(bez.Ax(), bez.Bx(), bez.Cx(), bez.Dx() - sub.P0().x, tcand);
Nr += cubic_roots(bez.Ay(), bez.By(), bez.Cy(), bez.Dy() - sub.P0().y, tcand + Nr);
for(i=0;i<Nr;i++)
{
p = bez.Eval(tcand[i]);
if(fabs(p.x - sub.P0().x) < epsilon && fabs(p.y - sub.P0().y) < epsilon)
{
ts[Ns++] = tcand[i];
}
}
/* same thing of sub curve end point */
Nr = cubic_roots(bez.Ax(), bez.Bx(), bez.Cx(), bez.Dx() - sub.P3().x, tcand);
Nr += cubic_roots(bez.Ay(), bez.By(), bez.Cy(), bez.Dy() - sub.P3().y, tcand + Nr);
for(i=0;i<Nr;i++)
{
p = bez.Eval(tcand[i]);
if(fabs(p.x - sub.P3().x) < epsilon && fabs(p.y - sub.P3().y) < epsilon)
{
te[Ne++] = tcand[i];
}
}
/* do an all by all to get matches (Ns, Ne will be small, but if
we have a degenerate, i.e. a loop, the loop intersection point is
where the mother curve is quite likely to be cut, so test everything*/
for(i = 0; i < Ns; i++)
{
double s,d;
double Ax, Bx, Cx, Dx;
double Ay, By, Cy, Dy;
for(ii=0;ii<Ne;ii++)
{
s = (te[ii] - ts[i]);
d = ts[i];
/* now substitute back */
Ax = bez.Ax() *s*s*s;
Bx = bez.Ax() *2*s*s*d + bez.Ax()*s*s*d + bez.Bx()*s*s;
Cx = bez.Ax()*s*d*d + bez.Ax()*2*s*d*d + bez.Bx()*2*s*d + bez.Cx() * s;
Dx = bez.Ax() *d*d*d + bez.Bx()*d*d + bez.Cx()*d + bez.Dx();
Ay = bez.Ay() *s*s*s;
By = bez.Ay() *2*s*s*d + bez.Ay()*s*s*d + bez.By()*s*s;
Cy = bez.Ay()*s*d*d + bez.Ay()*2*s*d*d + bez.By()*2*s*d + bez.Cy() * s;
Dy = bez.Ay() *d*d*d + bez.By()*d*d + bez.Cy()*d + bez.Dy();
if(fabs(Ax - sub.Ax()) < epsilon && fabs(Bx - sub.Bx()) < epsilon &&
fabs(Cx - sub.Cx()) < epsilon && fabs(Dx - sub.Dx()) < epsilon &&
fabs(Ay - sub.Ay()) < epsilon && fabs(By - sub.By()) < epsilon &&
fabs(Cy - sub.Cy()) < epsilon && fabs(Dy - sub.Dy()) < epsilon)
{
if(t)
{
t[0] = ts[i];
t[1] = te[ii];
}
return true;
}
}
}
return false;
}