多邊形分割算法在農(nóng)村土地確權(quán)中的應(yīng)用

劉蘇，植江瑜

(1.廣東瑞圖萬方科技股份有限公司廣州分公司，廣東 廣州 510655； 2.環(huán)境保護部華南環(huán)境科學(xué)研究所，廣東 廣州 510655)

多邊形分割算法在農(nóng)村土地確權(quán)中的應(yīng)用

劉蘇1*，植江瑜2

(1.廣東瑞圖萬方科技股份有限公司廣州分公司，廣東 廣州 510655； 2.環(huán)境保護部華南環(huán)境科學(xué)研究所，廣東 廣州 510655)

針對農(nóng)村土地確權(quán)工作中遇到的按比例劃分地塊的現(xiàn)實需求，提出了一種按照面積比例分割簡單多邊形的算法。該算法通過求取多邊形最小外接矩形(MABR)判定多邊形總體走勢，據(jù)此生成初始分割直線，再根據(jù)目標子多邊形面積與目標面積的差值調(diào)整分割線的位置，最終將多邊形分割為兩個邊界合理的子多邊形。利用該算法對實測農(nóng)田地塊進行一分為二的面積等分實驗，實驗結(jié)果表明：該算法適用于常見形狀的農(nóng)田地塊，分割結(jié)果合理，效率高、誤差小。

農(nóng)村土地確權(quán)；面積比例；最小凸包；最小外接矩形；分割線調(diào)整

1 引 言

農(nóng)村土地承包經(jīng)營權(quán)確權(quán)，是根據(jù)法律規(guī)定和要求，將農(nóng)民的土地承包經(jīng)營權(quán)以法律文本、證書的形式進行登記，進一步鞏固和完善農(nóng)村土地承包關(guān)系，賦予農(nóng)民更有保障的承包地財產(chǎn)權(quán)[1]。土地確權(quán)是農(nóng)村土地產(chǎn)權(quán)制度建設(shè)中的基礎(chǔ)一環(huán)，它提供的信息對公平、有效地進行土地管理十分必要，對科學(xué)制定相關(guān)政策、提高政策的針對性與有效性也大有助益。

當(dāng)前，農(nóng)村土地確權(quán)工作中常會遇到這樣一個問題：一塊耕地同時有多個承包經(jīng)營者，而這些經(jīng)營者通常是有血緣關(guān)系的，一般是兄弟。在這種情況下，若要進一步明確不同承包者對這一地塊的具體經(jīng)營權(quán)，則需根據(jù)實際情況，確定每個經(jīng)營者所占的面積比例，再依據(jù)這一比例將地塊分割成與經(jīng)營者數(shù)目相同的若干子地塊，從而確立每個經(jīng)營者與每個子地塊的一一對應(yīng)關(guān)系。

實際工作中，面積比例的確定通常根據(jù)等分原則來進行：若N個經(jīng)營者共同承包面積為S的地塊，則將地塊分割成N個面積為S/N的子地塊，當(dāng)然，有時也會遇到面積不等分的情況。若將地塊視為平面多邊形，則上述問題可抽象為基于目標面積的多邊形分割問題，即已知多邊形邊界、需要分割的子多邊形數(shù)目以及各子多邊形的面積比例，確定分割線的精確位置的過程。

對此，目前農(nóng)村土地確權(quán)工作者所采取的應(yīng)對方式主要有兩種：一是多次嘗試，根據(jù)結(jié)果再進行手動調(diào)整，直到分割結(jié)果滿意為止；二是直接忽略，不再進一步細分每個承包人的屬地。前者效率低下，誤差較大，且當(dāng)數(shù)據(jù)量較大時無從應(yīng)對；后者則不利于土地管理。而到目前為止，尚未有人就此提出更好的解決方法。

