土壤侵蝕形態演化數字攝影觀測系統設計與實驗

姜艷敏 郭明航 趙 軍 溫仲明 林 奇 史海靜

(1.中國科學院水利部水土保持研究所， 陜西楊凌 712100； 2.中國科學院大學， 北京 100049；3.西北農林科技大學水土保持研究所， 陜西楊凌 712100； 4.西北農林科技大學草業與草原學院， 陜西楊凌 712100；5.西安敦瑞測量技術有限公司， 西安 710065)

0 引言

中國是世界上水土流失最嚴重的國家之一，全國土壤侵蝕面積為356萬km2，約占國土面積的40%[1]，尤其是我國的黃土高原地區，水土流失情況十分嚴重，對生態環境及社會經濟的可持續發展帶來嚴重影響[2-3]。土壤侵蝕是侵蝕力與下墊面相互作用的過程，下墊面既是侵蝕作用的對象又是侵蝕作用的結果[4]。在侵蝕外力的作用下，土壤或其他地面組成物質被剝蝕、搬運以及沉積，從而產生侵蝕，而侵蝕作用的結果使得原來的下墊面形態及地表條件發生改變，從而引起新的侵蝕，如此反復，使得侵蝕作用不斷發展。因此，土壤侵蝕過程是一個逐漸發展演化的過程。

長期以來，各國學者圍繞土壤侵蝕量調查、侵蝕過程與機理，研發了很多的觀測方法和技術手段[5-9]。從傳統的人工手動測量法[10]、插釬法[11]、示蹤法[12]、徑流泥沙采樣法[13]，到目前基于測繪、測量技術以及信息技術為基礎的新的觀測方法，如三維激光掃描法[14]、高精度全球定位系統(Global positioning system，GPS)法[15]、攝影測量法[16]等，為土壤侵蝕觀測研究提供了多樣化的選擇。盡管在坡面侵蝕發生發展過程方法和手段的研究方面取得了一定進展，但無論哪種方法，尚不能解決連續降雨條件下土壤侵蝕形態的觀測問題，目前大多數研究集中在一次性降雨條件下侵蝕形態發育過程的主觀描述。如激光掃描觀測法雖有較高的空間定位精度，但不能在降雨過程中觀測，且在溝蝕階段會在溝道出現漏測[17-18]。高精度GPS法雖然測量精度高，但在地形復雜情況下易受衛星信號的影響，穩定性差，且不能在連續降雨條件下觀測。因此，目前土壤侵蝕坡面觀測方法仍存在測量時空關系不一致、精度不高、實時性差等問題。研究更為先進的技術和手段對土壤坡面侵蝕演化過程進行觀測，對于土壤侵蝕機理、土壤侵蝕演化過程等研究具有重要意義。

近年來，隨著攝影觀測技術的發展，數字化近景攝影觀測技術逐漸應用于土壤侵蝕測量[19-21]。文獻[22-26]利用數字攝影近景觀測技術，在室內模擬降雨條件下對土壤侵蝕的演變過程進行觀測，但無法實現在連續降雨過程中對土壤坡面的觀測。GUO等[27]研發了一種手持式數字化近景攝影觀測系統，實現了連續降雨條件下的觀測。該系統采用高幀率的電荷耦合裝置(Charge coupled device，CCD)工業相機，快速獲取高度重疊的影像，在獲取地面信息的整個過程中利用手持式的垂直掃描采集方式，巧妙避開了降雨過程中的大多數雨滴。但并未對影像中的雨滴噪聲做直接處理，另外，手持式的掃描使得該系統影像獲取的瞬時性低，限制了系統觀測的時間分辨率，從而降低了系統的實用性。因此，在連續降雨條件下去除雨滴干擾，瞬時獲取土壤侵蝕坡面形態變化的信息，是土壤侵蝕過程研究亟需解決的新問題。

