湖泊科学   2022, Vol. 34 Issue (1): 307-319.  DOI: 10.18307/2022.0125
0

研究论文——流域水文与水资源安全

引用本文 [复制中英文]

石希, 夏军强, 孙健, 基于热红外遥感影像的河流水温反演方法比较——以长江上游流域为例. 湖泊科学, 2022, 34(1): 307-319. DOI: 10.18307/2022.0125
[复制中文]
Shi Xi, Xia Junqiang, Sun Jian. Comparison of methods to derive river water temperature using thermal infrared imagery: A case study of the upper Yangtze River catchment. Journal of Lake Sciences, 2022, 34(1): 307-319. DOI: 10.18307/2022.0125
[复制英文]

基金项目

国家自然科学基金项目(51725902,U2040215)资助

通信作者

夏军强, E-mail: xiajq@whu.edu.cn

文章历史

2021-02-25 收稿
2021-05-17 收修改稿

码上扫一扫

基于热红外遥感影像的河流水温反演方法比较——以长江上游流域为例
石希1 , 夏军强1 , 孙健2     
(1: 武汉大学水资源与水电工程科学国家重点实验室, 武汉 430072)
(2: 清华大学水沙科学与水利水电工程国家重点实验室, 北京 100084)
摘要:水温是影响河流生态环境的重要因素.传统水温观测手段受限于测站的地理位置,数据分布稀疏.近20年来随着热红外卫星传感器在空间分辨率上的革新,热红外遥感影像逐步被应用于反演河流水温,但目前相关研究聚焦于利用单颗卫星数据并使用单一水温反演方法而忽略与其他卫星以及处理方法做对比.本文以长江上游流域为例,分别评估了水文反演流程中3个步骤的差异对水温精度的影响:遥感数据(MODIS和Landsat)、大气校正方法(辐射传递模型法与JM & S单通道算法)和大气校正参数(Atmcorr和GEOS-5 FP-IT).研究结果表明:(ⅰ)MODIS数据因为空间分辨率较低并不能满足反演长江上游流域水温的需求;(ⅱ)由于长江上游流域大气水汽含量较高,JM & S单通道算法反演水温误差过大;(ⅲ)当流域拥有实测水温资料时,最佳反演流程是使用辐射传递模型法与Atmcorr大气校正参数反演Landsat水温并进一步使用逻辑回归模型校正数据,最终所得反演水温RMSE为1.0~1.1℃;(ⅳ)对于无测站流域,较适合的流程为使用辐射传递模型法和GEOS-5 FP-IT大气校正参数并设定大气校正参数阈值,最终反演水温RMSE为1.2℃.同时研究还表明大气透过率以及大气上行辐射强度是影响反演水温精度的主要因素.因此本文建议在对不同流域进行水温反演时应尽可能对上述两种参数进行敏感性分析.
关键词河流水温    反演方法    遥感影像    Landsat    MODIS    长江上游    
Comparison of methods to derive river water temperature using thermal infrared imagery: A case study of the upper Yangtze River catchment
Shi Xi1 , Xia Junqiang1 , Sun Jian2     
(1: State Key Laboratory of Water Resources and Hydropower Engineering Science, Wuhan University, Wuhan 430072, P. R. China)
(2: State Key Laboratory of Hydroscience and Engineering, Tsinghua University, Beijing 100084, P. R. China)
Abstract: Water temperature is an essential factor that affects the ecology and environment of the river. Conventional river water temperature (RWT) measuring methods are limited by the geographic location of the hydrological station; thus, RWT data is under a sparse spatial resolution. With the improvement of the spatial resolution of thermal infrared sensors mounted on the satellites within the past twenty years, thermal infrared imagery has gradually been used to retrieve RWT. Currently, most research focuses on using sole sensor and specified processing methods without comparing with other sensors and retrieval approaches. Taking the upper Yangtze River catchment as a case study, this paper evaluates the impact of the following three steps within the entire process on the accuracy of RWT: thermal infrared data (MODIS and Landsat); atmospheric correction methods (the radiation transfer model and JM & S single-channel algorithm); atmospheric correction parameters (Atmcorr and GEOS-5 FP-IT). The results indicated that (ⅰ) MODIS data cannot meet the requirements of the study area due to its low spatial resolution; (ⅱ) As the atmospheric water vapor content is relatively high in the study catchment, the JM & S single-channel algorithm performs poorly on the retrieval; (ⅲ) When there is measured RWT within the study catchment, the proper retrieval process will be using the RTM combined with Atmcorr atmospheric correction parameter and further corrected by using the logistic regression while the final RMSE is about 1.0-1.1℃; (ⅳ) For the basin without measured RWT, the more suitable way is to use the GEOS-5 FP-IT instead and set an atmospheric correction parameter threshold. Its RMSE is approximately 1.2℃. Meanwhile, the results also present that atmospheric transmittance and upward radiation are the key parameters affecting the accuracy of RWT. Therefore, we highly recommended that a similar sensitivity analysis should be applied to other catchments for further studies.
Keywords: River water temperature    retrieval method    remote sensing images    Landsat    MODIS    upper Yangtze River    

水体温度是河流生态系统中的重要参数,其调控着大量物理、化学和生物过程,例如:水生生物的迁徙和营养物质的循环等[1-3]. 长期观测河流水温对于监管河流热动态和维护河流生态系统的健康起着至关重要的作用[4]. 同时,分析长时间下水温变化趋势也有利于探知全球变暖和人类活动对河流水温的影响[5]. 目前,传统水温观测手段主要为利用温度测量仪器实地探测测站附近的河流水温[6]. 然而由于水文测站分布较稀疏,大部分河段缺乏实测水温资料. 对于这些无测站地区,获取满足研究需求的水温数据往往需要花费大量人力物力建设水文观测站网并周期性地开展现场工作,常理上难以实现[7].

相较于传统观测手段,利用热红外遥感技术可以获取空间上连续的水温数据[8]. 近20年来,随着卫星热红外传感器在空间分辨率和辐射探测精度等技术方面上的不断进步,热红外遥感影像不再单单用于反演海洋、湖泊等大型水体水温[9-10],而逐步被应用于河流、河道型水库等中小尺度水体[11-12]. 但是受限于传感器硬件性能和卫星轨道设计等问题,卫星遥感数据普遍只能满足高时间分辨率或高空间分辨率两者中的一项[13]. 因此对于河流和河道型水库而言,使用契合研究需求的卫星至关重要[14-15].

