Method for rasterizing and designing rainstorm intensity based on time decomposition theory

Abstract

The invention discloses a method for rasterizing and designing rainstorm intensity based on a time decomposition theory, and the method comprises the steps: constructing a time decomposition model, decomposing the day-by-day rainfall raster data over the years into five-minute-by-minute rainfall raster data, deducing the frequency distribution of each historical rainstorm amount of each raster unit, and calculating the rainstorm intensity according to the frequency distribution. The method comprises the following steps: constructing a rainstorm intensity formula parameter optimization model based on a system differential response theory, completing rainstorm intensity formula compilation of a grid unit scale, and finally interpolating each grid unit parameter by adopting a spatial difference technology to obtain a grid data set of the rainstorm intensity formula parameters in a drainage basin range. The method can effectively solve the difficult problem of designing the rainstorm intensity in a region lacking minute-level high-resolution rainfall data, gridding of rainstorm intensity formula parameters is realized, the spatial difference of the designed rainstorm intensity is fully embodied, and in addition, the method has the advantages that the design efficiency is improved. The coupling nonparametric kernel density estimation method and the rainstorm intensity formula parameter optimization model based on the system differential response theory can significantly improve the rainstorm intensity calculation precision.

Description

Technical Field

[0001] The invention belongs to the technical field of rainstorm formula design, and in particular relates to a method for rasterizing rainstorm intensity design based on time decomposition theory. Background Art

[0002] The calculation of rainstorm intensity is of great significance to urban planning and flood control projects. Designing a rainstorm intensity calculation formula helps to understand the rainstorm characteristics and rainfall patterns in specific areas, so as to scientifically plan and design urban drainage pipelines and facilities, improve urban drainage and flood prevention capabilities, and reduce the risk of urban waterlogging. At present, the derivation of the design rainstorm intensity has formed a main approach of compiling a rainstorm intensity formula based on rainstorm frequency analysis.

[0003] The current method for compiling a rainstorm intensity formula based on rainstorm frequency analysis still has the following defects: First, this method requires a high temporal resolution of rainstorm data (minute level), but the current rainfall data in most areas is only on a daily or even monthly time scale. Occasionally, some special stations have hourly or minute-by-minute rainfall data, but the years are not long, which cannot meet the needs of compiling a rainstorm intensity formula; second, due to data limitations, the traditional rainstorm intensity formula is still compiled based on the rainstorm data of one or several stations, and other regions directly use the formula of the station, which fails to take into account the spatial heterogeneity of rainstorms and reduces the intensity of rainstorms. The reliability of the rain intensity calculation results; third, when the traditional rainstorm intensity formula is compiled, the rainstorm frequency analysis mainly estimates the frequency distribution characteristics of rainstorms based on certain frequency distribution models (such as P-Ⅲ distribution, exponential distribution and Gumbel distribution, etc.), which has a certain subjective blindness in the selection of distribution models, and the different parameter combinations of different distribution models also bring great uncertainty to the frequency distribution results; fourth, the determination of traditional rainstorm intensity parameters is usually optimized through intelligent algorithms, which has the problems of multiple local optimal solutions and different parameters with the same effect. It cannot be guaranteed that the obtained parameter values ​​are their true values, resulting in great uncertainty in the results.

[0004] Therefore, it is urgently necessary to adopt scientific and efficient methods to obtain spatially distributed rainfall data with high temporal resolution at the minute level, and then derive the rasterized rainstorm intensity formula that takes spatial differentiation into account; at the same time, for the design rainstorm intensity analysis in the new era, it is necessary to improve the accuracy of the design results while improving the temporal resolution and spatial representativeness of the data. Summary of the invention

[0005] The purpose of the present invention is to solve the problems in the prior art and to propose a method for rasterizing and designing rainstorm intensity based on time decomposition theory.

