Low-density electrical reconstruction data volume mixed weight constraint inversion imaging method

Abstract

The application discloses a low-density electrical property reconstruction data body mixed weight constraint inversion imaging method, which comprises the following steps: deploying an electrical monitoring system at a downhole working face; establishing an inversion grid for an observation area; determining observation point weights; determining depth weights; fusing the observation point weights and the depth weights to form mixed weights; constructing a diagonal matrix as a mixed weight matrix; applying prior constraint information to a target function; linearizing and minimizing the target function to which the prior constraint is applied to obtain a model update formula to which mixed weight constraints are applied, which is an inversion problem to which mixed weight constraints are applied; and solving the inversion problem to which the mixed weight constraints are applied by using an inversion algorithm to realize three-dimensional inversion imaging of a low-density electrical property reconstruction data body in a coal mine downhole.

Description

Technical Field

[0001] This invention relates to an imaging method for mine geophysical exploration data, belonging to the field of mine geophysical exploration, and specifically to a hybrid weighted constraint inversion imaging method for low-density electrical reconstruction data volumes. Background Technology

[0002] The "Detailed Rules for Water Prevention and Control in Coal Mines" stipulate that mine hydrogeological exploration should employ a combination of geophysical methods, with electrical resistivity methods being the primary approach for hydrogeological exploration. As an important means of making hydrogeological information transparent, electrical resistivity methods detect and identify anomalous areas based on the electrical differences between abnormal geological bodies and surrounding rocks. They are particularly sensitive to low-resistivity anomalies caused by changes in rock water content. The DC resistivity method, in particular, allows for flexible adjustment of the electrode array layout based on different detection targets, and therefore has been widely used in the prevention and control of water hazards in coal mines.

[0003] The DC resistivity method, based on the principle of geometric sounding, detects electrical anomalies at different depths by changing the transmitter-receiver distance. In its application in underground coal mines, this method has led to the development of numerous observation methods, among which electrical profiling and electrical fluoroscopy are two of the most common. The electrical profiling method typically involves arranging several electrodes in a roadway on one side of the working face. Each electrode is powered sequentially; when one electrode transmits a signal, the others act as receiving electrodes, observing and recording the corresponding voltage data. This voltage data can be further converted into apparent resistivity through calculation. A pseudo-section of the apparent resistivity can be obtained through specific data arrangement. Two-dimensional inversion of the voltage data can then yield the resistivity profile below the measuring line. Electro-optic transillumination (EOT) is typically conducted in the transport and return airways of the working face. Several transmitting electrodes are placed in one roadway, and several receiving electrodes are placed in the opposite roadway. The transmitting electrodes are powered sequentially, and the voltage data penetrating the working face is observed and recorded using the receiving electrodes. The electrodes in the two roadways are then interchanged, and the transmission and reception are repeated. The received voltage data is then subjected to CT imaging or resistivity inversion imaging to obtain the resistivity distribution within the detection area. Both electrical profiling and EOT analyze and interpret anomalous geological structures within the detection area by delineating resistivity anomalies. For underground coal mine working faces, DC resistivity data acquisition using roadways, regardless of the observation method, can only utilize survey lines arranged in the roadways surrounding the working face. Compared to high-density electrical resistivity methods at the surface, the resulting electrical reconstruction data volume density is extremely low.

[0004] With the increasing demand for transparency in underground geological information in the intelligent construction of coal mines, resistivity inversion imaging has gradually become an important means of processing and interpreting DC resistivity data in mines. However, due to the limited observation space in coal mines, the volume density of electrical reconstruction data obtained from DC resistivity exploration is low, leading to the problem of multiple solutions in the inversion results, which is particularly prominent in the application of coal mines. Lu Jingjin (2016) pointed out in his paper "Three-dimensional resistivity inversion imaging technology for water-bearing / conducting structures in coal mines" that the low-resistivity area obtained by the dual-lane perspective observation system is strip-shaped, and the imaging resolution can be improved by establishing a multi-lane observation system. However, the multi-lane observation system is generally difficult to implement in engineering practice due to the limited observation space in coal mines. Li Maofei et al. (2022) pointed out in their paper "Mine DC electro-perspective floor detection and three-dimensional inversion interpretation" that the electro-perspective detection results can only give the horizontal position of the water hazard in the working face, and cannot accurately obtain the abnormal top and bottom interfaces. Because the transmission and reception distance of the electro-optical imaging method is limited by the width of the working face due to cross-tunnel observation, it is difficult to distinguish the depth changes of the anomaly. Therefore, the imaging results have poor resolution for the spatial location and distribution range of the anomaly in the depth direction.

