一种基于海洋区域再分析资料的卫星红外载荷定标方法与流程

1.本发明属于遥感技术领域，尤其涉及一种基于海洋区域再分析资料的卫星红外载荷定标方法。


背景技术：

2.卫星载荷传感器的绝对辐射定标是遥感数据定量化应用的前提和关键，经过辐射定标后的遥感数据能够反映地物的特征及变化规律。绝对辐射定标指建立传感器观测影像像元亮度值与入瞳处表观辐亮度之间定量关系的过程。在卫星上天后，由于外太空辐射以及元器件发生老化，仪器的观测性能可能发生变化，影响遥感数据的真实性和准确性，所以需要开展卫星传感器的在轨绝对辐射定标。
3.对于红外载荷，常用的在轨辐射定标方法有星上定标、场地定标、交叉定标等。其中星上定标是指利用卫星载荷自身携带的辐射定标装置开展的定标工作，红外传感器星上定标一般以参考黑体和深冷空间为辐射定标参考。场地定标是利用地球表面大面积均匀的地物，例如红外定标常用水体，在卫星过境时，测量获取地表以及大气参数，利用辐射传输模型获得传感器入瞳等效辐亮度进行定标。交叉定标是以一个公认精度较高的传感器为参考，当待定标卫星和参考卫星观测同一均匀目标时，获取待定标传感器等效入瞳辐亮度进行定标系数计算。
4.上述三种红外载荷绝对辐射定标方法在实际定标中有着广泛应用，但也都存在自身的优势与局限性。星上定标方法不受时间、环境等客观因素的影响，能实时获取传感器定标系数，但星上定标系统受到外太空辐射影响，自身性能可能随时间发生衰减，影响定标精度。场地定标能够利用实测地表和大气参数，提供精确的绝对辐射定标系数，但该方法易受地理空间的限制，需要在卫星过境时开展同步观测，消耗大量的人力物力，定标成本较高，这也导致了定标周期相对较长、频次低。交叉定标成本低，无需建立地面校正场便可实现绝对辐射定标，但该方法定标精度取决于参考传感器自身观测精度，同时对于只能覆盖观测全球部分区域的卫星，与参考星交叉定标次数相对较少，不能达到较高的定标频次。
5.因此，如何在提高红外载荷定标频次的同时保证定标精度，是现有技术中需要解决的问题。


技术实现要素：