基于此，本文耦合數字近景攝影觀測技術和無線組網技術，設計一套能夠在連續降雨條件下對土壤侵蝕下墊面形態演變過程進行觀測的系統。通過組網技術，并行拍攝、解算降雨過程中拍攝的下墊面的數字影像，提取具有高時間分辨率和高空間分辨率的下墊面土壤侵蝕形態演化信息，從而為土壤侵蝕過程研究提供新的途徑和技術手段。

1 工作原理與系統組成

1.1 工作原理

系統基于無線網絡技術對若干相機進行組網。相機基于無線網絡命令，并行采集數據，數據采集時將各組傳感器單次采集的數字影像按時間排序，逐像素按其灰度值做二分類處理，進而實現雨滴噪聲的剔除。系統基于攝影測量技術完成下墊面對象的高精度、高密度三維點云重建。一場降雨可以獲得多個時間點的三維場景數據，以達到動態的觀測效果。

1.2 系統組成

圖1 數字近景攝影觀測系統的邏輯結構設計Fig.1 Structural design of digital close-range photogrammetry system

數字近景攝影觀測系統由影像采集、影像傳輸和影像解算3個功能子系統組成，且每個子系統由不同的軟硬件單元組成(圖1)。系統的各功能子系統都在一臺運行環境為Windows 7的高容量PICO(Participant intervention comparison outcome)計算機控制下運行。 并針對各功能子系統開發了對應的軟件系統，以z-map命名，其包括相機工作狀況診斷、影像采集、影像解算及數字高程模型(Digital elevation model，DEM)生成等功能界面。

1.2.1影像采集系統

影像采集系統負責土壤侵蝕下墊面數字影像的采集、觸發信號的接收、雨滴去除等工作。該系統的硬件部分主要由12臺索尼CMOS(Complementary metal oxide semiconductor)相機和工控機組成的相機組、直流電源、標靶、防水轉接件等部件構成。軟件部分由總控制PICO計算機z-map軟件的影像采集單元組成。影像采集信號觸發后，12臺相機組并行采集下墊面數字影像，通過工控機對各相機單次采集的多幅數字影像，按同一位置像素單元灰度大小排序，逐像素依據其灰度運用K-means算法進行聚類處理，去除雨滴在數字影像上所形成的噪聲，獲得去除雨滴后的下墊面的數字影像。

1.2.2影像解算系統

影像解算系統主要負責影像數字點云的匹配、三維重建、DEM生成以及土壤侵蝕量計算等。由超高容量的匹配機來實現數據的存儲、匹配、三維重建等解算工作，與影像采集系統的軟件部分一樣，只需通過設置數據解算后存放的路徑即可完成影像數據的解算。數字影像解算系統的軟件有3個模塊：并行計算管理模塊、點云匹配和編輯模塊、DEM生成和土壤侵蝕量計算模塊。在各模塊算法開發過程中采用python語言配合NumPy(Numerical python)計算庫來做原型的研發，之后再采用C++語言重新實現。這樣的流程減少了調試過程中的時間消耗，又能保證最終執行代碼的效率。

1.2.3影像傳輸系統

影像傳輸系統在采集系統和解算系統之間起連接作用，主要負責控制命令的發出、信號接收、影像數據的傳輸。無線路由器、網絡協議(Transmission control protocol/Internet protocol，TCP/IP)、千兆網硬件接口是影像傳輸系統的主要硬件單元。各子系統之間通過無線路由器組成一個局域網絡，控制和計算單元通過無線網絡發布并發采集命令，影像采集系統采集影像后并發作業，再通過影像傳輸系統把采集到的影像傳輸給影像解算系統。

2 關鍵技術與實現方法

2.1 多影像采集器的組網及空間坐標系的建立

2.1.1采集器的選取與組建

