Earth-rock dam leakage positioning method based on joint inversion and data fusion

Abstract

The invention relates to the technical field of data processing, discloses an earth and rockfill dam leakage positioning method based on joint inversion and data fusion, and aims to solve the problems of limitation and non-uniqueness of an existing single geophysical prospecting method in leakage hidden danger positioning. Converting the first apparent resistivity into a first apparent resistivity; second apparent resistivity distribution of the earth-rock dam is obtained through a high-density electrical method; acquiring electromagnetic wave velocity distribution by adopting a ground penetrating radar method, and converting the electromagnetic wave velocity distribution into third apparent resistivity; performing data standardization processing on the apparent resistivity to obtain standardized apparent resistivity data; performing dimension reduction fusion on the standardized apparent resistivity data corresponding to each apparent resistivity by adopting a principal component analysis method to obtain simulated apparent resistivity distribution; and determining the position of the leakage hidden danger channel in the earth-rock dam based on the pseudo-apparent resistivity distribution. The method improves objectivity, accuracy and reliability of leakage positioning, and is particularly suitable for danger removal and reinforcement of earth and rockfill dams.

Description

Technical Field

[0001] This invention relates to the field of data processing technology, specifically to a method for locating seepage in earth-rock dams based on joint inversion and data fusion. Background Technology

[0002] Earth-rock dams, as important water conservancy projects, play a vital role in my country's national economic development. However, over time, many earth-rock dams inevitably develop hidden dangers such as seepage after long-term service, which can lead to dam failure and endanger downstream safety in severe cases. Therefore, it is necessary to identify the locations of seepage inlets and outlets and channels in order to carry out subsequent reinforcement and safety measures to ensure the safety of the dams.

[0003] Traditional methods for detecting seepage in earth-rock dams mainly rely on manual inspections, borehole sampling, and water pressure tests. These methods are experience-based, highly subjective, and difficult to detect hidden dangers inside the dam. They also cannot achieve rapid, comprehensive, and accurate location of seepage channels in the dam, and therefore cannot meet the needs of modern dam safety monitoring.

[0004] With the development of geophysical technology, non-destructive testing methods such as surface wave method, high-density electrical resistivity tomography (EDT), and ground-penetrating radar (GPR) have been widely used in dam seepage detection. However, a single geophysical method has obvious limitations: surface wave method, based on the elastic differences of the medium, is sensitive to changes in wave velocity, but its resolution decreases with increasing depth and it is difficult to distinguish water-bearing fissures from other low-velocity anomalies; high-density EDT, based on the electrical differences of the medium, is sensitive to water-bearing areas, but it is greatly affected by topography and noise interference, and suffers from volume effects, resulting in limited resolution; GPR has high resolution, but its penetration depth in highly conductive media is limited, and its signal attenuation is severe in water-rich areas. More importantly, geophysical inversion problems have an inherent "non-uniqueness," meaning that different subsurface medium models may produce the same or similar observational data. This non-uniqueness leads to multiple interpretations of the results from a single method, reducing the reliability of seepage location.

[0005] To address the limitations and multiple solutions inherent in single-method approaches, integrated geophysical exploration techniques have been proposed. However, current integrated geophysical exploration primarily focuses on method combinations, employing multiple methods sequentially, interpreting each independently, and then conducting comparative analysis. While this approach improves detection reliability to some extent, it suffers from the following problems: first, the data formats, resolutions, and detection depths of different methods are inconsistent, making direct comparison difficult; second, manual comparative analysis is highly subjective and lacks quantitative fusion methods; and third, information utilization is low, failing to leverage the advantages of multi-parameter joint constraints. Summary of the Invention

[0006] This invention aims to address the limitations and non-uniqueness of existing single geophysical exploration methods in locating leakage hazards, and proposes a leakage location method for earth-rock dams based on joint inversion and data fusion.

[0007] The technical solution adopted by the present invention to solve the above-mentioned technical problems is as follows:

[0008] A method for locating seepage in earth-rock dams based on joint inversion and data fusion, the method comprising:

[0009] Step 1: Obtain the shear wave velocity distribution of the earth-rock dam using the surface wave method, and convert the shear wave velocity into the first apparent resistivity based on the Faust formula;

[0010] Step 2: Obtain the second apparent resistivity distribution of the earth-rock embankment using high-density electrical resistivity tomography.

[0011] Step 3: Use ground-penetrating radar to obtain the electromagnetic wave velocity distribution of the earth-rock dam, and convert the electromagnetic wave velocity into a third apparent resistivity based on the relationship between resistivity and relative permittivity.

[0012] Step 4: Perform data standardization processing on the first apparent resistivity, the second apparent resistivity and the third apparent resistivity respectively, including scale standardization and numerical standardization, to obtain standardized apparent resistivity data with consistent resolution and numerical range.

[0013] Step 5: Use principal component analysis to perform dimensionality reduction and fusion on the standardized apparent resistivity data corresponding to the first, second, and third apparent resistivity to obtain the pseudo-apparent resistivity distribution.

[0014] Step 6: Based on the apparent resistivity distribution, determine the location of potential leakage channels in the earth-rock embankment.

[0015] Further, in step 1, the method for converting the first apparent resistivity specifically includes:

[0016] Step 11: Use the damped least squares method to invert the surface wave detection data to obtain the shear wave velocity distribution of the earth-rock dam;

[0017] Step 12: Substitute the Garden formula into the standard Faust formula for conversion, establish the improved Faust formula, and use known geological data to fit the stratigraphic parameters in the improved Faust formula using the least squares method.

