一种库容监测方法

摘要

本发明涉及一种库容监测方法，包括以下步骤：对待测水体上方的遥感影像进行正射纠正，然后根据遥感影像提取待测水体边界，得到水体边界多边形；将水体边界多边形与地形DEM叠加，提取水体边界多边形与地形DEM的相交点作为三维水边界离散点；对三维水边界离散点进行插值，获得水面DEM；求取水面DEM与地形DEM之间包容的容积，即为待测水体的库容监测结果。本发明利用遥感和地理信息系统等空间信息技术，能准确监测水面为曲面，包括动库容在内的总库容。并且提出了基于河流分段的优化方式，创建缓冲区将计算控制在有效区域内，精度高效率高，具有业务化推广应用的前景。

权利要求

1.一种库容监测方法，其特征在于：包括以下步骤，

步骤1，对待测水体上方的遥感影像进行正射纠正，然后根据遥感影像提取待测水体边界，得到水体边界多边形；

步骤2，将水体边界多边形与地形DEM叠加，提取水体边界多边形与地形DEM的相交点作为三维水边界离散点；

步骤3，对三维水边界离散点进行插值，获得水面DEM；

步骤4，求取水面DEM与地形DEM之间包容的容积，即为待测水体的库容监测结果；所述水面DEM与地形DEM之间包容的最大容积求取方式为，水面DEM的空间参考、像素尺寸和分辨率属性与地形DEM保持一致，将水面DEM上每个像素分别与地形DEM对应位置上像素之间的空间作为正四棱柱处理，然后对所有正四棱柱体积求和；按以下公式对所有正四棱柱体积求和，式中，HWij是水面DEM在位置(i，j)处像素的高程值，HTij是地形DEM在位置(i，j)处像素的高程值，Scell是单个像素的面积，m、n分别为水面DEM长和宽像素个数；

对正四棱柱体积求和时进行分段简化处理，包括以下步骤，

步骤4.1，利用地形DEM进行流域分割，提取待测水体的水系线；

步骤4.2，参照待测水体的水面比降，对待测水体分段；

步骤4.3，根据流域分割所得待测水体宽度，以水系线为轴线，在每段待测水体两侧创建缓冲区；

步骤4.4，对每段待测水体的缓冲区，提取落入其内的三维水边界离散点所对应的DEM栅格单元高程值，并求平均值作为该段待测水体的平均水位，根据平均水位确定该段待测水体的正四棱柱体高度；

步骤4.5，按以下公式对所有正四棱柱体积求和，

式中，Hwk是水面DEM在第k段处的平均水位，HTij是地形DEM在位置(i，j)处像素的高程值，Scell是单个像素的面积，m、n分别为水面DEM长和宽像素个数，s为待测水体的分段数。

2.根据权利要求1所述的库容监测方法，其特征在于：缓冲区宽度大于待测水体宽度的1/2。

3.根据权利要求1或2所述的库容监测方法，其特征在于：所述遥感影像采用SPOT-5影像资料。

说明书

一种库容监测方法

技术领域

[0001] 本发明属于空间信息技术在水利工程学领域应用技术，特别涉及一种库容监测方法。

背景技术

[0002] 一直以来湖泊和水库蓄水量的研究一直是陆地水循环研究中的一个重要方面，然而就全球范围而言，给予陆面监测网络的传统方法对其研究和把握能力十分有限，早已不能满足实际应用中数据精度和更新速度的要求。随着空间信息技术的日益成熟和完善，国内外研究人员越来越多的将这项技术应用于水利水电建设中的各项领域。

[0003] 遥感信息具有周期性、宏观性和实用性等特点，因此，运用遥感技术进行水库水资源监测，有着快速、同步对比和实时性等优势。特别是在汛期，利用遥感手段能够非常迅速地对水量进行风险评估，以做出果断准确的决策。而且对于提取水体方面的研究开展得比较早，方法很多，技术也十分成熟。利用遥感图像，进行水库动库容的定量测算，是水库水资源信息研究的新途径，也是水库本身动态演化过程监测新方法。

[0004] 在库容监测方面已有大量研究成果。国外近期比较有代表性的有，美国学者Jordan和Bamey利用水文学方法来测定水库容量，并作长时间序列分析水量损失的原因。Panayiotis提出了利用三次样条插值生成等高线的方法来估算库容量。另外，还有不少学者研究了水库库容在水资源配置、水利发电中的应用。相关文献有：

