任中杰,萬永智,王 偉,王繼龍,宋銀燕,張方方
(江蘇省水文水資源勘測(cè)局徐州分局,江蘇徐州 221000)
水深以及水下地形測(cè)量是水文研究及生產(chǎn)活動(dòng)中的重要領(lǐng)域之一,傳統(tǒng)測(cè)量方法主要以人工方式進(jìn)行長(zhǎng)序列大范圍觀測(cè),耗費(fèi)大量人力物力、成本較高。面對(duì)大范圍水域,作業(yè)時(shí)間跨度較高,河道湖泊的水量調(diào)度及水位實(shí)時(shí)變化均會(huì)影響水深數(shù)據(jù)的真實(shí)性及精度,不能反映水域真實(shí)情況。受限于自然條件,部分區(qū)域船舶難以進(jìn)入,也形成了監(jiān)測(cè)空白地帶。
隨著遙感技術(shù)的快速發(fā)展,遙感影像成本低、覆蓋面積廣且具有較高空間分辨率、時(shí)間分辨率,逐漸被廣泛應(yīng)用于資源調(diào)查、環(huán)境監(jiān)測(cè)、災(zāi)害預(yù)警、區(qū)域時(shí)空演變等領(lǐng)域[1]。其中多光譜影像是通過多光譜電磁波,獲取目標(biāo)地物豐富的光譜信息,因此被廣泛應(yīng)用于水深反演領(lǐng)域。遙感反演水深是指利用遙感數(shù)據(jù),依據(jù)可測(cè)參數(shù)值去反推水深值的水深量測(cè)方法[2],主要分為3 類:理論解析法、半理論半經(jīng)驗(yàn)法和統(tǒng)計(jì)相關(guān)法。較為常用的是統(tǒng)計(jì)相關(guān)法,需要大量實(shí)測(cè)數(shù)據(jù)支撐,但不需要分析和研究水體內(nèi)部的光學(xué)參數(shù),只考慮數(shù)據(jù)之間的相關(guān)映射關(guān)系[3]。李經(jīng)緯等[4]通過對(duì)比不同大氣校正方法對(duì)水深因子相關(guān)性的影響,判斷研究區(qū)域適宜的大氣校正模型;趙順利等[5]基于OLI 影像,利用最佳水深因子OLI-4 波段建立了水深模型,對(duì)錯(cuò)戳龍錯(cuò)鹽湖進(jìn)行水深反演;吳忠強(qiáng)等[6]利用G-S 變換對(duì)landsat-8 影像的多光譜波段及全色波段進(jìn)行融合,提高了空間分辨率,基于融合后的影像進(jìn)行水深反演。
影響水深反演精度的因素主要以水體懸浮物質(zhì)、水體底質(zhì)、水深反演因子為主。而水深反演因子的空間尺度取決于遙感影像的空間分辨率。影像融合能夠有效提高影像的空間分辨率,但目前尚未有針對(duì)影像融合方法與水深因子相關(guān)性的研究,同時(shí)基于最新的Landsat-9 影像進(jìn)行水深反演研究,尚未在云龍湖地區(qū)開展過。本文基于Landsat-9 影像,通過影像融合提高了空間分辨率,構(gòu)建水深反演模型,確定了該區(qū)域的最佳水深反演模型。
云龍湖位于江蘇省徐州市泉山區(qū),是著名的國(guó)家5A 級(jí)景區(qū)。其東靠云龍山,西依韓山、天齊山。南偎泉山、珠山。三面環(huán)山,一面臨城。云龍湖以湖中路為界,東湖周長(zhǎng)約8.1 km,西湖長(zhǎng)約7 km,全湖周長(zhǎng)約12 km,水面面積約6 km2[7]。
云龍湖水庫是奎河源頭,集水面積60 km2,總庫容3 323萬m3,是一座以城市防洪為主,兼有灌溉、養(yǎng)殖、景觀、旅游開發(fā)等綜合利用功能的中型水庫。對(duì)云龍湖水庫進(jìn)行水深反演研究具有重要的現(xiàn)實(shí)意義。
本文以云龍湖東湖水上世界以南水面區(qū)域作為研究區(qū)域(拐點(diǎn)坐標(biāo)分別為1 號(hào)點(diǎn)117.1549°、34.2412° ,2 號(hào)點(diǎn)117.1666° 、34.2394° ,3 號(hào)點(diǎn)117.1582°、34.2298°,4號(hào)點(diǎn)117.1497°、34.2324°,5號(hào)點(diǎn)117.1502°、34.2349°,6號(hào)點(diǎn)117.1479°、34.2384°),面積約1.3 km2。
Landsat-9 是美國(guó)陸地衛(wèi)星計(jì)劃(Landsat)的第九顆衛(wèi)星,2021年9月發(fā)射,每隔16d對(duì)地球進(jìn)行一次成像,與Landsat-8 存在8d 的偏移。相較于Landsat8 的12 位量化及可識(shí)別4 096 種色調(diào),Landsat-9 具有更高的輻射分辨率(14 位量化),可以識(shí)別16 384 種色調(diào),同時(shí)具有更優(yōu)的信噪比。Landsat-9 共有11個(gè)波段,波段簡(jiǎn)介見表1。
為避免水位變化的影響,本文選取了與實(shí)測(cè)作業(yè)時(shí)間基本一致的Landsat-9 影像數(shù)據(jù),產(chǎn)品號(hào)為L(zhǎng)C09_L1TP_121036_20221107_20221129_02_T1。數(shù)據(jù)來源于空天院遙感數(shù)據(jù)服務(wù)系統(tǒng)(http://eds.ceode.ac.cn/nuds/test)。研究區(qū)域影像清晰,無云。對(duì)影像進(jìn)行預(yù)處理,主要包括幾何校正、影像裁切。
本文根據(jù)研究區(qū)域形狀,大致東西方向布設(shè)測(cè)量航線56條,各航線間隔25 m。將布設(shè)的航線導(dǎo)入到測(cè)深儀導(dǎo)航軟件中。將高頻測(cè)深儀中海達(dá)HD-27T、華測(cè)GPS X90 固定在船邊沿,保持測(cè)桿垂直水面。將記錄間距設(shè)置為10 m,按照導(dǎo)航軟件的指示駕駛船只按航線方向行進(jìn),軟件自動(dòng)記錄所在位置和水深。水位數(shù)據(jù)來源于云龍湖水文閣水位站現(xiàn)有水尺觀測(cè)的每日水位數(shù)據(jù)。采集數(shù)據(jù)期間,云龍湖水位穩(wěn)定在32.65 m,風(fēng)力為東南風(fēng)2級(jí)。
對(duì)測(cè)得的水深數(shù)據(jù)進(jìn)行檢查,剔除異常值。由于Landsat-9 多光譜影像分辨率為30 m,岸邊區(qū)域存在混合像元,影響水深精度,因此對(duì)岸邊實(shí)測(cè)點(diǎn)進(jìn)行剔除。從8 392 個(gè)點(diǎn)位中挑選出765 個(gè)點(diǎn)位作為訓(xùn)練樣本數(shù)據(jù),另挑出661 個(gè)點(diǎn)位作為檢驗(yàn)樣本數(shù)據(jù),挑選原則為點(diǎn)位均勻分布,具有足夠代表性,訓(xùn)練樣本和檢驗(yàn)樣本無重疊點(diǎn)位。
2.1.1 單波段相關(guān)性
選取Landsat-9 影像中多光譜數(shù)據(jù)B1,B2,B3,B4,B5,B6,B7 共7 個(gè)波段,將像元灰度值提取至點(diǎn)位與訓(xùn)練樣本水深進(jìn)行相關(guān)性分析。由于影像空間分辨率為30 m,同一個(gè)像元對(duì)應(yīng)多個(gè)實(shí)測(cè)點(diǎn),像元中心值不能反映像元區(qū)域內(nèi)部水下地形的差異性,因此采用雙線性內(nèi)插法計(jì)算點(diǎn)位對(duì)應(yīng)的灰度值,結(jié)果顯示,相較于像元中心值,水深相關(guān)性有一定提升。
從表2 可以看出,相關(guān)系數(shù)絕對(duì)值較高的波段為B5、B6、B7。理論上,透水波段應(yīng)為藍(lán)、綠波段,相關(guān)系數(shù)最高,但僅建立在水質(zhì)較清且水底反射率較高的前提下。由于不同水體的泥沙含量、葉綠素濃度各有不同,所以與水深相關(guān)系數(shù)較高的波段也會(huì)有所差異,這種現(xiàn)象被稱為“紅移”[8]。
表2 波段灰度值與水深相關(guān)系數(shù)
2.1.2 影像融合
Landsat-9 多光譜影像每個(gè)像元對(duì)應(yīng)900 m2的水域面積,將多光譜及全色影像進(jìn)行融合,能夠獲得高空間分辨率的多光譜影像,能有效降低像元尺度,能細(xì)致反映水下地形變化。本文選用多種融合方法,對(duì)融合影像進(jìn)行相關(guān)性分析,選擇水深相關(guān)系數(shù)最高的融合方法。
常用的融合方法有:Brovey、Gram-Schmidt、HSV、NNDiffuse、Pansharp、PCA、SFIM 融合方法等[9]。其中,Brovey、HSV、SFIM 對(duì)波段數(shù)存在限制,因此,僅選擇相關(guān)系數(shù)較高的B5、B6、B7 波段參與融合。對(duì)融合結(jié)果進(jìn)行相關(guān)性分析,結(jié)果如表3所示。
表3 不同融合方法水深相關(guān)系數(shù)
從表3 可知,相較于原始影像,SFIM 融合方法相關(guān)性保持較好,各波段沒有明顯降低,雖然有一定的光譜信息損失,但有效提高了影像空間分辨率;其次是NNDiffuse融合方法,僅B5波段相關(guān)性保持較好,其他融合方法均明顯降低了水深相關(guān)性。故本文采用SFIM 融合方法對(duì)B5、B6、B7 以及全色波段B8進(jìn)行融合。
2.1.3 水深因子確定
研究表明,單波段模型水深反演精度遠(yuǎn)不如多波段模型[10],且用對(duì)數(shù)因子做比值可以進(jìn)一步提高數(shù)據(jù)相關(guān)性。本文嘗試用對(duì)數(shù)因子、對(duì)數(shù)比值因子、比值因子等多種組合進(jìn)行相關(guān)性分析,尋求最佳水深因子,分析結(jié)果如表3所示,僅列出相關(guān)系數(shù)較高的4項(xiàng),如表4所示。相關(guān)系數(shù)較高的4個(gè)水深因子均包含B5 波段,其中B5+B6 相關(guān)系數(shù)最高,相較于B5波段有一定的提升。
基于最佳水深因子B5+B6進(jìn)行回歸性分析,分為指數(shù)回歸、傅里葉回歸、高斯回歸、多項(xiàng)式回歸(1次、2次)。
從表5可以看出,反演模型決定系數(shù)均較差,一方面是由于湖泊水體懸浮物相較于清澈海域較多,導(dǎo)致波譜紅移,影響反演精度,另一方面是由于影像空間分辨率仍然很低,通過影像融合將像元對(duì)應(yīng)水域面積900 m2降低至225 m2,仍然不能真實(shí)反映水下地形變化起伏情況,實(shí)測(cè)水深數(shù)據(jù)平均250 m2采集一個(gè)點(diǎn)位,但很難保證點(diǎn)位處于像元中心值,因此出現(xiàn)較多異常值。常規(guī)多項(xiàng)式擬合一般采用最小二乘法(LSM),使殘差平方和極小,計(jì)算簡(jiǎn)便,然而穩(wěn)健性(抗粗差性)較差。為了避免異常值對(duì)反演模型的影響,采用最小絕對(duì)殘差(LAR)對(duì)數(shù)據(jù)進(jìn)行擬合,最優(yōu)擬合結(jié)果見圖1。
圖1 二次多項(xiàng)式LAR擬合
表5 水深反演模型
表6 反演模型精度
RMSE 僅為0.0253 m,且決定系數(shù)顯著提升。反演模型公式為
式中,p1為-1.564e-06;p2為0.03792;p3為-227.1。
通過水深反演模型進(jìn)行波段運(yùn)算得到水深反演灰度圖,見圖2。
利用未參與水深反演模型計(jì)算的661個(gè)實(shí)測(cè)水深點(diǎn)作為驗(yàn)證組,對(duì)反演模型進(jìn)行精度評(píng)價(jià)。反演組(FYZ)和實(shí)測(cè)組(SCZ)回歸,如圖3所示。
從圖3及表7可以看出,基于最優(yōu)水深因子B5+B6 的水深反演模型反演組與實(shí)測(cè)組具有相關(guān)性。RMSE 小于0.1 m,得益于湖底地形起伏波動(dòng)不大,降低了反演誤差。從定性角度看,雖然總體精度不高,但在一定程度上真實(shí)反映了水下地形,湖南路一側(cè)及荷風(fēng)島區(qū)域水深較淺,東湖中部區(qū)域水深較深,符合實(shí)際情況。
表7 水深因子相關(guān)系數(shù)
本文基于Landsat-9影像的多光譜及全色波段,結(jié)合實(shí)測(cè)水深數(shù)據(jù)進(jìn)行對(duì)比分析,確定了適用于云龍區(qū)東湖南部區(qū)域水深反演的最佳影像融合方法為SFIM融合算法,該算法能夠在有效降低水深反演像元尺度的情況下較好保持影像光譜信息。確定了Landsat-9 影像水深相關(guān)性最高的波段為B5,最佳水深因子為B5+B6,并建立了多種水深反演模型,效果最好的反演模型為基于最小絕對(duì)殘差(LAR)的二次多項(xiàng)式模型。
采用實(shí)測(cè)水深數(shù)據(jù)進(jìn)行反演精度評(píng)價(jià),結(jié)果表明反演結(jié)果能應(yīng)用于實(shí)際,影像獲取難度低、計(jì)算簡(jiǎn)單,彌補(bǔ)了人工觀測(cè)的不足。但從評(píng)價(jià)結(jié)果看,反演精度有限,空間分辨率較低,反演模型受湖泊水質(zhì)因素及水位變化影響較大,最佳水深因子及反演模型隨時(shí)可能變化,后續(xù)有必要針對(duì)反演模型適用性開展進(jìn)一步研究。