CN113834547B - 一种河流虚拟站水位时序重建方法及系统 - Google Patents

一种河流虚拟站水位时序重建方法及系统 Download PDF

Info

Publication number
CN113834547B
CN113834547B CN202110838483.2A CN202110838483A CN113834547B CN 113834547 B CN113834547 B CN 113834547B CN 202110838483 A CN202110838483 A CN 202110838483A CN 113834547 B CN113834547 B CN 113834547B
Authority
CN
China
Prior art keywords
elevation
river
observation
waveform
water level
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
CN202110838483.2A
Other languages
English (en)
Other versions
CN113834547A (zh
Inventor
柯灵红
宋春桥
雍斌
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Hohai University HHU
Nanjing Institute of Geography and Limnology of CAS
Original Assignee
Hohai University HHU
Nanjing Institute of Geography and Limnology of CAS
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Hohai University HHU, Nanjing Institute of Geography and Limnology of CAS filed Critical Hohai University HHU
Priority to CN202110838483.2A priority Critical patent/CN113834547B/zh
Publication of CN113834547A publication Critical patent/CN113834547A/zh
Application granted granted Critical
Publication of CN113834547B publication Critical patent/CN113834547B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • GPHYSICS
    • G01MEASURING; TESTING
    • G01FMEASURING VOLUME, VOLUME FLOW, MASS FLOW OR LIQUID LEVEL; METERING BY VOLUME
    • G01F23/00Indicating or measuring liquid level or level of fluent solid material, e.g. indicating in terms of volume or indicating by means of an alarm
    • G01F23/22Indicating or measuring liquid level or level of fluent solid material, e.g. indicating in terms of volume or indicating by means of an alarm by measuring physical variables, other than linear dimensions, pressure or weight, dependent on the level to be measured, e.g. by difference of heat transfer of steam or water
    • G01F23/28Indicating or measuring liquid level or level of fluent solid material, e.g. indicating in terms of volume or indicating by means of an alarm by measuring physical variables, other than linear dimensions, pressure or weight, dependent on the level to be measured, e.g. by difference of heat transfer of steam or water by measuring the variations of parameters of electromagnetic or acoustic waves applied directly to the liquid or fluent solid material
    • G01F23/284Electromagnetic waves
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01SRADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S13/00Systems using the reflection or reradiation of radio waves, e.g. radar systems; Analogous systems using reflection or reradiation of waves whose nature or wavelength is irrelevant or unspecified
    • G01S13/88Radar or analogous systems specially adapted for specific applications
    • G01S13/882Radar or analogous systems specially adapted for specific applications for altimeters
    • YGENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
    • Y02TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02ATECHNOLOGIES FOR ADAPTATION TO CLIMATE CHANGE
    • Y02A90/00Technologies having an indirect contribution to adaptation to climate change
    • Y02A90/30Assessment of water resources

Landscapes

  • Physics & Mathematics (AREA)
  • Engineering & Computer Science (AREA)
  • Electromagnetism (AREA)
  • Radar, Positioning & Navigation (AREA)
  • Remote Sensing (AREA)
  • General Physics & Mathematics (AREA)
  • Computer Networks & Wireless Communication (AREA)
  • Thermal Sciences (AREA)
  • Fluid Mechanics (AREA)
  • Position Fixing By Use Of Radio Waves (AREA)
  • Testing Or Calibration Of Command Recording Devices (AREA)

Abstract

本发明公开了一种河流虚拟站水位时序重建方法及系统,在改进自适应子波分析的基础上,以中心脚点及临近点的平均波形为参照,提取多个满足参数条件的脉冲峰值并将其与沿轨波形的子波配对,形成多组纬向高程剖面,以波峰能量和高程标准误差构建多目标优化函数,确定目标水面高程观测组。最后通过时间序列噪声检测和高程重选,优化水位高程时间序列。本方法基于公开获取的卫星测高数据和辅助资料,可实现自动化的河流虚拟站点水位时序提取;并且能自动适应高噪声水平、多峰分布情形下的主波峰探测,提高复杂地表条件下的河流水位提取精度,可为河流水情监测、全球变化和水循环监测等应用提供方法支撑。

Description

一种河流虚拟站水位时序重建方法及系统
技术领域
本发明涉及水文信息学及遥感科学技术领域,特别涉及一种河流虚拟站水位时序重建方法及系统。
背景技术
河流承载着生态健康,经济可持续发展和人类福祉。水位是关键的河流水情要素之一,与气候变化、淡水资源供应、洪涝灾害等紧密联系。在气候变化和人类活动增强的影响之下,监测河流水位变化对于认识和应对水安全、能源风险、粮食安全等现存或潜在问题至关重要。
传统水位信息获取方法是设立水文监测站点进行观测,但是这种方式在偏远和气候恶劣区域受限,并且由于维护成本和人力限制,近年来全球水文观测站点数量显著下降。卫星测高技术具有大范围、长时序、自动化观测的优势,可以填补大量无观测资料区域的空缺,加密水文监测网络。在卫星轨道与河流网络的交叉点上,定期提取河流水情信息,形成类似无人值守的水文观测站,被称为虚拟观测站(virtual station)技术。目前此项技术在某些大河如亚马逊河流的某些河段取得了较成功的应用,但是成功监测的河段相比卫星实际覆盖度非常有限。以河流网系发达的中国东南部为例,基于Jason、Sentinel和CryoSat-2等测高卫星的河流虚拟站个数在个数级别,并且部分精度验证结果表明水位精度低至2-3米,难以满足实际应用需求。
当前卫星高度计观测河流的潜能很大受限于雷达高度计时空分辨率矛盾的约束。由于雷达高度计脚印点大小在公里或者几十公里级别,在观测复杂地形条件下或异质性地表的内陆水体时,雷达高度计返回的波形数据受地面异质性反射源污染呈现复杂的多峰信号和噪声,导致传统水位提取算法容易失效。目前针对内陆水体(湖泊、水库、河流)水位监测的方法主要分为两大类,一类是基于传统波形重订算法(如OCOG,ICE-1)的GDRs(Geographical Data Records)数据,主要对沿轨高程值进行噪声过滤确定水位观测值(Coss et al.,2019;Okeowo et al.,2016;Silva et al.,2010)。另一类方法基于原始波形数据SGDRs(Sensor Geophysical Data Record),通过发展新的波形重订算法,进行子波分析(多峰波形识别和分离)和波形重订来确定陆面水位信号,包括NPPR(Narrow PrimaryPeak Retracker)(Jain et al.,2015),MWaPP(Multiple Waveform Persistent Peak)(Villadsen et al.,2016),SMWTR(Multi-subwaveform Multi-weight ThresholdRetracker)(Yuan et al.,2017),TIC(50%Threshold and Ice-1Combined algorithm)(Huang et al.,2018),AMPDR(Automatic Multiscale-based Peak DetectionRetracker)(Chen et al.,2020)等。第一类算法在大型湖泊和水库应用中简单有效,而在河流监测中因为传统波形踪算法的失效而难以成功。后一类方法能够处理更复杂的波形情况,其难点在于如何确定多峰波形中与水体相关的主波峰信号,目前算法的主要策略包括,1)将最大返回能量的峰值为参照(Jain et al.,2015;Jiang et al.,2020),2)将第一个峰值能量超过最高能量20%作为主信号(Villadsen et al.,2016),3)对同时段观测高程值统计累计密度函数CDF确定参考高程(Chen et al.,2020),4)采用人工选取或者先验知识来确定(SMWTR和TIC)。第二类算法中大部分是针对湖泊水库研制,而对于更精细的地表水体河流而言,其水面范围更加有限,单纯基于统计分析的算法因为对于大样本量的需求仅适用于超大型河流。
综上所述,河流水位是关系到水安全、灾害评估和预防和气候变化的关键信息,基于卫星测高的河流虚拟站技术是目前水文信息自动化的前沿和方向。现有的水位提取算法需要人工干预或依赖较精细的先验知识,或基于大样本的统计算法,极大限制了河流虚拟站技术的推广和应用。在全球变化加剧的背景下,提出一种针对河流等精细地表水体、能自动适应高噪声水平和多峰分布的波形重跟踪方法,可为河流水情监测、区域水循环和水资源评估、灾害预警和评估等应用提供方法支撑,具有重要的科学意义和应用潜力。
发明内容
本发明的目的在于克服现有技术中的不足,提供一种河流虚拟站水位时序重建方法及系统,能自动适应高噪声水平、多峰分布情形下的主波峰探测,提高复杂地表条件下的河流水位提取精度。
为达到上述目的,本发明是采用下述技术方案实现的:
第一方面,本发明提供了一种河流虚拟站水位时序重建方法,包括以下步骤:
获取卫星观测数据;
根据所述卫星观测数据确定河流虚拟站位置和水面范围,根据所述河流虚拟站位置和水面范围定义缓冲区,根据所述缓冲区的范围提取出该范围内的所有观测脚点的雷达测高卫星观测波形和辅助数据,观测脚点即卫星星下点位置;
对缓冲区内的所有雷达测高卫星观测波形和辅助数据按照观测时段分组,对每组数据进行波形分析和子波分解,获得多组纬向高程剖面;
对所述多组纬向高程剖面,进行星下点偏移效应改正后,构建多目标优化函数,以损失函数最小的组纬向高程剖面的高程值作为初始水面高程序列,同时记录其它组纬向高程剖面的高程值作为备用高程组;
以初始水面高程序列为基础,进行噪声识别获取噪声点,并从所述备用高程组中选择高程值替换噪声点的高程初始值,获得优化后的水位高程时间序列。
进一步的,所述卫星观测数据包括测高卫星轨迹线、全球河网中心线和水体频率图以及卫星测高原始波形数据SGDRs;
根据所述卫星观测数据确定河流虚拟站位置和水面范围,根据所述河流虚拟站位置和水面范围定义缓冲区,根据所述缓冲区的范围提取出该范围内的所有观测脚点的雷达测高卫星观测波形和辅助数据的方法包括以下步骤:
将测高卫星轨迹线向外扩充2千米形成轨迹放大区;
根据所述轨迹放大区和全球河网中心线相交得到河流虚拟站位置,根据水体频率图提取河流水面范围并计算河段平均宽度,并根据河段平均宽度筛选出平均河宽大于阈值的河流虚拟站位置以确定可重建水位的目标地点;
基于筛选出的河流虚拟站位置,提取出该河段水面范围,并将河流范围向外扩充1.5千米形成缓冲区,在所述卫星测高原始波形数据SGDRs中提取出所述缓冲区内的所有雷达测高卫星观测点波形和辅助数据。
进一步的,对缓冲区内的所有雷达测高卫星观测波形和辅助数据按照观测时间分组,对每组数据分别进行波形分析和子波分解,获得多组纬向高程剖面的方法包括以下步骤:
选取组中每个观测时段所有观测脚点的中心点及其邻域观测点的波形合成平均波形,采用自适应峰值探测方法和逐步细致的参数策略提取所述平均波形中的脉冲峰值;
对该组观测中所有的波形数据,以最大能量、噪声和峰值水平为依据筛选单点波形数据,将筛选后的单点波形数据归一化后利用自适应脉冲峰值探测方法识别子波,并剔除微小的噪声子波,采用ICE-1跟踪算法确定每个子波对应的高程值;
对筛选后的每个单点波形,将其子波与平均波形中的脉冲峰值进行一对一匹配,并剔除脉冲峰值对应的高程差值的绝对值大于2米及重复选取的配对,获得多组纬向高程剖面。
进一步的,选取组中每个观测时段所有观测脚点的中心点及其邻域观测点的波形合成平均波形的方法包括以下步骤:
对每一观测时段的卫星波形,按照河段宽度调整参与合成的波形个数n,n为奇数且不小于5,以观测脚点中心点及其左右邻域各(n-1)/2个点进行波形合成;
所述波形合成的方法为:计算参与合成的波形数据采样门的高程范围,将每个波形数据插值到该范围内,然后按照中值合成,获得平均波形;
按照河段宽度调整参与合成的波形个数的方法为:参与合成的波形个数为河段宽度除以卫星采样距离的最近的奇数,并且最小为5个波形合成。
进一步的,采用自适应峰值探测方法和逐步细致的参数策略提取所述平均波形中的脉冲峰值的方法包括以下步骤:
对于平均波形,利用改进的自适应峰值自动探测算法识别脉冲峰值;
所述过滤参数为:脉冲峰值探测窗口大小为3个采样门宽度,脉冲峰值能量大于25%的最大能量,脉冲峰值上升幅度大于10%的最大能量,并且上升幅度大于下降幅度的35%;
如果提取出的脉冲峰值点对应的高程值的最小值小于前一观测时间的高程值,则减小窗口大小,降低脉冲峰值能量阈值,降低上升沿能量幅度阈值,再次检测脉冲峰值;
采用ICE-1算法确定跟踪门并计算脉冲峰值对应地面目标的高程。
进一步的,将每个单点波形的子波与平均波形中的脉冲峰值进行一对一匹配的方法包括:采用脉冲峰值和跟踪高程两种标准配对每个单点波形的子波与平均波形中的脉冲峰值,并以最小标准差为标准选取最终的配对;通过领域距离分析保证每个子波只对应一组脉冲峰值。
进一步的,以初始水面高程序列为基础,进行噪声识别获取噪声点,并从所述备用高程组中选择高程值替换噪声点的高程初始值的方法包括:以初始水面高程序列为基础,采用滑动窗口平均拟合和时间序列拟合进行噪声识别,从所述备用高程组中重新确定噪声点的目标高程;迭代停止的方式包括限制总迭代次数,是否有噪声和高程序列是否改变。
进一步的,所述方法还包括:基于附近的实测水位数据,实现对所构建的虚拟站水位序列的精度分析。
第二方面,本发明提供了一种河流虚拟站水位时序重建系统,所述系统包括:
获取模块:用于获取卫星观测数据;
提取模块:用于根据所述卫星观测数据确定河流虚拟站位置和水面范围,根据所述河流虚拟站位置和水面范围定义缓冲区,根据所述缓冲区的范围提取出该范围内的所有观测脚点的雷达测高卫星观测波形和辅助数据,观测脚点即卫星星下点位置;
高程剖面模块:用于对缓冲区内的所有雷达测高卫星观测波形和辅助数据按照观测时段分别进行波形分析和子波分解,获得多组纬向高程剖面;
水面高程模块:用于对所述多组纬向高程剖面,进行星下点偏移效应改正后,构建多目标优化函数,以损失函数最小的组纬向高程剖面的高程值作为初始水面高程序列,同时记录其它组纬向高程剖面的高程值作为备用高程组;
噪声剔除模块:用于以初始水面高程序列为基础,进行噪声识别获取噪声点,并从所述备用高程组中选择高程值替换噪声点的高程初始值,获得优化后的水位高程时间序列。
第三方面,本发明提供了一种河流虚拟站水位时序重建系统,其特征在于,包括处理器及存储介质;
所述存储介质用于存储指令;
所述处理器用于根据所述指令进行操作以执行根据第一方面所述方法的步骤。
与现有技术相比,本发明所达到的有益效果:
1、本发明以中心脚点及临近点的平均波形为参照,提取多个满足参数条件的脉冲峰值并将其与沿轨波形的子波配对,形成多组纬向高程剖面,以波峰能量和高程标准误差构建多目标优化函数,确定目标水面高程观测组。最后通过时间序列噪声检测和高程重选,优化水位高程时间序列,本方法可基于公开获取的卫星测高数据和辅助资料,可实现自动化的河流虚拟站点水位时序提取;并且能自动适应高噪声水平、多峰分布情形下的主波峰探测,提高复杂地表条件下的河流水位提取精度,可为河流水情监测、全球变化和水循环监测等应用提供方法支撑;
2、本发明提出的河流水位序列重建方法不依赖实测数据和先验知识,基于可公开获取的卫星数据和河网基础数据可完成;
3、本发明算法对雷达测高波形噪声水平不一、波形特征变化较大的情况有良好的适用性,能够提高精细尺度、复杂地表环境下的河流水位提取精度,扩展目前的河流虚拟站监测潜能,为全球变化背景下水文水资源研究提供了有效工具。
附图说明
图1本发明实施例样本,所选择的河流虚拟站位置及卫星观测覆盖脚点选取。黑色圆点为某一观测时段的纬向脚点,白色为参与平均波形计算的脚点;
图2本发明算法流程图;
图3本发明实施例中针对不同情形下的自适应峰值探测和子波分解结果图,以存在两组脉冲峰值情形为例(观测时间2009年3月10日),a)为纬向观测的波形序列,b)中间某一波形的子波分解和及其与平均波形的峰值对应结果,c)平均波形及其峰值检测,d)纬向高程剖面和初始选定高程;
图4本发明实施例中针对不同情形下的自适应峰值探测和子波分解结果图,以高噪声水平的多组峰值情形为例(观测时间2017年3月13日);
图5本发明实施例中星下点偏离效应改正方法及结果:a)存在明显星下点偏离效应,用二次多项式拟合得到定点高程值为该次测量目标高程;b)星下点偏离效应不显著,以高程剖面中值为该次测量目标高程;
图6本发明实施例中的时间序列去噪效果图,a)为去噪前后水位时间序列及噪声点,b)为水位时间序列的多年均值及标准茶;
图7本发明实施例中的河流虚拟站水位精度评价,a)本算法提取水位与实测水位、传统算法结果比较,b)本算法提取结果与实测数据散点图。
具体实施方式
下面结合附图对本发明作进一步描述。以下实施例仅用于更加清楚地说明本发明的技术方案,而不能以此来限制本发明的保护范围。
实施例一:
本实施例提供一种新型的智能化、自动化水位序列重建算法,以扩展河流监测潜能的河流虚拟站水位时序重建方法。该方法采用可公开获取的测高卫星原始波形数据、水系网络数据和地表水体频率分布图,在改进自适应子波分析的基础上,融合多波形多峰识别、多目标优化主波探测、噪声迭代识别等主要策略重建高精度的水位长时间序列,包括以下过程:
步骤1、将测高卫星轨迹线和全球河网中心线相交得到河流虚拟站位置,根据水体频率图提取河流水面范围并计算河段平均宽度,并根据平均河宽筛选河宽大于一定阈值的流虚拟站位置;并不是所有河段都能顺利提取出水位时间序列,只有河段河宽大于一定范围才有可能提取出水位值;基于河段宽度的筛选,意思是平均河宽大于1000米或者500米,才确定在这个位置上提取卫星观测值(具体阈值取决于卫星参数等,一般需要河宽大于百米级);确定河流虚拟站位置,才能确定在河流的哪些位置上得到水位,才能确定在哪里建缓冲区,从而提取缓冲区内的卫星观测数据;
步骤2、针对筛选出的河流虚拟站,将河流水面范围向外扩充1.5千米建成缓冲区,提取此缓冲区内的所有雷达测高卫星观测波形和辅助数据;
步骤3、将以上缓冲区内的卫星观测数据按照观测时段分组分别进行处理,选取每个时段所有观测脚点的中心点及其邻域的若干波形合成平均波形,采用自适应峰值探测方法和逐步细致的参数策略提取平均波形中的有效脉冲峰值及其参数;
步骤4、以最大能量、噪声和峰值水平为依据筛选单点波形数据,对筛选后的波形分别进行归一化后利用自适应脉冲峰值探测方法识别子波,并剔除微小的噪声子波,采用ICE-1跟踪算法确定每个子波对应的高程值;
步骤5、将每个单点波形的子波与平均波形中的脉冲峰值以最近距离进行一对一匹配,并剔除对应高程差值太大(>2米)的配对及重复选取情况;
步骤6、对每个脉冲峰值所形成的纬向高程剖面,进行星下点偏移效应改正后,以最大脉冲能量、高程组标准误差为变量异构建多目标优化函数,以损失函数最小的值作为该时段水面高程初始值,同时记录其它组高程值作为备用;
步骤7、以初始水面高程序列为基础,采用迭代方法进行噪声识别并从步骤6的备用高程组中重新目标高程,直到目标条件达到停止噪声过滤;
步骤8、基于附近的实测水位数据,实现对所构建的虚拟站水位序列的精度分析。
本发明所述的方法适用于雷达测高数据包括Jason系列卫星,Sentinel-3和CryoSat-2;全球河网中心线选用Global River Widths from Landsat(GRWL)Database;水体频率图为基于长时间光学卫星进行水体制图后形成的陆面水体分布频率图(Globalsurface water frequency)。
步骤1中,对于虚拟站位置的确定,将轨道或者河流中心线均做缓冲区后再求交,以扩大可能的虚拟站搜索范围;提取河流水面范围,可以根据水体频率图或者现有光学卫星数据提取。如采用光学卫星数据的方式提取水面范围,可以选取平均水位时期的光学影像如Landsat/Sentinel,利用水体分类方法直接提取出水面范围。
步骤2中,可以根据卫星脚点大小相应增加或者减少缓冲区大小以改变输入观测脚点个数,单时段落入缓冲区内的脚点个数应大于10。一个脚印点指雷达观测时刻的星下点空间位置;根据卫星观测提供的经纬度能得知该观测所代表的空间位置。
步骤3中,平均波形是以轨迹中心点附近的波形按照高程重采样后以中值合成的,参与合成的波形个数按照河段宽度调整(个数为奇数,并最小为5);平均波形中的峰值探测采用不同的参数策略,不仅探测出最大能量波峰,且探测出所有满足参数条件的前6个最大能量的脉冲峰值;也可以将平均波形做平滑去噪处理后再识别峰值。
步骤4中,采用更加宽松的策略过滤子波,尽可能保留更多的波形;自适应脉冲峰值探测设定为最小尺度(scale=3),利用峰值,脉冲上升幅度、回波识别等方案过滤噪声子波。
步骤5中,采用脉冲峰值和跟踪高程两种标准配对平均高程峰值和子波,并以最小标准差为标准选取最终的配对;通过领域距离分析保证每个子波只对应一组脉冲峰值。
步骤6中,采用多目标优化函数从多组高程剖面中选取水面目标高程,该目标函数以峰值能量最大、高程剖面标准误差最小、与前一时段参考高程绝对差异最小为综合考量选取目标水位测量值。
步骤7中,采用Savitzky-Golay拟合和滑动中值拟合值为约束,迭代进行噪声选取,并从步骤6的备用高程组中重新选取目标高程;迭代停止的方式包括限制总迭代次数,是否有噪声和高程序列是否改变。
本发明具有以下两点优势:
本发明提出的河流水位序列重建方法不依赖实测数据和先验知识,基于可公开获取的卫星数据和河网基础数据可完成;本发明算法对雷达测高波形噪声水平不一、波形特征变化较大的情况有良好的适用性,能够提高精细尺度、复杂地表环境下的河流水位提取精度,扩展目前的河流虚拟站监测潜能,为全球变化背景下水文水资源研究提供了有效工具。
实施例二:
本申请的实施例为长江流域黄石港站附近的河流虚拟站点,基于Jason2和Jason3卫星观测数据建立该站2008-2019年的水位时间序列。在该站3千米范围内有实测水文站(黄石港口水文站)数据,便于对算法精度的验证。如图1所示,该段河流平均宽度为1.15千米,河流周边包含城市建筑、自然植被、湖泊水库、人工水域等地表覆盖,导致波形形状、振幅复杂多变,在河流水位提取中具有很大挑战性。
如图2所示,是实施例的流程图,本实施例包括以下步骤:
步骤一,确定河流虚拟站位置和水面范围。此步骤初步确定可能建立虚拟站的位置和河段水面范围。卫星轨迹大致确定了卫星每次过境位置,而每次实际过境路线(卫星观测的星下点位置的连线)与预先确定的轨迹偏差有百米至公里的差距,因此将卫星轨迹做2千米缓冲区后再与河流中心线相交;以交点为中心,建立5千米缓冲区范围剪陆面水体分布频率图,利用10%频率进行栅格分割得到水体范围。图1显示了实例虚拟站位置,水体频率图以及提取出的河流水面范围。
步骤二,定义缓冲区和提取卫星观测数据。将河流水面范围向外扩充1.5千米成缓冲区,此大小的缓冲区包含所有可能观测到水面信号的卫星观测脚点。从原始SGDRs数据中解译出需要的20赫兹和1赫兹参数,分别包括星下点经纬度,卫星高度,未重订跟踪距离,以及1赫兹原始波形数据以及各种地球物理改正参数。
步骤三,平均波形合成和峰值探测。对每一观测时段的卫星波形,以卫星观测位置的中心点位置及其左右邻域观测点的波形进行波形合成,如图1所示。参与合成的波形个数为河段宽度除以卫星采样距离的最近的奇数,并且最小为5个波形合成。初始化的合成方法为:计算参与合成的波形数据采样门的高程范围,按照5厘米的分辨率将每个波形数据插值到该范围内,然后按照中值合成;
对于平均波形,利用改进的自适应峰值自动探测算法(Automatic Multi-scalePeak Detection)识别峰值,因为噪声水平的不同,采用不同的探测和过滤参数识别主峰。初始化参数为:峰值探测窗口大小为3个采样门宽度,峰值能量大于25%的最大能量,峰值上升幅度大于10%的最大能量,并且上升幅度大于下降幅度的35%(非回波)。如果提取出的峰值点对应的高程值的最小值小于前一观测时间的高程值,则减小窗口大小,降低峰值能量阈值,降低上升沿能量幅度阈值,再次检测峰值;峰值起点和终点通过跟踪方法确定,通过向前(向后)搜索到第一个下降拐点且能量降幅超过最大幅度的50%作为上升沿起点(终点);采用ICE-1算法确定跟踪门并计算对应地面目标的高程。图3c)和图4c)显示了不同时间、不同特征的子波合成的平均波形特征及识别出的峰值结果。
步骤四,波形过滤、子波分解和过滤以及波形重跟踪。首先计算每个波形的尖度指数(peak pulse)和噪声水平(前5个采样门和后5个采样门的能量均值,将最大能量小于10且脉冲峰度小于0.5且噪声水平大于30%振幅的波形作为噪声剔除;将波形能量按照最大最小值归一化到0至1的范围,再利用步骤三中的自适应峰值探测算法,选取窗口大小为3个采样门的尺度探测脉冲峰值,利用上升沿能量增加幅度大于10%振幅,且能量上升幅度大于下降幅度的35%的子波过滤纯噪声子波。利用步骤四中同样的搜索方法,确定子波起点和终点。波形跟踪点的计算方法同样采用ICE-1算法。该算法通过搜索第一个能量大于30%子波振幅的采样门,通过线性插值确定跟踪门。基于跟踪门,采用以下公式计算水面高程:
Hr=Hsensor-Rr-(Rnorm+Ciono+Cdry+Cwet+Cset+Cpt+hgeoid)
其中Hr为重跟踪点相对于EGM96的正高,Hsensor为传感器距离地面高度,Rnorm为预设重跟踪点采样门对应的距离,Ciono为电离层改正,Cdry为干改正,Cwet为是电离层改正项,Cpt为极潮改正,hgeoid为EMG96水准高,以上改正项均从卫星数据集中提供;Rr代表重跟踪点对应的改正距离,计算公式为:
Rr=S×F×(G-Gnorm)/2
其中S为光的传播速度速,G为ICE-1算法计算出的重跟踪点采样门;F和Gnorm分别为卫星阈门时间和预设跟踪点采样门,对应Jason卫星分别为3.125ns和32。
步骤五,搜索步骤四中的子波,确定与步骤三中的各峰值对应的子波序列,形成单组或多组纬向高程组。分别以峰值点和跟踪点的高程为参照,在绝对差值2米的范围内搜索最近距离的子波,在两组高程组中选取标准差最小的作为最匹配的高程组;反向搜索每个子波,确保其只对应某一脉冲峰值组,对于不同峰值组重复选取同一子波的情况,通过计算其与最临近脚点的高程距离确定属于哪一高程组。图3显示在存在两个有效主波峰的识别和配对情形,a)纬向波形序列,b)单波形的子波分解及重跟踪点位置,c)平均波形的脉冲峰值识别结果,d)平均波形中每个脉冲峰值对应的纬向(沿轨方向)高程剖面;图4显示了噪声水平较高、存在多组主波峰的对应结果,a),b),c)和d)分别显示波形序列,单波形的子波分解,平均波形的脉冲峰值识别和提取出的多组纬向高程剖面。
步骤六,星下点偏离效应改正和初始高程值确定(确定主波峰)。判断每组纬向高程剖面是否存在明显的星下点偏离效应效应并进行改正:利用二次多项式拟合高程序列和纬度的关系,判断拟合是否合理(拟合误差在0.5米以内且在观测纬度范围内存在最大值),对于合理的情况以二次多项式顶点高程值作为高程期望值,以拟合标准误差RMSE作为误差,对于不合理的情况,以中值作为高程期望值,以平均绝对误差MAE作为误差。图5显示了a)显著的偏离星下点效应,拟合的二次多项式和目标高程值确定,b)不存在明显的偏离星下点效应,利用中值和平均绝对误差作为高程组的均值和误差。对于每组高程值,建立多目标优化函数,按照以下公式对每组的平均峰值能量、标准误差计分,峰值能量最大、标准误差最小的组为初始目标主波峰结果。多目标优化函数的表现形式,以输入峰值能量为例:
Spower=Wpower×Floor(Bpower×Npower)
其中Spower为平均波形中的峰值能量得分,Ceil为向下取整函数,Wpower为能量值的权重(设为3),Bpower为对能量的分层个数(设为4),Npower为峰值能量按最大最小值归一化到[0,1]范围后的值,其中最大能量值为1,最小为0。类似地,以标准误差为输入变量,其得分计算公式为:
Sse=Wse×Ceil(Bse×Nse)
其中:
Figure BDA0003178027060000171
Error代表误差,count为该组观测值个数,Nse标准误差的倍数,Wse和Bse设为5。损失函数的计算公式为得分之和:
Stotal=Spower+Sse
按照此目标函数,能量值越高,标准误差越小的高程组为目标水面信号,图3d)和图4d)中实心圆点显示了多组高程里面确定的初始水面高程。
步骤七,对初始水位时间序列进行噪声识别和重选。首先进行多年均值统计,将每个观测时段的高程值基于3倍标准差过滤后计算中值和标准差,得到多年均值及其波动幅度。然后参照波动幅度水平,利用多种噪声识别和约束条件进行迭代噪声识别和重选。过程为:利用Savitzky-Golay对高程序列进行全局拟合,以拟合值为参照值,计算高程序列与参照值的绝对差值,超过3倍标准差范围的识别为潜在噪声点,搜索该时段内的高程组,以参照值为目标从备用组中选取最接近的高程值为新值,更新多年均值和标准差后再循环进行噪声探测和重选。第二步,以3个窗口大小对水位序列进行滑动平均拟合,观测值与拟合值的绝对差值超过3倍标准差判断为噪声,以拟合均值为参照值选取该时段的高程组中最接近的值,更新多年均值和标准差后再循环进行噪声探测和重选。两个步骤中满足以下三个条件之一则迭代停止,总迭代次数超过10,没有检测出潜在噪声,重选后高程序列没有改变。图6a)显示了去噪前后的水面高程序列、探测出的噪声点,b)为水面高程的多年均值信号及其标准差,其中标准差用于界定噪声。
步骤八,为评估算法精度,将实施案例的反演水位序列与周边的实测水位序列去除系统偏差后进行比对,如图7所示,a)显示本算法提取的水位序列与实测数据、传统算法比较,本算法提取的水位序列与实测值的决定系数(R2)为0.98,均方根误差(RMSE)为0.5米。而传统波形跟踪算法提取结果基本无法反应真实水位变化,特别在冬季低水位期,由于波形噪声水平高,提取值与实测数据相关性不显著。b)显示了本算法结果与实测结果的散点分布图。在雷达波形具有显著不同的特征情形下(图3,图4),本算法均能有效分离出水域信号,大部分情况下(95%概率)目标损失函数能识别出与水面信号相关的主波峰(图7b),其余通过噪声探测和高程重选能找出与实际水位最接近的峰值信号。
实施例三:
本实施例提供一种河流虚拟站水位时序重建系统,所述系统包括:
获取模块:用于获取卫星观测数据;
提取模块:用于根据所述卫星观测数据确定河流虚拟站位置和水面范围,根据所述河流虚拟站位置和水面范围定义缓冲区,根据所述缓冲区的范围提取出所述缓冲区内的所有雷达测高卫星观测脚印点及其波形和辅助数据;
高程剖面模块:用于对缓冲区内的所有雷达测高卫星观测脚印点及其波形和辅助数据进行波形分析和子波分解,获得多组纬向高程剖面;
水面高程模块:用于对所述多组纬向高程剖面,进行星下点偏移效应改正后,构建多目标优化函数,以损失函数最小的组纬向高程剖面的高程值作为初始水面高程序列,同时记录其它组纬向高程剖面的高程值作为备用高程组;
噪声剔除模块:用于以初始水面高程序列为基础,进行噪声识别获取噪声点,并从所述备用高程组中选择高程值替换噪声点的高程初始值,获得优化后的水位高程时间序列。
实施例四:
本发明实施例还提供了一种河流虚拟站水位时序重建系统,其特征在于,包括处理器及存储介质;
所述存储介质用于存储指令;
所述处理器用于根据所述指令进行操作以执行根据实施例一所述方法的步骤。
本领域内的技术人员应明白,本申请的实施例可提供为方法、系统、或计算机程序产品。因此,本申请可采用完全硬件实施例、完全软件实施例、或结合软件和硬件方面的实施例的形式。而且,本申请可采用在一个或多个其中包含有计算机可用程序代码的计算机可用存储介质(包括但不限于磁盘存储器、CD-ROM、光学存储器等)上实施的计算机程序产品的形式。
本申请是参照根据本申请实施例的方法、设备(系统)、和计算机程序产品的流程图和/或方框图来描述的。应理解可由计算机程序指令实现流程图和/或方框图中的每一流程和/或方框、以及流程图和/或方框图中的流程和/或方框的结合。可提供这些计算机程序指令到通用计算机、专用计算机、嵌入式处理机或其他可编程数据处理设备的处理器以产生一个机器,使得通过计算机或其他可编程数据处理设备的处理器执行的指令产生用于实现在流程图一个流程或多个流程和/或方框图一个方框或多个方框中指定的功能的系统。
这些计算机程序指令也可存储在能引导计算机或其他可编程数据处理设备以特定方式工作的计算机可读存储器中,使得存储在该计算机可读存储器中的指令产生包括指令装置的制造品,该指令装置实现在流程图一个流程或多个流程和/或方框图一个方框或多个方框中指定的功能。
这些计算机程序指令也可装载到计算机或其他可编程数据处理设备上,使得在计算机或其他可编程设备上执行一系列操作步骤以产生计算机实现的处理,从而在计算机或其他可编程设备上执行的指令提供用于实现在流程图一个流程或多个流程和/或方框图一个方框或多个方框中指定的功能的步骤。
以上所述仅是本发明的优选实施方式,应当指出,对于本技术领域的普通技术人员来说,在不脱离本发明技术原理的前提下,还可以做出若干改进和变形,这些改进和变形也应视为本发明的保护范围。