[0018] Step 13: Convert the shear wave velocity of the earth-rock dam into the first apparent resistivity based on the improved Faust formula, which is as follows:

[0019] ;

[0020] in, Indicates the first apparent resistivity. Indicates the shear wave velocity. Represents Poisson's ratio. , , This represents the stratigraphic parameters obtained by fitting known geological data.

[0021] Further, in step 2, the method for obtaining the second apparent resistivity distribution includes:

[0022] Step 21: Construct the geoelectric model of the earth-rock dam, whose model parameter vector is... ,in, Indicates the first True resistivity of each grid cell Indicates the total number of model parameters. Indicates matrix transpose;

[0023] Step 22: Perform forward modeling on the geoelectric model to obtain theoretical apparent resistivity data. ,in, Indicates the index of the observation point. Indicates the total number of observation data;

[0024] Step 23: Calculate theoretical apparent resistivity data Compared with measured apparent resistivity data residuals: ;

[0025] Step 24: Construct the error function To measure the overall deviation between theoretical and measured data, this error function is defined as the sum of the squares of the residuals at all observation points:

[0026] ;

[0027] Step 25: Iteratively adjust the model parameter vector using the least squares method. Solve for the correction amount of the model parameter vector. and in accordance with Update the geoelectric model, in which, and They represent the first The model parameter vector and correction value corresponding to the next iteration;

[0028] Step 26: Repeat steps 22 to 25 until the error function values ​​of two adjacent iterations meet the convergence condition. At this point, the true resistivity distribution corresponding to the final model parameter vector is the second apparent resistivity distribution.

[0029] Furthermore, in step 3, the method for converting the third apparent resistivity specifically includes:

[0030] Step 31: Based on the B-scan data obtained by ground penetrating radar, identify the hyperbola formed by the reflected wave of the underground target;

[0031] Step 32: Extract the electromagnetic wave travel time information of the vertex position and each observation point position of the hyperbola using the hyperbola fitting method, and calculate the propagation speed of electromagnetic waves in the underground medium. The calculation formula is as follows:

[0032] ;

[0033] in, This indicates the propagation speed of electromagnetic waves in underground media. Indicates the antenna position. This indicates the horizontal position corresponding to the vertex of the hyperbola. Indicates position When the reflected wave from the target location travels, Represents the vertex of the hyperbola Travel time of the reflected wave from the target at the location;

[0034] Step 33: Based on the formula for the propagation speed of electromagnetic waves in a medium and the assumption that the relative permeability of soil and rock is equal to 1, the relationship between the electromagnetic wave speed and the relative permittivity of the medium is derived.

[0035] The formula for the propagation speed of electromagnetic waves in a medium is as follows:

[0036] ;

[0037] The relationship between the propagation speed of electromagnetic waves and the relative permittivity of the medium is as follows:

[0038] ;

[0039] in, This represents the speed of light in a vacuum. Represents the relative permittivity of the medium. Indicates the relative magnetic permeability of soil and rock mass;

[0040] Step 34: Using the empirical relationship between true resistivity and relative permittivity obtained by fitting experimental data of soil and rock, convert the relative permittivity into true resistivity according to the empirical relationship.

[0041] The empirical relationship between true resistivity and relative permittivity is as follows:

[0042] ;

[0043] in, Indicates true resistivity;

[0044] Step 35: Substitute the relationship between electromagnetic wave velocity and relative permittivity of the medium obtained in Step 33 into the empirical relationship between true resistivity and relative permittivity obtained in Step 34, and eliminate intermediate variables. The propagation speed of electromagnetic waves is directly converted into the final expression of the third apparent resistivity. The propagation speed of electromagnetic waves obtained in step 32 is converted into the third apparent resistivity through this final expression.

[0045] The propagation speed of electromagnetic waves can be directly converted into the final expression for the third apparent resistivity as follows:

[0046] ;

[0047] in, This represents the third apparent resistivity.

[0048] Furthermore, in step 4, the scale is uniformly implemented using a linear sampling method, specifically including:

[0049] Step 41: Resample the data of the first apparent resistivity, the second apparent resistivity and the third apparent resistivity to a grid with the same preset spatial resolution.

[0050] Step 42: For each grid cell, its scale-unified apparent resistivity is determined by a weighted average of the apparent resistivity of all sampling points within the grid cell and their weights.

[0051] ;

[0052] in, Indicates the first Apparent resistivity after unifying the scale of each grid cell Indicates the first cell within the grid. Apparent resistivity at each sampling point Indicates the first cell within the grid. The weight of each sampling point This indicates the number of sampling points within a grid cell.

[0053] Furthermore, in step 4, the numerical values ​​are uniformly transformed using a linear transformation method to unify the apparent resistivity of each grid cell to the range of [0, 100], as shown in the following formula:

[0054] ;

[0055] in, Indicates the first Standardized apparent resistivity after unifying the values ​​of each grid cell and These represent the maximum and minimum values ​​of apparent resistivity after all grid cells have been scaled uniformly.

[0056] Furthermore, in step 5, the specific process of the principal component analysis method includes:

[0057] Step 51: Use the standardized apparent resistivity data as the original feature vector. The feature vector of each grid cell is The first feature The first characteristic represents the normalized apparent resistivity corresponding to the first apparent resistivity, and the second characteristic... The third characteristic represents the normalized apparent resistivity corresponding to the second apparent resistivity. This represents the standardized apparent resistivity corresponding to the third apparent resistivity.

[0058] Step 53: Decentralize the original feature matrix using the following formula:

