徐豐,李恒凱
(江西理工大學 土木與測繪工程學院,江西 贛州 341000)
由于光學傳感器在獲取影像數據時,容易受到許多不確定因素的影響,如氣候條件、成像時間、成像角度等,因此在長時間序列對地觀測研究中,經常難以獲得理想的影像,造成某一時間段的影像數據缺失,必須用其他傳感器的影像來彌補[1-2],尤其在多云多雨地區,由于回訪周期、天氣等因素的限制,單一來源的影像數據有較大局限性[3]。此外,不同年代的傳感器相互補充構建的長時序對地觀測數據集,也已被廣泛應用于環境演變監測。因此,針對不同傳感器的交互定標一直是遙感領域的研究熱點,檢查不同衛星信息之間的一致性和相互關系極為重要[4-5]。
針對傳感器的交互定標,國內外學者開展了大量的研究。胡昌苗[6]在輻射歸一化研究方面,以Landsat TM/ETM+及 HJ-1 A/B CCD數據為例,提出了相對和絕對輻射校正一體化處理方案,將 DN 值轉換到地表反射率。不少學者為了能讓波段間光譜特征盡可能保持一致,根據每個波段中的 DN 值計算各種類別像素的轉換系數[7-8],或通過最小二乘法對大量多源時序影像進行回歸,給出了表觀反射率、地表反射率的回歸系數[9-10]。目前,有研究針對不同的傳感器獲取的主要植被指數如歸一化植被指數(normalized difference vegetation index,NDVI)、土壤調整植被指數(soil-adjusted vegetation index,SAVI)、增強植被指數(enhanced vegetation index,EVI)[11-13]進行交互對比,查明之間的定量關系,求出其轉換方程。這些研究充分表明:不同傳感器數據可以通過構建數學方法,進行交叉定標。但是,不同數據傳感器數據之間有較大差異,需要開展有針對性的研究[14-15]。Landsat系列影像具有較為豐富的歷史存檔數據,是長時序對地監測的重要數據來源。但Landsat-5傳感器和Landsat-8傳感器在時間上并不能無縫銜接,有幾年空檔,另外兩種傳感器本身性能也有較大差異。HJ-1B CCD傳感器與Landsat系列影像具有相同空間分辨率,回訪周期僅兩天,且其數據跨越TM和OLI傳感器,利用其銜接和補充TM和OLI傳感器具有可行性[16]。基于此,本文結合NDVI植被指數和地表反射率,通過對同日過空的三對(六幅)Landsat TM/OLI和HJ-1B CCD影像對進行比較,分析不同波段的觀測能力,并建立其回歸關系。本研究有利于三者觀測數據的互相補充,為構建多源長時序遙感數據集提供重要的數據源支撐。
本文選取的Landsat TM/OLI和HJ-1B CCD傳感器在許多方面都較為相似,如HJ-1B CCD的四個波段與Landsat-5 TM的1~4波段以及Landsat-8 OLI的2~5波段的光譜范圍基本相同,相應為藍、綠、紅、近紅外波段;三者過境時間較為接近;HJ-1B CCD、Landsat TM/OLI的傳感器參數基本一致,空間分辨率都是30 m;軌道傾角也頗為接近,分別為98.2°、98.2°、98°[17-18]。因此,HJ-1B 與Landsat系列是滿足影像對比分析的遙感儀器。然而,Landsat-5 TM傳感器有七個多光譜波段,Landsat-8 OLI有九個多光譜波段,HJ-1B CCD只有四個多光譜波段(藍、綠、紅、近紅外波段),光譜響應函數曲線不一致,為了準確對比三種衛星傳感器的對地觀測能力,重要的前提條件就是獲取二者同日過境的遙感影像對[19-20]。此外,為了避免實驗結果的偶然性,不能僅用一組影像對來進行交互對比。但由于南方多云天氣較多,以及運行周期的不一致,導致部分數據只能選取日期最為接近的影像對[21]。因此,本研究選用了滿足以上條件的三對同日過境、相同區域的Landsat TM/OLI和HJ-1B CCD影像對。同時,選取2016年3月1日的影像數據進行驗證對比。選取的研究區域為江西省贛州市定南縣,地理坐標為:

