基于無人機(jī)可見光影像的華山松人工林計(jì)測(cè)參數(shù)提取

袁梓健,舒清態(tài),趙洪瑩,王柯人,譚德宏

(西南林業(yè)大學(xué) 林學(xué)院,昆明 650224)

0 引言

在陸地生態(tài)系統(tǒng)的功能運(yùn)轉(zhuǎn)中,森林是極為關(guān)鍵的組成要素,會(huì)對(duì)碳水資源的順利循環(huán)以及氣溫氣候的調(diào)節(jié)等起到重要的影響作用。單木是森林空間結(jié)構(gòu)中的重要組成,有效地獲取詳細(xì)的單木結(jié)構(gòu)參數(shù)能夠幫助我們更好地掌握森林資源信息。傳統(tǒng)林業(yè)信息統(tǒng)計(jì)與調(diào)查工作需要經(jīng)由實(shí)地測(cè)樹以及抽樣等方式[1-3],明確相當(dāng)范圍內(nèi)的樹木生長(zhǎng)相關(guān)的分布以及質(zhì)量等信息[4-9]。此種調(diào)查模式下,測(cè)量人員自身的工作情況將會(huì)影響到調(diào)查結(jié)果,呈現(xiàn)出很強(qiáng)的主觀性特點(diǎn),而且大范圍的統(tǒng)計(jì)需要耗費(fèi)很大的人力與物力資源,整體的效率水平偏低[10-11]。伴隨著科技的持續(xù)發(fā)展與進(jìn)步,人們嘗試將遙感技術(shù)引入到森林調(diào)查工作中,并且憑借著自身數(shù)據(jù)處理信息快及能夠大范圍進(jìn)行測(cè)量等優(yōu)勢(shì)而被廣泛運(yùn)用到林業(yè)發(fā)展等諸多生產(chǎn)領(lǐng)域[12-13]。20世紀(jì)七八十年代,國(guó)內(nèi)就在著手利用遙感技術(shù)來實(shí)現(xiàn)森林資源相關(guān)的調(diào)查統(tǒng)計(jì)工作,尤其是使用航空像片完成的森林蓄積量分析成果顯著。而且,所得像片的正射影像圖能夠更加清晰地將實(shí)地森林情況展現(xiàn)出來,為更好地實(shí)現(xiàn)現(xiàn)代化的林業(yè)經(jīng)營(yíng)及生產(chǎn)發(fā)展等工作提供信息參考[14]?；诖?宮鵬[15]使用數(shù)字?jǐn)z影測(cè)量相關(guān)技術(shù)與理論完成橡樹草原相關(guān)生長(zhǎng)情況的探測(cè)研究,針對(duì)其樹高以及郁閉度情況進(jìn)行記錄；王二麗等[16]基于涼山自然保護(hù)區(qū)0.2m分辨率的航空像片,選取子像元分類的處理方式來獲得其云冷杉、落葉松以及紅松相關(guān)的分布情況,構(gòu)建得到提取冠幅及其實(shí)測(cè)冠幅相關(guān)的回歸模型,各區(qū)目視解釋精度都超出80%,成效顯著。從近些年的發(fā)展來看，Popescu等[17]嘗試采取基于局部最大值法來對(duì)樹頂點(diǎn)進(jìn)行提取的方式以獲得單木數(shù)量,針對(duì)闊葉樹及松樹相關(guān)的樣地林分密度進(jìn)行預(yù)估,所得估測(cè)值與其實(shí)測(cè)值呈現(xiàn)的相關(guān)系數(shù)都要超過0.66；Koukoulas等[18]選擇應(yīng)用冠層高度模型獲得等高線的方式來分割單木樹冠,在落葉闊葉林的分析中,使用該方法獲得的樹頂精度顯著超出LiDAR數(shù)據(jù),精度能夠達(dá)到91%；Kwak等[19]經(jīng)由數(shù)據(jù)估測(cè)的方式來研究日本落葉松以及韓國(guó)紅松相關(guān)的單木樹高情況,并且使用形態(tài)學(xué)影像分析及局部最大值來提取數(shù)字冠層模型內(nèi)的樹頂點(diǎn),對(duì)比實(shí)測(cè)值明確其評(píng)估精度維持范圍為67.4%～86.7%。