[0059] ;

[0060] in, Indicates the first The first grid cell of the first grid cell The value of each feature after decentralization Indicates the first The first grid cell of the first grid cell One characteristic, Indicates the number of grid cells;

[0061] Step 54: Based on the decentralized values, calculate the 3×3 covariance matrix and solve the characteristic equation of the covariance matrix to obtain the eigenvalues ​​of the covariance matrix and its corresponding eigenvectors.

[0062] Step 54: Select the eigenvector corresponding to the largest eigenvalue to form a projection space, and project the original eigenvector into the projection space to obtain the apparent resistivity.

[0063] Furthermore, the apparent resistivity is obtained by projecting the original eigenvectors onto this projection space, and the projection formula is as follows:

[0064] ;

[0065] in, Indicates the first Apparent resistivity of each grid cell This represents the eigenvector corresponding to the largest eigenvalue. This indicates the matrix transpose.

[0066] Furthermore, the principal component analysis method also includes:

[0067] Calculate the contribution rate and cumulative contribution rate of the principal components corresponding to each feature, and evaluate the amount of original data information retained after dimensionality reduction based on the contribution rate and cumulative contribution rate;

[0068] The formula for calculating the contribution rate is as follows:

[0069] ;

[0070] in, Indicates the first Each feature corresponds to the contribution rate of the principal component. Indicates the first Each feature corresponds to an eigenvalue of a principal component. Indicates the summation index;

[0071] The formula for calculating the cumulative contribution rate is as follows:

[0072] ;

[0073] in, Indicates the preceding The cumulative contribution rate of each feature corresponds to the principal component.

[0074] Furthermore, the location of potential seepage channels within the earth-rock embankment was determined, specifically including:

[0075] Step 61: Based on the apparent resistivity distribution, identify areas with apparent resistivity below the abnormal threshold as low-resistivity abnormal areas, which correspond to areas with potential leakage risks due to high moisture content.

[0076] Step 62: Based on the continuous distribution of potential leakage areas, reconstruct the three-dimensional spatial morphology of the leakage channels using spatial interpolation methods to determine the direction, depth, and range of the leakage channels.

[0077] Step 63: Based on the structural characteristics of the dam, identify the inlet and outlet positions of the seepage channels in the reconstructed three-dimensional morphology of the seepage channels.

[0078] The beneficial effects of this invention are as follows: The seepage location method for earth-rock dams based on joint inversion and data fusion provided by this invention achieves multi-parameter joint constraints of elastic, electrical, and electromagnetic parameters through the organic combination of surface wave method, high-density electrical resistivity tomography (EDS), and ground-penetrating radar (GPR) method, effectively reducing the non-uniqueness of geophysical inversion; the use of principal component analysis for data fusion retains the effective information of each method while highlighting common anomaly characteristics, improving the ability to distinguish seepage channels; and the use of scale and numerical standardization solves the problem of inconsistent data formats and resolutions among different methods, laying the foundation for data fusion. This invention establishes a complete quantitative process from data acquisition to hazard identification, reduces human interference, improves the objectivity, accuracy, and reliability of seepage location, and provides a scientific basis for the reinforcement and repair of earth-rock dams. Attached Figure Description

[0079] Figure 1 A flowchart illustrating the method for locating seepage in earth-rock dams based on joint inversion and data fusion, provided as an example.

[0080] Figure 2 A schematic diagram of the hyperbola in the B-scan provided for the embodiment;

[0081] Figure 3 A schematic cross-sectional view of the arrangement of the preset leakage channel 1 provided in the embodiment;

[0082] Figure 4 A schematic cross-sectional view of the arrangement of the preset leakage channel 2 provided in the embodiment;

[0083] Figure 5 A schematic cross-sectional view of the arrangement of the preset leakage channel 3 provided in the embodiment;

[0084] Figure 6 A schematic elevation view of the layout of the preset leakage channels provided in the embodiment;

[0085] Figure 7 A schematic diagram showing the location distribution of the preset leakage channels and their corresponding potential fusion areas provided in the embodiment;

[0086] Figure 8 A schematic diagram showing the location distribution of the preset leakage channel 1 and its corresponding potential hazard area provided in the embodiment;

[0087] Figure 9 A schematic diagram showing the location distribution of the preset leakage channel 2 and its corresponding potential fusion area provided in the embodiment;

[0088] Figure 10 A schematic diagram showing the location distribution of the preset leakage channel 3 and its corresponding potential hazard area provided in the embodiment;

[0089] Explanation of reference numerals in the attached figures:

[0090] - Antenna location - The horizontal position corresponding to the vertex of the hyperbola -Location When the reflected wave from the target location travels, - Hyperbola vertex The travel time of the reflected wave from the target at that location. Detailed Implementation

[0091] To address the limitations and non-uniqueness of existing single geophysical exploration methods in locating leakage hazards, and to achieve joint inversion and data fusion of multiple geophysical data, as well as to solve the problems of comparison difficulties, strong subjectivity, and low information utilization in existing integrated geophysical exploration technologies, this invention proposes a technical solution.