[0006] In order to achieve the above object, the present invention adopts the following technical solutions: A method for designing rainstorm intensity based on gridding of time decomposition theory includes the following steps: S1. Collection and preprocessing of rainfall data: including daily rainfall grid data of the study area over the years and five-minute rainfall extract data of hydrological and meteorological stations, and data preprocessing; S2. Construction of time decomposition model and time decomposition of rainfall data: First, the five-minute rainfall excerpt data of the station measured in step 1 is used to calibrate and verify the time decomposition model, and then the calibrated time decomposition model is used to decompose the daily rainfall grid data into five-minute rainfall grid data; S3, sampling of rainstorm samples for each duration of grid unit: combining the five-minute rainfall grid data of step S2 with GIS pixel statistics technology, fishing net generation technology and annual maximum sampling method to obtain the maximum value of rainstorm amount for each grid unit for multiple durations over the years; S4. Derivation of the frequency distribution of heavy rain of each duration in grid cells: using the non-parametric kernel density estimation method to deduce the frequency distribution of the heavy rain of each duration in each grid cell in step S3, and then deduce the relationship between the intensity-duration-recurrence period of the heavy rain of the grid cell; S5. Optimization of grid unit rainstorm intensity formula parameters and determination of rainstorm intensity formula: Based on the obtained rainstorm intensity-duration-return period relationship, a rainstorm intensity formula parameter optimization model based on system differential response theory is constructed, and the parameters of each grid unit rainstorm intensity formula are selected to complete the compilation of the grid unit scale rainstorm intensity formula; S6. Generation of grid data atlas of designed rainstorm intensity in the watershed: Using spatial interpolation technology, the parameters of each grid cell are interpolated to obtain a grid data set of parameters of the rainstorm intensity formula within the watershed.

[0007] Preferably, in step S3, the maximum heavy rainfall in each grid unit over the years, which lasts for 15 periods of 5 min, 10 min, 15 min, 20 min, 30 min, 45 min, 60 min, 90 min, 120 min, 180 min, 240 min, 360 min, 540 min, 720 min and 1440 min, is obtained.

[0008] Preferably, the preprocessing in step S1 includes performing missing value check, consistency check and white noise processing on the rainfall excerpt data; for the daily rainfall grid data of previous years, the grid data whose rainfall values ​​of all grid cells are 0 are eliminated.

[0009] Preferably, the model in step S2 is a time decomposition model based on random cascade.

[0010] Preferably, in step S3, ARCGIS fishing net technology is used to obtain the rainfall of each grid unit. On the basis of the site scale weight, the grid scale weight is assigned based on the Thiessen polygon theory. The range area controlled by each site is determined by the Thiessen polygon method. According to the hydrological transplantation theory, it is believed that the rainfall time distribution process of all grid units within the range area controlled by the site is similar to that of the site. Therefore, the rainfall time division pattern and its probability of occurrence in the grid unit can be obtained. Finally, the corresponding weight is determined and time decomposition is performed according to the same method to obtain the rainfall raster data of each duration.

[0011] Preferably, in step S4, a non-parametric kernel density estimation method is used to deduce the frequency distribution of the duration of heavy rainfall in each grid unit. First, a kernel function is used to deduce the probability density function, and then the marginal distribution function of the duration of heavy rainfall samples in the grid unit is derived based on the probability density function.

[0012] When the traditional rainstorm intensity formula is compiled, multiple distribution models such as P-III distribution, exponential distribution and Gumbel distribution are used to estimate the frequency distribution characteristics of rainstorms. Different models may give different results. The present invention uses the kernel density estimation method to avoid the uncertainty brought by multiple different models.

[0013] Based on the frequency distribution results of heavy rainfall, the series recurrence period p is set, which can be 2 years, 3 years, 5 years, 10 years, 20 years, 30 years, 50 years, and 100 years. Then, based on the given recurrence period combined with the frequency distribution model, the heavy rainfall intensity i of each duration t under a specific recurrence period is deduced, and then the relationship between heavy rainfall intensity, duration, and recurrence period of the grid unit, p~i~t, is deduced.

[0014] Preferably, in step S5, based on the obtained rainstorm intensity-duration-return period relationship, a rainstorm intensity formula is determined, which is expressed as:

[0015] In the formula, i is the design rainstorm intensity, in mm / min, t is the rainstorm duration, in min; n is the rainstorm attenuation coefficient; b is the regional parameter; A is the rainfall force, in mm / min.

[0016] Preferably, in step S6, based on the grid data set of the rainstorm intensity formula parameters, the design rainstorm intensity grid data under specific design conditions is further calculated, and the specific conditions are a certain duration and a certain recurrence period.

