基于遥感数据的三维温度场参数化分析方法研究

邢霄波，徐永生*，贾永君，黄超

(1.中国科学院海洋研究所，山东 青岛 266071；2.中国科学院大学，北京 100049；3.青岛海洋科学与技术试点国家实验室海洋动力过程与气候功能实验室，山东 青岛 266237；4.国家卫星海洋应用中心，北京 100081)

1 引言

从古至今，海洋对于人类来说，都扮演着极其重要的角色。海洋中蕴藏着丰富的资源，具有广阔的开发前景。海洋在我国的经济发展、科学研究以及军事等方面，都发挥着越来越重要的作用，例如远洋货物运输、海上科研活动、军事探测等日益增多。一系列的变化使得海洋科学越来越被人们所重视，人们对水下环境信息的精度要求也随之提高[1]。海水的温度、盐度是海洋水文的基本要素，温度和盐度的时空分布和变化，是海洋科学研究的基本内容的组成部分，海洋中发生的一切现象几乎都与海水温度和盐度密切相关。

邻近我国的海域具有极度复杂的海洋环境要素分布情况，这些要素包括海水的温度、盐度、密度及其跃层。中国邻近海域以及西北太平洋区域，同时被西边界流（黑潮及黑潮延伸体）和季风影响。该海域的地形复杂，存在陆架海区（黄海、南海北部、渤海、东海）和深海区域（南海的海盆、西北太平洋），因此对于这一海洋区域的温度、盐度及其跃层的分布和变化机理的研究具有重要并且典型的科学意义。温跃层的分布情况不仅影响水团的垂直边界，而且温跃层的形成机理也与海洋环流、内波等密切相关[2]。因此，关于温度变化情况的研究对于水下物质扩散、湍流扩散和海洋其他的理论研究都有重大的意义。

海水温度、盐度作为最基本的海洋环境参数，是描述海水性质的重要物理量。然而在海洋的现场观测中，Argo 浮标具有随波自由漂动的特点，导致了观测剖面在时间和位置上的不定性，以及不太均匀的空间分布。现有的CTD、XBT 等剖面实时观测资料不足以构成实时的、大范围、满足精度要求的三维温盐场[3-4]。

由于现场观测的种种不足，卫星遥感观测反演海洋动力要素顺势而生。利用海面信息重构温盐剖面的方法在20 世纪已被提出。经过数十年的发展，形成了3 种较为成熟的方法：动力学方法[5]、变分法[6]和统计分析法[7-9]。还有基于经验正交函数（EOF）分解从而回归分析的海表面温度（SST）和海表面高度（SSH）的地下热结构重建的方法[10]。基于美国海军的模块化海洋数据同化系统（MODAS）的静态气候学（隐含平均季节信号）数据与历史水文资料，构建“残余最大的经验模式（GEM）”，据此构建三维温度场[11]。我国在这方面的研究工作开展的时间较晚，但发展迅速。国内的海洋科学家基于EOF 分析方法，生成了中国海域(此处指中国海洋领土)及其邻近海域历史温度的主要特征向量场，在多次温度观测的基础上重建了良好的温度剖面，并利用该模型对温跃层进行了估算[12]；基于人工神经网络方法，利用海面温度、海面风速以及海面气压数据反演了南海区域近海面气温[13]；基于历史观测的温盐剖面数据资料，采用回归分析，研究了海面温度异常、海面动力高度异常与温度剖面异常的相关关系，接着利用卫星遥感海表面温度和卫星观测的海面高度构建了三维温盐场[3]；基于卫星测高和测温数据结合实测数据，运用统计分析和数据同化方法重构了南海海域三维温盐场[1]；基于海面位势高度和海表温度数据，构建了西太平洋三维温度场[14]。由于海洋内部三维温盐场的复杂性，利用卫星遥感观测反演温盐场的技术仍然需要不断改进，以便为海洋科学研究提供更高精度的温盐数据。

西北太平洋黑潮是著名的西边界强流，携带着太平洋海域巨大的热量进入到中国海域，从而影响了中国近海的温盐分布情况[15]。台湾暖流、对马暖流和黄海暖流等是黑潮的主要分支，并且携带着高温高盐水，它们的变化对于黄海、渤海、东海的温盐分布、水团消长都有很重要的影响。黑潮及黑潮延伸体水系的各层水团分布变化情况对东海、东南海域温盐的垂直分布有十分重要的影响。本文利用经过严格的质量控制的西北太平洋海域的历史实测数据，应用参数化的方法对该海域的三维温度结构进行拟合。