上述問題與基于目標面積修改多邊形邊界[2，3]的問題較為類似，不同的是前者不能改變原多邊形邊界，且要將多邊形分成幾份；后者則是通過修改多邊部分頂點的位置，使多邊形的面積達到目標面積[4]。本文受多邊面積修改的算法思想啟發(fā)[5]，提出了一種按面積比例分割簡單多邊形的全自動算法。

2 基本概念

為了下面敘述的方便，先給出幾個相關(guān)的定義[6]。

定義1：設(shè)pi=(xi，yi)，i=1，2，3，…，n，pn+1=p1是給定多邊形的n個頂點，若對任意i，j(i≠j)，i，j=1，2，3，…，n，線段pipi+1與pjpj+1或是相鄰且相交于一端點或不相關(guān)，則稱該多邊形為簡單多邊形。

定義2：設(shè)p1，…，pn，pn+1=p1是一個簡單多邊形。若線段pi-1pi與線段pipi+1所成的內(nèi)角(即由該多邊形所圍成有界區(qū)域內(nèi)所形成的角)是一個不超過180°的角，則稱pi頂點是凸的，否則稱是pi凹的。

由此定義可知，對任意一個多邊形，其每個頂點或是凸的，或是凹的。

定義3：設(shè)p1，…，pn，pn+1=p1是一個簡單多邊形。若沿p1→p2…pn→pn+1方向走，該簡單多邊形所圍成的有界區(qū)域總在左邊，則稱該多邊形的走向是逆時針的；反之，稱其為順時針的。

分割原則

實際工作表明，對地塊進行分割時，需要分割的份數(shù)越多，遇到的概率就越小，因此，將地塊一分為二的情況是最常見的。就算法實現(xiàn)的難度而言，當(dāng)需要分割的數(shù)目越大，分割線求取的難度也就越大，因此將多邊形一分為二的分割線的求取是最容易實現(xiàn)的。而實際上，將多邊形分成N(N>2)份，可以轉(zhuǎn)化成多次一分為二的運算，如要將一個多邊形分割成面積相等的3個子多邊形，可以先將多邊形分成原來面積的1/3和2/3兩份，再將2/3的部分平均分成兩份。因此，本文提出的算法針對的是將多邊形一分為二的情況。

根據(jù)既定的面積比例將多邊形一分為二可以有無數(shù)種方案，而本文算法的目的是確定其中最合理且操作性最強的一種作為最終的解決方案。要從無數(shù)種方案中選出一種作為最終方案，則必須要有一定的準則可對不同的方案進行評價，因此，基于土地管理的實用性，文章提出了以下兩大多邊形分割的技術(shù)原則：

原則一：分割線上的控制點盡可能少；

原則二：分割線段盡可能短。

圖1 將多邊形一分為二的分割方案

由原則一可知，將多邊形一分為二時，應(yīng)采用直線分割法，因其只需要最基本的兩個控制點即可，而折線或圓弧分割則需要較多的控制點，徒然增加實際工作量，因此上圖中的方案a不可取；方案b、c、d均為直線分割法，控制點均只有2個，符合分割原則一，三者的區(qū)別在于分割線的長短，根據(jù)分割盡可能短的原則，方案d是最優(yōu)的分割方案。分割線長短實質(zhì)上反映了子多邊形的面積緊湊度，在面積相同的情況下，周長越短，多邊形面積越緊湊[7]，從土地利用的角度而言就具有更高的實用性。

上述例子中的目標多邊形均為凸多邊形，現(xiàn)實中會經(jīng)常碰到形狀為凹多邊形的農(nóng)田，對此類多邊形進行一分為二的分割操作時，在確定采用直線分割法的前提下，還必須加上一個約束條件：直線分割操作能且只能將多邊形分割成兩個子多邊形，以防出現(xiàn)如圖2所示的將之邊形分成了3個子多邊形的情況。

圖2 分割線將凹多邊形分成3個子多邊形

3 算法描述