[0005] Jordan Furnans，and Barney Austin.Hydrographic Survey Method forDetermining Reservoir Volume[J].Environmental Modelling&Software，2008，23:139-146.

[0006] Panayiotis Foteinopoulos.Cubic Spline Interpolation to Develop Contours ofLarge Reservoirs and Evaluate area and Volume[J].Advances in EngineeringSoftware，2008，doi：10.1016/j.advengsoft.2008.03.005.[0007] Nien-Sheng Hsu，Wei-Chen Cheng，Wen-Ming Cheng，et al.Optimization andCapacity Expansion of a Water Distribution System[J].Advances in Water Resources，2008，31:776-786.

[0008] S.M.H.Hosseini，F.Forouzbakhsh，M.Fotouhi，et al.Determination ofInstallation Capacity in Reservoir Hydro-power Plants Considering Technical，Economical and Reliability indices[J].Electrical Power and Energy Systems，2008，30：393-402.

[0009] Thomas A.Mcmahon，Geoffrey G.S.Pegram，Richard M.Vogel，et al.Revisiting Reservoir Storage-yield Relationships Using a Global StreamflowDatabase[J].Advances in Water Resources，2007，30：1858-1872.

[0010] 国内学者对如何监测库容也做了大量研究，这些研究大致可以概括为三方面：(1)将静库容作为水体的总库容。一般有横断面法、等高线法、三角网格法以及方格网法等，这些库容监测方法都是按照微积分学原理，将河道水库进行切割、分块，然后求和得到静库容。它们都需要给定水位的监测值，并将水面看成水平，计算的总库容就是水平水面以下包容的最大体积；(2)利用水文水力学模型来监测包括动库容、瞬时库容在内的总库容。当考虑动库容时，水库的库容定义为：最远回水断面至坝前的总蓄水量。以水库为例，如图1所示，①标识水库大坝；②标识河道底坡；③标识坝前水深H时水面线。根据库容定义，坝前水深H所对应的库容在纵向上表示为封闭多边形ANDEFBC对应的体积。于是，总库容包括静库容Ws(ABC)和动库容Wd(NDEFB)两部分之和，而动库容又可以分为河道槽蓄动库容Wr(EFBD)与库区回水引起的附加楔蓄动库容Wu(BDN)两部分，其中，静库容Ws是坝前水深H的函数，而动库容Wd需要复杂的水文水动力模型求解；(3)也有探讨了空间信息技术运用于库容监测的。

[0011] 相关文献有：

[0012] 米鸿燕，宰建，蒋兴华.静库容计算方法的比较分析[J].地矿测绘，2007，23(2)：1-4.

[0013] 许海军，陈守煜.水库动库容调洪计算的数值-解析解法[J].水利学报，2002，3：69-73.

[0014] 秦惠承.动态库容调洪的数值计算[J].水利学报，1983，1：42-50.[0015] 韩爱民，武淑华，高军，等.用数字地图计算洪泽湖库容等特征参数的方法初探[J].水文，2001，21(5):35-37.

[0016] 袁勇.基于DEM库容计算及可视化研究[D].武汉:武汉大学硕士学位论文，2004.[0017] 吕杰堂，王治华，周成虎.西藏易贡滑坡堰塞湖的卫星遥感监测方法初探[J].地球学报，2002，23(4):363-368.

[0018] 田雨，梁勇，林宗坚，等.遥感技术在雪野水库库容测算中的应用[J].人民黄河，2007，29(4):76-77.

[0019] 曹波.基于遥感图像和DEM测定水库动库容的方法研究[D].武汉：华中科技大学硕士学位论文，2006.

[0020] 长期以来，动库容的研究方法大多依赖于简化复杂的实际河槽形态，即简化成典型的河槽断面处理，然后，依据水力学方法来计算。但是，由于所依赖的许多边界条件都有经验性、实验性和近似性，并且计算所需要的各种水文、地质资料，往往要耗费大量的人力物力去搜集，在遇到大灾时更加难于获取，致使计算结果常常难以令人满意。另外，水力学水文学模型的求解需要长时间的迭代运算，这些模型往往都是计算密集型的，对计算资源要求较高，对快速业务化库容应急监测是一个较大的挑战。引入空间信息技术来改变传统繁琐的资料搜集方式，以较为简捷的方式快速精确的计算库容曲线，是解决动库容监测问题的一种新思路。