2 数据与方法

2.1 数据资料

2.1.1 Argo 剖面资料

Argo 温度剖面数据来自于全球海洋Argo 浮标温盐数据集（http://www.argo.ucsd.edu/About_Argo.html）。Argo 建立在其他海洋上层观测网络的基础之上，在时间空间上扩大了它的深度范围、覆盖范围以及提高了精度，并增加了盐度测量的数据。Argo 浮标通过质量控制，将对上层海洋变化情况和海洋气候变化，包括热量输送等提供较为定量的描述[16]。本文实验中使用的是太平洋海域某一区域的Argo 数据集。实验区域如图1 所示。由于该区域2002 年之前的数据资料较少，资料参差不齐，所以选取2002 年至2019年5 月的数据资料。首先，对数据进行质量控制，去掉了测量中明显不合乎该海域水文特性的点以及坏点。然后对温深、温盐深剖面曲线进行质量控制，得到实验剖面。图2 是研究区域2002 年至2019 年5 月的Argo 数据年分布。

图1 西北太平洋海域及研究区域地形分布Fig.1 Topography of the northwest Pacific Ocean and the study area

2.1.2 海表温度

美国国家海洋和大气管理局（NOAA）每日（1/4）°最优插值海表面温度（OISST）是通过将来自不同观测平台的观测结果结合在一个常规全球网格上构建的一种数据集。通过插值方法来填补空白区域，得到一个较为完整的海洋温度图。这里包括了卫星观测和船舶观测（相对于参照浮标）的偏差调整，补偿平台之间差异和传感器产生的偏差。在1991 年皮娜图博火山爆发期间，有研究证明上述补充操作是至关重要的，其中火山气溶胶的广泛存在致使红外卫星观测温度比实际温度低很多[17]。

2.1.3 验证数据

我们选取了两种以Argo 为基础的网格化的产品，使用中国Argo 实时资料中心的BOA-Argo 三维网格温盐资料集以及来自英国气象局的EN4.2.1 数据集来验证文中重建的温度剖面。其中，BOA-Argo 是通过Argo 数据衍生的数据产品。中国Argo 实时资料中心利用一种简单、有效，并且易于操作的逐步订正的方法[18]结合一种混合层模型（用来反推对应的Argo 表层温度和盐度）[19]构建了全球海洋2004 年1 月至2018 年12 月的三维网格温、盐度资料集及相应的衍生产品。该产品包括了逐年、逐月以及年、月平均等多种时间分辨率的资料。空间分辨率为1°× 1°的BOA-Argo 资料不仅保留了原始数据中的信号，也消除了一些噪声。因此BOA-Argo 数据集在捕捉中尺度信号方面更优质。

组织实施阶段首先要明确各项目中不同阶段的学习任务和相应的学习内容，即明确在项目中“先什么”“后做什么”的问题，也就是将项目各阶段的学习内容用一个合乎逻辑的顺序排列出来，将每个子项目的重难点呈现，内容设计是根据旅游专业发展需求，结合学校实际实训条件，对普通话所有教学内容进行加工整合再放入所设计的项目中。项目教学法的应用在确保普通话的教学质量的同时，也全面提高学生的普通话水平及应用能力。具体实施内容举例：

图2 2002 年至2019 年5 月Argo 数据年分布Fig.2 Annual distribution of Argo data from 2002 to May,2019

EN 数据提供逐月的温度和盐度（如果有盐度数据的话）资料，数据集有规则的1°的水平网格和42 层垂直层。EN4 数据接收了来自各种海洋剖面仪器的数据，该数据集的主要数据源是WOD09，以此为基础添加了其他数据：首先，是北极气候盆地海洋学（the Arctic Synoptic Basin Wide Oceanography，ASBO）项目期间编译的数据，以提高北极地区数据覆盖率；第二，来自全球的温度和盐度剖面计划（the Global Temperature and Salinity Profile Program，GTSPP）的数据；第三，来自于Argo 数据全球编译中心（the Argo global data assembly centres，GDACs）的2000 年以后的Argo数据。本文采用的EN4.2.1 最新数据集相较于EN4.2.0数据集有两个方面的改善：EN4.2.1 包含了源数据的最新下载；拆除了6 个泊系浮标，在1988-1900 年间这些浮标造成拉布拉多海异常的温度。除此之外，还增加了不确定性估计。因此，基于上述良好性能，两种数据集适用于验证算法的有效性。