華山松(PinusarmandiiFranch.)是云南省廣泛栽植的用材樹種。近年來，華山松人工林和天然林的經(jīng)營(yíng)規(guī)模和種植面積不斷增加,在經(jīng)濟(jì)、社會(huì)發(fā)展和生態(tài)安全中占據(jù)著重要地位[20]。本文以石林縣國(guó)有林場(chǎng)華山松人工林為研究對(duì)象,以樣地實(shí)際調(diào)查數(shù)據(jù)結(jié)合無人機(jī)影像數(shù)據(jù),通過兩種不同的單木分割方法,得到單木計(jì)測(cè)參數(shù),以期為華山松人工林的動(dòng)態(tài)監(jiān)測(cè)提供參考。

1 研究區(qū)概況及數(shù)據(jù)獲取

1.1 研究區(qū)概況

石林縣處在滇東高原腹心,地面起伏和緩,切割較輕。山脈統(tǒng)屬云嶺余脈,沿縣境東西兩側(cè)向東北展布延伸,整個(gè)地勢(shì)自東向西呈階梯狀傾斜。平均海拔1 737m,屬低緯高原山地季風(fēng)氣候區(qū),雨量充沛。研究區(qū)位于云南省昆明市石林縣國(guó)有林場(chǎng)(北緯 24°30′～25°03′,東經(jīng)103°10′～103°41′)。此國(guó)有林場(chǎng)的總面積為4.46萬hm2,其中林業(yè)用地面積為3.45萬hm2,約占總面積的77.35%。

1.2 無人機(jī)數(shù)據(jù)獲取

2019年11月19日中午進(jìn)行無人機(jī)數(shù)據(jù)的采集工作：此日晴空萬里無云,并且風(fēng)速相對(duì)較小,能夠很好地滿足無人機(jī)工作的要求；選擇中午,能夠減小樹冠的陰影面積對(duì)后期數(shù)據(jù)處理的誤差。

實(shí)驗(yàn)選用的無人機(jī)為大疆-悟系列(Inspire),相機(jī)使用專業(yè)級(jí)三軸增穩(wěn)云臺(tái)相機(jī)(型號(hào):禪思 Zenmuse X4S)。飛行前對(duì)無人機(jī)傳感器進(jìn)行校準(zhǔn),在飛行區(qū)域內(nèi)布設(shè)3個(gè)像控點(diǎn),并進(jìn)行參數(shù)設(shè)定(航向重疊率設(shè)置為75%,旁向重疊率設(shè)置為65%)。飛行速度7.1M/S,飛行高度為100m。

1.3 樣地?cái)?shù)據(jù)獲取

選取石林縣國(guó)有林場(chǎng)內(nèi)301林班12小班內(nèi)共5塊邊長(zhǎng)為28.28m方形標(biāo)準(zhǔn)樣地(圖1)。使用南方測(cè)繪(T66 RTK)進(jìn)行樣地四角拐點(diǎn)的單木坐標(biāo)獲取,使用全站儀采集樣地內(nèi)華山松單木坐標(biāo),使用激光測(cè)距儀測(cè)量單木樹高,使用胸徑尺測(cè)量樣地單木胸徑,使用皮尺測(cè)量單木東西冠幅以及南北冠幅。樣地調(diào)查數(shù)據(jù)(胸徑、樹高、冠幅)概況如表1所示。

圖1 研究區(qū)的樣地及樣木Fig.1 Sample land and sample wood in the study area

表1 樣地?cái)?shù)據(jù)概況Tab.1 Overview of sample site data

2 研究方法與數(shù)據(jù)處理

2.1 單木分割方法

2.1.1基于冠層髙度模型的單木提取方法

