石 希，夏軍強(qiáng)**，孫 健

(1：武漢大學(xué)水資源與水電工程科學(xué)國家重點(diǎn)實(shí)驗(yàn)室，武漢 430072)

(2：清華大學(xué)水沙科學(xué)與水利水電工程國家重點(diǎn)實(shí)驗(yàn)室，北京 100084)

水體溫度是河流生態(tài)系統(tǒng)中的重要參數(shù)，其調(diào)控著大量物理、化學(xué)和生物過程，例如：水生生物的遷徙和營養(yǎng)物質(zhì)的循環(huán)等[1-3]. 長期觀測河流水溫對于監(jiān)管河流熱動態(tài)和維護(hù)河流生態(tài)系統(tǒng)的健康起著至關(guān)重要的作用[4]. 同時(shí)，分析長時(shí)間下水溫變化趨勢也有利于探知全球變暖和人類活動對河流水溫的影響[5]. 目前，傳統(tǒng)水溫觀測手段主要為利用溫度測量儀器實(shí)地探測測站附近的河流水溫[6]. 然而由于水文測站分布較稀疏，大部分河段缺乏實(shí)測水溫資料. 對于這些無測站地區(qū)，獲取滿足研究需求的水溫?cái)?shù)據(jù)往往需要花費(fèi)大量人力物力建設(shè)水文觀測站網(wǎng)并周期性地開展現(xiàn)場工作，常理上難以實(shí)現(xiàn)[7].

相較于傳統(tǒng)觀測手段，利用熱紅外遙感技術(shù)可以獲取空間上連續(xù)的水溫?cái)?shù)據(jù)[8]. 近20年來，隨著衛(wèi)星熱紅外傳感器在空間分辨率和輻射探測精度等技術(shù)方面上的不斷進(jìn)步，熱紅外遙感影像不再單單用于反演海洋、湖泊等大型水體水溫[9-10]，而逐步被應(yīng)用于河流、河道型水庫等中小尺度水體[11-12]. 但是受限于傳感器硬件性能和衛(wèi)星軌道設(shè)計(jì)等問題，衛(wèi)星遙感數(shù)據(jù)普遍只能滿足高時(shí)間分辨率或高空間分辨率兩者中的一項(xiàng)[13]. 因此對于河流和河道型水庫而言，使用契合研究需求的衛(wèi)星至關(guān)重要[14-15].

在探測河流及河道型水庫水溫方面，Landsat系列衛(wèi)星搭載的Landsat 7 ETM+、Landsat 8 OIRS和已經(jīng)失效的Landsat 5 TM傳感器憑借著其60、100和120 m的熱紅外波段中等空間分辨率和近20年穩(wěn)定的歷史數(shù)據(jù)庫這兩點(diǎn)優(yōu)勢使用最為廣泛[16]. 其反演所得水溫已被應(yīng)用于模擬河流日水溫[17]、評估河流熱污染狀況[18]、分析建造水庫對河流上下游水溫影響[19]等各項(xiàng)研究中. 而對于另一個(gè)常用的衛(wèi)星傳感器MODIS而言，由于其熱紅外波段僅具有1000 m的低空間分辨率，目前只被應(yīng)用于河道型水庫、湖泊和長江中下游干流等河道較寬的水體[20-22]. 但是相較于Landsat系列衛(wèi)星16 d的時(shí)間分辨率，MODIS的重訪周期僅為1 d，且MODIS具有較高的輻射/溫度探測精度NEΔT=0.05℃(在27℃時(shí)的等效溫度探測精度，Landsat NEΔT=0.4℃)[23]. 因此對于不同的研究區(qū)域、研究尺度選擇合適的衛(wèi)星尤為關(guān)鍵. 同時(shí)，目前相關(guān)研究都僅聚焦于評估單顆傳感器的反演精度，但研究結(jié)果表明單顆傳感器往往具有一定的局限性[24]. 因此，利用不同衛(wèi)星傳感器獲取相同河段的水溫并橫向?qū)Ρ绕渚?，可以更直觀地展示傳感器之間的差異，也有助于評估不同傳感器在相應(yīng)流域中的實(shí)用性.

除此之外，使用不同的大氣校正方法和大氣校正參數(shù)對反演水溫的精度也有一定影響. 通常情況下，具有單一熱紅外波段的衛(wèi)星傳感器，例如Landsat 5 TM 與Landsat 7 ETM+, 在大氣頂端探測到的熱紅外輻射可以通過兩種大氣校正方法：輻射傳遞模型法和JM&S (Jimenez-Munoz & Sobrino)單通道算法轉(zhuǎn)化為水體溫度[25]. 以往研究表明，由于兩種方法對大氣水汽含量的敏感性有所差異，使用兩者在具有不同氣候條件的研究區(qū)域內(nèi)反演所得的水溫精度各有高低[26-27]. 因此在大范圍反演任意流域水溫之前，應(yīng)盡可能地使用兩種方法得到的水溫與實(shí)測水溫比對并分析兩者在該研究區(qū)域內(nèi)的適用性. 另一方面，兩種反演方法均需要輸入額外的大氣校正參數(shù)，這些大氣校正參數(shù)通常由不同的大氣模型輔以實(shí)測氣象資料模擬所得[25]. 查閱相關(guān)資料發(fā)現(xiàn)當(dāng)大氣校正參數(shù)處于異常范圍時(shí)會嚴(yán)重影響校正精度[28]. 因此，在評估大氣校正方法精度的同時(shí)，進(jìn)行大氣校正參數(shù)敏感性分析也尤為重要.