[0092] In this invention, the apparent resistivity of earth-rock embankments is obtained using surface wave method, high-density electrical resistivity method, and ground-penetrating radar method, respectively. Surface wave method is insensitive to saturation but can effectively identify material softening; high-density electrical resistivity method is sensitive to moisture but cannot distinguish between water-bearing fissures and clay interlayers; ground-penetrating radar has high resolution but limited detection depth. Although these three physical fields have different detection principles, they all respond to the same geological target, namely seepage channels, forming a multi-physics coupled response. Based on this, this invention converts the different physical quantities obtained by the three methods into a unified apparent resistivity parameter, establishing a common feature expression. Then, scale and numerical unification are performed on each apparent resistivity, normalizing data from different support volumes to the same geological grid unit and eliminating the influence of dimensions while maintaining the relative relationship of data distribution. Finally, principal component analysis is used to find new features that can best reflect the variance of the original data through linear transformation. While ensuring minimal information loss, the three-dimensional apparent resistivity feature is reduced to a one-dimensional pseudo-apparent resistivity, and the pseudo-apparent resistivity distribution is used to achieve accurate identification and location of seepage channels. This invention establishes a complete quantitative process from data collection to hazard identification, reduces human interference, and improves the objectivity, accuracy, and reliability of leakage location, providing a scientific basis for the reinforcement and repair of earth-rock embankments.

[0093] The technical solutions in this embodiment will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only some embodiments of the present invention, and not all embodiments.

[0094] Figure 1 A flowchart illustrating a method for locating seepage in earth-rock dams based on joint inversion and data fusion is shown. Please refer to [link / reference]. Figure 1 The method includes the following steps:

[0095] Step 1: Obtain the shear wave velocity distribution of the earth-rock dam using the surface wave method, and convert the shear wave velocity into the first apparent resistivity based on the Faust formula.

[0096] It is understandable that the surface wave method is based on the dispersion characteristics of Rayleigh surface waves in layered media. The propagation speed of Rayleigh waves is closely related to the shear wave velocity of the medium, which in turn reflects the mechanical properties of the medium. When there are seepage channels in the dam body, the filling of these channels with water or weak materials will cause the wave velocity to decrease, thus forming a low-velocity anomaly.

[0097] In this embodiment, the method for converting the first apparent resistivity specifically includes:

[0098] Step 11: Use the damped least squares method to invert the surface wave detection data to obtain the shear wave velocity distribution of the earth-rock dam;

[0099] Step 12: Substitute the Garden formula into the standard Faust formula for conversion, establish the improved Faust formula, and use known geological data to fit the stratigraphic parameters in the improved Faust formula using the least squares method.

[0100] Step 13: Convert the shear wave velocity of the earth-rock dam into the first apparent resistivity based on the improved Faust formula, which is as follows:

[0101] ;

[0102] in, Indicates the first apparent resistivity. Indicates the shear wave velocity. Represents Poisson's ratio. , , This represents the stratigraphic parameters obtained by fitting known geological data.

[0103] The standard Faust formula is an empirical formula describing the relationship between formation resistivity and acoustic velocity, and its basic form is: ,in, Indicates the speed of sound. Indicates depth, Indicates true resistivity. , , This represents formation parameters. Considering that sonic logging data is typically lacking for earth-rock dams, this embodiment uses a modified Faust formula for conversion, i.e., the Garden formula. Substituting into the standard Faust formula, where Represent density, and use constants The value is set to 0.25, and the scaling factor is... The value is set to 0.31, and a transformation is performed. Further consideration is that the surface wave method yields the shear wave velocity. It needs to be converted into longitudinal wave velocity. : This leads to the improved Faust formula used for the first apparent resistivity conversion.

[0104] By obtaining the shear wave velocity of the dam body using the surface wave method and converting it into apparent resistivity, electrical parameters related to the mechanical properties of the medium are provided. This can effectively identify weak interlayers and loose areas in the dam body, providing an important basis for identifying seepage channels.

[0105] Step 2: Obtain the second apparent resistivity distribution of the earth-rock embankment using high-density electrical resistivity method.

[0106] It is understandable that high-density electrical resistivity tomography (EDT) is based on the difference in conductivity between soil and rock masses. In earth-rock dams, seepage channels are usually filled with soft materials with high water content, whose resistivity is significantly lower than that of the surrounding intact dam body, thus forming a low-resistivity anomaly.

[0107] In this embodiment, the method for obtaining the second apparent resistivity distribution includes:

[0108] Step 21: Construct the geoelectric model of the earth-rock dam, whose model parameter vector is... ,in, Indicates the first True resistivity of each grid cell Indicates the total number of model parameters. Indicates matrix transpose;

[0109] Step 22: Perform forward modeling on the geoelectric model to obtain theoretical apparent resistivity data. ,in, Indicates the index of the observation point. Indicates the total number of observation data;

[0110] Step 23: Calculate theoretical apparent resistivity data Compared with measured apparent resistivity data residuals: ;

[0111] Step 24: Construct the error function To measure the overall deviation between theoretical and measured data, this error function is defined as the sum of the squares of the residuals at all observation points:

[0112] ;

[0113] Step 25: Iteratively adjust the model parameter vector using the least squares method. Solve for the correction amount of the model parameter vector. and in accordance with Update the geoelectric model, in which, and They represent the first The model parameter vector and correction value corresponding to the next iteration;

[0114] Step 26: Repeat steps 22 to 25 until the error function values ​​of two adjacent iterations meet the convergence condition. At this point, the true resistivity distribution corresponding to the final model parameter vector is the second apparent resistivity distribution.

[0115] High-density electrical resistivity analysis directly obtains the apparent resistivity distribution, is sensitive to water-bearing areas, can effectively identify water-rich areas in leakage channels, and provides the most direct evidence of electrical anomalies.

[0116] Step 3: Use ground-penetrating radar to obtain the electromagnetic wave velocity distribution of the earth-rock dam, and convert the electromagnetic wave velocity into a third apparent resistivity based on the relationship between resistivity and relative permittivity.

