周小成 王鋒克 黃洪宇 馮芝清 肖祥希 李 媛

(1.福州大學(xué)空間數(shù)據(jù)挖掘與信息共享教育部重點(diǎn)實(shí)驗(yàn)室， 福州 350108;2.福建金森林業(yè)股份有限公司， 三明 353300； 3.福建省林業(yè)科學(xué)研究院， 福州 350012)

0 引言

造林挖穴是植樹(shù)造林的第一步，也是造林過(guò)程中極其重要的一環(huán)，坑穴的數(shù)量和質(zhì)量關(guān)乎植樹(shù)造林的成活率，影響最終的造林成果質(zhì)量。當(dāng)前對(duì)造林坑穴的參數(shù)提取和質(zhì)量評(píng)價(jià)主要依賴工作人員現(xiàn)場(chǎng)抽檢，一方面工作量大、效率低，另一方面客觀準(zhǔn)確性難以得到保證。無(wú)人機(jī)遙感技術(shù)應(yīng)用于植被信息提取和林業(yè)資源調(diào)查方面已經(jīng)進(jìn)行了大量研究，各種方法被廣泛應(yīng)用于植被高度、植被覆蓋度、造林規(guī)劃設(shè)計(jì)、林分參數(shù)提取等方面[1-6]。隨著圖像處理技術(shù)的發(fā)展，基于無(wú)人機(jī)影像的人工智能和模板匹配技術(shù)在目標(biāo)探測(cè)[7-8]和圖像自動(dòng)解譯[9-11]等方面得到廣泛應(yīng)用。已有文獻(xiàn)證明無(wú)人機(jī)遙感在各類目標(biāo)探測(cè)和識(shí)別方面具有較好的實(shí)用價(jià)值[12-15]。造林坑穴目標(biāo)具有特定的形狀、光譜和空間分布特征，相互之間具有較大的相似度。模板匹配技術(shù)被廣泛應(yīng)用于相似光譜特征目標(biāo)提取方面。文獻(xiàn)[16]利用基于模板的特征和基于目標(biāo)的多光譜圖像特征構(gòu)建出蟻穴自動(dòng)檢測(cè)的分類器。文獻(xiàn)[17]采用面向?qū)ο蟮哪０迤ヅ浞ê椭С窒蛄繖C(jī)法對(duì)草地鼠洞進(jìn)行自動(dòng)識(shí)別，總體精度較高?；舴蜃儞Q在提取特定形狀目標(biāo)方面有較多應(yīng)用。文獻(xiàn)[9]基于無(wú)人機(jī)遙感影像和數(shù)字表面模型(Digital surface model,DSM)數(shù)據(jù)綜合運(yùn)用圓形霍夫變換的方法實(shí)現(xiàn)柑橘樹(shù)自動(dòng)化提取。文獻(xiàn)[18]基于圓形霍夫變換提取圖像中圓形對(duì)象的方法提取高分辨率遙感影像中的圓形建筑物。

本文采用大疆精靈4 Pro型無(wú)人機(jī)，實(shí)地航拍獲取0.01 m空間分辨率的可見(jiàn)光影像數(shù)據(jù)，進(jìn)行伐區(qū)造林坑穴的數(shù)量和寬深參數(shù)提取研究，以期為造林坑穴質(zhì)量評(píng)價(jià)提供自動(dòng)化的解決方案。

1 研究區(qū)和數(shù)據(jù)獲取

研究區(qū)位于福建省將樂(lè)縣某伐區(qū)(26°51′13″N，117°34′50″E)，該林場(chǎng)位于福建省西北部，地處武夷山脈東南部、金溪河畔，多為中、低山地貌，以低山丘陵為主，平均海拔400～800 m，最高海拔為1 403 m，最低海拔為140 m。研究區(qū)整體坡度約為25°，面積約3 712 m2，該研究區(qū)的具體位置如圖1所示。