在探测河流及河道型水库水温方面,Landsat系列卫星搭载的Landsat 7 ETM+、Landsat 8 OIRS和已经失效的Landsat 5 TM传感器凭借着其60、100和120 m的热红外波段中等空间分辨率和近20年稳定的历史数据库这两点优势使用最为广泛[16]. 其反演所得水温已被应用于模拟河流日水温[17]、评估河流热污染状况[18]、分析建造水库对河流上下游水温影响[19]等各项研究中. 而对于另一个常用的卫星传感器MODIS而言,由于其热红外波段仅具有1000 m的低空间分辨率,目前只被应用于河道型水库、湖泊和长江中下游干流等河道较宽的水体[20-22]. 但是相较于Landsat系列卫星16 d的时间分辨率,MODIS的重访周期仅为1 d,且MODIS具有较高的辐射/温度探测精度NEΔT=0.05℃(在27℃时的等效温度探测精度,Landsat NEΔT=0.4℃)[23]. 因此对于不同的研究区域、研究尺度选择合适的卫星尤为关键. 同时,目前相关研究都仅聚焦于评估单颗传感器的反演精度,但研究结果表明单颗传感器往往具有一定的局限性[24]. 因此,利用不同卫星传感器获取相同河段的水温并横向对比其精度,可以更直观地展示传感器之间的差异,也有助于评估不同传感器在相应流域中的实用性.

除此之外,使用不同的大气校正方法和大气校正参数对反演水温的精度也有一定影响. 通常情况下,具有单一热红外波段的卫星传感器,例如Landsat 5 TM与Landsat 7 ETM+, 在大气顶端探测到的热红外辐射可以通过两种大气校正方法:辐射传递模型法和JM & S (Jimenez-Munoz & Sobrino)单通道算法转化为水体温度[25]. 以往研究表明,由于两种方法对大气水汽含量的敏感性有所差异,使用两者在具有不同气候条件的研究区域内反演所得的水温精度各有高低[26-27]. 因此在大范围反演任意流域水温之前,应尽可能地使用两种方法得到的水温与实测水温比对并分析两者在该研究区域内的适用性. 另一方面,两种反演方法均需要输入额外的大气校正参数,这些大气校正参数通常由不同的大气模型辅以实测气象资料模拟所得[25]. 查阅相关资料发现当大气校正参数处于异常范围时会严重影响校正精度[28]. 因此,在评估大气校正方法精度的同时,进行大气校正参数敏感性分析也尤为重要.

本文选取长江上游流域内具有不同河道宽度以及气候环境条件的3个研究区域:向家坝坝前、三峡坝前和三峡库位寸滩,分析使用不同的遥感数据、大气校正方法以及大气校正参数对水温精度的影响:(1)遥感数据:MODIS和Landsat;(2)大气校正方法:辐射传递模型法和JM & S单通道算法;(3)大气校正参数:Atmcorr(atmosphere correction/大气校正参数数据库)和GEOS-5 FP-IT(goddard earth observing system-version 5-forward processing for instrument teams/戈达德地球观测系统模型-版本5-仪器团队版本),并对大气校正参数进行了敏感性分析. 此外,本文通过对比各类反演水温精度以及适用条件,为今后利用热红外遥感数据反演长江上游流域其他河段水温提供了可供参考的方案.

1 研究区域与数据 1.1 研究区域

本研究聚焦于长江上游流域,包括金沙江及长江至三峡段,共选取3个研究区域:(a)向家坝坝前、(b)三峡坝前和(c)寸滩(图 1). 其中向家坝水库位于云南省和四川省交界处,坐落于金沙江下游河段,距宜宾市33 km. 选取坝上测站作为研究区域内基点,坝前河道宽度在800~1000 m之间. 其属中亚热带性季风气候,气候湿润,常年气温处于0~ 30℃. 三峡大坝坐落于湖北省宜昌市,是长江上游和中游的分界线. 本文选取大坝下游测站黄陵庙站作为研究区域内基点,坝前河宽为2~3 km,黄陵庙站附近河宽在550~750 m. 其属亚热带季风气候,秋季多雨,常年气温徘徊于5~35℃,夏季大气水汽含量约为4 g/cm2. 寸滩位于重庆市,属于三峡水库变动回水区. 选取寸滩水文站作为研究区域内基点,寸滩站距离上游嘉陵江与长江交汇处约7 km,距上游向家坝约406 km,距下游三峡大坝约570 km. 河道宽度处于600~850 m. 其气候类型、气温与大气水汽含量与三峡坝前相似.

图 1 研究区域示意图 Fig.1 Sketch map of the study area
1.2 实测水温数据

本研究共选用3个水文测站的实测水温数据:向家坝坝上测站、黄陵庙站和寸滩站,分别对应上述3个研究区域(图 1). 其中向家坝坝上测站实测水温数据可用年份为2015年,寸滩站与黄陵庙站实测水温数据范围为2004-2016年. 3个测站的水温皆测量于每日上午8:00,所有实测数据皆来源于水文年鉴.

1.3 卫星遥感数据

本文选取Landsat 5 TM与Landsat 7 ETM+热红外遥感数据以及MODIS温度数据产品MOD11. 而对于Landsat 8卫星,由于其热红外传感器收到杂散光影响精度较低,因此并未使用. 3种卫星遥感数据的基本信息如下.

1.3.1 Landsat热红外数据

