杜蒙蒙，李瀚遠，金 鑫，ROSHNIANFARD Ali，周 浩

（1.河南科技大學(xué)農(nóng)業(yè)裝備工程學(xué)院，河南 洛陽 471003；2.莫哈格達阿比德利大學(xué)農(nóng)業(yè)與自然資源系，伊朗 阿爾達比勒 566199）

0 引言

近年來，極端暴雨天氣與自然災(zāi)害頻發(fā)，導(dǎo)致農(nóng)田受損，影響農(nóng)作物生長，降低了產(chǎn)量和品質(zhì)[1-2]。在雨水沖刷或其他自然災(zāi)害等外界因素作用下，農(nóng)田中往往出現(xiàn)尺寸不等的凸起與洼地等微地形特征，不利于進行農(nóng)田耕作[3]。凸起特征在形態(tài)上表現(xiàn)為中間高、四周低的封閉地勢，與經(jīng)典地貌學(xué)中的山脊特征相似，指包含單個或多個高程極大值點的局部區(qū)域；而洼地特征則為四周高、中間低的封閉地勢，與經(jīng)典地貌學(xué)中的山谷相似，指包含一個或多個高程極小值點的局部區(qū)域[4]。

傳統(tǒng)土地平整忽視了農(nóng)田微地形特征的結(jié)構(gòu)尺寸和位置分布，易造成平整地機具過載[5]。準確識別農(nóng)田微地形特征，并在此基礎(chǔ)上對其進行精準平整，可避免平整地作業(yè)過載率偏高的問題。數(shù)字地形模型（Digital Terrain Model，DTM）是分析地形特征的重要工具[6]。仲騰等[7]利用反轉(zhuǎn)數(shù)字地形模型提取山頂特征，該方法需人為確定閾值，導(dǎo)致提取精度易受人員主觀判斷影響。易煒[8]采用圖像處理算法劃定最大標高值區(qū)域，根據(jù)人為設(shè)定的閾值提取山頂特征，該方法受人員主觀判斷的影響較大。DENG Y 等[9]使用多個地形因子綜合評價用于提取山頂特征的閾值，在模型中有效剔除了不符合實際狀況的山頂特征點。

本研究利用低空無人機遙感圖像，建立高精度農(nóng)田數(shù)字地形模型（Farmland Digital Terrain Model，F(xiàn)DTM），實現(xiàn)農(nóng)田微地形特征的定位與輪廓提取，可為精準平整地作業(yè)提供方法支持，具有迫切的現(xiàn)實需求和重要的理論意義。

1 數(shù)據(jù)獲取

DTM 空間分辨率大小直接決定了表達地形曲面的精細化程度，影響地形分析的精度。因此，建立高精度FDTM，可有效提高農(nóng)田微地形特征的提取精度。研 究 區(qū) 域 位 于 陜 西 省 寶 雞 市（107°9′14.893″E，34°39′182″N），其尺寸約為50 m×100 m。采用極飛C2000 型智能測繪無人機，搭載機械全局快門2 000 萬像素的航拍相機，以距地面相對高度80 m 進行航拍測繪作業(yè)，共獲取42 張航片，航片橫向重疊率85%、縱向重疊率80%。基于SfM 技術(shù)，獲取空間分辨率3 cm的FDTM，如圖1 所示。試驗田高程范圍為853.527～854.821 m，高程極差為1.294 m。

圖1 試驗地塊的FDTMFig.1 FDTM of experimental farmland

2 研究方法

2.1 確定局部窗口尺度

統(tǒng)計不同分析窗口尺度下FDTM 的高程方差，分析其隨窗口尺度的變化規(guī)律，確定分析窗口的尺度范圍，如圖2 所示。

圖2 高程方差變化率Fig.2 Change rate of elevational variance