[0021] 综合以上列举的文献资料，已有不少学者针对空间信息技术在库容监测中的应用进行了一系列研究，有些学者利用数字地图或DEM估算库容，这些方法只考虑了水面水平的情况，监测到的是静库容；有学者利用遥感影像提取水体面积，研究水面面积与库容之间的关系，但是没有通过遥感影像得到库容具体数值；也有的将动库容的槽蓄动库容和楔蓄动库容分别近似成平行四边形和三角形空间体来计算，该方法依赖于河床形状。因此，水利工程学领域缺乏能够实时监测库容的有效技术。

发明内容

[0022] 本发明目的在于提供一种库容监测方法，以高效监测山区具有较高水位比降的湖泊河流库容，为防洪泻流提供技术支持。

[0023] 本发明的技术方案包括以下步骤：

[0024] 步骤1，对待测水体上方的遥感影像进行正射纠正，然后根据遥感影像提取待测水体边界，得到水体边界多边形；

[0025] 步骤2，将水体边界多边形与地形DEM叠加，提取水体边界多边形与地形DEM的相交点作为三维水边界离散点；

[0026] 步骤3，对三维水边界离散点进行插值，获得水面DEM；

[0027] 步骤4，求取水面DEM与地形DEM之间包容的容积，即为待测水体的库容监测结果。

[0028] 而且，所述水面DEM与地形DEM之间包容的最大容积求取方式为，水面DEM的空间参考、像素尺寸和分辨率属性与地形DEM保持一致，将水面DEM上每个像素分别与地形DEM对应位置上像素之间的空间作为正四棱柱处理，然后对所有正四棱柱体积求和。

[0029] 而且，按以下公式对所有正四棱柱体积求和，

[0030]

[0031] 式中，HWij是水面DEM在位置(i，j)处像素的高程值，HTij是地形DEM在位置(i，j)处像素的高程值，Scell是单个像素的面积，m、n分别为水面DEM长和宽像素个数。

[0032] 而且，对正四棱柱体积求和时进行分段简化处理，包括以下步骤，[0033] 步骤4.1，利用地形DEM进行流域分割，提取待测水体的水系线；

[0034] 步骤4.2，参照待测水体的水面比降，对待测水体分段；

[0035] 步骤4.3，根据流域分割所得待测水体宽度，以水系线为轴线，在每段待测水体两侧创建缓冲区；

[0036] 步骤4.4，对每段待测水体的缓冲区，提取落入其内的三维水边界离散点所对应的DEM栅格单元高程值，并求平均值作为该段待测水体的平均水位，根据平均水位确定该段待测水体的正四棱柱体高度；

[0037] 步骤4.5，按以下公式对所有正四棱柱体积求和，

[0038]

[0039] 式中，Hwk是水面DEM在第k段处的平均水位，HTij是地形DEM在位置(i，j)处像素的高程值，Scell是单个像素的面积，m、n分别为水面DEM长和宽像素个数，s为河流的分段数。

[0040] 而且，缓冲区宽度大于待测水体宽度的1/2。

[0041] 而且，所述遥感影像采用SPOT-5影像资料。

[0042] 本发明提出利用遥感和地理信息技术相结合，不需要监测的水位信息，直接在遥感影像上提取水体边界，通过对应地形DEM提取的高程值内插水表面的方法，建立了水面DEM，监测结果既考虑了动库容对总库容的贡献，又不受水体底部的地形约束，实现了库容动态实时监测。并且提出了基于河流分段的优化方式，创建缓冲区将计算控制在有效区域内，获得了提高库容监测的速度和精度的技术效果。

附图说明

[0043] 图1为河流水库动库容剖面示意图；

[0044] 图2为本发明实施例的库容监测流程示意图；

[0045] 图3为本发明实施例的分段简化处理示意图。

具体实施方式

[0046] 以下基于具体实施例详述本发明技术方案：

[0047] 步骤1，对待测水体上方的遥感影像进行正射纠正，然后根据遥感影像提取待测水体边界，得到水体边界多边形。