Landsat 5和Landsat 7热红外数据均具有30 m空间分辨率(原始分辨率为120 m和60 m)和16 d时间分辨率,属于Level-1级遥感影像. 对于Landsat热红外数据,结合图 1中各研究区域所处景号以及实测水温可用年份,共选定5组数据,如表 1所示. 其中由于Landsat 5从2011年起不再提供有效数据,因此所有Landsat 5 TM热红外数据的截止年份为2011年. 此外,不同于Landsat 5 TM热红外波段B6,Landsat 7 ETM+的热红外波段B6具有两个通道:High gain和Low gain(高/低增益模式). 前者具有更高的探测精度,后者具有更大的温度/辐射探测范围[29]. 查阅Landsat 7用户使用手册后发现,由于正常河流水温不会导致传感器饱和,因此选取具有更高探测精度的High gain通道. 针对遥感数据质量,选择高质量的T1级别数据,同时对于因扫描线校正器失效的Landsat 7热红外影像,本文选择具有空白条带的SLC-off数据以避免修正带来的误差. 所有满足条件的数据皆从美国地质调查局地球探索数据库(https://earthexplorer.usgs.gov/)中免费下载. 下载的Landsat遥感影像均为Geotiff文件,并采用了WGS84椭球坐标基准以及UTM投影,数据获取时间约为每日上午11:00.

表 1 所用Landsat热红外数据基本信息 Tab. 1 Basic information of Landsat TIR data
1.3.2 MODIS温度数据

MODIS数据库中具有两个地表温度数据集:MOD11和MOD28,分别为陆地表面温度与海洋和大型水体温度[26]. 经查阅后发现,MOD28数据集并不涵盖长江上游流域. 因此本次研究只选用MOD11温度数据集. MOD11A1数据为MOD11数据集中的原始数据,其具有1 d的时间分辨率和1 km的空间分辨率(精确的网格大小为0.928 km),属于Level-3级别遥感影像. 其主要采用劈窗法对MODIS的31、32热红外波段原始数据进行了大气校正,并结合土地类型和云层覆盖数据生成了晴天状态下的地表温度. 本文选用覆盖3个研究区的MODIS温度产品,时间跨度选取2004-2016年,所有满足条件的数据皆下载于LAADS DAAC(大气和一级遥感数据存储和查询系统分布式存档中心)数据库(https://ladsweb.modaps.eosdis.nasa.gov/search/). 下载的MODIS温度数据格式为HDF-EOS,数据获取时间约为每日上午10:00.

2 研究方法

根据所用数据不同,研究方法分为两部分(图 2):(1)Landsat热红外数据处理;(2)MODIS温度数据处理.

图 2 研究方法流程图 Fig.2 Flow chart of the methodology implemented in this study
2.1 Landsat热红外数据处理

Landsat所提供的原始热红外数据为0~255的亮度值,并不能直接反映辐射强度. 所以,首先需将原始数据通过辐射定标转化为传感器在大气顶端所探测到的光谱辐射亮度,公式为:

$ {L_{{\rm{sensor}}}} = {\rm{gain}} \cdot {\rm{DN + bias}} $ (1)

式中,Lsensor为Landsat卫星传感器在大气顶端接收到的热红外光谱辐射亮度/辐射值,W/(sr ·m2 ·nm);DN为原始热红外数据中的亮度值;gain和bias为定标参数,Landsat 5和Landsat 7热红外波段的定标参数如表 2所示.

表 2 Landsat热红外波段定标参数与普朗克公式常数 Tab. 2 Radiometric calibration parameters and Plank's equation parameters of Landsat TIR band

其次,因为大气中的水汽和其他气体会一定程度上吸收并向外发射热红外波段辐射,大气顶端探测到的光谱辐射值有别于水体自身向外发出的辐射值. 因此,为了精确反演水体温度,需剔除大气对辐射产生的影响,即进行大气校正. 通常情况下,具有两个热红外探测波段的卫星,可以直接结合两个波段的数据并利用劈窗法实现大气校正[30]. 而对于Landsat 5与Landsat 7两颗只具有单一热红外波段的卫星而言,需额外辅以来自其他数据库的大气校正参数方可实现校正. 根据所使用大气校正参数的不同,校正方法可以分为两种:辐射传递模型法与Jimenez-Munoz和Sobrino(JM & S)单通道算法.

2.1.1 辐射传递模型法

辐射传递模型法(RTM)的主要思路是利用卫星过空时的大气校正参数来消除大气对辐射传递过程的影响. 通常来说,卫星传感器在大气顶端探测到的辐射主要由3部分组成:水体热辐射、大气上行辐射以及经由地面反射的大气下行辐射. 如果将大气视为一个整体,大气对水体辐射的影响主要取决于3个大气校正参数:大气透过率、大气总下行辐射和大气总上行辐射,公式为:

$ {L_{{\rm{sensor}}}} = \left[ {\varepsilon \cdot B + \left( {1 - \varepsilon } \right){L_{{\rm{down}}}}} \right]\tau + {L_{{\rm{up}}}} $ (2)

式中,ε为地表/水体比辐射率,不同类型的地表具有不同比辐射率,对于水体本文参考以往文章选取0.9885[31]B为与水体同温的黑体辐亮度,W/(sr ·m2 ·nm);Ldown为大气下行辐射强度,W/(sr ·m2 ·nm);Lup为大气上行辐射强度,W/(sr ·m2 ·nm);τ为大气透过率.

由于Landsat自身不具有可供使用的大气校正参数,因此额外选取了两个专门用于校正Landsat系列卫星热红外波段的大气参数数据库:(a)Atmcorr;(b) GEOS-5 FP-IT. Atmcorr是美国国家航天航空局以美国国家环境预报中心收集的卫星和地面观测大气资料为原始数据,使用MODTRAN4.0大气辐射传输模型模拟生成的大气校正参数数据库(https://atmcorr.gsfc.nasa.gov/). Atmcorr数据库中包含用于校正Landsat 5 TM、Landsat 7 ETM+以及Landsat 8 OIRS热红外波段的大气透过率、上行与下行辐射强度以及一些大气基本信息. 本文主要使用了Landsat 5与7热红外波段下的3个大气校正参数,其具有1°×1°的空间分辨率,时间分辨率为6 h[32-33]. 第2种大气参数数据库GEOS-5 FP-IT是美国国家航天航空局基于ESMF(地球系统建模框架)和大量实测气象数据构建的大气数据同化系统,其中包含了大量经重分析后的大气信息[34]. 本文提取了所需的3个大气校正参数,其空间分辨率为0.625°×0.5°,时间分辨率为3 h.

两种大气校正参数的根本区别在于其所使用的大气模型. Atmcorr使用的为专门模拟辐射传输过程的MODTRAN 4.0,属于大气校正中常用模型. 而GEOS-5 FP-IT所使用的是模拟整个地球系统,包括大气、水体,真实物理情况的ESMF框架,其可用于气象模拟、气溶胶与微量气体模拟分析等众多领域. 另一方面,Atmcorr数据库经过大量的地面校正,其数据精度较为可信,而相比之下GEOS-5 FP-IT依然处于测试阶段,因此需要后续更多的校准和优化.

两种大气校正参数经时空三维插值以后,生成相应时刻、相应空间位置上的数据,并代入公式(2)中获得与水体同温黑体辐射值B. 然后进一步使用变形后的普朗克公式得到水体温度:

$ RWT = \frac{{{K_2}}}{{\ln \left( {\frac{{{K_1}}}{B} + 1} \right)}} $ (3)

式中,RWT为河流水温,K;K1K2为普朗克公式等效常数,与传感器热红外波段的有效波长相关,Landsat 5 TM和Landsat 7 ETM+的K1K2表 2所示.

2.1.2 JM & S单通道算法

JM & S单通道算法由Jiménez-Mu${\rm{\tilde n}}$o和Sobrino于2003年提出[35]. 其简化了辐射传递模型法中传统的三参数公式. 使用大气水汽含量与3个大气校正参数之间的经验关系替代三者,从而得到只需大气水汽含量的单通道算法,其主要公式为:

$ RWT = \gamma \left[ {\frac{1}{\varepsilon }\left( {{\psi _1} \cdot {L_{{\rm{sensor}}}} + {\psi _2}} \right) + {\psi _3}} \right] + \delta $ (4)

式中,RWT为反演所得河流水温,K;γδ为由普朗克公式决定的中间变量;ε为比辐射率,与辐射传递模型法中一致选取0.9885作为水体比辐射率;Lsensor为Landsat卫星传感器在大气顶端接收到的热红外光谱辐射亮度,由公式(1)所得,W/(sr·m2·nm);ψ1ψ2ψ3为大气函数,公式为:

$ {\psi _1} = \frac{1}{\tau };{\psi _2} = - {L_{{\rm{down}}}} - \frac{{{L_{{\rm{up}}}}}}{\tau };{\psi _3} = {L_{{\rm{down}}}} $ (5)

其中,ψ1ψ2ψ3为3个大气函数取代了辐射传递模型中3个大气校正参数:大气透过率τ,大气上行辐射Lup以及大气下行辐射Ldown. JM & S单通道算法根据二项式拟合,将公式(5)中定义的大气函数与大气水汽含量进行近似计算,这种近似可由矩阵表示为:

$ \left[ {\begin{array}{*{20}{c}} {{\psi _1}}\\ {{\psi _2}}\\ {{\psi _3}} \end{array}} \right] = \left[ {\begin{array}{*{20}{c}} {{c_{11}}}&{{c_{12}}}&{{c_{13}}}\\ {{c_{21}}}&{{c_{22}}}&{{c_{23}}}\\ {{c_{31}}}&{{c_{32}}}&{{c_{33}}} \end{array}} \right]\left[ {\begin{array}{*{20}{c}} {{w^2}}\\ w\\ 1 \end{array}} \right] $ (6)

式中,w为大气水汽含量,g/cm2cij(ij=1、2、3)为大气函数的参数,与获取大气水汽含量的观测手段和Landsat热红外波段等效波长相关[36]. 本文中使用的cij表 3.

表 3 Landsat 5与7热红外波段大气函数cij参数 Tab. 3 Values of cij for the Landsat 5 and 7 TIR bands

除大气函数之外,JM & S单通道算法还将辐射传递模型法中所需的普朗克公式估算成为公式(4)中两个中间变量γδ

$ \gamma = {\left\{ {\frac{{{c_2} \cdot {L_{{\rm{sensor}}}}}}{{T_{{\rm{sensor}}}^2}}\left[ {\frac{{{\lambda ^4} \cdot {L_{{\rm{sensor}}}}}}{{{c_1}}} + \frac{1}{\lambda }} \right]} \right\}^{ - 1}} $ (7)
$ \delta = - \gamma \cdot {L_{{\rm{sensor}}}} + {T_{{\rm{sensor}}}} $ (8)
$ {T_{{\rm{sensor}}}} = \frac{{{c_2}}}{{\lambda \ln \left( {\frac{{{c_1}}}{{{\lambda ^5} \cdot {L_{{\rm{sensor}}}}}} + 1} \right)}} $ (9)

式中,λ为热红外波段的等效波长,μm;c1c2为普朗克辐射公式中的常数,c1=1.19104×108(W ·μm4)/(m2 ·sr),c2=14387.7 μm ·K;Tsensor为传感器探测到的亮温,K.

对于JM & S单通道算法所需要的大气水汽含量,本研究选取Uwyoming大学收集处理的中国高空气象站无线电探空仪探测数据(http://weather.uwyo.edu/upperair/sounding.html). 由于高空气象站分布较为稀疏(全中国范围内仅有89个测站),大部分研究区域内并没有实测数据. 但大气水汽含量在小尺度范围内变化较小,因此本文选取距研究区域(b)和(c)基点最近的测站:宜宾站、重庆站的数据. 而对于研究区域(a)向家坝,由于其与最近测站西昌站海拔相差过多,因此未使用JM & S法反演研究区域(a)内水温. 至此,已对Landsat热红外数据进行了辐射定标和大气校正两步处理,如图 2所示.

2.1.3 Landsat反演水温处理

无论是经过辐射传递模型法还是JM & S单通道算法校正获得的Landsat水温影像都还存在以下两个问题:(1)河流水温尚与地面温度掺杂在一起;(2)选取何种数据采样窗口才能精确反映研究区域内的真实水温. 针对上述问题,本文采用以往研究中提出的TSM-S2法对Landsat水温影像进行进一步处理[11, 19](图 2). 具体步骤包括:水体提取、针对Landsat热红外波段水体-陆地边界效应的水体修正以及用于将栅格数据转化为点数据的矩形选取. 同时,对使用辐射传递模型法和Atmcorr大气校正参数反演所得的Landsat 7水温应用了TSM-S2法中基于实测数据与反演数据的逻辑回归经验公式,以修正其在高温以及低温区域出现的系统误差. 由于本小节内使用的4种数据处理方法和其对反演水温精度带来的改善, 这部分内容不是本文核心所在,详细讨论见参考文献[19]. 处理后的Landsat水温数据按照数据来源与处理方法不同共分为7种,如表 4所示.

表 4 Landsat反演所得7种类型水温数据 Tab. 4 Seven types of Landsat-derived RWT (river water temperature)
2.2 处理MODIS温度数据

本文主要对MODIS温度数据进行如下两步处理:(1)提取每个研究区域内水体占比最多的5个栅格;(2)利用质量管控图层剔除质量较差的数据,如图 2所示. 具体处理原因和方法为:由于MOD11A1空间分辨率仅为1 km,在寸滩以及向家坝坝前部分河宽小于1 km的区域内并不具有纯水栅格,即包含水体的栅格中参杂有一定的陆地. 因此本文首先使用Landsat生成的30 m水体掩膜代入至MODIS 1 km栅格中,确定MODIS数据中单个栅格的水体占比,并分别选出3个研究区域内水体占比最高的5个栅格,提取研究时段内这些点的MODIS水温. 第2步则是通过MOD11A1数据中的质量管控图层,剔除掉质量较差(QC_day > 1)的数据,获得有效的MODIS水温数据.

3 结果分析

经数据整合,将表 4中同一颗卫星下不同种类数据归为一组,2004-2016年期间,寸滩共有27组Landsat 5有效水温数据,68组Landsat 7数据和807幅MODIS数据,三峡坝前共有37组Landsat 5水温数据、77组Landsat 7数据和1370幅MODIS数据. 而对于向家坝坝前区域,2015年全年共有6组Landsat 7水温数据和64幅MODIS数据. 从可用数据量的角度,MODIS水温数据远远多于Landsat系列卫星. 同时,以往研究结果[11]表明3个研究区域之间的Landsat反演水温误差差异较小,因此本文将3个研究区域视为一个整体,在不考虑卫星过空时间和实测时间差的情况下,分析了表 4中7类Landsat温度数据的绝对误差、均方根误差(RMSE)和偏差(Bias),并加入了3个研究区域内的MODIS水温误差,如表 5所示.

表 5 MODIS水温产品和Landsat反演水温的绝对误差、RMSE和偏差 Tab. 5 Error, RMSE and Bias of MODIS and Landsat-derived RWT
3.1 MODIS水温产品与Landsat反演水温比较

总体而言,MODIS水温产品相较于Landsat反演水温误差偏大. 尤其是在寸滩以及向家坝坝前两个研究区域,MODIS水温RMSE分别高达6.24和4.86℃. 经初步分析,可能原因为两个区域内河道宽度较窄,水体温度受陆地影响严重. 因此针对3个不同的研究区域,分别选取研究区域内水体占比最高的5个栅格及其反演水温误差,汇总于图 3中.

图 3 MODIS数据水体占比与均方根误差的关系 Fig.3 Relationship between MODIS waterbody ratio and RWT RMSE

图 3所示,寸滩站MODIS栅格水体占比处于0.49~0.61之间,RMSE为6.24~8.25℃;向家坝坝前水体占比为0.79~0.99,RMSE为4.86~5.92℃;三峡坝前水体占比为0.99~1.00,RMSE为2.80~4.17℃. 对比不同研究区域水体占比与反演水温误差可以发现,随着水体占比的增加,无论是3个研究区域之间还是每个研究区域内不同栅格之间,反演水温误差均呈现下降趋势. 因此使用MODIS水温数据之前,应需考察研究区域有内无纯水栅格或河道宽度是否满足MODIS要求. 一般而言,由于MODIS数据使用的是地理位置固定的网格,单个网格的长度为1 km,考虑到河流流向问题,理论上当河道具有一段长度和河宽均大于2 2 km的河段,即存在至少一个纯水栅格. 然而比较三峡坝前多个纯水栅格可以很明显地发现,除去水体占比这一影响因素,MODIS反演河道水温数据精度可能也与其他条件相关,例如:水体栅格与岸边陆地的距离、与水文测站的位置关系以及河水表面有无覆盖物等. 现阶段尚无办法判断并剔除这些因素对水温数据的影响. 同时,对于长江上游流域而言,三峡坝前属于其中河道较为宽广的区域,然而三峡坝前水温RMSE远远高于以往研究中MODIS数据在大型湖泊上的误差(RMSE为1~1.5℃)[26]. 因此,MODIS水温数据暂不具备精确获取长江上游河道水温的能力. 而相比之下,陆地-水体边缘效应对Landsat卫星反演水温的影响则更容易被识别和剔除. 以往研究表明[11, 19],利用TSM-S2法除去水体边缘效应后寸滩站的L7-RTM-Atmcorr水温均方根误差从1.93℃降低到了1.66℃,降幅约14 %. 因此目前阶段相较于MODIS水温数据,Landsat反演所得水温受陆地影响较小,多种方法反演所得水温RMSE处于1~2℃. 综上所述,在不考虑有效水温数据量的前提下,在长江上游流域使用Landsat热红外影像反演所得水温相较于MODIS水温产品精度更高,可用性更强.

3.2 JM & S单通道算法与辐射传递模型法比较分析

在衡量过遥感数据差异对反演水温的影响之后,本文继续比较分析了不同反演方法的适用性. 如表 5所示,在3个研究区域内,利用辐射传递模型法反演的Landsat 5、Landsat 7水温RMSE处于1.00~1.62℃之间;然而使用JM & S单通道算法获取的同研究区域、时段内的水温数据RMSE却处于2.47~3.28℃之间,远远低于以往研究中JM & S法反演精度[36],因此对JM & S法在长江上游流域可能存在的问题进行了研究,分别生成了Landsat 5与Landsat 7反演水温偏差与大气水汽含量的关系(图 4).

图 4 JM & S反演水温偏差与大气水汽含量的关系 Fig.4 Relationship between water vapor content and bias of JM & S retrieved RWT

分析图 4,可以清楚地观察到偏差与大气水汽含量间存在一定的关系:当大气水汽含量w < 4 g/cm2时,大部分水温偏差都处于[-4℃, +4℃]的灰色区间内,且大多数聚集于0℃线两边;而当w>4 g/cm2时,随着w的增长大部分点的水温偏差游离于灰色区域以外,且0℃中心线附近数据点较稀疏. 当剔除掉大气水汽含量过大(w>4 g/cm2)反演所得的水温数据后,Landsat 5反演水温RMSE从3.28℃下降到1.96℃,Landsat 7反演水温RMSE从2.47℃下降到1.76℃. 在Jiménez-Mu${\rm{\bar n}}$o对该方法的分析中曾指出,当大气水汽含量小于0.5 g/cm2或大于3 g/cm2时,利用经验公式(5)、(6)模拟生成的大气校正参数会与实际大气参数间存在较大的差异,导致反演水温误差过大[32]. 因此,本研究收集了长江上游流域内4个高空气象站探测的多年大气水汽含量数据,数据表明长江上游流域一年365 d中平均有186.5 d大气水汽含量大于3 g/cm2,且主要分布在6-9月. 所以,当剔除这些因大气水汽含量过高而出现问题的数据,并忽略季节性云层对Landsat卫星有效数据的影响,使用JM & S单通道算法得到的温度数据量约为辐射传递模型法的1/2. 除此之外,由于JM & S单通道算法缺少夏季高温数据,结合以往文献的结论[16],其反演所得水温相较辐射传递模型法会更难精确模拟高温区域日水温. 另一方面,在长江上游流域内高空气象站分布较为稀疏,能精确获取大气水汽含量的区域相对较少. 综合上述结论,目前利用JM & S单通道算法反演长江上游流域内河道水温存在两方面的局限性:(1)研究区域受高空气象站地理位置限制;(2)反演水温精度受大气水汽含量影响. 相较于JM & S单通道算法,辐射传递模型法的适用区域更广,反演精度更高.

3.3 不同大气参数下辐射传递模型法反演水温的精度分析

在比较过不同大气校正方法的反演精度之后,本文继续对采用辐射传递模型法并利用不同大气校正参数和修正方法所得的5种数据进行了精度分析. 其中包括(1)L7-RTM-GEOS-5;(2)L7-RTM-Atmcorr; (3)L7-RTM-Atmcorr-Log; (4)L5-RTM-GEOS-5;(5)L5-RTM- Atmcorr. 由于5种数据均使用了辐射传递模型法RTM,因此后续章节及图中省略所有RTM.

表 5所示,对于Landsat 5热红外波段数据,使用Atmcorr和GEOS-5大气校正参数反演所得水温RMSE分别为1.13和1.33℃. 而在相同研究区域内,L7-Atmcorr反演水温RMSE为1.62℃,大于L7-GEOS-5的1.52℃. 因此,通过比较均方根误差无法得到两种大气校正参数之间的区别. 因此针对上述5种数据,各自生成了基于数据偏差的小提琴图和箱形图,如图 5所示. 从图 5中可以观察到,使用Atmcorr大气校正参数反演所得水温相较GEOS-5具有更矮的箱形与更宽的密度峰值,但其偏差密度峰值和偏差的中位数距离0℃更远,例如L7-Atmcorr(2)的密度峰值出现在大约-0.3℃,L7-Atmcorr的密度峰值(4)位于0.4~0.5℃. 而对于GEOS-5数据,其箱形更长且小提琴图图宽更狭长,并具有更大和更小的数据偏差极值. 以上两种现象说明,使用GEOS-5反演所得的水温数据相较Atmcorr有更高的概率存在数据异常值,从而影响整体精度. 而Atmcorr反演水温,相较于GEOS-5可能存在更严重的系统性偏差.

图 5 5种反演水温数据偏差的箱形和小提琴复合图 Fig.5 Box plot and violin plot of five types of Landsat-derived river water temperature

对于Atmcorr存在的系统性偏差,往往可以进一步使用经验公式进行修正. 本文对L7-Atmcorr数据进行了逻辑回归修正并得到了L7-Atmcorr-Log数据. 如图 5(3)所示,相较于L7-Atmcorr,修正后的数据具有更小的系统性偏差,即数据中位数以及密度峰值更趋近于0℃,同时具有误差极值也回归到正常区间. 而针对GEOS-5中存在的数据异常值问题,本文分别生成了Landsat 7和Landsat 5数据下的GEOS-5、Atmcorr反演水温与实测水温散点图(图 6).

图 6 Landsat 5与7热红外数据的GEOS-5、Atmcorr反演水温与实测水温散点图 Fig.6 Scatter plot of measured RWT with Landsat-derived RWT by using GEOS-5 and Atmcorr

图 6所示,尽管Atmcorr与GEOS-5之间存在大量重合或相似数据组(红色空心圈包含蓝色×的数据组), 仍有4组GEOS-5数据出现了较大偏差. Landsat 7数据中1组,Landsat 5数据中共有3组,且4组数据所在温度区间相对分散,说明数据异常与温度之间不存在一定的关系,无法通过经验公式修正. 因此需进一步分析其成因,并尽可能地在数据处理过程中找到方法进行剔除.

3.4 辐射传递模型法中各参数敏感性分析

首先,通过解析辐射传递模型法中的核心公式,公式(2),可知由于水体在热红外波段的比辐射率较高,大部分大气下行辐射会被水体吸收,因此在已知大气下行辐射与大气上行辐射为同量级数据的情况下,影响反演水温大小以及精度的主要因子为Lup/τ,即大气上行辐射与大气透过率的比值. 且当Lupτ的比值越大时,反演水温会对大气校正参数越敏感,即更容易产生较大水温误差,因此为了尝试分析和解决GEOS-5反演水温存在的异常值问题,生成了GEOS-5反演水温偏差与大气上行辐射,大气透过率以及两者比值的相关关系(图 7).

图 7 Landsat 5与7 GEOS-5反演水温与大气参数间的关系 Fig.7 Relationship between atmospheric parameters with Landsat-derived RWTs by using GEOS-5

图 7中3张图的形状特征相似,GEOS-5反演水温偏差数据上下包络线呈喇叭口型:当大气上行辐射相对较小或大气透过率校高时,GEOS-5反演水温偏差聚集在[-4℃, 4℃]之间,但随着大气上行辐射的变大和大气透过率的变小,水温数据的偏差相对0℃刻度线越来越分散,且有部分数据超出偏差正常区间,这部分数据即为3.3节中所提到的异常数据点. 根据上述特性,本文通过对大气校正参数设定相应阈值来减少异常数据点的出现概率. 设定Lup < 4.5,τ>0.4,Lup/τ < 11.5 3条阈值分界线进行L7-GEOS-5数据划分,3种阈值分别剔除了28、23、25组数据,皆约占数据总量的1/6,整体水温数据误差RMSE由1.52℃降低到了1.20℃,低于L7-Atmcorr反演水温数据的均方根误差(1.62℃),但高于L7-Atmcorr-log修正数据的均方根误差(1.00℃). 但相较于L7-Atmcorr-log修正数据,这种处理方法的优势是其不需要实测水温资料,可以应用于无测站地区. 综上所述,对于GEOS-5反演水温而言,由于其具有与大气校正参数相关的数据异常值,建议分析并设定校正参数阈值以提高反演精度. 而对于Atmcorr反演水温所具有的系统性误差,在有实测水温资料地区可以使用逻辑回归等方法进行校正从而提高水温精度.

4 结论

本文利用不同数据处理流程获取了长江上游河段MODIS和Landsat水温数据,结果表明在长江上游河段各种方法和数据都具有各自的优缺点:

1) MODIS温度产品的优点在于其可以提供日水温数据,但由于自身空间分辨率无法满足长江上游流域河道要求,导致数据精度过低,因此不建议使用.

2) 对于JM & S单通道算法,其主要优点为只需大气水汽含量这一个大气参数. 但是由于长江上游流域在6-9月大部分日期大气水汽含量过高,会严重影响反演水温精度,因此同样不建议使用.

3) 对于辐射传递模型法RTM以及大气参数数据库GEOS-5 FP-IT反演所得数据,其主要缺点是存在数据异常值,本文提出通过设置大气校正参数的阈值(Lup < 4.5,τ>0.4,Lup/τ < 11.5)剔除部分数据来减小数据异常点带来的误差. 但阈值的普适性有待未来在更多研究区域内进一步的研究.

4) 使用辐射传递模型法RTM以及大气参数数据库Atmcorr为众多方法中最为稳定的流程,其数据精度略位于方法(3)所获得的原始水温数据精度方法(3)修正后水温数据精度之间. 其缺点在于存在一定的系统性偏差,如果研究区域内有适量的实测水温数据,可以通过逻辑回归等经验公式进行修正,修正后的数据精度为所有方法中最高.