[0117] It is understandable that ground-penetrating radar (GPR) detects electromagnetic waves based on their propagation characteristics in a medium. The propagation speed of electromagnetic waves in a medium is closely related to the dielectric constant, which in turn is related to the water content. Increased water content in a leakage channel leads to a decrease in the speed of electromagnetic waves.

[0118] In this embodiment, the method for converting the third apparent resistivity specifically includes:

[0119] Step 31: Based on the B-scan data obtained by ground-penetrating radar, identify the hyperbola formed by the reflected waves from underground targets. Please refer to [link to relevant documentation]. Figure 2 ;

[0120] Step 32: Extract the electromagnetic wave travel time information of the vertex position and each observation point position of the hyperbola using the hyperbola fitting method, and calculate the propagation speed of electromagnetic waves in the underground medium. The calculation formula is as follows:

[0121] ;

[0122] in, This indicates the propagation speed of electromagnetic waves in underground media. Indicates the antenna position. This indicates the horizontal position corresponding to the vertex of the hyperbola. Indicates position When the reflected wave from the target location travels, Represents the vertex of the hyperbola Travel time of the reflected wave from the target at the location;

[0123] Step 33: Based on the formula for the propagation speed of electromagnetic waves in a medium and the assumption that the relative permeability of soil and rock is equal to 1, the relationship between the electromagnetic wave speed and the relative permittivity of the medium is derived.

[0124] The formula for the propagation speed of electromagnetic waves in a medium is as follows:

[0125] ;

[0126] The relationship between the propagation speed of electromagnetic waves and the relative permittivity of the medium is as follows:

[0127] ;

[0128] in, This represents the speed of light in a vacuum. Represents the relative permittivity of the medium. Indicates the relative magnetic permeability of soil and rock mass;

[0129] Step 34: Using the empirical relationship between true resistivity and relative permittivity obtained by fitting experimental data of soil and rock, convert the relative permittivity into true resistivity according to the empirical relationship.

[0130] The empirical relationship between true resistivity and relative permittivity is as follows:

[0131] ;

[0132] in, Indicates true resistivity;

[0133] Step 35: Substitute the relationship between electromagnetic wave velocity and relative permittivity of the medium obtained in Step 33 into the empirical relationship between true resistivity and relative permittivity obtained in Step 34, and eliminate intermediate variables. The propagation speed of electromagnetic waves is directly converted into the final expression of the third apparent resistivity. The propagation speed of electromagnetic waves obtained in step 32 is converted into the third apparent resistivity through this final expression.

[0134] The propagation speed of electromagnetic waves can be directly converted into the final expression for the third apparent resistivity as follows:

[0135] ;

[0136] in, This represents the third apparent resistivity.

[0137] Ground penetrating radar (GPR) is characterized by high resolution, which can finely characterize the morphological features of shallow seepage channels, thus compensating for the shortcomings of other methods in shallow resolution.

[0138] Step 4: Perform data standardization processing on the first apparent resistivity, the second apparent resistivity and the third apparent resistivity respectively, including scale standardization and numerical standardization, to obtain standardized apparent resistivity data with consistent resolution and numerical range.

[0139] It is understandable that, since the three methods have different detection volumes, detection depths, and resolutions, it is necessary to unify the apparent resistivity data onto a grid of the same scale to achieve scale uniformity. In this embodiment, scale uniformity specifically includes:

[0140] Step 41: Resample the data of the first apparent resistivity, the second apparent resistivity and the third apparent resistivity to a grid with the same preset spatial resolution.

[0141] Step 42: For each grid cell, its scale-unified apparent resistivity is determined by a weighted average of the apparent resistivity of all sampling points within the grid cell and their weights.

[0142] ;

[0143] in, Indicates the first Apparent resistivity after unifying the scale of each grid cell Indicates the first cell within the grid. Apparent resistivity at each sampling point Indicates the first cell within the grid. The weight of each sampling point This indicates the number of sampling points within a grid cell.

[0144] In this embodiment, the numerical values ​​are uniformly transformed using a linear transformation method to unify the apparent resistivity of each grid cell to the range of [0, 100]. The formula is as follows:

[0145] ;

[0146] in, Indicates the first Standardized apparent resistivity after unifying the values ​​of each grid cell and These represent the maximum and minimum values ​​of apparent resistivity after all grid cells have been scaled uniformly.

[0147] By standardizing the scale, apparent resistivity data from different support volumes were normalized to the same geological grid unit. Numerical standardization eliminated the influence of dimensions, making different physical quantities comparable and maintaining the relative relationships of data distribution without affecting subsequent statistical feature extraction. This standardization process of scale and numerical standardization laid the foundation for subsequent multi-source data fusion.

[0148] Step 5: Use principal component analysis to perform dimensionality reduction and fusion on the standardized apparent resistivity data corresponding to the first, second, and third apparent resistivity to obtain the pseudo-apparent resistivity distribution.

[0149] Principal component analysis is a data dimensionality reduction method that transforms multiple correlated variables into a few uncorrelated variables (principal components) through linear transformation, while preserving as much information as possible from the original data.

[0150] In this embodiment, the specific process of principal component analysis includes:

[0151] Step 51: Use the standardized apparent resistivity data as the original feature vector. The feature vector of each grid cell is The first feature The first characteristic represents the normalized apparent resistivity corresponding to the first apparent resistivity, and the second characteristic... The third characteristic represents the normalized apparent resistivity corresponding to the second apparent resistivity. This represents the standardized apparent resistivity corresponding to the third apparent resistivity.

