作为一个兼职项目,我正在研究一些几何实用程序,遇到了一个相对简单的问题,似乎有一个不那么简单的解决方案。
问题在于EPSILON对于这个问题来说太小了。为了观察两个三角形是否相似,我为每个三角形计算出余弦形式的3个内角,然后对它们进行排序。然后我测试Math.abs(t1[0]-t2[0]) < EPSILON
,其中t1是一个三角形,t2是另一个,每个三角形都包含三个角。
在我知道相似的三角形上,我的失败率大约为20%-80%。当我将EPSILON设置为一个更大的值时,例如,仍然是一个非常小的0.0000001,没有出现故障(在我让测试运行的时间内)。
下面是提取的相关三角形函数,我还将测试代码作为演示包含在其中。单击按钮,它将运行测试并显示结果。三角形是随机生成的。每隔一段时间就会创建两个类似的三角形,其中大约一半是精确的副本,其余的是副本,但缩放、镜像、旋转和矢量顺序打乱,同时仍然保持相似性
我想知道如何计算一个合理的EPSILON,以减少不正确的结果,但保持系统尽可能准确?
尽管也有可能在测试代码中出现另一个错误,我将继续检查。
const EPSILON = Number.EPSILON
function Triangle(p1,p2,p3){
this.p1 = p1;
this.p2 = p2;
this.p3 = p3;
}
Triangle.prototype.p1 = undefined;
Triangle.prototype.p2 = undefined;
Triangle.prototype.p3 = undefined;
Triangle.prototype.isSimilar = function(triangle){
var a1,b1,c1,a2,b2,c2,aa1,bb1,cc1,aa2,bb2,cc2; //
var t1 = [];
var t2 = [];
var sortF = function(a,b){ return a-b };
// get the length squared and length of each side
a1 = Math.sqrt(aa1 = Math.pow(this.p1.x - this.p2.x, 2) + Math.pow(this.p1.y - this.p2.y, 2));
b1 = Math.sqrt(bb1 = Math.pow(this.p2.x - this.p3.x, 2) + Math.pow(this.p2.y - this.p3.y, 2));
c1 = Math.sqrt(cc1 = Math.pow(this.p3.x - this.p1.x, 2) + Math.pow(this.p3.y - this.p1.y, 2));
a2 = Math.sqrt(aa2 = Math.pow(triangle.p1.x - triangle.p2.x, 2) + Math.pow(triangle.p1.y - triangle.p2.y, 2));
b2 = Math.sqrt(bb2 = Math.pow(triangle.p2.x - triangle.p3.x, 2) + Math.pow(triangle.p2.y - triangle.p3.y, 2));
c2 = Math.sqrt(cc2 = Math.pow(triangle.p3.x - triangle.p1.x, 2) + Math.pow(triangle.p3.y - triangle.p1.y, 2));
// get the cosin of each angle for both triangle
t1[0] = (cc1 - (aa1 + bb1)) / (-2 * a1 * b1);
t1[1] = (aa1 - (cc1 + bb1)) / (-2 * c1 * b1);
t1[2] = (bb1 - (cc1 + aa1)) / (-2 * c1 * a1);
t2[0] = (cc2 - (aa2 + bb2)) / (-2 * a2 * b2);
t2[1] = (aa2 - (cc2 + bb2)) / (-2 * c2 * b2);
t2[2] = (bb2 - (cc2 + aa2)) / (-2 * c2 * a2);
t1.sort(sortF);
t2.sort(sortF);
if(Math.abs(t1[0] - t2[0]) < EPSILON && Math.abs(t1[1] - t2[1]) < EPSILON && Math.abs(t1[2] - t2[2]) < EPSILON){
return true;
}
return false;
}
function Vec(x,y){
this.x = x;
this.y = y;
}
Vec.prototype.x = undefined;
Vec.prototype.y = undefined;
更新
更多信息。
使用角度余弦的类似三角形失败EPSILON:2.220446049250313e-16
Failed Triangles ID : 94
Method : compare cosine of angles
Both Compare T1 to T2 and T2 to T1 failed
Both Triangles known to be similare
Triangle 1
p1.x = -149241116087155.97;
p1.y = -1510074922190599.8;
p2.x = -2065214078816255.8;
p2.y = 6756872141691895;
p3.x = -7125027429739231;
p3.y = -5622578541875555;
Triangle 2
p1.x = -307440480802857.2;
p1.y = -404929352172871.56;
p2.x = -3020163594243123;
p2.y = -355583557775981.75;
p3.x = 595422457974710.8;
p3.y = 2291176238828451.5;
Compare T1 to T2 Result : false
Computed values
Triangle 1 length of side and square length
length a : 8486068945686473 squared : 7.201336615094433e+31
length b : 13373575078230092 squared : 1.78852510373057e+32
length c : 8097794805726894 squared : 6.557428071565746e+31
Unsorted cosines C is angle opposite side c
cosine C : 0.8163410767815653
cosine A : 0.7960251614312384
cosine B : -0.30024590551189423
ratio a : undefined
ratio b : undefined
ratio c : undefined
Triangle2
length a : 2713171888697380.5 squared : 7.36130169761771e+30
length b : 4480825808030667.5 squared : 2.0077799921913682e+31
length c : 2843263414467020.5 squared : 8.08414684404666e+30
Unsorted cosines C is angle opposite side c
cosine C : 0.7960251614312384
cosine A : 0.8163410767815651
cosine B : -0.3002459055118942
Compare T2 to T1 Result : false
Triangle1
Computed values
Triangle 1 length of side and square length
length a : 2713171888697380.5 squared : 7.36130169761771e+30
length b : 4480825808030667.5 squared : 2.0077799921913682e+31
length c : 2843263414467020.5 squared : 8.08414684404666e+30
Unsorted cosines C is angle opposite side c
cosine a : 0.7960251614312384
cosine b : 0.8163410767815651
cosine c : -0.3002459055118942
ratio a : undefined
ratio b : undefined
ratio c : undefined
Triangle2
length a : 8486068945686473 squared : 7.201336615094433e+31
length b : 13373575078230092 squared : 1.78852510373057e+32
length c : 8097794805726894 squared : 6.557428071565746e+31
cosine a : 0.8163410767815653
cosine b : 0.7960251614312384
cosine c : -0.30024590551189423
更新2
结果输出和错误修复(抱歉@lhf我没有sqrt epsilon我仍然使用原始常量)
这显示了对同一组三角形的测试结果。不一致意味着将三角形1与三角形2进行比较与将三角形2与三角形1进行比较的结果不同。不正确表示两个已知的相似三角形失败,不正确不一致表示两个未知的相似三角形在一次测试中失败,并通过了另一次测试。
使用长度比给出了最差的结果,但使用余弦并没有好到哪里去,除了使用长度比比较t1与t2和t2与t1之间具有非常高的不一致率的不正确不一致性类似三角形。但这是有道理的,因为根据测试的顺序,比率的大小会有很大的变化。
正如您所看到的,使用EPSILON的平方根完全消除了两种方法的错误。
如果lhf希望将sqrt(epsilon)注释作为答案,我将接受它作为答案。感谢大家的投入,感谢Salix ,我还有更多的阅读
======================================
Default EPSILON : 2.220446049250313e-16
======================================
Via cosine of angles
All Inconsistency failed : 0 of 10000
Similar Incorrect failed : 1924 of 5032
Similar Incorrect Inconsistency failed : 0 of 5032
======================================
Via ratio of lengths
All Inconsistency failed : 1532 of 10000
Similar Incorrect failed : 2082 of 5032
Similar Incorrect Inconsistency failed : 1532 of 5032
======================================
Squaring EPSILON : 1.4901161193847656e-8
======================================
Via cosine of angles
All Inconsistency failed : 0 of 10000
Similar Incorrect failed : 0 of 5032
Similar Incorrect Inconsistency failed : 0 of 5032
======================================
Via ratio of lengths
All Inconsistency failed : 0 of 10000
Similar Incorrect failed : 0 of 5032
Similar Incorrect Inconsistency failed : 0 of 5032
const EPSILON = Number.EPSILON
function Triangle(p1,p2,p3){
this.p1 = p1;
this.p2 = p2;
this.p3 = p3;
}
Triangle.prototype.p1 = undefined;
Triangle.prototype.p2 = undefined;
Triangle.prototype.p3 = undefined;
Triangle.prototype.isSimilar = function(triangle){
var a1,b1,c1,a2,b2,c2,aa1,bb1,cc1,aa2,bb2,cc2; //
var t1 = [];
var t2 = [];
var sortF = function(a,b){ return a-b };
// get the length squared and length of each side
a1 = Math.sqrt(aa1 = Math.pow(this.p1.x - this.p2.x, 2) + Math.pow(this.p1.y - this.p2.y, 2));
b1 = Math.sqrt(bb1 = Math.pow(this.p2.x - this.p3.x, 2) + Math.pow(this.p2.y - this.p3.y, 2));
c1 = Math.sqrt(cc1 = Math.pow(this.p3.x - this.p1.x, 2) + Math.pow(this.p3.y - this.p1.y, 2));
a2 = Math.sqrt(aa2 = Math.pow(triangle.p1.x - triangle.p2.x, 2) + Math.pow(triangle.p1.y - triangle.p2.y, 2));
b2 = Math.sqrt(bb2 = Math.pow(triangle.p2.x - triangle.p3.x, 2) + Math.pow(triangle.p2.y - triangle.p3.y, 2));
c2 = Math.sqrt(cc2 = Math.pow(triangle.p3.x - triangle.p1.x, 2) + Math.pow(triangle.p3.y - triangle.p1.y, 2));
// get the cosin of each angle for both triangle
t1[0] = (cc1 - (aa1 + bb1)) / (-2 * a1 * b1);
t1[1] = (aa1 - (cc1 + bb1)) / (-2 * c1 * b1);
t1[2] = (bb1 - (cc1 + aa1)) / (-2 * c1 * a1);
t2[0] = (cc2 - (aa2 + bb2)) / (-2 * a2 * b2);
t2[1] = (aa2 - (cc2 + bb2)) / (-2 * c2 * b2);
t2[2] = (bb2 - (cc2 + aa2)) / (-2 * c2 * a2);
t1.sort(sortF);
t2.sort(sortF);
if(Math.abs(t1[0] - t2[0]) < EPSILON && Math.abs(t1[1] - t2[1]) < EPSILON && Math.abs(t1[2] - t2[2]) < EPSILON){
return true;
}
return false;
}
function Vec(x,y){
this.x = x;
this.y = y;
}
Vec.prototype.x = undefined;
Vec.prototype.y = undefined;
var iterations = 1000; // number of tests
var presentSimilar = 1/2; // odds of similar triangle
var presentSuperSimilar = 1/2; // odds of triangles being identical
var presentInfinity = 0;//1/20; // odds of a divide by zero
var presentDegenerate = 0;//1/100; // odds of a degenerate triangle can be colinear or degenerate right triangle
var v; // temp for swap
var xdx,xdy,ydx,ydy; // vars for rotation
var copyVec = [["p1","p2","p3"],["p2","p3","p1"],["p3","p1","p2"]]; // pick a vec for selecting vecs
// the triangles for testing;
var tri1 = new Triangle(new Vec(0,0), new Vec(0,0), new Vec(0,0));
var tri2 = new Triangle(new Vec(0,0), new Vec(0,0), new Vec(0,0));
// max Random
function rMax(){
return ((Math.random()*2)-1) * Number.MAX_SAFE_INTEGER;
}
// rotate function
function rotate(vec){
var x = vec.x;
var y = vec.y;
vec.x = x * xdx + y * ydx;
vec.y = x * xdy + y * ydy;
};
function translateVec(vec,x,y){
vec.x += x;
vec.y += y;
}
function translateTriangle(tri,x,y){
translateVec(tri.p1);
translateVec(tri.p2);
translateVec(tri.p3);
}
// make infinite vec to simulate the result of a divide by zero
function doInfinity(vec){
if(Math.random() < presentInfinity){
if(Math.random() < 0.5){
vec.x = Infinity;
vec.y = Infinity;
}else{
vec.x = -Infinity;
vec.y = -Infinity;
}
}
}
// create a random vector;
function randomVec(vec){
vec.x = rMax();
vec.y = rMax();
doInfinity(vec);
}
// create a random triangle
function randomTriangle(tri){
var p,r;
randomVec(tri.p1);
randomVec(tri.p2);
randomVec(tri.p3);
if(Math.random() < presentDegenerate){
r = Math.random();
if(r < 1/3){ // Degenerate right triangle
p = copyVec[Math.floor(Math.random()*3)]; // get two vec to be at the same location
tri[p[0]].x = tri[p[1]].x;
tri[p[0]].y = tri[p[1]].y;
}else // Degenerate colinear triangle
if(r < 2/3){
p = copyVec[Math.floor(Math.random()*3)]; // get two vec to be at the same location
r = Math.random();
tri[p[0]].x = (tri[p[1]].x - tri[p[2]].x) * r + tri[p[2]].x;
tri[p[0]].y = (tri[p[1]].y - tri[p[2]].y) * r + tri[p[2]].y;
}else{ // degenerate undimentioned triangle. Has not area
tri.p1.x = tri.p2.x = tri.p3.x;
tri.p1.y = tri.p2.y = tri.p3.y;
}
}
}
function runTest(){
var result1,result2,mustBeSimilar;
var countSimilar = 0;
var countNorm = 0;
var error1 = 0;
var error2 = 0;
for(var i = 0; i < iterations; i ++){
randomTriangle(tri1);
if(Math.random() < presentSimilar){
mustBeSimilar = true;
countSimilar += 1;
tri2.p1.x = tri1.p1.x;
tri2.p1.y = tri1.p1.y;
tri2.p2.x = tri1.p2.x;
tri2.p2.y = tri1.p2.y;
tri2.p3.x = tri1.p3.x;
tri2.p3.y = tri1.p3.y;
if(Math.random() >= presentSuperSimilar){
if(Math.random() < 0.5){ // swap two
v = tri2.p1;
tri2.p1 = tri2.p2;
tri2.p2 = v;
}
if(Math.random() < 0.5){ // swap two
v = tri2.p2;
tri2.p2 = tri2.p3;
tri2.p3 = v;
}
if(Math.random() < 0.5){ // swap two
v = tri2.p1;
tri2.p1 = tri2.p3;
tri2.p3 = v;
}
// scale and or mirror the second triangle
v = Math.random() * 2 - 1;
tri2.p1.x *= v;
tri2.p1.y *= v;
tri2.p2.x *= v;
tri2.p2.y *= v;
tri2.p3.x *= v;
tri2.p3.y *= v;
// rotate the triangle
v = (Math.random()- 0.5) * Math.PI * 4;
ydy = xdx = Math.cos(v);
ydx = -(xdy = Math.sin(v));
rotate(tri2.p1);
rotate(tri2.p2);
rotate(tri2.p3);
}
}else{
randomTriangle(tri2);
mustBeSimilar = false;
}
countNorm += 1;
result1 = tri1.isSimilar(tri2);
result2 = tri2.isSimilar(tri1);
if(result1 !== result2){
error1 += 1;
}
if(mustBeSimilar && (!result1 || !result2)){
error2 += 1;
}
for(var j = 0; j < 10; j++){
translateTriangle(tri1,Math.random(),Math.random());
translateTriangle(tri2,Math.random(),Math.random());
if(mustBeSimilar){
countSimilar += 1;
}
countNorm += 1;
result1a = tri1.isSimilar(tri2);
result2a = tri2.isSimilar(tri1);
if(result1a !== result2a || result1 !== result1a || result2 !== result2a){
error1 += 1;
}
if(mustBeSimilar && (!result1a || !result2a)){
error2 += 1;
}
}
}
divResult1.textContent = "Inconsistancy result failed : "+error1 + " of "+ countNorm;
divResult2.textContent = "Incorrect result failed : "+error2 + " of "+ countSimilar
}
var button = document.createElement("input");
button.type = "button"
button.value = "Run test"
button.onclick = runTest;
var divResult1 = document.createElement("div");
var divResult2 = document.createElement("div");
document.body.appendChild(button);
document.body.appendChild(divResult1);
document.body.appendChild(divResult2);
遵循willyokka的评论。你可以看看是否有一个单一的比例因子。
// get the length squared and length of each side
a1 = Math.sqrt(...);
....
// Sort the values so a1 < b1 < c1, a2 < b2 < c2
if(b1 < a1) { tmp = b1; b1 = a1; a1 = tmp }
if(c1 < a1) { tmp = c1; c1 = a1; a1 = tmp }
if(c1 < b1) { tmp = c1; c1 = b1; b1 = tmp }
if(b2 < a2) { tmp = b2; b2 = a2; a2 = tmp }
if(c2 < a2) { tmp = c2; c2 = a2; a2 = tmp }
if(c2 < b2) { tmp = c2; c2 = b2; b2 = tmp }
// work out the scale factors
ka = a2 / a1;
kb = b2 / b1;
kc = c2 / c1;
if( abs(ka - kb) < epsilon && abs(kc - ka) < epsilon && abs(kc - kb) < epsilon )
// similar
else
// not similar
与其使用绝对ε,您可能需要一个在该值的x%以内的ε。使得如果ka-x%<kb<ka+x%,即(1-x/100)ka<kb<(1+x/100)ka。或者(1-x/100)<kb/ka<(1+x/100)或abs(kb/ka)<x/100.
对于这个问题,有一种更严格的统计方法。这将涉及到对我们所说的相似的含义进行相当精确的定义,并检查三角形的分布。统计形状分析是一个糟糕的起点。大卫·乔治·肯德尔研究了三角形的形状。
如果真的只是角度,你可以简单地计算两个三角形的三个内角。只需使用余弦,
cos(angle) = dot (normalized(edge[i]),normalized(edge[(i+1)%3])).
with edge[i] = p[(i+1)%3] - p[i].
每个三角形一和二都有三个cos(角度)。然后检查每个排列。三角形只有六种排列。(http://mathworld.wolfram.com/Permutation.html)
besterr = max;
for i=1..6 perm(i) in tri1
for j=1..6 perm(j) in tri2
err = 0
for k=1..3 angle
err += abs(angletri1[perm[i,k]] - angletri2[perm[j,k]])
if (err<besterr) besterr = err;
return besterr;
这是你预期的结果吗?我们当然可以做得更有效率。但这是蛮力测试算法。需要注意的一点是,它只适用于三角形——三角形中的任何顶点排列都是相同的三角形轮廓。对于更大的多边形来说,情况并非如此。
一旦成功,你就可以开始试验了。角度和cos(角度)的结果相同吗?对于err+=abs(d)和err+=d*d?你能通过角度排序来检查2个排列吗?记住三角形的角度总和(https://en.wikipedia.org/wiki/Sum_of_angles_of_a_triangle)。哪些计算是多余的?
最后:这真的是你想要的指标吗?两个绕组相反的三角形真的相似吗?一个大的和一个小的?