借助于無線路由器通過TCP/IP網絡將若干組數字影像采集器進行組網，實現影像采集器的并發作業，獲取同一時間節點下的土壤下墊面信息，每組數字影像采集器包括一個數碼照相機和一個工業控制級別的計算機。本系統采用的相機是索尼CMOS相機，相機分辨率為3 264像素×2 448像素，配有12 mm鏡頭，實用光圈為F 1.2，相機幀率為15 f/s，為保證更大的拍攝視角，共選取了12臺相機。12臺相機共同組建在距離地面高18 m的鋼筋板架上，且相機之間呈均勻排列，與地面土槽呈垂直方向布設。與每個相機相匹配工作的硬件單元是電源和工控機，電源負責給相機和工控機供電，工控機控制相機的影像采集、雨滴去除等工作。

2.1.2空間坐標系的建立

空間坐標系的建立是將所有影像采集器的坐標進行統一。影像采集器所記錄的是空間物體信息的二維圖像，為了獲取實際空間物體表面某點的三維幾何位置，必須建立物體的三維空間坐標和對應的二維圖像坐標之間的對應關系，影像采集器的幾何成像模型決定了目標物體表面點的坐標與其在二維圖像中的像素坐標的對應關系，而解算相機參數是建立幾何成像模型的前提。參數的求解通過相機標定來完成[28]，相機的參數包括內部參數和外部參數。內部參數包括相機的焦距f，圖像主點的x、y坐標(cx,cy)，畸變參數K1、K2和K3；外部參數是拍攝圖像時相機的投影中心點坐標(x0，y0，z0)和3個旋轉角(ψ，ω，κ)[29]。本文以針孔模型為相機標定的理論基礎，借助棋盤格和控制點作為相機標定參照物，以Microsoft Visual Studio為開發平臺，采用開源計算機視覺庫(Open source computer vision library，OpenCV)編譯相機標定程序，求解相機的參數。

圖2 相機標定Fig.2 Calibration of camera

標定實驗在黃土高原土壤侵蝕與旱地農業國家重點實驗室人工模擬降雨實驗大廳進行，實驗場地為10 m(長)×4 m(寬)的液壓式可調坡度鋼制土槽。相機內部參數依據張正友棋盤格標定算法獲取[30]。首先，將小網格長和寬均為0.05 m的棋盤格面板平整地放置于標定實驗場土槽平面(圖2a)。采用所有相機拍攝多幅圖像，改變棋盤格面板的方位及傾斜度，再次拍攝多幅圖像。將圖像導入算法中求解內部參數。相機的外姿態通過在標定實驗場內布設控制點來獲取。選取102個小方塊標志作為標定控制點，將所有控制點均勻地布設于土槽表面(圖2b)。采用所有相機拍攝多幅圖像，調整控制點的位置和距離，再次拍攝多幅圖像。在相機內部參數已知的基礎上，通過編制解算程序代碼求解相機外部參數，從而得到每個控制點的x、y、z坐標。

經過標定算法及編譯程序的迭代運算，得到各參數的近似值，若近似值在一定的容許范圍內收斂，則結束迭代運算，得到最終的參數值。實驗中各相機內部參數的標定結果如表1所示。

表1 各相機的標定參數Tab.1 Calibration parameters of each camera

2.2 數字影像中雨滴的去除

坡面數字影像的采集是在降雨條件下進行的，降雨過程中，雨滴在空間的場分布近似于隨機場，相機拍攝得到的影像混合了雨滴和下墊面對象的兩類信息。雨滴的去除是獲取坡面物點精確信息的前提和必要工作。在短暫的時間段內，比如幾秒的時間段內，下墊面對象可以認為是一個穩定的空間對象，主要變化的是隨機性很高的雨場數據。根據這一思路，對各組傳感器單次采集的數字影像按時間排序，逐像素按其灰度值做二分類處理，并通過K-means[31-32]算法去除雨滴在數字影像上所形成的噪聲。K-means算法是一種基于形心劃分的聚類算法[33]，它以數據到形心的距離作為目標函數，并以誤差平方和準則函數作為聚類質量的度量函數，不斷進行迭代計算求極值優化聚類結果[34]。具體算法過程如下：