表1 實驗影像對參數
114°47′49″E~115°22′48″E,24°33′37″N~25°03′21″N,行政面積為 1 318.72 km2,境內屬低山丘陵地形,植被覆蓋率在82%以上,具有耕地、林地、礦區、城鎮等多種土地利用類型。
在進行不同傳感器數據的交互比較時,需要對其進行輻射校正,以使對比影像之間的輻射信號能相互一致。不同傳感器之間由于光譜波段和波長之間設計的不同,使得如果僅將DN(digital number)值反演成輻射值還可能會給對比結果帶來不確定性,因此還必須進一步將其反演成反射率。對Landsat-5 TM、HJ-1B遙感影像進行輻射定標,定標參數從中國資源衛星應用中心獲取(http://www.cresda.com/CN/),將影像的DN值轉換成輻亮度,進一步將輻亮度轉換為反射率,通過對不同影像之間的日-地距離和太陽天頂角的差異進行正規化,以減少它們之間因日照和地形條件的不同所造成的影響(式(1))。
(1)
式中:ρλ為像元在傳感器處的反射率;DN·Gainλ+biasλ為傳感器處的輻亮度值;Gainλ和biasλ是傳感器的定標增益值和偏置值;ESUNλ為大氣頂部平均太陽輻照度,單位為 W/(m2·μm);d為日-地天文單位距離;θs為太陽天頂角。這些參數都可以從影像的頭文件中獲得。
由于Landsat-8衛星改進后,在表觀反射率反演方面與以往Landsat系列有較大的差異,減少了d、ESUNλ等參數的計算,可以直接采用式(2)得到傳感器處的表觀反射率[22]。
ρλ=(Mρ·Qcal+Aρ)∕cosθs
(2)
式中:Qcal為影像的DN 值;Mρ為波段λ的反射率調整因子;Aρ為波段λ反射率調整參數,可以從影像頭文件中獲得。對Landsat level1和HJ-1B level 2數據進行大氣校正,以Landsat TM/OLI為基準影像對HJ-1B影像進行幾何校正,校正后的均方根誤差(RMSE)小于0.5個像元。
在Landsat TM/OLI和HJ-1B CCD遙感圖像中采用國際慣用的樣區比較法[23-24],即在二者影像上選取范圍相同的感興趣區(region of interest,ROI)作為樣區,然后以各樣區的均值來進行對比。樣區要盡量選擇均勻的地物,避免光譜差異大的地區,可有效避免配準精度可能產生的問題。據此,本文在每對實驗影像上共選了66個均質的植被樣區(林地:25;耕地:10;建設用地:10;礦區:10;裸土:7;水體:4),每個樣區的像元數在數十個到數百個之間,共計16 442個。圖1以2013年Landsat-8影像(543真彩色合成)為例列出部分樣區。

圖1 樣區例圖
通過計算兩兩傳感器各對應波段之間的NDVI指數,地表反射率的最大值、最小值、平均值和動態范圍,標準誤差(standard deviation,Std.Dev),決定系數(R2),均方根誤差(root mean square error,RMSE),相對偏差率(mean error,ME),來衡量兩組數據之間的偏差程度。決定系數(R2)是代表各對應波段之間線性擬合程度,以Landsat TM/OLI為自變量,HJ-1B CCD為因變量進行擬合,擬合程度越高,自變量對因變量的解釋程度越高,觀察點在回歸直線附近越密集。計算如式(3)至式(4)所示。

(3)
(4)

植被指數是用來反映植被覆蓋度及其生長狀況的度量參數,有益于增強遙感影像的解譯能力。計算Landsat TM/OLI和HJ-1B CCD影像的植被指數NDVI,并進一步計算出研究區域各樣區 NDVI的統計特征值以及二者之間的ME和RMSE,然后將各樣區的NDVI均值散點投影到二維光譜特征空間上進行回歸分析,以查明它們之間的定量關系。
從表2和圖2來看,Landsat NDVI均與HJ-1B NDVI數據顯示出良好的相關性。當NDVI小于0.2時,HJ-1B數據比Landsat稍大,但是,隨著NDVI的增加,HJ-1B數據開始偏小于Landsat。散點總體均勻分布在 1∶1 線附近,它們的相對偏差率(ME)不大于10%,RMSE小于0.1,三組對比中產生的回歸方程的R2都大于0.93。ME都為正數,表明HJ-1B的NDVI數據均略大于Landsat,二者的NDVI數據十分接近。進一步分析發現,二者的NDVI數據仍存在一定的差距。
1)HJ-1B CCD的NDVI最小值、最大值、均值和動態范圍都大于Landsat TM/OLI的相應值,而標準差小于Landsat TM/OLI的相應值,這說明HJ-1B CCD獲取的植被信息量、動態范圍以及植被的可分性會大于Landsat TM/OLI的相應值,即使用HJ-1B CCD數據計算的NDVI敏感性更高,獲取的植被信息也更豐富且獲取植被信息的穩定性較好。
2)Landsat獲取的植被指數信號會強于HJ-1B,因為與HJ-1B相比,Landsat TM/OL的NDVI 表現為低估,其均值小于HJ-1B系列。在植被覆蓋率高的地方,HJ-1B影像獲得的植被指數值低于Landsat影像,而在植被覆蓋率較低的地方,情況相反。表明在植被覆蓋區HJ-1B植被指數獲得的信號弱于Landsat影像。
通過NDVI指數表明,HJ-1B影像與Landsat在植被應用方面有著相當,甚至更高的優勢。Landsat TM/OLI和HJ-1B CCD三者的植被觀測能力很接近,但在不同的獲取信息上表現各不相同。魏宏偉等[25]的研究證實發現,HJ-1B植被信號要弱于后者,主要原因是由于兩種傳感器在構成植被指數的紅光與近紅外波段的光譜響應函數不同造成的。如圖4可見,HJ-1B的光譜范圍明顯小于其他傳感器。

