评估微波辐射计亮温数据长时间序列稳定性的方法

1.本发明属于被动微波遥感技术领域，特别涉及评估微波辐射计亮温数据长时间序列稳定性的方法。


背景技术：

2.微波辐射计亮温数据可以用来反演海表面温度、风速、大气水汽含量、云液水含量等地球物理参数。其中海表面温度是影响全球气候变化与长期天气过程的主要因素之一，对全球水循环和全球表面的能量收支平衡具有重要作用。获取长时间序列的海表面温度变化对于研究全球气候变化具有重要意义。而获取正确的海表面温度变化的前提是用来反演海表面温度的亮温数据是稳定无漂移的，这样才能保证海表面温度的变化是实际的变化而不是由于亮温的漂移带来的变化。
3.微波辐射计是一个获取地表辐射特性的高灵敏度的被动遥感载荷，其成像过程复杂，任何一个环节的变化都有可能引起辐射计测量亮温的漂移或者错误。比如传感器的老化、卫星姿态的变化以及卫星轨道引起的入射角的变化等。还有一种正常的变化或者漂移是由实际地球物理参数的变化引起的。所以辐射计亮温的稳定性评估就是获取非地球物理参数变化引起的亮温漂移。
4.由于微波辐射计的亮温漂移往往比较复杂，是多种因素的共同作用引起的变化，传统的分析和计算手段，很难定量每一个引起亮温漂移的因素的贡献量，最终也很难评估出亮温的漂移量。


技术实现要素：

