CN116758049A - 一种基于主被动卫星遥感的城市洪涝立体监测方法 - Google Patents

Abstract

本发明涉及洪涝监测技术领域，具体公开了一种基于主被动卫星遥感的城市洪涝立体监测方法，包括获取Sentinel‑2遥感影像和ICESat‑2卫星测高数据，对Sentinel‑2遥感影像进行云去除，对ICESat‑2卫星测高数据进行光子去噪，结合洪水轮廓数据及洪水水位数据分析洪水体积实现洪水的立体检测。本发明技术方案通过获取Sentinel‑2遥感影像和ICESat‑2卫星测高数据，并且基于Otsu算法提取水体边界并分析洪水淹没范围、基于ICESat‑2卫星测高数据估算了洪水水位高度，并结合DEM数据分析了洪水淹没水深，最后结合洪水范围数据及洪水水位数据分析洪水体积，从而实现洪水的立体检测，解决现有技术中的洪水监测手段时间长、范围小的问题。

Description

一种基于主被动卫星遥感的城市洪涝立体监测方法

技术领域

本发明涉及洪涝监测技术领域，具体为一种基于主被动卫星遥感的城市洪涝立体监测方法。

背景技术

目前洪水立体监测是指通过洪水淹没范围监测和洪水水位监测来获取全面的洪水信息。卫星遥感技术由于其广阔的覆盖范围、短响应时间和强大的时效性，在洪水淹没范围监测中得到了广泛应用。提取水体信息是洪水范围监测的关键。目前常用的方法包括单波段法、水指数法和光谱关系法。单波段法依赖于遥感图像在近红外波段的强吸收特性，并通过设定阈值来提取水体信息。然而，这种方法只适用于提取特征简单的地区的水体信息，在复杂地区的结果通常包含混乱的信息，提取精度较低。研究者发现，使用水指数法和光谱关系法等多波段组合可以更准确地区分水体和背景，从而获得更好的水体提取结果。其中，水指数法是效率高且精度高的水体提取方法。

归一化差异植被指数(NDVI)的基础上提出了归一化差异水指数(NDWI)，该指数通过比较绿色波段和近红外波段来强调水体信息。他设置了一个0的阈值，因此当像素值大于0时被认为是水体。然而，不同传感器在不同地区和不同时间拍摄的多光谱图像具有不同的特征，仅使用固定阈值无法有效提取许多地区的水体信息。

传统的基于模型的方法可以模拟洪水淹没过程并提供洪水水位信息。通过结合CLUE-S、SCS和等体积淹没算法成功模拟了深圳茅洲河流域的暴雨入渗-径流汇流-淹没过程。然而，基于模型的洪水测量方法依赖于有限的测量点分布，很难快速获得大规模、高精度的洪水水位测量信息。光学卫星可以提供大范围的洪水水位测量结果，但由于天气和水质条件的限制，精度较低，因此，我们提出一种基于主被动卫星遥感的城市洪涝立体监测方法。

发明内容

本发明的目的在于提供一种基于主被动卫星遥感的城市洪涝立体监测方法，解决了背景技术中所提出的问题。

为实现上述目的，本发明提供如下技术方案：一种基于主被动卫星遥感的城市洪涝立体监测方法，包括以下步骤：

S1:获取Sentinel-2遥感影像和ICESat-2卫星测高数据；

S2:对Sentinel-2遥感影像进行云去除，对ICESat-2卫星测高数据进行光子去噪；

S3:对云去除后的Sentinel-2遥感影像进行归一化差分水面指数计算，归一化差分水面指数为NDWSI，应用Otsu算法确定的最优分割阈值T，应用Otsu算法提取洪水灾害后研究区域的水体信息，并与土地利用分类数据中的原始水体边界进行比较，确定洪水轮廓数据；

S4：进一步对光子去噪后的ICESat-2卫星测高数据提取水面光子信息，利用粗略提取和高斯函数拟合后的保留光子的高程分布的直方图，得到了水面光子的高程分布范围，提取洪水灾害前后的水面信号光子，可以获得研究区域的水位数据，进一步分析水位变化和洪水水位信息，利用Stumpf经验模型，利用蓝绿波段反射率的差值来预测水深，获得洪水水位数据；

S5：结合洪水轮廓数据及洪水水位数据分析洪水体积，实现洪水的立体检测。

作为本申请技术方案的一种优选实施方式，步骤S5中，将洪水轮廓数据、洪水水位数据与DEM数据相结合，计算研究区域的蓄洪量，研究区域蓄洪量的计算方法采用以下步骤：根据单个湖区的DEM数据，计算出洪水淹没深度，计算洪水前后的水位变化的公式为：

ΔH＝H_n-H