(1)

(2)

min(Pix(i,j)1,Pix(i,j)2,…,Pix(i,j)n))/2

(3)

式中i、j——像素行、列數

n——迭代次數

(2)以像素的灰度距離聚類，并且構造選擇集

(4)

(5)

(6)

Pix|Pix=min(dis0,dis1,dis2)

(7)

(3)以三分類選擇集內元素聚合，平均得到新的類中心

(8)

(9)

(10)

(11)

對每個工控機采集的60幅原始影像像素灰度進行方差分析。結果顯示，60幅原始影像的方差在28.96～29.95之間，而去除雨滴后的影像方差為24.14，由此可知K-means算法能夠較好地去除雨滴噪聲。

2.3 點云匹配

數字點云匹配是將多幅影像進行匹配得到同名點的過程。在匹配過程中首先提取影像上的SIFT(Scale invariant feature transform)特征[35]，SIFT特征的提取主要是依據多幅影像上具有明顯特征的點，包括邊緣點，具有明顯特征的點、孤立的點[36]，以SIFT算子作為特征提取的工具，獲取多幅影像的SIFT特征。通過SIFT特征匹配影像之間的同名點，同名點匹配是基于特征點及其描述子的相似性來進行的。獲得同名點后，以光束法區域網平差原理[37]進行平差處理，以單幅影像的光線束為平差單元，以中心投影的共線方程作為平差的基礎方程，將相機攝影點、相機成像點及其相應地面物點坐標作為一個整體，組建共線誤差方程，通過迭代計算求出該匹配點的坐標。再通過空間前方交繪計算同名點的三維坐標，所有匹配成功的同名點都計算完畢后，可獲得土壤侵蝕坡面的三維點云。

2.4 點云的修補和DEM生成

當土壤侵蝕出現溝道后，徑流就會沿溝道匯集，數字影像難以拍攝到水流下方物點的影像，即難以得到水流下方物點對應的像點信息，這給水流下方數字點云的匹配帶來了困難。通過對溝道存在水流情形下數字影像的分析以及匹配解算發現，由于溝底凸凹不平，總能匹配出一些稀疏的點云。基于這種客觀存在，根據地形變化的連續性，利用地學普遍采用的反距離權重法[38](Inverse distance weighted，IDW)，依據水流周邊和水流區域稀疏的點云坐標，便可內插得到水淹區域空缺物點的三維坐標，從而擬合出其數字點云(圖3)。修補后的點云，利用同樣的距離平方反比法，在系統軟件模塊下插值便可生成下墊面DEM。反距離權重法計算式為

(12)

式中a——參與計算的像素點數

m——像素點總數

Va——控制點的屬性值

da——控制點與當前計算點間的距離

V——計算所得當前點的屬性值

圖3 溝道底部數字點云修補Fig.3 Digital point cloud reparation at bottom of flow channel

3 實例分析

3.1 實驗布設

應用實驗在黃土高原土壤侵蝕與旱地農業國家重點實驗室人工模擬降雨實驗大廳進行。實驗小區為可調坡度鋼制土槽，小區規格：長×寬×深為10 m×1.0 m×0.5 m，土槽坡度為20°(圖4a)，土槽下端設置集流裝置，用于收集徑流泥沙。供試土壤為黃綿土，裝填土容重為1.3 g/cm3。小區布設后，正式降雨前在實驗土槽表面均勻撒水，再靜置24 h，使土槽內部土壤水分的再分配達到應力均勻、土壤結構穩定。降雨強度設置為120 mm/h，降雨歷時為150 min。至坡面產流后開始收集徑流含沙量全樣，每隔5 min采樣1次(圖4b)。降雨結束后，采用便攜式徑流泥沙儀測量泥沙含量。在收集徑流泥沙的同時，對土壤侵蝕坡面進行全覆蓋的數字影像采集，降雨開始前采集第1次坡面的數字影像，降雨開始后每隔5 min采集1次，直至降雨結束。所采集圖像的重疊度至少為4°，采集幀率不少于15 f/s。