5.针对微波辐射计的亮温漂移的影响因素多，传统方法难以较为准确地计算出亮温的漂移量的上述问题，本发明提供一种评估微波辐射计亮温数据长时间序列稳定性的方法，根据亮温自身特性来确定亮温的非正常变化，从而计算实际的亮温漂移量，避免了从原因出发去计算漂移量的困难和复杂，能够高效和准确地评估微波辐射计亮温的长期稳定性。
6.为了实现上述目的，本发明提供的评估微波辐射计亮温数据长时间序列稳定性的方法，包括如下步骤：
7.步骤(1)确定一个时间尺度，以所述时间尺度读取微波辐射计亮温观测数据，选取无陆地、无海冰的海上亮温观测数据；
8.步骤(2)使用步骤(1)中选取的时间尺度的微波辐射计亮温数据组成频率直方图，根据该直方图确定冷亮温集搜索的亮温上限阈值；
9.步骤(3)根据步骤(2)得到的亮温上限阈值，在步骤(1)中得到的整月亮温中，搜索得到各通道的冷亮温集；
10.步骤(4)根据步骤(3)中得到的各通道冷亮温集，计算直方图中小于n％对应的亮温平均值，n为区间[3,10]，区间内相邻数值的间隔为0.2，得到包含36个元素的向量，每个
元素对应一个百分比；所述小于n％对应的亮温平均值的计算方法是，先对所有亮温进行排序，再用亮温总数乘n％得到待平均的亮温个数，在排序好的亮温中从小到大取这些个数的亮温进行平均，即为小于n％对应的亮温平均值；
[0011]
步骤(5)利用步骤(4)中得到的向量进行多项式拟合，外推得到0％对应的亮温，即为每个月对应的最冷亮温；
[0012]
步骤(6)根据步骤(5)得到的每个时间尺度的最冷亮温，进行线性拟合，得到直线斜率，即可等效计算出亮温漂移速度。
[0013]
进一步的，所述步骤(1)时间尺度优选为月。
[0014]
进一步的，微波辐射计优选多频微波辐射计。
[0015]
进一步的，所述步骤(2)确定搜索冷亮温集的阈值上限的具体方法是根据直方图第一个波峰的1/3对应的亮温确定，公式如下：
[0016]
tb
上限阈值
＝1/3
·
tb
第一个波峰
[0017]
tb
第一个波峰
为直方图第一个波峰对应的亮温；tb
上限阈值
为搜索冷亮温集的上限阈值。
[0018]
进一步的，所述步骤(4)计算小于n％对应的亮温平均值之前，先对冷亮温集进行
±
3σ控制，排除可能出现的异常值。
[0019]
进一步的，所述步骤(5)可以进行一次多项式拟合、二次多项式拟合，优选三次多项式拟合；利用步骤(4)得到的包含36个元素的向量进行3次拟合外推计算，公式如下：
[0020]
tb(f)＝a0 a1f a2f2 a3f3[0021]
f代表3％到10％以0.2％为间隔的36个百分比，tb(f)为对应的36个平均亮温，利用最小二乘可求得系数a0～a3，根据得到的a0～a3，将f＝0％代入上式，即可得到0％对应的亮温，即待求的最冷亮温。
[0022]
进一步的，所述步骤(6)对每个月的最冷亮温进行线性拟合，公式如下：
[0023]
y＝slope
·
x bias
[0024]
χ表示长时间序列的时间，χ的值大于等于2，y表示每个时间对应的最冷亮温，slope表示长时间序列亮温拟合直线的斜率，bias表示长时间序列亮温拟合直线的偏置，拟合出slope和bias值后，其中根据slope计算等效的年漂移速度drift
year
：drift
year
＝slope.12。
[0025]
由于微波辐射计的亮温漂移往往比较复杂，是多种因素的共同作用引起的变化，如果从引起漂移的原因去分析和计算漂移量，是非常困难的，很难定量每一个引起亮温漂移的因素的贡献量，最终也很难评估出亮温的漂移量。本发明从亮温本身出发，根据亮温自身特性来确定亮温的非正常变化，这样不但可以计算实际的亮温漂移量，还巧妙地避免了从原因出发去计算的困难。本发明就是基于这样一种思路，根据一定时间尺度内的全球海上亮温数据的直方图分布特性，找到这一时间尺度内的最冷亮温，最冷亮温往往存在于低海温、低风速、低大气水汽的海洋区域，理论上最冷亮温并不存在，所以通过搜索一个冷亮温集，利用外推的方法拟合得到最冷亮温，再根据最冷亮温不变的特性，可得到长时间序列的亮温漂移特性。
[0026]
本发明的有益效果如下：
[0027]
1.本发明从微波辐射计长时间序列的亮温自身特性出发去评估亮温的漂移状况，避免了从引起亮温漂移的原因出发的困难和复杂；
[0028]
2.本发明提出的最冷亮温集的上限阈值确定方法，对后续步骤中不同百分比平均亮温计算和外推最冷亮温的结果影响较小，即提高了该方法整体的鲁棒性；
[0029]
3.本发明可进行任意时长的微波辐射计亮温数据的稳定性评估，并可以给出定量化的年度漂移速度。
附图说明
[0030]
图1为本发明实施例评估微波辐射计亮温数据长时间序列稳定性的方法的流程图；
[0031]
图2为本发明实施例smr各通道全部亮温的直方图以及各通道小于阈值上限的亮温直方图；
[0032]
图3为本发明实施例计算不同百分比的平均亮温的示意图；
[0033]
图4为本发明实施例10.7ghz-h通道亮温根据3％-10％的对应平均亮温外推到0％的示意图；
[0034]
图5为本发明实施例smr的10.7ghz-h通道在2018年11月到2021年10月的亮温漂移速度评估结果。
具体实施方式
[0035]
下面结合具体实施例对本发明做进一步详细说明。
[0036]
实施例
[0037]
实施例中选用了hy-2b卫星的扫描微波辐射计(smr)多频亮温数据，数据可由国家海洋卫星数据分发中心获取(https://osdds.nsoas.org.cn/#/)。hy-2b卫星于2018年10月发射，其上搭载有包括smr在内的多个载荷。本实施例使用本发明所述的技术方法评估了hy-2b卫星发射以来三年(2018年11月-2021年10月)的smr亮温稳定性状况。
[0038]
本技术技术方法的理论依据是：辐射计的地球观测本身提供了一种自然的、可靠的定标标准，据此可以判断仪器漂移，由星载微波辐射计观测到的最冷亮度温度(tb)通常是很容易建模的，而且随着时间的推移是可重复的。也就是说在一定时间尺度内，全球海洋亮温组成的频率直方图具有稳定的下界，也就是低亮温端是稳定的，选择一定时长的海洋亮温为单位，获得每个时间单位内的直方图亮温下界，就可以得到这些单位在长时间序列的亮温趋势，即总体亮温的漂移信息。
[0039]
如图1所示，以月为单位，基于数据质量标识、陆地标识、海冰标识，获取smr的各通道海面观测亮温数据；将每个通道的全部亮温组成频率直方图；根据频率直方图确定一个亮温上限阈值，选取所有小于该阈值的亮温，形成冷亮温集；利用冷亮温集中小于n％(3≤n≤10，间隔为0.2)的亮温的平均值组成亮温序列，利用这个亮温序列使用3次多项式拟合外推0％对应的亮温，即最冷亮温。得到每个月的最冷亮温就可以利用直线拟合，得到长时间序列的亮温漂移速度。
[0040]
下面对本技术技术方法在评估smr亮温稳定性的实施例中的具体工作过程做详细说明。
[0041]
步骤(1)以月为单位，读取smr亮温观测数据，根据数据质量标识、陆地标识、海冰标识等，选取无陆地、无海冰的亮温观测数据。
[0042]
步骤(2)使用步骤(1)中选取的整月的smr亮温数据组成频率直方图，根据该直方图确定下一步冷亮温集搜索的上限阈值。具体方法是根据直方图第一个波峰的1/3对应的亮温确实，公式如下：
[0043]
tb
上限阈值
＝1/3
·
tb
第一个波峰
[0044]
tb
第一个波峰
为直方图第一个波峰对应的亮温；tb
上限阈值
为搜索冷亮温集的上限阈值。
[0045]
步骤(3)根据步骤(2)中的各通道冷亮温集搜索的上限阈值，在步骤(1)中得到的smr整月亮温中，搜索得到各通道的冷亮温集。
[0046]
图2所示为各通道全部亮温的直方图以及各通道小于阈值上限的部分(粗线表示)，也就是步骤(2)中得到的冷亮温集部分。
[0047]
步骤(4)根据步骤(3)中得到的各通道冷亮温集，计算小于n％的亮温对应的平均值，n为3到10，间隔为0.2，即可得到包含36个元素的向量。
[0048]
首先，对步骤3中得到的冷亮温集进行3-σ控制，排除可能出现的异常值。再计算直方图中小于n％的所有亮温的平均值。例如n为10，即先对所有亮温进行排序，再用亮温总数乘10％得到待平均的亮温个数k，在排序好的亮温中从小到大取k个亮温进行平均，即为小于10％对应的平均亮温。以0.2％为间隔，得到小于3％到10％的平均亮温，即包含36个元素的向量，每个元素对应一个百分比。如图3所示，计算不同百分比的平均亮温，如10％即计算黑色虚线内包含的全部亮温的平均值。表1展示了以10.7ghz-h通道为例，不同百分比对应的平均亮温的计算结果。
[0049]
表1 不同百分比及其对应的平均亮温(以10.7ghz-h为例)
[0050][0051][0052]
步骤(5)利用步骤(4)中得到的包含36个元素(如表1所示)的向量进行3次多项式拟合，外推得到0％对应的亮温，即为smr每个月对应的最冷亮温。
[0053]
拟合公式如下：
[0054]
tb(f)＝a0 a1f a2f2 a3f3[0055]
f代表小于3％到10％的36个百分比，tb(f)为对应的36个平均亮温，利用最小二乘可求得系数a0～a3，根据得到的a0～a3，将f＝0％代入上式，即可得到0％对应的亮温，即待求的最冷亮温。
[0056]
图4所示为外推到0％的示意图。
[0057]
步骤(6)根据步骤(5)得到的smr每个月的最冷亮温，进行线性拟合，得到直线斜率，即可等效计算出smr在长时间序列的亮温漂移速度。
[0058]
线性拟合公式如下：
[0059]
y＝slope
·
x bias
[0060]
χ表示长时间序列的时间，χ的值大于等于2，y表示每个时间对应的最冷亮温，slope表示长时间序列亮温拟合直线的斜率，bias表示长时间序列亮温拟合直线的偏置；拟合出slope和bias值后，其中根据slope，可以计算等效的年度漂移速度：
[0061]
drift
year
＝slope
·
12
[0062]
drift
year
即为年漂移速度，单位为k/year。
[0063]
图5为smr的10.7ghz-h通道在2018年11月到2021年10月的亮温漂移速度评估结果。
[0064]
以上所述的实施例仅是对本发明的优选实施方式进行描述，并非对本发明的范围进行限定，在不脱离本发明设计精神的前提下，本领域普通技术人员对本发明的技术方案作出的各种变形和改进，均应落入本发明权利要求书确定的保护范围内。