本文選取長江上游流域內(nèi)具有不同河道寬度以及氣候環(huán)境條件的3個(gè)研究區(qū)域：向家壩壩前、三峽壩前和三峽庫位寸灘，分析使用不同的遙感數(shù)據(jù)、大氣校正方法以及大氣校正參數(shù)對水溫精度的影響：(1)遙感數(shù)據(jù)：MODIS和Landsat；(2)大氣校正方法：輻射傳遞模型法和JM&S單通道算法；(3)大氣校正參數(shù)：Atmcorr(atmosphere correction/大氣校正參數(shù)數(shù)據(jù)庫)和GEOS-5 FP-IT(goddard earth observing system-version 5-forward processing for instrument teams/戈達(dá)德地球觀測系統(tǒng)模型-版本5-儀器團(tuán)隊(duì)版本)，并對大氣校正參數(shù)進(jìn)行了敏感性分析. 此外，本文通過對比各類反演水溫精度以及適用條件，為今后利用熱紅外遙感數(shù)據(jù)反演長江上游流域其他河段水溫提供了可供參考的方案.

1 研究區(qū)域與數(shù)據(jù)

1.1 研究區(qū)域

本研究聚焦于長江上游流域，包括金沙江及長江至三峽段，共選取3個(gè)研究區(qū)域：(a)向家壩壩前、(b)三峽壩前和(c)寸灘(圖1). 其中向家壩水庫位于云南省和四川省交界處，坐落于金沙江下游河段，距宜賓市33 km. 選取壩上測站作為研究區(qū)域內(nèi)基點(diǎn)，壩前河道寬度在800～1000 m之間. 其屬中亞熱帶性季風(fēng)氣候，氣候濕潤，常年氣溫處于0～ 30℃. 三峽大壩坐落于湖北省宜昌市，是長江上游和中游的分界線. 本文選取大壩下游測站黃陵廟站作為研究區(qū)域內(nèi)基點(diǎn)，壩前河寬為2～3 km，黃陵廟站附近河寬在550～750 m. 其屬亞熱帶季風(fēng)氣候，秋季多雨，常年氣溫徘徊于5～35℃，夏季大氣水汽含量約為4 g/cm2. 寸灘位于重慶市，屬于三峽水庫變動回水區(qū). 選取寸灘水文站作為研究區(qū)域內(nèi)基點(diǎn)，寸灘站距離上游嘉陵江與長江交匯處約7 km，距上游向家壩約406 km，距下游三峽大壩約570 km. 河道寬度處于600～850 m. 其氣候類型、氣溫與大氣水汽含量與三峽壩前相似.

圖1 研究區(qū)域示意圖

1.2 實(shí)測水溫?cái)?shù)據(jù)

本研究共選用3個(gè)水文測站的實(shí)測水溫?cái)?shù)據(jù)：向家壩壩上測站、黃陵廟站和寸灘站，分別對應(yīng)上述3個(gè)研究區(qū)域(圖1). 其中向家壩壩上測站實(shí)測水溫?cái)?shù)據(jù)可用年份為2015年，寸灘站與黃陵廟站實(shí)測水溫?cái)?shù)據(jù)范圍為2004-2016年. 3個(gè)測站的水溫皆測量于每日上午8：00，所有實(shí)測數(shù)據(jù)皆來源于水文年鑒.

1.3 衛(wèi)星遙感數(shù)據(jù)

本文選取Landsat 5 TM與Landsat 7 ETM+熱紅外遙感數(shù)據(jù)以及MODIS溫度數(shù)據(jù)產(chǎn)品 MOD11. 而對于Landsat 8衛(wèi)星，由于其熱紅外傳感器收到雜散光影響精度較低，因此并未使用. 3種衛(wèi)星遙感數(shù)據(jù)的基本信息如下.

