Method for estimating regional flood volume based on fusion of multi-source satellite remote sensing images and DEM

Abstract

This invention proposes a regional flood volume estimation method based on the fusion of multi-source satellite remote sensing images and DEM (Digital Elevation Model). Utilizing flood extent extraction results based on multi-source SAR and optical satellite data, the water surface elevation is adaptively determined within a local neighborhood, avoiding systematic biases introduced by the global assumption. A high-resolution water depth distribution map is constructed by traversing the regional flood extent through a sliding window, and then the total flood volume is obtained by summing these maps. This method requires no hydrological station network calibration and only two types of basic data: multi-source satellite remote sensing images and DEM, enabling rapid estimation of flood discharge from remote sensing observation data. This invention overcomes the systematic biases of the global constant water surface elevation assumption in complex terrain, achieving rapid and high-precision estimation of flood volume in areas without a hydrological station network through local adaptive water depth inversion.

Description

Technical Field

[0001] This invention relates to the field of remote sensing and hydrology, and in particular to a method for estimating regional flood volume based on the fusion of multi-source satellite remote sensing images and DEM. Background Technology

[0002] Flood disasters have a severe impact on the natural environment and power grid facilities. Traditional flood monitoring methods mainly rely on local hydrological parameters from surface hydrological stations, which cannot effectively address the following issues: 1) There are a large number of data-free areas globally, and the spatial distribution of hydrological stations is sparse, with insufficient station density to support detailed regional flood mapping; 2) Extreme flood events can easily damage hydrological station infrastructure, resulting in data loss during critical periods; 3) Hydrological station monitoring can only obtain cross-sectional flow information and cannot directly reflect the spatial distribution of floods and total water storage at the regional scale.

[0003] The rapid development of multi-source satellite remote sensing technology has provided a technological path to overcome the aforementioned bottlenecks. Among them, various synthetic aperture radars (SAR) and optical remote sensing satellites have become core data sources for flood disaster monitoring. However, existing multi-source satellite remote sensing data is mainly used for two-dimensional flood inundation range identification, lacking the ability to quantitatively invert flood volume.

[0004] The core technical challenge of the aforementioned problem lies in comprehensively utilizing two-dimensional flood inundation range products and three-dimensional DEM data to determine the spatial distribution of water surface elevation and water depth, and then estimating the water depth at each pixel. The global constant water surface elevation method assumes a uniform water surface elevation throughout the inundation area, then the water depth at any point is the difference between this constant elevation and the DEM elevation at that point. However, in areas with significant elevation differences, the global constant water surface elevation method tends to underestimate the water surface elevation at higher elevations upstream of the inundation area and overestimate the water surface elevation at lower elevations downstream. The river cross-section hydraulic extrapolation method calculates the water surface line based on hydraulic equations combined with river cross-sectional morphology. However, this method has high data requirements and complex calculations, requiring detailed river cross-section data, roughness parameters, and boundary conditions, making it difficult to apply to rapid estimation over large-scale watersheds. Summary of the Invention

[0005] The present invention aims to at least partially solve one of the technical problems in the related art.

[0006] Therefore, the first objective of this invention is to propose a method for estimating regional flood volume based on the fusion of multi-source satellite remote sensing images and DEM.

[0007] To achieve the above objectives, the first aspect of the present invention proposes a method for estimating regional flood volume based on the fusion of multi-source satellite remote sensing images and DEM, comprising: acquiring multi-source satellite remote sensing flood range sequence data and digital elevation model data of the area to be estimated, and performing spatial registration and binarization processing on the cleaned digital elevation model data and the flood range sequence data to obtain registered flood inundation area mask data and cleaned terrain elevation data;

[0008] Based on the registered flood inundation area mask data, a sliding window is constructed within the flood inundation area, and the topographic elevation difference index of the flood area within the sliding window is calculated. The sliding window size at the current pixel position is adaptively determined according to the topographic elevation difference index. Using a sliding window of a determined size, each submerged pixel in the flood inundation area mask data is traversed. The maximum topographic elevation value corresponding to all submerged pixels in the window is extracted as the local water surface elevation. The difference between the local water surface elevation and the topographic elevation data at the submerged pixel is calculated to obtain the water depth data of each submerged pixel to construct a water depth distribution map. The area per unit pixel is calculated based on the geographic coordinate parameters of the registered flood inundation area mask data. The water depth data of each inundated pixel in the water depth distribution map is multiplied by the area per unit pixel to obtain the volume of a single pixel water body. All single pixel water body volumes are then summed to obtain the total volume of the regional flood.