[0017] Compared with the prior art, the present invention has the following beneficial effects: The present invention can effectively solve the problem of designing rainstorm intensity in areas where minute-level high-resolution rainfall data are missing, and realizes the gridding of rainstorm intensity formula parameters, which fully reflects the spatial differences of design rainstorm intensity. In addition, the coupled nonparametric kernel density estimation method and the rainstorm intensity formula parameter optimization model based on system differential response theory can significantly improve the calculation accuracy of rainstorm intensity. BRIEF DESCRIPTION OF THE DRAWINGS

[0018] Figure 1 It is a flowchart of the method of the present invention; Figure 2 A schematic diagram of the accuracy evaluation results of the random cascade-based time decomposition model applied to the Getang station; Figure 3 A schematic diagram of the accuracy evaluation results of the random cascade-based time decomposition model applied to Xiaoqiao Station; Figure 4 A schematic diagram of the accuracy evaluation results of the random cascade-based time decomposition model applied to Nanjing Station; Figure 5 A schematic diagram of the accuracy evaluation results of the random cascade-based time decomposition model applied to Dongshan Station; Figure 6 This is a schematic diagram of the calculation results of the 5-60 minute duration heavy rainfall frequency distribution of the 1# grid cell based on the non-parametric kernel density estimation method; Figure 7 This is a schematic diagram of the calculation results of the 5-60min duration heavy rainfall frequency distribution of the 2# grid cell based on the non-parametric kernel density estimation method; Figure 8 This is a schematic diagram of the calculation results of the frequency distribution of heavy rainfall in the 3# grid cell with a duration of 5-60 minutes based on the non-parametric kernel density estimation method; Fig. 9 This is a schematic diagram of the calculation results of the frequency distribution of heavy rainfall in the 1# grid cell from 90 to 1440 minutes based on the non-parametric kernel density estimation method; Fig.10 This is a schematic diagram of the calculation results of the frequency distribution of heavy rainfall in the 2# grid cell from 90 to 1440 minutes based on the non-parametric kernel density estimation method; Fig.11 This is a schematic diagram of the calculation results of the frequency distribution of heavy rainfall in the 3# grid cell from 90 to 1440 minutes based on the non-parametric kernel density estimation method; Fig.12 This is a schematic diagram of the relationship between rainstorm intensity, duration and return period of grid unit 1# in the embodiment; Fig.13 This is a schematic diagram of the relationship between rainstorm intensity, duration and return period for grid unit 2# in the embodiment; Fig.14This is a schematic diagram of the relationship between rainstorm intensity, duration and return period for grid unit 3# in the embodiment; Fig.15 It is a grid data result diagram of parameter A1 of the Nanjing rainstorm intensity formula in the embodiment; Fig.16 It is a grid data result diagram of parameter C of the rainstorm intensity formula in Nanjing in the embodiment; Fig.17 It is a grid data result diagram of parameter b of the rainstorm intensity formula in Nanjing in the embodiment; Fig.18 It is a grid data result diagram of parameter n of the rainstorm intensity formula in Nanjing in the embodiment; Fig.19 This is a schematic diagram of a 24h-1h rainfall time decomposition model used in an embodiment of the present invention; Fig. 20 Schematic diagram of the 60min-5min rainfall time decomposition model used in the embodiment of the present invention. DETAILED DESCRIPTION

[0019] The technical solutions in the embodiments of the present invention will be clearly and completely described below in conjunction with the drawings in the embodiments of the present invention. Obviously, the described embodiments are only part of the embodiments of the present invention, rather than all of the embodiments. Based on the embodiments of the present invention, all other embodiments obtained by ordinary technicians in this field without making creative work are within the scope of protection of the present invention.

[0020] like Figure 1 , a method for designing rainstorm intensity based on time decomposition theory, including the following steps: S1. Collection and preprocessing of rainfall data: including daily rainfall raster data of previous years and five-minute rainfall excerpt data of hydrological and meteorological stations, and data preprocessing.

[0021] S2. Construction of time decomposition model and time decomposition of rainfall data: First, the five-minute rainfall excerpt data measured at the site in step 1 is used to calibrate and verify the time decomposition model based on random cascade. Finally, the calibrated time decomposition model is used to decompose the daily rainfall grid data into five-minute rainfall grid data.

[0022] The temporal decomposition model structure based on random cascade is as follows Fig.19 , 20 As shown in the figure, each level represents a time series with a specific resolution. For example, level 0 represents 24 hours, i.e., daily resolution, and level 4 represents 1 hour resolution. When the daily rainfall time decomposition model is used, it is only necessary to calibrate the determined branch b0 and the assigned weight wi.