1.3.1 Landsat熱紅外數(shù)據(jù) Landsat 5和Landsat 7熱紅外數(shù)據(jù)均具有30 m空間分辨率(原始分辨率為120 m和60 m)和16 d時(shí)間分辨率，屬于Level-1級遙感影像. 對于Landsat熱紅外數(shù)據(jù)，結(jié)合圖1中各研究區(qū)域所處景號以及實(shí)測水溫可用年份，共選定5組數(shù)據(jù)，如表1所示. 其中由于Landsat 5從2011年起不再提供有效數(shù)據(jù)，因此所有Landsat 5 TM熱紅外數(shù)據(jù)的截止年份為2011年. 此外，不同于Landsat 5 TM熱紅外波段B6，Landsat 7 ETM+的熱紅外波段B6具有兩個(gè)通道：High gain和Low gain(高/低增益模式). 前者具有更高的探測精度，后者具有更大的溫度/輻射探測范圍[29]. 查閱Landsat 7用戶使用手冊后發(fā)現(xiàn)，由于正常河流水溫不會導(dǎo)致傳感器飽和，因此選取具有更高探測精度的High gain通道. 針對遙感數(shù)據(jù)質(zhì)量，選擇高質(zhì)量的T1級別數(shù)據(jù)，同時(shí)對于因掃描線校正器失效的Landsat 7熱紅外影像，本文選擇具有空白條帶的SLC-off數(shù)據(jù)以避免修正帶來的誤差. 所有滿足條件的數(shù)據(jù)皆從美國地質(zhì)調(diào)查局地球探索數(shù)據(jù)庫(https://earthexplorer.usgs.gov/)中免費(fèi)下載. 下載的Landsat遙感影像均為Geotiff文件，并采用了WGS84橢球坐標(biāo)基準(zhǔn)以及UTM投影，數(shù)據(jù)獲取時(shí)間約為每日上午11：00.

表1 所用Landsat熱紅外數(shù)據(jù)基本信息

1.3.2 MODIS溫度數(shù)據(jù) MODIS數(shù)據(jù)庫中具有兩個(gè)地表溫度數(shù)據(jù)集：MOD11和MOD28，分別為陸地表面溫度與海洋和大型水體溫度[26]. 經(jīng)查閱后發(fā)現(xiàn)，MOD28數(shù)據(jù)集并不涵蓋長江上游流域. 因此本次研究只選用MOD11溫度數(shù)據(jù)集. MOD11A1數(shù)據(jù)為MOD11數(shù)據(jù)集中的原始數(shù)據(jù)，其具有1 d的時(shí)間分辨率和1 km的空間分辨率(精確的網(wǎng)格大小為0.928 km)，屬于Level-3級別遙感影像. 其主要采用劈窗法對MODIS的31、32熱紅外波段原始數(shù)據(jù)進(jìn)行了大氣校正，并結(jié)合土地類型和云層覆蓋數(shù)據(jù)生成了晴天狀態(tài)下的地表溫度. 本文選用覆蓋3個(gè)研究區(qū)的MODIS溫度產(chǎn)品，時(shí)間跨度選取2004-2016年，所有滿足條件的數(shù)據(jù)皆下載于LAADS DAAC(大氣和一級遙感數(shù)據(jù)存儲和查詢系統(tǒng)分布式存檔中心)數(shù)據(jù)庫(https://ladsweb.modaps.eosdis.nasa.gov/search/). 下載的MODIS溫度數(shù)據(jù)格式為HDF-EOS，數(shù)據(jù)獲取時(shí)間約為每日上午10：00.

2 研究方法

根據(jù)所用數(shù)據(jù)不同，研究方法分為兩部分(圖2)：(1)Landsat熱紅外數(shù)據(jù)處理；(2)MODIS溫度數(shù)據(jù)處理.

圖2 研究方法流程圖

2.1 Landsat熱紅外數(shù)據(jù)處理

Landsat所提供的原始熱紅外數(shù)據(jù)為0～255的亮度值，并不能直接反映輻射強(qiáng)度. 所以，首先需將原始數(shù)據(jù)通過輻射定標(biāo)轉(zhuǎn)化為傳感器在大氣頂端所探測到的光譜輻射亮度，公式為：

Lsensor=gain·DN+bias

(1)

式中，Lsensor為Landsat衛(wèi)星傳感器在大氣頂端接收到的熱紅外光譜輻射亮度/輻射值，W/(sr·m2·nm)；DN為原始熱紅外數(shù)據(jù)中的亮度值；gain和bias為定標(biāo)參數(shù)，Landsat 5和Landsat 7熱紅外波段的定標(biāo)參數(shù)如表2所示.

表2 Landsat熱紅外波段定標(biāo)參數(shù)與普朗克公式常數(shù)

其次，因?yàn)榇髿庵械乃推渌麣怏w會一定程度上吸收并向外發(fā)射熱紅外波段輻射，大氣頂端探測到的光譜輻射值有別于水體自身向外發(fā)出的輻射值. 因此，為了精確反演水體溫度，需剔除大氣對輻射產(chǎn)生的影響，即進(jìn)行大氣校正. 通常情況下，具有兩個(gè)熱紅外探測波段的衛(wèi)星，可以直接結(jié)合兩個(gè)波段的數(shù)據(jù)并利用劈窗法實(shí)現(xiàn)大氣校正[30]. 而對于Landsat 5與Landsat 7兩顆只具有單一熱紅外波段的衛(wèi)星而言，需額外輔以來自其他數(shù)據(jù)庫的大氣校正參數(shù)方可實(shí)現(xiàn)校正. 根據(jù)所使用大氣校正參數(shù)的不同，校正方法可以分為兩種：輻射傳遞模型法與Jimenez-Munoz和Sobrino(JM&S)單通道算法.

2.1.1 輻射傳遞模型法 輻射傳遞模型法(RTM)的主要思路是利用衛(wèi)星過空時(shí)的大氣校正參數(shù)來消除大氣對輻射傳遞過程的影響. 通常來說，衛(wèi)星傳感器在大氣頂端探測到的輻射主要由3部分組成：水體熱輻射、大氣上行輻射以及經(jīng)由地面反射的大氣下行輻射. 如果將大氣視為一個(gè)整體，大氣對水體輻射的影響主要取決于3個(gè)大氣校正參數(shù)：大氣透過率、大氣總下行輻射和大氣總上行輻射，公式為：

Lsensor=[ε·B+(1-ε)Ldown]τ+Lup

(2)

式中，ε為地表/水體比輻射率，不同類型的地表具有不同比輻射率，對于水體本文參考以往文章選取0.9885[31]；B為與水體同溫的黑體輻亮度，W/(sr·m2·nm)；Ldown為大氣下行輻射強(qiáng)度，W/(sr·m2·nm)；Lup為大氣上行輻射強(qiáng)度，W/(sr·m2·nm)；τ為大氣透過率.

由于Landsat自身不具有可供使用的大氣校正參數(shù)，因此額外選取了兩個(gè)專門用于校正Landsat系列衛(wèi)星熱紅外波段的大氣參數(shù)數(shù)據(jù)庫：(a)Atmcorr；(b) GEOS-5 FP-IT. Atmcorr是美國國家航天航空局以美國國家環(huán)境預(yù)報(bào)中心收集的衛(wèi)星和地面觀測大氣資料為原始數(shù)據(jù)，使用MODTRAN4.0大氣輻射傳輸模型模擬生成的大氣校正參數(shù)數(shù)據(jù)庫(https://atmcorr.gsfc.nasa.gov/). Atmcorr數(shù)據(jù)庫中包含用于校正Landsat 5 TM、Landsat 7 ETM+以及 Landsat 8 OIRS熱紅外波段的大氣透過率、上行與下行輻射強(qiáng)度以及一些大氣基本信息. 本文主要使用了Landsat 5與7熱紅外波段下的3個(gè)大氣校正參數(shù)，其具有1°×1°的空間分辨率，時(shí)間分辨率為6 h[32-33]. 第2種大氣參數(shù)數(shù)據(jù)庫GEOS-5 FP-IT是美國國家航天航空局基于ESMF(地球系統(tǒng)建?？蚣?和大量實(shí)測氣象數(shù)據(jù)構(gòu)建的大氣數(shù)據(jù)同化系統(tǒng)，其中包含了大量經(jīng)重分析后的大氣信息[34]. 本文提取了所需的3個(gè)大氣校正參數(shù)，其空間分辨率為0.625°×0.5°，時(shí)間分辨率為3 h.

兩種大氣校正參數(shù)的根本區(qū)別在于其所使用的大氣模型. Atmcorr使用的為專門模擬輻射傳輸過程的MODTRAN 4.0，屬于大氣校正中常用模型. 而GEOS-5 FP-IT所使用的是模擬整個(gè)地球系統(tǒng)，包括大氣、水體，真實(shí)物理情況的ESMF框架，其可用于氣象模擬、氣溶膠與微量氣體模擬分析等眾多領(lǐng)域. 另一方面，Atmcorr數(shù)據(jù)庫經(jīng)過大量的地面校正，其數(shù)據(jù)精度較為可信，而相比之下GEOS-5 FP-IT依然處于測試階段，因此需要后續(xù)更多的校準(zhǔn)和優(yōu)化.

兩種大氣校正參數(shù)經(jīng)時(shí)空三維插值以后，生成相應(yīng)時(shí)刻、相應(yīng)空間位置上的數(shù)據(jù)，并代入公式(2)中獲得與水體同溫黑體輻射值B. 然后進(jìn)一步使用變形后的普朗克公式得到水體溫度：

(3)

式中，RWT為河流水溫，K；K1、K2為普朗克公式等效常數(shù)，與傳感器熱紅外波段的有效波長相關(guān)，Landsat 5 TM和Landsat 7 ETM+的K1和K2如表2所示.

2.1.2 JM&S單通道算法 JM&S單通道算法由Jiménez-Muo和Sobrino于2003年提出[35]. 其簡化了輻射傳遞模型法中傳統(tǒng)的三參數(shù)公式. 使用大氣水汽含量與3個(gè)大氣校正參數(shù)之間的經(jīng)驗(yàn)關(guān)系替代三者，從而得到只需大氣水汽含量的單通道算法，其主要公式為：

(4)

式中，RWT為反演所得河流水溫，K；γ和δ為由普朗克公式?jīng)Q定的中間變量；ε為比輻射率，與輻射傳遞模型法中一致選取0.9885作為水體比輻射率；Lsensor為Landsat衛(wèi)星傳感器在大氣頂端接收到的熱紅外光譜輻射亮度，由公式(1)所得，W/(sr·m2·nm)；ψ1、ψ2、ψ3為大氣函數(shù)，公式為：

(5)

其中，ψ1、ψ2、ψ3為3個(gè)大氣函數(shù)取代了輻射傳遞模型中3個(gè)大氣校正參數(shù)：大氣透過率τ，大氣上行輻射Lup以及大氣下行輻射Ldown. JM&S單通道算法根據(jù)二項(xiàng)式擬合，將公式(5)中定義的大氣函數(shù)與大氣水汽含量進(jìn)行近似計(jì)算，這種近似可由矩陣表示為：

(6)

式中，w為大氣水汽含量，g/cm2；cij(i，j=1、2、3)為大氣函數(shù)的參數(shù)，與獲取大氣水汽含量的觀測手段和Landsat熱紅外波段等效波長相關(guān)[36]. 本文中使用的cij見表3.

表3 Landsat 5與7熱紅外波段大氣函數(shù)cij參數(shù)

除大氣函數(shù)之外，JM&S單通道算法還將輻射傳遞模型法中所需的普朗克公式估算成為公式(4)中兩個(gè)中間變量γ和δ：

(7)

δ=-γ·Lsensor+Tsensor

(8)

(9)

式中，λ為熱紅外波段的等效波長，μm；c1和c2為普朗克輻射公式中的常數(shù)，c1=1.19104×108(W·μm4)/(m2·sr)，c2=14387.7 μm·K；Tsensor為傳感器探測到的亮溫，K.

對于JM&S單通道算法所需要的大氣水汽含量，本研究選取Uwyoming大學(xué)收集處理的中國高空氣象站無線電探空儀探測數(shù)據(jù)(http://weather.uwyo.edu/upperair/sounding.html). 由于高空氣象站分布較為稀疏(全中國范圍內(nèi)僅有89個(gè)測站)，大部分研究區(qū)域內(nèi)并沒有實(shí)測數(shù)據(jù). 但大氣水汽含量在小尺度范圍內(nèi)變化較小，因此本文選取距研究區(qū)域(b)和(c)基點(diǎn)最近的測站：宜賓站、重慶站的數(shù)據(jù). 而對于研究區(qū)域(a)向家壩，由于其與最近測站西昌站海拔相差過多，因此未使用JM&S法反演研究區(qū)域(a)內(nèi)水溫. 至此，已對Landsat熱紅外數(shù)據(jù)進(jìn)行了輻射定標(biāo)和大氣校正兩步處理，如圖2所示.

2.1.3 Landsat反演水溫處理 無論是經(jīng)過輻射傳遞模型法還是JM&S單通道算法校正獲得的Landsat水溫影像都還存在以下兩個(gè)問題：(1)河流水溫尚與地面溫度摻雜在一起；(2)選取何種數(shù)據(jù)采樣窗口才能精確反映研究區(qū)域內(nèi)的真實(shí)水溫. 針對上述問題，本文采用以往研究中提出的TSM-S2法對Landsat水溫影像進(jìn)行進(jìn)一步處理[11,19](圖2). 具體步驟包括：水體提取、針對Landsat熱紅外波段水體-陸地邊界效應(yīng)的水體修正以及用于將柵格數(shù)據(jù)轉(zhuǎn)化為點(diǎn)數(shù)據(jù)的矩形選取. 同時(shí)，對使用輻射傳遞模型法和Atmcorr大氣校正參數(shù)反演所得的Landsat 7水溫應(yīng)用了TSM-S2法中基于實(shí)測數(shù)據(jù)與反演數(shù)據(jù)的邏輯回歸經(jīng)驗(yàn)公式，以修正其在高溫以及低溫區(qū)域出現(xiàn)的系統(tǒng)誤差. 由于本小節(jié)內(nèi)使用的4種數(shù)據(jù)處理方法和其對反演水溫精度帶來的改善,這部分內(nèi)容不是本文核心所在，詳細(xì)討論見參考文獻(xiàn)[19]. 處理后的Landsat水溫?cái)?shù)據(jù)按照數(shù)據(jù)來源與處理方法不同共分為7種，如表4所示.