[0152] Step 53: Decentralize the original feature matrix using the following formula:

[0153] ;

[0154] in, Indicates the first The first grid cell of the first grid cell The value of each feature after decentralization Indicates the first The first grid cell of the first grid cell One characteristic, Indicates the number of grid cells;

[0155] Step 54: Based on the decentralized values, calculate the 3×3 covariance matrix. And solve for the covariance matrix. Characteristic equation ,in, Represent the identity matrix and obtain the eigenvalues ​​of the covariance matrix. and its corresponding eigenvectors ;

[0156] Step 54: Select the eigenvector corresponding to the largest eigenvalue to construct the projection space, and project the original eigenvector into this projection space to obtain the apparent resistivity. The projection formula is as follows:

[0157] ;

[0158] in, Indicates the first Apparent resistivity of each grid cell This represents the eigenvector corresponding to the largest eigenvalue. This indicates the matrix transpose.

[0159] Principal component analysis (PCA) essentially treats the three apparent resistivity values ​​of each grid cell as a point in three-dimensional space. By finding a new coordinate system, the variance of the projection of the data onto the first axis of this new coordinate system is maximized. This projected value is the pseudo-apparent resistivity. The pseudo-apparent resistivity is a fusion parameter obtained by fusing and refining the three different sources of apparent resistivity data. It represents the direction of maximum variance in the original three-dimensional data, i.e., the component that best reflects the common trend of the three apparent resistivity data. The pseudo-apparent resistivity obtained through dimensionality reduction and fusion using PCA retains the effective information from each method while highlighting common anomaly characteristics, thus achieving effective fusion of multi-source data. This solves the problem of non-uniqueness in geophysical exploration and significantly improves the ability to distinguish leakage channels.

[0160] In this embodiment, the principal component analysis method further includes:

[0161] Step 55: Calculate the contribution rate and cumulative contribution rate of the principal components corresponding to each feature, and evaluate the amount of original data information retained after dimensionality reduction based on the contribution rate and cumulative contribution rate.

[0162] The formula for calculating the contribution rate is as follows:

[0163] ;

[0164] in, Indicates the first Each feature corresponds to the contribution rate of the principal component. Indicates the first Each feature corresponds to an eigenvalue of a principal component. Indicates the summation index;

[0165] The formula for calculating the cumulative contribution rate is as follows:

[0166] ;

[0167] in, Indicates the preceding The cumulative contribution rate of each feature corresponds to the principal component.

[0168] It is understandable that eigenvalues The magnitude of the eigenvalue directly represents the amount of original data information carried by the corresponding principal component. The larger the eigenvalue, the greater the data fluctuation in the direction of the principal component, and the richer the information it contains. Based on this, this embodiment accurately calculates the information retention ratio through the contribution rate to ensure that information loss is within an acceptable range.

[0169] Step 6: Based on the apparent resistivity distribution, determine the location of potential leakage channels in the earth-rock embankment.

[0170] In this embodiment, determining the location of potential seepage channels in an earth-rock dam specifically includes:

[0171] Step 61: Based on the apparent resistivity distribution, identify areas with apparent resistivity below the abnormal threshold as low-resistivity abnormal areas, which correspond to areas with potential leakage risks due to high moisture content.

[0172] Step 62: Based on the continuous distribution of potential leakage areas, reconstruct the three-dimensional spatial morphology of the leakage channels using spatial interpolation methods to determine the direction, depth, and range of the leakage channels.

[0173] Step 63: Based on the structural characteristics of the dam, identify the inlet and outlet positions of the seepage channels in the reconstructed three-dimensional morphology of the seepage channels.

[0174] In practical applications, statistical methods can first be used to determine the anomaly threshold. If the apparent resistivity value is lower than the anomaly threshold, it indicates that the corresponding area is rich in moisture and is identified as a potential leakage area. Then, Kriging interpolation can be used to reconstruct the low-resistivity anomaly area in three dimensions to determine the direction, depth, and range of the leakage channel. Finally, the inlet and outlet locations of the leakage channel are identified in the reconstructed three-dimensional model. Through apparent resistivity anomaly analysis, the precise location and quantitative evaluation of leakage channels are achieved, providing a scientific basis for dam reinforcement and safety.

[0175] To verify the effectiveness and accuracy of the method proposed in this embodiment, a physical model test of an earth-rock dam with pre-set leakage risks was conducted.

[0176] To simulate the potential leakage risks of earth-rock dams under the influence of long-term complex natural environments, and based on the characteristics of natural earth-rock dam hazards, three pre-designed leakage channels were set up in the dam body, continuously and uniformly running through the upstream and downstream. Two of the pre-designed leakage channels were simulated by wrapping sand and gravel in gauze, while the natural channel was simulated by filling it with fine sand.

[0177] Please see Figures 3 to 6 The following are the designations for seepage channels: 1) **Leakage Channel 1:** 10cm diameter, arranged in a straight line below the upstream water level. The inlet and outlet sections are 3m from the right bank wall of the dam and 1.0m from the bottom of the dam. The gravel is wrapped with mesh. 2 **Leakage Channel 2:** 15cm diameter, arranged in an arc shape below the upstream water level. The inlet and outlet sections are 3m from the left bank wall of the dam and 1.5m from the bottom of the dam. The arc section bends 1m towards the wall. The gravel is wrapped with mesh. 3 **Leakage Channel 3 (Natural Seepage Channel):** 10cm diameter, arranged in a straight line below the upstream water level. The inlet and outlet sections are 5m from the left bank wall of the dam and 1.25m from the bottom of the dam. The channel is filled with fine sand.