冠層高度模型(CHM)能夠表達(dá)林分中每棵樹木的枝葉到地面的垂直距離。由于華山松的冠層大體分布呈現(xiàn)傘狀,邊緣位置相較地面偏近且中心位置相較地面偏遠(yuǎn),故而繪制得到的CHM灰度圖像內(nèi),單木表現(xiàn)是近圓形狀,邊緣值偏小,中心值高,圖像中冠狀最大值就是單木的最高點(diǎn)位置,基于局部最大值探測(cè)來實(shí)現(xiàn)具體的單木提取工作。將可變圓形窗口設(shè)置為探測(cè)窗口(圖2),將所得局部最大值標(biāo)作樹頂來完成統(tǒng)計(jì)。

圖2 圓形探測(cè)窗口Fig.2 Circular probe window

2.1.2基于點(diǎn)云分割的單木提取方法

利用原始航片和地面控制點(diǎn)生成調(diào)查樣地的三維點(diǎn)云數(shù)據(jù)模型。點(diǎn)云數(shù)據(jù)中的高程值不僅包含了樣地內(nèi)樹木、灌木叢等高程,而且還包含了地面點(diǎn)高程。因?yàn)檠芯繀^(qū)所在地地勢(shì)高低起伏,由于地面高程值疊加的影響作用存在,高大樹木所得高程值存在小于低矮樹木的情況。如此，不能夠順利經(jīng)由點(diǎn)云數(shù)據(jù)提取的方式來展現(xiàn)單木的高度信息。因此,在對(duì)華山松展開提取研究期間,需要對(duì)點(diǎn)云數(shù)據(jù)進(jìn)行歸一化操作,將地面點(diǎn)高程值由所得高程值內(nèi)減掉,所得的數(shù)據(jù)就是地物的實(shí)際高度。使得點(diǎn)云數(shù)據(jù)與沒有歸一化之前形成對(duì)比,呈現(xiàn)華山松人工林的真實(shí)形態(tài)。

該算法從種子點(diǎn)A(即全局最大值)開始(圖3)。主要根據(jù)間距臨界值和最小間距規(guī)則,通過對(duì)更低的點(diǎn)進(jìn)行估計(jì),將種子點(diǎn)A發(fā)展為一個(gè)樹聚類。例如:A點(diǎn)是最高點(diǎn),因此將點(diǎn)A作為一號(hào)樹木的樹頂(目標(biāo))。此時(shí)把低于A的點(diǎn)相繼進(jìn)行分類,首先,點(diǎn)B被分類成二號(hào)樹木,因?yàn)殚g距dAB大于一個(gè)設(shè)定的臨界值；然后設(shè)置點(diǎn)C,點(diǎn)C相對(duì)于B點(diǎn)又可定義為二號(hào)待分割點(diǎn)的二號(hào)樹木。因?yàn)閐AC小于dBC,通過與點(diǎn)B和點(diǎn)C比較,點(diǎn)D相對(duì)于點(diǎn)C又被分類成二號(hào)樹木；通過與點(diǎn)C和點(diǎn)D的比較,點(diǎn)E相對(duì)于點(diǎn)D又被分類成二號(hào)樹木。臨界值應(yīng)當(dāng)與冠層半徑相等,當(dāng)設(shè)置臨界值太大或太小就會(huì)出現(xiàn)分割不足和過度分割的情況。

圖3 點(diǎn)云分割算法原理Fig.3 Principle of point cloud segmentation algorithm

2.2 華山松人工林蓄積量模型建立及精度驗(yàn)證

2.2.1蓄積量估測(cè)模型

選取5塊標(biāo)準(zhǔn)地的華山松純林為研究對(duì)象,對(duì)兩種方法所提取的冠幅和胸徑與外業(yè)調(diào)查的數(shù)據(jù)進(jìn)行一次多項(xiàng)式函數(shù)和二次多項(xiàng)式函數(shù)擬合,并建立華山松二元材積方程(表2)。

表2 二元材積模型估測(cè)方程Tab.2 Estimation equation of binary volume model

2.2.2精度驗(yàn)證

本研究建立的模型由決定系數(shù)(R2)、均方根誤差(RMSE)進(jìn)行模型精度驗(yàn)證和評(píng)價(jià)。

1)決定系數(shù)(R2)。