式中，H_n为洪水淹没深度，H为地形高度，ΔH为洪水前后的水位变化；

研究区域蓄洪量的计算方法为：研究区域蓄洪量＝洪水面积*各湖区的洪水蓄洪深度。

作为本申请技术方案的一种优选实施方式，步骤S3中计算NDWSI的步骤为：通过使用绿色和近红外波段计算NDWSI，可以突出显示图像中的水体信息，从而揭示了水与陆地的分离，NDWSI的计算公式如下：

式中，Green指遥感影像中的绿波段，Nir指遥感影像中的近红外波段，应用Otsu算法提取水体信息，NDWSI图像像素值范围为0到L-1，L为灰度值的个数，所使用的类间方差公式如下：

式中，i和j表示NDWSI图像的行列数；Sum_land和Sum_water分别表示整个NDWSI图像中陆地像素和水体像素的比例；μ_land和μ_water分别表示目标水体和背景像素的平均灰度值；σ(NDWSI)为类间方差最大的分割阈值；

应用Otsu算法确定最优分割阈值T的计算公式如下：

作为本申请技术方案的一种优选实施方式，步骤S4中，利用粗略提取和高斯函数拟合后的保留光子的高程分布的直方图并且得到水面光子的高程分布范围的详细执行步骤如下：

首先，对原始光子数据进行分割，并利用各段光子的直方图统计，得高程分布直方图。计算所有高程切片中的平均光子数，以确定在水覆盖区域中光子的高程分布范围；

其次，对水覆盖区域的光子数据进行直方图统计，得到最大直方图HN_max及其相应的高程值H_maxN，从高程值HN_max/2较小的一侧搜索直方图的最大位置，直到找到该位置对应的高程值H_opt为止，采用下列公式计算高斯函数拟合的初始参数，

公式为：

并且绘制高斯分布拟合曲线进行直方图拟合，得到均值μ_fit和标准差σ_fit，最后，将间隔[μ_fit-3σ_fit，μ_fit+3σ_fit]中的光子视为水面信号光子；

通过提取洪水灾害前后的水面信号光子，可以获得研究区域的水位数据，进一步分析水位变化和洪水水位信息，利用Stumpf经验模型，利用蓝绿波段反射率的差值来预测水深，获得水位数据，结合实际水位数据，可以得到更准确的水位测量结果，Stumpf经验模型为:

水位高度用H表示，来自于云去除处理后的Sentinel-2图像,利用经验参数m₁和m₀将计算结果转化为水位高度,R_rsi和R_rsj分别表示蓝色和绿色波段的遥感反射率。

与现有技术相比，本发明的有益效果如下：

本申请技术方案通过获取Sentinel-2遥感影像和ICESat-2卫星测高数据，并且基于Otsu算法提取水体边界并分析洪水淹没范围、基于ICESat-2卫星测高数据估算了洪水水位高度，并结合DEM数据分析了洪水淹没水深，最后结合洪水范围数据及洪水水位数据分析洪水体积，从而实现洪水的立体检测，解决现有技术中的洪水监测手段时间长、范围小的问题。

附图说明

通过阅读参照以下附图对非限制性实施例所作的详细描述，本发明的其它特征、目的和优点将会变得更明显：

图1为本发明一种基于主被动卫星遥感的城市洪涝立体监测方法的整体技术路线图；

图2为本发明一种基于主被动卫星遥感的城市洪涝立体监测方法的局部技术路线图；

图3为本发明水体边界提取方法的流程示意图。

具体实施方式

下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本发明一部分实施例，而不是全部的实施例。

实施例1，如图1-3所示，本发明提供一种技术方案：一种基于主被动卫星遥感的城市洪涝立体监测方法，包括以下步骤：

S1:获取Sentinel-2遥感影像和ICESat-2卫星测高数据；

S2:对Sentinel-2遥感影像进行云去除，对ICESat-2卫星测高数据进行光子去噪；

S3:对云去除后的Sentinel-2遥感影像进行归一化差分水面指数计算，归一化差分水面指数为NDWSI，应用Otsu算法确定的最优分割阈值T，应用Otsu算法提取洪水灾害后研究区域的水体信息，并与土地利用分类数据中的原始水体边界进行比较，确定洪水轮廓数据；

S4：进一步对光子去噪后的ICESat-2卫星测高数据提取水面光子信息，利用粗略提取和高斯函数拟合后的保留光子的高程分布的直方图，得到了水面光子的高程分布范围，提取洪水灾害前后的水面信号光子，可以获得研究区域的水位数据，进一步分析水位变化和洪水水位信息，利用Stumpf经验模型，利用蓝绿波段反射率的差值来预测水深，获得洪水水位数据；

S5：结合洪水轮廓数据及洪水水位数据分析洪水体积，实现洪水的立体检测。