算法的核心思想是繪制初始分割直線將多邊形一分為二，選取其中一個子多邊形作為計算對象，計算該子多邊形與目標面積的差值(目標面積=對應(yīng)面積比例*原多邊形面積)，將面積差值除以分割線長度，得出分割線的調(diào)整位移，再根據(jù)調(diào)整位移作原分割線的平行線，得出新的分割直線，直到分割出來的子多邊形的面積與目標面積的差值小于設(shè)定的誤差限值為止，其運算過程如圖3所示。

圖3 算法簡單圖示

輸入：多邊形頂點序列；用戶指定的其中一個子多邊形的面積S目標(通過面積比例和原多邊形面積求得)；用戶根據(jù)分割原則繪制的分割直線。

輸出：分割直線與多邊形的兩個交點的坐標。

步驟①根據(jù)分割原則繪制初始分割直線AB，如圖3所示；

步驟②計算分割直線與多邊形的交點A和B的坐標；

步驟③以圖3為例，選擇分割線左側(cè)的子多邊形(即圖中陰影部分)作為調(diào)整對象，計算該子多邊形的面積S子，多邊形的面積計算公式為[8]：

步驟④計算子多邊形的面積與目標面積的差值△S=S子-S目標；

步驟⑤若△S

步驟⑥計算分割線的調(diào)整位移offset=△S/|AB|；

步驟⑦將直線AB沿其垂線方向平移offset個單位，得到調(diào)整后分割直線CD；

步驟⑧計算調(diào)整后的分割線與多邊形的新交點C和D的坐標，將C→A；D→B

重復(fù)步驟④到步驟⑧。

完整的算法流程如圖4所示；算法運行的效果如圖5所示，圖5中黑色粗線為初始分割線，顏色各異的平行細線則是每次調(diào)整后的分割線，從其軌跡可以看出，算法通過多次迭代讓分割線逐步趨近最終的精確位置。

圖4算法流程圖

圖5 算法計算效果圖

在多邊形的形狀和誤差限值給定的前提下，初始分割線距目標位置越遠，其初次迭代時調(diào)整位移offset的值就越大，但隨著迭代次數(shù)的增加，該值會越來越小，最終固定在滿足精度要求的目標位置上。而算法的精度只與設(shè)定的誤差限值相關(guān)，與多邊形的形狀無關(guān)，只要運算結(jié)果的誤差值小于設(shè)定誤差就會終止迭代，其具體數(shù)值則與多邊的形狀有關(guān)，具有一定的隨機性，但均會滿足精度要求。如圖5十邊形面積等分結(jié)果與圖6喇叭狀多邊形的面積等分結(jié)果均為誤差限值設(shè)定為0.1個單位的運算結(jié)果，前者的實際誤差值為0.093個單位，后者則為0.089個單位，就統(tǒng)計意義上而言，兩者精度等級是一樣的。理論上，設(shè)定的誤差限值越小，運算結(jié)果的精度就越高，但迭代次數(shù)會明顯增加，嚴重降低運算效率。

圖6 喇叭狀多邊形面積等分運算結(jié)果

4 算法進階

上述多邊形分割算法只是半自動算法，用戶仍需要手動輸入分割直線，若分割線輸入不符合最優(yōu)化原則，得出的分割結(jié)果也會不理想，如圖7所示。

圖7 手動輸入不合理分割線情況

用戶在繪制分割線前，實際上還包含了多邊形總體走勢的判斷，然后用戶根據(jù)判斷結(jié)果繪制其認為最合理的分割線。若多邊形的邊界十分復(fù)雜，人工判斷多邊形的總體走勢較為準確，但當(dāng)數(shù)據(jù)量較大時，依靠人工判斷多邊形總體走勢然后據(jù)此畫出分割線效率太低，不能進行批處理，若能將多邊形走勢判斷及分割線繪制都交給計算機處理，則能應(yīng)對數(shù)據(jù)量較大的情況。