采用決定系數(shù)(R2)對(duì)估測(cè)模型進(jìn)行精度檢驗(yàn),R2的變化范圍在0～1之間,值越大,表明自變量與因變量之間的相關(guān)性較好,模型精度高。表達(dá)式為:

2)均方根誤差(RMSE)。

均方根誤差(RMSE)用于量化模型精度,RMSE值越小,表明模型預(yù)測(cè)效果越好、精度越高。表達(dá)式為:

2.3 數(shù)據(jù)處理

本研究5塊樣地共采集無人機(jī)遙感影像875余張,利用無人機(jī)數(shù)據(jù)處理軟件Pix4D mapper 生成樣地的數(shù)字表面模型(DSM)和數(shù)字正射影像(DOM),用DSM與已有的DEM高程模型相減獲取調(diào)查地的冠層高度模型(CHM),并利用調(diào)查樣地時(shí)設(shè)置的控制點(diǎn)與原始航拍數(shù)據(jù)生成調(diào)查地點(diǎn)云數(shù)據(jù)。

3 結(jié)果與分析

3.1 株數(shù)的提取結(jié)果與驗(yàn)證

實(shí)現(xiàn)單木分割的關(guān)鍵所在就是順利由CHM冠層高度模型內(nèi)獲得單木樹冠最大值,進(jìn)而完成單木標(biāo)記。考慮到樹冠最大值所在位置就是樹頂點(diǎn),故而選取局部最大值來實(shí)現(xiàn)樹頂點(diǎn)相關(guān)的識(shí)別工作。由于所調(diào)查標(biāo)準(zhǔn)地內(nèi)的華山松基本是同一年種下的人工林且為單層林,冠幅長(zhǎng)勢(shì)相對(duì)較為平均,所以采用固定窗口對(duì)最大值進(jìn)行識(shí)別。由于華山松在生長(zhǎng)時(shí),冠幅越大的單木對(duì)應(yīng)的樹高和胸徑也是在逐漸變大,所以窗口的大小通過樹高與胸徑和冠幅之間的關(guān)系確定,所以根據(jù)5塊標(biāo)準(zhǔn)地內(nèi)實(shí)測(cè)數(shù)據(jù)建立樹高與胸徑和冠幅之間的聯(lián)系。

通常局部最大值的識(shí)別窗口為正方形。但是本文的研究對(duì)象為華山松的樹冠區(qū)域,樹冠的投影近似為圓形,若使用正方形窗口則會(huì)對(duì)我們的識(shí)別有一定影響,所以,選擇使用圓形窗口對(duì)冠層高度模型進(jìn)行局部最大值的識(shí)別誤差就會(huì)減小。搜索窗口通過樹高與冠幅、胸徑進(jìn)行擬合得出結(jié)果。比如：若圓心處像元值為15,通過樹高-冠幅、胸徑擬合結(jié)果求得冠幅為9.4m即為圓形搜索窗口直徑；樹冠假定為圓形,則樹冠半徑即為4.7m；所使用的CHM空間分辨率為0.4m,對(duì)應(yīng)識(shí)別窗口像元值為11。

在5塊調(diào)查樣地內(nèi),地面實(shí)測(cè)共計(jì)獲得443棵單木,通過兩種單木提取方法提取單木株數(shù)與實(shí)測(cè)地面單木位置進(jìn)行對(duì)比(表3),未分離株數(shù)是指一定范圍內(nèi)沒有提取出任何單木,過度分離株數(shù)是指一株單木被分離成了多株單木。

表3 株數(shù)提取對(duì)比分析Tab.3 Comparison of plant number extraction and analysis

兩種提取單木屬性的方法所獲得的單木數(shù)量都小于樣地實(shí)測(cè)的單木數(shù)量,5塊樣地中基于生成的點(diǎn)云數(shù)據(jù)分離了406棵樹,基于CHM分離了413棵樹。由于地勢(shì)、冠幅遮擋等影響,提取到的單木與實(shí)際單木數(shù)據(jù)對(duì)應(yīng)會(huì)存在一定誤差。利用實(shí)測(cè)的單木坐標(biāo)位置與提取到的單木數(shù)據(jù)進(jìn)行對(duì)比。在對(duì)比的過程中,主要發(fā)現(xiàn)單木分割會(huì)出現(xiàn)未分離,多分離的錯(cuò)誤分割。未分割指由于冠幅過于密集,在大范圍內(nèi)并沒有分割出單木；多分離是指在大范圍內(nèi)一棵單木被分割成了許多單木。