表4 Landsat反演所得7種類型水溫?cái)?shù)據(jù)

2.2 處理MODIS溫度數(shù)據(jù)

本文主要對MODIS溫度數(shù)據(jù)進(jìn)行如下兩步處理：(1)提取每個(gè)研究區(qū)域內(nèi)水體占比最多的5個(gè)柵格；(2)利用質(zhì)量管控圖層剔除質(zhì)量較差的數(shù)據(jù)，如圖2所示. 具體處理原因和方法為：由于MOD11A1空間分辨率僅為1 km，在寸灘以及向家壩壩前部分河寬小于1 km的區(qū)域內(nèi)并不具有純水柵格，即包含水體的柵格中參雜有一定的陸地. 因此本文首先使用Landsat生成的30 m水體掩膜代入至MODIS 1 km柵格中，確定MODIS數(shù)據(jù)中單個(gè)柵格的水體占比，并分別選出3個(gè)研究區(qū)域內(nèi)水體占比最高的5個(gè)柵格，提取研究時(shí)段內(nèi)這些點(diǎn)的MODIS水溫. 第2步則是通過MOD11A1數(shù)據(jù)中的質(zhì)量管控圖層，剔除掉質(zhì)量較差(QC_day > 1)的數(shù)據(jù)，獲得有效的MODIS水溫?cái)?shù)據(jù).

