セグメントセグメントの交差点コードの精度問題

Question

複雑な自己交差ポリゴンをすべての交差点を識別して外側の周りを歩くことによって単純なポリゴンに変換するアルゴリズムをデバッグしています。

ストレステストを生成する一連のランダムデータを書きましたが、何千回も正しく作業した後、私のアルゴリズムが失敗するという面白い状況に遭遇しました。

double fRand(double fMin, double fMax) 
{ 
    double f = (double)rand()/RAND_MAX; // On my machine RAND_MAX = 2147483647 
    return fMin + f * (fMax - fMin); 
} 
// ... testing code below: 
srand(1); 
for(int j=3;j<5000;j+=5) { 
    std::vector<E_Point> geometry; 
    for(int i=0;i<j;i++) { 
     double radius = fRand(0.6,1.0); 
     double angle = fRand(0,2*3.1415926535); 
     E_Point pt(cos(angle),sin(angle)); 
     pt *= radius; 
     geometry.push_back(pt); // sending in a pile of shit 
    } 
    // run algorithm on this geometry 

これは、これがどのように関係しているか、私が今どこにいるかの概要です。私が放棄している詳細はたくさんあります。私は何をすることができたことは、私が使用しているセグメントのセグメントの交差点コードへのダウン狭い問題である

：起こっ

bool intersect(const E_Point& a0, const E_Point& a1, 
       const E_Point& b0, const E_Point& b1, 
       E_Point& intersectionPoint) { 

    if (a0 == b0 || a0 == b1 || a1 == b0 || a1 == b1) return false; 
    double x1 = a0.x; double y1 = a0.y; 
    double x2 = a1.x; double y2 = a1.y; 
    double x3 = b0.x; double y3 = b0.y; 
    double x4 = b1.x; double y4 = b1.y; 

    //AABB early exit 
    if (b2Max(x1,x2) < b2Min(x3,x4) || b2Max(x3,x4) < b2Min(x1,x2)) return false; 
    if (b2Max(y1,y2) < b2Min(y3,y4) || b2Max(y3,y4) < b2Min(y1,y2)) return false; 

    float ua = ((x4 - x3) * (y1 - y3) - (y4 - y3) * (x1 - x3)); 
    float ub = ((x2 - x1) * (y1 - y3) - (y2 - y1) * (x1 - x3)); 
    float denom = (y4 - y3) * (x2 - x1) - (x4 - x3) * (y2 - y1); 
    // check against epsilon (lowest normalized double value) 
    if (fabs(denom) < DBL_EPSILON) { 
     //Lines are too close to parallel to call 
     return false; 
    } 
    ua /= denom; 
    ub /= denom; 

    if ((0 < ua) && (ua < 1) && (0 < ub) && (ub < 1)) { 
     intersectionPoint.x = (x1 + ua * (x2 - x1)); 
     intersectionPoint.y = (y1 + ua * (y2 - y1)); 
     return true; 
    } 
    return false; 
} 

何が、私はこの機能はまったく同じ交差点を返す2つの交点を有することですのポイント値。垂直線は点 (0.3871953044519425, -0.91857980824611341), (0.36139704793723609, 0.91605957361605106)

点によって

緑色ほぼ水平ライン(0.8208980020500205, 0.52853407296583088), (0.36178501611208552, 0.88880385168617226)

とによって規定される

：ここに非常に関連するジオメトリを考慮してズームされます白い行by (0.36178501611208552, 0.88880385168617226), (-0.43211245441046209, 0.68034202227710472)

最後の2行は共通点を共有しています。

私の関数は、これらの交点の両方について(0.36178033094571277, 0.88880245640159794) の解を与えます（そのうちの1つは赤い点として表示されます）。

これが大きな問題である理由は、私の周辺アルゴリズムは各エッジの交点をソートすることに依存するからです。これらの交点は両方とも同じ正確な値を持つように計算されていたため、ソートはそれらを間違った方向に置きました。周囲のパスは上から来て、緑の線が残っているのではなく、左の白い線に沿っています。これは、もはや私の多角形の外周をたどっていないことを意味します。

この問題を解決するには、できることは多分ありますが、交差点リストをすべて検索して、他の点と比較して位置が等しいかどうかを確認する必要はありません。よりよい解決策は、私の交差点機能の精度を上げることです。

私は何を求めているのですか。なぜソリューションポイントが不正確なのですか？それはラインの1つがほぼ垂直なのでですか？最初に何らかの変換を実行する必要がありますか？両方の場合において、ほぼ垂直な線はa0およびa1として渡された。

更新：ねえ、これを見て：

TEST(intersection_precision_test) { 
    E_Point problem[] = { 

    {0.3871953044519425, -0.91857980824611341}, // 1559 
    {0.36139704793723609, 0.91605957361605106}, // 1560 
    {-0.8208980020500205, 0.52853407296583088}, // 1798 
    {0.36178501611208552, 0.88880385168617226}, // 1799 
    {-0.43211245441046209, 0.6803420222771047} // 1800 
    }; 
    std::cout.precision(16); 
    E_Point result; 
    intersect(problem[0],problem[1],problem[2],problem[3],result); 
    std::cout << "1: " << result << std::endl; 
    intersect(problem[0],problem[1],problem[3],problem[2],result); 
    std::cout << "2: " << result << std::endl; 
    intersect(problem[1],problem[0],problem[2],problem[3],result); 
    std::cout << "3: " << result << std::endl; 
    intersect(problem[1],problem[0],problem[3],problem[2],result); 
    std::cout << "4: " << result << std::endl; 

    intersect(problem[2],problem[3],problem[0],problem[1],result); 
    std::cout << "rev: " << result << std::endl; 
    intersect(problem[3],problem[2],problem[0],problem[1],result); 
    std::cout << "revf1: " << result << std::endl; 
    intersect(problem[2],problem[3],problem[1],problem[0],result); 
    std::cout << "revf2: " << result << std::endl; 
    intersect(problem[3],problem[2],problem[1],problem[0],result); 
    std::cout << "revfboth: " << result << std::endl; 
} 

出力：

Starting Test intersection_precision_test, at Polygon.cpp:1830 
1: <0.3617803309457128,0.8888024564015979> 
2: <0.3617803309457128,0.8888024564015979> 
3: <0.3617803314022162,0.8888024239374175> 
4: <0.3617803314022162,0.8888024239374175> 
rev: <0.3617803635476076,0.8888024344185281> 
revf1: <0.3617803313928456,0.8888024246235207> 
revf2: <0.3617803635476076,0.8888024344185281> 
revfboth: <0.3617803313928456,0.8888024246235207> 

は、私が実際に仮数ビットのうち、実行しているか、私は賢くアルゴリズムで有意に優れて行うことができますか？

ここで問題となるのは、頂点が本当に本当に本当に別の線の近くに設定されているかを簡単に判断できないということです。私はそれを動かすことも、完全にそれを噛んだりすることもないでしょう。

Answer 1

あなたがダブルスにあなたのフロート中間体（UA、UB、denom）を変更し、（分割後）UA値を印刷する場合は、あなたがこれらを取得します：

0x1.f864ab6b36458p-1 in the first case 
0x1.f864af01f2037p-1 in the second case 

私は六角でそれらを印刷しましたビットを見やすくする。これらの2つの値は、最初の22ビット（1.f864a + bとfの上位ビット）に一致します。浮動小数点数は23ビットの仮数しか持っていません！浮動小数点で中間体を計算すると、同じ答えに丸められることは驚くことではありません。

この場合、おそらく浮動小数点の代わりに倍精度を使って中間体を計算することで問題を回避できます。 （浮動小数点を使用している場合、double型の中間体を計算すると助けになるかどうかは分かりません）。しかし、垂直線分がさらに近づく場合2つの水平セグメントの交点までの距離は、2倍が提供できるよりもさらに高い精度を必要とする可能性がある。

垂直セグメントが水平セグメントの共有エンドポイントを正確に通過した場合はどうしますか？私はあなたがこのケースを正しく扱っていると思います。

あなたは倍精度で計算UB値（分割後）、見れば、あなたはこれらを取得する：

0.9999960389052315 
5.904388076838819e-06 

この交差点は、水平セグメントの共有エンドポイントに非常に、非常に近いことを意味しています。

私はあなたができると思います。交差点を計算するたびに、ubを見てください。十分に1に近い場合は、の交点に移動します。その後、最初の場所で正確なパススルーケースであるかのように交差点を処理します。実際にデータを変更する必要はありませんが、エンドポイントを共有している両方のセグメント（現在の交差テストと次のラインセグメント上のテスト）で、エンドポイントを移動したものとして扱う必要があります。