由圖2 可知，高程方差在[0，30] 像素和[60，100]像素區(qū)間內(nèi)的變化率較大，呈發(fā)散狀態(tài)；而在[30，60] 像素區(qū)間內(nèi)的變化率較小，收斂于9.77×10?5m2。另外，相關(guān)文獻表明，在高程方差收斂區(qū)間內(nèi)，地形特征的識別精度隨局部窗口尺度的增大而增加[10]。因此，本研究選擇51 像素×51 像素窗口確定農(nóng)田高程極大值點，共計315 個；結(jié)合31 像素×31像素、41 像素×41 像素和51 像素×51 像素3 個窗口，從農(nóng)田高程極大值點集中定位凸起特征的頂點。

2.2 分析地形因子綜合隸屬度

農(nóng)田現(xiàn)實環(huán)境中存在大量細碎土塊和秸稈殘茬，嚴重影響農(nóng)田凸起特征的識別效率及精度。地形因子是數(shù)字地形分析的基礎(chǔ)與核心內(nèi)容，其主要分為宏觀和微觀地形因子[11]。本研究根據(jù)農(nóng)田微地形結(jié)構(gòu)，確定高程、坡度和地形起伏度作為農(nóng)田高程極大值點的評價指標。高程表征某點沿鉛垂線到大地水準面的距離；坡度計算如式（1）所示，指坡面豎直高度和水平距離的比值；地形起伏度計算如式（2）所示，指局部范圍內(nèi)地面高程最大值與最小值之差[11]。

式中S?坡度，（°）

fx?高程在水平方向的變化率

fy?高程在垂直方向的變化

RF?地形起伏度，m

Hmax、Hmin?局部窗口下的最大高程、最小高程，m

周侗等[12]研究表明，在單一窗口下受干擾因素影響，無法精準識別試驗地塊的微地形特征。因此，本研究首先利用式（3）對51 像素×51 像素窗口下的高程、地形起伏度和坡度數(shù)據(jù)進行歸一化處理，使用熵權(quán)法對其進行權(quán)重賦值，如式（4）～式（6）所示，分別確定高程、地形起伏度和坡度的權(quán)重為0.384、0.487和0.129。于是，根據(jù)式（7）計算高程極值點在51 像素×51 像素窗口下的地形因子綜合隸屬度。進一步，分別計算31 像素×31 像素與41 像素×41 像素窗口下高程極值點的地形因子綜合隸屬度，并使用線性權(quán)重法確定31 像素×31 像素、41 像素×41 像素和51 像素×51 像素窗口的權(quán)重為0.252、0.333 和0.415，如式（8）所示。最終，利用式（9）計算多窗口地形因子綜合隸屬度，如表1 所示。

表1 多窗口地形因子綜合隸屬度Tab.1 Comprehensive membership degree of multi window terrain factors

式中i?地形因子序號，i=1，2，3，分別為高程、地形起伏度和坡度

j?高程極值點序號，j=1，2，…，315

k?窗口尺度大小，k=31，41，51

fijk?地形因子的歸一化值

min（xijk）、max（xijk）?在窗口尺度k下第i個地形因子的最小值和最大值

xijk?在窗口尺度k下第j個高程極值點第i個

地形因子的值

Pijk?在窗口尺度k下第j個高程極值點第i個地形因子的概率

Hik?在窗口尺度k下各地形因子的信息熵

ωik?在窗口尺度k下各地形因子的權(quán)重

φjk?在窗口尺度k下3 個地形因子對各農(nóng)田高程極值點的地形因子綜合隸屬度

γk?局部窗口的權(quán)重

ψjk?多窗口綜合地形因子綜合隸屬度

3 結(jié)果與分析

3.1 提取農(nóng)田凸起特征