[0048] 将遥感影像进行正射纠正，就是利用控制点改正原始图像的几何形变，产生一幅符合某种地图投影或图件表达要求的新图像。它将图上存在的各种误差和畸变作为一个整体，利用一定数量的控制点建立一个数学模型，用该模型来建立待校正图像与标准图像空间的某种影射变换或对应变换关系，再依据这种变换函数将待校正图像空间的所有像元一一变换到标准空间中去。通过正射纠正后，影像与地形DEM可以无缝叠加在一起，根据遥感影像提取湖区水体边界，得到水体边界的多边形；本发明实施例采用武汉大学研制的星载光学和SAR处理软件对研究区域灾后SPOT-5影像资料进行基于RPC模型的正射纠正；地形DEM直接采用国家基础地理信息数据即可，这种数字高程模型的建立通常由国家花费大量人类物力量测构建，是可靠的数据来源。所需的地形DEM也可先进行流域分割，提取出待测水体所在地理位置的对应部分，包括水体和部分边岸，本发明简称为河床DEM。待测水体可能是河流湖泊或水库，而本发明的目标主要是监测山区具有较高水位比降的湖泊河流库容。

[0049] 本发明所提供的实施例所选区域为汶川地震后的唐家山堰塞湖区，所用的遥感影像为2008年5月16日获取的SPOT-5影像，1：5万的地形DEM，分辨率为25米。在计算库容之前需对数据进行预处理，首先对研究区域的SPOT-5影像进行正射纠正，然后，根据提取堰塞湖区回水线的水体边界多边形。回水线提取的精度取决于影像的空间分辨率和几何纠正精度，因此，遥感影像几何纠正是堰塞湖库容量测的基础，其精度对结果的准确性至关重要。本发明所采用地形DEM的空间分辨率为25米，几何纠正精度在DEM一个格网以内，能够满足数据精度要求。

[0050] 步骤2，将水体边界多边形与地形DEM叠加，提取水体边界多边形与地形DEM的相交点作为三维水边界离散点。

[0051] 三维水边界离散点上附着该点对应的地形DEM像素的高程值，因此这些点具有高程属性，也即是水体边界的水位值。在求交之前需要对提取的水边界多边形与地形DEM进行精配准，否则叠加获取的离散点无法正确表达真实水边界水位高低；

[0052] 步骤3，对三维水边界离散点进行插值，获得水面DEM。

[0053] 具体实施时可以采用普通Kriging插值方法，根据三维水边界离散点的高程值生成三维水面，即水面DEM。该方法对于稀疏点内插具有很好的效果。在插值过程中可以设置结果的像素大小，以及剪切边界，以便将插值结果的栅格图像限制在水体范围之内，来减少后面计算量。插值时可以设定参数，使得插值所得水面DEM与地形DEM具有相同的空间参考和分辨率等属性

[0054] 步骤4，求取水面DEM与地形DEM之间包容的容积，即为待测水体的库容监测结果。

[0055] 水面DEM与地形DEM叠加将形成一个封闭的几何体，该几何体所包容的最大体积即使所求的总库容。本发明提供了进一步技术方案，基于像素实现库容监测，以便能够提高求取效率：根据积分原理，将该几何体中的水面DEM像素与地形DEM中对应位置上的像素分别作为上下底面，形成一个正四棱柱(底面为正方形的直棱柱)，水柱高度为水面DEM像素值与地形DEM对应位置像素值之差。于是，整个几何体被划分为若干底面积相等，水柱高度不同的四棱柱。库容监测过程中，首先以地形DEM为基准，针对每一个像素找到水面DEM上与之对应的位置上的像素，形成一个正四棱柱，正四棱柱的高即为水面DEM像素的值与地形DEM对应像素的值之差，依此过程循环，计算各正四棱柱体积之和即为堰塞湖总库容。插值处理时，需要将地形DEM的分辨率和像素大小传递给插值的水面DEM，否则，计算库容时，水面DEM与地形DEM上的像素无法对应形成正四棱柱。山区堰塞湖的上下游水位相差很大，必须考虑动库容才能准确计算其总库容，这种基于像素的库容监测方法更能满足实际需求。本发明研究了其数学模型，在理论情况下，对于连续曲面：

[0056] 式中，Hw是水面高程，HT是地形高程。