表2 基于 ROI 的Landsat TM/OLI和HJ-1B CCD的NDVI統計特征

圖2 不同年份Landsat TM/OLI、HJ-1B影像對的NDVI散點圖
圖3展現的是Landsat TM/OLI和HJ-1B CCD對應波段像元地表反射率散點分布(紅色實線)以及過原點的擬合直線1∶1(黑色虛線)。表3列出了各實驗影像對所有樣區所計算得出的統計特征。
由表3的統計特征可以得出,Landsat TM/OLI和HJ-1B CCD的各對應波段回歸方程的決定系數R2均大于0.8,最大值可達到0.93,表明三個傳感器的地表反射率具有較好的一致性,HJ-1B CCD在各波段的地表反射率均值基本都大于Landsat TM/OLI數據,且兩個影像的均值都隨著波段范圍的提升而提高,大小順序為:藍光<綠光<紅光<近紅外,同時,其RMSE也隨之提高。由于波段波譜差異,紅光、近紅外波段的地表反射率明顯增強,遠大于藍綠波段,近紅外的平均值遠超其他三個波段。
從圖3可以看出,Landsat TM/OLI和HJ-1B CCD傳感器在所有頻段之間都具有明顯的相關性,四個頻段的決定系數R2在0.821 0至0.932 3范圍內。各波段的變化規律基本一致,藍、綠光波段散點大部分在1∶1線下方;紅光波段散點分布較均勻,分布范圍較大,從0.02~0.4都有所分布,但都基于1∶1線對稱分布;近紅外波段的散點主要在0.2以上,地面反射率均大于其他光譜,說明該區域地表吸收的近紅外較少。由于這兩個波段具有非常不同的波長范圍,HJ-1B CCD 近紅外波段的光譜響應函數與 Landsat-8 OLI 近紅外波段光譜響應函數僅有很小的一部分相交。

圖3 不同年份Landsat、HJ-1B影像對的地表反射率散點圖
通過提取66個樣區中的16 442個像元的地表反射率值,對每一對影像進行比較,可獲得的轉換方程如圖3所示(黑色虛線為1∶1線,紅色實線為回歸方程)。