[0005] The problem of ambiguity is an inherent challenge in resistivity inversion. Generally, based on model regularization constraints, prior information constraints are applied to minimize ambiguity. Prior information is typically based on the physical variation patterns of electrical parameters or derived from anomaly information obtained through geological inference or other geophysical exploration methods. To reduce the ambiguity of resistivity inversion using low-density electrical reconstruction data from underground coal mines, it is crucial to apply appropriate prior information constraints to the resistivity inversion process for these data volumes, thereby improving the imaging resolution of the inversion results regarding the spatial location and distribution range of anomalies. Summary of the Invention

[0006] This invention primarily addresses the technical problem of poor imaging resolution for the spatial location and distribution range of anomalies in inversion results of low-density electrical reconstruction data volumes in underground coal mines, a problem inherent in existing technologies. It provides a hybrid weighted constraint inversion imaging method for low-density electrical reconstruction data volumes. This method applies prior information constraints to the inversion process using hybrid weights based on observation point weights and depth weights. This effectively improves the imaging resolution along the working face dip and depth directions, suppresses banded anomalies in the imaging results, and enhances the ability of the imaging results to resolve deep anomalies.

[0007] The technical solution adopted by this invention to solve its technical problem is:

[0008] A hybrid weighted inversion imaging method for low-density electrical reconstruction data volumes includes the following steps:

[0009] Step 1: Deploy an electrical resistivity monitoring system at the downhole working face and collect monitoring data;

[0010] Step 2: Process the monitoring data in real time. Establish a rectangular coordinate system with the center point of the working face as the origin. The x-direction is consistent with the direction of the working face, the y-direction is consistent with the dip of the working face, and the z-direction is perpendicular to the working face and downwards. Establish an inversion grid for the observation area under this coordinate system.

[0011] The inversion mesh is established as follows: the number of meshes in the x-direction is N. x The number of grid cells in the y-direction is N. y The number of grid cells in the z-direction is N. z The total number of grid cells is N = N x ×N y ×N z Arbitrary grid cell numbers are denoted as (i, j, k), where 1, i ≤ N. x , 1, j≤N y 1, k≤N z The coordinates of the center point of any grid cell are defined as (x... i y j , z k ), 1, i≤N x , 1, j≤N y 1, k≤N z ;

[0012] Step 3: Determine the observation point weight w for each grid cell. c (i,j,k), i = 1, ...,N x j = 1, ..., N y k = 1, ..., N z ;

[0013] Step 4: Determine the depth weight w for each grid cell. z (i,j,k), i = 1, ...,N x j = 1, ..., N y k = 1, ..., N z ;

[0014] Step 5: Combine the observation point weights and depth weights of each grid cell to form a hybrid weight w. m (i,j,k);

[0015] Step 6, using the mixed weight w for each grid cell m (i,j,k) are the diagonal elements, and an N×N diagonal matrix is ​​constructed as a mixed weight matrix containing observation point information and depth information;

[0016] Step 7: Apply prior constraint information to the objective function of three-dimensional resistivity inversion using the hybrid weight matrix obtained in step 6 to obtain the objective function with prior constraints.

[0017] Step 8: Linearize and minimize the objective function with prior constraints obtained in Step 7 to obtain the model update formula with mixed weight constraints. This model update formula is the inversion problem with mixed weight constraints to be solved.

[0018] Step 9: Solve the inversion problem with mixed weight constraints using the inversion algorithm to reconstruct the resistivity value of each grid cell in the observation area, thereby realizing the three-dimensional inversion imaging of the low-density electrical reconstruction data volume in the coal mine.

[0019] Furthermore, in step 1, the electrical monitoring system includes a monitoring substation, monitoring cables, and monitoring electrodes, and the monitoring electrodes have both transmitting and receiving functions.

[0020] Furthermore, step 3 includes the following sub-steps:

[0021] Step 31, denote the observation point weight of any grid cell as w. c (i,j,k), initialize the observation point weights of all grid cells to 1, that is, let w c (i, j, k) = 1, i = 1, ..., N x j = 1, ..., N y k = 1, ..., N z The grid cell containing the observation point is numbered (i0,j0,k0), the grid cells adjacent to (i0,j0,k0) are numbered (i1,j1,k1), the grid cells adjacent to (i1,j1,k1) are numbered (i2,j2,k2), the grid cells adjacent to (i2,j2,k2) are numbered (i3,j3,k3), the grid cells adjacent to (i3,j3,k3) are numbered (i4,j4,k4), and so on.