[0057] 本发明求取基于像素的离散曲面，结果逼近理想的连续曲面情况：

[0058]

[0059] 式中，HWij是水面DEM在位置(i，j)处像素的高程值，HTij是地形DEM在位置(i，j)处像素的高程值，Scell是单个像素的面积，m、n分别为水面DEM长和宽像素个数。参见附图3所示。具体实施时可以直接采用该公式，由计算机自动运行计算。

[0060] 具体实施时，可以采用计算机手段自动化完成上述步骤：在程序设计时创建一个空的点集合，利用函数将水体边界多边形叠加在地形DEM之上，然后每隔若干个像素对水体边界多边形与地形DEM像素求交，水体边界多边形覆盖的像素点即为交点，并获取该像素的高程赋值给交点，并将带高程的交点存入点集合中，通过得到的三维点集合生成要素类，利用该要素类进行Kriging插值，并设置属性，如插值的搜索半径，输出栅格图像像素尺寸，以及空间参考；然后将输出的栅格图像通过剪裁函数，利用水体边界剪裁有效区域；输入水面DEM和地形DEM到库容计算函数，以地形DEM为基准，从第一行第一列的像素开始，搜索水面DEM上对应为之的像素，若为空，则移动地形DEM像素，直到可以在水面DEM上找到对应的像素，根据两像素间构成的正四棱柱计算其体积，依此循环所有像素，并累加所计算的所有体积，直到地形DEM上所有像素遍历完，才停止计算，最后所得的体积和即为库容量。

[0061] 为了进一步提高监测效率，本发明还提供了对正四棱柱体积求和时进行分段简化处理，基于水体分段的优化监测方法，包括以下步骤：利用地形DEM进行流域分割，提取待测水体的水系线；参照待测水体的水面比降，对待测水体分段；根据流域分割所得待测水体宽度，以水系线为轴线，在每段待测水体两侧创建缓冲区；对每段待测水体的缓冲区，提取落入其内的三维水边界离散点所对应的DEM栅格单元高程值，并求平均值作为该段待测水体的平均水位，根据平均水位确定该段待测水体的正四棱柱体高度。

[0062] 本发明实施例首先对1：5万DEM填洼和削峰来消除地形DEM生产中的系统误差，然后给定提取阈值，进行流域分割。因为地形DEM实时性不强，获取的是待测水体大概范围。然后根据流域提取水系线，也被称为河网水系线。从地形DEM提取流域界线，首先要对地形中的洼地和尖峰进行处理，洼地和尖峰的存在使得在计算水流方向时会出现水流逆流的情况，给水流线的跟踪和流域界线的确定带来困难。然后确定水流在地形DEM每个栅格单元格内的流动方向，根据相邻的八个单元格来确定。在通过流域汇流分析来确定河流网络，进而在河流网络的基础上确定流域边界即分水线。可以借用ESRI公司开发的ArcHydro工具实现河网水系提取。

[0063] 在程序计算前，考虑水面比降，在水面比降较大的地方将水体分段，于是整个水体被分割成若干段；根据河流湖泊宽度，通过对每段水网线创建缓冲区，以保证生成的缓冲区包括了所有水体范围，将计算限制在缓冲区之内以提高计算效率；对落入每个缓冲区内的三维水体边界离散点高程值求平均值，作为该段河流或湖面的平均水位，然后按像素构建正四棱柱的方法动态监测库容。在对于同一缓冲区内的水面，不需要进行水位高程值的搜索，而直接利用平均值来计算，因此可以适当减少处理时间，提高效率。缓冲区就是在线两边一定距离的范围，宽度应略大于河流湖泊宽度的一半。求平均值的过程在像素级库容循环过程之前实施。求取库容的数学模型可以表示为：

[0064]

[0065] 式中，Hwk是水面DEM在第k段处的平均水位，HTij是地形DEM在位置(i，j)处像素的高程值，Scell是单个像素的面积，m、n分别为水面DEM长和宽像素个数，s为待测水体的分段数。

[0066] 为了说明本发明效果，选择ESRI公司集成在桌面系统ArcGIS中的监测方法，与本发明所提出的基于像素的监测方式与基于水体分段的监测方式的效率进行实验比较。实验数据包括1：5万25米分辨率的DEM，大小为207×568个像素，以及2008年5月16日地震后获取的SPOT-5影像。实验数据覆盖了唐家山堰塞湖回水尾部到堰塞坝部分。实施例监测的库容结果及比较如表1所示。

