一种晨昏轨道卫星微光通道杂散光去除方法及装置与流程

1.本发明涉及卫星气象学技术领域，具体涉及一种晨昏轨道卫星微光通道杂散 光去除方法及其装置。


背景技术：

2.微光遥感对于观测与理解气候系统有着重要意义，它可用于天气云图分析， 探测海冰的范围、海水的浑浊程度，检测火点、火山灰烟羽和火山碎屑流等。相 比于其它卫星数据，微光数据的独特之处在于可以直接反映人类活动，比如识别 船舶以打击非法捕鱼，识别天然气燃除以评估温室气体排放，以及通过夜间灯光 研究城市相关问题，因此微光遥感逐渐成为遥感领域的研究热点。
3.目前国际上主要的微光探测数据来自美国snpp和noaa-20卫星搭载的可 见红外成像辐射计套件(visible infrared imaging radiometer suite，viirs)的白 天/夜间通道(day/night band，dnb)通道。viirs/dnb有16个探元，空间分 辨率750m，幅宽约3000km，具有低、中、高三增益阶段，分别在白天太阳光照 条件、晨昏暮曙光条件和夜间低照度条件下成像，从而实现7个量级的大动态 范围辐亮度观测。低增益阶段可进行星上定标，通过定标传递实现中、高增益阶 段的定标。snpp和noaa-20都是下午星，降轨过境地方时约凌晨1：30，在 南北极跨越晨昏线。viirs/dnb在晨昏线附近受杂散光污染严重，使图像质量 大幅下降。所谓杂散光，泛指卫星接收到的非目标物发射或反射的辐射，其产生 原因有以下几种可能：
4.1、由于镜面粗糙或不干净造成光学镜面散射；
5.2、遮光板散射的光进入探测器的光学路径；
6.3、仪器内部结构的孔洞或缝隙造成的漏光。
7.以viirs/dnb 2012年6月10日的一段数据为例(如图1a所示)，杂散光 造成图像上部一层灰白的蒙版。viirs图像产品研发团队的研究发现dnb杂散 光与以下几个因素相关：探元序号、扫描角度(即像元序号)、太阳光相对于星 下法向量的角度(简称太阳入射角)、太阳相对于星下法向量的方位角(简称太 阳方位角)和半角镜的不同面(半角镜有a、b两面)。基于此，他们杂散光去 除方案如下：
8.(1)筛选取黑暗像元：去除太阳天顶角或月球天顶角小于105
°
的像元，另 外使用2003年世界夜间灯光数据集去除存在灯光的像元。
9.(2)将黑暗像元分组：每个探元的每个半角镜面的数据都按扫描角度分组， 一条扫描线上4064个像元分成127组，每组32个像元，为减小噪声取，这32 个像元辐亮度的中位数，同时记录平均太阳光入射角。
10.(3)拟合太阳光入射角与辐亮度的关系：以2021年6月18日11时半角 镜b面第108组像元(在扫描线的右端)为例，图1b展示了16个探元所接收 的辐亮度与太阳入射角余弦(cosθ)的关系。cosθ＝0.27对应晨昏区，随cosθ增 加，太阳光直射入探测器导致辐亮度增强，cosθ＝0.45太阳光进入地球半影区， 辐亮度迅速降低到一个基本稳定的背景值，辐亮
度的背景值由气辉造成，因此将 卫星接收辐亮度减背景值即得到杂散光的强度。对杂散光进行分段拟合，在辐亮 度上升区间用二次多项式拟合辐亮度与cosθ的关系，地球半影区用四次多项式 拟合，在最左端辐亮度减弱区间需要通过外插得到。
11.(4)以上三步得到每个探元在不同半角镜面状态不同扫描角情景下杂散光 与太阳光入射角的关系，由此每个像元都可以通过探元序号、半角镜面状态、像 元序号和太阳光入射角估算出杂散光大小，之后使用原始辐亮度减去杂散光进 行去除。使用各月份的数据重复步骤(1)～(3)用于该月的杂散光去除，由此 解决太阳方位角改变造成影响。
12.2021年7月5日我国第二代极轨气象卫星风云三号e星(fy-3e)成功发 射，它是风云卫星家族里首颗晨昏轨道卫星，也是世界上首科晨昏轨道业务气象 卫星，其降轨过境地方时约为5：40，升轨约为17：40，与上午、下午卫星组网 实现了每6小时全球资料完全覆盖，满足了数值天气预报6小时同化时间窗的 要求。fy-3e的有效载荷中分辨率成像仪(mersi)配有微光通道(low lightband,llb)，星下点分辨率1km，幅宽2500km。llb具有低、中、高三增益阶 段，可观测3.5
×
10-5～50w/m2/sr的辐亮度，使风云系列卫星首次具有夜间微光 的探测能力。fy-3e始终处于晨昏线附近，mersi在设计时考虑晨昏轨道的特 殊性，为了避免杂散光影响增加了遮光板面积，减小了扫描镜向白天一侧的刈幅， 然而出乎意料的是llb在夜间一侧的数据受杂散光污染严重。初步分析发现 mersi/llb杂散光特征与viirs/dnb的差异很大，具体体现在：同一条扫描线 上杂散光随扫描角的变化大；杂散光随太阳入射角的变化复杂；杂散光强度存在 周期性的帧间变化。将viirs/dnb的杂散光去除方案施用于mersi/llb要求 步骤(2)中减小每组像元个数，以及在步骤(3)中做更多的分段来拟合辐亮度 和太阳入射角的关系，这极大的增加了计算量，而且最重要的是不能去除 mersi/llb帧间变化的杂散光。