年份對應波段Landsat TM/OLIHJ-1B CCD最大值最小值均值最大值最小值均值RMSER2轉換方程藍光 0.1890.0200.0570.1920.0160.0590.014 20.825 6y=0.842 1x+0.008 42009綠光 0.2960.0290.0920.2930.0180.0750.025 40.892 5y=0.793 7x+0.000 5紅光 0.3970.0190.0990.3860.0210.0910.026 30.884 9y=0.831 5x+0.005 9近紅外0.4970.0370.2730.4990.0500.2760.032 60.932 3y=0.942 3x+0.012 8藍光 0.1670.0120.0540.1750.0360.0660.017 30.841 2y=0.995 8x+0.013 02013綠光 0.2400.0240.0790.2350.0410.0990.021 60.891 1y=1.039 8x+0.013 5紅光 0.2860.0190.0820.2960.0530.1020.025 60.845 2y=0.843 0x+0.027 0近紅外0.4400.0330.2380.3630.1010.2370.038 80.897 6y=0.851 7x+0.025 5藍光 0.2180.0150.0600.2320.0480.0840.025 30.821 0y=0.978 5x+0.012 52014綠光 0.2950.0270.0830.2810.0650.1090.028 70.860 0y=1.062 0x+0.017 7紅光 0.2890.0170.0810.2910.0430.0980.027 50.827 7y=0.904 9x+0.020 3近紅外0.5690.0240.2970.5840.0770.3760.041 20.926 9y=1.152 7x+0.024 9
每兩個傳感器的回歸線y越靠近1∶1線時,則表明二者的反射率差異越小;回歸線位于1∶1線上方時,說明HJ-1B的反射率高于Landsat的反射率;回歸線位于其下方時則反之。從圖3中的回歸方程可看出,Landsat-5 TM與HJ-1B CCD的回歸線都在1∶1線下方,說明Landsat-5 TM比HJ-1B CCD接收的地表反射率更高;而Landsat-8 OLI與HJ-1B CCD的回歸線在1∶1線上方,說明HJ-1B CCD比Landsat-8 OLI接收的地表反射率更高。三個傳感器在藍、綠和紅光波段的RMSE和ME差異程度較一致,在近紅外波段則表現出較大差異,其地表反射率和R2均大于藍、綠和紅波段,表明三種傳感器信號差異在不同波段的表現并不相同。近紅外波段地表反射率差異為0.2左右,R2的差異為0.3~0.5,說明三個傳感器的地表反射率散點在近紅外波段的協調性優于其他三個波段。對于Landsat TM/OLI,NIR波段反射率與HJ-1B CCD產生極為顯著的相關性。
上述對比結果表明,在植被觀測能力上,Landsat TM/OLI傳感器相較于 HJ1B-CCD 具有較強的植被探測能力;在對地觀測能力上,HJ-1B CCD獲取的地表反射率總體高于Landsat系列數據,但這兩個系列數據在不同波段表現有所差異,差異的原因可能有以下方面。
1)光譜響應函數差異。如圖4所示,Landsat TM/OLI和HJ-1B CCD在紅光、近紅外波段之間的光譜響應函數存在一定的差異。HJ-1B的光譜響應在紅光、近紅外的波長范圍均大于Landsat TM/OLI,特別是近紅外,峰值、范圍都有較大差異,這是Landsat TM/OLI和HJ-1B CCD近紅外波段在對地觀測能力上RMSE和R2差值最大的原因之一。但是,與傳感器的天頂角相比,光譜響應函數對Landsat TM/OLI和HJ-1B CCD之間的差異影響相對較小。

圖4 Landsat TM/OLI與HJ-1B CCD光譜響應函數曲線
2)傳感器天頂角的影響。隨著傳感器天頂角度的增加,各波段反射率和NDVI數據隨之變化。當傳感器天頂角上升到大約34°然后下降時,三個可見波段的反射率會變大。但是,當傳感器的天頂角小于45°時,NIR波段反射率始終會增加,NDVI數據變小,直到傳感器天頂角上升到34°,NDVI變大,這意味著超過34°的傳感器天頂角對NDVI、反射率值產生很大影響。在0°和17°之間的傳感器天頂角差將帶來9.54%的反射率誤差,而光譜誤差僅導致0.53%的誤差。與每個頻段的反射率相比,NDVI受傳感器天頂角度的影響很小。
3)傳感器的定標精度差異。研究表明,Landsat TM/OLI和HJ-1B CCD自發射以來,其傳感器的性能都發生了衰減。特別是Landsat-5發射時間過長,各種觀測能力下降較為嚴重。HJ-1B的定標文件不確定,在中國資源衛星應用中心官網,影像頭文件XML里面的定標系數各有一套定標文件,定標系數不一致。
為了檢驗基于上述方法所得出的Landsat TM/OLI與HJ-1B的關系方程的精確度,本文選取了2016年3月1日同日過境的Landsat-8和HJ-1B的影像對進行驗證。
1)NDVI植被指數驗證。將研究區域2009、2013、2014年的66個樣本合并起來進行回歸分析,以此建立Landsat TM/OLI和HJ-1B CCD 植被指數NDVI的統一轉換,表達如式(5)所示。
HJ-1BNDVI=1.020 2LandsatNDVI+0.007
(5)
首先,對驗證影像對進行預處理,保證與另三組影像進行了相同的預處理;其次,在Landsat-8的驗證影像上選取了同質的66個樣本,將轉換方程代入Landsat的驗證影像中,計算模擬HJ-1B的NDVI值;然后,與實際的HJ-1B的NDVI值進行回歸分析(圖5),并計算出它們的ME與RMSE。對照圖6可看出,在經過轉換后,Landsat和HJ-1B植被指數NDVI 的散點擬合曲線更加靠近 1∶1 線,差距在變小;R2也較之前提高了3~4個百分點,達到0.991 8;ME降低為0.36%,取三年的ME平均值為0.65,降幅為80.5%;RMSE降低為0.02,取三年的RMSE平均值為0.046,降幅為130%。顯然,經過轉換的HJ-1BNDVI數據已經達到了較高的精度,接近于實測數據,可以通過轉換方程促進這兩種傳感器NDVI數據的協同使用,可以用于兩種傳感器的相互替代。