[0178] Please see Figure 7The preset leakage channels are plotted on the fusion hazard location distribution map obtained according to the method of this embodiment. It can be seen that the fusion hazard locations obtained in this embodiment are quite close to the preset leakage channel locations, which can better reflect the leakage hazard information. To better highlight the fusion effect of this embodiment, please refer to... Figures 8 to 10 The coordinates of the fused hidden danger area and the preset hidden danger area obtained by the method of this embodiment are extracted and drawn. It can be seen that the method of this embodiment has high accuracy, with a maximum error of 0.2m and a minimum error of only 0.075m, which is highly accurate.

[0179] In summary, this embodiment organically combines surface wave method, high-density electrical resistivity tomography (EDS), and ground-penetrating radar (GPR) method to achieve multi-parameter joint constraints on elastic, electrical, and electromagnetic parameters. This effectively reduces the non-uniqueness of geophysical inversion and compensates for the multiple solutions and limitations inherent in the individual application of each geophysical method. The detection data source is a multi-physics field, and through scale and numerical unification processing, the problem of inconsistent data formats and resolutions among different methods is resolved, laying the foundation for data fusion. Principal component analysis is used for data fusion, which not only preserves the effective information of each method but also highlights common anomaly characteristics, meeting the requirements for data richness, accuracy, and reliability. Furthermore, this embodiment establishes a complete quantitative process from data acquisition to hazard identification, reducing human interference and improving the objectivity, accuracy, and reliability of seepage location, providing a scientific basis for the reinforcement and repair of earth-rock dams.

Claims

1. A method for locating seepage of earth-rock dam based on joint inversion and data fusion, characterized in that, The method comprises: Step 1, obtaining the shear wave velocity distribution of the earth-rock dam by using the surface wave method, and converting the shear wave velocity into the first apparent resistivity based on the Faust formula; Step 2, obtaining the second apparent resistivity distribution of the earth-rock dam by using the high-density electrical method; Step 3, obtaining the electromagnetic wave velocity distribution of the earth-rock dam by using the ground penetrating radar method, and converting the electromagnetic wave velocity into the third apparent resistivity based on the relationship between the resistivity and the relative dielectric constant; Step 4, performing data standardization processing on the first apparent resistivity, the second apparent resistivity and the third apparent resistivity, including scale unification and numerical unification, to obtain standardized apparent resistivity data with consistent resolution and unified numerical range; Step 5, performing dimension reduction fusion on the standardized apparent resistivity data corresponding to the first apparent resistivity, the second apparent resistivity and the third apparent resistivity by using the principal component analysis method, to obtain a pseudo apparent resistivity distribution; Step 6, determining the location of the leakage hidden danger channel in the earth-rock dam based on the pseudo apparent resistivity distribution.

2. The earth and rock dam leakage positioning method based on joint inversion and data fusion according to claim 1, characterized in that, In step 1, the conversion method of the first apparent resistivity specifically comprises: Step 11, obtaining the shear wave velocity distribution of the earth-rock dam by using the damped least squares method to invert the surface wave method detection data; Step 12, substituting the Garden formula into the standard Faust formula to perform conversion, establishing an improved Faust formula, and fitting the stratum parameters in the improved Faust formula by using the known geological data through the least squares method; Step 13, converting the shear wave velocity of the earth-rock dam into the first apparent resistivity based on the improved Faust formula, and the improved Faust formula is as follows: ； wherein, represents the first apparent resistivity, represents the shear wave velocity, represents the Poisson's ratio, , , represents the formation parameters fitted by known geological data.

3. The earth and rock dam leakage positioning method based on joint inversion and data fusion according to claim 1, characterized in that, In step 2, the obtaining method of the second apparent resistivity distribution comprises: Step 21, build the geoelectric model of the earth-rock dam, the model parameter vector of which is wherein, represents the true resistivity of the th grid unit, represents the total number of model parameters, represents the matrix transpose; Step 22, forward calculation is performed on the geoelectric model to obtain theoretical apparent resistivity data wherein, denotes the observation point index, denotes the total number of observation data; Step 23, calculate theoretical apparent resistivity data with the measured apparent resistivity data residuals: ; Step 24, construct error function The error function is defined as the sum of the squares of the residuals for all observation points, to measure the overall deviation of the theoretical data from the measured data: ； Step 25, iteratively adjust the model parameter vector using least square method , solve the correction of the model parameter vector , and update the geoelectric model according to , wherein, and respectively represent the corresponding model parameter vector and correction at the first iteration. Step 26, repeating steps 22 to 25 until the error function values of adjacent two iterations meet the convergence condition, and the true resistivity distribution corresponding to the final model parameter vector is the second apparent resistivity distribution.