[0022] Step 32, define the observation point weights: Let w c (i0,j0,k0)=C, where C is a constant greater than 1, then the observation point weights of the grid cells near the observation point are defined as follows: w c (i1,j1,k1)=α1w c (i0,j0,k0), w c (i2,j2,k2)=α2w c (i1,j1,k1), w c (i3,j3,k3)=α3w c (i2,j2,k2), w c(i4,j4,k4)=α4w c (i3,j3,k3), and so on; α1, α2, α3, α4, etc. are the attenuation coefficients of the observation point weights, which are constants greater than 0 and less than 1;

[0023] Step 33: Starting from the grid cell numbered (1, 1, 1), the observation point weights of the grid cells described in step 31 are assigned values ​​according to the weight definition described in step 32.

[0024] Furthermore, step 33 includes the following operations:

[0025] The weight value w to be assigned to the current grid cell is calculated according to the weight definition described in step 32. It is then determined whether the weight value w to be assigned to the current grid cell is greater than the initial weight value w of the observation point of that grid cell. c (i,j,k) is greater than w: if w>w c (i,j,k), let w c (i,j,k)=w; if w <w c (i,j,k), w c The values ​​of (i,j,k) remain unchanged from their initial weight values.

[0026] Furthermore, in step 4, the depth weight takes the following form: w z (i,j,k)=(z k / z1) -a i = 1, ..., N x j = 1, ..., N y k = 1, ..., N z ; a is a constant greater than 0.

[0027] Furthermore, in step 5, the mixed weights take the following form: w m (i,j,k)=w c (i,j,k)w z (i,j,k), i = 1, ...,N x j = 1, ..., N y k = 1, ..., N z .

[0028] Furthermore, in step 6, the diagonal matrix is ​​as follows:

[0029] W m (n,n)=w m (i,j,k),

[0030] In the formula, 1≤n≤N, 1≤i≤N x , 1≤j≤N y , 1≤k≤N zn = (k-1)N x N y +(j-1)N x +i.

[0031] Furthermore, in step 7, the objective function with applied prior constraints is as follows:

[0032]

[0033]

[0034] In the formula, Δd is the residual between the measured data and the simulated observation data d; ρ is the resistivity model generated by the inversion fitting; ρ ref It is a reference model; W d W is the data weight matrix. ρ W is the model weight matrix; β is the regularization parameter. m I is the mixed weight matrix; G is the identity matrix; x G y G z The gradient operators are defined in the x, y, and z directions, respectively; λ is the model matrix W. ρ The weights of the diagonal elements; μ, γ, and η are the gradient weights in the x, y, and z directions, respectively.

[0035] Compared with the prior art, the beneficial effects of the present invention are as follows:

[0036] This invention makes the imaging results more focused by applying observation point weights, effectively improves the resolution of the imaging results in the depth direction by applying depth weights, and simultaneously improves the imaging resolution of the inversion results along the working face dip and depth direction by applying a hybrid weight based on observation point weights and depth weights, thereby enhancing the imaging results' ability to distinguish deep anomalies. The hybrid weight constraint method of this invention is applicable to the inversion imaging of DC resistivity detection data obtained in underground coal mines using different observation devices and systems, and is also applicable to the inversion imaging of DC resistivity detection data obtained on the surface or under combined well-to-surface conditions using different observation devices and systems. Attached Figure Description

[0037] Figure 1 Convergence curves of the inversion algorithm when different weights are applied for h=0m;

[0038] Figure 2 The inversion imaging results when h = 0m without applying prior information constraints;

[0039] Figure 3 The inversion imaging result is when observation point weights are applied at h = 0m;

[0040] Figure 4The inversion imaging result is when depth weights are applied at h = 0m;

[0041] Figure 5 The inversion imaging result when mixed weights are applied at h=0m;

[0042] Figure 6 The inversion imaging results when h = 20m without applying prior information constraints;

[0043] Figure 7 The inversion imaging results are obtained when observation point weights are applied at h = 20m.

[0044] Figure 8 The inversion imaging results are when depth weighting is applied at h = 20m.

[0045] Figure 9 The inversion imaging results when a mixed weight is applied at h=20m. Detailed Implementation

[0046] The specific embodiments of the present invention will be further described below with reference to the accompanying drawings.

[0047] This invention provides a hybrid weighted inversion imaging method for low-density electrical reconstruction data volumes, which specifically includes the following steps:

[0048] Step 1: Design a monitoring scheme based on the actual working conditions of the working face, ensuring effective signal strength. Select an electrical monitoring system suitable for the working face. Deploy the electrical monitoring system on the downhole working face according to the monitoring scheme, and collect and process the monitoring data in real time. The electrical monitoring system includes monitoring substations, monitoring cables, and monitoring electrodes, among which the monitoring electrodes have both transmitting and receiving functions.