[0009] In one embodiment of the present invention, the step of spatially registering and binarizing the cleaned digital elevation model data with the flood range sequence data to obtain registered flood inundation area mask data and cleaned topographic elevation data includes: The flood range data of each frame in the multi-source satellite remote sensing flood range sequence data is converted into a binary matrix, wherein the pixels determined to be flooded are assigned a value of 1, and the pixels of dry land or other land types are assigned a value of 0, and a set of pixel coordinates of the flooded area is constructed based on the pixels assigned a value of 1. The values ​​without data markers and abnormal elevation values ​​less than or equal to -999 in the digital elevation model data are uniformly set to zero, while the remaining elevation values ​​remain unchanged to obtain the cleaned terrain elevation data. The geographic coordinates corresponding to each pixel are calculated based on the affine transformation parameters of the flood range sequence data. The geographic coordinates are then used to resample and crop the cleaned terrain elevation data so that it has the same matrix dimension as the binary matrix. Spatial registration is then completed and registered flood inundation area mask data is generated.

[0010] In one embodiment of the present invention, the step of converting each frame of flood extent data in the multi-source satellite remote sensing flood extent sequence data into a binary matrix, wherein pixels determined to be flooded are assigned a value of 1, pixels of dry land or other land types are assigned a value of 0, and constructing a set of pixel coordinates of the flooded area based on the pixels assigned a value of 1, includes: Obtain the total number of multi-source satellite remote sensing images of the current area and the number of rows and columns of a single frame image, and construct an initial flood range sequence data matrix. The pixel value is defined as 255 if the pixel is determined to be flooded, and 0 if the pixel is dry land or other land types. Iterate through each pixel in the initial flood range sequence data matrix. When the pixel value is 255, set the value of its corresponding position in the binary matrix to 1. When the pixel value is zero, set the value of its corresponding position in the binary matrix to 0, thus generating a binary matrix in a mathematical sense. Extract the row and column indices of all pixels with a value of 1 in the binary matrix to form a set of pixel coordinates of the flood-inundated area, and count the total number of elements in the set as the total number of pixels in the flood-inundated area.

[0011] In one embodiment of the present invention, the step of calculating the geographic coordinates corresponding to each pixel based on the affine transformation parameters of the flood range sequence data, resampling and cropping the cleaned terrain elevation data using the geographic coordinates to make it have the same matrix dimension as the binary matrix, completing spatial registration, and generating registered flood inundation area mask data includes: Read the affine transformation parameter six-tuple of the flood range sequence data. The six-tuple includes the east coordinate of the top left corner of the top left pixel, the pixel width, the row rotation parameter, the north coordinate of the top left corner of the top left pixel, the pixel height, and the column rotation parameter. The mapping relationship between pixel row and column indices and geographic coordinates is constructed using the affine transformation parameter six-tuple, and the east and north coordinates corresponding to the center point of each pixel in the binary matrix are calculated. Based on the calculated geographic coordinates, the cleaned terrain elevation data is subjected to geographic transformation and resampling processing. The overlapping area of ​​the two geographic ranges is extracted, and the resampled terrain elevation data and the binary matrix are aligned to the same number of rows and columns to obtain the registered flood inundation area mask data and the cleaned terrain elevation data with consistent dimensions.

[0012] In one embodiment of the present invention, the step of adaptively determining the sliding window size of the current pixel position based on the terrain elevation difference index includes: Centered on the currently traversed submerged pixel, extract all pixels belonging to the flood-inundated area within the preset initial sliding window range, and obtain the set of elevation values ​​corresponding to these pixels in the cleaned terrain elevation data. Calculate the difference between the maximum and minimum elevation values ​​in the set of elevation values, and define the difference as the terrain elevation difference index of the neighborhood of the current pixel location, which is used to characterize the degree of undulation of the local terrain. Based on the preset maximum and minimum window size pixel values, and combined with the statistical maximum and minimum values ​​of the current area elevation difference, the initial window size is calculated using an inverse proportional function relationship. The smaller the terrain elevation difference index, the larger the calculated window size; the larger the terrain elevation difference index, the smaller the calculated window size. The calculated initial window size is subjected to boundary constraints to ensure that it is between the minimum and maximum window size pixel values. The constrained values ​​are then rounded down and oddized to obtain the final sliding window size applied to the current pixel position.

[0013] In one embodiment of the present invention, the calculation of the initial window size using an inverse proportional function relationship, wherein a smaller terrain elevation difference index results in a larger calculated window size, and a larger terrain elevation difference index results in a smaller calculated window size, includes: The maximum window size pixel value is initialized to 200, the minimum window size pixel value is initialized to 10, and a very small positive number is set as a correction term to prevent division by zero errors; A formula for calculating the window size is constructed, in which the maximum window size pixel value is subtracted from a correction amount as the numerator. The correction amount is equal to the difference between the maximum window size pixel value and the minimum window size pixel value multiplied by the terrain elevation difference index of the neighborhood of the current pixel position. The quotient is obtained by dividing the numerator by a denominator, and the denominator is equal to the difference between the statistical maximum and statistical minimum values ​​of the current area elevation difference plus the correction term to prevent division by zero error. The absolute value of the quotient is taken and rounded down to obtain the integer value of the initial window size. This ensures that a large window is used in flat terrain areas to suppress noise, and a small window is used in areas with drastic terrain undulations to avoid crossing terrain watersheds.

[0014] In one embodiment of the present invention, the step of performing boundary constraint processing on the calculated initial window size to ensure that it is between the minimum window size pixel value and the maximum window size pixel value, and performing floor-down and odd-number processing on the constrained value to obtain the final sliding window size applied to the current pixel position includes: First, a maximum value constraint operation is performed on the initial window size, which is compared with the minimum window size pixel value and the larger one is retained. Then, a minimum value constraint operation is performed on the comparison result, which is compared with the maximum window size pixel value and the smaller one is retained. This gives the window size after boundary constraints. Divide the window size after boundary constraints by two and perform a floor operation, then multiply the result by two and add one to obtain the final odd-numbered sliding window size; Based on the final odd-numbered sliding window size, determine the pixel index range covered by the square sliding window centered on the current submerged pixel, ensuring that the window does not exceed the effective image domain boundary, and extract all pixel subsets belonging to the flood-inundated area within the window for subsequent local water surface elevation calculation.

[0015] In one embodiment of the present invention, the step of extracting the maximum topographic elevation value corresponding to all submerged pixels within the extraction window as the local water surface elevation, and calculating the difference between the local water surface elevation and the topographic elevation data at the submerged pixel to obtain the water depth data of each submerged pixel to construct a water depth distribution map includes: For the currently traversed submerged pixels, within the coverage area of ​​the sliding window after the size is determined, filter out all pixels marked as submerged in the registered flood submerged area mask data; Extract the elevation values ​​of the selected pixels in the cleaned terrain elevation data, and determine the maximum value among these elevation values ​​as the local water surface elevation of the current submerged pixel location; Calculate the difference between the local water surface elevation and the topographic elevation value of the current submerged pixel in the cleaned topographic elevation data, and perform non-negative processing on the difference. If the difference is greater than zero, retain the original value; if the difference is less than or equal to zero, force it to be set to zero to obtain the water depth data of the current submerged pixel. Iterate through all submerged pixels in the flood-inundated area mask data and repeat the above steps for determining local water surface elevation and calculating water depth data. Fill the calculated water depth data into the corresponding pixel positions, and set the water depth data of non-submerged area pixels to zero to generate a complete water depth distribution map.

[0016] In one embodiment of the present invention, the step of calculating the unit pixel area based on the geographic coordinate parameters of the registered flood inundation area mask data, multiplying the water depth data of each inundated pixel in the water depth distribution map by the unit pixel area to obtain the single pixel water volume, and accumulating all single pixel water volumes to obtain the total flood volume of the region, includes: Pixel width and pixel height parameters are extracted from the affine transformation parameters of the registered flood inundation area mask data. The absolute value of the pixel width parameter and the absolute value of the pixel height parameter are calculated, and the product is used as the unit pixel area. Traverse each pixel in the water depth distribution map. When a pixel is located in the flood inundation area, read the water depth data of the pixel and multiply it by the unit pixel area to obtain the single pixel water volume corresponding to the pixel. When a pixel is located in the non-flooded area, set its single pixel water volume to zero. The individual pixel water volume of all flood-inundated pixels is summed to obtain the total flood volume of the region in cubic meters, which is used as the flood storage capacity estimation result for the region to be estimated.

[0017] In one embodiment of the present invention, the method further includes a step of evaluating the validity of the flood volume estimation results: The total area of ​​the flood zone at the current moment is calculated by multiplying the total number of flood zone pixels in the registered flood inundation area mask data with the area of ​​the unit pixel. The total volume of regional floods and the corresponding total area of ​​flood range were estimated at multiple different times. The two sets of data sequences were scaled to the same scale for trend analysis. Calculate the Pearson correlation coefficient between the estimated total flood volume sequence and the calculated total flood area sequence. If the correlation coefficient is higher than a preset threshold, the regional flood volume inversion is considered successful; otherwise, the inversion result is considered to deviate from the actual value.

[0018] The method of this invention adaptively determines water surface elevation through a local sliding window, effectively overcoming the systematic bias caused by the global constant assumption in areas with undulating terrain, and significantly improving the accuracy of flood volume estimation. This method does not rely on hydrological station network data, enabling rapid and quantitative inversion of regional flood storage in areas without available data.

[0019] Additional aspects and advantages of the invention will be set forth in part in the description which follows, and in part will be obvious from the description, or may be learned by practice of the invention. Attached Figure Description

[0020] The above and / or additional aspects and advantages of the present invention will become apparent and readily understood from the following description of the embodiments taken in conjunction with the accompanying drawings, wherein: Figure 1 A flowchart illustrating a regional flood volume estimation method based on the fusion of multi-source satellite remote sensing images and DEMs, provided in an embodiment of the present invention; Figure 2 An example diagram illustrating the correlation analysis between the flood volume estimation results and the actual flood area provided in this embodiment of the invention. Detailed Implementation

[0021] It should be noted that, unless otherwise specified, the embodiments and features described in the present invention can be combined with each other. The present invention will now be described in detail with reference to the accompanying drawings and embodiments.

[0022] To enable those skilled in the art to better understand the present invention, the technical solutions of the present invention will be clearly and completely described below with reference to the accompanying drawings of the embodiments of the present invention. Obviously, the described embodiments are only some embodiments of the present invention, and not all embodiments. Based on the embodiments of the present invention, all other embodiments obtained by those skilled in the art without creative effort should fall within the scope of protection of the present invention.

[0023] The following describes, with reference to the accompanying drawings, a method for estimating regional flood volume based on the fusion of multi-source satellite remote sensing images and DEM, according to an embodiment of the present invention.

[0024] This invention utilizes flood extent extraction results based on multi-source SAR and optical satellite data to adaptively determine water surface elevation within a local neighborhood, avoiding systematic biases caused by global assumptions. By traversing the regional flood extent through a sliding window, a high-resolution water depth distribution map is constructed, and then the total flood volume is obtained by summing the data. This method requires no hydrological station network calibration and only two types of basic data: multi-source satellite remote sensing images and DEM (Digital Emission Model), enabling rapid estimation of flood discharge from remote sensing observation data.

[0025] Figure 1 This is a flowchart of a regional flood volume estimation method based on the fusion of multi-source satellite remote sensing images and DEM, according to an embodiment of the present invention. Figure 1 As shown, the method includes the following steps: S1. Acquire multi-source satellite remote sensing flood range sequence data and digital elevation model data of the area to be estimated, and perform spatial registration and binarization processing on the cleaned digital elevation model data and the flood range sequence data to obtain the registered flood inundation area mask data and the cleaned terrain elevation data. S2, Based on the registered flood inundation area mask data, a sliding window is constructed within the flood inundation area, and the topographic elevation difference index of the flood area within the sliding window is calculated. The sliding window size of the current pixel position is adaptively determined according to the topographic elevation difference index. S3. Using a sliding window of a determined size, traverse each submerged pixel in the flood inundation area mask data, extract the maximum topographic elevation corresponding to all submerged pixels in the window as the local water surface elevation, and calculate the difference between the local water surface elevation and the topographic elevation data at the submerged pixel to obtain the water depth data of each submerged pixel to construct a water depth distribution map. S4. Calculate the area per unit pixel based on the geographic coordinate parameters of the registered flood inundation area mask data. Multiply the water depth data of each inundated pixel in the water depth distribution map by the area per unit pixel to obtain the volume of a single pixel water body. Then, sum up all the volumes of single pixels water bodies to obtain the total volume of the regional flood.

[0026] Specifically, the regional flood volume estimation method based on the fusion of multi-source satellite remote sensing images and DEMs of the present invention includes data preprocessing, local sliding window water depth inversion, and flood volume summation calculation. The mathematical principles of each step are explained in detail below.

[0027] In one embodiment of the present invention, the I data preprocessing step is as follows: (1) Input data 1) Flood extent sequence data extracted from multi-source satellite remote sensing images: The flood extent sequence data is represented as a matrix. .in, , This indicates the total number of multi-source satellite remote sensing images of the current area. This represents the current image row number. This represents the current column number of the image. Given a binary matrix, the pixel values ​​are defined as follows:

[0028] in, , For pixel row and column indices.

[0029] 2) Digital Elevation Model (DEM) The digital elevation model (DEM) is represented as a matrix. The unit is meters (m). It contains terrain and elevation information, but there may be areas with no data or outliers, so preprocessing is required.

[0030] (2) Binarization of flood extent Will Convert to a binary matrix in the mathematical sense. :

[0031] Define the set of pixel coordinates of the flood-inundated area as:

[0032] The total number of pixels in the flooded area is denoted as .

[0033] (3) DEM data cleaning The DEM data contains missing data values ​​and abnormal elevation values, requiring data cleaning. Let the original missing data marker value of the DEM be... The cleaned DEM is denoted as :

[0034] (4) Spatial registration because and The spatial resolution, map projection, and geographic extent may differ. This step achieves accurate registration between the two through geographic coordinate transformation.

[0035] set up The affine transformation parameters (GeoTransform) are six-tuples. ,in: The eastward (X) coordinate of the top-left pixel; : Pixel width (East-West resolution, usually a positive value); : Row rotation parameter (usually 0); The north (Y) coordinate of the top-left pixel; Pixel height (north-south resolution, usually a negative value); : Column rotation parameter (usually 0); then the pixel Corresponding geographic coordinates for:

[0036] Resampling and cropping based on geographic coordinates, and Having the same matrix dimensions ,in, .

[0037] (5) Pixel area calculation For those already registered and Data, area per pixel (Unit: square meters, m²) Calculated from affine transformation parameters:

[0038] in, This represents absolute value operations. The physical meaning of this formula is: the distance of a pixel from the ground in the east-west direction. Length of the ground in the north-south direction The product of. For example, for DEM data with standard UTM projection, if If the resolution is 30 meters, then .

[0039] In one embodiment of the present invention, the II local sliding window water depth inversion step is as follows: This step aims to address the systematic errors caused by the assumption of a globally constant water surface elevation.

[0040] (1) Method Principle Local horizontal plane assumption: Within any local neighborhood of the flood-inundated area, the water surface can be approximated as a horizontal plane, meaning that all inundated pixels within that neighborhood have the same water surface elevation. This water surface elevation is determined by the highest point of the terrain within the neighborhood, and the highest point is usually located at the edge of the flood-inundated area.

[0041] Physical basis: Floodwaters tend to form a horizontal surface under the influence of gravity. In local areas with relatively small topographic relief, the water surface elevation is limited by the highest point in the vicinity—this point is usually the boundary between the flooded area and the dry surface, and its DEM elevation is the water surface elevation of that local water area.

[0042] (2) Definition of sliding window Set the half width of the square sliding window to (Unit: pixels), then the total width of the window is For pixels centered in the flooded area. The window covers a range of pixel indices:

[0043] in, For the effective image domain, This represents the intersection operation, ensuring the window does not exceed the image boundaries. The submerged pixel subset within the window is:

[0044] The corresponding DEM elevation subset is:

[0045] (3) Determination of local water surface elevation For the current center pixel Its local water surface elevation Defined as the maximum DEM elevation of all submerged pixels within the neighborhood:

[0046] (4) Pixel water depth calculation The water depth of a pixel is calculated based on the difference between the local water surface elevation and the current pixel terrain elevation. :

[0047] in, The difference between the water surface elevation and the ground elevation is theoretically the water depth. The calculation ensures that the water depth is non-negative. When anomalies appear in the local terrain that are higher than the water surface elevation, the depth is forcibly set to 0 to avoid negative water depths.

[0048] (5) Construction of complete water depth map Traverse all submerged pixels Perform the sliding window calculation described above to construct a complete water depth distribution map. :

[0049] In the non-submerged zone, the water depth is forcibly set to 0.

[0050] (6) Window size parameters Selection strategy Window size This is a key adjustment parameter in this method, directly affecting the inversion accuracy. Traditional fixed windows are difficult to adapt to complex and changing terrain conditions. This invention proposes a method for adaptively determining the window size based on the terrain elevation difference of the regional DEM.

[0051] 1) Definition of terrain complexity index For the current center pixel Define its neighborhood topographic elevation index The range of DEM elevations within the flood zone during the sliding window:

[0052] This index characterizes the degree of undulation in local terrain. This indicates that the terrain in a localized area is flat and the water surface is approximately horizontal. A larger value indicates significant terrain undulations, requiring the window to be narrowed to ensure the assumption of a local horizontal plane holds true.

[0053] 2) Adaptive window size calculation The window size is set to have a negative correlation with the terrain elevation difference, using an inverse proportional function:

[0054] in, This represents the maximum window size in pixels, initialized to 200. This represents the minimum window size, initialized to 10. and These represent the maximum and minimum statistical values ​​of the current area's elevation difference, respectively. Represents extremely small positive numbers (e.g.) ), to prevent division by zero errors; Rounding down ensures the window size is an odd number of pixels. In other words, when there is a difference in elevation... near hour, tending to A large window is used to fully utilize the statistical sample and suppress DEM noise; when the terrain elevation difference... near hour, tending to Small windows are used to ensure that the local horizontal plane assumption holds and to avoid crossing the topographic watershed.

