New reservoir automatic identification method based on two-step remote sensing mapping

Abstract

The invention provides a new reservoir automatic identification method based on two-step remote sensing mapping, and the method comprises the steps: determining the spatial position and the maximum water area range of a potential new reservoir based on the multi-year water body frequency synthesized by medium and high resolution remote sensing images, and dividing and comparing the accumulated water body frequency difference in different time periods according to different years; for a potential newly-built reservoir object, a year-by-year average water body area time sequence under the maximum water area range constraint is extracted to carry out trend analysis and breakpoint detection, and reservoir building / water storage activities are positioned according to trend characteristics before and after breakpoints disclosed by a time sequence analysis result. Spatial position extraction, water area spatial range matching and reservoir water storage time identification of the newly-built reservoir / artificial lake are effectively realized. According to the method, medium-high spatial resolution images can be obtained freely, the application range of the method is expanded, and an important technical method support is provided for establishing efficient large-scale newly-built reservoir cataloguing and watershed reservoir water storage and drainage information parameterization.

Description

Technical Field

[0001] The invention relates to the field of remote sensing technology, and in particular to a new reservoir automatic identification method based on two-step remote sensing mapping. Background Art

[0002] Reservoirs are water bodies formed by artificial dams. They play a very important role in the development and utilization of water resources. They are of great significance to water storage and flood control, agricultural irrigation, hydropower generation, transportation and shipping, tourism and entertainment. With the continuous growth of economy and energy demand, the number of reservoirs is still increasing significantly, but the data describing the geographical characteristics and distribution of reservoirs is incomplete. For example, the basic information such as the location, water storage time and water storage capacity of most small and medium-sized reservoirs in the world is difficult to obtain accurately. In particular, newly built reservoirs have not been updated in time, which greatly affects the accuracy of hydrological assessment and the scientific nature of basin regulation. Therefore, designing a method for automatic monitoring of spatiotemporal information of newly built reservoirs is of great significance for scientific research and government decision-making.

[0003] Automatic identification of newly built reservoirs refers to overcoming obstacles such as geographical location, determining the spatial location, spatial range and water storage time of newly built reservoirs within a larger spatial range, and solving the problem of difficulty in obtaining spatiotemporal information such as reservoir construction and operation. It mainly serves the decision-making of water conservancy, environment and other departments from regional to global, as well as research in the fields of hydrology and environment. Accurately obtaining the water storage time of newly built reservoirs can be used to estimate the changes in the area and storage capacity of global reservoirs, update the spatiotemporal information of global water resource distribution, and is of great significance for evaluating the impact of humans on river systems and global water cycles.

[0004] Existing technologies mainly focus on reservoir construction engineering technology and the estimation of water area, water depth, storage capacity and other attributes after the completion of the reservoir. For the information update of newly built reservoirs, except for the relevant information publicly released by relevant departments and units, there is currently no method that can realize large-scale identification of newly built reservoirs; and the relevant departments or units often disclose information on reservoir construction activities relatively late, and there is currently no similar method that can achieve low-cost, high-precision, large-scale real-time monitoring of reservoir construction activities, which makes it difficult for research in the fields of hydrology and environment to obtain new reservoir information, and related research is difficult to carry out in a timely and smooth manner, forming unnecessary obstacles. Summary of the invention

[0005] The purpose of the present invention is to overcome the defects existing in the above-mentioned prior art, including the lack of spatiotemporal information of newly built reservoirs, the delay in the disclosure of information by relevant departments, the lack of similar methods for realizing low-cost, high-precision, large-scale real-time monitoring of reservoir construction activities, and the like, and to provide a method for automatic identification of newly built reservoirs based on two-step remote sensing mapping.

[0006] To achieve the above object, the technical solution adopted by the present invention is as follows: A method for automatic identification of newly built reservoirs based on two-step remote sensing mapping, comprising: Obtain multi-year remote sensing image data of the study area; Dividing the multi-year remote sensing image data into several time periods according to time sequence, extracting water body pixels from the remote sensing image data, and calculating the pixel water body cumulative frequency of each time period respectively; Compare the pixel water body accumulation frequency of adjacent time periods. If the difference between the pixel water body accumulation frequency of the latter time period and the previous time period exceeds the preset threshold, the pixel is a newly added reservoir pixel. Match the selected new reservoir pixels with the water body dataset to obtain the newly built reservoirs and construct the annual average water area time series of the newly built reservoirs; The water storage time of the newly built reservoir is determined based on the area change trend and area mutation point of the newly built reservoir.

[0007] As a preferred implementation, the method of extracting water body pixels from remote sensing image data is: calculating the normalized water index NDWI for the remote sensing image, and extracting water body pixels in combination with a threshold segmentation method.

[0008] As a preferred implementation, the pixel water body cumulative frequency is the annual water body frequency mean of the pixel in the corresponding time period.

[0009] As a preferred implementation, the time scale of the period is ≥10 years.

[0010] As a preferred implementation, the annual average water area time series of the newly built reservoir is constructed as follows: Based on the water body dataset, the maximum water body boundary range of the newly built reservoir is obtained, and the annual maximum water body area of ​​each reservoir is determined within the annual maximum water body boundary range; Calculate the annual water body frequency mean of each pixel within the boundary range corresponding to the annual maximum water body area to obtain the annual average water body frequency; Calculate the product of the annual average water body frequency and the annual maximum water body area to obtain the annual average water body area; Calculate the annual average water area year by year and obtain the annual average water area time series.

[0011] As a preferred implementation, the annual maximum water area is the sum of the permanent water area and the seasonal water area within the annual maximum water boundary.

[0012] As a preferred implementation, the method further includes processing missing values ​​and outliers on the annual average water area time series, including: Traverse the value of each observation in the annual average water body time series, and for all missing values ​​or outliers, take the average of the two observations before the missing value or outlier and the two observations after the missing value or outlier to fill the missing values ​​or outliers; If missing values ​​or outliers appear in the first or last observation, the mean of several observations before or after the missing value or outlier is used to replace the missing value or outlier.

[0013] As a preferred implementation, the MK mutation test method is used to analyze the reservoir area change trend and find potential mutation points in the annual average water area time series.

[0014] As a preferred implementation, the area changes of the newly built reservoirs before and after the mutation point are compared. If the water area after the mutation point shows an increasing trend, the year of the mutation point is the year of water storage for the newly built reservoir. Otherwise, the mutation point is not the year of water storage for the newly built reservoir. Since the breakpoints detected by the MK mutation test algorithm are not all caused by reservoir construction, a significant reduction in water area caused by other reasons will also cause the MK mutation test algorithm to detect breakpoints. The UF before and after the mutation point can be analyzed. k The size of the change, if the breakpoint after UF k If the value of is greater than zero and continues to rise, it means that the water area after the breakpoint continues to expand. The year of the mutation point can represent the water storage year of the newly built reservoir. Otherwise, it is not the water storage year of the newly built reservoir.

[0015] As a preferred implementation, the spatial resolution of the remote sensing image data is ≤30m.

[0016] The present invention is based on free access to medium and high resolution remote sensing images, such as Landsat series satellite TM, ETM+ and OLI data, and uses two-step mapping (the first step: determine the spatial location of potential new reservoirs and obtain their water area: the second step: calculate the water storage time of the new reservoir) to achieve automatic identification of new reservoirs. The method of the present invention accurately extracts the spatial location, spatial range and construction time of new reservoirs, accurately grasps the time information of reservoir water storage activities, improves observation accuracy, fills the gap in related methods in this field, and provides important scientific and technological support for hydrological data updating, hydrological analysis and research, river water and sediment transport research, river runoff and its change research, regional climate change, government decision-making, environmental monitoring and other fields. BRIEF DESCRIPTION OF THE DRAWINGS

[0017] The drawings are not intended to be drawn to scale. In the drawings, each identical or nearly identical component shown in various figures may be represented by the same reference numeral. For clarity, not every component is labeled in every figure. Embodiments of various aspects of the present invention will now be described by way of example and with reference to the accompanying drawings, in which: Figure 1 This is a sample area diagram provided by an embodiment of the present invention.

[0018] Figure 2 It is the flow chart of the algorithm of the present invention.

[0019] Figure 3 It is the maximum range distribution map of the newly built reservoir extracted by the present invention based on the water body frequency difference.

[0020] Figure 4 This is the trend information and potential mutation point situation obtained by the embodiment of the present invention based on the MK test.

[0021] Figure 5 This is the extraction result of the newly built reservoir in the research sample area provided by the embodiment of the present invention.

[0022] In the above-mentioned Figures 1-5, the coordinates, symbols or other expressions expressed in English are all well known in the art and will not be described in detail in this example. DETAILED DESCRIPTION

[0023] The specific implementation of the present invention is further described in detail below in conjunction with the accompanying drawings and examples. The following implementation cases are used to illustrate the present invention, but are not used to limit the scope of the present invention.

[0024] The embodiment of the present invention takes Yunnan, Guizhou and Sichuan provinces in the southwest of my country as examples to further illustrate the method of the present invention.

[0025] like Figure 1 As shown, this study area is located in Yunnan, Guizhou and Sichuan, with a total area of ​​more than 500,000 square kilometers. It is home to the Yangtze River system, the Pearl River system, the Mekong River system, etc. It has abundant water resources and numerous rivers.

[0026] like Figure 2 FIG. 1 is a flow chart of the present invention. This embodiment includes the following steps: Step 1: Obtain multi-year remote sensing image data of the study area. The downloaded image data should contain all available continuous images within the longest possible time span.

[0027] In this example, Landsat5 TM data from 1984 to 2011, Landsat7 ETM+ data from 2011 to 2013, and Landsat8 OLI data from 2013 to 2020 were obtained in the study area.

[0028] Step 2: Based on the normalized water index NDWI algorithm proposed by McFeeters, the NDWI is calculated pixel by pixel using the green light band and the near-infrared band, and the appropriate threshold is automatically determined to distinguish water and land pixels in combination with the edge Otsu algorithm. The multi-temporal water frequency is extracted pixel by pixel. Taking 1984-2020 as an example, it is divided into two statistical decades, 1984-1999 and 2000-2020, and the cumulative frequency of water in the two periods is calculated respectively: the number of times each pixel is detected as a water body in a year (Water Detection Count, WDC) and the number of times it is effectively monitored in a year (ValidObservation Count, VOC) are calculated, and the annual water frequency (Water Occurrence, WO) of each pixel is obtained by dividing the two: (1) The mean of the annual water body frequency in each period is calculated to obtain the cumulative frequency of water bodies in the two periods. The changes in the cumulative frequency of water bodies in different periods are compared and compared with the threshold to determine whether the pixel is a new reservoir pixel after 2000, screen out new reservoirs that may exist, and determine the maximum boundary range. For the setting of the threshold, in order to exclude receding water bodies and increased seasonal water bodies, while minimizing the omission of new permanent water bodies that may be caused by setting too high a threshold, retain as many new reservoir pixels as possible, and set the threshold to 50%, that is, compared with 1984-1999, the cumulative frequency of water bodies in this pixel increased by 50% between 2000 and 2020.

[0029] like Figure 3 As shown in the figure, the changes in water body frequency between the two periods were screened according to the threshold value. Pixels with a water body frequency change of less than 50% were determined to be non-new permanent water bodies. They were filtered to obtain new reservoir pixels. After matching them with the existing reservoir dataset, a total of 248 new reservoirs were obtained in the study area. The maximum water body boundary range of the new reservoirs was obtained based on the reservoir dataset.

[0030] Step 3. The annual maximum water area (MWA) refers to the maximum coverage of each pixel identified as a water body in a year, that is, the sum of the areas of pixels that have been classified as water bodies at least once in a year. Based on the water body frequency of each pixel, water body determination is performed in months. If the pixel is determined to be a water body in a month, the pixel is determined to be a water body, otherwise it is a non-water body. If it is determined to be a water body every month in a year, the pixel is a permanent water body (Permanent Water, PW). If the months determined to be a water body in a year are between 1 and 11 months, the pixel is a seasonal water body (Seasonal Water, SW). Within the maximum water body boundary of the new reservoir obtained in step 2, the permanent water body and seasonal water body areas of each year are added to obtain the annual maximum water body area: (2) Furthermore, the boundary range corresponding to the annual maximum water area of ​​each reservoir is obtained, and the annual water frequency of each pixel within the boundary range of each reservoir is calculated ( WO ) to obtain the annual average water body frequency ( ). Multiply the annual average water body frequency by the annual maximum water body area to get the annual average water area (Water Area, WA): (3) The annual average water area of ​​the newly built reservoirs is calculated, and the start and end times (2000-2020) and the observation interval (1 year) are set to complete the time series construction.

[0031] Step 4: Traverse the values ​​of each observation in the time series of the newly built reservoir. For all missing values, take the average of the two observations before and after the missing value to fill the missing value; for abnormal values, the replacement method is similar. If the missing value / abnormal value appears in the first or last observation, take the average of the three observations before and after the missing value / abnormal value to replace the missing value / abnormal value.

[0032] Step 5: Run the MK mutation test algorithm. When UF k The curve fluctuates within the confidence interval of 0.95, indicating that the trend and mutation of the change curve are not obvious; UF k If the value of is greater than zero, it indicates that the sequence is in an upward trend, otherwise it is in a downward trend; when it exceeds the critical line, it indicates that the upward or downward trend is significant. k and UB k If the two curves intersect between the critical lines, it means that a mutation has occurred. The moment corresponding to the intersection is the time when the mutation begins, and the intersection is the potential mutation point.

[0033] like Figure 4 As shown in the figure, three large reservoirs in the study area were randomly selected as samples of newly built reservoirs: Hongjiadu Reservoir, Xiangjiaba Reservoir and Xiluodu Reservoir for MK mutation test analysis. The average water area of ​​Hongjiadu Reservoir showed a downward trend from 2000 to 2004, and then the area continued to expand; the area of ​​Xiangjiaba Reservoir expanded slightly from 2000 to 2002, slightly decreased in 2003, and slightly expanded from 2003 to 2006. After that, the water area fluctuated and decreased from 2006 to 2012, and the area continued to expand after 2012; the water area of ​​Xiluodu Reservoir has been in a state of growth except for a slight decrease in 2002. In the MK mutation test, by observing the UF k and UB k The intersection of the curves determines the mutation point, and the UF of the three reservoirs k and UB k The two curves intersect in 2005, 2014 and 2014 respectively, so the water storage time of these three reservoirs is 2005, 2014 and 2014.

[0034] Step 6: The breakpoints detected by the MK mutation test algorithm are not all caused by reservoir construction. The substantial reduction of water area caused by other reasons will also cause the MK mutation test algorithm to detect breakpoints. If it is not filtered, its accuracy will be unacceptable. k The size of the change, if the breakpoint after UF k If the value of is greater than zero and continues to rise, it means that the water area after the breakpoint continues to expand. The year of the mutation point can represent the water storage year of the newly built reservoir. Otherwise, it is not the water storage year of the newly built reservoir.