[0023] Determine the number of branches \(b_0\). The number of branches \(b_0\) is mainly determined by the time resolution of the initial level and the time resolution of the final level, and comprehensively consider the complexity of time decomposition. The initial level rainfall is daily rainfall (i.e., 24 hours). When decomposed on the hourly scale, it can be finally decomposed to 1h. Therefore, when decomposed on the hourly scale, the number of branches is taken as (2, 2, 2, 3) respectively, that is, four-time time decomposition is carried out, and the number of branches for each decomposition is 2, 2, 2, and 3 respectively, as Fig.19 shown. On the minute scale, the initial level rainfall is 1h (60min), and the final decomposition is to 5min. Therefore, when decomposed on the minute scale, the number of branches is taken as (2, 2, 3) respectively, that is, three-time time decomposition is carried out, and the number of branches for each decomposition is 2, 2, and 3 respectively, as Fig. 20 shown.

[0024] Calibrate and determine the weight \(w_i\). The weight \(w_i\) reflects the distribution characteristics of rainfall in time, and the distribution characteristics of rainfall in time have certain statistical laws. Therefore, by statistically analyzing the site data of measured rainfall at the minute level, the time distribution characteristic pattern of the site and the probability of its occurrence are statistically analyzed. The weight \(w_i\) of the rainfall time decomposition model of a certain site can be determined as follows: When the number of branches \(b_0 = 2\), there are three possible cases of decomposition weights, that is, there are three partitioning modes: the first is \(w_1 = 1\) and \(w_2 = 0\), which is called 1 / 0 partitioning; the second is the opposite, that is, \(w_1 = 0\) and \(w_2 = 1\), called 0 / 1 partitioning; the third possibility is \(w_1 = v\) and \(w_2 = 1 - v\), which is called \(v / 1 - v\) partitioning, where \(0 < v < 1\). Let the probabilities of the 1 / 0, 0 / 1, or \(v / 1 - v\) partitioning modes be \(P(1, 0)\), \(P(0, 1)\), and \(P(v, 1 - v)\) respectively. Generate a uniformly distributed random number \(r\) in the interval \((0, 1)\). If \(r < P(1 / 0)\), then perform 1 / 0 partitioning; if \(P(1 / 0) < r < \{P(1 / 0)+P(0 / 1)\}\), then perform 0 / 1 partitioning; if \(r > P(1 / 0)+P(0 / 1)\), then perform \(v / 1 - v\) partitioning.

[0025] When the number of branches \(b_0 = 3\), there are seven possible cases of decomposition weights, that is, there are seven partitioning modes: the first is \(w_1 = 1\), \(w_2 = 0\), \(w_3 = 0\), which is called 1 / 0 / 0 partitioning; the second is \(w_1 = 0\), \(w_2 = 1\), \(w_3 = 0\), which is called 0 / 1 / 0 partitioning; the third is \(w_1 = 0\), \(w_2 = 0\), \(w_3 = 1\), which is called 0 / 0 / 1 partitioning; the fourth is \(w_1 = v\), \(w_2 = 1 - v\), \(w_3 = 0\), which is called \(v / 1 - v / 0\) partitioning; the fifth is \(w_1 = v\), \(w_2 = 0\), \(w_3 = 1 - v\), which is called \(v / 0 / 1 - v\) partitioning; the sixth is \(w_1 = 0\), \(w_2 = v\), \(w_3 = 1 - v\), which is called \(0 / v / 1 - v\) partitioning; the seventh is \(w_1 = v\), \(w_2 = u\), \(w_3 = 1 - v - u\), which is called \(v / u / 1 - v - u\) partitioning.