Claims (10)

1.一种河流虚拟站水位时序重建方法,其特征在于,包括以下步骤:
获取卫星观测数据;
根据所述卫星观测数据确定河流虚拟站位置和水面范围,根据所述河流虚拟站位置和水面范围定义缓冲区,根据所述缓冲区的范围提取出该范围内的所有观测脚点的雷达测高卫星观测波形和辅助数据,观测脚点即卫星星下点位置;
对缓冲区内的所有雷达测高卫星观测波形和辅助数据按照观测时段分组,对每组数据进行波形分析和子波分解,获得多组纬向高程剖面;
对所述多组纬向高程剖面,进行星下点偏移效应改正后,构建多目标优化函数,以损失函数最小的组纬向高程剖面的高程值作为初始水面高程序列,同时记录其它组纬向高程剖面的高程值作为备用高程组;
以初始水面高程序列为基础,进行噪声识别获取噪声点,并从所述备用高程组中选择高程值替换噪声点的高程初始值,获得优化后的水位高程时间序列。
2.根据权利要求1所述的河流虚拟站水位时序重建方法,其特征在于,所述卫星观测数据包括测高卫星轨迹线、全球河网中心线和水体频率图以及卫星测高原始波形数据SGDRs;
根据所述卫星观测数据确定河流虚拟站位置和水面范围,根据所述河流虚拟站位置和水面范围定义缓冲区,根据所述缓冲区的范围提取出该范围内的所有观测脚点的雷达测高卫星观测波形和辅助数据的方法包括以下步骤:
将测高卫星轨迹线向外扩充2千米形成轨迹放大区;
根据所述轨迹放大区和全球河网中心线相交得到河流虚拟站位置,根据水体频率图提取河流水面范围并计算河段平均宽度,并根据河段平均宽度筛选出平均河宽大于阈值的河流虚拟站位置以确定可重建水位的目标地点;
基于筛选出的河流虚拟站位置,提取出该河段水面范围,并将河流范围向外扩充1.5千米形成缓冲区,在所述卫星测高原始波形数据SGDRs中提取出所述缓冲区内的所有雷达测高卫星观测点波形和辅助数据。
3.根据权利要求2所述的河流虚拟站水位时序重建方法,其特征在于,对缓冲区内的所有雷达测高卫星观测波形和辅助数据按照观测时间分组,对每组数据分别进行波形分析和子波分解,获得多组纬向高程剖面的方法包括以下步骤:
选取组中每个观测时段所有观测脚点的中心点及其邻域观测点的波形合成平均波形,采用自适应峰值探测方法和逐步细致的参数策略提取所述平均波形中的脉冲峰值;
对该组观测中所有的波形数据,以最大能量、噪声和峰值水平为依据筛选单点波形数据,将筛选后的单点波形数据归一化后利用自适应脉冲峰值探测方法识别子波,并剔除微小的噪声子波,采用ICE-1跟踪算法确定每个子波对应的高程值;
对筛选后的每个单点波形,将其子波与平均波形中的脉冲峰值进行一对一匹配,并剔除脉冲峰值对应的高程差值的绝对值大于2米及重复选取的配对,获得多组纬向高程剖面。
4.根据权利要求3所述的河流虚拟站水位时序重建方法,其特征在于,选取组中每个观测时段所有观测脚点的中心点及其邻域观测点的波形合成平均波形的方法包括以下步骤:
对每一观测时段的卫星波形,按照河段宽度调整参与合成的波形个数n,n为奇数且不小于5,以观测脚点中心点及其左右邻域各(n-1)/2个点进行波形合成;
所述波形合成的方法为:计算参与合成的波形数据采样门的高程范围,将每个波形数据插值到该范围内,然后按照中值合成,获得平均波形;
按照河段宽度调整参与合成的波形个数的方法为:参与合成的波形个数为河段宽度除以卫星采样距离的最近的奇数,并且最小为5个波形合成。
5.根据权利要求3所述的河流虚拟站水位时序重建方法,其特征在于,采用自适应峰值探测方法和逐步细致的参数策略提取所述平均波形中的脉冲峰值的方法包括以下步骤:
对于平均波形,利用改进的自适应峰值自动探测算法识别脉冲峰值;
改进的自适应峰值自动探测算法的过滤参数为:脉冲峰值探测窗口大小为3个采样门宽度,脉冲峰值能量大于25%的最大能量,脉冲峰值上升幅度大于10%的最大能量,并且上升幅度大于下降幅度的35%;
如果提取出的脉冲峰值点对应的高程值的最小值小于前一观测时间的高程值,则减小窗口大小,降低脉冲峰值能量阈值,降低上升沿能量幅度阈值,再次检测脉冲峰值;
采用ICE-1算法确定跟踪门并计算脉冲峰值对应地面目标的高程。
6.根据权利要求3所述的河流虚拟站水位时序重建方法,其特征在于,将每个单点波形的子波与平均波形中的脉冲峰值进行一对一匹配的方法包括:采用脉冲峰值和跟踪高程两种标准配对每个单点波形的子波与平均波形中的脉冲峰值,并以最小标准差为标准选取最终的配对;通过领域距离分析保证每个子波只对应一组脉冲峰值。
7.根据权利要求3所述的河流虚拟站水位时序重建方法,其特征在于,以初始水面高程序列为基础,进行噪声识别获取噪声点,并从所述备用高程组中选择高程值替换噪声点的高程初始值的方法包括:以初始水面高程序列为基础,采用滑动窗口平均拟合和时间序列拟合进行噪声识别,从所述备用高程组中重新确定噪声点的目标高程;迭代停止的方式包括限制总迭代次数,是否有噪声和高程序列是否改变。
8.根据权利要求3所述的河流虚拟站水位时序重建方法,其特征在于,所述方法还包括:基于附近的实测水位数据,实现对所构建的虚拟站水位序列的精度分析。
9.一种河流虚拟站水位时序重建系统,其特征在于,所述系统包括:
获取模块:用于获取卫星观测数据;
提取模块:用于根据所述卫星观测数据确定河流虚拟站位置和水面范围,根据所述河流虚拟站位置和水面范围定义缓冲区,根据所述缓冲区的范围提取出该范围内的所有观测脚点的雷达测高卫星观测波形和辅助数据,观测脚点即卫星星下点位置;
高程剖面模块:用于对缓冲区内的所有雷达测高卫星观测波形和辅助数据按照观测时段分别进行波形分析和子波分解,获得多组纬向高程剖面;
水面高程模块:用于对所述多组纬向高程剖面,进行星下点偏移效应改正后,构建多目标优化函数,以损失函数最小的组纬向高程剖面的高程值作为初始水面高程序列,同时记录其它组纬向高程剖面的高程值作为备用高程组;
噪声剔除模块:用于以初始水面高程序列为基础,进行噪声识别获取噪声点,并从所述备用高程组中选择高程值替换噪声点的高程初始值,获得优化后的水位高程时间序列。
10.一种河流虚拟站水位时序重建系统,其特征在于,包括处理器及存储介质;
所述存储介质用于存储指令;
所述处理器用于根据所述指令进行操作以执行根据权利要求1~8任一项所述方法的步骤。
CN202110838483.2A 2021-07-23 2021-07-23 一种河流虚拟站水位时序重建方法及系统 Active CN113834547B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202110838483.2A CN113834547B (zh) 2021-07-23 2021-07-23 一种河流虚拟站水位时序重建方法及系统

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202110838483.2A CN113834547B (zh) 2021-07-23 2021-07-23 一种河流虚拟站水位时序重建方法及系统