在本发明的一个具体实施例中，Sentinel-2是欧洲空间局开发的一组卫星，用于进行地球观测和环境监测，Sentinel-2卫星通过收集多光谱数据提供了许多有用的信息，ICESat-2是美国国家航空航天局发射的一颗激光雷达卫星，ICESat-2卫星通过激光高度计仪器提供测高数据；对Sentinel-2遥感影像进行云去除采用多光谱归一化差异植被指数阈值法，属于比较常见的预处理，处理步骤如下，首先计算NDVI图像，其中NIR代表近红外波段，Red代表红色波段。然后设置一个合适的阈值，将NDVI图像中低于该阈值的像素标记为可能为云的区域；ICESat-2卫星的测高数据是通过激光高度计仪器(ATLAS)获取的，光子去噪是为了减少测量误差和提高数据质量，为了降低数据处理的计算复杂度，有必要粗略地提取光子点云数据，以确定噪声和信号光子的近似区间。在ATL03数据中，每个光子都在官方网站上进行标记，通过删除置信度最低的可疑光子，可以粗略地提取点云数据，Otsu算法是一种用于图像分割的自适应阈值选择方法，该算法通过寻找图像灰度级的一个阈值，将图像分为前景和背景两个部分，Otsu算法是基于图像的灰度直方图来确定最佳阈值，它通过最小化类间方差或最大化类内方差的方式来选择阈值，类间方差反映了前景和背景之间的灰度差异，而类内方差衡量了同一区域内各像素值的灰度相似性。

在优选技术方案中，步骤S5中，将洪水轮廓数据、洪水水位数据与DEM数据相结合，计算研究区域的蓄洪量，洪水区域主要位于平原区，因此洪水区底部相对平坦，洪水截面相对垂直，因此，可以假定一个洪水的体积近似等于一个棱镜的体积，研究区域蓄洪量的计算方法采用以下步骤：根据单个湖区的DEM数据，计算出洪水淹没深度，计算洪水前后的水位变化的公式为：

ΔH＝H_n-H

式中，H_n为洪水淹没深度，H为地形高度，ΔH为洪水前后的水位变化；

研究区域蓄洪量的计算方法为：研究区域蓄洪量＝洪水面积*各湖区的洪水蓄洪深度。

在优选技术方案中，步骤S3中计算NDWSI的步骤为：通过使用绿色和近红外波段计算NDWSI，可以突出显示图像中的水体信息，从而揭示了水与陆地的分离，NDWSI的计算公式如下：

式中，Green指遥感影像中的绿波段，Nir指遥感影像中的近红外波段，应用Otsu算法提取水体信息，NDWSI图像像素值范围为0到L-1，L为灰度值的个数，所使用的类间方差公式如下：

式中，i和j表示NDWSI图像的行列数；Sum_land和Sum_water分别表示整个NDWSI图像中陆地像素和水体像素的比例；μ_land和μ_water分别表示目标水体和背景像素的平均灰度值；σ(NDWSI)为类间方差最大的分割阈值；

应用Otsu算法确定最优分割阈值T的计算公式如下：

在优选技术方案中，步骤S4中，利用粗略提取和高斯函数拟合后的保留光子的高程分布的直方图并且得到水面光子的高程分布范围的详细执行步骤如下：

首先，对原始光子数据进行分割，并利用各段光子的直方图统计，得高程分布直方图。计算所有高程切片中的平均光子数，以确定在水覆盖区域中光子的高程分布范围；

其次，对水覆盖区域的光子数据进行直方图统计，得到最大直方图HN_max及其相应的高程值H_maxN，从高程值HN_max/2较小的一侧搜索直方图的最大位置，直到找到该位置对应的高程值H_opt为止，采用下列公式计算高斯函数拟合的初始参数，

公式为：

并且绘制高斯分布拟合曲线进行直方图拟合，得到均值μ_fit和标准差σ_fit，最后，将间隔[μ_fit-3σ_fit+3σ_fit]中的光子视为水面信号光子；

通过提取洪水灾害前后的水面信号光子，可以获得研究区域的水位数据，进一步分析水位变化和洪水水位信息，利用Stumpf经验模型，利用蓝绿波段反射率的差值来预测水深，获得水位数据，结合实际水位数据，可以得到更准确的水位测量结果，Stumpf经验模型为:

水位高度用H表示，来自于云去除处理后的Sentinel-2图像,利用经验参数m₁和m₀将计算结果转化为水位高度,R_rsi和R_rsj分别表示蓝色和绿色波段的遥感反射率。