2.2 算法简介

大多数温度剖面都表现出混合层、温跃层和深层（图3），为了使模型的使用更具有普遍性，我们假设在温跃层下面有两个深层（图4）。如果在某些海域或者某些时期两个深层的垂直梯度相等，则仍然变为1 个深层。除此之外，我们在混合层和温跃层之间也添加1 个区域，称为夹层；在温跃层和深层之间添加1 个过渡层[20]。这样我们就把整个温度剖面划分为了6 层，那么就可以用6 层的参数模型来表示研究区域的温度曲线的垂直分布情况。这样每一个观测剖面都可以用1 组参数表示，其中包括了海表面温度（SST）、混合层深度（MLD）、温跃层梯度和深层梯度等。我们假设混合层梯度（约等于0）和温跃层梯度（很大）是恒定的，假设剩余其他各层的梯度随深度线性变化[21]，用平均梯度作为这些层的代表值，而夹层梯度是混合层梯度和温跃层梯度的平均值，表示为

图3 Argo 观测数据表征的温度剖面Fig.3 Temperature profile characterized by Argo observation data

我们在梯度空间分析剖面资料[22]，基于

式中，GT是温度梯度;z是重构温度对应深度;T则为该深度温度;b为校正的常量。我们化简得到下列参数模型[21]：

其中各变量含义如表1 所示。

因为Argo 现场观测的剖面数据某一部分会较为稀疏，可能会导致误差增大，所以我们先将Argo 观测资料插值设为Δz=5 m，则zi=zi+1-5 m,并且z0=0，得到高分辨率的温度剖面。在把这6 个参数方程应用于每一个观测剖面之前，我们需要在高分辨率的温度剖面的垂直方向上确定5 个深度猜想值（d1,d2,d3,d4,d5）作为第一猜想，并确定高分辨率的温度点，从而获得6 个温度梯度这种方法可以从观测数据中得到最重要的特征。

图4 温度剖面表征的模型特征Fig.4 Model characteristics represented by temperature profile

得到第一猜想值以及海表面温度值后，我们将每一个观测剖面拟合到参数方程中，得到零阶的拟合梯度结果。接着，我们通过使用迭代方法获得最佳的拟合剖面。每一次迭代只能调节1 个深度，那么则需要有242（=35-1）种不同的组合（除去每个深度都变化零的组合），需要迭代242 次，

本研究使用6 个参数方程迭代242 次，直到得到最佳拟合结果。我们使用两种方法来判断迭代终止：第一种判断依据是如果RMSE＜Rc（标准误差）则终止迭代，如果RMSE＞Rc 则继续迭代，直到k=kmax；第二种判断标准是如果迭代到kmax的时候RMSE 仍然大于Rc，那么就去除这个剖面数据。标准差的选取标准是通过观测数据来确定的。

遥感SST 数据不能真实地表示现场观测[23]。为了避免这一情况产生相关的误差，文中采用线性回归的方法建立遥感SST 数据与Argo 观测SST 数据之间的联系，回归分析得到SST 和Argo 观测数据表面温度的回归曲线，确定算法中的初始值，计算公式为

由SST 和Argo 观察数据得到回归系数k、b，就可以得到二者的线性关系，即可得到重构三维温度场的初始条件。

表1 参数式中各变量名及其含义Table 1 Variable names and their meanings in the parameter formula

3 结果与分析

3.1 算法实验

本文研究区域是西北太平洋黑潮及黑潮延伸体海域，我们在该海域应用上文的算法。首先插值得到高分辨率的Argo 剖面数据；然后对每个剖面进行分段参数拟合，并计算垂直梯度；最后构建了三维温盐场并进行了验证。

根据2.2 节回归分析得到遥感SST 数据与Argo观测的SST 数据的回归曲线，确定算法中的初始值，回归结果如图5 所示。

图5 基于回归分析的遥感SST 与Argo 剖面最浅测量值之间的线性关系Fig.5 Linear relationship between remote sensing SST and the shallowest measured value of Argo profile based on regression analysis