圖5 基于驗證影像建立的Landsat和HJ-1B NDVI轉換方程的擬合結果

圖6 基于驗證影像建立的Landsat和HJ-1B地表反射率值轉換方程的擬合結果
2)地表反射率(bottom of atmospheric reflectance,BOA)驗證。將研究區域2009、2013、2014年的66個樣本共16 442個像元值合并起來進行回歸分析,建立Landsat TM/OLI和HJ-1B CCD 各波段的地表反射率值的統一轉換,表達如式(6)所示。
藍 光:HJ-1BBOA=0.954 0LandsatBOA+0.013 2
綠 光:HJ-1BBOA=0.861 0LandsatBOA+0.013 8
紅 光:HJ-1BBOA=0.839 8LandsatBOA+0.018 3
近紅外:HJ-1BBOA=0.882 8LandsatBOA+0.021 5
(6)
與上節相同的方法,在驗證影像選取同質的66個樣本區共16 442個像元,將像元值代入相對應波段的轉換方程中,計算出HJ-1B的模擬地表反射率值,與驗證影像HJ-1B的實際值進行回歸分析,并計算它們之間的RMSE和R2,并將結果進行比較。

表4 驗證影像地表反射率值模擬轉換前后的統計特征值對比
對照圖3和圖6可以看出,經過模擬轉換后各波段的RMSE和R2較模擬前都有較顯著的優化。對比實驗影像,轉換后的地表反射率值散點擬合曲線與1∶1線幾乎重疊,R2提高了3~12個百分點,都在0.93以上,接近于1,說明轉換后具有極為顯著的線性關系。RMSE的降幅最大為60%,最小值可以達到0.007 5,這表明Landsat數據與HJ-1B數據能進行有效的轉換,在經過模擬轉換后,兩種傳感器之間的數據差異大大減小,能用于兩種傳感器的交替使用。
本研究通過對三對同日過境的影像采用基于樣區的比較法對Landsat TM/OLI和HJ-1B CCD的傳感器數據進行交互對比分析。通過分析發現,三種傳感器的NDVI數據和地表反射率值都具有較高的相關性,各對應波段的R2都處于0.82~0.95之間,RMSE也小于0.05,處于一個較高的水平。對比分析發現,三者的植被觀測能力以及對地觀測能力都存在著一定的差異,主要表現在植被觀測能力上,HJ-1B影像與Landsat在植被應用方面有著相當甚至更高的優勢,Landsat TM/OLI和HJ-1B CCD三者的植被觀測能力很接近,但在不同的獲取信息上表現各不相同。在對地觀測能力上,HJ-1B獲取的地表反射率總體高于Landsat系列數據。分析表明,Landsat系列數據與HJ-1B數據之間的差異是由于光譜響應函數、定標精度以及傳感器天頂角的影響造成的。
因此,通過波段光譜轉換,將不同傳感器類型的數據以數值回歸校正法來消除傳感器之間的波段間差異,可以實現波段間光譜轉換,從而將多源數據轉換為連續的單一類型的時序數據,以提高研究數據的一致性。由光譜轉換結果可知,經過光譜定標轉換的HJ-1B影像數據明顯提高了與TM/OLI數據的擬合度,通過光譜定標轉換方法,在精度允許范圍內數據間可以相互轉化與模擬,為存檔數據利用提供新的思路與方法。
由于本文所獲取的HJ-1B數據量有限,僅以定南縣的Landsat系列和HJ-1B圖像對擬合得到了NDVI和地表反射率的轉換方程,該方程在其他區域的適用效果還有待進一步驗證。在今后的研究中,需要使用更多區域、更長時間序列的不同影像對來檢驗該方程的穩定性。