3 結(jié)果分析

經(jīng)數(shù)據(jù)整合，將表4中同一顆衛(wèi)星下不同種類數(shù)據(jù)歸為一組，2004-2016年期間，寸灘共有27組Landsat 5有效水溫?cái)?shù)據(jù)，68組Landsat 7數(shù)據(jù)和807幅MODIS數(shù)據(jù)，三峽壩前共有37組Landsat 5水溫?cái)?shù)據(jù)、77組Landsat 7數(shù)據(jù)和1370幅MODIS數(shù)據(jù). 而對于向家壩壩前區(qū)域，2015年全年共有6組Landsat 7水溫?cái)?shù)據(jù)和64幅MODIS數(shù)據(jù). 從可用數(shù)據(jù)量的角度，MODIS水溫?cái)?shù)據(jù)遠(yuǎn)遠(yuǎn)多于Landsat系列衛(wèi)星. 同時(shí)，以往研究結(jié)果[11]表明3個(gè)研究區(qū)域之間的Landsat反演水溫誤差差異較小，因此本文將3個(gè)研究區(qū)域視為一個(gè)整體，在不考慮衛(wèi)星過空時(shí)間和實(shí)測時(shí)間差的情況下，分析了表4中7類Landsat溫度數(shù)據(jù)的絕對誤差、均方根誤差(RMSE)和偏差(Bias)，并加入了3個(gè)研究區(qū)域內(nèi)的MODIS水溫誤差，如表5所示.