综上所述，本发明提供了一种基于主被动卫星遥感的城市洪涝立体监测方法，通过获取Sentinel-2遥感影像和ICESat-2卫星测高数据，并且基于Otsu算法提取水体边界并分析洪水淹没范围、基于ICESat-2卫星测高数据估算了洪水水位高度，并结合DEM数据分析了洪水淹没水深，最后结合洪水范围数据及洪水水位数据分析洪水体积，从而实现洪水的立体检测，解决现有技术中的洪水监测手段时间长、范围小的问题。

Claims (4)

1.一种基于主被动卫星遥感的城市洪涝立体监测方法，其特征在于：包括以下步骤：

S1:获取Sentinel-2遥感影像和ICESat-2卫星测高数据；

S2:对Sentinel-2遥感影像进行云去除，对ICESat-2卫星测高数据进行光子去噪；

S3:对云去除后的Sentinel-2遥感影像进行归一化差分水面指数计算，归一化差分水面指数为NDWSI，应用Otsu算法确定的最优分割阈值T，应用Otsu算法提取洪水灾害后研究区域的水体信息，并与土地利用分类数据中的原始水体边界进行比较，确定洪水轮廓数据；

S4：进一步对光子去噪后的ICESat-2卫星测高数据提取水面光子信息，利用粗略提取和高斯函数拟合后的保留光子的高程分布的直方图，得到了水面光子的高程分布范围，提取洪水灾害前后的水面信号光子，可以获得研究区域的水位数据，进一步分析水位变化和洪水水位信息，利用Stumpf经验模型，利用蓝绿波段反射率的差值来预测水深，获得洪水水位数据；

S5：结合洪水轮廓数据及洪水水位数据分析洪水体积，实现洪水的立体检测。

2.根据权利要求1所述的一种基于主被动卫星遥感的城市洪涝立体监测方法，其特征在于：所述步骤S5中，将洪水轮廓数据、洪水水位数据与DEM数据相结合，计算研究区域的蓄洪量，研究区域蓄洪量的计算方法采用以下步骤：根据单个湖区的DEM数据，计算出洪水淹没深度，计算洪水前后的水位变化的公式为：

ΔH＝H_n-H

式中，H_n为洪水淹没深度，H为地形高度，ΔH为洪水前后的水位变化；

研究区域蓄洪量的计算方法为：研究区域蓄洪量＝洪水面积*各湖区的洪水蓄洪深度。

3.根据权利要求1所述的一种基于主被动卫星遥感的城市洪涝立体监测方法，其特征在于：所述步骤S3中计算NDWSI的步骤为：通过使用绿色和近红外波段计算NDWSI，可以突出显示图像中的水体信息，从而揭示了水与陆地的分离，NDWSI的计算公式如下：

式中，Green指遥感影像中的绿波段，Nir指遥感影像中的近红外波段，应用Otsu算法提取水体信息，NDWSI图像像素值范围为0到L-1，L为灰度值的个数，所使用的类间方差公式如下：

式中，i和j表示NDWSI图像的行列数；Sum_land和Sum_water分别表示整个NDWSI图像中陆地像素和水体像素的比例；μ_land和μ_water分别表示目标水体和背景像素的平均灰度值；σ(NDWSI)为类间方差最大的分割阈值；

应用Otsu算法确定最优分割阈值T的计算公式如下：

4.根据权利要求1所述的一种基于主被动卫星遥感的城市洪涝立体监测方法，其特征在于：所述步骤S4中，利用粗略提取和高斯函数拟合后的保留光子的高程分布的直方图并且得到水面光子的高程分布范围的详细执行步骤如下：

首先，对原始光子数据进行分割，并利用各段光子的直方图统计，得高程分布直方图。计算所有高程切片中的平均光子数，以确定在水覆盖区域中光子的高程分布范围；

其次，对水覆盖区域的光子数据进行直方图统计，得到最大直方图HN_max及其相应的高程值H_maxN，从高程值HN_max/2较小的一侧搜索直方图的最大位置，直到找到该位置对应的高程值H_opt为止，采用下列公式计算高斯函数拟合的初始参数，

公式为：

并且绘制高斯分布拟合曲线进行直方图拟合，得到均值μ_fit和标准差σ_fit，最后，将间隔[μ_fit-3σ_fit，μ_fit+3σ_fit]中的光子视为水面信号光子；

通过提取洪水灾害前后的水面信号光子，可以获得研究区域的水位数据，进一步分析水位变化和洪水水位信息，利用Stumpf经验模型，利用蓝绿波段反射率的差值来预测水深，获得水位数据，结合实际水位数据，可以得到更准确的水位测量结果，Stumpf经验模型为:

水位高度用H表示，来自于云去除处理后的Sentinel-2图像,利用经验参数m₁和m₀将计算结果转化为水位高度,R_rsi和R_rsj分别表示蓝色和绿色波段的遥感反射率。