由于Argo 现场观测的剖面数据某一部分会较为稀疏，可能会导致误差增大。因此，根据2.2 节所述，首先将Argo 观测资料插值设为Δz=5 m，则zi=zi+1-5 m，并且z0=0。如果我们提前知道（d1，d2，d3，d4，d5）5 个深度值，那么我们就可以把整个剖面分成6 部分（混合层，夹层，温跃层，过渡层，第一深层，第二深层）并得到6 个垂直梯度来进行分段参数拟合。实际上，处理数据之前我们并不知道5 个深度值，文中通过迭代方法找到最佳的深度值。首先，我们以一组深度数据和梯度数据的第一猜想值开始：

在本文中，我们选择第一猜想值如下：

每一次迭代计算后都要观察误差值的情况，误差值计算公式如下：

然后，我们按照2.2 节介绍的迭代方法进行迭代计算，通过迭代的方法得到最佳的拟合剖面。我们采取两种方式来终止迭代过程，细节见2.2.2 节，本文kmax=300，Rc=0.5。这个标准误差的选择是由观测资料决定的。由于迭代次数一般设置为算法本身能够承受的次数，此外还要看可行解的空间大小，空间大的则需要更多的迭代次数，所以一般将迭代次数设置多一点，然后设置结束误差，达到循序的误差跳出即可。因此文中设置kmax=300，并设置标准误差Rc=0.5作为终止迭代的判据。

3.2 验证与结果分析

参考模型可以把任意观测剖面转换为一系列参数：5 个深度（d1，d2，d3，d4，d5）和6 个垂直梯度这些参数代表了温度垂直剖面中重要的物理特征。并且实验结果显示它们对第一猜测值不敏感，较好的第一猜测值选择只是会减少迭代次数，因此该算法对于预处理工作的要求会有所降低。

拟合的结果的误差值变化趋势如图6a 所示，从图中可以看出，随着拟合次数的增加，拟合结果的RMSE 值减小，拟合精度提高。图6b 很好地反映了垂直温度随深度的变化，随着深度的增加海水温度逐渐降低，温度递减速率（温度梯度）在不同深度也不尽相同。在海洋表层的混合层的反演误差相较于深层更大，这表明影响表层海水和次表层海水的因素较多致使反演误差增大。

如前所述，对每个Argo 剖面采用分段拟合方法，得到一系列垂直温度梯度。图7 中黑潮延伸体不同月份剖面的拟合结果与Argo 观测结果的对比中，RMSE 分别为0.48 和0.53。因为三维温度场重构计算的是垂直温度梯度，所以拟合过程引入的最终温度场误差有限。此外，温度趋势在1 000 m 以下往往表现良好，也会使误差降低。

图8 是重构的4 个季节剖面结果，图中春季（2-4 月）重构结果与观测结果的温跃层梯度的误差为3.70%；夏季（5-7 月）重构结果与观测结果的温跃层梯度误差为17.91%；秋季（8-10 月）重构结果与观测结果的温跃层梯度误差为1.39%；冬季（11 月至翌年1 月）重构结果与观测结果的温跃层梯度误差为15.14%[12]。导致夏季和冬季拟合误差较大的原因可能是黑潮延伸体海域在夏冬两季辐射强迫增强或变弱时，混合层都会发生较大变化。图9 显示研究海域的温度水平分布，不论是在南北向还是东西向都存在一定的温度差异，所以水平温度的空间不均匀性可能也是导致夏季冬季误差增大的原因之一。从图8 中可以看出，黑潮延伸体海域夏季存在很强的跃层，这是因为表层增温，从而产生强跃层。冬季的混合层加深，则是因为表层温度降低，导致混合层向下扩展。由此可见拟合结果与观测结果得到的结论一致。

图6 拟合结果误差分析（a），拟合结果与观测数据反应的垂直温度随深度的变化（b）Fig.6 Error analysis of the fitting results (a),and the vertical temperature of the fitting result and the observed data changes with depth (b)

图7 通过参数方程拟合得到的黑潮延伸体不同月份（a.2018 年2 月; b.2018 年3 月）拟合剖面与Argo 观测剖面的对比Fig.7 Comparison of the Kuroshio extension body in different months (a.February 2018; b.March 2018) between the fitting profile and the Argo observation profile obtained through parametric equation fitting

图8 重构的4 个季节剖面结果Fig.8 Profile results reconstructed in 4 seasons

