天宫二号InIRA数据系统延迟改正估算方法①
2020-11-06陈嘉明廖静娟
陈嘉明 廖静娟
(*中国科学院空天信息创新研究院数字地球重点实验室 北京 100094)(**中国科学院大学 北京 100049)
0 引 言
据估计,大约2.1%的陆地表面为面积大于1 hm2的湖泊和水塘所覆盖[1]。湖泊是区域和全球范围内环境变化的关键因素[2],探测湖泊的水位变化可以为水资源短缺或洪水预测提供早期预警[3]。虽然传统的水文站观测方式能够获取连续、高精度的水位数据,但是运行成本高且站点数量有限,许多地方的水文站数据受地方政策限制不能公开使用[4,5];全球水文站的数量在近年来也呈下降趋势[6,7],因此湖泊水位变化的长时间探测显得极为迫切。从20世纪90年代开始,卫星测高技术在湖泊水位变化监测方面取得了广泛的应用[8-12],但是对于中小型湖泊,雷达回波的波形受到周围陆地表面污染,反演水位精度较低,Cryosat-2的出现在一定程度上解决了中小型湖泊水位的反演问题[12-15],但是传统星载雷达高度计受到地面分辨率和刈幅宽度的影响,观测效率仍然较低。
2016年9月15日,中国天宫二号空间实验站中搭载的3维成像微波高度计(interferometric image radar altimeter,InIRA)[16]成功升空,打破了传统星载雷达高度计只能沿星下点航迹方向进行高度测量,以及测量刈幅宽度小而无法进行大范围观测的局限,同时能够提供地表水体高程变化的2维观测[17]。与美国宇航局(National Aeronautics and Space Administration,NASA)2021年即将发射的SWOT(surface water ocean topography)搭载的Ka波段雷达干涉仪(Ka-band radar interferometers,KaRIn)[18]类似,InIRA利用双天线干涉测量获取高精度的地面回波干涉相位,结合波形跟踪及合成孔径使得InIRA 可在35~40 km的幅宽范围内以10 km×10 km的分辨率进行观测[19],从而有效避免刈幅内采样时间的不一致,获取高质量水体2维观测数据。迄今为止,InIRA已经在轨运行900余天,成功获取了陆地和海洋的观测数据,为进一步提高大型湖泊的水位反演精度和实现对中小湖泊的水位变化监测提供了重要支撑。
目前关于InIRA数据用于内陆水体的研究还很少[20],且由于InIRA仅为实验数据,数据的系统计时延迟、系统信号延迟均未进行改正,这对水位反演的精度影响极大,达到了米级的偏差[21],因此必须进行系统延迟改正来提高数据的测量精度。为此,本文基于有限的InIRA数据,提出了一种大气路径延迟改正和系统延迟改正的方法。利用该方法提高InIRA数据的测量精度,并通过青海湖的实测水位,验证及评价本文所提出方法的可靠性。
1 研究区与数据
1.1 研究区
为了使所求解的InIRA系统延迟改正量具有足够的代表性,本文选择了青藏高原地区16个具有不同特征的湖泊用以求解系统延迟改正量。所选湖泊位于中国青藏高原腹地,该区域气候寒冷干燥,年平均降水量低至247.3 mm[22],由于降水稀少,该地区的冰川为许多湖泊提供了重要的水源,湖泊水位在短时间内不会出现较大的变化。
图1显示了所选的16个湖泊的地理位置、形态以及周围的环境。表1列出了所选湖泊的中心坐标、面积。所选取的湖泊具有如下特征:湖泊形态各不相同,周边环境多为荒原、草地和高山,具有一定的地形起伏;所有湖泊在冬季都有结冰现象,冰期时间一般为当年11月至次年4月;湖泊受人类干扰影响较小;Cryosat-2过境时获取数据的工作模式均为SARin(synthetic aperture interferometric)模式。
图1 研究区地理位置、形态及周边地形
1.2 湖泊掩模数据
本文以2014年青藏高原湖泊子数据集[23]为掩模数据。该数据集全面、准确地展现了1 171个青藏高原1 km2以上湖泊的最新状态,包括名称、位置、面积和形态等特征,通过解译2014年获取的136幅中国高分1号宽视场(wide field of view, WFV)影像和11幅Landsat8 陆地成像仪(operational land imager, OLI)影像获得。
表1 InIRA过境所选16个湖泊数据统计
1.3 天宫二号 InIRA数据
InIRA是国际上第1个实现宽刈幅3维成像的微波高度计,它是为了评估近星下点成像和干涉测量技术在测量海面高度、水体范围和陆地高度的性能而设计的。仪器的主要参数如表2所示,InIRA突破了传统星载高度计只能沿星下点方向测量、观测刈幅只有数千米的局限,与传统星载雷达高度计的高度测量所不同的是,InIRA采用了偏离星下点1~8 °的小入射角、短基线的干涉测量方式,以获取具有高相干的干涉相位,并通过模型[24]进行3维地形的反演;干涉波束同时也为干涉基线倾角的测量提供了帮助[17],以弥补卫星姿态测量能力难以满足宽刈幅范围内高精度测高要求的不足,这对湖泊水位高精度反演具有重要意义。
表2 InIRA主要参数列表
InIRA数据可在中国载人航天工程空间应用数据推广服务平台(http://www.msadc.cn/sy/)中申请下载,有4种级别的数据可供下载,分别为0级、1级、2级和4级。天宫二号空间实验站常有其他观测任务,因此InIRA不能持续工作[24],故在内陆水体观测的数据覆盖区域有限,无法进行水位时间序列的研究。本文选取了InIRA过境的16个青藏高原湖泊1级数据进行系统延迟改正,同时为保证所求系统延迟量的稳定性,所选过境的1级数据均为同一天的数据;并选取青海湖2017-2018年的数据与实测数据进行对比验证。
1.4 Cryosat-2 SARin数据
Cryosat-2于2010年4月8日发射,主要用于监测海洋冰盖和大陆冰盖,所搭载的干涉合成孔径雷达高度计(SAR interferometer radar altimeter, SIRAL)是世界上第1台采用延迟多普勒技术的星载高度计[25],根据不同的地理模式可在3种不同的测量模式上工作[26],分别为LRM (low rate mode)、SAR (synthetic aperture radar)和SARin模式,本文中所采用的是青藏高原地区的SARin Baseline C版本的1b级与2级数据(ftp://science-pds.cryosat.esa.int)。Level 1b数据中提供了SARin 20 Hz波形,与SAR波形类似。SARin波形的前缘更加陡峭、后缘衰减也更快,这使得多波峰波形里包含多个可区分的表面高度[25],同时提供了通过各类算法得到的各星下点数据的对流层、电离层、固体潮、极潮等改正;Level 2级数据中提供了由Wingham/Wallis算法计算的高度,同时通过干涉处理重新定位了观测点位置[27],这对于偏离星下点距离改正[28]和筛选湖泊过境数据点具有重要意义。图2展示了部分本文所选取的青藏高原地区湖泊的InIRA及Cryosat-2数据的覆盖情况及周边环境图。
图2 研究区部分湖泊InIRA及Cryosat-2数据覆盖图
1.5 实测数据
青海湖实测水位数据通过青海省水文水资源勘测局获得,是由青海湖东南部下社站(36°35′16.0″N /100°29′28.4″E)水位计采集的非冰期水位数据[11]。该实测水位数据以1985中国国家黄海高程基准为参考,而本文反演的湖泊水位数据是以大地水准面为参考基准,故本文采用翟振和等人[29]提出的方法将黄海高程参考面转换至EGM2008大地水准面。
2 基于Cryosat-2数据的InIRA系统延迟改正估算
2.1 地球物理改正
InIRA数据在用于监测海洋或是内陆水体时,必须要对湿对流层、干对流层、电离层、固体潮、极潮、海潮及大地水准面等进行延迟改正。由于所发布的InIRA数据中没有进行地球物理改正,必须通过延迟改正模型进行地球物理改正,采用的各延迟改正模型以及相应的延迟改正范围见表3。采用的大地水准面模型为地球引力模型2008(EGM2008)[30]。
为进一步提高InIRA数据反演湖泊水位的精度,传感器本身的延迟误差需进行改正,如系统计时延迟、系统信号延迟等均会对高度计测高精度带来影响,而InIRA的L1b级数据并未进行相应的系统延迟改正,因此为获取更精确的湖泊水位需要进行InIRA系统延迟改正的估测,在2.2节中提出了基于Cryosat-2 SARin参考水位的InIRA系统延迟改正估算方法。
表3 地球物理改正模型统计
2.2 InIRA系统延迟改正
InIRA在内陆水体获取的观测数据有限,过境的内陆水体大都在青藏高原腹地,没有实测水位数据可以用于InIRA系统的延迟改正。Cryosat-2 SARin数据在青藏高原地区已经有大量的研究,即使对中小型湖泊的水位反演精度也较高[12,14]。故在无实测水位的情况下可采用Cryosat-2 SARin数据反演水位作为参考水位来替代实测水位对InIRA进行系统延迟改正。Cryosat-2 SARin数据采用双天线干涉测量方式获取,并可通过干涉处理获取观测点的水位信息,这与InIRA数据获取方式类似。此外由于Cryosat-2漂移轨道及92 °的近极地轨道的特性,使得Cryosat-2数据可以过境大部分InIRA数据过境的青藏高原区域湖泊。为此,本文选择Cryosat-2 SARin数据反演的高精度水位作为参考水位来估算InIRA系统延迟改正。主要选取InIRA数据过境30 d范围内的Cryosat-2数据进行处理获取高精度湖泊水位,并将此水位视为对应湖泊的精确水位对InIRA数据进行标定,从而估测InIRA系统的延迟。选取的16个湖泊的Cryosat-2数据以及InIRA数据的详情如表4所示。
2.2.1 Cryosat-2湖泊水位提取
对于Cryosat-2 SARin数据水位的提取,主要是通过处理1b级数据的波形数据获取精确的水面回波位置,其湖泊水位的计算公式如下:
-Ngeoid
(1)
其中,Halt为卫星质心至参考椭球的高度;Ngeoid为大地水准面至参考椭球的改正距离,本文选用1 min的EGM2008格网数据;c为光速;WD为窗口延迟;Rgeo为地球物理改正,包括电离层、干/湿对流层、固体潮、极潮和海潮改正,均可由 1b数据读取;Rretrack为重跟踪改正,本文采用ImpMWapp[12]重跟踪算法。
由于研究区过境湖泊面积较小,为最大程度利用Cryosat-2 SARin数据过境点数,本文选取湖泊过境点时未采用星下点坐标,而是采用Cryosat-2 SARin 2级数据干涉重定后的坐标点,这样可以更精确地反映高度计实际观测点的地理位置。此外,为精确获取Cryosat-2 湖泊水位,还需进行如下处理步骤:(1) 利用Cryosat-2 SARin 2级数据进行偏离星下点距离改正[28],消除实际观测点与星下点不在同一位置而产生的距离误差; (2) 对过境轨迹采用1.5倍中误差判断方法去除明显异常的水位值;(3) 利用高斯柯西误差混合模型[36]描述单点水位噪声项,采用最大似然估计法求解湖泊过境沿轨均水位。
2.2.2 InIRA湖泊水位提取
本节为InIRA未经系统延迟改正的湖泊水位反演方法,具体计算方法如下。
(1) 为避免Cryosat-2 SARin数据与InIRA数据的观测点坐标相差过大,本文只选取Cryosat-2 SARin 2级数据干涉重定后地面观测点的150 m范围内的所有InIRA过境数据值(沿轨数据)。
(2) 由于InIRA数据在内陆水体观测的过程中受数据质量的制约会出现误差,因此需对所获取的数据按照如下规则进行筛选:InIRA过境湖泊水位数据点大于10 000;过境数据点需离岸大于300 m;后向散射系数大于0;像素点幅度值大于4 500。
(3) 对于每景InIRA数据,沿轨数据的所有像素值按式(2)计算,得到未经系统延迟改正的湖泊沿轨单点水位。
H=Halt-(Hrange+Rgeo)-Hsystem-Ngeoid
(2)
式中所有参数均与式(1)中含义相同,其中Hsystem为InIRA系统延迟估计量,在本文中设为0;Hrange为传感器与星下点之间的距离;InIRA数据的地球物理改正量Rgeo,需根据2.1节中模型计算得到。
(4) 对每个湖泊对应的沿轨单点水位数据进行高斯均值滤波,减弱陆地回波信号对InIRA数据造成的干扰。
(5) 对过境沿轨数据采用1.5倍中误差异常水位值去除法,迭代去除超过1.5倍中误差范围的湖泊单点水位,直到得到没有超过1.5倍中误差的单点水位。
(6) 利用高斯柯西误差混合模型描述单点水位噪声项,采用最大似然估计法求解InIRA湖泊过境沿轨均水位。
表4 InIRA数据及Cryosat-2数据过境湖泊的详细信息
2.2.3 InIRA系统延迟改正量估测
通过上述计算得到Cryosat-2 SARin数据参考水位及InIRA数据沿轨水位,从而可以估算InIRA系统延迟改正。
Cryosat-2 SARin数据在经过地球物理改正、重跟踪改正、偏离星下点距离改正后可以获取湖泊的高精度水位,在无实测水位的情况下可以采用Cryosat-2 SARin数据为InIRA数据进行系统延迟改正。此外,研究区域气候干燥,降水量非常低,在短时间内湖泊水位变化较小。因此,InIRA数据在进行地球物理改正、高斯滤波、异常水位去除等处理后,所反演的湖泊沿轨水位与Cryosat-2参考水位在理论上的差值即为系统延迟改正。系统延迟改正估算的流程图如图3所示,具体处理过程如下。
(1) 获取研究区16个湖泊过境时间差为30 d之内的Cryosat-2参考水位以及未经系统延迟改正的InIRA数据的沿轨水位。
(2) 计算这16个湖泊范围内的Cryosat-2参考水位与InIRA沿轨均水位的差值。
(3) 为减弱不同湖泊计算得到的差值水位之间的误差,采用高斯柯西误差混合模型来描述差值水位的噪声项,可以有效避免异常差值水位对系统延迟改正的干扰,采用最大似然估计法求解InIRA系统延迟改正值。
图3 InIRA湖泊水位反演流程图
3 结果与讨论
3.1 InIRA系统延迟改正
按照2.2节中所述方法分别利用过境的InIRA沿轨数据及Cryosat-2数据对所选16个湖泊的绝对水位进行估算,得到InIRA和Cryosat-2 SARin数据反演的16个湖泊的绝对水位差值(表5)。
由表5可见,InIRA和Cryosat-2在同一湖泊过境范围内反演的湖泊绝对水位差值主要集中在4.5 m附近。但对于琼浆湖、美菊湖、玉琳湖等湖泊,该绝对水位差值表现为负值,这主要是由于上述出现水位差值异常的湖泊所覆盖的InIRA数据受陆地区域的影响,陆地表面对InIRA信号造成较大污染,产生大量异常观测,如InIRA过境银波湖区域存在大量的浅滩和低缓岗丘,使回波无法获取正确的水体信息,导致无法获取有效的湖泊水位数据。
表5 InIRA和Cryosat-2湖泊绝对水位差值
为获取稳定的InIRA系统延迟改正,本文进行InIRA数据的系统偏差估算时,剔除上述绝对水位差值小于0的湖泊数据,将剩下的绝对水位差值采用高斯柯西误差混合分布模型进行水位差的稳态均值计算,所得结果为4.6287 m,此值即为InIRA系统延迟改正的估算量。
3.2 青海湖湖泊水位验证
为验证本文提出的InIRA系统延迟改正估算方法的有效性,将加入系统延迟改正量后反演的2017-2018年青海湖过境水位与下社水文站实测水位进行比较,以此来间接验证本文方法的有效性。按照2.2.2节中水位提取方法计算得到2017-2018 年内InIRA过境青海湖沿轨均水位,并与下社站实测水位进行比较。不同时间段青海湖InIRA反演水位与实测水位比较结果如图4所示,沿轨水位数据已经过高斯滤波、中误差去除等方法进行了噪声去除,其中用于计算InIRA沿轨均水位数据的水位标准差及计算点数统计如表6所示。
从图4中可见在青海湖加入系统延迟改正后的InIRA沿轨均水位与实测水位的拟合程度较好,除图中的4个异常水位外,InIRA反演的青海湖沿轨均水位与实测水位的差值的绝对值在50 cm以下。同时水位差值的绝对值较小的数据相应的单点水位标准差均在1 m以下,沿轨均水位参与计算的点数也越多,局部单点水位的波动较小,能够获取较稳定的湖泊水位。为更直观地反映本文系统延迟改正计算的有效性,对InIRA数据反演青海湖湖泊水位进行精度评估,计算InIRA沿轨均水位及实测水位的无偏均方根误差ubRMSE,计算方法见式(3)。InIRA数据反演青海湖湖泊水位的无偏均方根误差为0.3327 m。结果表明加入本文的系统延迟改正后,InIRA数据反演湖泊水位的精度有显著提高,间接表明本文提出的系统延迟改正方法的有效性。
图4 青海湖InIRA沿轨均水位与实测水位比较图
表6 InIRA过境青海湖水位与实测水位差值
ubRMSE=
(3)
其中,E[*]为期望算子,WL为通过InIRA反演的湖泊水位,insitu为湖泊实测水位。
3.3 InIRA反演湖泊水位性能分析
为进一步分析InIRA数据在湖泊水位反演方面的性能,选取20170803[E]、20171006[NE_2]、20171216[SW]、20171216[SE_3]4景数据进行分析,其中前2景数据得到的青海湖沿轨均水位与实测水位的差值在米级,后2景数据的沿轨均水位与实测水位吻合程度较高,图5及图6分别展示了上述4景InIRA过境数据覆盖区域及沿轨单点水位的分布情况。
图5 部分InIRA过境数据覆盖区域
从中可以观察到,InIRA数据过境点数越多,单点湖泊水位越稳定;处于湖岸边InIRA数据,如20170803 [E]、20171006 [NE_2]局部数据非常不稳定,有效数据点数较少,所获取的单点水位最大波动可达10 m,这表明受到复杂地形的陆地回波的影响,可能会导致水面信号未进入系统跟踪窗口内或被噪声信号所淹没,造成局部单点水位的测量出现严重错误。而对于20171216 [SW]、20171216 [SE_3]2景数据,其大量数据处于湖中心区域,参与沿轨水位计算的点数非常丰富,单点水位标准偏差较小,数据观测质量稳定,有利于获取准确的相对水位。因此与传统高度计需选择优质波形反演水位相类似,InIRA反演水位时需选取过境点数丰富、单点水位标准偏差小的数据。此外,上述分析也在一定程度上反映了InIRA数据的稳定性存在问题。
4 结 论
本文针对InIRA数据特点,提出了基于Cryosat-2 SARin参考水位的InIRA系统延迟改正量估算方法。该方法的核心在于通过ImpMWapp算法获取得到一种稳定的高精度湖泊水位,同时筛选出Cryosat-2数据的干涉重定后的观测点150 m范围的InIRA水位数据,再通过两者的水位差来反映InIRA系统延迟改正量。通过有实测水位的青海湖来间接验证所获取的InIRA系统延迟改正的准确性。
图6 单点水位分布图
验证结果显示,2017-2018年间InIRA反演青海湖湖泊水位精度为0.3327 m,表明本文获取的InIRA系统延迟改正精度良好,这对于InIRA用于湖泊水位监测有着重要的意义。由于InIRA数据具有3维测高的能力,能够获取覆盖湖区范围的水位信息的空间分布图,对单点水位反演的精度更加关键。为进一步分析InIRA数据的精度,本文分析了若干湖泊的InIRA数据,发现同一湖泊不同区域的InIRA影像数据反演得到的水位有时存在明显差异,部分影像中不同像素点对应的单点水位的标准偏差较大,这从一定程度上反映了InIRA数据的稳定性存在一定问题。通过本文方法进行了InIRA数据的系统延迟改正,在一定程度上可提高InIRA数据反演湖泊水位的精度,但要从根本上提升3维成像高度计单点水位反演的精度,还需要依赖高度计本身设备性能的稳定,因此,未来的3维成像高度计需要全面提升其3维测高的性能。
同时,由于系统计时误差会随着时间的推移以及地域的不同而产生一定的变化,因此需要针对不同时间段、不同区域对系统进行时延校正,从而提高水位反演的精度,这将在进一步研究中开展相关工作。