[0049] Step 2: When processing the monitoring data in real time, establish a rectangular coordinate system with the center point of the working face as the origin. The x-direction is consistent with the direction of the working face, the y-direction is consistent with the dip of the working face, and the z-direction is perpendicular to the working face and downwards. Invert the grid for the observation area under this coordinate system.

[0050] Specifically, the inversion mesh is established as follows: the number of meshes in the x-direction is N. x The number of grid cells in the y-direction is N. y The number of grid cells in the z-direction is N. z The total number of grid cells is N = N x ×N y ×N z Arbitrary grid cell numbers are denoted as (i, j, k) (i ≤ N). x 1 means j≤N y 1 means k≤N zThe coordinates of the center point of any grid cell are defined as (x...). i y j , z k (1 heart i≤N) x 1 heart j≤N y 1, k≤N z ).

[0051] Step 3: Determine the weight of the observation point for each grid cell.

[0052] Specifically, the reliability of the resistivity information of each grid cell obtained through inversion reconstruction gradually decreases with increasing distance from the grid cell to the observation point. Based on this characteristic, observation point information constraints are applied to form observation point weights. Observation point weights are generally constants greater than 1, and the observation point weight values ​​of different grid cells gradually decrease with increasing distance from the grid cell to the observation point. Based on this, this step specifically includes the following sub-steps:

[0053] Step 31, denote the observation point weight of any grid cell as w. c (i,j,k), the observation point weights are defined through the following steps:

[0054] The initial weights of all observation points in the grid cells are set to 1, i.e., w is set to 1. c (i,j,k)=1(i=1,...,N x j = 1, ..., N y k = 1, ..., N z The grid cell containing the observation point is denoted as (i0,j0,k0); the grid cells adjacent to (i0,j0,k0) are denoted as (i1,j1,k1) (excluding the grid cells marked as (i0,j0,k0); the grid cells adjacent to (i1,j1,k1) are denoted as (i2,j2,k2) (excluding the grid cells marked as (i0,j0,k0) and (i1,j1,k1); and the grid cells adjacent to (i2,j2,k1) are denoted as (i2,j2,k2). The grid cells adjacent to (i3,j3,k3) are denoted as (i3,j3,k3) (excluding grid cells labeled (i0,j0,k0), (i1,j1,k1), (i2,j2,k2)). The grid cells adjacent to (i3,j3,k3) are denoted as (i4,j4,k4) (excluding grid cells labeled (i0,j0,k0), (i1,j1,k1), (i2,j2,k2), (i3,j3,k3)), and so on.

[0055] Step 32, define the weights: Let w c (i0,j0,k0)=C (C is a constant greater than 1), then the observation point weights of the grid cells near the observation point are defined as follows: w c (i1,j1,k1)=α1wc (i0,j0,k0), w c (i2,j2,k2)=α2w c (i1,j1,k1), w c (i3,j3,k3)=α3w c (i2,j2,k2), w c (i4,j4,k4)=α4w c (i3,j3,k3), and so on; α1, α2, α3, α4, etc. are the attenuation coefficients of the observation point weights, which are generally constants greater than 0 and less than 1;

[0056] Step 33: Starting from the grid cell numbered (1, 1, 1), the observation point weights of the grid cells (i0, j0, k0), (i1, j1, k1), (i2, j2, k2), (i3, j3, k3), and (i4, j4, k4) mentioned in Step 31 are reassigned according to the weight definition described in Step 32. Specifically, the weight value w to be assigned to the current grid cell is calculated according to the weight definition described in Step 32, and it is determined whether the weight value w to be assigned to the current grid cell is greater than the initial value w of the observation point weight of the grid cell. c (i,j,k) is greater than w: if w>w c (i,j,k), let w c (i,j,k)=w; if w <w c (i,j,k), w c The values ​​of (i,j,k) remain unchanged from their initial weight values. Through the above process, it can be ensured that the weight of the observation point in each grid cell satisfies the law that the weight value of different grid cells gradually decreases with increasing distance from the grid cell to the observation point, while also satisfying w c (i,j,k)≥1(1≤i≤N x , 1≤j≤N y 1k≤N z ).

[0057] Step 4: Determine the depth weight for each grid cell.

[0058] Specifically, for the DC resistivity method, the electric field intensity exhibits a characteristic of decreasing with depth, and the resolution for anomalies also gradually decreases with depth. Based on this characteristic, depth information constraints are applied to form depth weights. The depth weights are related to the depth z of the current mesh cell. k (1≤k≤N z Related to this, the depth weight of different grid cells gradually decreases as the depth of the grid cell increases. Let w be the depth weight of any grid cell. z (i,j,k), in the following form:

[0059] w z (i,j,k)=(z k / z1) -a (i = 1, ..., N) x j = 1, ..., N y k = 1, ..., N z )

[0060] In the formula, a is a constant greater than 0.

[0061] Step 5: The observation point weights and depth weights of each grid cell are fused to form a hybrid weight. This step enables the simultaneous application of observation point information constraints and depth information constraints to each grid cell.

[0062] Specifically, the mixed weight of any grid cell is denoted as w. m (i,j,k) is a mixed weight formed by the product of depth weight and observation point weight, and its form is as follows:

[0063] w m (i,j,k)=w c (i,j,k)w z (i,j,k)(i=1,……,N x j = 1, ..., N y k = 1, ..., N z ).

[0064] Step 6: Construct a diagonal matrix using the mixed weights of each grid cell as diagonal elements, and use it as a mixed weight matrix containing observation point information and depth information.

[0065] Specifically, with w m (i,j,k)(i=1,……,N x j = 1, ..., N y k = 1, ..., N z Using the diagonal elements as the basis, construct an N×N diagonal matrix W. m :

[0066] W m (n,n)=w m (i,j,k)(1≤n≤N,1≤i≤N x , 1≤j≤N y , 1≤k≤N z )

[0067] In the formula, n = (k-1)N x N y +(j-1)N x +i, W m It is a mixed weight matrix.

[0068] Step 7: Apply prior constraint information to the objective function of three-dimensional resistivity inversion using the hybrid weight matrix obtained in step 6 to obtain the objective function with prior constraints.

[0069] Specifically, for the three-dimensional resistivity regularization inversion objective function of the following form:

[0070]

[0071] In the formula: Δd is the residual between the measured data and the simulated observation data d; ρ is the resistivity model generated by the inversion fitting; ρ ref It is a reference model; W d W is the data weight matrix. ρ W is the model weight matrix; β is the regularization parameter, which is expressed in this invention through the model weight matrix W. ρ Imposing prior information constraints on the objective function:

[0072]

[0073] In the formula: W m The weighted matrix is ​​a mixture of observation point information and depth information; I is the identity matrix; G is the weighted matrix. x G y G z The gradient operators are defined in the x, y, and z directions, respectively; λ is the model matrix W. ρ The weights of the diagonal elements; μ, γ, and η are the gradient weights in the x, y, and z directions, respectively. (Gradient weights can be user-defined.)

[0074] Step 8: Linearize and minimize the objective function obtained in Step 7 to obtain the model update formula with mixed weight constraints. This model update formula is the inversion problem to be solved with mixed weight constraints.

[0075] Specifically, the objective function for the three-dimensional resistivity regularization inversion described in step 7 is linearized and minimized to obtain the following model update formula with mixed weight constraints:

[0076]

[0077] Where α is the linear search parameter, This is the sensitivity matrix.

[0078] Step 9: Solve the inversion problem with mixed weight constraints using the inversion algorithm to reconstruct the resistivity value of each grid cell in the observation area, thereby realizing the three-dimensional inversion imaging of the low-density electrical reconstruction data volume in the coal mine.

[0079] The beneficial effects of this invention are that applying observation point weights can make the imaging results more focused, applying depth weights can effectively improve the resolution of the imaging results in the depth direction, and applying a hybrid weight based on observation point weights and depth weights can simultaneously improve the imaging resolution of the inversion results along the working face dip and the depth direction, thereby enhancing the imaging results' ability to distinguish deep anomalies. The hybrid weight constraint method is applicable to the inversion imaging of DC resistivity detection data obtained by different observation devices and systems in underground coal mines, and is also applicable to the inversion imaging of DC resistivity detection data obtained by different observation devices and systems under surface or well-to-surface combined conditions.

[0080] To verify the feasibility and effectiveness of this invention, the following experiment was designed for verification:

[0081] Each of the following conditions is not subject to prior information constraints (w) m =1) Apply observation point weights (w) m =w c Apply depth weights (w) m =w z ) and applying mixed weights (w m =w c w z Four weighting methods were used for calculation. The working face has a strike length of 200m and a dip width of 100m. Electrodes were placed in the working face transport roadway and return airway, with an adjacent electrode spacing of 10m. The survey line length was 200m, and 21 electrodes were placed in each roadway. Electrofluoroscopic data were collected using a dipole-dipole observation device. A low-resistivity spherical anomaly exists below the working face floor, with a radius of 30m and a resistivity of 20Ω·m, while the surrounding rock resistivity is 500Ω·m. The sphere burial depth h was compared with 0m and 20m. When analyzing the imaging results, a resistivity anomaly was identified when the change in resistivity value relative to the background resistivity (500Ω·m) reached 10% or more. The threshold for low-resistivity anomalies was 450Ω·m.

[0082] Figure 1The convergence curves of the inversion algorithm with different weights applied when h=0m are shown. 1 represents the convergence curve without prior information constraints, 2 represents the convergence curve with observation point weights, 3 represents the convergence curve with depth weights, and 4 represents the convergence curve with mixed weights. Without prior information constraints and with observation point weights, the convergence curve becomes almost flat after about 5 iterations, and the relative residual cannot be further reduced after more than 10 iterations, leading to termination of the inversion. With depth weights and mixed weights, the convergence curve decreases faster. After more than 10 iterations, the convergence curve approaches a straight line parallel to the horizontal axis, indicating that the decay of the data residual is extremely slow, and further iterations cannot further improve the inversion results. The comparison shows that the algorithm converges better with depth weights and mixed weights, with the algorithm with mixed weights showing slightly better convergence than the algorithm with depth weights. To avoid invalid iterative calculations, an inversion stopping criterion is established based on the rate of decrease of the relative residual. When χ²... i / χ i-1 ≥0.99(χ i Let χ be the relative residual after the i-th inversion iteration. i-1 Let be the relative residual before the i-th inversion iteration (i≥1). If i ≥ 1, it is considered that further iterations cannot improve the inversion result, and the inversion is stopped at this point. The model examples below all perform inversion calculations based on this criterion.

[0083] When h = 0m, the spatial distribution range of the sphere is -30m, -30m, 0m, and 60m. The depth of the sphere's center is z = 30m. Figures 2-5 The inversion results of electrofluoroscopic data at h=0m under different model weights are given, and the extreme point of the low-resistivity region in the figure is defined as the center of the low-resistivity region. Figure 2 The image shows the imaging results when h = 0m without applying prior information constraints. As can be seen from the figure, the spatial distribution range of the low-resistivity region is as shown in the figure for -20m, -70m, and 0m. The center depth of the low-resistivity region is z = 0m. Both the spatial distribution range and depth deviate significantly from the model. Figure 3 The image shows the imaging results when the observation point weights are applied at h=0m. As can be seen from the figure, the spatial distribution range of the low-resistivity region (the result when the weights are -20m, -40m, and 0m) ​​is close to the model range, but the center depth of the low-resistivity region z=0m deviates significantly from the model depth. Figure 4This is the imaging result when depth weight is applied at h=0m. As can be seen from the figure, the distribution range of the low-resistivity region along the working surface (-20m≤x≤20m) is close to the model range, while the distribution range along the working surface (-60m≤y≤60m) is larger than the model range. The distribution range in the depth direction (0m≤z≤60m) basically matches the model range. The center depth of the low-resistivity region z=30m, which basically matches the model depth. Figure 5 The image shows the imaging results when h=0m and mixed weights are applied. As can be seen from the figure, the distribution range of the low-resistivity region along the working surface (image with a weight of -30m) and in the depth direction (0m depth) basically matches the model range. The distribution range along the working surface dip (basically matches the -50m range) is larger than the model range. The center depth of the low-resistivity region, z=30m, basically matches the model depth. In comparison, for the spherical model with h=0m, the imaging results deviate significantly from the model without prior information constraints. Applying observation point weights effectively improves the resolution along the working surface dip, but the improvement in depth direction resolution is poor. Applying depth weights effectively improves the resolution in the depth direction, but the resolution along the working surface dip is basically unimproved. Applying mixed weights simultaneously improves the resolution along both the working surface dip and depth directions. The strip-shaped anomalies along the working surface dip are improved to some extent. The improvement in depth direction resolution is slightly better than when applying depth weights, but the improvement in working surface dip resolution is slightly worse than when applying observation point weights.

[0084] When h = 20m, the spatial distribution range of the sphere is -30m ≤ x ≤ 30m, -30m ≤ y ≤ 30m, 20m ≤ z ≤ 80m, and the depth of the sphere's center is z = 50m. Figures 6-9 The inversion results of electrofluoroscopy data with different model weights are given when h=20m. Figure 6 The image shows the imaging results when h=20m without prior information constraints. As can be seen from the figure, the low-resistivity anomaly intensity is weak, the low-resistivity region is not focused, and obvious false low-resistivity anomalies appear below the survey line. The distribution range of the low-resistivity region in the depth direction (0m plus prior information) deviates significantly from the model range. Figure 7 The image shows the imaging results when the observation point weights are applied at h=20m. As can be seen from the figure, the low-resistivity anomaly intensity is relatively weak. The distribution range of the low-resistivity region along the working face (with the weight at -20m) and along the working face dip (with the working face dip at -40m) are close to the model range. The distribution range in the depth direction (closer to the model range at 0m) is smaller than the model range. The center depth of the low-resistivity region is located at z=20m, which is shallower than the model depth. Figure 8The image shows the imaging results when depth weighting is applied at h=20m. As can be seen from the figure, the low-resistivity anomaly is significantly enhanced. The distribution range of the low-resistivity region along the working face (-20m≤x≤20m) is close to the model range, while the distribution range along the working face dip (-70m≤y≤70m) is larger than the model range. The distribution range in the depth direction (10m≤z≤50m) is smaller than the model range. The center depth of the low-resistivity region is located at z=30m, which is shallower than the model depth. Figure 9 The image shows the imaging results when mixed weights are applied at h=20m. As can be seen from the figure, the low-resistivity anomaly is significantly enhanced. The distribution range of the low-resistivity region along the working face (-30m≤x≤30m) is basically consistent with the model range. The distribution range along the working face dip (-50m≤y≤50m) is larger than the model range. The distribution range in the depth direction (10m≤z≤70m) is close to the model range. The center depth of the low-resistivity region is located at z=40m, which is close to the model depth. The comparison shows that for the spherical model with h = 20m, the imaging results deviate significantly from the model without prior information constraints. Applying observation point weights improves the resolution of the imaging results along the working surface, but the improvement in depth resolution is poor. Applying depth weights improves the resolution of the imaging results in the depth direction, but does not improve the resolution along the working surface. Applying mixed weights significantly improves the resolution of the imaging results in the depth direction, with a better improvement effect than applying depth weights alone. The strip-shaped anomalies along the working surface are improved to some extent, but low-resistivity false anomalies appear below the survey line. In addition, low-resistivity anomalies are significantly enhanced when depth weights and mixed weights are applied, with the intensity of low-resistivity anomalies being the greatest when mixed weights are applied, which is beneficial for improving the detection capability of deep low-resistivity anomaly structures.

[0085] A comprehensive comparison shows that, for electro-optical imaging data obtained using a dipole-dipole observation device, applying depth weights can effectively improve the resolution of the imaging results in the depth direction; applying observation point weights can make the distribution range of the imaging results more focused; applying mixed weights can significantly improve the imaging results' ability to distinguish deep anomalies, and the strip-shaped anomalies along the working surface are improved to a certain extent; when the burial depth of the sphere changes, both the imaging results obtained by applying depth weights and those obtained by applying mixed weights can identify the depth changes of the sphere, among which the imaging results obtained by applying mixed weights reflect the depth changes more closely to the actual changes in the model.

[0086] The low-density electrical reconstruction data volume hybrid weight constraint inversion imaging method of the present invention is applicable to the inversion imaging of DC resistivity detection data obtained in underground coal mines using different observation devices and systems, and is also applicable to the inversion imaging of DC resistivity detection data obtained on the ground or under combined well-ground conditions using different observation devices and systems.

Claims

1. A hybrid weighted inversion imaging method for low-density electrical reconstruction data volumes, characterized in that, Specifically, the steps include the following: Step 1: Deploy an electrical resistivity monitoring system at the downhole working face and collect monitoring data; Step 2: Process the monitoring data in real time. Establish a rectangular coordinate system with the center point of the working face as the origin. The x-direction is consistent with the direction of the working face, the y-direction is consistent with the dip of the working face, and the z-direction is perpendicular to the working face and downwards. Establish an inversion grid for the observation area under this coordinate system. The inversion mesh is established as follows: the number of meshes in the x-direction is N. x The number of grid cells in the y-direction is N. y The number of grid cells in the z-direction is N. z The total number of grid cells is N = N x ×N y ×N z Arbitrary grid cell numbers are denoted as (i, j, k), 1 ≤ i ≤ N. x , 1≤j≤N y , 1≤k≤N z The coordinates of the center point of any grid cell are defined as (x... i y j , z k ), 1≤i≤N x , 1≤j≤N y , 1≤k≤N z ; Step 3: Determine the observation point weight w for each grid cell. c (i,j,k), i = 1, ...,N x j = 1, ..., N y k = 1, ..., N z ; Step 4: Determine the depth weight w for each grid cell. z (i,j,k), i = 1, ...,N x j = 1, ..., N y k = 1, ..., N z ; Step 5: Combine the observation point weights and depth weights of each grid cell to form a hybrid weight w. m (i,j,k); Step 6, using the mixed weight w for each grid cell m (i,j,k) are the diagonal elements, and an N×N diagonal matrix is ​​constructed as a mixed weight matrix containing observation point information and depth information; Step 7: Apply prior constraint information to the objective function of three-dimensional resistivity inversion using the hybrid weight matrix obtained in step 6 to obtain the objective function with prior constraints. Step 8: Linearize and minimize the objective function with prior constraints obtained in Step 7 to obtain the model update formula with mixed weight constraints. This model update formula is the inversion problem with mixed weight constraints to be solved. Step 9: Solve the inversion problem with mixed weight constraints using the inversion algorithm to reconstruct the resistivity value of each grid cell in the observation area, thereby realizing the three-dimensional inversion imaging of the low-density electrical reconstruction data volume in the coal mine.

2. The low-density electrical reconstruction data volume hybrid weighted constrained inversion imaging method as described in claim 1, characterized in that, In step 1, the electrical monitoring system includes a monitoring substation, monitoring cables, and monitoring electrodes, and the monitoring electrodes have both transmitting and receiving functions.

3. The low-density electrical reconstruction data volume hybrid weighted constrained inversion imaging method as described in claim 1, characterized in that, Step 3 includes the following sub-steps: Step 31, denote the observation point weight of any grid cell as w. c (i,j,k), initialize the observation point weights of all grid cells to 1, that is, let w c (i, j, k) = 1, i = 1, ..., N x j = 1, ..., N y k = 1, ..., N z The grid cell containing the observation point is numbered (i0,j0,k0), the grid cells adjacent to (i0,j0,k0) are numbered (i1,j1,k1), the grid cells adjacent to (i1,j1,k1) are numbered (i2,j2,k2), the grid cells adjacent to (i2,j2,k2) are numbered (i3,j3,k3), the grid cells adjacent to (i3,j3,k3) are numbered (i4,j4,k4), and so on. Step 32, define the observation point weights: Let w c (i0,j0,k0)=C, where C is a constant greater than 1, then the observation point weights of the grid cells near the observation point are defined as follows: w c (i1,j1,k1)=α1w c (i0,j0,k0), w c (i2,j2,k2)=α2w c (i1,j1,k1), w c (i3,j3,k3)=α3w c (i2,j2,k2), w c (i4,j4,k4)=α4w c (i3,j3,k3), and so on; α1, α2, α3, α4, etc. are the attenuation coefficients of the observation point weights, which are constants greater than 0 and less than 1; Step 33: Starting from the grid cell numbered (1, 1, 1), the observation point weights of the grid cells described in step 31 are assigned values ​​according to the weight definition described in step 32.

4. The low-density electrical reconstruction data volume hybrid weighted constrained inversion imaging method as described in claim 3, characterized in that, Step 33 includes the following operations: The weight value w to be assigned to the current grid cell is calculated according to the weight definition described in step 32. It is then determined whether the weight value w to be assigned to the current grid cell is greater than the initial weight value w of the observation point of that grid cell. c (i,j,k) is greater than w: if w>w c (i,j,k), let w c (i,j,k)=w; if w <w c (i,j,k), w c The values ​​of (i,j,k) remain unchanged from their initial weight values.

5. The low-density electrical reconstruction data volume hybrid weighted constrained inversion imaging method as described in claim 1, characterized in that, In step 4, the depth weight is in the following form: w z (i,j,k)=(z k / z1) -a i = 1, ..., N x j = 1, ..., N y k = 1, ..., N z ; a is a constant greater than 0.

6. The low-density electrical reconstruction data volume hybrid weighted constrained inversion imaging method as described in claim 1, characterized in that, In step 5, the hybrid weighting takes the following form: w m (i,j,k)=w c (i,j,k)w z (i,j,k), i = 1, ...,N x j = 1, ..., N y k = 1, ..., N z .

7. The low-density electrical reconstruction data volume hybrid weighted constrained inversion imaging method as described in claim 1, characterized in that, In step 6, the diagonal matrix is ​​as follows: W m (n,n)=w m (i,j,k), In the formula, 1≤n≤N, 1≤i≤N x , 1≤j≤N y , 1≤k≤N z n = (k-1)N x N y +(j-1)N x +i.

8. The low-density electrical reconstruction data volume hybrid weighted constrained inversion imaging method as described in claim 1, characterized in that, In step 7, the objective function with applied prior constraints is as follows: In the formula, Δd is the residual between the measured data and the simulated observation data d; ρ is the resistivity model generated by the inversion fitting; ρ ref It is a reference model; W d W is the data weight matrix. ρ It is the model weight matrix; β is the regularization parameter: W m I is the mixed weight matrix; G is the identity matrix; x G y G z The gradient operators are defined in the x, y, and z directions, respectively; λ is the model matrix W. ρ The weights of the diagonal elements; μ, γ, and η are the gradient weights in the x, y, and z directions, respectively.

Images