本文中用BOA-Argo 数据集以及EN4 数据集验证算法的有效性。将初始条件与垂直梯度相结合，利用参数化方程拟合的剖面，记重构剖面为TFR(Temperature fitting result)。BOA-Argo 数据库的垂直层为0～2 000×104Pa，BOA-Argo 数据库的性能证明了中尺度特征的准确性和持久性。EN4 月数据也是可靠的，因为许多类型的海洋剖面仪已经被使用。因此，我们使用BOA-Argo 和EN4 网格化产品来验证得到的TFR结果，其目标分析也较为可靠。

为了验证算法的有效性，本文分别计算分析了TFR 结果与BOA-Argo、EN4 数据的RMSE 值。图10显示了TFR 结果与BOA-Argo 和EN4 数据库的相关数据随深度的分布，与两种不同数据集的相关系数均为0.84 以上，RMSE 则分别为0.46 和0.37。为了更进一步确定算法的准确性，我们以BOA-Argo 为参考标准计算BOA-Argo 和EN4 两数据集的RMSE，得到的RMSE 为0.55，与TFR 结果与BOA-Argo 的误差相近，间接证明了算法的有效性，表明TFR 结果与两个现有的数据库具有较高的一致性。

图9 来自TFR、BOA-Argo 和EN4 数据集的表层、100×104 Pa 和500×104 Pa 深度的温度分布Fig.9 Temperature distributions at the surface,100×104 Pa,and 500×104 Pa depths from the TFR,BOA-Argo,and EN4 datasets

经过实验的对比分析，重构得到的剖面分辨率为0.25°× 0.25°（根据现有的观测资料，经过反复计算试验，确定的模型是具有价值的），因此更多的信息被反映在重建中。我们比较分析了表层，100×104Pa、500×104Pa 的水层的温度分布情况（图9）。图中可以看到，拟合结果很好地反映了海洋温度的分布特征，分辨率有所提高，能够捕获到尺度较小的信息，如黑潮延伸体主体部分水平温度存在南高北低、东低西高的特征。如图9 所示，BOA-Argo 数据集的分辨率最低，消除噪声时可能会把原始信息消除。

4 结论

本文提出了一种基于Argo 剖面资料以及卫星海表面温度数据重构的三维温度场的有效算法，通过分段拟合方法得到连续的垂直温度变化，然后将该算法以海表面作为开始，向海洋深层应用。

图10 TFR 结果与BOA-Argo（a），TFR 结果与EN4 产品（b），TFR 结果与BOA-Argo、EN4 产品的温度-水深散点图（c）Fig.10 Scatter plots between TFR results and BOA-Argo (a),TFR results and EN4 products (b),TFR results and BOA-Argo,EN4 products (c)

为了验证该算法的有效性，本文在黑潮区域进行了三维温度场重建的相关实验并进行验证。与EN4和BOA-Argo 数据相比较而言，重构后温度场的精度较高，分辨率也较高。对于现有的Argo 衍生产品而言，是把Argo 测量值直接融合，那么原始信息则有可能会被剔除；而BOA-Argo 产品是合并了Argo 数据来改变初始条件，提高了数据精度。然而海洋温度是一个非常复杂的变量，其变化情况及其分布规律也同样是一个复杂的问题。那么在缺乏现场观测数据的情况下，使用插值方法得到连续的温度剖面从而估计温度是不太准确的，在一定程度上也影响了反演的准确性。然而温度梯度相对于温度来说具有一定的滞后性，于是对温度梯度插值的方法就成为了较好的一种方法。此外，我们的算法只使用了原始的Argo 温度剖面数据，并直接应用于我们的算法中，没有与其他任何分析进行交互，因此我们的重建在获得较小尺度信号方面得到了很大的改善。

本文的算法仍有不足之处，需要在后续工作中加以改进。首先文中的参数量较大，参数较多会增加拟合误差，所以接下来我们需要减少参数量，提高拟合精度。其次，对于夏冬两季拟合误差较大这一情况，需要我们在后续工作中去改善我们的方法，以提高精度。例如，在算法中增加辐射强迫、风、热量以及动力因素等。

目前，三维温度重建算法已经在黑潮等相对较小的海域进行了验证，并且结果较好，未来该算法可以在更大的海域甚至全球海域得到验证并广泛应用。