[0067] 表1监测所得库容量比较

[0068]

[0069]

[0070] ArcGIS中集成的库容监测方法是将静库容作为监测的总库容，仅考虑了水位，并将DEM上低于该水位所有空间都包括在库容中。可以从实验结果看出，ArcGIS集成的监测方法与基于像素的三维水面的库容监测方法得到的库容量减少了228.35万立方米，该结果实际上也是堰塞湖动库容与静库容之差，导致该结果的原因主要是由于山区堰塞湖上下游水位相差约60米所致。另外，由于地形DEM生产的误差，在计算范围内非河道处存在的某些高程极低或极高的情况，在该方法中没有排除其参与计算，使得计算结果与真值偏差较大，因此，该方法计算速度虽然比较快，但是，只是适合于水面水平的静态湖泊水库库容的监测。

[0071] 基于像素的库容监测方法充分考虑山区堰塞湖的地形特征，而且在计算中，参与计算的像素限制在河道附近的水面DEM上，在空白区域(该区域像素值为NoData)只做空循环不参与计算，因此，最大限度的排除了DEM生产误差的影响，可以认为该方法得到的结果接近真值。基于河流分段的优化监测方法所得库容与之仅相差1.35％，效率大幅度提高。此外，这种方式不依赖于河槽形态和河流边界形状，在计算过程中不需要各种水文地质资料，节省了大量的人力物力。相对于水力学计算，省去了过多的简化假设条件对精度的影响。

[0072] 表2运行效率比较

[0073]

[0074] 表2比较了三种监测方法的运行效率，ArcGIS集成的监测方法虽然具有很可观的效率，但由上面的分析可知，该方法不适合用于崎岖山区上下游水位比降过大的堰塞湖库容监测。基于像素的监测方法，从监测库容的原理上具有合理性，但是，由于需要大量搜索对应像素而导致其效率较低。基于水体分段的优化监测方法从量级上提高了计算效率，而监测的库容精度没有过多损失，从而使该方法更适合于日常监测。

[0075] 可见，基于水体分段的优化监测方法充分考虑了山区地形中总库容动态监测的理论合理性，同时也考虑了时间消耗，具有普遍适用性。因此，具体实施时采用该方法对河流水库库容进行日常业务化监测，效果最佳。参见附图2，本发明提供了实施例的库容业务化监测流程以供实施参考：

[0076] 1、利用经过精正射纠正的高分辨率SPOT—5影像进行水体提取，得到待测水体边界，预先从地形DEM提取待测水体处的河床DEM；

[0077] 2、将提取的水体边界与河床DEM叠加求交，求解水体边界线与河床DEM上相交的像素点，获得三维水边界离散点；

[0078] 3、通过对离散点进行插值，并剪裁出实施区域，获得湖泊河流的三维水面，即水面DEM，插值时通过参数设置，使水面DEM的分辨率、像素尺寸和空间参考等属性与地形DEM保持一致；

[0079] 4、水面DEM与河床DEM之间所包容的空间即为包括动库容在内的总库容。以河床DEM为基准，针对其上的每一个像素，搜索水面DEM上与之对应的像素。或者，在求解库容之前，对地形DEM进行河网提取，提取待测水体的河网水系线，在数据预处理时，针对水位比降大的地方进行分段并创建缓冲区。对于日常业务化运行，在水位没有剧烈增长时，没有必要每次计算都做一次河流分段，在一次分段后可多次使用，因此在下次计算时可以跳过。落入每个缓冲区内的三维水体边界离散点高程值求平均值，作为该段河流或湖面的平均水位，视为该段水体具有水平水面。缓冲区内水平水面与河床DEM间通过像素构造正四棱柱。

[0080] 5、两个像素之间的空间组成一个正四棱柱，在此将DEM栅格图像上的像素看成正方形。正四棱柱的高即为水面DEM像素的高程值与河床DEM对应像素的高程值之差，或者水平水面的高度值与河床DEM对应像素的高程值之差。然后利用体积计算公式求算该正四棱柱的体积。针对河床DEM每个像素循环，计算每个正四棱柱的体积并累积，即可得到待测水体区域库容量。