使用大疆精靈4 Pro型無(wú)人機(jī)進(jìn)行遙感影像數(shù)據(jù)采集，該型無(wú)人機(jī)搭載1英寸2 000萬(wàn)像素CMOS攝像頭，等效35 mm焦距，攜帶方便，在復(fù)雜環(huán)境下也可快速執(zhí)行影像采集任務(wù)。因測(cè)區(qū)地形起伏大且無(wú)人機(jī)不具備仿地形飛行能力，結(jié)合實(shí)際將測(cè)區(qū)劃分為兩部分分別設(shè)計(jì)航線來(lái)獲取數(shù)據(jù)，其目的是防止地面分辨率變化過(guò)大影響圖像拼接效果。航線設(shè)計(jì)為航向重疊率90%，旁向重疊率70%，飛行相對(duì)高度約為36 m，空間分辨率為0.01 m。為了保證成果的平面精度和高程精度以及將成果統(tǒng)一至CGCS2000(China geodetic coordinate system 2000)坐標(biāo)系下,現(xiàn)場(chǎng)采用華測(cè)i70型GNSS接收機(jī)獲取11個(gè)高精度地面控制點(diǎn)(Ground control point,GCP)，實(shí)時(shí)差分(Real-time kinematic，RTK)模式下的平面精度為±(8 mm+1×10-6D)，高程精度為±(15 mm+1×10-6D),其中D表示移動(dòng)站RTK與基準(zhǔn)站間距離。

使用主流無(wú)人機(jī)攝影測(cè)量處理軟件Pix4D進(jìn)行數(shù)據(jù)預(yù)處理，得到研究區(qū)正射影像(Digital orthophoto map, DOM)、數(shù)字表面模型(DSM)及無(wú)人機(jī)攝影測(cè)量點(diǎn)云產(chǎn)品。對(duì)使用GCP和不使用GCP生產(chǎn)的DSM分別進(jìn)行精度驗(yàn)證，得到的標(biāo)準(zhǔn)誤差(Root mean square error, RMSE)分別為0.056 m與0.158 m。根據(jù)林業(yè)生產(chǎn)要求，標(biāo)準(zhǔn)坑穴深度0.400 m，誤差小于15%，使用GCP生產(chǎn)的DSM滿足使用精度要求。

2 研究方法

根據(jù)無(wú)人機(jī)可見(jiàn)光影像進(jìn)行造林坑穴參數(shù)提取，整體思路為：基于無(wú)人機(jī)可見(jiàn)光影像生產(chǎn)的DOM、DSM和影像點(diǎn)云數(shù)據(jù)，采用目視解譯、模板匹配和圓形霍夫變換方法對(duì)造林坑穴的數(shù)量、寬度和深度參數(shù)進(jìn)行提取，整體流程如圖2所示。

2.1 目視解譯法提取坑穴數(shù)量

目視解譯主要基于人工判讀，在DOM上借助坑穴的影子信息容易人工辨別坑穴與非坑穴對(duì)象，在坑穴上量取南北和東西兩方向長(zhǎng)度并計(jì)算均值獲取寬度信息，在無(wú)人機(jī)影像點(diǎn)云上手工量取深度信息。根據(jù)人工解譯結(jié)果，最終在研究區(qū)得到坑穴總計(jì)1 021個(gè)。

2.2 基于DOM模板匹配法提取坑穴數(shù)量

模板匹配是在目標(biāo)檢測(cè)領(lǐng)域最早被提出以及相對(duì)來(lái)說(shuō)較為簡(jiǎn)單的檢測(cè)方法，它的原理是依據(jù)相關(guān)策略,根據(jù)已知模塊在搜索圖像中尋找逼近模塊的匹配過(guò)程[19]。造林坑穴有標(biāo)準(zhǔn)尺寸，相互之間存在很大的相似度，因此選擇模板匹配技術(shù)作為坑穴數(shù)量檢測(cè)方法之一。

2.2.1模板創(chuàng)建與測(cè)試