隨著多窗口地形因子綜合隸屬度的增大，相應(yīng)的高程極大值點成為凸起特征頂點的概率隨之增大。因此，多窗口地形因子綜合隸屬度的閾值取值過小，作物殘茬或細小土塊易被誤判為農(nóng)田凸起特征；反之，識別農(nóng)田凸起特征的準確率降低。本研究利用斯特吉斯公式對多窗口地形因子綜合隸屬度進行抽樣處理，以確定最佳的地形因子綜合隸屬度閾值[13]。通過式（10）和式（11）分別計算分級組數(shù)n=9 及分級組距D=0.075，確定多窗口地形因子綜合隸屬度與對應(yīng)的高程極值點數(shù)量，繪制高程極值點數(shù)量的直方圖，如圖3 所示。

圖3 多窗口地形因子綜合隸屬度頻率Fig.3 Comprehensive membership degree frequency of multiwindow terrain factors

式中n?分級組數(shù)

C?農(nóng)田高程極值點數(shù)量（315 個）

D?分級組距，m

φmax、φmin?多窗口地形因子綜合隸屬度的最大值和最小值

由圖3 可知，在0.176～0.627 區(qū)間范圍內(nèi)，隨著地形因子綜合隸屬度的增大，高程極值點的數(shù)量未見明顯變化，表征農(nóng)田中大量不規(guī)則分布的細碎土塊與作物殘茬；而在綜合隸屬度0.627 處高程極值點數(shù)量急劇減小，即出現(xiàn)尺寸顯著增大的農(nóng)田凸起特征。因此，確定0.627 作為地形因子綜合隸屬度的閾值，并提取表1 中地形因子綜合隸屬度大于該閾值的農(nóng)田高程極值點，即農(nóng)田凸起特征的頂點，共計16 個，實現(xiàn)農(nóng)田凸起特征的定位，其空間分布如圖4a 所示。

圖4 提取農(nóng)田凸起特征Fig.4 Extracting hump features of farmland

等高線無法定量描述坡度和地形起伏度等地形信息，但能夠直觀地表達地形結(jié)構(gòu)特征。研究表明，等高距約等于地形圖空間分辨率時可較為完整地表達地形信息[14]。因此，本研究利用3 cm 等高距的等高線圖，并提取上述農(nóng)田凸起特征頂點的等高線閉合區(qū)域，識別農(nóng)田凸起特征的外形輪廓，如圖4b 所示。

3.2 提取農(nóng)田洼地特征

RDTM 是一種非客觀存在的地形，以DTM 的平均高程作為水平面基準，經(jīng)反向翻轉(zhuǎn)構(gòu)成的數(shù)字地形模型。RDTM 與DTM 存在反轉(zhuǎn)對稱關(guān)系，兩者坡度相差180°，原始DTM 的高程值極小點表現(xiàn)為RDTM 的高程值極大點。因此，原始FDTM 的洼地特征在形態(tài)上表現(xiàn)為RFDTM 的凸起特征。建立研究區(qū)域的RFDTM，重復(fù)前文提取農(nóng)田凸起特征的方法，獲取9 個農(nóng)田洼地特征，其空間分布如圖5 所示。

圖5 提取農(nóng)田洼地特征Fig.5 Extracting concave features of farmland

4 結(jié)束語

以陜西省寶雞市某試驗區(qū)域為研究對象，采用極飛C2000 型智能測繪無人機進行航拍測繪作業(yè)，通過分析農(nóng)田高程方差隨局部窗口尺度的變化規(guī)律，確定適應(yīng)于試驗田微地形分析的局部窗口尺度范圍為31～51 像素。利用相關(guān)系數(shù)確定高程、地形起伏度和坡度3 個地形因子，并計算農(nóng)田高程極大值點在31 像素×31 像素、41 像素×41 像素和51 像素×51 像素窗口下的多窗口地形因子綜合隸屬度，根據(jù)閾值法定位16個凸起特征與9 個洼地特征，最終結(jié)合農(nóng)田等高線圖提取上述微地形特征的外形輪廓，實現(xiàn)農(nóng)田微地形特征的定位與輪廓提取，為精準平整地作業(yè)提供支持。