基于CHM的單木提取，沒有被分割的單木與實(shí)測(cè)數(shù)據(jù)相比，均為樹高較小的單木。在CHM上的像元值較低,就使得一部分低矮單木無法通過像元大小被識(shí)別出來,造成了漏分的情況。

基于CHM的單木提取共分離出 413棵樹,通過與5塊樣地實(shí)測(cè)單木坐標(biāo)數(shù)據(jù)對(duì)比,刪除誤差較大的單木,還有360棵樹能夠正確被識(shí)別,分離正確率為87.17%；基于點(diǎn)云數(shù)據(jù)的單木提取共分離出 406棵樹,通過與5塊樣地單木坐標(biāo)實(shí)測(cè)數(shù)據(jù)對(duì)比,刪除誤差較大的單木,還有328棵樹能夠正確被識(shí)別,分離正確率為80.79%。

3.2 樹高的提取結(jié)果與驗(yàn)證

分別將兩種單木分割方法所得樹高與實(shí)測(cè)值進(jìn)行分析比較,兩類方法所得樹高與實(shí)測(cè)樹高均表現(xiàn)出很強(qiáng)的線性相關(guān)。其中：基于CHM所得樹高h(yuǎn)1與實(shí)測(cè)值h兩者的R2對(duì)應(yīng)是0.71,相應(yīng)的回歸方程式為h=0.7183h1+3.5542(圖4)；基于點(diǎn)云分割所得樹高h(yuǎn)2與實(shí)測(cè)值h兩者的R2對(duì)應(yīng)是0.69,相應(yīng)的回歸方程式為h=0.7854h2+2.9254(圖5)。

圖4 基于CHM提取樹高相關(guān)性分析Fig.4 High correlation analysis of tree based on CHM

圖5 基于點(diǎn)云提取樹高相關(guān)性分析Fig.5 High correlation analysis of tree based on Point Cloud

3.3 蓄積量的估測(cè)與驗(yàn)證

由兩種方法精度對(duì)比可以得出,基于CHM比點(diǎn)云的單木計(jì)測(cè)參數(shù)的提取精度高。通過對(duì)CHM單木分割獲取的單木參數(shù)進(jìn)行D=a+bCW+cH的模型估測(cè)(D為單木胸徑；CW為單木冠幅；H為單木樹高；a,b,c為常數(shù))。所得胸徑—冠幅、樹高模型為D=1.76+4.66CW+0.27H。

用統(tǒng)計(jì)分析方法建立的模型,其計(jì)算結(jié)果影響到林業(yè)數(shù)表應(yīng)用中的可靠性和森林經(jīng)營(yíng)決策的正確性。建立單木材積與胸徑、樹高模型,使用最小二乘法進(jìn)行回歸模型的建立。采用最小二乘法建立回歸模型的誤差應(yīng)該滿足“獨(dú)立、正態(tài)、等方差”條件。進(jìn)行多元非線性回歸計(jì)算參數(shù)值；査云南省華山松二元材積表,得出實(shí)測(cè)單株材積；經(jīng)由模型相關(guān)的擬合研究可得詳細(xì)的單木二元材積模型參數(shù)值。石林縣國(guó)有林場(chǎng)華山松人工林二元材積模型擬合結(jié)果如表4所示。

表4 石林縣國(guó)有林場(chǎng)華山松人工林二元材積模型擬合結(jié)果Tab.4 Fitting results of binary volume model of P.armandi in State-owned forest farm of Shilin county