4. The earth and rock dam leakage positioning method based on joint inversion and data fusion according to claim 1, characterized in that, In step 3, the conversion method of the third apparent resistivity specifically comprises: Step 31, identifying hyperbolas formed by underground target reflection waves based on the B-scan data obtained by the ground penetrating radar method; Step 32, extracting the vertex position of the hyperbola and the electromagnetic wave travel time information of each observation point position by using the hyperbolic fitting method, and calculating the propagation velocity of the electromagnetic wave in the underground medium, and the calculation formula is as follows: ； wherein, denotes the propagation speed of electromagnetic waves in the subsurface medium, denotes the antenna position, denotes the horizontal position of the hyperbolic vertex, denotes the position of the target reflected wave, denotes the propagation speed of electromagnetic waves in the subsurface medium, denotes the propagation speed of electromagnetic waves in the subsurface medium. Step 33, based on the propagation velocity formula of the electromagnetic wave in the medium and the assumption that the relative magnetic permeability of the rock-soil body is equal to 1, the relationship between the electromagnetic wave velocity and the relative dielectric constant of the medium is derived; The propagation velocity formula of the electromagnetic wave in the medium is as follows: ； The relationship between the electromagnetic wave velocity and the relative dielectric constant of the medium is as follows: ； wherein, denotes the speed of light in vacuum, denotes the relative permittivity of the medium, denotes the relative permeability of the rock-soil mass; Step 34, converting the relative dielectric constant into the true resistivity according to the empirical relationship between the true resistivity and the relative dielectric constant fitted by the rock-soil body experimental data; The empirical relationship between the true resistivity and the relative dielectric constant is as follows: ； wherein represents the true resistivity; Step 35, the electromagnetic wave velocity obtained in step 33 is substituted into the empirical relationship between true resistivity and relative dielectric constant in step 34 to eliminate the intermediate variable , to obtain a final expression for converting the propagation velocity of electromagnetic waves directly into a third apparent resistivity, and the propagation velocity of electromagnetic waves obtained in step 32 is converted into the third apparent resistivity through the final expression; The final expression of the electromagnetic wave velocity directly converted into the third apparent resistivity is as follows: ； wherein represents the third apparent resistivity.

5. The earth and rock dam leakage positioning method based on joint inversion and data fusion according to claim 1, characterized in that, In step 4, the scale unification adopts a linear sampling method, specifically comprising: Step 41, resample the data of the first apparent resistivity, the second apparent resistivity and the third apparent resistivity to a preset same spatial resolution grid respectively; Step 42, for each grid cell, the size of the apparent resistivity after uniform is determined by the weighted average of the apparent resistivity and its weight of all sampling points in the grid cell: ； wherein, represents the apparent resistivity of the grid cell, represents the apparent resistivity of the grid cell, represents the apparent resistivity of the sample point, represents the apparent resistivity of the sample point, represents the weight of the sample point, represents the weight of the sample point, represents the number of sample points in the grid cell.

6. The earth and rockfill dam leakage locating method based on joint inversion and data fusion according to claim 5, characterized in that, In step 4, the numerical uniform adopts a linear transformation method to uniformly the apparent resistivity of each grid cell after size uniform to the range of [0, 100], and the formula is as follows: ； wherein, represents the normalized apparent resistivity of the i-th grid cell, and respectively represent the maximum and minimum values of the apparent resistivity of all grid cells.​ 7. The earth and rock dam leakage positioning method based on joint inversion and data fusion according to claim 1, characterized in that, In step 5, the specific process of the principal component analysis method includes: Step 51: Use the standardized apparent resistivity data as the original feature vector. The feature vector of each grid cell is The first feature The first characteristic represents the normalized apparent resistivity corresponding to the first apparent resistivity, and the second characteristic... The third characteristic represents the normalized apparent resistivity corresponding to the second apparent resistivity. This represents the standardized apparent resistivity corresponding to the third apparent resistivity. Step 53, the original feature matrix is processed by decentralization, and the decentralization formula is as follows: ； wherein, represents the value of the decentering of the first feature of the first grid cell, represents the value of the decentering of the first feature of the first grid cell, represents the number of grid cells; Step 54, based on the values after decentralization, a 3*3 covariance matrix is calculated, and the characteristic equation of the covariance matrix is solved to obtain the eigenvalue and the corresponding eigenvector of the covariance matrix; Step 54, the eigenvector corresponding to the maximum eigenvalue is selected to form a projection space, and the original feature vector is projected into the projection space to obtain the pseudo apparent resistivity.

8. The earth and rock dam leakage positioning method based on joint inversion and data fusion according to claim 7, characterized in that, The projection formula of the original feature vector into the projection space to obtain the pseudo apparent resistivity is as follows: ； wherein, represents the pseudo resistivity of the th grid cell, represents the eigenvector corresponding to the largest eigenvalue, represents the matrix transpose.

9. The earth and rockfill dam leakage locating method based on joint inversion and data fusion according to claim 7, characterized in that, The principal component analysis method further includes: The contribution rate and the cumulative contribution rate of each feature corresponding to the principal component are calculated, and the amount of information of the original data retained after dimension reduction is evaluated according to the contribution rate and the cumulative contribution rate; The calculation formula of the contribution rate is as follows: ； wherein, denotes the contribution rate of the th feature to the principal component, denotes the eigenvalue of the th feature to the principal component, denotes the summation index; The calculation formula of the cumulative contribution rate is as follows: ； wherein indicates the cumulative contribution rate of the first characteristic to the principal component.

10. The joint inversion and data fusion based earth-rock dam leakage positioning method according to claim 1, characterized in that, The position of the leakage hidden danger channel in the earth-rock dam is determined, and specifically includes: Step 61, based on the pseudo apparent resistivity distribution, the area with pseudo apparent resistivity lower than the abnormal threshold is identified as a low-resistance abnormal area, and the low-resistance abnormal area corresponds to a water-rich leakage hidden danger area; Step 62, based on the continuous leakage hidden danger area distribution, the three-dimensional spatial form of the leakage channel is reconstructed by a spatial interpolation method to determine the trend, depth and range of the leakage channel; Step 63, in combination with the structural characteristics of the dam, the inlet position and the outlet position of the leakage channel are identified in the reconstructed three-dimensional form of the leakage channel.

Images