圖4 實驗布設Fig.4 Parallel experiments and layout

3.2 精度檢測

精度是指對同一對象多次測量值的穩定程度，選用標準差來衡量。為檢測數字近景攝影觀測系統的測量精度，將具有標準尺寸的標尺均勻地布設于土槽表面，并在土槽表面任意位置布設2個標靶以進行長度約束。采用數字近景攝影觀測法，在相同的光照和紋理條件下重復拍照60次，單獨對每次的圖像集合進行匹配計算，測量每把標尺的尺寸，對60次測量的結果進行統計分析并選取中誤差作為衡量測量精度的指標。

通過SPSS 18軟件對標尺的60次測量結果統計分析可知，標尺測量的平均長度為309.270 3 mm，測量的標準差為1.711 3 mm，單次最小測量誤差為0.006 2 mm，說明該觀測的精度達到毫米級。經過K-S(Kolmogorov-Smirnov)檢驗，得出該標尺尺寸測量結果的Z值為0.392，P值為0.999，大于0.05，由此可知該標尺尺寸的測量結果均服從正態分布，測量數據的分布如圖5所示。

圖5 標尺測量數據分布圖Fig.5 Histogram of measured results of ruler

3.3 準確度檢測

3.3.1凹槽尺寸觀測法

準確度是測量值與實際值之間的偏差，以相對誤差來衡量。在人工模擬降雨條件下，通過不同雨強(30、60、90、120 mm/h)、土槽坡度(0°、5°、10°、15°)共20種組合條件下，采用數字近景攝影觀測系統獲取土槽坡面3個已知凹槽長、寬、深的尺寸，并計算其與實際值之間的相對誤差，從而評估數字近景攝影觀測系統對土壤侵蝕坡面幾何尺寸的觀測準確度。

(1)同雨強不同坡度觀測

對數字近景攝影觀測系統在雨強60 mm/h，0°、5°、10°、15°共4個不同坡度條件下的觀測準確度進行檢測(表2)。測量值與實際值兩者之間的最大相對誤差為-2.556 2%，最高精度可達到99.996 8%；4個不同坡度下的平均相對誤差分別為0.005 0%、-0.251 3%、-0.353 9%、-0.396 5%。對相對誤差進行頻率分布分析(圖6a)，發現大部分觀測的相對誤差都較小且在0附近分布，相對誤差為-0.5%～1%約占85%。以上結果表明數字化攝影觀測系統對土壤侵蝕坡面的幾何尺寸的觀測是準確的，且坡度對該系統的觀測準確度無顯著影響。

表2 同雨強不同坡度觀測結果Tab.2 Observation results of same rain intensity and different slopes

圖6 測量值與實際值相對誤差分布Fig.6 Distributions of relative error between measured and actual values

(2)同坡度不同雨強觀測

對數字近景攝影觀測系統在坡度為10°，雨強為30、60、90、120 mm/h條件下的觀測準確度進行檢測(表3)。數字近景攝影觀測系統得測量值與實際值兩者之間的最大相對誤差為-2.968 3%，最高精度可達到99.990 1%；不同雨強條件下觀測的平均相對誤差分別為-0.495 8%、-0.353 9%、-0.475 1%、-0.637 6%。對相對誤差進行頻率分布分析(圖6b)，發現大部分觀測的相對誤差在0附近分布，相對誤差為-0.5%～1%約占75%。以上結果表明數字化攝影觀測系統對土壤侵蝕坡面的幾何尺寸的觀測是準確的，且雨強對該系統的觀測準確度無顯著影響。

表3 同坡度不同雨強觀測結果Tab.3 Observation results of same slope and different rain intensities

3.3.2三維激光掃描法

