轉角法判斷點和多邊形的關系大家都知道，原理比較簡單，在多邊形內掃過的轉角一定是360度，在邊界上和外面則不一定。
實現起來也比較麻煩，浮點誤差比較大，而且還要考慮些特殊情況。
在網上找到一種叫做改進弧長法的算法，原理和轉角法類似，但是做了很多重要的改進。比如，計算轉角改成了計算叉積，根據叉積決定
旋轉方向，還要根據計算下一個點的象限決定偏轉多少，每次偏轉的都是90度的倍數。
該算法可以方便判斷出點在多邊形內，還是邊界上，還是在多邊形外面。
摘自別人對該算法的描述如下：
首先從該書中摘抄一段弧長法的介紹：“弧長法要求多邊形是有向多邊形，一般規定沿多邊形的正向，邊的左側為多邊形的內側域。
以被測點為圓心作單位圓，將全部有向邊向單位圓作徑向投影，並計算其中單位圓上弧長的代數和。若代數和為0，則點在多邊形外部；
若代數和為2π則點在多邊形內部；若代數和為π，則點在多邊形上。”
按書上的這個介紹，其實弧長法就是轉角法。但它的改進方法比較厲害：將坐標原點平移到被測點P，這個新坐標系將平面劃分為4個
象限，對每個多邊形頂點P ，只考慮其所在的象限，然後按鄰接順序訪問多邊形的各個頂點P，分析P和P[i+1]，有下列三種情況：
（1）P[i+1]在P的下一象限。此時弧長和加π/2；
（2）P[i+1]在P的上一象限。此時弧長和減π/2；
（3）P[i+1]在Pi的相對象限。首先計算f=y[i+1]*x-x[i+1]*y（叉積），若f=0，則點在多邊形上；若f<0，弧長和減π；若f>0，弧長和加π。
最後對算出的代數和和上述的情況一樣判斷即可。
實現的時候還有兩點要注意，第一個是若P的某個坐標為0時，一律當正號處理；第二點是若被測點和多邊形的頂點重合時要特殊處理。
以上就是書上講解的內容，其實還存在一個問題。那就是當多邊形的某條邊在坐標軸上而且兩個頂點分別在原點的兩側時會出錯。
如邊(3,0)-(-3,0)，按以上的處理，象限分別是第一和第二，這樣會使代數和加π/2，有可能導致最後結果是被測點在多邊形外。而實際上
被測點是在多邊形上（該邊穿過該點）。
對於這點，我的處理辦法是：每次算P和P[i+1]時，就計算叉積和點積，判斷該點是否在該邊上，是則判斷結束，否則繼續上述過程。
這樣犧牲了時間，但保證了正確性。
具體實現的時候，由於只需知道當前點和上一點的象限位置，所以附加空間只需O(1)。實現的時候可以把上述的“π/2”改成1，“π”改成2，
這樣便可以完全使用整數進行計算。不必考慮頂點的順序，逆時針和順時針都可以處理，只是最後的代數和符號不同而已。整個算法編寫
起來非常容易。
代碼如下：
#include <stdio.h>
#include <math.h>
#define MAX (100 + 10)
struct Point
{
double x,y;
};
Point pts[MAX];
const int OUT = 0;
const int IN = 1;
const int EDGE = 2;
const double fPre = 1e-8;
int DblCmp(double fD)
{
if (fabs(fD) < fPre)
{
return 0;
}
else
{
return fD > 0 ? 1 : -1;
}
}
int GetQuadrant(Point p)
{
return DblCmp(p.x) >= 0 ? (DblCmp(p.y) >= 0 ? 0 : 3) :
(DblCmp(p.y) >= 0 ? 1 : 2);
}
double Det(double fX1, double fY1, double fX2, double fY2)
{
return fX1 * fY2 - fX2 * fY1;
}
int PtInPolygon(Point* pts, int nN, Point p)
{
int i, j, k;
for (j = 0; j < nN; ++j)
{
pts[j].x -= p.x;
pts[j].y -= p.y;
}
int nA1, nA2;
int nSum = 0;
nA1 = GetQuadrant(pts[0]);
for (i = 0; i < nN; ++i)
{
k = (i + 1) % nN;
if (DblCmp(pts[k].x) == 0 && DblCmp(pts[k].y) == 0)
{
break;//與頂點重合
}
int nC = DblCmp(Det(pts[i].x, pts[i].y,
pts[k].x, pts[k].y));
if (!nC && DblCmp(pts[i].x * pts[k].x) <= 0
&& DblCmp(pts[i].y * pts[k].y) <= 0)
{
break;//邊上
}
nA2 = GetQuadrant(pts[k]);
if ((nA1 + 1) % 4 == nA2)
{
nSum += 1;
}
else if ((nA1 + 2) % 4 == nA2)
{
if (nC > 0)
{
nSum += 2;
}
else
{
nSum -= 2;
}
}
else if ((nA1 + 3) % 4 == nA2)
{
nSum -= 1;
}
nA1 = nA2;
}
for (j = 0; j < nN; ++j)
{
pts[j].x += p.x;
pts[j].y += p.y;
}
if (i < nN)
{
return EDGE;
}
else if (nSum)//逆時針nSum == 4, 順時針nSum == -4
{
return IN;
}
else
{
return OUT;
}
}
int main()
{
int nN, nM;
int nCase = 1;
while (scanf("%d%d", &nN, &nM), nN)
{
if (nCase > 1)
{
printf("\n");
}
for (int i = 0; i < nN; ++i)
{
scanf("%lf%lf", &pts[i].x, &pts[i].y);
}
printf("Problem %d:\n", nCase++);
for (int i = 0; i < nM; ++i)
{
Point p;
scanf("%lf%lf", &p.x, &p.y);
if (PtInPolygon(pts, nN, p))
{
printf("Within\n");
}
else
{
printf("Outside\n");
}
}
}
return 0;
}