通過上述材積模型的建立,由表4可知,模型3的擬合程度最好。雖然二元材積模型擬合不能只看精度,但模型3在國(guó)內(nèi)外二元材積模型表達(dá)式中應(yīng)用較廣,并且也符合我國(guó)的二元材積表編制公式。模型3決定系數(shù)(R2=0.937843),均方根誤差(RMSE=0.0338)具有比較好的擬合度和使用性。選用模型3作為最優(yōu)材積模型,其表達(dá)示為V=0.000073D2.001569H0.788883；將已求得胸徑-冠幅、樹高模型帶入已知參數(shù)中,得到最優(yōu)單木二元材積模型為V=0.000073(1.76+4.66CW+0.27H)2.001569×H0.788883。根據(jù)最優(yōu)單木二元材積模型公式,分別計(jì)算5塊樣地的對(duì)應(yīng)蓄積量。其計(jì)算結(jié)果如表5所示。

4 結(jié)果與討論

由提取單木得到的冠幅以及樹高等相關(guān)參數(shù)完成評(píng)價(jià)分析,明確使用CHM能夠在樣地內(nèi)分離獲得413株單木,與實(shí)測(cè)得到單木株數(shù)的正確匹配量是360棵,準(zhǔn)確率高達(dá)87.17%；使用點(diǎn)云分割能夠在樣地內(nèi)分離獲得406株單木,與實(shí)測(cè)得到的正確匹配量是328株,準(zhǔn)確率高達(dá)80.79%。表現(xiàn)出顯著的欠分離問題。提取得到的樹高與其地面實(shí)測(cè)所得樹高的R2相比較,使用CHM方法,R2為0.71；而使用點(diǎn)云算法,R2為0.69。提取所得數(shù)據(jù)相較實(shí)際所得數(shù)據(jù)呈現(xiàn)出很大的偏差現(xiàn)象。針對(duì)各項(xiàng)參數(shù)評(píng)價(jià)情況進(jìn)行綜合判斷可知,基于CHM實(shí)現(xiàn)的單木提取效率更加突出,可以更加充分地展現(xiàn)森林中單木具體的分布情況及其參數(shù)信息。

無人機(jī)數(shù)據(jù)具有便攜性、采集數(shù)據(jù)快等特點(diǎn)。本文由無人機(jī)原始航片生成樣地點(diǎn)云數(shù)據(jù),以及通過樣地模型數(shù)據(jù)生成冠層高度模型。兩種方法具有一定的便攜性和可行性,但是對(duì)本研究可言,也存在一定誤差。1)在對(duì)CHM使用最大像元值窗口搜索時(shí),會(huì)存在少部分的無效像元值,這對(duì)后續(xù)的冠幅分割會(huì)產(chǎn)生一定的影響；2)使用的點(diǎn)云數(shù)據(jù)為無人機(jī)原始航片生成,地形因子在點(diǎn)云歸一化會(huì)后存在一定的影響；3)無人機(jī)數(shù)據(jù)源也是重要的誤差來源,無人機(jī)在采集數(shù)據(jù)時(shí),像控點(diǎn)的選取對(duì)生成正射影像及生成點(diǎn)云數(shù)據(jù)會(huì)有一定的誤差；4)樣地采集數(shù)據(jù)過程使用南方測(cè)繪(T66 RTK)測(cè)量單木坐標(biāo),由于林場(chǎng)內(nèi)信號(hào)較弱,RTK的固定解坐標(biāo)會(huì)存在一定誤差,對(duì)后續(xù)單木坐標(biāo)的定位存在影響；5)在測(cè)量單木計(jì)測(cè)參數(shù)時(shí),由于所在樣地為人工林,林分結(jié)構(gòu)完整,對(duì)于樹高及胸徑的測(cè)量能夠減少后續(xù)數(shù)據(jù)處理的誤差,但是測(cè)量冠幅時(shí)，存在冠幅被遮擋、生成彎曲等情況，導(dǎo)致冠幅的測(cè)量結(jié)果存在一定的誤差。在今后的研究中,為了提高無人機(jī)估測(cè)森林計(jì)測(cè)參數(shù),應(yīng)該減少數(shù)據(jù)獲取過程中帶來的誤差,以達(dá)到更高精度的單木計(jì)測(cè)參數(shù)的提取。