[0026] Let the probabilities of the above partition patterns be P(1, 0, 0), P(0, 1, 0), P(0, 0, 1), P(v, 1 - v, 0), P(v, 0, 1 - v), P(0, v, 1 - v), and P(v, u, 1 - v - u), respectively. Similarly, generate a uniformly distributed random number r in the interval (0, 1). If r < P(1, 0, 0), then perform a 1 / 0 / 0 partition; if P(1, 0, 0) < r < P(1, 0, 0) + P(0, 1, 0), then perform a 0 / 1 / 0 partition; if P(1, 0, 0) + P(0, 1, 0) < r < P(1, 0, 0) + P(0, 1, 0) + P(0, 0, 1), then perform a 0 / 0 / 1 partition; if P(1, 0, 0) + P(0, 1, 0) + P(0, 0, 1) < r < P(1, 0, 0) + P(0, 1, 0) + P(0, 0, 1) + P(v, 1 - v, 0), then perform a v / 1 - v / 0 partition; if P(1, 0, 0) + P(0, 1, 0) + P(0, 0, 1) + P(v, 1 - v, 0) < r < P(1, 0, 0) + P(0, 1, 0) + P(0, 0, 1) + P(v, 1 - v, 0) + P(v, 0, 1 - v), then perform a v / 0 / 1 - v partition; if P(1, 0, 0) + P(0, 1, 0) + P(0, 0, 1) + P(v, 1 - v, 0) + P(v, 0, 1 - v) + P(v, 0, 1 - v) < r < (1, 0, 0) + P(0, 1, 0) + P(0, 0, 1) + P(v, 1 - v, 0) + P(v, 0, 1 - v) + P(v, 0, 1 - v) + P(0, v, 1 - v), then perform a 0 / v / 1 - v partition; if r > P(1, 0, 0) + P(0, 1, 0) + P(0, 0, 1) + P(v, 1 - v, 0) + P(v, 0, 1 - v) + P(v, 0, 1 - v) + P(0, v, 1 - v), then perform a v / u / 1 - v - u partition. Here, v and u are determined as random numbers of the empirical distribution. v can take any value in the range (0, 1), but u is conditional on v, so u can be assumed to be any number between (0, 1 - v).

[0027] S3. Sampling of rainstorm samples for each time period of grid cells: Combine the five - minute rainfall grid data in step S2 with GIS pixel statistics technology, fishing net generation technology, and the annual maximum sampling method to obtain the rainstorm amounts for 15 time periods of each grid cell over the years, namely the maximum 5min, 10min, 15min, 20min, 30min, 45min, 60min, 90min, 120min, 180min, 240min, 360min, 540min, 720min, and 1440min.

[0028] Based on the daily rainfall grid data, the ARCGIS fishing net technology is used to obtain the daily rainfall of each grid unit. On the basis of the station scale weight, the grid scale weight is assigned based on the Thiessen polygon theory. Due to the lack of measured minute-level rainfall data on the grid unit, the Thiessen polygon method can be used to determine the range area controlled by each station. According to the hydrological transplantation theory, it is believed that the rainfall time distribution process of all grid units in the station control area is similar to that of the station. Therefore, the rainfall time division pattern in the grid unit and its probability of occurrence P can be obtained. Finally, the corresponding weights are determined and time decomposition is performed according to the same method to obtain 5-minute rainfall grid data.

[0029] S4. Derivation of the frequency distribution of heavy rain of each duration in grid cells: Based on the sampling results of heavy rain samples, the non-parametric kernel density estimation method is used to deduce the frequency distribution of heavy rain of each duration in each grid cell, and then the relationship between the intensity, duration and recurrence period of heavy rain in the grid cell is derived.

[0030] First, the kernel function is used to deduce the probability density function, and then the marginal distribution function of the heavy rain samples of each duration in the grid unit is derived based on the probability density function; assuming m Heavy rain sample points x 1 ，x 2 ，....x m It comes from a univariate continuous population distribution, the population probability density function f(x) , its kernel density estimation function can be expressed as:

[0031] In the formula, yes f(x) The kernel density estimation function is K ( x ) is the kernel function, and , h is called the window width or bandwidth (h>0), j’ Represents the number of rainstorm sample points.

[0032] Based on the frequency distribution results of heavy rainfall, the series recurrence period p is set, which can be 2 years, 3 years, 5 years, 10 years, 20 years, 30 years, 50 years, and 100 years. Then, based on the given recurrence period combined with the frequency distribution model, the heavy rainfall intensity i of each duration t under a specific recurrence period is deduced, and then the relationship between heavy rainfall intensity, duration, and recurrence period of the grid unit, p~i~t, is deduced.

[0033] S5. Optimization of grid unit rainstorm intensity formula parameters and determination of rainstorm intensity formula: Based on the obtained rainstorm intensity-duration-return period relationship, a rainstorm intensity formula parameter optimization model based on system differential response theory is constructed to optimize the parameters of each grid unit rainstorm intensity formula and complete the compilation of the grid unit scale rainstorm intensity formula.