實際工作中，農(nóng)田地塊大多是較為規(guī)則的四邊形，或者總體呈四邊形走勢。進行多邊形總體走勢判斷，在此可以歸納為多邊形最小面積外接矩形(Minimum Area Bounding Rectangle，簡稱MABR)的求取。除MABR外，多邊形還有一種常用的外接矩形，稱為最小綁定矩形(Minimum Bounding Rectangle，簡稱MBR)，即以多邊形頂點中最大、最小坐標確定的矩形，如圖8所示；其計算非常簡單[9]，在此不再贅述。

圖8 多邊形的MBR與MABR

本文采用文獻[10]提出的MABR計算方法，其算法計算步驟如下：

第1步計算多邊形最小凸包。計算多邊形最小凸包的算法有很多，如Graham算法、卷包算法、分治算法等[11，12]。其中由于Graham算法的時間復(fù)雜性最低[13，14]，故此處選用該算法。多邊形最小凸包計算涉及多邊形的正反向判斷[15]、多邊形頂點的凹凸性判斷等算法[6]，在此不再累贅。

第2步選取所得凸包中一條邊作為起始邊并對該凸包以選中邊的左端點為中心旋轉(zhuǎn)使平行坐標橫軸，計算并保存其MBR的坐標、該邊的編號、旋轉(zhuǎn)角度及面積。

第3步按順序選擇其他邊，并按照第2步的方法計算并保存其MBR的坐標、該邊的編號、旋轉(zhuǎn)角度及面積。

第4步比較所得的MBR的面積，其中面積最小者按其記錄的旋轉(zhuǎn)角度以該邊的左端點為圓心逆向旋轉(zhuǎn)即為所求的MABR。

根據(jù)上述算法求得多邊形的MABR后，比較其MABR的兩相鄰邊的長度，取較長對邊的中點作為分割線的兩個端點，然后執(zhí)行上述步驟2到步驟9的計算過程，就完成了基于目標面積的多邊形全自動分割，完整的算法流程圖如圖9所示。

圖9全自動算法流程圖

5 算例分析

為驗證本文算法的實用性，取廣東省某行政村XXX村的農(nóng)田數(shù)據(jù)作為實算基礎(chǔ)數(shù)據(jù)。如上文所述，項目組在對該村進行土地確權(quán)時，就碰到了十幾例需要對農(nóng)田地塊進行等分切割的情況。為了更好地驗證算法是否適用于各種形狀的農(nóng)田地塊，現(xiàn)假設(shè)該村所有490個農(nóng)田多邊形均需要進行等面積分割，分割的誤差限值設(shè)為 0.1 m2，分割前后的情況對比如圖10所示。

圖10 算例效果圖

XXX村農(nóng)田地塊多邊形信息統(tǒng)計 表1

算例運算信息統(tǒng)計 表2

圖11 各多邊形的分割誤差

表1顯示，XXX村的農(nóng)田多為四邊形，其占比接近50%，而隨著頂點的增加，多邊形的數(shù)量急劇下降，由此可見，該村農(nóng)用地塊多為形狀較為規(guī)整的四邊形。由圖10的分割結(jié)果可以看出，本算法的多邊形面分割算法結(jié)果合理；表2顯示，算法分割誤差要小于 0.1 m2，平均迭代次數(shù)不超過3次，而分割的平均誤差則不到 0.03 m2。

6 結(jié) 論

本文針對農(nóng)村土地確權(quán)中普遍存在的“一地多承包人”的問題，提出了根據(jù)面積比例分割多邊形的算法，該算法綜合了多邊形最小凸包計算、多邊形最小外接矩形計算、分割直線調(diào)整等子算法，實現(xiàn)了全自動分割多邊形的功能。利用本文算法對實測農(nóng)田地塊數(shù)據(jù)進行了面積等分實驗，結(jié)果表明，對常見形狀的農(nóng)用地塊，算法適用性強、分割結(jié)果合理、誤差較小、效率較高：完成490個多邊形的分割時間只需29ms，能勝任多數(shù)情況下農(nóng)田多邊形分割運算，具備較強的批處理功能，能切實提高實際工作中解決此類問題的精度及效率。