Publications (2)

Publication Number Publication Date
CN113834547A CN113834547A (zh) 2021-12-24
CN113834547B true CN113834547B (zh) 2022-08-16

Family

ID=78962905

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202110838483.2A Active CN113834547B (zh) 2021-07-23 2021-07-23 一种河流虚拟站水位时序重建方法及系统

Country Status (1)

Country Link
CN (1) CN113834547B (zh)

Families Citing this family (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN114969664B (zh) * 2022-06-01 2023-03-21 广州市城市规划勘测设计研究院 一种水位改正方法、装置、设备及介质
CN116719025B (zh) * 2023-08-09 2023-10-13 中国科学院地理科学与资源研究所 一种内陆水面高程测量方法

Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN107167786A (zh) * 2017-06-05 2017-09-15 中国测绘科学研究院 卫星激光测高数据辅助提取高程控制点方法
CN111192282A (zh) * 2019-12-19 2020-05-22 中国科学院南京地理与湖泊研究所 湖滨带虚拟站的湖库时序水位重建方法
CN112414510A (zh) * 2020-10-10 2021-02-26 武汉大学 卫星测高河流水位监测方法及系统
CN112697232A (zh) * 2020-12-15 2021-04-23 自然资源部国土卫星遥感应用中心 基于多源卫星测高数据的水位测量及变化监测方法及装置

Family Cites Families (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN109489637B (zh) * 2018-11-08 2019-10-18 清华大学 水量变化监测方法、装置、计算机设备和存储介质

Patent Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN107167786A (zh) * 2017-06-05 2017-09-15 中国测绘科学研究院 卫星激光测高数据辅助提取高程控制点方法
CN111192282A (zh) * 2019-12-19 2020-05-22 中国科学院南京地理与湖泊研究所 湖滨带虚拟站的湖库时序水位重建方法
CN112414510A (zh) * 2020-10-10 2021-02-26 武汉大学 卫星测高河流水位监测方法及系统
CN112697232A (zh) * 2020-12-15 2021-04-23 自然资源部国土卫星遥感应用中心 基于多源卫星测高数据的水位测量及变化监测方法及装置

Non-Patent Citations (1)

* Cited by examiner, † Cited by third party
Title
《基于ICESat/CryoSat-2卫星测高及站点观测的纳木错湖水位趋势变化监测》;宋春桥等;《科学通报》;20151231;第60卷(第21期);1 *

Also Published As

Publication number Publication date
CN113834547A (zh) 2021-12-24

Similar Documents

Publication Publication Date Title
CN110427857B (zh) 一种基于遥感数据融合的输电线路地质灾害分析方法
CN113834547B (zh) 一种河流虚拟站水位时序重建方法及系统
Trevisani et al. Surface texture analysis of a high-resolution DTM: Interpreting an alpine basin
CN103363962A (zh) 一种基于多光谱影像的湖泊水储量遥感估算方法
Miao et al. Developing efficient procedures for automated sinkhole extraction from Lidar DEMs
CN107247927B (zh) 一种基于缨帽变换的遥感图像海岸线信息提取方法及系统
CN104239884A (zh) 一种基于遥感植被指数时间序列的异常淹没区域检测方法
CN103236063A (zh) 基于多尺度谱聚类及决策级融合的sar图像溢油检测方法
CN103106658A (zh) 一种海岛、礁岸线快速提取方法
Ogashawara et al. The use of optical remote sensing for mapping flooded areas
CN109992635B (zh) 一种震后泥石流早期识别方法
CN112418506B (zh) 基于机器学习的海岸带湿地生态安全格局优化方法、装置
CN108537116B (zh) 一种基于多尺度特征的海岸线二级类型提取方法及系统
Ghaderi et al. Detecting shoreline change employing remote sensing images (Case study: Beris Port-east of Chabahar, Iran)
CN110569733B (zh) 基于遥感大数据平台的湖泊长时序连续水域变化重建方法
Thissen Automating Surface Water Detection for Rivers
Athanasiou et al. Global Coastal Characteristics (GCC): A global dataset of geophysical, hydrodynamic, and socioeconomic coastal indicators
Maze et al. Profile classification models
CN114782211B (zh) 一种海山分布范围信息的获取方法及系统
CN116882731A (zh) 一种基于斜坡单元的地质灾害危险性评估方法及系统
Yu et al. Digital terrain model extraction from airborne LiDAR data in complex mining area
Wang et al. Deriving natural coastlines using multiple satellite remote sensing images
Schmidt et al. Monitoring concepts for coastal areas using lidar data
Novák et al. The Potential and Implications of Automated Pre-Processing of LiDAR-Based Digital Elevation Models for Large-Scale Archaeological Landscape Analysis
Li et al. On the capacity of ICESat-2 laser altimetry for river level retrieval: An investigation in the Ohio River basin

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant