How can I determine whether a 2D Point is within a Polygon? Ask Question

How can I determine whether a 2D Point is within a Polygon? Ask Question

I'm trying to create a fast 2D point inside polygon algorithm, for use in hit-testing (e.g. Polygon.contains(p:Point)). Suggestions for effective techniques would be appreciated.

ベストアンサー1

For graphics, I'd rather not prefer integers. Many systems use integers for UI painting (pixels are ints after all), but macOS, for example, uses float for everything. macOS only knows points and a point can translate to one pixel, but depending on monitor resolution, it might translate to something else. On retina screens half a point (0.5/0.5) is pixel. Still, I never noticed that macOS UIs are significantly slower than other UIs. After all, 3D APIs (OpenGL or Direct3D) also work with floats and modern graphics libraries very often take advantage of GPU acceleration.

Now you said speed is your main concern, okay, let's go for speed. Before you run any sophisticated algorithm, first do a simple test. Create an axis aligned bounding box around your polygon. This is very easy, fast and can already save you a lot of calculations. How does that work? Iterate over all points of the polygon and find the min/max values of X and Y.

E.g. you have the points (9/1), (4/3), (2/7), (8/2), (3/6). This means Xmin is 2, Xmax is 9, Ymin is 1 and Ymax is 7. A point outside of the rectangle with the two edges (2/1) and (9/7) cannot be within the polygon.

// p is your point, p.x is the x coord, p.y is the y coord
if (p.x < Xmin || p.x > Xmax || p.y < Ymin || p.y > Ymax) {
    // Definitely not within the polygon!
}

This is the first test to run for any point. As you can see, this test is ultra fast but it's also very coarse. To handle points that are within the bounding rectangle, we need a more sophisticated algorithm. There are a couple of ways how this can be calculated. Which method works also depends on whether the polygon can have holes or will always be solid. Here are examples of solid ones (one convex, one concave):

穴のない多角形

And here's one with a hole:

穴のある多角形

The green one has a hole in the middle!

The easiest algorithm, that can handle all three cases above and is still pretty fast is named ray casting. The idea of the algorithm is pretty simple: Draw a virtual ray from anywhere outside the polygon to your point and count how often it hits a side of the polygon. If the number of hits is even, it's outside of the polygon, if it's odd, it's inside.

光線がポリゴンを切断する方法のデモンストレーション

The winding number algorithm would be an alternative, it is more accurate for points being very close to a polygon line but it's also much slower. Ray casting may fail for points too close to a polygon side because of limited floating point precision and rounding issues, but in reality that is hardly a problem, as if a point lies that close to a side, it's often visually not even possible for a viewer to recognize if it is already inside or still outside.

上記の境界ボックスはまだ残っていますね。覚えていますか? 境界ボックスの外側の点を選択し、それをレイの開始点として使用します。たとえば、点は(Xmin - e/p.y)確実にポリゴンの外側にあります。

しかし、 とは何でしょうかe? そうですね、e(実際にはイプシロン) は境界ボックスにいくらかのパディングを与えます。 前に述べたように、ポリゴンの線に近すぎるとレイ トレーシングは失敗します。 境界ボックスがポリゴンと等しくなる可能性があるので (ポリゴンが軸に沿った長方形の場合、境界ボックスはポリゴン自体と等しくなります!)、これを安全に行うにはいくらかのパディングが必要です。それだけです。 はどのくらいの大きさを選択すればよいでしょうかe? 大きすぎてはいけません。 描画に使用する座標系のスケールによって異なります。 ピクセル ステップ幅が 1.0 の場合は、1.0 を選択します (0.1 でも機能します)。

開始座標と終了座標を持つ光線が得られたので、問題は「ポリゴン内の点」から「光線がポリゴンの辺と交差する頻度」に移ります。したがって、以前のようにポリゴンの点だけで作業することはできず、実際の辺が必要になります。辺は常に 2 つの点によって定義されます。

side 1: (X1/Y1)-(X2/Y2)
side 2: (X2/Y2)-(X3/Y3)
side 3: (X3/Y3)-(X4/Y4)
:

光線をすべての面に対してテストする必要があります。光線をベクトル、すべての面をベクトルとみなします。光線は各面に 1 回だけ当たるか、まったく当たらないかのどちらかです。同じ面に 2 回当たることはできません。2D 空間の 2 本の線は、平行でない限り、常に 1 回だけ交差します。平行の場合は、交差することはありません。ただし、ベクトルの長さには制限があるため、2 つのベクトルは平行ではなく、短すぎて互いに出会うことがないため、交差しない可能性があります。

// Test the ray against all sides
int intersections = 0;
for (side = 0; side < numberOfSides; side++) {
    // Test if current side intersects with ray.
    // If yes, intersections++;
}
if ((intersections & 1) == 1) {
    // Inside of polygon
} else {
    // Outside of polygon
}

ここまでは順調ですが、2 つのベクトルが交差するかどうかをどのようにテストするのでしょうか? 次は、目的の C コード (テストされていません) です。

#define NO 0
#define YES 1
#define COLLINEAR 2