表5 MODIS水溫產(chǎn)品和Landsat反演水溫的絕對誤差、RMSE和偏差

3.1 MODIS水溫產(chǎn)品與Landsat反演水溫比較

總體而言，MODIS水溫產(chǎn)品相較于Landsat反演水溫誤差偏大. 尤其是在寸灘以及向家壩壩前兩個(gè)研究區(qū)域，MODIS水溫RMSE分別高達(dá)6.24和4.86℃. 經(jīng)初步分析，可能原因?yàn)閮蓚€(gè)區(qū)域內(nèi)河道寬度較窄，水體溫度受陸地影響嚴(yán)重. 因此針對3個(gè)不同的研究區(qū)域，分別選取研究區(qū)域內(nèi)水體占比最高的5個(gè)柵格及其反演水溫誤差，匯總于圖3中.

圖3 MODIS數(shù)據(jù)水體占比與均方根誤差的關(guān)系

3.2 JM&S單通道算法與輻射傳遞模型法比較分析

在衡量過遙感數(shù)據(jù)差異對反演水溫的影響之后，本文繼續(xù)比較分析了不同反演方法的適用性. 如表5所示，在3個(gè)研究區(qū)域內(nèi)，利用輻射傳遞模型法反演的Landsat 5、Landsat 7水溫RMSE處于1.00～1.62℃之間；然而使用JM&S單通道算法獲取的同研究區(qū)域、時(shí)段內(nèi)的水溫?cái)?shù)據(jù)RMSE卻處于2.47～3.28℃之間，遠(yuǎn)遠(yuǎn)低于以往研究中JM&S法反演精度[36]，因此對JM&S法在長江上游流域可能存在的問題進(jìn)行了研究，分別生成了Landsat 5與Landsat 7反演水溫偏差與大氣水汽含量的關(guān)系(圖4).

圖4 JM&S反演水溫偏差與大氣水汽含量的關(guān)系

分析圖4，可以清楚地觀察到偏差與大氣水汽含量間存在一定的關(guān)系：當(dāng)大氣水汽含量w<4 g/cm2時(shí)，大部分水溫偏差都處于[-4℃,+4℃]的灰色區(qū)間內(nèi)，且大多數(shù)聚集于0℃線兩邊；而當(dāng)w>4 g/cm2時(shí)，隨著w的增長大部分點(diǎn)的水溫偏差游離于灰色區(qū)域以外，且0℃中心線附近數(shù)據(jù)點(diǎn)較稀疏. 當(dāng)剔除掉大氣水汽含量過大(w>4 g/cm2)反演所得的水溫?cái)?shù)據(jù)后，Landsat 5反演水溫RMSE從3.28℃下降到1.96℃，Landsat 7反演水溫RMSE從2.47℃下降到1.76℃. 在Jiménez-Muo對該方法的分析中曾指出，當(dāng)大氣水汽含量小于0.5 g/cm2或大于3 g/cm2時(shí)，利用經(jīng)驗(yàn)公式(5)、(6)模擬生成的大氣校正參數(shù)會與實(shí)際大氣參數(shù)間存在較大的差異，導(dǎo)致反演水溫誤差過大[32]. 因此，本研究收集了長江上游流域內(nèi)4個(gè)高空氣象站探測的多年大氣水汽含量數(shù)據(jù)，數(shù)據(jù)表明長江上游流域一年365 d中平均有186.5 d大氣水汽含量大于3 g/cm2，且主要分布在6-9月. 所以，當(dāng)剔除這些因大氣水汽含量過高而出現(xiàn)問題的數(shù)據(jù)，并忽略季節(jié)性云層對Landsat衛(wèi)星有效數(shù)據(jù)的影響，使用JM&S單通道算法得到的溫度數(shù)據(jù)量約為輻射傳遞模型法的1/2. 除此之外，由于JM&S單通道算法缺少夏季高溫?cái)?shù)據(jù)，結(jié)合以往文獻(xiàn)的結(jié)論[16]，其反演所得水溫相較輻射傳遞模型法會更難精確模擬高溫區(qū)域日水溫. 另一方面，在長江上游流域內(nèi)高空氣象站分布較為稀疏，能精確獲取大氣水汽含量的區(qū)域相對較少. 綜合上述結(jié)論，目前利用JM&S單通道算法反演長江上游流域內(nèi)河道水溫存在兩方面的局限性：(1)研究區(qū)域受高空氣象站地理位置限制；(2)反演水溫精度受大氣水汽含量影響. 相較于JM&S單通道算法，輻射傳遞模型法的適用區(qū)域更廣，反演精度更高.