[0035] Figure 5 The method of the present invention is used to delineate the spatial location and spatial range of the newly built reservoirs in the study sample area and extract the water storage time. Taking Hongjiadu Reservoir and Xiangjiaba Reservoir as examples, the extracted reservoir water storage time is compared with the water storage time revealed by the Google Earth historical images of the reservoir dam, showing that the method proposed by the present invention has a better extraction result for the spatial location, spatial range and water storage time of the newly built reservoir.

[0036] The above method can be used to accurately extract the spatial scope and spatial position of the newly built reservoir and accurately determine the water storage time, correctly grasp the temporal and spatial information of the newly built reservoir, improve the accuracy of hydrological analysis and environmental monitoring, and provide important scientific and technological support for the study of hydrological and environmental changes.

[0037] Although the present invention has been disclosed as above with preferred embodiments, it is not intended to limit the present invention. A person with ordinary knowledge in the technical field to which the present invention belongs may make various changes and modifications without departing from the spirit and scope of the present invention. Therefore, the protection scope of the present invention shall be determined by the definition of the claims.

Claims

1. A method for automatic identification of newly built reservoirs based on two-step remote sensing mapping, characterized in that: include: Obtain multi-year remote sensing image data of the study area; Dividing the multi-year remote sensing image data into several time periods according to time sequence, extracting water body pixels from the remote sensing image data, and calculating the pixel water body cumulative frequency of each time period respectively; Compare the pixel water body accumulation frequency of adjacent time periods. If the difference between the pixel water body accumulation frequency of the latter time period and the previous time period exceeds the preset threshold, the pixel is a newly added reservoir pixel. Match the selected new reservoir pixels with the water body dataset to obtain the newly built reservoirs and construct the annual average water area time series of the newly built reservoirs; The water storage time of the newly built reservoir is determined based on the area change trend and area mutation point of the newly built reservoir.

2. The method according to claim 1, characterized in that: The method for extracting water body pixels from remote sensing image data is: calculating the normalized water index NDWI for the remote sensing image, and extracting water body pixels in combination with a threshold segmentation method.

3. The method according to claim 1, characterized in that: The pixel water body cumulative frequency is the annual water body frequency mean of the pixel in the corresponding time period.

4. The method according to claim 1, characterized in that: The time scale of the period is ≥ 10 years.

5. The method according to claim 1, characterized in that The construction method of the annual average water area time series of the newly built reservoir is: Based on the water body dataset, the maximum water body boundary range of the newly built reservoir is obtained, and the annual maximum water body area of ​​each reservoir is determined within the annual maximum water body boundary range; Calculate the annual water body frequency mean of each pixel within the boundary range corresponding to the annual maximum water body area to obtain the annual average water body frequency; Calculate the product of the annual average water body frequency and the annual maximum water body area to obtain the annual average water body area; Calculate the annual average water area year by year and obtain the annual average water area time series.

6. The method according to claim 1, characterized in that The annual maximum water area is the sum of the areas of permanent water bodies and seasonal water bodies within the boundary of the annual maximum water body.

7. The method according to claim 1, characterized in that It also includes the processing of missing values ​​and outliers for the annual average water area time series, including: Traverse the value of each observation in the annual average water body time series, and for all missing values ​​or outliers, take the average of the two observations before the missing value or outlier and the two observations after the missing value or outlier to fill the missing values ​​or outliers; If missing values ​​or outliers appear in the first or last observation, the mean of several observations before or after the missing value or outlier is used to replace the missing value or outlier.

8. The method according to claim 1, characterized in that The MK mutation test method was used to analyze the reservoir area change trend and find potential mutation points in the annual average water area time series.

9. The method according to claim 1, characterized in that: Compare the changes in the area of ​​the newly built reservoir before and after the mutation point. If the water area shows an increasing trend after the mutation point, the year of the mutation point is the water storage year of the newly built reservoir; otherwise, the mutation point is not the water storage year of the newly built reservoir.

10. The method according to claim 1, characterized in that The spatial resolution of the remote sensing image data is ≤30m.

Images