int areIntersecting(
    float v1x1, float v1y1, float v1x2, float v1y2,
    float v2x1, float v2y1, float v2x2, float v2y2
) {
    float d1, d2;
    float a1, a2, b1, b2, c1, c2;

    // Convert vector 1 to a line (line 1) of infinite length.
    // We want the line in linear equation standard form: A*x + B*y + C = 0
    // See: http://en.wikipedia.org/wiki/Linear_equation
    a1 = v1y2 - v1y1;
    b1 = v1x1 - v1x2;
    c1 = (v1x2 * v1y1) - (v1x1 * v1y2);

    // Every point (x,y), that solves the equation above, is on the line,
    // every point that does not solve it, is not. The equation will have a
    // positive result if it is on one side of the line and a negative one 
    // if is on the other side of it. We insert (x1,y1) and (x2,y2) of vector
    // 2 into the equation above.
    d1 = (a1 * v2x1) + (b1 * v2y1) + c1;
    d2 = (a1 * v2x2) + (b1 * v2y2) + c1;

    // If d1 and d2 both have the same sign, they are both on the same side
    // of our line 1 and in that case no intersection is possible. Careful, 
    // 0 is a special case, that's why we don't test ">=" and "<=", 
    // but "<" and ">".
    if (d1 > 0 && d2 > 0) return NO;
    if (d1 < 0 && d2 < 0) return NO;

    // The fact that vector 2 intersected the infinite line 1 above doesn't 
    // mean it also intersects the vector 1. Vector 1 is only a subset of that
    // infinite line 1, so it may have intersected that line before the vector
    // started or after it ended. To know for sure, we have to repeat the
    // the same test the other way round. We start by calculating the 
    // infinite line 2 in linear equation standard form.
    a2 = v2y2 - v2y1;
    b2 = v2x1 - v2x2;
    c2 = (v2x2 * v2y1) - (v2x1 * v2y2);

    // Calculate d1 and d2 again, this time using points of vector 1.
    d1 = (a2 * v1x1) + (b2 * v1y1) + c2;
    d2 = (a2 * v1x2) + (b2 * v1y2) + c2;

    // Again, if both have the same sign (and neither one is 0),
    // no intersection is possible.
    if (d1 > 0 && d2 > 0) return NO;
    if (d1 < 0 && d2 < 0) return NO;

    // If we get here, only two possibilities are left. Either the two
    // vectors intersect in exactly one point or they are collinear, which
    // means they intersect in any number of points from zero to infinite.
    if ((a1 * b2) - (a2 * b1) == 0.0f) return COLLINEAR;

    // If they are not collinear, they must intersect in exactly one point.
    return YES;
}

入力値は、ベクトル 1 (および) とベクトル 2 (および) の2 つの端点です。したがって、2 つのベクトル、4 つの点、8 つの座標があります。およびは明確です。交差を増やしますが、何もしません。v1x1/v1y1v1x2/v1y2v2x1/v2y1v2x2/v2y2YESNOYESNO

COLLINEAR はどうでしょうか? これは、両方のベクトルが同じ無限の直線上にあることを意味します。位置と長さに応じて、まったく交差しないか、無限の数の点で交差します。このケースをどのように処理するかは完全にはわかりませんが、どちらにしても交差としてカウントしません。まあ、このケースは実際には浮動小数点の丸め誤差のため、どちらにして== 0.0fもまれです。より良いコードでは、 をテストするのではなく、 のようなものをテストします< epsilon。ここで、イプシロンはかなり小さな数です。

より多くのポイントをテストする必要がある場合、ポリゴンの辺の線形方程式の標準形式をメモリに保存して、毎回再計算しなくても済むようにすることで、全体の処理速度を確かに少し上げることができます。これにより、ポリゴンの辺ごとに 3 つの浮動小数点値をメモリに保存する代わりに、テストごとに 2 回の浮動小数点乗算と 3 回の浮動小数点減算を節約できます。これは、メモリと計算時間の典型的なトレードオフです。

最後になりましたが、3D ハードウェアを使用してこの問題を解決できる場合は、興味深い代替手段があります。GPU にすべての作業を任せましょう。画面外のペイント サーフェスを作成します。それを完全に黒色で塗りつぶします。次に、OpenGL または Direct3D でポリゴン (ポイントがポリゴンのいずれかにあるかどうかをテストしたいだけで、どのポリゴンかは気にしない場合は、すべてのポリゴン) をペイントし、ポリゴンを別の色 (たとえば白) で塗りつぶします。ポイントがポリゴン内にあるかどうかを確認するには、描画サーフェスからそのポイントの色を取得します。これは、O(1) のメモリ フェッチで済みます。

もちろん、この方法は描画面が巨大である必要がない場合にのみ使用できます。描画面が GPU メモリに収まらない場合、この方法は CPU で実行するよりも遅くなります。描画面が巨大になる必要があり、GPU が最新のシェーダーをサポートしている場合は、上記のレイ キャスティングを GPU シェーダーとして実装することで GPU を使用できます。これは絶対に可能です。テストするポリゴンやポイントの数が多い場合は、この方法が効果的です。一部の GPU では、64 ~ 256 ポイントを並行してテストできます。ただし、CPU から GPU へのデータ転送とその逆の転送は常にコストがかかることに注意してください。そのため、ポイントまたはポリゴンのいずれかが動的で頻繁に変更される、いくつかのポイントをいくつかの単純なポリゴンに対してテストするだけでは、GPU アプローチが効果的になることはほとんどありません。

おすすめ記事