3.3 不同大氣參數(shù)下輻射傳遞模型法反演水溫的精度分析

在比較過不同大氣校正方法的反演精度之后，本文繼續(xù)對采用輻射傳遞模型法并利用不同大氣校正參數(shù)和修正方法所得的5種數(shù)據(jù)進(jìn)行了精度分析. 其中包括(1)L7-RTM-GEOS-5;(2)L7-RTM-Atmcorr;(3)L7-RTM-Atmcorr-Log;(4)L5-RTM-GEOS-5;(5)L5-RTM- Atmcorr. 由于5種數(shù)據(jù)均使用了輻射傳遞模型法RTM，因此后續(xù)章節(jié)及圖中省略所有RTM.

如表5所示，對于Landsat 5熱紅外波段數(shù)據(jù)，使用Atmcorr和GEOS-5大氣校正參數(shù)反演所得水溫RMSE分別為1.13和1.33℃. 而在相同研究區(qū)域內(nèi)，L7-Atmcorr反演水溫RMSE為1.62℃，大于L7-GEOS-5的1.52℃. 因此，通過比較均方根誤差無法得到兩種大氣校正參數(shù)之間的區(qū)別. 因此針對上述5種數(shù)據(jù)，各自生成了基于數(shù)據(jù)偏差的小提琴圖和箱形圖，如圖5所示. 從圖5中可以觀察到，使用Atmcorr大氣校正參數(shù)反演所得水溫相較GEOS-5具有更矮的箱形與更寬的密度峰值，但其偏差密度峰值和偏差的中位數(shù)距離0℃更遠(yuǎn)，例如L7-Atmcorr(2)的密度峰值出現(xiàn)在大約-0.3℃，L7-Atmcorr的密度峰值(4)位于0.4～0.5℃. 而對于GEOS-5數(shù)據(jù)，其箱形更長且小提琴圖圖寬更狹長，并具有更大和更小的數(shù)據(jù)偏差極值. 以上兩種現(xiàn)象說明，使用GEOS-5反演所得的水溫?cái)?shù)據(jù)相較Atmcorr有更高的概率存在數(shù)據(jù)異常值，從而影響整體精度. 而Atmcorr反演水溫，相較于GEOS-5可能存在更嚴(yán)重的系統(tǒng)性偏差.

圖5 5種反演水溫?cái)?shù)據(jù)偏差的箱形和小提琴復(fù)合圖

對于Atmcorr存在的系統(tǒng)性偏差，往往可以進(jìn)一步使用經(jīng)驗(yàn)公式進(jìn)行修正. 本文對L7-Atmcorr數(shù)據(jù)進(jìn)行了邏輯回歸修正并得到了L7-Atmcorr-Log數(shù)據(jù). 如圖5(3)所示，相較于L7-Atmcorr，修正后的數(shù)據(jù)具有更小的系統(tǒng)性偏差，即數(shù)據(jù)中位數(shù)以及密度峰值更趨近于0℃，同時(shí)具有誤差極值也回歸到正常區(qū)間. 而針對GEOS-5中存在的數(shù)據(jù)異常值問題，本文分別生成了Landsat 7和Landsat 5數(shù)據(jù)下的GEOS-5、Atmcorr反演水溫與實(shí)測水溫散點(diǎn)圖(圖6).

圖6 Landsat 5與7熱紅外數(shù)據(jù)的GEOS-5、Atmcorr反演水溫與實(shí)測水溫散點(diǎn)圖

如圖6所示，盡管Atmcorr與GEOS-5之間存在大量重合或相似數(shù)據(jù)組(紅色空心圈包含藍(lán)色×的數(shù)據(jù)組), 仍有4組GEOS-5數(shù)據(jù)出現(xiàn)了較大偏差. Landsat 7數(shù)據(jù)中1組，Landsat 5數(shù)據(jù)中共有3組，且4組數(shù)據(jù)所在溫度區(qū)間相對分散，說明數(shù)據(jù)異常與溫度之間不存在一定的關(guān)系，無法通過經(jīng)驗(yàn)公式修正. 因此需進(jìn)一步分析其成因，并盡可能地在數(shù)據(jù)處理過程中找到方法進(jìn)行剔除.

3.4 輻射傳遞模型法中各參數(shù)敏感性分析

首先，通過解析輻射傳遞模型法中的核心公式，公式(2)，可知由于水體在熱紅外波段的比輻射率較高，大部分大氣下行輻射會被水體吸收，因此在已知大氣下行輻射與大氣上行輻射為同量級數(shù)據(jù)的情況下，影響反演水溫大小以及精度的主要因子為Lup/τ，即大氣上行輻射與大氣透過率的比值. 且當(dāng)Lup與τ的比值越大時(shí)，反演水溫會對大氣校正參數(shù)越敏感，即更容易產(chǎn)生較大水溫誤差，因此為了嘗試分析和解決GEOS-5反演水溫存在的異常值問題，生成了GEOS-5反演水溫偏差與大氣上行輻射，大氣透過率以及兩者比值的相關(guān)關(guān)系(圖7).