[0034] The rainstorm intensity formula is expressed as:

[0035] In the formula, i is the design rainstorm intensity, in mm / min, t is the rainstorm duration, in min; n is the rainstorm attenuation coefficient; b is the regional parameter; A It is the rain force in mm / min.

[0036] Rain Power A It varies with the return period p, and its relationship with p is often expressed by the following formula:

[0037] In the formula, A 1 , C is the regional parameter, p is the return period, so Can be rewritten as .

[0038] The principle of the rainstorm intensity formula parameter optimization model based on system differential response theory is as follows: Since the rainstorm intensity formula itself is continuous and differentiable in the parameter neighborhood, it can be obtained by sampling Function with parameters The intersection position obtains the true value information of the parameters, assuming that the nonlinear function model can be expressed as:

[0039] Where: θ is the model parameter, X is the model input variable, and:

[0040] Where: q is the total number of parameters, k is the number of independent variables.

[0041] The above nonlinear function model is expanded by a first-order Taylor series, and its expression is:

[0042] Where: j is the number of iterations; For parameters No. j Sub-rate setting value; The deviation between the observed value and the calculated value.

[0043] Assume there are m groups of rainstorm observation samples, that is , and the vector form of the function model is:

[0044]

[0045] And there is Furthermore, it can be deduced that: .

[0046] Through the above process, the optimal parameter value can be obtained by iterative calculation from the initial parameter value. Considering the influence of the high-order Taylor series expansion term, a correction coefficient c is introduced, where c > 1. Thus, we can get:

[0047] In the formula, 0 < B < 1, and it can be determined by the least squares estimation of the function calculated value and the measured value.

[0048] In summary, the solution steps of the rainstorm intensity formula parameter optimization model based on the system differential response theory are as follows: (1) Randomly generate a set of initial parameter values , that is, the initialization of the parameter A 1 , C, b, n ; (2) Combine the initial parameters with the samples, calculate the sensitivity matrix S and for each sample; (3) Determine the search direction and the parameter B ; (4) Calculate and obtain ; (5) Judge whether the convergence principle is satisfied. If it is satisfied, terminate the optimization. If not, let and go back to step (2);至此,确定出每个栅格单元的暴雨强度公式的参数,即完成了每个栅格单元的暴雨公式编制。

[0049] S6. Generation of the grid data atlas of the designed rainstorm intensity of the basin: Using spatial interpolation technology, interpolate the parameters of each grid cell to obtain a grid data set of the rainstorm intensity formula parameters within the basin range.

[0050] Based on the longitude and latitude coordinates where each grid is located, then using spatial interpolation technology, interpolate the rainstorm intensity formula A 1 , C, b, n parameters of each grid cell to obtain a grid data set of the rainstorm intensity formula parameters within the basin range.

[0051] Based on the raster dataset of the rainstorm intensity formula parameters, the rainstorm intensity raster data under specific design conditions (a certain duration and a certain return period) can be further calculated.

[0052] Example

[0053] Taking the calculation of design rainstorm intensity in Nanjing as an example, the method of the present invention is explained: S1. Daily rainfall raster data of China from 1978 to 2021 were collected with a data resolution of 0.1°. The daily rainfall raster data of Nanjing from 1978 to 2021 were obtained by clipping with ACRGIS software. For the five-minute rainfall data, the five-minute rainfall data of 11 rain gauges in Nanjing, including Liuhe Station, Hongshanyaozha Station, Getang Station, Xiaoqiao Station, Nanjing Station, Qinhuai Xinhezha Station, Dongshan Station, Qian'an Village Station, Tianshengqiao Station, Yangjiawan Station, and Gaochun Station, were collected from 1992 to 2021. The basic attributes of the rain gauges are shown in Table 1. The collected rainfall data were preprocessed.

[0054] Table 1 Basic attributes of rainfall stations

[0055] S2. Based on the pre-processed rainfall data of Nanjing City, a temporal decomposition model based on random cascade is constructed. Fig.19 , 20 As shown, the five-minute rainfall excerpt data measured at the above-mentioned 11 rainfall stations were used to calibrate and verify the model. Finally, the calibrated time decomposition model was used to decompose the daily rainfall raster data in Nanjing from 1978 to 2021 into five-minute rainfall raster data.