基于DOM創(chuàng)建初始模板時(shí)應(yīng)選擇具有代表性以及不同光照條件下的坑穴。標(biāo)準(zhǔn)坑穴穴面寬為0.5 m,因此設(shè)置固定窗口為50個(gè)像素的正方形來(lái)選取樣本。分別計(jì)算所有樣本在紅、綠、藍(lán)波段的平均值生成初始模板，并統(tǒng)計(jì)各樣本和初始模板之間的平均相關(guān)系數(shù)，樣本生成模板示意圖如圖3所示。實(shí)驗(yàn)得紅、綠、藍(lán)3個(gè)波段平均相關(guān)系數(shù)分別為0.75、0.73、0.70，相關(guān)系數(shù)越大說(shuō)明該波段越適合坑穴的目標(biāo)提取。因此本文選擇紅光波段數(shù)據(jù)作為源數(shù)據(jù)進(jìn)行坑穴提取。

初始模板生成后應(yīng)用于隨機(jī)選擇的某區(qū)域，正確識(shí)別率未達(dá)到80%以上則手動(dòng)糾正并將正確識(shí)別的坑穴加入樣本重新生成模板。重復(fù)上述過(guò)程，直至滿足要求，此時(shí)模板創(chuàng)建完成，紅光波段下樣本與模板的平均相關(guān)系數(shù)提升至0.77。

2.2.2模板應(yīng)用

將最終生成的模板應(yīng)用于正射影像的紅光波段圖像中，采用滑動(dòng)窗口的方式在每個(gè)位置上計(jì)算備選區(qū)與模板之間相關(guān)系數(shù)。設(shè)置合適的相關(guān)系數(shù)閾值，本實(shí)驗(yàn)設(shè)置閾值為0.7，即相關(guān)系數(shù)不小于0.7可判斷為坑穴。相關(guān)系數(shù)計(jì)算式為

(1)

式中 (x,y)——待匹配圖像像元坐標(biāo)

(x′,y′)——模板圖像像元坐標(biāo)

T——模板圖像

I——待匹配圖像(圖4)

w——模板寬h——模板高

R——相關(guān)系數(shù)

2.3 基于DSM圓形霍夫變換法提取坑穴數(shù)量、寬深參數(shù)

由正射影像圖觀察可知，基于DOM提取坑穴主要利用坑穴影子的光譜信息，實(shí)際提取得到的是坑穴的影子，因此使用DOM模板匹配法只能獲取坑穴數(shù)量信息而無(wú)法獲取真實(shí)的寬度和深度信息。此外，因山區(qū)地形和天氣變化影響，DOM容易出現(xiàn)色彩不均勻的情況，這對(duì)于模板匹配法坑穴提取較為不利。

人工造林坑穴坑面近似于圓形，運(yùn)用邊緣檢測(cè)技術(shù)對(duì)坑穴對(duì)象進(jìn)行邊界提取，之后基于邊界圖像進(jìn)行圓形霍夫變換找到其中潛在的圓形可有效提取坑穴對(duì)象。在本研究中，基于圓形霍夫變換對(duì)圓形坑穴對(duì)象進(jìn)行提取主要步驟包括：邊緣檢測(cè)、圓形霍夫變換找圓和坑穴參數(shù)計(jì)算3個(gè)步驟。

2.3.1Canny邊緣檢測(cè)

圖像邊緣分布在圖像中像素的灰度、紋理等特征劇烈變化的地方，圖像中的目標(biāo)與圖像背景的灰度在這些地方必然會(huì)存在階躍，表現(xiàn)為函數(shù)圖像呈現(xiàn)出劇烈的變化。常用的邊緣檢測(cè)算子有Sobel算子、Laplacian算子和Canny算子等，其中，Sobel算子適用于漸變圖像，對(duì)邊緣的定位精度不高；Laplacian算子是各向同性的，但是對(duì)噪聲比較敏感。對(duì)于坑穴DSM數(shù)據(jù)來(lái)說(shuō)，坑穴邊界高程發(fā)生跳躍式分布，坑穴輪廓清晰,地面情況復(fù)雜，本文選擇Canny算子進(jìn)行邊緣檢測(cè)。