圖7中3張圖的形狀特征相似，GEOS-5反演水溫偏差數(shù)據(jù)上下包絡(luò)線呈喇叭口型：當(dāng)大氣上行輻射相對較小或大氣透過率校高時(shí)，GEOS-5反演水溫偏差聚集在[-4℃, 4℃]之間，但隨著大氣上行輻射的變大和大氣透過率的變小，水溫?cái)?shù)據(jù)的偏差相對0℃刻度線越來越分散，且有部分?jǐn)?shù)據(jù)超出偏差正常區(qū)間，這部分?jǐn)?shù)據(jù)即為3.3節(jié)中所提到的異常數(shù)據(jù)點(diǎn). 根據(jù)上述特性，本文通過對大氣校正參數(shù)設(shè)定相應(yīng)閾值來減少異常數(shù)據(jù)點(diǎn)的出現(xiàn)概率. 設(shè)定Lup<4.5，τ>0.4，Lup/τ<11.5 3條閾值分界線進(jìn)行L7-GEOS-5數(shù)據(jù)劃分，3種閾值分別剔除了28、23、25組數(shù)據(jù)，皆約占數(shù)據(jù)總量的1/6，整體水溫?cái)?shù)據(jù)誤差RMSE由1.52℃降低到了1.20℃，低于L7-Atmcorr反演水溫?cái)?shù)據(jù)的均方根誤差(1.62℃)，但高于L7-Atmcorr-log修正數(shù)據(jù)的均方根誤差(1.00℃). 但相較于L7-Atmcorr-log修正數(shù)據(jù)，這種處理方法的優(yōu)勢是其不需要實(shí)測水溫資料，可以應(yīng)用于無測站地區(qū). 綜上所述，對于GEOS-5反演水溫而言，由于其具有與大氣校正參數(shù)相關(guān)的數(shù)據(jù)異常值，建議分析并設(shè)定校正參數(shù)閾值以提高反演精度. 而對于Atmcorr反演水溫所具有的系統(tǒng)性誤差，在有實(shí)測水溫資料地區(qū)可以使用邏輯回歸等方法進(jìn)行校正從而提高水溫精度.

圖7 Landsat 5與7 GEOS-5反演水溫與大氣參數(shù)間的關(guān)系

4 結(jié)論

本文利用不同數(shù)據(jù)處理流程獲取了長江上游河段MODIS和Landsat水溫?cái)?shù)據(jù)，結(jié)果表明在長江上游河段各種方法和數(shù)據(jù)都具有各自的優(yōu)缺點(diǎn)：

1)MODIS溫度產(chǎn)品的優(yōu)點(diǎn)在于其可以提供日水溫?cái)?shù)據(jù)，但由于自身空間分辨率無法滿足長江上游流域河道要求，導(dǎo)致數(shù)據(jù)精度過低，因此不建議使用.

2)對于JM&S單通道算法，其主要優(yōu)點(diǎn)為只需大氣水汽含量這一個(gè)大氣參數(shù). 但是由于長江上游流域在6-9月大部分日期大氣水汽含量過高，會嚴(yán)重影響反演水溫精度，因此同樣不建議使用.

3)對于輻射傳遞模型法RTM以及大氣參數(shù)數(shù)據(jù)庫GEOS-5 FP-IT反演所得數(shù)據(jù)，其主要缺點(diǎn)是存在數(shù)據(jù)異常值，本文提出通過設(shè)置大氣校正參數(shù)的閾值(Lup<4.5，τ>0.4，Lup/τ<11.5)剔除部分?jǐn)?shù)據(jù)來減小數(shù)據(jù)異常點(diǎn)帶來的誤差. 但閾值的普適性有待未來在更多研究區(qū)域內(nèi)進(jìn)一步的研究.

4)使用輻射傳遞模型法RTM以及大氣參數(shù)數(shù)據(jù)庫Atmcorr為眾多方法中最為穩(wěn)定的流程，其數(shù)據(jù)精度略位于方法(3)所獲得的原始水溫?cái)?shù)據(jù)精度方法(3)修正后水溫?cái)?shù)據(jù)精度之間. 其缺點(diǎn)在于存在一定的系統(tǒng)性偏差，如果研究區(qū)域內(nèi)有適量的實(shí)測水溫?cái)?shù)據(jù)，可以通過邏輯回歸等經(jīng)驗(yàn)公式進(jìn)行修正，修正后的數(shù)據(jù)精度為所有方法中最高.

綜上所述，在反演長江上游流域水溫且研究區(qū)域內(nèi)有適量實(shí)測水溫?cái)?shù)據(jù)時(shí)，優(yōu)先選用RTM-Atmorr-Log數(shù)據(jù). 當(dāng)不考慮有效數(shù)據(jù)量且研究區(qū)域內(nèi)無實(shí)測水溫時(shí)，建議使用RTM-GEOS-5數(shù)據(jù)并利用閾值法剔除部分?jǐn)?shù)據(jù)點(diǎn). 對于數(shù)據(jù)量有限且區(qū)域內(nèi)無實(shí)測水溫?cái)?shù)據(jù)時(shí)，應(yīng)選用RTM-Atmcorr原始水溫?cái)?shù)據(jù).