采用激光掃描儀觀測降雨前后侵蝕坡面的三維數字地形，解算觀測所得到的土壤侵蝕量(表4)。結果表明，數字近景攝影測量的土壤侵蝕總量為452 180 cm3，激光掃描觀測的土壤侵蝕總量為407 971.36 cm3，由此可得數字近景攝影測量和激光掃描相對于徑流泥沙法的土壤侵蝕總量觀測誤差分別為3.87%和6.28%，數字近景攝影測量系統可以更加精確地量化侵蝕坡面的土壤侵蝕量。激光掃描儀由于受掃描視角的限制，測量時存在掃描盲區、漏測的問題；而數字攝影觀測系統由于采用多影像采集器組網技術，增大了影像采集的視野范圍，可采集到足夠數量溝道底部、溝壁的數字影像，從而彌補了激光掃描法在溝道觀測時存在的數據缺失等缺陷。

表4 數字攝影測量與激光掃描觀測對比Tab.4 Comparison of digital photogrammetry and laser scanning observation

3.3.3徑流泥沙觀測法

采集徑流泥沙全樣，是觀測土壤流失量最為可靠的方法，在本次檢測中以徑流泥沙觀測法得到的結果作為實際值。采用數字攝影觀測系統將降雨過程中不同時間點的數字影像進行體積解算得到土壤流失量結果，并與相同時間段收集的徑流泥沙含量進行對比(表5)。結果表明，兩種觀測方法的平均相對誤差為-1.73%；從不同侵蝕階段兩種觀測方法的相對誤差來看，在降雨初期，兩種觀測方法相對誤差比較大，在降雨歷時達到50 min時，觀測精度開始變高，說明此時正是坡面溝道快速發育的明顯分界點，溝道快速發育前后兩種觀測方法觀測的相對誤差分別為20.85%和3.47%；隨著降雨歷時的延長，土壤坡面侵蝕溝發育形態的變化越來越明顯，數字近景攝影觀測系統的觀測精度逐漸提高，觀測精度最高可達99.26%。

數字點云密度代表著數字點云對坡面地表形態表達的準確程度，高密度的數字點云可以更加精確地表達坡面侵蝕溝的形態信息。通過不同時間節點侵蝕坡面數字點云的數量和坡面觀測面積計算數字點云密度。結果表明，數字點云的平均數量為1.335 1×106個，平均點云密度為0.134個/mm2。表5中列出了由不同時間節點侵蝕坡面的高密度數字點云轉換生成的DEM，其空間分辨率可達到2 mm，可準確表達侵蝕形態的空間變化，實現了對土壤侵蝕坡面形態變化過程的動態監測。

表5 不同時間點土壤流失量的觀測結果Tab.5 Observation results of soil loss at different times

4 結論

(1)提出了一種在連續降雨過程中對土壤侵蝕坡面動態監測的數字近景攝影觀測方法，通過無線組網技術，并行拍攝、解算降雨過程中下墊面的數字影像，提取坡面精細地貌動態變化的信息。該系統時間觀測分辨率可達到分鐘級別，空間分辨率達到2 mm。可從時間和空間尺度上更加準確地描述土壤侵蝕過程，解決了侵蝕觀測中時空不一致的問題。

(2)與傳統徑流泥沙法的平行觀測結果比較表明，數字化攝影觀測法在坡面土壤侵蝕過程的不同階段其準確度不同，隨降雨歷時的延長，數字近景攝影觀測系統的觀測精度逐漸提高，土壤流失量估算平均相對誤差為-1.73%，單次觀測精度最高可達99.26%。

(3)與激光掃描法平行觀測結果比較表明，數字近景攝影測量系統對坡面土壤侵蝕量的觀測精度高于激光掃描儀。數字近景攝影觀測法采用多影像采集器組網技術，增大了影像采集的視野范圍，克服了激光掃描儀觀測時溝道底部激光線不能投射到位而造成的漏測現象，實現了全覆蓋的數字影像采集。