[0056] The rate period of the random cascade time decomposition model of the present invention is from January 1, 1992 to December 31, 2010, and the verification period is from January 1, 2011 to December 31, 2015. During the rate period, the probability values ​​of three division modes appearing when the number of branches b0=2 are mainly counted at 11 rain gauges in Nanjing, namely, P(1,0), P(0,1) and P(v,1-v), and the probability values ​​of seven division modes appearing when b0=3 are P(1,0,0), P(0,1,0) and P(0,0,1), P(v,1-v,0), P(v,0,1-v), P(0,v,1-v) and P(v,u,1-vu). During the verification period, it is assumed that the probability of each group of division models statistically generated by the rate period remains unchanged, and the daily rainfall of each rain gauge is divided based on the determined division weights. Based on five statistical variables reflecting the nature and characteristics of rainfall, namely, percentage of zero values, volume of a single divided rain box, volume of rainfall events, duration of events, and drought period, the division results were compared with the measured results. The multiple average relative errors of the above five variables between the decomposition model results and the measured results were analyzed.

[0057] Figure 2-Figure 5 Taking Getang Station, Xiaoqiao Station, Nanjing Station and Dongshan Station as examples, the accuracy evaluation results of the time decomposition model are shown. Figure 2-Figure 5 It can be seen that the random cascade-based time decomposition model constructed by the present invention has high accuracy, can well maintain the nature and characteristics of rainfall, and the relative errors of various statistical variables are within ±15%, which can be used for subsequent rainfall time decomposition.

[0058] The constructed random cascade-based time decomposition model is used to decompose the daily rainfall raster data of Nanjing for a total of 44 years from 1978 to 2021 into five-minute rainfall data, where each daily rainfall raster data file can be finally decomposed into 288 five-minute rainfall raster data. The raster data is a tiff format file, and the data volume is large, so it will not be shown in this embodiment.

[0059] S3. After obtaining the five-minute rainfall grid data of Nanjing, GIS pixel statistics technology, fishing net generation technology and annual maximum sampling method are combined to obtain the maximum 15-hour heavy rainfall in each grid unit in the past years, and finally complete the sampling of heavy rainfall samples of each duration in the grid unit.

[0060] S4. After the sampling is completed, the frequency distribution of the duration of heavy rainfall in each grid unit is calculated using the non-parametric kernel density estimation method. In this embodiment, three grid units are randomly selected, which are marked as 1#, 2#, and 3#, respectively, to display the frequency distribution results, where Figure 6-Figure 8 The following are the calculation results of the frequency distribution of heavy rainfall in grid units 1#~3# for 5~60 minutes. Figure 9-11 This is the calculation result of the frequency distribution of heavy rainfall in grid cells 1#~3# from 90 to 1440 minutes. Figure 6-Figure 11 It can be seen that the frequency curves of each grid unit at different durations are generally coordinated, without any unreasonable phenomena such as intersection, and can be used for the analysis of rainstorm intensity formula.

[0061] Based on the frequency distribution of each grid unit in Nanjing, the relationship between rainstorm intensity, duration and return period of each grid unit in Nanjing is further deduced. Figure 12-14 This is the relationship diagram of rainstorm intensity, duration and return period for grid cells 1# to 3#. The rainstorm intensity and duration are logarithmically processed in the figure.

[0062] S5. The rainstorm intensity formula parameter optimization model based on the system differential response theory determines the rainstorm intensity parameters corresponding to each grid unit. Finally, the grid data of the four parameters A1, C, b and n of the rainstorm intensity formula in Nanjing can be obtained. The results are as follows: Figure 15-18 shown.

[0063] S6. Based on the grid data of the rainstorm intensity formula parameters, the design rainstorm intensity grid data under specific design conditions (a certain duration and a certain return period) is calculated. This embodiment shows the rainstorm intensity grid data results with a return period of 3 years, as shown in Table 2.

[0064] Table 2 Design rainstorm intensity of each grid unit when the return period is 3 years (unit: mm / min)

[0065] The above description is only a preferred specific implementation manner of the present invention, but the protection scope of the present invention is not limited thereto. Any technician familiar with the technical field can make equivalent replacements or changes according to the technical scheme and inventive concept of the present invention within the technical scope disclosed by the present invention, which should be covered by the protection scope of the present invention.