[0055] 3) Window size constraints and odd number handling Boundary constraints and odd-number processing are applied to the calculated window size:

[0056] In one embodiment of the present invention, the steps for calculating the summation of flood volume III are as follows: (1) Single pixel water volume For pixels in flooded areas The corresponding water column volume The product of water depth and pixel area:

[0057] (2) Total volume of water The total flood volume is obtained by summing the water depth values ​​of discrete pixels and then multiplying by the area per unit pixel, expressed in cubic meters (m³). The formula for summing the water column volumes of all flood zone pixels is: Total Flood Volume. It can be defined as:

[0058] In one embodiment of the present invention, the evaluation steps for the IV flood volume estimation results are as follows: The flood volume is estimated by following steps I to III above, and then compared with the actual flood area to evaluate the effectiveness of the flood volume estimation results.

[0059] (1) Calculation of the area of ​​the flood zone Total area of ​​flood for:

[0060] (2) Correlation analysis between flood volume estimation results and actual flood area The estimated flood volume and calculated flood area within the current flood zone are scaled to the same scale for trend analysis. A higher correlation coefficient indicates a more successful flood volume inversion, while a lower correlation coefficient indicates that the inverted flood volume deviates from the actual value. An example of the correlation analysis between the estimated flood volume and the actual flood area is shown below. Figure 2 As shown. Among them, Figure 2 In the figure, (a) represents the trend of the time series change. Figure 2 (b) in the figure represents the volume-area correlation.

[0061] This invention presents a regional flood volume estimation method based on the fusion of multi-source satellite remote sensing images and DEMs. Utilizing flood extent extraction results from multi-source SAR and optical satellite data, it adaptively determines water surface elevation within a local neighborhood, avoiding systematic biases introduced by global assumptions. By traversing the regional flood extent through a sliding window, a high-resolution water depth distribution map is constructed, and then the total flood volume is obtained through summation. This method requires no hydrological station network calibration and only two types of basic data—multi-source satellite remote sensing images and DEMs—to achieve rapid estimation of flood discharge from remote sensing observation data.

[0062] In the description of this specification, the references to terms such as "one embodiment," "some embodiments," "example," "specific example," or "some examples," etc., refer to specific features, structures, materials, or characteristics described in connection with that embodiment or example, which are included in at least one embodiment or example of the present invention. In this specification, the illustrative expressions of the above terms do not necessarily refer to the same embodiment or example. Furthermore, the specific features, structures, materials, or characteristics described may be combined in any suitable manner in one or more embodiments or examples. Moreover, without contradiction, those skilled in the art can combine and integrate the different embodiments or examples described in this specification, as well as the features of different embodiments or examples.

[0063] Furthermore, the terms "first" and "second" are used for descriptive purposes only and should not be construed as indicating or implying relative importance or implicitly specifying the number of technical features indicated. Thus, a feature defined as "first" or "second" may explicitly or implicitly include at least one of that feature. In the description of this invention, "a plurality of" means at least two, such as two, three, etc., unless otherwise explicitly specified.

Claims

1. A method for estimating regional flood volume based on the fusion of multi-source satellite remote sensing images and DEM, characterized in that, include: Acquire multi-source satellite remote sensing flood extent sequence data and digital elevation model data of the area to be estimated, and perform spatial registration and binarization processing on the cleaned digital elevation model data and the flood extent sequence data to obtain the registered flood inundation area mask data and the cleaned terrain elevation data. Based on the registered flood inundation area mask data, a sliding window is constructed within the flood inundation area, and the topographic elevation difference index of the flood area within the sliding window is calculated. The sliding window size at the current pixel position is adaptively determined according to the topographic elevation difference index. Using a sliding window of a determined size, each submerged pixel in the flood inundation area mask data is traversed. The maximum topographic elevation value corresponding to all submerged pixels in the window is extracted as the local water surface elevation. The difference between the local water surface elevation and the topographic elevation data at the submerged pixel is calculated to obtain the water depth data of each submerged pixel to construct a water depth distribution map. The area per unit pixel is calculated based on the geographic coordinate parameters of the registered flood inundation area mask data. The water depth data of each inundated pixel in the water depth distribution map is multiplied by the area per unit pixel to obtain the volume of a single pixel water body. All single pixel water body volumes are then summed to obtain the total volume of the regional flood.