[1] 李秀川. 山西省農(nóng)村土地承包經(jīng)營權(quán)確權(quán)技術(shù)模式與推進措施研究[D]. 西安：西北農(nóng)林科技大學(xué)，2012.

[2] DAVID A，DAVID P. Equal-area locus-based convex polygon decomposition[C]. Structrual Information and Communication Complexity，2008：141～155.

[3] SIVAKHOLUNDU K M，PRABAHARAN N. A program to compute the area of an irregular polygon on a spheroidal surface[J]. Computers & Geosciences，1998，24(8)：823～826.

[4] Mark K J，TZVETALIN S V. Algorithm for optimal area triangulations of a convex polygon[J]. Computational Geometry；Theory and Applications，2006，25(3)：173～187.

[5] 方雷，張華鑫，姚申君. 一種基于面積修改簡單多邊形的算法[J]. 華東師范大學(xué)學(xué)報·自然科學(xué)版，2011，2：77～88.

[6] 劉潤濤. 任意多邊形頂點凸、凹性判別的簡捷算法[J]. 軟件學(xué)報，2002，13(7)：1309～1312.

[7] 毛亮，李滿春，劉永學(xué)等. 一種基于面積緊湊度的二維空間形狀指數(shù)及其應(yīng)用[J]. 地理與地理信息科學(xué)，2005，21(5)：11～14.

[8] 羅志強，鐘爾杰. 任意多邊形面積公式的推導(dǎo)及其應(yīng)用[J]. 大學(xué)數(shù)學(xué)，2005，21(1).

[9] 閆浩文. 空間方向關(guān)系理論研究[M]. 成都：成都地圖出版社，2003：23～24.

[10] 程鵬飛，閆浩文，韓振輝. 一個求解多邊形最小面積外接矩形的算法[J]. 工程圖學(xué)學(xué)報，2008，1：122～126.

[11] 吳文周，李利番，王結(jié)臣. 平面點集凸包Graham算法的改進[J]. 測繪科學(xué)，2010，35(6)：123～125.

[12] 周培德，盧開澄. 計算幾何——算法分析與設(shè)計[M]. 北京：清華大學(xué)出版社，1999： 57～58.

[13] Graham R L. An Efficient Algorithm for Determine the Convex Hull of a Finite Planar Set[J]. Information Processing Letters，1972，1(1).

[14] Chand D R，Kaqur S S. An Algorithm for Convex Polytopes[J]. JACM，1970，17(1).

[15] Yao A C. A Lower Bound to Finding Convex Hulls[J]. JACM，1981，28(4).

ApplicationofAlgorithmforDividingPolygoninConfirmationofRuralLand

Liu Su1，Zhi Jiangyu2
(1.Guangdong Ruituwanfang Technology Co.，Ltd，Guangzhou 510655，China；2.South China Institute of Environmental Science.MEP，Guangzhou 510655，China)

In the confirmation of rural land contractual operation right，it often occurs that more than one person share a single piece of land，which asks a demand for a useful method of dividing a piece of land at specific area ratio. Aiming at this requirement，a new algorithm for dividing simple polygon according to target area was proposed in the paper. First，the algorithm identified the general trend of polygon by calculating its minimum area bounding rectangle (MABR)，then generated the initial dividing line based on its MABR. After that，it adjusted the position of dividing line according to the difference between the area of the divided sub-polygon and the target area until the value of area difference descended to zero. Finally，we got two sub-polygons with target area as well as reasonable borders. For testing the practicality of this method，an experiment was taken. The results show that the algorithm is efficient，accurate and practical.

rural land contractual operation right；area ratio；minimum convex hull；minimum area bounding rectangle；dividing line adjustment

1672-8262(2017)06-146-06

P209

B

2017—02—21

劉蘇(1989—)，女，碩士，工程師，現(xiàn)主要從事測繪地理信息相關(guān)工作。

植江瑜(1987—)，男，碩士，工程師，現(xiàn)主要從事環(huán)境信息系統(tǒng)方向研究。