6.为解决上述问题，本发明提供一种基于海洋区域再分析资料的卫星红外载荷定标方法，精度高、成本低，能够实现卫星红外载荷的高频次定标。
7.一种基于海洋区域再分析资料的卫星红外载荷定标方法，包括以下步骤：
8.s1：基于卫星运行及红外载荷观测的地表规律，确定感兴趣的时间段与空间区域条件；
9.s2：将感兴趣的时间段与空间区域条件作为检索条件，从卫星红外载荷相应的遥感影像数据服务系统筛选出待定标的红外影像数据集；
10.s3：逐景从待定标的红外影像数据集中筛选出包含无云、均匀海洋区域的影像作为初步影像；同时，分别对各景初步影像进行格网划分，其中，格网的大小与再分析资料中的地表数据分辨率和大气数据分辨率的较小值相同；
11.s4：对每一景初步影像，根据初步影像中各像元与各格网的经纬度信息确定像元与网格之间的包含关系，并根据所述包含关系筛选出感兴趣格网，其中，感兴趣格网包含的像元均为无云、均匀海洋像元；
12.s5：获取各景初步影像包含的感兴趣格网在卫星红外载荷的红外波段影像中对应的亮度值；
13.s6：分别根据各景初步影像包含的感兴趣格网的经纬度信息与卫星成像时间，从再分析数据库中进行空间匹配和时间匹配，提取出对应空间范围、时间范围的地表数据集与大气参数数据集，并分别对各景初步影像对应的地表数据集与大气参数数据集进行空间聚合以及时间插值，得到与所述待定标的红外影像数据集相对应的再分析资料数据集；
14.s7：基于所述再分析资料数据集，利用大气辐射传输模型计算卫星红外载荷入瞳处红外波段的等效辐亮度l
toa
；
15.s8：采用最小二乘法分别对所有感兴趣格网在红外波段影像中对应的亮度值与等效辐亮度l
toa
进行拟合，得到红外波段的定标增益系数gain与定标偏移系数bias如下：
16.l
toa
＝gain
·
dn bias
17.其中，dn为红外波段下所有感兴趣格网的亮度值的集合。
18.进一步地，步骤s3中所述的初步影像的获取方法如下：
19.s31：根据各待定标的红外影像中各像元的经纬度信息与全球海陆框架经纬度数据集，分别判断各像元的海陆属性，得到各景待定标的红外影像对应的二值化海陆区域掩膜矩阵；
20.s32：将各待定标的红外影像与各自对应的二值化海陆区域掩膜矩阵相乘，得到海陆影像，同时，利用卫星近红外波段的观测数据分别计算各海陆影像中各像元的大气层顶表观反射率ρ
toa
；
21.s33：根据ρ
toa
对海陆影像中的海洋区域进行云层和海洋的分割，其中，对于任一像元，若其ρ
toa
小于设定阈值，则该像元代表海洋，反之，该像元代表云层，从而得到二值化无云海洋区域掩膜；
22.s34：将各待定标的红外影像与各自对应的二值化无云海洋区域掩膜矩阵相乘，得到无云海洋影像，分别将无云海洋影像中无云海洋区域的各个像元作为滑窗中心，计算各滑窗内所有像元亮度值的离散系数，将离散系数低于设定阈值的像元记为均匀像元，否则为非均匀像元，得到二值化无云、均匀海洋区域掩膜矩阵；
23.s35：对二值化无云、均匀海洋区域掩膜矩阵进行形态学膨胀，得到无云、均匀海洋区域膨胀掩膜；
24.s36：将各待定标的红外影像与各自对应的无云、均匀海洋区域膨胀掩膜相乘，得到初步影像。
25.进一步地，每一景初步影像中，各像元与各网格之间的包含关系分别由各像元的中心点的经纬度与各格网四个角点之间的经纬度决定。
26.进一步地，各感兴趣格网在卫星红外载荷的红外波段下的亮度值为：各感兴趣格
网所包含的所有像元在卫星红外载荷的红外波段下的亮度值均值。
27.进一步地，步骤s7所述卫星红外载荷入瞳处红外波段的等效辐亮度l
toa
的获取方法包括以下步骤：
28.s71：基于所述再分析资料数据集中的大气廓线参数以及成像时卫星观测几何参数计算卫星成像时目标海域-卫星路径的大气上行总透过率、大气上行辐射以及大气下行辐射；
29.s72：基于所述再分析资料数据集中的海表温度数据，利用黑体辐射定律计算热红外谱段黑体发射辐射b(λ,t)，计算公式如下：
30.b(λ,t)＝2
·h·
c2/(λ5·ehc/λkt-1)
31.其中，h为普朗克常数，c为光速，λ为波长，e为自然对数的底，k为玻尔兹曼常数，t为黑体温度；
32.s73：利用热红外波段卫星观测地表接收的单色辐亮度计算公式，获取卫星红外载荷的入瞳辐亮度l
toa
(λ)，计算公式如下：
33.l
toa
(λ)＝ρ
ocean
·
b(λ,t)
·
τ l
up
(λ) l
down
(λ)
·
(1-ρ
ocean
)
34.其中，ρ
ocean
为海水发射率，τ为大气上行总透过率，l
up
(λ)为大气上行辐射，l
down
(λ)为大气下行辐射；
35.s74：根据所述入瞳辐亮度l
toa
(λ)及卫星红外载荷的光谱响应函数，计算卫星红外载荷入瞳处红外波段的等效辐亮度l
toa
，计算公式为：
[0036][0037]
其中，s(λ)为卫星红外载荷在红外波段下的光谱响应函数，λ1与λ2分别为光谱响应函数中红外波段覆盖光谱范围的最小值、最大值。
[0038]
有益效果：
[0039]
1、本发明提供一种基于海洋区域再分析资料的卫星红外载荷定标方法，根据筛选提取的待定标红外载荷影像数据集，获取无云、均匀感兴趣海域对应的感兴趣格网的经纬度信息，基于此提取地表数据集与大气参数数据集等再分析资料，并进行空间插值、时间聚合的预处理；利用大气辐射传输模型，结合再分析海表、大气参数及卫星观测辅助信息，计算卫星观测入瞳等效辐亮度，与感兴趣格网的亮度值线性拟合计算定标系数，实现绝对辐射定标；由此可见，本发明将再分析资料数据作为辐射定标基准源，不需要开展实地的地表、大气参数测量工作，数据获取简便且低成本；同时，基于无云、均匀海面开展辐射定标工作，能够有效提高定标频次，总体来说，本发明可操作性强，优势明显且易于实现。
[0040]
2、本发明提供一种基于海洋区域再分析资料的卫星红外载荷定标方法，在空间上基于卫星所经过海洋区域的遥感观测数据开展定标，地球海洋区域覆盖范围大，在轨运行卫星在良好运行状态下，能够获取覆盖不同纬度带海洋区域的遥感数据，相应海洋表面温度所覆盖的动态范围能够满足定标需求；在海洋区域观测可用数据量足够时，理论上可实现每日或数天的高频次定标，对于全球覆盖观测以及区域观测卫星，相较于其他定标方法，在定标频次上均有较大的提高。
[0041]
3、本发明提供一种基于海洋区域再分析资料的卫星红外载荷定标方法，定标数据源选用基于同化模型的再分析资料数据，额外投入的人力和资源较少；也就是说，本发明是
对已有观测资料的充分利用，同时利用大数据量的观测资料，在统计意义上降低所选用数据的不确定性，保证了辐射定标精度，弥补了场地定标单次观测样本点少、资源消耗大的缺点。
[0042]
4、本发明提供一种基于海洋区域再分析资料的卫星红外载荷定标方法，有效提高定标时遥感影像数据的利用率，在遥感数据选择空间范围和时间跨度上有了较大的扩展；同时，本发明还挖掘了再分析数据在遥感定标领域的应用潜力，同时可以实现载荷红外观测影像的历史数据在统一定标基准下的长时间序列的高频次定标。
附图说明
[0043]
图1为本发明提供的一种基于海洋区域再分析资料的卫星红外载荷定标方法的流程图；
[0044]
图2为本发明提供的感兴趣格网的获取方法流程图；
[0045]
图3为本发明提供的等效辐亮度的获取方法流程图。
具体实施方式
[0046]
为了使本技术领域的人员更好地理解本技术方案，下面将结合本技术实施例中的附图，对本技术实施例中的技术方案进行清楚、完整地描述。
[0047]
本发明的目的是提供一种基于海洋区域再分析资料的卫星红外载荷定标方法，它使用再分析资料数据，低成本且能够实现高频次且精度较高的红外谱段辐射定标，同时该方法可以针对传感器红外载荷历史观测数据开展在统一定标基准下的长时间序列绝对辐射定标。
[0048]
具体的，如图1所示，一种基于海洋区域再分析资料的卫星红外载荷定标方法，包括以下步骤：
[0049]
s1：基于卫星运行及红外载荷观测的地表规律，确定感兴趣的时间段与空间区域条件。
[0050]
需要说明的是，在利用卫星待定标影像开展辐射定标工作时，需覆盖较大的辐射亮度观测范围，以此保证所计算定标系数的精确性、可靠性，同时，考虑到基于海表面开展红外载荷辐射定标提供的待定标数据量较大，需合理确定时空范围。因此，本发明首先基于卫星运行及红外载荷观测地表规律，对过境海洋区域的空间范围、时间跨度进行分析，为后续步骤中筛选出待定标的红外载荷观测数据集时提供筛选条件，其中，筛选条件即为感兴趣的时间段与空间区域条件。
[0051]
具体的，不同卫星红外载荷成像原理及设计、运行规律、观测的全球覆盖范围、星下点空间分辨率均有一定差异，需基于所选红外载荷成像特点，设定合适的数据统计时空范围条件。例如，对于能够观测覆盖全球范围的红外载荷，可在全球海域中选择相对较为稳定，覆盖高、中、低纬度的海洋区域，而在时间跨度上可以较小，可选择数天或一周作为数据统计的时间段；对于只能观测覆盖全球部分区域的红外载荷，可将所覆盖的全部海域范围进行统计，但在时间跨度上的统计条件需要增大，基于红外载荷能够覆盖的海域纬度范围，可选择月或季度作为数据统计的时间段。
[0052]
s2：将感兴趣的时间段与空间区域条件作为检索条件，从卫星红外载荷相应的遥
感影像数据服务系统筛选出待定标的红外影像数据集。
[0053]
需要说明的是，筛选出待定标的红外影像数据集后，可以剔除其中不满足定标应用需求的红外影像数据，如基本只覆盖陆地范围、海域完全被云覆盖、影像数据存在异常等红外影像，获得无异常的待定标红外影像数据集后再进行后续步骤。
[0054]
s3：逐景从待定标的红外影像数据集中筛选出包含无云、均匀海洋区域的影像作为初步影像；同时，分别对各景初步影像进行格网划分，其中，格网的大小与地表数据和大气数据中分辨率的较大值相同。
[0055]
进一步地，初步影像的获取方法如下：
[0056]
s31：根据各待定标的红外影像中各像元的经纬度信息与全球海陆框架经纬度数据集，分别判断各像元的海陆属性，得到各景待定标的红外影像对应的二值化海陆区域掩膜矩阵，其中，1代表海洋，0代表陆地；
[0057]
s32：将各待定标的红外影像与各自对应的二值化海陆区域掩膜矩阵相乘，得到海陆影像，同时，利用卫星近红外波段的观测数据分别计算各海陆影像中各像元的大气层顶表观反射率ρ
toa
；
[0058]
进一步地，大气层顶表观反射率ρ
toa
的计算公式如下：
[0059]
ρ
toa
＝π
·
l
·
d2/(esun
·
cosθ)
[0060]
其中，ρ
toa
为红外波段的大气表观反射率，d为天文单位日地距离，θ为太阳天顶角，esun为红外波段大气层顶的平均太阳光谱辐照度，l为入瞳辐亮度。
[0061]
基于卫星近红外波段观测的大气表观反射率，确定阈值条件进行海洋上空云覆盖区域剔除，获取影像无云海洋范围区域掩膜，具体参见步骤s33；
[0062]
s33：根据ρ
toa
对海陆影像中的海洋区域进行云层和海洋的分割，其中，对于任一像元，若其ρ
toa
小于设定阈值，如ρ
toa
＜0.1，则该像元代表海洋，反之，该像元代表云层，从而得到二值化无云海洋区域掩膜，其中，1代表海洋，0代表陆地和云层；
[0063]
s34：将各待定标的红外影像与各自对应的二值化无云海洋区域掩膜矩阵相乘，得到无云海洋影像，分别将无云海洋影像中无云海洋区域的各个像元作为滑窗中心，计算各滑窗内所有像元亮度值的离散系数，将离散系数低于设定阈值的像元记为均匀像元，否则为非均匀像元，得到二值化无云、均匀海洋区域掩膜矩阵，其中，1代表均匀海洋，0代表陆地、云层以及非均匀海洋；
[0064]
其中，离散系数的计算公式如下：
[0065]
cv＝σ/μ
[0066]
其中，cv为区域范围的离散系数，σ为区域范围内亮度值标准差，μ为区域范围内亮度值平均值。
[0067]
基于离散系数计算结果与设定的阈值条件，即可获取影像无云、均匀海洋区域掩膜。此外，滑动窗口的大小，可基于影像及再分析资料空间分辨率确定，而变异系数阈值条件，可基于其变化范围及所代表的滑窗区域内观测数据的均匀程度来确定，本发明不一一列举。
[0068]
s35：对二值化无云、均匀海洋区域掩膜矩阵进行形态学膨胀，得到无云、均匀海洋区域膨胀掩膜；
[0069]
其中，膨胀运算公式如下：
[0070][0071]
其中，x
′
为膨胀后影像无云、均匀海洋区域的掩膜，x为原影像无云、均匀海域范围区域初步掩膜，s为菱形结构元素；
[0072]
菱形结构元素的大小通过下式计算：
[0073]
n＝w/r
[0074]
其中，n为菱形结构元素的长/宽，w为设定的腐蚀宽度，r为影像分辨率；
[0075]
通过形态学膨胀，获得影像无云、均匀海洋区域的膨胀后掩膜。
[0076]
s36：将各待定标的红外影像与各自对应的无云、均匀海洋区域膨胀掩膜相乘，得到初步影像。
[0077]
s4：如图2所示，对每一景初步影像，根据初步影像中各像元与各格网的经纬度信息确定像元与网格之间的包含关系，并根据所述包含关系筛选出感兴趣格网，其中，感兴趣格网包含的像元均为无云、均匀海洋像元。
[0078]
其中，每一景初步影像中，各像元与各网格之间的包含关系分别由各像元的中心点的经纬度与各格网四个角点之间的经纬度决定，若某一像元中心点落入某一格网四个角点之中，则认为该像元被包含在该格网内；同时，初步影像中的每个像元分为两种，一种为无云、均匀海洋像元，另一种为非无云、均匀海洋的像元，如陆地属性像元或者云属性像元等。
[0079]
s5：获取各景初步影像包含的感兴趣格网在卫星红外载荷的红外波段影像中对应的亮度值，其中，亮度值为各感兴趣格网所包含的所有像元在卫星红外载荷的红外波段下的亮度值均值。
[0080]
需要说明的是，由于卫星红外影像分辨率与再分析资料的空间分辨率不同，需将其统一于相同的空间尺度。例如，可基于再分析资料的空间分辨率，提取观测影像在该感兴趣格网大小区域内的像元亮度值，通过均值计算代表卫星载荷在此空间分辨率下的观测亮度值。
[0081]
s6：分别根据各景初步影像包含的感兴趣格网的经纬度信息与卫星成像时间，从再分析数据库中进行空间匹配和时间匹配，提取出对应空间范围、时间范围的地表数据集与大气参数数据集，并分别对各景初步影像对应的地表数据集与大气参数数据集进行空间聚合以及时间插值，得到与所述待定标的红外影像数据集相对应的再分析资料数据集。
[0082]
需要说明的是，本发明所提取的大气参数数据集包括气压、大气温度、湿度、臭氧含量等大气廓线数据，此外，还可以将所选感兴趣海洋区域内卫星过境当日各参数数据均进行提取。
[0083]
进一步地，空间聚合计算公式如下：
[0084][0085]
其中，para
re-agg
为空间聚合后的再分析参数数据，n为聚合空间内再分析数据的数量，para
re
为提取的再分析参数数据。
[0086]
时间插值计算公式如下：
[0087]
para
re-inter
(t)＝(t-t1)
·
para
re-agg
(t0)/(t
0-t1) (t-t0)
·
para
re-agg
(t1)/(t
1-t0)
[0088]
其中，t为影像成像时间，t0为成像前再分析数据中时间点，t1为成像后再分析数据中时间点。
[0089]
s7：基于所述再分析资料数据集，利用大气辐射传输模型计算卫星红外载荷入瞳处红外波段的等效辐亮度l
toa
。
[0090]
进一步地，如图3所示，等效辐亮度l
toa
的获取方法包括以下步骤：
[0091]
s71：基于所述再分析资料数据集中的大气廓线参数以及成像时卫星观测几何参数计算卫星成像时目标海域-卫星路径的大气上行总透过率、大气上行辐射以及大气下行辐射；
[0092]
其中，大气上行总透过率指从海表到卫星红外载荷入瞳处大气的总透过率，该参数与卫星观测几何、大气状况等输入相关。大气上行辐射、大气下行辐射指大气自身在红外谱段的发射辐射。例如，可利用大气辐射传输模型modtran，基于所输入的再分析大气廓线数据、卫星观测几何及相关辅助信息进行模拟计算，输出文件中可提取波段相对应的大气上行总透过率及大气上行辐射、大气下行辐射。
[0093]
s72：基于所述再分析资料数据集中的海表温度数据，利用黑体辐射定律计算热红外谱段黑体发射辐射b(λ,t)，其中，基于再分析海表温度数据计算黑体红外谱段发射辐射是为后续计算海表在红外谱段的发射辐射做准备，黑体发射辐射计算公式如下：
[0094]
b(λ,t)＝2
·h·
c2/(λ5·ehc/λkt-1)
[0095]
其中，h为普朗克常数，c为光速，λ为波长，e为自然对数的底，k为玻尔兹曼常数，t为黑体温度；
[0096]
s73：利用热红外波段卫星观测地表接收的单色辐亮度计算公式，获取卫星红外载荷的入瞳辐亮度l
toa
(λ)；
[0097]
其中，卫星红外载荷的入瞳观测的辐射来源可分别三部分：第一部分来源于海面的热红外发射辐射，与海表面温度、海面发射率以及海面与载荷之间的大气透过率相关；第二部分来源于大气自身的上行热红外发射辐射，即大气上行辐射，直接进入载荷传感器；第三部分来源于大气自身的下行热红外发射辐射，与大气下行辐射、海表面发射率大气透过率相关；本发明入瞳辐亮度l
toa
(λ)的计算公式如下：
[0098]
l
toa
(λ)＝ρ
ocean
·
b(λ,t)
·
τ l
up
(λ) l
down
(λ)
·
(1-ρ
ocean
)
[0099]
其中，ρ
ocean
为海水发射率，τ为大气上行总透过率，l
up
(λ)为大气上行辐射，l
down
(λ)为大气下行辐射；
[0100]
s74：根据所述入瞳辐亮度l
toa
(λ)及卫星红外载荷的光谱响应函数，计算卫星红外载荷入瞳处红外波段的等效辐亮度l
toa
；
[0101]
其中，由于传感器硬件的限制，载荷传感器在所设定波长范围内对每个波长并不是100％响应，传感器对不同波长的相对响应强度构成了光谱响应曲线，一般来说均为单峰函数。为将基于大气辐射传输模型计算的入瞳辐亮度与卫星观测值联系起来，需要使用光谱响应函数与步骤s73计算的入瞳辐亮度l
toa
(λ)进行卷积运算，计算整个热红外响应波段的入瞳等效辐亮度，具体计算公式为：
[0102][0103]
其中，s(λ)为卫星红外载荷在红外波段下的光谱响应函数，λ1与λ2分别为光谱响应
函数中红外波段覆盖光谱范围的最小值、最大值。
[0104]
s8：采用最小二乘法分别对所有感兴趣格网在红外波段影像中对应的亮度值与等效辐亮度l
toa
进行拟合，得到红外波段的定标增益系数gain与定标偏移系数bias如下：
[0105]
l
toa
＝gain
·
dn bias
[0106]
其中，dn为红外波段下所有感兴趣格网的亮度值的集合。
[0107]
也就是说，本发明结合所述规则格网经纬度信息，统计所设定时空范围内，卫星红外载荷观测各景影像数据感兴趣海域的区域亮度值，将其与所述卫星待定标影像的观测入瞳等效辐亮度一一对应，采用最小二乘法线性拟合出该对应关系，也即确定出红外载荷的定标系数gain与bias，最终实现卫星红外载荷的辐射定标。
[0108]
综上所述，本发明提供一种基于海洋区域再分析资料的卫星红外载荷定标方法，首先确定待定标红外载荷观测数据集统计的时空范围；获取待定标传感器海洋区域过境红外遥感影像；基于卫星遥感影像，筛选获得无云、均匀的感兴趣海域包含的格网经纬度信息；提取再分析数据库中与感兴趣格网的经纬度信息空间相匹配、与卫星成像时间相匹配的海表温度、大气廓线数据，并对其进行空间聚合及时间插值预处理；通过影像头文件中的观测几何参量、热红外通道光谱响应函数及预处理后的再分析资料数据，利用大气辐射传输模型计算待定标传感器热红外波段入瞳等效辐亮度；根据获取的热红外波段入瞳等效辐亮度，与待定标影像中感兴趣格网的亮度值进行线性拟合，计算传感器热红外波段定标系数，实现辐射定标。由此可见，本发明基于海洋区域再分析资料进行辐射定标，精度高、成本低，能够实现高频次定标，同时可以针对传感器红外载荷历史观测数据进行定标，它在遥感数据处理业务化应用领域具有推广前景。
[0109]
当然，本发明还可有其他多种实施例，在不背离本发明精神及其实质的情况下，熟悉本领域的技术人员当然可根据本发明作出各种相应的改变和变形，但这些相应的改变和变形都应属于本发明所附的权利要求的保护范围。