2. The method as described in claim 1, characterized in that, The step of spatially registering and binarizing the cleaned digital elevation model data with the flood range sequence data to obtain registered flood inundation area mask data and cleaned topographic elevation data includes: The flood range data of each frame in the multi-source satellite remote sensing flood range sequence data is converted into a binary matrix, wherein the pixels determined to be flooded are assigned a value of 1, and the pixels of dry land or other land types are assigned a value of 0, and a set of pixel coordinates of the flooded area is constructed based on the pixels assigned a value of 1. The values ​​without data markers and abnormal elevation values ​​less than or equal to -999 in the digital elevation model data are uniformly set to zero, while the remaining elevation values ​​remain unchanged to obtain the cleaned terrain elevation data. The geographic coordinates corresponding to each pixel are calculated based on the affine transformation parameters of the flood range sequence data. The geographic coordinates are then used to resample and crop the cleaned terrain elevation data so that it has the same matrix dimension as the binary matrix. Spatial registration is then completed and registered flood inundation area mask data is generated.

3. The method as described in claim 2, characterized in that, The step of converting each frame of flood extent data in the multi-source satellite remote sensing flood extent sequence data into a binary matrix, wherein pixels identified as flooded are assigned a value of 1, and pixels of dry land or other land types are assigned a value of 0, and constructing a set of pixel coordinates of the flooded area based on the pixels assigned a value of 1, includes: Obtain the total number of multi-source satellite remote sensing images of the current area and the number of rows and columns of a single frame image, and construct an initial flood range sequence data matrix. The pixel value is defined as 255 if the pixel is determined to be flooded, and 0 if the pixel is dry land or other land types. Iterate through each pixel in the initial flood range sequence data matrix. When the pixel value is 255, set the value of its corresponding position in the binary matrix to 1. When the pixel value is zero, set the value of its corresponding position in the binary matrix to 0, thus generating a binary matrix in a mathematical sense. Extract the row and column indices of all pixels with a value of 1 in the binary matrix to form a set of pixel coordinates of the flood-inundated area, and count the total number of elements in the set as the total number of pixels in the flood-inundated area.

4. The method as described in claim 2, characterized in that, The process involves calculating the geographic coordinates of each pixel based on the affine transformation parameters of the flood range sequence data, resampling and cropping the cleaned terrain elevation data using these geographic coordinates to ensure it has the same matrix dimension as the binary matrix, completing spatial registration, and generating registered flood-inundated area mask data. This includes: Read the affine transformation parameter six-tuple of the flood range sequence data. The six-tuple includes the east coordinate of the top left corner of the top left pixel, the pixel width, the row rotation parameter, the north coordinate of the top left corner of the top left pixel, the pixel height, and the column rotation parameter. The mapping relationship between pixel row and column indices and geographic coordinates is constructed using the affine transformation parameter six-tuple, and the east and north coordinates corresponding to the center point of each pixel in the binary matrix are calculated. Based on the calculated geographic coordinates, the cleaned terrain elevation data is subjected to geographic transformation and resampling processing. The overlapping area of ​​the two geographic ranges is extracted, and the resampled terrain elevation data and the binary matrix are aligned to the same number of rows and columns to obtain the registered flood inundation area mask data and the cleaned terrain elevation data with consistent dimensions.

5. The method as described in claim 1, characterized in that, The step of adaptively determining the sliding window size of the current pixel position based on the terrain elevation difference index includes: Centered on the currently traversed submerged pixel, extract all pixels belonging to the flood-inundated area within the preset initial sliding window range, and obtain the set of elevation values ​​corresponding to these pixels in the cleaned terrain elevation data. Calculate the difference between the maximum and minimum elevation values ​​in the set of elevation values, and define the difference as the terrain elevation difference index of the neighborhood of the current pixel location, which is used to characterize the degree of undulation of the local terrain. Based on the preset maximum and minimum window size pixel values, and combined with the statistical maximum and minimum values ​​of the current area elevation difference, the initial window size is calculated using an inverse proportional function relationship. The smaller the terrain elevation difference index, the larger the calculated window size; the larger the terrain elevation difference index, the smaller the calculated window size. The calculated initial window size is subjected to boundary constraints to ensure that it is between the minimum and maximum window size pixel values. The constrained values ​​are then rounded down and oddized to obtain the final sliding window size applied to the current pixel position.