综上所述,在反演长江上游流域水温且研究区域内有适量实测水温数据时,优先选用RTM-Atmorr-Log数据. 当不考虑有效数据量且研究区域内无实测水温时,建议使用RTM-GEOS-5数据并利用阈值法剔除部分数据点. 对于数据量有限且区域内无实测水温数据时,应选用RTM-Atmcorr原始水温数据.

5 参考文献

[1]
McCullough DA. A review and synthesis of effects of alterations to the water temperature regime on freshwater life stages of salmonids, with special reference to Chinook salmon. Washington: Environmental Protection Agency, 1999.
[2]
Gooding RA, Harley CDG, Tang E. Elevated water temperature and carbon dioxide concentration increase the growth of a keystone echinoderm. PNAS, 2009, 106(23): 9316-9321. DOI:10.1073/pnas.0811143106
[3]
Cole GA, Weihe PE eds. Textbook of limnology. Waveland: Waveland Press, 2015.
[4]
Guo WX, Wang HX, Xia ZQ et al. Effects of Three Gorges and Gezhouba reservoirs on river water temperature regimes. Journal of Hydroelectric Engineering, 2009, 28(6): 182-187. [郭文献, 王鸿翔, 夏自强等. 三峡-葛洲坝梯级水库水温影响研究. 水力发电学报, 2009, 28(6): 182-187.]
[5]
van Vliet MTH, Franssen WHP, Yearsley JR et al. Global river discharge and water temperature under climate change. Global Environmental Change, 2013, 23(2): 450-464. DOI:10.1016/j.gloenvcha.2012.11.002
[6]
Hu H. Impact of reservoirs on downstream water temperature variation research based on river temperature remote sensing[Dissertation]. Wuhan: Huazhong University of Science & Technology, 2018. [胡海畅. 基于河流温度遥感反演的水库蓄水对下游水温影响研究[学位论文]. 武汉: 华中科技大学, 2018.]
[7]
Webb BW, Hannah DM, Moore RD et al. Recent advances in stream and river temperature research. Hydrological Processes, 2008, 22(7): 902-918. DOI:10.1002/hyp.6994
[8]
Torgersen CE, Faux RN, McIntosh BA et al. Airborne thermal remote sensing for water temperature assessment in rivers and streams. Remote Sensing of Environment, 2001, 76(3): 386-398. DOI:10.1016/S0034-4257(01)00186-9
[9]
Minnett PJ, Alvera-Azcárate A, Chin TM et al. Half a century of satellite remote sensing of sea-surface temperature. Remote Sensing of Environment, 2019, 233: 111366. DOI:10.1016/j.rse.2019.111366
[10]
Syariz MA, Jaelani LM, Subehi L et al. Retrieval of sea surface temperature over Poteran Island water of Indonesia with Landsat 8 TIRS image: A preliminary algorithm. The International Archives of Photogrammetry, Remote Sensing and Spatial Information Sciences, 2015, 40: 87. DOI:10.5194/isprsarchives-xl-2-w4-87-2015
[11]
Shi X, Sun J, Shi LD. Derivation of river surface temperature from Landsat thermal infrared data. Journal of Hydroelectric Engineering, 2021, 40(2): 121-130. [石希, 孙健, 史立地. 基于Landsat卫星遥感资料的河流水温反演研究. 水力发电学报, 2021, 40(2): 121-130.]
[12]
Handcock RN, Torgersen CE, Cherkauer KA et al. Thermal infrared remote sensing of water temperature in riverine landscapes. Fluvial Remote Sensing for Science and Management, 2012, 1: 85-113. DOI:10.1002/9781119940791.ch5
[13]
Huang B, Zhao YQ. Research status and prospect of spatiotemporal fusion of multi-source satellite remote sensing imagery. Acta Geodaetica et Cartographica Sinica, 2017, 46(10): 1492-1499. [黄波, 赵涌泉. 多源卫星遥感影像时空融合研究的现状及展望. 测绘学报, 2017, 46(10): 1492-1499. DOI:10.11947/j.AGCS.2017.20170376]
[14]
Alcântara EH, Stech JL, Lorenzzetti JA et al. Remote sensing of water surface temperature and heat flux over a tropical hydroelectric reservoir. Remote Sensing of Environment, 2010, 114(11): 2651-2665. DOI:10.1016/j.rse.2010.06.002
[15]
Lieberherr G, Wunderle S. Lake surface water temperature derived from 35 years of AVHRR sensor data for European lakes. Remote Sensing, 2018, 10(7): 990. DOI:10.3390/rs10070990
[16]
Fu P, Weng QH. A time series analysis of urbanization induced land use and land cover change and its impact on land surfacetemperature with Landsat imagery. Remote Sensing of Environment, 2016, 175: 205-214. DOI:10.1016/j.rse.2015.12.040
[17]
Tavares MH, Cunha AHF, Motta-Marques D et al. Derivation of consistent, continuous daily river temperature data series by combining remote sensing and water temperature models. Remote Sensing of Environment, 2020, 241: 111721. DOI:10.1016/j.rse.2020.111721
[18]
Feng L, Giles MF, Hao D et al. Monitoring thermal pollution in rivers downstream of dams with Landsat ETM+ thermal infrared images. Remote Sensing, 2017, 9(11): 1175.
[19]
Shi X, Sun J, Xiao ZJ. Investigation on river thermal regime under dam influence by integrating remote sensing and water temperature model. Water, 2021, 13(2): 133. DOI:10.3390/w13020133
[20]
Wan W, Li H, Xie HJ et al. A comprehensive data set of lake surface water temperature over the Tibetan Plateau derived from MODIS LST products 2001-2015. Scientific Data, 2017, 4(1): 1-10. DOI:10.1038/sdata.2017.95
[21]
Wang YH, Luo Y, Tan WJ et al. Analysis of temporal and spatial variation process of Dianchi Lake surface water temperature based on MODIS remote sensing images. IOP Conference Series: Earth and Environmental Science, 2021, 658: 012005. DOI:10.1088/1755-1315/658/1/012005
[22]
Xiong YJ, Yin J, Kyaw Tha Paw U et al. How the three Gorges Dam affects the hydrological cycle in the mid-lower Yangtze River: A perspective based on decadal water temperature changes. Environmental Research Letters, 2020, 15(1): 014002. DOI:10.1088/1748-9326/ab5d9a
[23]
Barnes WL, Pagano TS, Salomonson VV. Prelaunch characteristics of the moderate resolution imaging spectroradiometer (MODIS) on EOS-AM1. IEEE Transactions on Geoscience and Remote Sensing, 1998, 36(4): 1088-1100. DOI:10.1109/36.700993
[24]
Ruiz-Verdú A, Jiménez JC, Lazzaro X et al. Comparison of MODIS and Landsat-8 retrievals of chlorophyll-a and water temperature over Lake Titicaca. Beijing: 2016 IEEE International Geoscience and Remote Sensing Symposium (IGARSS), 2016: 7643-7646. DOI: 10.1109/IGARSS.2016.7730993.
[25]
Li ZL, Tang BH, Wu H et al. Satellite-derived land surface temperature: Current status and perspectives. Remote Sensing of Environment, 2013, 131: 14-37. DOI:10.1016/j.rse.2012.12.008
[26]
Tavares M, Cunha A, Motta-Marques D et al. Comparison of methods to estimate lake-surface-water temperature using Landsat 7 ETM+ and MODIS imagery: Case study of a large shallow subtropical lake in southern Brazil. Water, 2019, 11(1): 168. DOI:10.3390/w11010168
[27]
Yu XL, Guo XL, Wu ZC. Land surface temperature retrieval from Landsat 8 TIRS-comparison between radiative transfer equation-based method, split window algorithm and single channel method. Remote Sensing, 2014, 6(10): 9829-9852. DOI:10.3390/rs6109829
[28]
Xu HQ, Lin ZL, Pan WH. Some issues in land surface temperature retrieval of Landsat thermal data with the single-channel algorithm. Geomatics and Information Science of Wuhan University, 2015, 40(4): 487-492. [徐涵秋, 林中立, 潘卫华. 单通道算法地表温度反演的若干问题讨论——以Landsat系列数据为例. 武汉大学学报·信息科学版, 2015, 40(4): 487-492. DOI:10.13203/j.whugis20130733]
[29]
Liu CG, Lu XF, Gao SF. Comparison of brightness temperatures inversed from high and low gain state data of Lansat-7 ETM+ thermal infrared band. Journal of Henan Polytechnic University: Natural Science, 2011, 30(5): 561-566. [刘春国, 卢晓峰, 高松峰. Lansat-7 ETM+热红外波段高低增益状态数据反演亮度温度比较研究. 河南理工大学学报: 自然科学版, 2011, 30(5): 561-566.]
[30]
Hulley GC, Hook SJ, Schneider P. Optimized split-window coefficients for deriving surface temperatures from inland water bodies. Remote Sensing of Environment, 2011, 115(12): 3758-3769. DOI:10.1016/j.rse.2011.09.014
[31]
Chander G, Markham BL, Helder DL. Summary of current radiometric calibration coefficients for Landsat MSS, TM, ETM+, and EO-1 ALI sensors. Remote Sensing of Environment, 2009, 113(5): 893-903. DOI:10.1016/j.rse.2009.01.007
[32]
Barsi JA, Barker JL, Schott JR. An atmospheric correction parameter calculator for a single thermal band earth-sensing instrument. Toulouse: 2003 IEEE International Geoscience and Remote Sensing Symposium. Proceedings (IEEE Cat. No. 03CH37477), 2003: 3014-3016. DOI: 10.1109/IGARSS.2003.1294665.
[33]
Barsi JA, Schott JR, Palluconi FD et al. Validation of a web-based atmospheric correction tool for single thermal band instruments. San Diego: Optics and Photonics 2005. Proc SPIE 5882, Earth Observing Systems X, 2005: 58820E. DOI: 10.1117/12.619990.
[34]
Lucchesi R. File Specification for GEOS-5 FP-IT (Forward Processing for Instrument Teams), 2013.
[35]
Jiménez-Muñoz JC, Sobrino JA. A generalized single-channel method for retrieving land surface temperature from remote sensing data. Journal of Geophysical Research: Atmospheres, 2003, 108(D22): 2003JD003480. DOI:10.1029/2003jd003480
[36]
Jimenez-Munoz JC, Cristobal J, Sobrino JA et al. Revision of the single-channel algorithm for land surface temperature retrieval from Landsat thermal-infrared data. IEEE Transactions on Geoscience and Remote Sensing, 2009, 47(1): 339-349. DOI:10.1109/TGRS.2008.2007125