技术实现要素：

13.由于mersi/llb的杂散光随扫描角和太阳入射角变化复杂且存在周期性 的帧间变化，viirs/dnb的杂散光去除方案对mersi/llb的施用效果欠佳。本 申请实施例的目的在于针对mersi/llb杂散光的特点，提供一种杂散光去除方 法，能够较好的去除mersi/llb的杂散光。
14.为实现上述目的，第一方面，本发明实施例提供了一种晨昏轨道卫星微光通 道杂散光去除方法，包括：
15.获取mersi/llb原始图像；所述原始图像仅包括纯夜间数据；
16.采用杂散光去除算法对所述待处理图像进行处理，以得到目标图像；
17.其中，所述杂散光去除算法包括如下步骤：
18.对所述原始图像进行逐扫描线滑动拟合，设某条扫描线的幅度亮度值为r1， 拟合辐亮度值为f1；
19.根据幅度亮度值r1和拟合幅度亮度f1的差异提取连续灯光像元的首尾像 元的像元号和辐亮度；
20.对连续灯光像元的辐亮度进行线性拟合，得到拟合结果，并用所述拟合结果 替代对应像元的幅度亮度值r1数值，得到新幅度亮度值r2；
21.对所述新幅度亮度值r2进行滑动拟合，得到拟合辐亮度f2；
22.根据幅度亮度值r1、拟合辐亮度f2以及预设的背景辐亮度数值，得到去 除杂散光后的目标图像。
23.在本技术的某些优选实施方式中，获取mersi/llb原始图像之前，所述方 法还包括：
24.读取mersi/llb的太阳天顶角数据，并从中选取纯夜间数据；所述纯夜间 数据为太阳天顶角大于101
°
的数据。
25.作为本技术的一种具体实施方式，对幅度亮度值r1进行拟合，得到拟合辐 亮度f1，具体为：
26.采用savitzky-golay滤波器，以61像元为窗口，三次多项式对幅度亮度值 r1做滑动拟合，得到拟合辐亮度f1；或
27.采用滑动平均、傅里叶滤波或经验模态分解对幅度亮度值r1做滑动拟合， 得到拟合辐亮度f1。
28.作为本技术的一种具体实施方式，根据幅度亮度值r1和拟合幅度亮度f1 提取连续灯光像元的首尾像元的像元号和辐亮度，具体为：
29.计算所述拟合幅度亮度f1与幅度亮度值r1的差值，取该差值的绝对值后 进行11点滑动平均，得到滑动平均结果；
30.采用设定阈值对所述滑动平均结果进行判断，对连续大于所述阈值的像元 判定为连续灯光像元，并提取首尾像元的像元号和辐亮度。
31.作为本技术的一种具体实施方式，对连续灯光像元的辐亮度进行线性拟合， 得到拟合结果，具体为：
32.以提取的首尾像元的像元号为自变量，其辐亮度为因变量做线性拟合，拟合 首尾像元之间像元的辐亮度，得到拟合结果。
33.作为本技术的一种具体实施方式，对所述新幅度亮度值r2进行拟合，得到 拟合辐亮度f2，具体为：
34.采用savitzky-golay滤波器对所述新幅度亮度值r2进行拟合，得到拟合辐 亮度f2；或
35.采用滑动平均、傅里叶滤波或经验模态分解对所述新幅度亮度值r2做滑动 拟合，得到拟合辐亮度f2。
36.作为本技术的一种具体实施方式，根据幅度亮度值r1、拟合辐亮度f2以 及预设的背景辐亮度数值，得到去除杂散光后的目标图像，具体为：
37.计算幅度亮度值r1与拟合辐亮度f2的差值，再加上气辉造成的背景辐亮 度数值，即得到去除杂散光的结果；
38.其中，气辉造成的背景辐亮度数值由mersi/llb在6～7月南极地区无杂散 的夜间观测而定。
39.第二方面，本技术实施例还提供了一种晨昏轨道卫星微光通道杂散光去除 装置，包括处理器、输入设备、输出设备和存储器，所述处理器、输入设备、输 出设备和存储器相互连接，其中，所述存储器用于存储计算机程序，所述计算机 程序包括程序指令，所述处理器被配置用于调用所述程序指令，执行上述第一方 面的方法。
40.第三方面，本技术实施例还提供了另一种晨昏轨道卫星微光通道杂散光去 除装
置，包括：
41.获取单元，用于获取mersi/llb原始图像；所述原始图像仅包括纯夜间数 据；
42.处理单元，用于采用杂散光去除算法对所述原始图像进行处理，以得到目标 图像；
43.其中，所述杂散光去除算法包括如下步骤：
44.对所述原始图像逐扫描线进行滑动拟合，设其中某条扫描线的幅度亮度值 为r1，拟合辐亮度为f1；
45.根据幅度亮度值r1和拟合幅度亮度f1的差异提取连续灯光像元的首尾像 元的像元号和辐亮度；
46.对连续灯光像元的辐亮度进行线性拟合，得到拟合结果，并用所述拟合结果 替代对应像元的幅度亮度值r1数值，得到新幅度亮度值r2；
47.对所述新幅度亮度值r2进行拟合，得到拟合辐亮度f2；
48.根据幅度亮度值r1、拟合辐亮度f2以及预设的背景辐亮度数值，得到去 除杂散光后的目标图像。
49.实施本技术实施例，具有如下优点：
50.1、本技术相比于原有的杂散光去除方案更简单易行，不需要提前计算用于 估算杂散光的一系列参数。
51.2、原有杂散光去除方案每个月更新一次估算杂散光的参数，由于太阳方位 角在一个月内是不断变化的，因此所估算的杂散光的准确性会有变化。而本技术 的方案不受季节变化的影响，理论上去除杂散光的准确性是不变的。
52.3、本技术的方案能够去除帧间变化的杂散光，而原有方案不能实现。
附图说明
53.为了更清楚地说明本发明具体实施方式或现有技术中的技术方案，下面将对 具体实施方式或现有技术描述中所需要使用的附图作简单地介绍。
54.图1a是2012年6月10日viirs/dnb的数据图像；
55.图1b是2021年6月18日11时半角镜b面第108组像元不同探元(不同 颜色细线表示)辐射强度与太阳入射角余弦的关系；
56.图2a是2022年1月2日22：15utc经过北京的一条扫描线的辐亮度(散 点)与使用savitzky-golay滤波器对辐亮度的滑动拟合值(红线)的示意图；
57.图2b是图2a中辐亮度与拟合值的差的绝对值(散点)与其十一点滑动平均 值(红线)的示意图；
58.图2c是去除灯光像元后的辐亮度(散点)与其滑动拟合值(红线)的示意 图；
59.图2d是去除基础辐亮度后的辐亮度的示意图；
60.图3是本技术实施例提供的晨昏轨道卫星微光通道杂散光去除方法的流程 图；
61.图4a是2022年1月2日22：15utc mersi/llb的原始图像，等值线代 表太阳高度角；
62.图4b是图4a去除杂散光之后的图像，只展示纯夜间区域。
63.图5是本技术实施例提供的晨昏轨道卫星微光通道杂散光去除装置的结构 示意图。
具体实施方式
64.下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清 楚、完整地描述，显然，所描述的实施例是本发明一部分实施例，而不是全部的 实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前 提下所获得的所有其他实施例，都属于本发明保护的范围。
65.应当理解，当在本说明书和所附权利要求书中使用时，术语“包括”和“包 含”指示所描述特征、整体、步骤、操作、元素和/或组件的存在，但并不排除 一个或多个其它特征、整体、步骤、操作、元素、组件和/或其集合的存在或添 加。
66.为更好地理解本发明实施例，下面以2022年1月2日22：15utc经过北 京的一条扫描线的纯夜间数据为例解释我们的方法原理：
67.图2a中的横坐标为该扫描线上像元的序号，太阳天顶角随序号增加而减小， 蓝色散点为该扫描线上辐亮度，单位w/m2/sr，可见辐亮度存在随太阳天顶角增 加而增加的趋势。这一增长趋势主要是受杂散光影响，是本实施例方法中所要要 去除的对象。
68.另外，图2a中第200～270、630～700、770～840像元受城市灯光影响，辐亮 度值显著高于周围像元。据此可以将辐亮度拆分为两部分，一是基础辐亮度，二 是城市灯光的辐亮度。其中基础辐亮度还可以拆分成两部分，一部分是存在趋势 变化的杂散光，另一部分是由气辉造成的恒定背景值。
69.为了得到基础辐亮度的变化趋势，我们使用savitzky-golay滤波器做滑动拟 合(图2a中红线)，然而在灯光像元处拟合结果会明显偏离基础辐射值，这就需 要去除灯光再做拟合。灯光像元的辐亮度与拟合值的差异较大(图2b中散点)， 可以作为判定灯光像元的依据。为了消除城市内个别暗像元的影响，我们对二者 的差异的绝对值做滑动平均(图2b中曲线)，大于阈值的像元判断为灯光像元。
70.选出灯光像元后，对其起始和末尾的像元以像元序号为自变量，基础辐亮度 值为因变量做线性拟合，用拟合的辐亮度值替代原辐亮度值(图2c中散点)从 而去除灯光的辐亮度。
71.之后再次使用savitzky-golay滤波器对去除灯光的辐亮度做滑动拟合(图2c 中曲线)，得到平滑的基础辐亮度，再将原始辐亮度值与基础辐亮度相减(图2d)， 如此去除了杂散光和背景辐亮度。
72.最后，统一加一个代表气辉辐亮度的定常值，就得到了去除杂散光的辐亮度。 其中，气辉的辐亮度值由mersi/llb在6～7月南极地区无杂散的夜间观测而 定。
73.基于上述算法原理，本技术提供一种晨昏轨道卫星微光通道杂散光去除方 法，主要针对llb在夜间一侧的数据严重受到杂散光污染的问题进行处理。如 图3所示，该杂散光去除方法主要包括如下步骤：
74.获取mersi/llb原始图像；所述原始图像仅包括纯夜间数据；
75.采用杂散光去除算法对所述原始图像进行处理，以得到目标图像；
76.其中，所述杂散光去除算法包括如下步骤：
77.对所述原始图像逐扫描线进行滑动拟合，设其中某条扫描线的幅度亮度值 为r1，拟合辐亮度为f1；即对每一条扫描线进行滑动拟合；
78.根据幅度亮度值r1和拟合幅度亮度f1的差异提取连续灯光像元的首尾像 元的像
元号和辐亮度；
79.对连续灯光像元的辐亮度进行线性拟合，得到拟合结果，并用所述拟合结果 替代对应像元的幅度亮度值r1数值，得到新幅度亮度值r2；
80.对所述新幅度亮度值r2进行拟合，得到拟合辐亮度f2；
81.根据幅度亮度值r1、拟合辐亮度f2以及预设的背景辐亮度数值，得到去 除杂散光后的目标图像。
82.作为一种具体的实施方式，杂散光去除算法的具体步骤如下：
83.(1)读取mersi/llb辐亮度数据和太阳天顶角，选取纯夜间数据(太阳 天顶角》101
°
)，设其中一条扫描线的原始辐亮度值为r1，使用savitzky-golay 滤波器，以61像元为窗口，三次多项式对r1做滑动拟合，得到拟合辐亮度f1。
84.(2)f1与r1做差并取绝对值后做11点滑动平均，对滑动平均结果设定 阈值进行判断，对连续大于阈值的像元提取首尾像元号。
85.(3)以提取首尾像元的像元号为自变量，其辐亮度为因变量做线性拟合， 拟合首尾像元之间像元的辐亮度，并用拟合结果代替对应像元的r1数值得到 r2。
86.(4)使用savitzky-golay滤波器对r2进行拟合，得到拟合辐亮度f2。
87.(5)r1减去f2，再加气辉造成的背景辐亮度数值，即得到去除杂散光的 结果。
88.需要说明的是，由于滤波的方式有多种，第1步和第4步拟合辐亮度的方 法并不唯一，比如使用滑动平均、傅里叶滤波、经验模态分解等，但其本质都是 提取辐亮度的低频变化趋势。第2步提取灯光像元的方法也有多种，比如通过 辐亮度变化的剧烈程度判断。第3步获取灯光像元基础辐亮度的方法也不唯一， 可以使用多项式拟合或直接赋值为首尾像元的辐亮度值。
89.再请参考图4a及图4b，是本技术算法的效果展示图。以2022年1月2日 22：15utc fy-3e在我国上空时mersi/llb的一个5分钟段数据为例验证杂 散光的去除算法(图4a)，由太阳高度角可知图像右侧为晨昏区，左侧为黑夜， 其左侧的白亮区为受杂散光影响所致。经我们的方案去除杂散后如图4b所示 (只显示纯夜间数据)，可见去除了左侧的白亮区，一些灯光像元得以展示出来， 效果基本令人满意。
90.从以上描述可以得知，实施本技术实施例提供的晨昏轨道卫星微光通道杂 散光去除方法，具有如下优点：
91.1、本技术相比于原有的杂散光去除方案更简单易行，不需要提前计算估算 杂散光的一系列参数。
92.2、原有杂散光去除方案每个月更新一次估算杂散光的参数，由于太阳方位 角在一个月内是不断变化的，因此所估算的杂散光的准确性会有变化。而本技术 的方案不受季节变化的影响，理论上去除杂散光的准确性是不变的。
93.3、本技术的方案能够去除帧间变化的杂散光，而原有方案不能实现。
94.基于相同的发明构思，本技术实施例提供了一种晨昏轨道卫星微光通道杂 散光去除装置，包括：
95.获取单元，用于获取mersi/llb原始图像；所述原始图像仅包括纯夜间数 据；
96.处理单元，用于采用杂散光去除算法对所述原始图像进行处理，以得到目标 图像；
97.其中，所述杂散光去除算法包括如下步骤：
98.(1)对所述原始图像逐扫描线进行滑动拟合，设其中某条扫描线的幅度亮 度值为r1，拟合辐亮度为f1；
99.(2)根据所述原始幅度亮度值r1和拟合幅度亮度f1提取连续灯光像元的 首尾像元的像元号和辐亮度；
100.(3)对连续灯光像元的辐亮度进行线性拟合，得到拟合结果，并用所述拟 合结果替代对应像元的原始幅度亮度值r1数值，得到新幅度亮度值r2；
101.(4)对所述新幅度亮度值r2进行拟合，得到拟合辐亮度f2；
102.(5)根据所述原始幅度亮度值r1、拟合辐亮度f2以及预设的背景辐亮度 数值，得到去除杂散光后的目标图像。
103.其中，步骤(1)具体为：
104.采用savitzky-golay滤波器，以61像元为窗口，三次多项式对原始幅度亮 度值r1做滑动拟合，得到拟合辐亮度f1；或
105.采用滑动平均、傅里叶滤波或经验模态分解对所述原始幅度亮度值r1做滑 动拟合，得到拟合辐亮度f1。
106.其中，步骤(2)具体为：
107.计算所述拟合幅度亮度f1与所述原始幅度亮度值r1的差值，取该差值的 绝对值后进行11点滑动平均，得到滑动平均结果；
108.采用设定阈值对所述滑动平均结果进行判断，对连续大于所述阈值的像元 进行首尾像元号的提取。
109.其中，步骤(3)具体为：
110.以提取的首尾像元的像元号为自变量，其辐亮度为因变量做线性拟合，拟合 首尾像元之间像元的辐亮度，得到拟合结果。
111.其中，步骤(4)具体为：
112.采用savitzky-golay滤波器对所述新幅度亮度值r2进行拟合，得到拟合辐 亮度f2；或
113.采用滑动平均、傅里叶滤波或经验模态分解对所述新幅度亮度值r2做滑动 拟合，得到拟合辐亮度f2。
114.其中，步骤(5)具体为：
115.计算所述原始幅度亮度值r1与拟合辐亮度f2的差值，再加上气辉造成的 背景辐亮度数值，即得到去除杂散光的结果；
116.其中，气辉造成的背景辐亮度数值由mersi/dnb在6～7月南极地区无杂 散的夜间观测而定。
117.可选地，本技术实施例还提供了另一种晨昏轨道卫星微光通道杂散光去除 装置。如图5所示，该杂散光去除装置可以包括：一个或多个处理器101、一个 或多个输入设备102、一个或多个输出设备103和存储器104，上述处理器101、 输入设备102、输出设备103和存储器104通过总线105相互连接。存储器104 用于存储计算机程序，所述计算机程序包括程序指令，所述处理器101被配置 用于调用所述程序指令执行上述方法实施例部分的方法。
118.应当理解，在本发明实施例中，所称处理器101可以是中央处理单元 (central processing unit，cpu)，该处理器还可以是其他通用处理器、数字信号 处理器(digital signal processor，dsp)、专用集成电路(application specificintegrated circuit，asic)、现成可编程门阵列(field-programmable gate array， fpga)或者其他可编程逻辑器件、分立门或者晶体管逻辑器件、分立硬件组件 等。通用处理器可以是微处理器或者该处理器也可以是任何常规的处理器等。
119.输入设备102可以包括键盘等，输出设备103可以包括显示器(lcd等)、 扬声器等。该存储器104可以包括只读存储器和随机存取存储器，并向处理器 101提供指令和数据。存储器104的一部分还可以包括非易失性随机存取存储 器。例如，存储器104还可以存储设备类型的信息。
120.具体实现中，本发明实施例中所描述的处理器101、输入设备102、输出设 备103可执行本发明实施例提供的杂散光去除方法的实施例中所描述的实现方 式，在此不再赘述。
121.需要说明的是，关于杂散光去除装置的具体实现方式，请参考前述方法实施 例部分，在此不再赘述。
122.以上所述，仅为本发明的具体实施方式，但本发明的保护范围并不局限于此， 任何熟悉本技术领域的技术人员在本发明揭露的技术范围内，可轻易想到各种 等效的修改或替换，这些修改或替换都应涵盖在本发明的保护范围之内。因此， 本发明的保护范围应以权利要求的保护范围为准。