6. The method as described in claim 5, characterized in that, The calculation of the initial window size using an inverse proportional function relationship, wherein a smaller terrain elevation difference index results in a larger calculated window size, and a larger terrain elevation difference index results in a smaller calculated window size, includes: The maximum window size pixel value is initialized to 200, the minimum window size pixel value is initialized to 10, and a very small positive number is set as a correction term to prevent division by zero errors; A formula for calculating the window size is constructed, in which the maximum window size pixel value is subtracted from a correction amount as the numerator. The correction amount is equal to the difference between the maximum window size pixel value and the minimum window size pixel value multiplied by the terrain elevation difference index of the neighborhood of the current pixel position. The quotient is obtained by dividing the numerator by a denominator, and the denominator is equal to the difference between the statistical maximum and statistical minimum values ​​of the current area elevation difference plus the correction term to prevent division by zero error. The absolute value of the quotient is taken and rounded down to obtain the integer value of the initial window size. This ensures that a large window is used in flat terrain areas to suppress noise, and a small window is used in areas with drastic terrain undulations to avoid crossing terrain watersheds.

7. The method as described in claim 5, characterized in that, The process of performing boundary constraint processing on the calculated initial window size to ensure that it lies between the minimum and maximum window size pixel values, and then performing floor-rounding and odd-numbering on the constrained value to obtain the final sliding window size applied to the current pixel position, includes: First, a maximum value constraint operation is performed on the initial window size, which is compared with the minimum window size pixel value and the larger one is retained. Then, a minimum value constraint operation is performed on the comparison result, which is compared with the maximum window size pixel value and the smaller one is retained. This gives the window size after boundary constraints. Divide the window size after boundary constraints by two and perform a floor operation, then multiply the result by two and add one to obtain the final odd-numbered sliding window size; Based on the final odd-numbered sliding window size, determine the pixel index range covered by the square sliding window centered on the current submerged pixel, ensuring that the window does not exceed the effective image domain boundary, and extract all pixel subsets belonging to the flood-inundated area within the window for subsequent local water surface elevation calculation.

8. The method as described in claim 1, characterized in that, The maximum topographic elevation corresponding to all submerged pixels within the extraction window is used as the local water surface elevation. The difference between the local water surface elevation and the topographic elevation data at the submerged pixel is calculated to obtain the water depth data of each submerged pixel to construct a water depth distribution map, including: For the currently traversed submerged pixels, within the coverage area of ​​the sliding window after the size is determined, filter out all pixels marked as submerged in the registered flood submerged area mask data; Extract the elevation values ​​of the selected pixels in the cleaned terrain elevation data, and determine the maximum value among these elevation values ​​as the local water surface elevation of the current submerged pixel location; Calculate the difference between the local water surface elevation and the topographic elevation value of the current submerged pixel in the cleaned topographic elevation data, and perform non-negative processing on the difference. If the difference is greater than zero, retain the original value; if the difference is less than or equal to zero, force it to be set to zero to obtain the water depth data of the current submerged pixel. Iterate through all submerged pixels in the flood-inundated area mask data and repeat the above steps for determining local water surface elevation and calculating water depth data. Fill the calculated water depth data into the corresponding pixel positions, and set the water depth data of non-submerged area pixels to zero to generate a complete water depth distribution map.

9. The method as described in claim 1, characterized in that, The method involves calculating the unit pixel area based on the geographic coordinate parameters of the registered flood inundation area mask data, multiplying the water depth data of each inundated pixel in the water depth distribution map by the unit pixel area to obtain the single pixel water volume, and then summing all the single pixel water volumes to obtain the total flood volume of the region, including: Pixel width and pixel height parameters are extracted from the affine transformation parameters of the registered flood inundation area mask data. The absolute value of the pixel width parameter and the absolute value of the pixel height parameter are calculated, and the product is used as the unit pixel area. Traverse each pixel in the water depth distribution map. When a pixel is located in the flood inundation area, read the water depth data of the pixel and multiply it by the unit pixel area to obtain the single pixel water volume corresponding to the pixel. When a pixel is located in the non-flooded area, set its single pixel water volume to zero. The individual pixel water volume of all flood-inundated pixels is summed to obtain the total flood volume of the region in cubic meters, which is used as the flood storage capacity estimation result for the region to be estimated.

10. The method as described in claim 9, characterized in that, The method also includes a step for evaluating the validity of flood volume estimation results: The total area of ​​the flood zone at the current moment is calculated by multiplying the total number of flood zone pixels in the registered flood inundation area mask data with the area of ​​the unit pixel. The total volume of regional floods and the corresponding total area of ​​flood range were estimated at multiple different times. The two sets of data sequences were scaled to the same scale for trend analysis. Calculate the Pearson correlation coefficient between the estimated total flood volume sequence and the calculated total flood area sequence. If the correlation coefficient is higher than a preset threshold, the regional flood volume inversion is considered successful; otherwise, the inversion result is considered to deviate from the actual value.

Images