Claims

1. A method for designing rainstorm intensity based on rasterization of time decomposition theory, characterized in that: The steps include: S1. Collection and preprocessing of rainfall data: including daily rainfall grid data of the study area over the years and five-minute rainfall extract data of hydrological and meteorological stations, and data preprocessing; S2. Construction of time decomposition model and time decomposition of rainfall data: First, the five-minute rainfall excerpt data of the station measured in step 1 is used to calibrate and verify the time decomposition model, and then the calibrated time decomposition model is used to decompose the daily rainfall grid data into five-minute rainfall grid data; S3, sampling of rainstorm samples for each duration of grid unit: combining the five-minute rainfall grid data of step S2 with GIS pixel statistics technology, fishing net generation technology and annual maximum sampling method to obtain the maximum value of rainstorm amount for each grid unit for multiple durations over the years; S4. Derivation of the frequency distribution of heavy rain of each duration in grid cells: using the non-parametric kernel density estimation method to deduce the frequency distribution of the heavy rain of each duration in each grid cell in step S3, and then deduce the relationship between the intensity-duration-recurrence period of the heavy rain of the grid cell; S5. Optimization of grid unit rainstorm intensity formula parameters and determination of rainstorm intensity formula: Based on the obtained rainstorm intensity-duration-return period relationship, a rainstorm intensity formula parameter optimization model based on system differential response theory is constructed, and the parameters of each grid unit rainstorm intensity formula are selected to complete the compilation of the grid unit scale rainstorm intensity formula; S6. Generation of grid data atlas of designed rainstorm intensity in the watershed: Using spatial interpolation technology, the parameters of each grid cell are interpolated to obtain a grid data set of parameters of the rainstorm intensity formula within the watershed.

2. The method for designing rainstorm intensity based on time decomposition theory according to claim 1, characterized in that: In step S3, the maximum rainfall in each grid cell over the years, which lasts for 5 min, 10 min, 15 min, 20 min, 30 min, 45 min, 60 min, 90 min, 120 min, 180 min, 240 min, 360 min, 540 min, 720 min and 1440 min, a total of 15 periods, is obtained.

3. The method for designing rainstorm intensity based on time decomposition theory according to claim 1, characterized in that: The preprocessing in step S1 includes performing missing value test, consistency test and white noise processing on the rainfall excerpt data; for the daily rainfall grid data of previous years, the grid data with the rainfall of all grid cells being 0 are eliminated.

4. The method for designing rainstorm intensity based on time decomposition theory according to claim 1, characterized in that: The model in step S2 is a time decomposition model based on random cascade.

5. The method for designing rainstorm intensity based on time decomposition theory according to claim 1 is characterized by: In step S3, ARCGIS fishing net technology is used to obtain the rainfall of each grid unit. On the basis of the site scale weight, the grid scale weight is assigned based on the Thiessen polygon theory. The range area controlled by each site is determined by the Thiessen polygon method. According to the hydrological transplantation theory, it is believed that the rainfall time distribution process of all grid units within the range controlled by the site is similar to that of the site. Therefore, the rainfall time division pattern and its probability of occurrence in the grid unit are deduced. Finally, the corresponding weight is determined and time decomposition is performed according to the same method to obtain the rainfall raster data of each duration.

6. The method for designing rainstorm intensity based on time decomposition theory according to claim 5, characterized in that: In step S4, a non-parametric kernel density estimation method is used to deduce the frequency distribution of the heavy rainfall of each grid unit. First, the probability density function is derived using the kernel function, and then the marginal distribution function of the heavy rainfall samples of each grid unit is derived based on the probability density function.

7. The method for designing rainstorm intensity based on time decomposition theory according to claim 6, characterized in that: In step S5, based on the obtained rainstorm intensity-duration-return period relationship, a rainstorm intensity formula is determined, which is expressed as: In the formula, i is the design rainstorm intensity, in mm / min, t is the rainstorm duration, in min; n is the rainstorm attenuation coefficient; b is the regional parameter; A is the rainfall intensity, the unit is mm / min.

8. The method for designing rainstorm intensity based on time decomposition theory according to claim 7, characterized in that: In step S6, based on the grid data set of the rainstorm intensity formula parameters, the design rainstorm intensity grid data under specific design conditions is calculated, and the specific conditions are a certain duration and a certain recurrence period.

Images