• <ins id="pjuwb"></ins>
    <blockquote id="pjuwb"><pre id="pjuwb"></pre></blockquote>
    <noscript id="pjuwb"></noscript>
          <sup id="pjuwb"><pre id="pjuwb"></pre></sup>
            <dd id="pjuwb"></dd>
            <abbr id="pjuwb"></abbr>

            zoj 1081 Points Within - 改進弧長法判斷點和簡單多邊形的關系

               轉角法判斷點和多邊形的關系大家都知道,原理比較簡單,在多邊形內掃過的轉角一定是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;
            }

            posted on 2012-07-20 10:30 yx 閱讀(1543) 評論(0)  編輯 收藏 引用 所屬分類: 計算幾何

            <2012年7月>
            24252627282930
            1234567
            891011121314
            15161718192021
            22232425262728
            2930311234

            導航

            統計

            公告

            常用鏈接

            留言簿(3)

            隨筆分類

            隨筆檔案

            me

            好友

            同學

            網友

            搜索

            最新評論

            閱讀排行榜

            評論排行榜

            亚洲国产婷婷香蕉久久久久久| 女人高潮久久久叫人喷水| 久久精品夜夜夜夜夜久久| 97热久久免费频精品99| 久久99精品久久久久久齐齐| 久久免费99精品国产自在现线| 日韩人妻无码一区二区三区久久99| 久久久无码精品亚洲日韩按摩 | 亚洲国产精品无码久久久蜜芽| 亚洲va中文字幕无码久久| 国产香蕉97碰碰久久人人| 久久久久亚洲av综合波多野结衣| 久久精品国产亚洲av麻豆小说| 精品一久久香蕉国产线看播放| 性欧美大战久久久久久久久| 久久嫩草影院免费看夜色| 精品无码久久久久国产| 国产精品美女久久福利网站| 久久国产精品波多野结衣AV| 日本欧美久久久久免费播放网| 久久久久亚洲精品男人的天堂 | 久久综合狠狠色综合伊人| 欧美亚洲国产精品久久高清| 久久九九久精品国产免费直播| 国产精品无码久久综合| 色欲综合久久中文字幕网| 一本一本久久a久久精品综合麻豆| 国内精品久久久久影院免费| 中文字幕无码免费久久| 午夜视频久久久久一区| 久久伊人影视| 日韩精品久久久久久久电影| 欧美一区二区久久精品| 久久久久九国产精品| 久久久久一本毛久久久| 久久无码国产| 久久久久亚洲AV片无码下载蜜桃| 久久九九久精品国产免费直播| 久久久久高潮综合影院| 日产精品99久久久久久| 国产精品久久波多野结衣|