Canny算子作為一種標(biāo)準(zhǔn)的邊緣檢測(cè)算子，由CANNY在1986年提出，被廣泛應(yīng)用于各類圖像邊緣提取[20]。Canny算子濾波步驟包括：高斯濾波、梯度幅度和方向計(jì)算、非極大值抑制以及雙閾值算法檢測(cè)邊緣與連接邊緣。經(jīng)過(guò)Canny算子的處理，坑穴邊界像素大于0，其他像素值全部歸零，從而得到準(zhǔn)確的坑穴邊界信息，如圖5b所示。

2.3.2圓形霍夫變換找圓

霍夫變換現(xiàn)在被廣泛用于自動(dòng)化檢測(cè)圖像中的直線、圓或橢圓，其最初的目的是從黑白圖像中提取直線。圓形霍夫變換則主要用于自動(dòng)化檢測(cè)圖像中的圓形，其檢測(cè)圓的基本原理為

(x-a)2+(y-b)2=r2

(2)

其中

x=a+rcosθ

(3)

y=b+rsinθ

(4)

式中a、b——圓心橫、縱坐標(biāo)r——圓半徑

θ——直線法線與x軸的夾角

根據(jù)極坐標(biāo)，圓上任意一點(diǎn)坐標(biāo)(x,y)(式(2))可以用式(3)、(4)來(lái)表示。通過(guò)在0～2π之間改變夾角θ，就可以確定任意點(diǎn)是否位于半徑為r、圓心為(a，b)的圓上。

基于Canny算子得到的邊界圖像應(yīng)用圓形霍夫變換找圓，將物理空間長(zhǎng)度參數(shù)(m)轉(zhuǎn)換為圖像空間中像素?cái)?shù)，通過(guò)考慮坑穴的尺寸與坑穴之間距離來(lái)確定霍夫變換參數(shù)。圓形霍夫變換使用的參數(shù)為(rmin、rmax、mdist、p1、p2、dp)，其中rmin為坑穴最小半徑、rmax為最大半徑、p1為梯度，p2為累加器閾值，dp是累加器分辨率和圖像分辨率的反比[21]，mdist為霍夫變換檢測(cè)到圓心之間的最小距離。本文采用Python-OpenCV編程的方式實(shí)現(xiàn)該方法，得到坑穴提取結(jié)果，該過(guò)程如圖5所示。

2.3.3坑穴寬深參數(shù)計(jì)算

人工檢測(cè)坑穴寬度實(shí)際上是量測(cè)坑穴面寬的均值，該值約等于圓形霍夫變換得到的圓形直徑d。

由單個(gè)坑穴的三維影像點(diǎn)云(圖6)可知，坑穴深度為坑面點(diǎn)與坑底點(diǎn)高程之差。分別統(tǒng)計(jì)每個(gè)坑穴對(duì)象內(nèi)部像素點(diǎn)高程的最大值和最小值，兩者的差值即為坑穴的深度參數(shù)Dp。Z(x,y)表示該坐標(biāo)像素高程，pv表示坑穴對(duì)象內(nèi)包含的所有像素點(diǎn)。提取公式為

(5)

3 結(jié)果與分析

以目視解譯數(shù)量和實(shí)際測(cè)量的坑穴參數(shù)作為真值驗(yàn)證提取參數(shù)值，引入正確率(Correctness)、漏檢率(Omission)和誤檢率(Commission)3個(gè)檢驗(yàn)指標(biāo)，對(duì)3種提取方法進(jìn)行精度檢驗(yàn)。

3.1 坑穴數(shù)量提取精度分析

選取坑穴局部提取結(jié)果(圖7)，將目視解譯結(jié)果作為正確的識(shí)別結(jié)果，對(duì)比分析其他2種提取結(jié)果?？梢钥闯?，基于DOM模板匹配法(圖7c)相比于基于DSM圓形霍夫變換法(圖7b)有較多的坑穴誤提和漏提現(xiàn)象；基于DSM圓形霍夫變換法則與目視解譯結(jié)果最為接近。

利用檢驗(yàn)指標(biāo)，分別計(jì)算2種方法提取結(jié)果的識(shí)別精度，如表1所示。結(jié)果表明，基于DSM圓形霍夫變換坑穴提取法精度最高，其正確率為95.15%，誤檢率為4.85%，漏檢率為6.56%。該法相較于基于DOM模板匹配法提取的結(jié)果，正確率提升明顯，漏檢率得到進(jìn)一步的控制，更加符合質(zhì)檢需要。

