黄 路，邵 浩，张彩云，夏文杰，杨星晨

(1.国家海洋技术中心漳州基地筹建办公室, 福建 厦门 361007；2.厦门大学海洋与地球学院, 福建 厦门 361102；3.91039部队, 北京 102401；4.声场声信息国家重点实验室、中国科学院声学研究所, 北京 100190)

高分辨率海表温度(Sea Surface Temperature,SST)在亚中尺度和小尺度海洋锋结构[1]、河流水动力过程[2]、海气相互作用[3]、水产养殖[4]、海冰热异常[5]、海岸带水质监测[6]等研究中发挥着重要作用。LandSAT 8卫星于2013年2月11日发射，有10、11两个热红外波段，中心波长分别为10.9 μm和12.0 μm，经重采样后空间分辨率均为30 m，其长时间稳定运行为获取高分辨率SST提供了大量可靠数据源。目前，LandSAT8SST相关反演算法主要有单窗算法[7-8]、普适性单通道算法[9-10]、劈窗算法[8,10]以及辐射传输方程算法等。其中，基于美国国家环境预报中心(NECP)再分析数据集、中等光谱分辨大气透过率和辐射计算软件包(Moderate Spectral Resolution Atmospheric Transmitance, MODTRAN)等大气辐射传输模型建立的辐射传输方程算法[11-13](以下简称辐射方程法)不依赖于地面观测站点资料，因而非常便于业务化遥感反演SST，美国航空航天局(NASA)也对其各个输入参数的计算过程进行集成，并开发出大气校正参数计算器(https://atmcorr.gsfc.nasa.gov/)，内有专门针对LandSAT8第10波段的大气校正模块。

虽然NECP再分析数据集提供了地表到30 km高度的大气廓线数据作为输入参数，较大程度提高了整层大气廓线的估算精度，但其较低的时间和空间分辨率、相对简单的插值计算方式以及数据集本身的准确性等都会影响辐射方程法输入参数的计算精度，同时也使得该算法最终反演的SST结果存疑。有学者针对辐射方程法在反演环境一号B卫星(HJ-1B)SST存在较大偏差的情况，从修正大气透过率参数入手，提出普适性的辐射传输方程改进算法[13](以下简称改进算法)，并获得了较高HJ-1BSST反演精度，但该改进算法还未推广应用在LandSAT8卫星上，其实际反演效果还有待进一步检验。此外，改进算法使用研究区中心位置附近一点处的大气透过率来近似代表整个研究海域的大气透过率，但是即使在小范围研究海域，大气透过率也可能会出现较大变化，尤其是在水汽含量较高的夏季，或有较多的不易剔除的薄云存在等极端情况下，该算法的反演精度和适用条件也需要进行进一步评估和分析。

因此，本研究选择冬、夏各2景LandSAT8图像，其中1景夏季图像研究区内云覆盖较多，在修正算法部分参数的基础上，分别反演出算法改进前后的LandSAT8SST，并使用MODISSST和浮标SST对结果进行验证和比较。最后以电厂温排水遥感监测作为改进算法的应用案例，分析了其时空分布情况。

1 数据与方法

1.1 研究区

以厦门湾及其邻近海域作为研究区，如图1所示。取金门岛南侧海域一点(24.38°N，118.41°E)作为研究海域的中心点，在中心点位置附近分别平行岸线方向和垂直岸线方向划设2条断面，在断面上随机均匀选取20个点作为算法验证点。

图1 厦门湾及其邻近海域研究区示意图

1.2 数据源及预处理

4景LandSAT8第10波段热红外图像来源于对地观测中心(http://ids.ceode.ac.cn/)，成像日期分别为2013-08-04、2014-01-27、2019-08-21、2019-12-11，除2019-08-21图像外，其他3景图像质量较好，云覆盖少，水汽含量相对均匀；与LandSAT8图像同日1 km的Aqua MODISSST产品[14-15]下载于美国NASA水色卫星网站(http://oceancolor.gsfc.nasa.gov/)。此外，还使用了2014年1月27日的5个浮标验证数据，来源于福建省海洋预报台海洋观测网(http://www.fjhyyb.cn/Ocean863Web_MAIN/default.aspx?r=1#hygc)，浮标数据采集时间与当日卫星成像时间相差小于5 min。

对4景LandSAT8热红外图像进行几何精校正、图像配准，控制点误差小于0.5个像元；裁剪出研究区，并利用水体像元和非水体像元的光谱特征差异，人机交互判读提取出研究海域纯水体像元。

1.3 基于改进算法反演LandSAT8 SST

大气表观辐亮度LTOA可用下式近似表示[16-17]：

LTOA=[ε·LT+(1-ε)Ld]·τ1+Lu

(1)

式(1)中：ε是海水发射率，这里取0.988 8；Lu和Ld分别表示中心点处的大气上行和下行辐亮度[W/(m2·sr·μm)]，由大气校正参数计算器求得；LT是海表辐亮度[W/(m2·sr·μm)]；τ1是中心点处修正的大气透过率。

LandSAT8LTOA通过以下定标公式计算：

LTOA=ML·Qcal+AL

(2)

式(2)中：系数ML和AL可在MTL.txt元数据文件中查找得到，Qcal为第10波段DN值。

对普朗克公式进行简化计算[18]，将海区SST分成10～20 ℃和21～33 ℃两个区间，分别进行线性拟合得到LandSAT8LT和LandSAT8SST转换关系式：

SST10～20= 7.771 5LT-47.316(R2= 0.999 8)

(3)

SST21～33= 6.992 3LT-40.481(R2= 0.999 8)

(4)

τ1的计算过程参考前人做法[13]，即将中心点处的MODISSST近似看作真值并作为该点处的LandSAT8SST，分不同的SST区间将其代入公式(3)或(4)，得到中心点处的LandSAT8LT，再将LandSAT8LT连同中心点处的Lu、Ld和LTOA一起代入式(1)即可得到修正的τ1。

此外，为了便于比较，本研究还将未修正的τ，即基于大气校正参数计算器计算的中心点处的τ代入式(1)，计算得到基于辐射方程法反演的LandSAT8 SST。

1.4 算法精度评价

2种算法计算的大气透过率偏差用均方根误差rmse来表示，分别为：

(5)

(6)

式(5)、(6)中：τ2i为使用1.3节计算τ1的方法求取的第i个MODIS验证点处的大气透过率τ2，下同。

2种算法反演的LandSAT8SST偏差主要用平均偏差bias和均方根误差rmse来表示，分别为：

(7)

(8)

式(7)、(8)中：SSTi为MODISSST或浮标SST。

使用大气透过率标准差σ(τ2)定性表示研究海域大气透过率分布的均匀程度，计算公式为：

(9)

2 结果与讨论

2.1 改进算法精度验证及适用性探讨

MODISSST精度验证结果如表1所示，辐射方程法计算的τ及其反演的LandSAT8SST存在不同程度的偏差，偏差大小与季节无关，SST偏差(包括平均偏差bias和均方根误差rmse)最大超过2.5 ℃，最小仅为0.5 ℃左右。相比辐射方程法，改进算法则对大气透过率偏差和SST偏差均有明显改善，除2019-08-21图像外，其他3景图像τ1偏差在0.005左右，SST偏差在0.5 ℃左右；而对于2019-08-21图像，改进算法对其SST反演精度也有较大提高，平均偏差由2.52 ℃降为1.09 ℃。浮标SST精度验证结果如表2所示，改进算法在近岸海域对SST偏差改善效果同样显著，辐射方程法反演的SST偏差超过了1.5 ℃，而算法改进后平均偏差bias和均方根误差rmse分别减少到0.14、0.18 ℃。

表1 MODIS SST精度验证结果

表2 浮标SST精度验证结果

进一步对改进算法的适用条件进行分析，如图2所示，σ(τ2)与τ1偏差和SST偏差相关系数分别为R2[rmse(τ1)]=0.648 5、R2[bias(SST)]=0.920 7、R2[rmse(SST)]=0.934 0，可看出改进算法SST反演精度与研究海域大气透过率分布均匀程度呈显著的正相关关系，而研究海域大气透过率分布不均匀则可能是该算法产生较大偏差的主要原因，即对于大气透过率分布均匀的研究海域，改进算法对其SST偏差具有显著改善，而当水汽含量不均匀，有云覆盖，尤其是不易剔除的薄云存在时，往往会造成局部大气透过率分布存在较大突变，在此情况下，基于中心点计算的τ1可能就难以代表整片研究海域大气透过率的平均状况，如2019-08-21图像σ(τ2)在4景图像中最大，相应的SST偏差也最大。或可以通过缩小研究海域范围、增加τ1的计算数量，以及进行更严格细致的薄云剔除等措施来提高改进算法的反演精度。

图2 线性相关性分析结果

2.2 改进算法在电厂温排水遥感监测上的应用

大中型核电或火力发电厂循环冷却水排放入海，会引起一定范围的海水增温，进而影响浮游动、植物的多样性，增加海洋生物病害和赤潮暴发频率[19]。厦门湾及其邻近海域主要有后石、嵩屿、晋江3座火力发电厂，如图1所示，后石电厂总装机容量可达420万千瓦，最后一号机组于2009年建成投产；嵩屿电厂总装机容量120万千瓦，二期两台机组分别于2006年3月、8月转入商业运行；晋江电厂总装机容量140万千瓦，4台燃-蒸气联合循环机组于2010年10月竣工投产。

基于改进算法反演的各电厂表层温排水分布，如图3所示。后石电厂表层温升现象最显著，呈射流状，沿岸线向东南方向流动，主流轴越过了鸟嘴尾角，羽流尾部形态不固定，差异较大；嵩屿电厂未形成较为明显的温升羽流痕迹，其排水口附近SST分布特征与相邻海岸海水较为相似，温排水增温效应基本可忽略不计；晋江电厂温升羽流呈扇形分布，主流轴近似垂直岸线向东，流轴短，影响范围小，部分羽流尾部存在向南扩散的特征。

进一步定量分析后石、晋江电厂的温排水特征。从排水口出发，沿主流轴方向做温升羽流的扩散断面，如图3黑线曲线所示，查找断面无明显降温变化的初始位置，将该位置处的温度作为基准温度，分别统计温升高于基准温度0、3、5 ℃的热污染半径和面积，如表3所示。后石电厂排水口附近存在较为稳定的高温水体(温升>3 ℃的水体，下同)，其温升幅度和扩散影响范围最大，各图像最大温升和热污染面积均超过了5 ℃和1.5 km2；晋江电厂热污染半径和面积较小，主要集中在排水口附近，但不排除部分年份排水口附近出现高温水体的可能，如在2019年8月21日图像中，高温水体面积占比达到了20.69%。

表3 羽流温升结果统计

图3 厦门湾及其邻近海域各电厂表层温排水分布图

从应用效果来看，改进算法对于遥感监测电厂表层温排水具有很好的优势，但也需注意遥感手段仅能描述和说明特定时刻电厂表层温排水情况，对于其实际排放通量和表层以下温排水情况，还需结合实测数据和数值模拟手段做进一步分析。

3 结论

本研究使用MODIS SST和浮标SST分别对算法改进前后反演的LandSAT8SST进行验证和研究比较，并初步分析了厦门湾及其邻近海域的电厂温排水时空分布情况，得出的主要结论如下：

(1)辐射方程法反演的LandSAT8SST精度不稳定，可能存在较大偏差。而改进算法则对τ1和SST偏差均有明显改善，MODISSST验证表明除2019-08-21图像外，其他3景图像SST偏差在0.5 ℃左右；浮标SST验证表明其平均偏差bias和均方根误差rmse分别减少到0.14、0.18 ℃，基本可以忽略不计。

(2)改进算法SST反演精度与研究海域大气透过率分布均匀程度呈显著的正相关关系，相关系数分别为R2[bias(SST)]=0.920 7、R2[rmse(SST)]=0.934 0，大气透过率分布不均匀可能是该算法产生较大偏差的主要原因。

(3)改进算法对于遥感监测电厂表层温排水具有很好的优势，从反演的SST图像来看，嵩屿电厂表层海水温升羽流不明显；后石电厂表层温升现象最显著，温升幅度和扩散影响范围最大，其排水口附近存在较为稳定的高温水体；晋江电厂温排水呈扇形分布，影响范围小，主要集中在排水口附近。