表1 2種方法坑穴數(shù)量提取結(jié)果Tab.1 Results of two methods in number extraction

3.2 基于DSM圓形霍夫變換法坑穴寬、深參數(shù)提取精度分析

實(shí)際人工采集20個(gè)坑穴的寬、深參數(shù)對(duì)提取結(jié)果進(jìn)行精度驗(yàn)證，由實(shí)測(cè)寬度與提取寬度、實(shí)測(cè)深度與提取深度線性擬合結(jié)果(圖8、9)分析，二者的R2分別達(dá)到0.93和0.92，RMSE分別為1.02 cm和1.67 cm，提取參數(shù)和實(shí)測(cè)值之間呈現(xiàn)明顯的相關(guān)性，提取結(jié)果較為可靠。

3.3 提取方法評(píng)價(jià)

對(duì)比3種方法在提取坑穴數(shù)量和寬度、深度參數(shù)的耗時(shí)、精度、工作量以評(píng)估處理速度和適用性，由表2可以看出，相比于人工目視解譯方法的耗時(shí)，

表2 3種方法評(píng)價(jià)Tab.2 Validation results of three methods

其余2種方法在處理速度上顯然更具優(yōu)勢(shì)；在考量精度、工作量和自動(dòng)化程度基礎(chǔ)上，基于DSM圓形霍夫變換法既可提取數(shù)量又可提取寬度和深度信息，因而具有最高的適用性。

對(duì)于山區(qū)地形來(lái)說(shuō)，因?yàn)榈匦纹鸱?、太?yáng)光線變化以及周邊高大樹(shù)木遮擋的影響，最終生產(chǎn)得到的DOM數(shù)據(jù)可能會(huì)存在色彩不均勻的現(xiàn)象?；谀０迤ヅ浞ㄌ崛】友▽?duì)象要求DOM色彩分布均勻，當(dāng)坑穴目標(biāo)位于陰影下或者曝光過(guò)度區(qū)域則容易導(dǎo)致漏檢和誤檢，從而使得整體漏檢率和誤檢率偏高。此外，模板質(zhì)量對(duì)模板匹配法提取坑穴的精度影響較大，樣本的選取要求較高，這在一定程度上降低了處理效率和精度。

基于DSM圓形霍夫變換法一方面利用簡(jiǎn)單幾何知識(shí)和常識(shí)快速確定算法所使用的參數(shù)，且不需要人工選擇樣本，人工干預(yù)較少，處理效率和精度相比于基于DOM模板匹配法提升明顯；另一方面它在獲取坑穴數(shù)量信息的同時(shí)也可得到高精度的寬度和深度信息，適用于現(xiàn)階段造林坑穴數(shù)量和寬度、深度參數(shù)的自動(dòng)化提取。

4 結(jié)束語(yǔ)

無(wú)人機(jī)遙感應(yīng)用于造林坑穴數(shù)量和坑穴寬度、深度參數(shù)提取切實(shí)可行。結(jié)合無(wú)人機(jī)影像，相較于基于DOM模板匹配法92.60%的正確率，基于DSM數(shù)據(jù)的圓形霍夫變換法坑穴數(shù)量提取正確率達(dá)到95.15%，誤檢率和漏檢率均處于較低水平，坑穴寬度和深度提取結(jié)果與實(shí)測(cè)驗(yàn)證值R2分別達(dá)到0.93和0.92，RMSE分別為1.02 cm和1.67 cm，滿足質(zhì)檢要求。該方法相比于現(xiàn)場(chǎng)人員抽查檢測(cè)更加客觀、科學(xué)和高效，是坑穴信息提取的最佳解決方案。本研究在數(shù)據(jù)采集過(guò)程中，地面控制點(diǎn)坐標(biāo)采集降低了工作效率，建議在實(shí)際生產(chǎn)中使用自帶RTK型無(wú)人機(jī)來(lái)獲取影像數(shù)據(jù)，生產(chǎn)得到高精度DSM。