# 3D Stratum Interpolation Algorithm of Metro Tunnel Based on BIM and Data Fusion

Authors: Haidong Jiang; Jian Yin; Yonggui Chen; Guofa Zhang; Xiaoyong Peng
Journal: Security and Communication Networks
 (2022)

---

## Abstract

The three-dimensional stratum space interpolation algorithm plays a vital role in promoting the safety, management, and decision-making of subway tunnels and is a current research hotspot. The purpose of full text is to study the 3D stratum space interpolation algorithm of subway tunnels and to study the 3D stratum space interpolation algorithm based on building information modeling (BIM) and data fusion technology. This article first uses computer three-dimensional visualization technology to simulate the surface of the ground-fall tunnel based on known data. It is to do the interpolation processing on the two-dimensional space of the subway tunnel to create an accurate surface model of the subway tunnel. It then uses the block model method to divide the internal data of the tunnel according to actual needs. In order to create a model of the internal data of the tunnel, three-dimensional spatial interpolation is performed. After the three-dimensional interpolation is completed, it is convenient to quickly estimate the taste and reserves according to the interpolation results. The experiments in the text show that the average error value of the algorithm in this paper is 38.98 and the performance is improved by 23.4%. The algorithm in the text has obvious advantages in interpolation effect, and the interpolation efficiency is faster. A bird’s eye view of the current study is shown in Figure 1.


---

## Body

## 1. Introduction

At the beginning of geoscience research, most of the data interpolation algorithms were for two-dimensional, but relatively few for three-dimensional space research. But with the advent of the information age, the development of computer three-dimensional visualization technology has become more and more in-depth in our lives and has become an indispensable part of us. With its increasingly powerful development and increasingly perfect functions, interpolation of three-dimensional spatial data has become more and more realistic. With its unparalleled advantages and unstoppable attraction, it has injected incomparably fresh blood into geoscience researchers and opened up a new sky for geoscience research. It frees people from the predicament and shackles of manual handling of taste and reserve forecasts in the past and leads people into a whole new world [1, 2]. Under the guidance of 3D visualization technology, not only can people quickly complete the estimation of taste and reserves, but the precision and accuracy are also guaranteed.Through the computer three-dimensional visualization technology to reproduce the real three-dimensional structure of the geological body, people can reproduce the three-dimensional real structure of the underground tunnel. It can also be used to analyze and judge the spatial distribution of the tunnel and its spatial structure. Essentially, spatial interpolation algorithms are all methods of moving weighted average based on linear theory. Therefore, it is easy to cause the loss of information caused by smoothing the data locally, and there is a problem of inefficient calculation of massive three-dimensional data. The geological data itself is three-dimensional and dynamic. When a three-dimensional voxel model is used to express a geological body, the amount of data will increase significantly, and the amount of calculation will also increase significantly. In addition, its description of the spatial distribution of tunnel grades with nonlinear characteristics also has certain limitations in accuracy, and it is necessary to improve the accuracy and efficiency of its algorithm. Therefore, it is particularly urgent and important to seek a method that has powerful processing capabilities for nonlinear systems and does not require sampling data distribution and assumptions. It is of great practical significance to use the improved interpolation algorithm to verify and apply the actual exploration data, and to combine the three-dimensional visualization technology with the interpolation algorithm to visually display the processing results of the interpolation algorithm. Experts can analyze its internal structure, causes, and laws based on the spatial characteristics displayed by the three-dimensional visualization, which is of great significance for further in-depth study of the properties of the tunnel.The innovation of full text is based on BIM and data fusion technology to study the interpolation algorithm of subway tunnel three-dimensional stratum space. Through BIM and data fusion technology, research and build a three-dimensional geological model, follow the design principles and design patterns, analyze the data structure of the design model, the functional structure of the design model, and the tasks to be performed by each module, and describe how they collaborate work, and then optimize the spatial interpolation algorithm based on the actual needs of modeling, and propose a modeling method based on feature data interpolation. This method first divides the geological blocks near the surface and performs spatial interpolation on the internally completed feature data. Modeling avoids the difficulty of constructing models with layers close to the surface by drawing triangular patches in complex areas. This realizes the local variable range of the variogram and the filling of the neighborhood point set. It can be applied to the interpolation of the three-dimensional attribute data volume, and the spatial position distribution relationship and structural characteristics can be obtained. Through comparative analysis with the traditional inverse distance and the kriging algorithm, the actual performance of the algorithm in the text is obtained (Figure1).Figure 1
Aerial view of the study.


## 2. Related Work

At present, most parallel spatial interpolation algorithms only use one computing unit to speed up the calculation, which leads to a waste of parallel resources. To solve this problem, Wang et al. proposed a collaborative parallel thin plate spline interpolation algorithm to accelerate the DEM generation of massive LiDAR point clouds. In this collaborative parallel algorithm, the input point cloud is first decomposed into a collection of discrete blocks and encapsulated as a general task object to shield the heterogeneous execution models of different processing units. Experimental results show that the collaborative parallel algorithm they proposed can achieve a maximum acceleration time of about 19.6. Compared with pure CPU and GPU parallel algorithms, the performance improvement rates are 54% and 44%, respectively. But their algorithm is still slightly insufficient in actual calculations [3]. Lee et al. proposed a received signal interpolation method to improve the performance of the multiple signal classification (MUSIC) algorithm. The simulation results show that the angular resolution of this method is better than that of the traditional MUSIC algorithm. In addition, they applied the proposed scheme to the actual data measured at the test site, providing a more enhanced DOA estimation result. However, their algorithm is not obvious in actual effect [4]. Mei et al. focused on the design and implementation of parallel adaptive inverse distance weighting (AIDW) interpolation algorithm using graphics processing unit (GPU). AIDW is an improved version of the standard IDW, which can adaptively determine power parameters according to the spatial distribution pattern of data points to achieve more accurate predictions than IDW predictions. Experimental results show that there is no significant difference in the calculation efficiency when using different data layouts; the tiled version is always slightly faster than the plain version; the acceleration achieved in single-precision is up to 763 (on GPU M5000), and the highest acceleration obtained in double-precision is 197 (on GPU K40c). However, the algorithm they proposed has not improved much in computational efficiency [5]. Zhang et al. used machine learning methods to design an interpolation algorithm based on Gaussian process regression. Their method uses a multiscale kernel function to process two-dimensional space meteorological ocean processes and introduces multiscale physical feature information (sea surface wind stress, sea surface heat flux, and ocean current velocity). This greatly improves the spatial resolution and interpolation accuracy of ocean features. They verified the effectiveness of the algorithm through interpolation experiments related to sea surface temperature (SST). The root mean square error (RMSE) of the interpolation algorithm is 38.9%, 43.7%, and 62.4% lower than bilinear interpolation, bicubic interpolation, and nearest neighbor interpolation, respectively. Their algorithm has acceptable running time cost and good time and space generalization. However, the calculation of the algorithm is slightly complicated and needs to be simplified [6]. Liu et al. proposed a non-uniform spatiotemporal kriging interpolation method. It breaks through the limitation of Euclidean distance in the space dimension and at the same time breaks through the linear relationship in the time dimension. They used the time-space optimal weight combination to construct the space-time deformation field model and then optimized it through the particle swarm optimization algorithm. They generalized ordinary kriging interpolation to non-uniform spatiotemporal kriging interpolation under space-time constraints. Their method has been successfully applied to the interpolation of landslide displacement monitoring data, which provides better data for the study of landslide disasters, and has important practical significance for the prevention and prediction of landslide disasters. It is just that the algorithm is not yet fully perfected, and further research is needed [7]. The air quality index (AQI) monitoring stations are sparsely distributed, and the spatial interpolation is less accurate than the existing methods. Therefore, Shi et al. proposed a new algorithm based on the extended field strength model. Experimental results prove that the AQI interpolation accuracy of the algorithm designed by Shi et al. is higher, and the two-parameter model obtains the highest accuracy. The algorithm is suitable for spatial interpolation of sparse data with a fixed number and location and can be used for spatial data with other types and dimensions. However, the applicability of this algorithm in spatial interpolation is not enough to meet practical applications [8].


## 3. Implementation Method of 3D Stratum Space Interpolation Algorithm for Subway Tunnels Based on BIM and Data Fusion

### 3.1. Data Fusion

Data fusion refers to the process of using computer technology to process information under specific standards [3, 9–12]. The process can automatically analyze, optimize, and integrate sensor observations that require sufficient time to complete the required decision-making and evaluation tasks. This process can automatically analyze, optimize, and synthesize time series multisensor observations to complete the required set of decision-making and evaluation tasks [6, 13]. With the continuous development of science and technology, the data environment in future life will become more complex, and the amount of data will increase exponentially. At the same time, there is still a lot of uncertain and false information [14]. In this case, the continuous development of data fusion technology will improve the efficiency of information processing and provide accurate, timely, and effective information support for data decision-making [15].Figure2 is an illustration of a global data management strategy for data fusion.Figure 2
Data fusion strategy management.Special tools are provided internally to solve the problem of data specification adjustment, without the need to propose an overall specification in advance, and only perform dynamic integration when business occurs [16–20]. While realizing standardized central data management, it realizes centralized data management and supports multiple synchronous or asynchronous distributed data applications [21, 22]. The architecture diagram of the data fusion center is shown in Figure 3.Figure 3
Data fusion center architecture diagram.


### 3.2. BIM

BIM (building information modeling) is an engineering data model that is based on three-dimensional digital technology and inherits various related information of construction projects. BIM is a digital expression of the physical and functional characteristics of engineering project facilities. Building information modeling is also a digital method applied to design, construction, and management. This method supports the integrated management environment of the construction project, which can significantly improve the efficiency and greatly reduce the risk of the construction project in its entire process.As a technology and tool for integrated construction project design, construction, and operation life cycle information management, BIM is relatively mature in international applications. With the rapid development of the domestic tunnel engineering industry, the tunnel design technology is getting higher and higher, and the design of the subway section tunnel is taken as an example. Considering the existing problems in its current design and the development trend of the industry, the reform of the design concept is the general trend and imperative. The advanced design concepts and thinking of BIM can bring good solutions to design [23, 24]. The future development trend of computer-aided architectural design is bound to be from 2D to 3D, from CAD to 3D modeling with BIM technology as the core [25]. Based on the fact that the BIM model of subway section tunnels is very rare, it has great research value. This also brings new ways and methods to the design of subway section tunnels.


### 3.3. Spatial Interpolation Algorithm

Among the interpolation algorithms, the research of interpolation algorithms in two-dimensional space started relatively early. Therefore, it is more complete than the difference of three-dimensional spatial data, which is roughly divided into functional methods, spatial statistical methods, and other categories [26, 27]. But if it wants to study and analyze the internal space of the tunnel, it must be interpolated in the three-dimensional space [28]. Although the research content of three-dimensional space interpolation is relatively limited at the beginning, however, with the rapid development and improvement of computer three-dimensional visualization technology, more and more scholars have begun to invest in research. This makes the three-dimensional spatial data interpolation of the tunnel get more and more attention, and even develop into an independent subject [29, 30]. There are many spatial interpolation methods, and there are many different classification methods based on different standards. According to the mathematical principle of its realization, it can be divided into two types, namely, the global interpolation method and the local interpolation method. Local interpolation includes inverse distance weighting method and radial basis interpolation method. The radial basis interpolation method is suitable for interpolating a large number of point data, while requiring a smooth surface. However, there are great uncertainties for short distances, so this article will not introduce them.


#### 3.3.1. Inverse Distance Weighting Method

The inverse distance weighting method is the most commonly used spatial interpolation method. It is based on the principle of similarity; that is, the similarity between objects depends on the distance between the two. The closer the distance is, the more similar it is, and on the contrary, the farther the distance is, the smaller the similarity will be [31, 32]. That is to say, the estimation of the interpolation point value is greatly affected by the closer sample points in the research range and less affected by the distant sample points. Therefore, the distance as the only influencing factor directly determines the data value of the point to be inserted. The basic idea of this method is that the smaller the distance between the insertion point and the sample point, the greater the weight, and the greater the distance, the smaller the weight. In other words, the relationship between the two is inversely proportional. In feature data interpolation, the feature value at the insertion point is calculated by calculating the weighted average of the features of all known points in the search area. If it uses the interpolation function R(a, b, c) to represent the attribute value of the point to be interpolated, the attribute value of each known point is Ri. Then, the attribute value R(a, b, c) at the point to be inserted can be expressed by the following formula:(1)Ra,b,c=∑i=1nsia,b,cη,(2)Ra,b,c=∑i=1nRi×Qia,b,c.Among them,Qia,b,c can be expressed as(3)Qia,b,c=∑i=1n∏j=1sia,b,cη.Among them,sia,b,c can be expressed as(4)sia,b,c=a−ai2+b−bi2+c−ci2.The previous formula represents the distance from pointai,bi,ci to unknown point a,b,c, and η is the power exponent, and the default value is 2. Under normal circumstances, the most reasonable interpolation result can be obtained by controlling the value within 0.5–3. The algorithm diagram is shown in Figure 4.Figure 4
Schematic diagram of inverse distance weighted interpolation algorithm.In spatial interpolation, the longer the distance of the spatial bright spot, the smaller the similarity. When the distance is far enough, the similarity can be approximated to zero. Therefore, one point can be removed from the calculation of another point. If the sample points used in the inverse distance weighting interpolation are unevenly distributed, the accuracy of the interpolation result will decrease. There are two ways to deal with this problem. One is to sort a specified number of points according to the distance from small to large, and the other is to determine the search neighborhood. This article implements the inverse distance weighting method; the specific steps are as follows:(5)R∗a0=∑i=1nωiRai.Here,a0 is the estimated value, n is the number of sample points used for interpolation, Rai is the measured value at the sample point ai, and ωi is the contribution weight of the i-th sample point to the estimated point. Its formula is expressed as(6)ωi=ti−s∑i=1nti−s,∑i=1nωi=1.Here,ti is the distance between the estimated point and the sample point and S is the power of the distance, which controls the change of the weighting factor in the process of increasing the distance and the minimum average error when selecting the standard.(7)Kijs=∑i=0nKiDjs,0≤s≤1.Among them,Kijs represents the curve section of the K section, and Djs in the above formula is the basis function of the D section, and its expression is(8)Djs=∑j=0n−i−1js+n−i−jn.Whenn = 1, the basis function is(9)D0,1s=1−s,D1,1s=s,0≤s≤1.Whenn = 2, its basis function is(10)sD0,2s=121−s2,D1,2s=12−2s2+2s+1,D2,2s=12s2.0≤s≤1,Whenn = 3, its basis function is(11)sD0,3s=161−s3,D1,3s=163s3+6s2+4,D2,3s=16−3s3+3s2+3s+1,D3,3s=16s3.0≤s≤1,The cubic spline curve is determined by four adjacent vertices, and the curve formula of the s-th segment is(12)Wi,js=D0,3sWk+D1,3sWk+1+D2,3sWk+2+D3,3sWk+3.The hierarchical interpolation implemented in this paper adopts a uniform spline function, which has a continuous nature at the connection.In order to make the difference result more accurate, a multi-level interpolation algorithm is added on this basis. The characteristic of this algorithm is that it only produces a fitted surface within a certain error range and cannot guarantee that the interpolated surface completely passes through all discrete point sets.Letψ be a data volume existing in the coordinate system, there is a set of discrete points R=ai,bi,ci,di in the three-dimensional space, and ai,bi,ci,di is any point located in the data volume ψ. In order to better represent the discrete point set R, it can be expressed as(13)sa,b,c=∑r=03∑i=03∑j=03TraTibTjc,where T is the sequence control points in the grid, and these control points define the unknown value in the grid.


#### 3.3.2. Spatial Self-Oblique Variance Interpolation Method

The spatial self-oblique variance interpolation method is also called the kriging method, which is the best interpolation method applied to the mineral grade in geostatistics. Through regionalized variables, the spatial change of the spatial self-oblique variance interpolation method can be divided into three parts: trend, correlation structure, and error. The advantage of kriging interpolation is based on geostatistics as its solid theoretical foundation. It can overcome the problem that the error in interpolation is difficult to analyze, and can make a point-by-point theoretical estimation of the error.Its formula is expressed as(14)R=∑i=1mωiRi,∑i=1Mωi=1,(15)∑i=1mωi=10≤m≤1.Here, the estimated value of the point to be estimated is expressed asR, the measured value of the sample point is expressed as Ri, the number of sample points involved in the interpolation calculation is m, and ωi is the weight coefficient of the sample point.The determination of the weight coefficient must satisfy the unbiased estimation ofR, and the estimated variance is smaller than the variance produced by the linear combination of the actual measured values at the interpolation point. Through calculation, the minimum variance of R can be obtained.(16)εe2=∑i=1mωiσηi,η+γ,(17)∑i=1mωiσηi,ηj+γ=σηj,η.Here,σηi,ηj is the semivariance of the two samples i and j, and σηj,η is the semivariance between the jth sample and the sample point to be tested.(18)σk=12Nk∑i=1NkRsi−Rsi+k2.Here, the semivariance function isσk and the step size is k.The scope of application of kriging interpolation is the spatial correlation of regionalized variables, that is, if the results of the variogram and structural analysis show that there is a spatial correlation between the regionalization variables. Kriging method can be used for interpolation or extrapolation; otherwise it is not feasible.


### 3.4. Interpolation Accuracy Evaluation Index

There are many methods for evaluating the accuracy of the interpolation results of the above two methods. Currently, cross-validation methods are usually used for testing, that is, the method of removing the data of some known points and using the data of other points to estimate these points to test the accuracy of the interpolation. The calculation formula is(19)MAE=1N∑i=1N|R∗a0−ra0,(20)RMSE=∑i=1NR∗a0−ra02N.To judge whether the interpolation method is good or bad, it is necessary to evaluate them. That is, the accuracy of each interpolation method is quantitatively compared. Among the above-mentioned interpolation methods, except for the ordinary kriging method, the error can be theoretically estimated point by point without the aid of the error test model. No other interpolation method can use its own model to theoretically estimate the interpolation error.


## 3.1. Data Fusion

Data fusion refers to the process of using computer technology to process information under specific standards [3, 9–12]. The process can automatically analyze, optimize, and integrate sensor observations that require sufficient time to complete the required decision-making and evaluation tasks. This process can automatically analyze, optimize, and synthesize time series multisensor observations to complete the required set of decision-making and evaluation tasks [6, 13]. With the continuous development of science and technology, the data environment in future life will become more complex, and the amount of data will increase exponentially. At the same time, there is still a lot of uncertain and false information [14]. In this case, the continuous development of data fusion technology will improve the efficiency of information processing and provide accurate, timely, and effective information support for data decision-making [15].Figure2 is an illustration of a global data management strategy for data fusion.Figure 2
Data fusion strategy management.Special tools are provided internally to solve the problem of data specification adjustment, without the need to propose an overall specification in advance, and only perform dynamic integration when business occurs [16–20]. While realizing standardized central data management, it realizes centralized data management and supports multiple synchronous or asynchronous distributed data applications [21, 22]. The architecture diagram of the data fusion center is shown in Figure 3.Figure 3
Data fusion center architecture diagram.


## 3.2. BIM

BIM (building information modeling) is an engineering data model that is based on three-dimensional digital technology and inherits various related information of construction projects. BIM is a digital expression of the physical and functional characteristics of engineering project facilities. Building information modeling is also a digital method applied to design, construction, and management. This method supports the integrated management environment of the construction project, which can significantly improve the efficiency and greatly reduce the risk of the construction project in its entire process.As a technology and tool for integrated construction project design, construction, and operation life cycle information management, BIM is relatively mature in international applications. With the rapid development of the domestic tunnel engineering industry, the tunnel design technology is getting higher and higher, and the design of the subway section tunnel is taken as an example. Considering the existing problems in its current design and the development trend of the industry, the reform of the design concept is the general trend and imperative. The advanced design concepts and thinking of BIM can bring good solutions to design [23, 24]. The future development trend of computer-aided architectural design is bound to be from 2D to 3D, from CAD to 3D modeling with BIM technology as the core [25]. Based on the fact that the BIM model of subway section tunnels is very rare, it has great research value. This also brings new ways and methods to the design of subway section tunnels.


## 3.3. Spatial Interpolation Algorithm

Among the interpolation algorithms, the research of interpolation algorithms in two-dimensional space started relatively early. Therefore, it is more complete than the difference of three-dimensional spatial data, which is roughly divided into functional methods, spatial statistical methods, and other categories [26, 27]. But if it wants to study and analyze the internal space of the tunnel, it must be interpolated in the three-dimensional space [28]. Although the research content of three-dimensional space interpolation is relatively limited at the beginning, however, with the rapid development and improvement of computer three-dimensional visualization technology, more and more scholars have begun to invest in research. This makes the three-dimensional spatial data interpolation of the tunnel get more and more attention, and even develop into an independent subject [29, 30]. There are many spatial interpolation methods, and there are many different classification methods based on different standards. According to the mathematical principle of its realization, it can be divided into two types, namely, the global interpolation method and the local interpolation method. Local interpolation includes inverse distance weighting method and radial basis interpolation method. The radial basis interpolation method is suitable for interpolating a large number of point data, while requiring a smooth surface. However, there are great uncertainties for short distances, so this article will not introduce them.


### 3.3.1. Inverse Distance Weighting Method

The inverse distance weighting method is the most commonly used spatial interpolation method. It is based on the principle of similarity; that is, the similarity between objects depends on the distance between the two. The closer the distance is, the more similar it is, and on the contrary, the farther the distance is, the smaller the similarity will be [31, 32]. That is to say, the estimation of the interpolation point value is greatly affected by the closer sample points in the research range and less affected by the distant sample points. Therefore, the distance as the only influencing factor directly determines the data value of the point to be inserted. The basic idea of this method is that the smaller the distance between the insertion point and the sample point, the greater the weight, and the greater the distance, the smaller the weight. In other words, the relationship between the two is inversely proportional. In feature data interpolation, the feature value at the insertion point is calculated by calculating the weighted average of the features of all known points in the search area. If it uses the interpolation function R(a, b, c) to represent the attribute value of the point to be interpolated, the attribute value of each known point is Ri. Then, the attribute value R(a, b, c) at the point to be inserted can be expressed by the following formula:(1)Ra,b,c=∑i=1nsia,b,cη,(2)Ra,b,c=∑i=1nRi×Qia,b,c.Among them,Qia,b,c can be expressed as(3)Qia,b,c=∑i=1n∏j=1sia,b,cη.Among them,sia,b,c can be expressed as(4)sia,b,c=a−ai2+b−bi2+c−ci2.The previous formula represents the distance from pointai,bi,ci to unknown point a,b,c, and η is the power exponent, and the default value is 2. Under normal circumstances, the most reasonable interpolation result can be obtained by controlling the value within 0.5–3. The algorithm diagram is shown in Figure 4.Figure 4
Schematic diagram of inverse distance weighted interpolation algorithm.In spatial interpolation, the longer the distance of the spatial bright spot, the smaller the similarity. When the distance is far enough, the similarity can be approximated to zero. Therefore, one point can be removed from the calculation of another point. If the sample points used in the inverse distance weighting interpolation are unevenly distributed, the accuracy of the interpolation result will decrease. There are two ways to deal with this problem. One is to sort a specified number of points according to the distance from small to large, and the other is to determine the search neighborhood. This article implements the inverse distance weighting method; the specific steps are as follows:(5)R∗a0=∑i=1nωiRai.Here,a0 is the estimated value, n is the number of sample points used for interpolation, Rai is the measured value at the sample point ai, and ωi is the contribution weight of the i-th sample point to the estimated point. Its formula is expressed as(6)ωi=ti−s∑i=1nti−s,∑i=1nωi=1.Here,ti is the distance between the estimated point and the sample point and S is the power of the distance, which controls the change of the weighting factor in the process of increasing the distance and the minimum average error when selecting the standard.(7)Kijs=∑i=0nKiDjs,0≤s≤1.Among them,Kijs represents the curve section of the K section, and Djs in the above formula is the basis function of the D section, and its expression is(8)Djs=∑j=0n−i−1js+n−i−jn.Whenn = 1, the basis function is(9)D0,1s=1−s,D1,1s=s,0≤s≤1.Whenn = 2, its basis function is(10)sD0,2s=121−s2,D1,2s=12−2s2+2s+1,D2,2s=12s2.0≤s≤1,Whenn = 3, its basis function is(11)sD0,3s=161−s3,D1,3s=163s3+6s2+4,D2,3s=16−3s3+3s2+3s+1,D3,3s=16s3.0≤s≤1,The cubic spline curve is determined by four adjacent vertices, and the curve formula of the s-th segment is(12)Wi,js=D0,3sWk+D1,3sWk+1+D2,3sWk+2+D3,3sWk+3.The hierarchical interpolation implemented in this paper adopts a uniform spline function, which has a continuous nature at the connection.In order to make the difference result more accurate, a multi-level interpolation algorithm is added on this basis. The characteristic of this algorithm is that it only produces a fitted surface within a certain error range and cannot guarantee that the interpolated surface completely passes through all discrete point sets.Letψ be a data volume existing in the coordinate system, there is a set of discrete points R=ai,bi,ci,di in the three-dimensional space, and ai,bi,ci,di is any point located in the data volume ψ. In order to better represent the discrete point set R, it can be expressed as(13)sa,b,c=∑r=03∑i=03∑j=03TraTibTjc,where T is the sequence control points in the grid, and these control points define the unknown value in the grid.


### 3.3.2. Spatial Self-Oblique Variance Interpolation Method

The spatial self-oblique variance interpolation method is also called the kriging method, which is the best interpolation method applied to the mineral grade in geostatistics. Through regionalized variables, the spatial change of the spatial self-oblique variance interpolation method can be divided into three parts: trend, correlation structure, and error. The advantage of kriging interpolation is based on geostatistics as its solid theoretical foundation. It can overcome the problem that the error in interpolation is difficult to analyze, and can make a point-by-point theoretical estimation of the error.Its formula is expressed as(14)R=∑i=1mωiRi,∑i=1Mωi=1,(15)∑i=1mωi=10≤m≤1.Here, the estimated value of the point to be estimated is expressed asR, the measured value of the sample point is expressed as Ri, the number of sample points involved in the interpolation calculation is m, and ωi is the weight coefficient of the sample point.The determination of the weight coefficient must satisfy the unbiased estimation ofR, and the estimated variance is smaller than the variance produced by the linear combination of the actual measured values at the interpolation point. Through calculation, the minimum variance of R can be obtained.(16)εe2=∑i=1mωiσηi,η+γ,(17)∑i=1mωiσηi,ηj+γ=σηj,η.Here,σηi,ηj is the semivariance of the two samples i and j, and σηj,η is the semivariance between the jth sample and the sample point to be tested.(18)σk=12Nk∑i=1NkRsi−Rsi+k2.Here, the semivariance function isσk and the step size is k.The scope of application of kriging interpolation is the spatial correlation of regionalized variables, that is, if the results of the variogram and structural analysis show that there is a spatial correlation between the regionalization variables. Kriging method can be used for interpolation or extrapolation; otherwise it is not feasible.


## 3.3.1. Inverse Distance Weighting Method

The inverse distance weighting method is the most commonly used spatial interpolation method. It is based on the principle of similarity; that is, the similarity between objects depends on the distance between the two. The closer the distance is, the more similar it is, and on the contrary, the farther the distance is, the smaller the similarity will be [31, 32]. That is to say, the estimation of the interpolation point value is greatly affected by the closer sample points in the research range and less affected by the distant sample points. Therefore, the distance as the only influencing factor directly determines the data value of the point to be inserted. The basic idea of this method is that the smaller the distance between the insertion point and the sample point, the greater the weight, and the greater the distance, the smaller the weight. In other words, the relationship between the two is inversely proportional. In feature data interpolation, the feature value at the insertion point is calculated by calculating the weighted average of the features of all known points in the search area. If it uses the interpolation function R(a, b, c) to represent the attribute value of the point to be interpolated, the attribute value of each known point is Ri. Then, the attribute value R(a, b, c) at the point to be inserted can be expressed by the following formula:(1)Ra,b,c=∑i=1nsia,b,cη,(2)Ra,b,c=∑i=1nRi×Qia,b,c.Among them,Qia,b,c can be expressed as(3)Qia,b,c=∑i=1n∏j=1sia,b,cη.Among them,sia,b,c can be expressed as(4)sia,b,c=a−ai2+b−bi2+c−ci2.The previous formula represents the distance from pointai,bi,ci to unknown point a,b,c, and η is the power exponent, and the default value is 2. Under normal circumstances, the most reasonable interpolation result can be obtained by controlling the value within 0.5–3. The algorithm diagram is shown in Figure 4.Figure 4
Schematic diagram of inverse distance weighted interpolation algorithm.In spatial interpolation, the longer the distance of the spatial bright spot, the smaller the similarity. When the distance is far enough, the similarity can be approximated to zero. Therefore, one point can be removed from the calculation of another point. If the sample points used in the inverse distance weighting interpolation are unevenly distributed, the accuracy of the interpolation result will decrease. There are two ways to deal with this problem. One is to sort a specified number of points according to the distance from small to large, and the other is to determine the search neighborhood. This article implements the inverse distance weighting method; the specific steps are as follows:(5)R∗a0=∑i=1nωiRai.Here,a0 is the estimated value, n is the number of sample points used for interpolation, Rai is the measured value at the sample point ai, and ωi is the contribution weight of the i-th sample point to the estimated point. Its formula is expressed as(6)ωi=ti−s∑i=1nti−s,∑i=1nωi=1.Here,ti is the distance between the estimated point and the sample point and S is the power of the distance, which controls the change of the weighting factor in the process of increasing the distance and the minimum average error when selecting the standard.(7)Kijs=∑i=0nKiDjs,0≤s≤1.Among them,Kijs represents the curve section of the K section, and Djs in the above formula is the basis function of the D section, and its expression is(8)Djs=∑j=0n−i−1js+n−i−jn.Whenn = 1, the basis function is(9)D0,1s=1−s,D1,1s=s,0≤s≤1.Whenn = 2, its basis function is(10)sD0,2s=121−s2,D1,2s=12−2s2+2s+1,D2,2s=12s2.0≤s≤1,Whenn = 3, its basis function is(11)sD0,3s=161−s3,D1,3s=163s3+6s2+4,D2,3s=16−3s3+3s2+3s+1,D3,3s=16s3.0≤s≤1,The cubic spline curve is determined by four adjacent vertices, and the curve formula of the s-th segment is(12)Wi,js=D0,3sWk+D1,3sWk+1+D2,3sWk+2+D3,3sWk+3.The hierarchical interpolation implemented in this paper adopts a uniform spline function, which has a continuous nature at the connection.In order to make the difference result more accurate, a multi-level interpolation algorithm is added on this basis. The characteristic of this algorithm is that it only produces a fitted surface within a certain error range and cannot guarantee that the interpolated surface completely passes through all discrete point sets.Letψ be a data volume existing in the coordinate system, there is a set of discrete points R=ai,bi,ci,di in the three-dimensional space, and ai,bi,ci,di is any point located in the data volume ψ. In order to better represent the discrete point set R, it can be expressed as(13)sa,b,c=∑r=03∑i=03∑j=03TraTibTjc,where T is the sequence control points in the grid, and these control points define the unknown value in the grid.


## 3.3.2. Spatial Self-Oblique Variance Interpolation Method

The spatial self-oblique variance interpolation method is also called the kriging method, which is the best interpolation method applied to the mineral grade in geostatistics. Through regionalized variables, the spatial change of the spatial self-oblique variance interpolation method can be divided into three parts: trend, correlation structure, and error. The advantage of kriging interpolation is based on geostatistics as its solid theoretical foundation. It can overcome the problem that the error in interpolation is difficult to analyze, and can make a point-by-point theoretical estimation of the error.Its formula is expressed as(14)R=∑i=1mωiRi,∑i=1Mωi=1,(15)∑i=1mωi=10≤m≤1.Here, the estimated value of the point to be estimated is expressed asR, the measured value of the sample point is expressed as Ri, the number of sample points involved in the interpolation calculation is m, and ωi is the weight coefficient of the sample point.The determination of the weight coefficient must satisfy the unbiased estimation ofR, and the estimated variance is smaller than the variance produced by the linear combination of the actual measured values at the interpolation point. Through calculation, the minimum variance of R can be obtained.(16)εe2=∑i=1mωiσηi,η+γ,(17)∑i=1mωiσηi,ηj+γ=σηj,η.Here,σηi,ηj is the semivariance of the two samples i and j, and σηj,η is the semivariance between the jth sample and the sample point to be tested.(18)σk=12Nk∑i=1NkRsi−Rsi+k2.Here, the semivariance function isσk and the step size is k.The scope of application of kriging interpolation is the spatial correlation of regionalized variables, that is, if the results of the variogram and structural analysis show that there is a spatial correlation between the regionalization variables. Kriging method can be used for interpolation or extrapolation; otherwise it is not feasible.


## 3.4. Interpolation Accuracy Evaluation Index

There are many methods for evaluating the accuracy of the interpolation results of the above two methods. Currently, cross-validation methods are usually used for testing, that is, the method of removing the data of some known points and using the data of other points to estimate these points to test the accuracy of the interpolation. The calculation formula is(19)MAE=1N∑i=1N|R∗a0−ra0,(20)RMSE=∑i=1NR∗a0−ra02N.To judge whether the interpolation method is good or bad, it is necessary to evaluate them. That is, the accuracy of each interpolation method is quantitatively compared. Among the above-mentioned interpolation methods, except for the ordinary kriging method, the error can be theoretically estimated point by point without the aid of the error test model. No other interpolation method can use its own model to theoretically estimate the interpolation error.


## 4. Experiment of 3D Stratum Space Interpolation Algorithm for Subway Tunnels Based on BIM and Data Fusion

### 4.1. Modeling Visualization Process

The three-dimensional visualization of a geological model generally includes several steps such as data collection and processing, model construction, and mapping. Drilling holes to obtain geological information is the most important method of field data collection, and drilling data is also the most common type of data in the process of geological modeling. Drilling is the process of arranging sampling points in the sampling area according to a certain rule. Because of its high cost, the number of sample points for drilling holes in a certain geological area is usually limited, and the fewest points should be used to maximize the effect. Data processing is mainly to properly eliminate or encrypt sample points obtained from data collection. It first needs to delete invalid or overlapping data. In addition, when there are few data points, the grid distribution of the data is generally achieved by spatial interpolation to make the data dense and continuous. The model structure is mainly divided into the structure of the surface model and the structure of the solid model. The surface model is the use of geological section or discrete point data to generate a three-dimensional geological model represented by a curved surface. The solid model uses a series of geophysical data to obtain the three-dimensional attribute model of the geological body through trend surface fitting or spatial data interpolation. Drawing mapping realizes the visualization of geological model drawing on the computer screen, which is a typical application of computer graphics in the field of geosciences. In order to make the model imaging more realistic, it is necessary to set attributes such as material, lighting, and color.


### 4.2. 3D Stratum Hierarchical Structure Model of Subway Tunnels

In the actual data collection process, it is difficult to obtain points regularly, so the spatial data obtained by sampling is generally discretely distributed. The values of other unknown points in the same area are usually estimated by spatial data interpolation to achieve the purpose of data encryption and regular grid distribution. The main goal of the algorithm in this paper is to use a 3D solid model suitable for the subway tunnel layer to correctly represent the robustness of the 3D layer. It also designs a modeling algorithm for the model to achieve seamless integration of each sub-body model of the complex stratigraphic body. Therefore, it is very important to find a three-dimensional solid model that can accurately mark the complex strata.As BIM becomes more mature in the construction industry, the use of BIM in the tunnel industry is still in the exploratory stage. Thanks to the powerful computer-aided design function, the design efficiency of tunnel engineering has been greatly improved. In order to complete the modeling and seamless integration of the complex stratigraphic body itself and various internal geological sub-bodies, a hierarchical structure model suitable for the three-dimensional strata of the subway tunnel was constructed. The model is shown in Figure5.Figure 5
Three-dimensional hierarchical structure model of subway tunnels.The three-dimensional hierarchical structure model has a clear hierarchical structure, which is divided into two layers, namely, the logic layer and the realization layer. The logical layer conceptually expresses the stratigraphic body, and the polyhedron belongs to VGTP. It is used to express special geological layers such as lenses, intrusions, and folds. According to the geological description of the VGTP geological sub-body, the realization layer separates the protons of various places from the VGTP in the form of polyhedrons through the cutting of the surface. The realization layer uses a unified TEN voxel to accurately express the protons of different regions. According to the geological description of VGTP, it uses a designed transformation algorithm to cut the surface through several steps. It separates the protons located in it in the form of middle polyhedrons and then transforms these middle polyhedrons into a series of seamlessly connected tetrahedrons. When the logic layer model is transformed into the realization layer TENs, the spatial analysis, virtual space manipulation, and 3D visualization functions will become practical, reasonable, and feasible. This is also more in line with the natural characteristics of the three-dimensional solid of the geological body.


### 4.3. The Business Process of the Three-Dimensional Stratigraphic Space Hierarchical Structure Model

In order to form the triangulation layer, different geological layers must be managed, the points with the same attribute value are placed in the same layer, and they are divided into three exemptions. Due to the limited number of points with the same flexural value, interpolation and fitting of spatial discrete data are also needed. This article intends to use the kriging interpolation method to achieve. After all the strata are drawn, these strata need to be stitched together at once, and the stitching between strata also needs to be connected by triangles. After filling in the colors, the geological body model has been initially realized. Finally, there are related interfaces for the scalability of the model under the interactive function module. The business process of the 3D geological model is shown in Figure6.Figure 6
Business flow chart of the 3D hierarchical structure model of subway tunnel.


### 4.4. Interpolation Analysis of Three-Dimensional Stratigraphic Space

It is necessary to find the 10 nearest sample sites of the strata to be interpolated, obtain the actual measured value of each stratigraphic sample site closest tot, and repeat the steps until all stratigraphic sites are estimated by the expansion method. Table 1 shows the cross-validation results of the reduction and expansion methods using IDW spatial interpolation.Table 1
Cross-validation results of reduction and expansion using IDW spatial interpolation.
Interpolation methodInspection standardsReductionExpansion methodIDWMAE1.5361.562RMSIE2.16332.1655IDWMAE197.6204.6RMSIE275.6286.8It takes all the original points as the detection points, and their inclination angle value as the standard value. It compares the inclination angle values at these points with the trends generated by the above methods and calculates the inclination angle value with the occurrence point. Figure7 is the result of intersecting the ground plane of the subway tunnel.Figure 7
Intersection result of subway tunnel ground plane.After the above-mentioned interpolation, the triangulation network at each level and the triangulation network at partial fault levels have elevation information. This kind of situation can also occur when the formation is missing or the formation is pinched. Such an intersection can be corrected by adjusting the elevation of the formation. As shown in Figure7, after processing, in some areas, the elevation values of stratum 2 and 3 are equal to the elevation value of stratum 4. In essence, the formations 2 and 3 are missing, and after cross-processing, some formations can be treated as pinching.Table2 shows the relative error indicators of the inclination angle of different methods.Table 2
Various indicators of relative error of inclination angle of different methods.
MethodMaximumMinimumSumMeanSDTIN262.42.011193.3760.636572.048592Higher-order surface function58.334.88488.1526.171413.285446Tension spline function55.01561.04262255.164213.1916213.24285It can be seen that the relative error of the inclination angle of the tension spline function has the smallest average value and standard deviation compared with other methods. Therefore, we believe that the trend surface generated by the spatial interpolation method of the tension spline function is more accurate. When this method is used to interpolate the occurrence of the strata, the interpolation accuracy can be effectively improved.


### 4.5. Comparison and Analysis of Interpolation Algorithms

We will test the two interpolation methods introduced in this article and compare the interpolation results with the test values at each checkpoint. The accuracy of various interpolation methods is determined by the calculation error, the average absolute error, and the root mean square. For attack indicators, the interpolation accuracy of the IDW interpolation method is generally higher than that of the conventional kriging method, and the average absolute error and root mean square error are relatively small and superior. The result of the verification point test value is shown in Figure8.Figure 8
Verification point test value result.The degree of density deviation from the test value of the verification point is shown in Figure9.Figure 9
The degree of density deviation from the test value of the verification point.In the interpolation calculation, the inverse distance weighting method and the ordinary kriging method are two extremely widely used methods. The above two methods can simulate the spatial distribution characteristics of the erosion factor attributes in the buffer, but the accuracy of the interpolation data needs to be considered.The interpolation results vary with the density of the mesh. When the control point mesh is sparsely divided, the interpolation result has good smoothness, but the accuracy is low. Conversely, as the control mesh is finer, the interpolation result will be closer to the known discrete sampling points in its neighborhood, the interpolation accuracy is higher, and the original surface and the interpolation surface have good fusion. However, unrestricted fine control meshing also has certain drawbacks. First, as the meshing density increases, the amount of calculation will increase. On the other hand, after the mesh is subdivided to a certain extent, there will be few data points in each divided area, and the interpolation may cause local peaks or even serious distortion of the interpolation results. Figure10 is a comparison between the original data and the spatial self-oblique variance interpolation data.Figure 10
Original data and spatial self-oblique variance interpolation algorithm.The contrast space is self-oblique variance interpolation method, although the inverse distance weighting method has the advantages of low computational cost and small storage space. However, when the data distribution is extremely uneven, the interpolation results obtained by using this method often cannot meet the accuracy requirements in actual work. Therefore, the use of inverse distance-weighted interpolation is generally considered when data collection is dense and evenly distributed. In addition, it can be known from the incompleteness of the sampled data and the realization principle of the interpolation method. The interpolation data obtained by the inverse distance weighting method has the characteristics of being greater than the minimum value of the data and less than the maximum value of the data. Using the spatial self-oblique variance method for interpolation avoids the “bull eye” phenomenon and the “bulge” phenomenon. In the same slope area, the effect of using spatial self-oblique variance interpolation is smoother, while the surface obtained by inverse distance weighted interpolation appears uneven. Table3 shows the comparison results between the original data and the interpolated data.Table 3
Comparison of the original data and the data obtained by each interpolation.
IDOriginal coordinatesOriginal elevationOrdinary krigingImprove krigingABCCC15120−3032.55−3035.19−3038.4720522.64−3045.63−3045.22−3056.133015542.84−3051.17−3059.85−3088.29402821.26−3060.46−3071.63−3132.68503830.96−3170.27−3211.24−3255.60………………812762.4425882.1424.65418.42419.068131523.1726262.4477.92478.10479.338142260.832677.6524.51520.51532.448153744.222893.5543.16589.27609.53The realization principle of the spatial self-oblique variance interpolation algorithm is to use a polynomial fitting method to generate a smooth interpolation surface and reduce the error between the interpolation surface and the original discrete points through multiple fittings. Therefore, the interpolation result obtained by this method has higher accuracy. In places with less data, the general trend of the surface can be maintained, and in places with a large number, surfaces closer to scattered data can be generated. In terms of interpolation efficiency, spatial self-oblique variance interpolation is suitable for large-scale data interpolation. Inverse distance-weighted interpolation is generally used for medium- and small-scale data interpolation, and spatial self-oblique variance interpolation has obvious speed advantages. The statistical information of various interpolation transformations is shown in Table4.Table 4
Data transformation statistics.
StatisticsNo transformationLogarithmic transformationPower exponent transformationMean3974.68.19639077240Std. dev.1543.20.43156752000Kurtosis2.64612.87414.0166Skewness0.4316−0.524031.1702Among them, kurtosis and skewness can reflect whether the data obey normal distribution. Kurtosis is used to describe the index of the height of data distribution, skewness is used to describe the symmetry of data distribution, and the kurtosis and skewness of normal distribution should be equal to 0.


### 4.6. Comparison Results of Various Interpolation Algorithms

Kriging interpolation method is to find the optimal, linear, and fair geometric features and spatial structure after considering the sample point distribution, density, size, and estimated spatial distribution position of the sample points. The normal distribution and quality of the data are the key to the success of the kriging interpolation method. Figure11 is a comparative analysis of the performance of the elevation interpolation results obtained by different interpolation algorithms on the same stratum surface.Figure 11
Performance comparison of elevation interpolation.Table5 is a comparison of the results of the interpolation of the three-dimensional attribute data volume for each difference algorithm.Table 5
Comparison of the results of the difference algorithm in the interpolation of the three-dimensional attribute data volume.
IDWKrigingSmooth krigingImprove algorithmMaximum error137.2583.2139.3743.23Minimum error1.343.013.8227.87Average error96.6262.4529.1638.98Among them, the sliding kriging algorithm is an interpolation method that finds the global variable range and performs the calculation of the local kriging equations according to the known variable range. The above data shows that the overall error of the improved interpolation algorithm is smaller than the inverse distance weighted interpolation. Compared with the inverse distance weighted interpolation method, it also estimates the unknown points by calculating the weighted average of the known points. The kriging method is based on the spatial variogram, and the variogram is determined by the existing spatial sampling points. It makes full use of the spatial structure and randomness of data points, and the interpolation effect is more realistic. The inverse distance weighting method is an interpolation method that constructs an approximation function to fit the original point set. This method has a certain error, and the interpolation result is too smooth due to the influence of the grid density, and some details are lost. Compared with the inverse distance weighting method and the ordinary kriging method, the improved interpolation algorithm in this paper improves the performance by 23.4%. And the average error value of the algorithm is 38.98, which has been greatly improved compared with other algorithms.


## 4.1. Modeling Visualization Process

The three-dimensional visualization of a geological model generally includes several steps such as data collection and processing, model construction, and mapping. Drilling holes to obtain geological information is the most important method of field data collection, and drilling data is also the most common type of data in the process of geological modeling. Drilling is the process of arranging sampling points in the sampling area according to a certain rule. Because of its high cost, the number of sample points for drilling holes in a certain geological area is usually limited, and the fewest points should be used to maximize the effect. Data processing is mainly to properly eliminate or encrypt sample points obtained from data collection. It first needs to delete invalid or overlapping data. In addition, when there are few data points, the grid distribution of the data is generally achieved by spatial interpolation to make the data dense and continuous. The model structure is mainly divided into the structure of the surface model and the structure of the solid model. The surface model is the use of geological section or discrete point data to generate a three-dimensional geological model represented by a curved surface. The solid model uses a series of geophysical data to obtain the three-dimensional attribute model of the geological body through trend surface fitting or spatial data interpolation. Drawing mapping realizes the visualization of geological model drawing on the computer screen, which is a typical application of computer graphics in the field of geosciences. In order to make the model imaging more realistic, it is necessary to set attributes such as material, lighting, and color.


## 4.2. 3D Stratum Hierarchical Structure Model of Subway Tunnels

In the actual data collection process, it is difficult to obtain points regularly, so the spatial data obtained by sampling is generally discretely distributed. The values of other unknown points in the same area are usually estimated by spatial data interpolation to achieve the purpose of data encryption and regular grid distribution. The main goal of the algorithm in this paper is to use a 3D solid model suitable for the subway tunnel layer to correctly represent the robustness of the 3D layer. It also designs a modeling algorithm for the model to achieve seamless integration of each sub-body model of the complex stratigraphic body. Therefore, it is very important to find a three-dimensional solid model that can accurately mark the complex strata.As BIM becomes more mature in the construction industry, the use of BIM in the tunnel industry is still in the exploratory stage. Thanks to the powerful computer-aided design function, the design efficiency of tunnel engineering has been greatly improved. In order to complete the modeling and seamless integration of the complex stratigraphic body itself and various internal geological sub-bodies, a hierarchical structure model suitable for the three-dimensional strata of the subway tunnel was constructed. The model is shown in Figure5.Figure 5
Three-dimensional hierarchical structure model of subway tunnels.The three-dimensional hierarchical structure model has a clear hierarchical structure, which is divided into two layers, namely, the logic layer and the realization layer. The logical layer conceptually expresses the stratigraphic body, and the polyhedron belongs to VGTP. It is used to express special geological layers such as lenses, intrusions, and folds. According to the geological description of the VGTP geological sub-body, the realization layer separates the protons of various places from the VGTP in the form of polyhedrons through the cutting of the surface. The realization layer uses a unified TEN voxel to accurately express the protons of different regions. According to the geological description of VGTP, it uses a designed transformation algorithm to cut the surface through several steps. It separates the protons located in it in the form of middle polyhedrons and then transforms these middle polyhedrons into a series of seamlessly connected tetrahedrons. When the logic layer model is transformed into the realization layer TENs, the spatial analysis, virtual space manipulation, and 3D visualization functions will become practical, reasonable, and feasible. This is also more in line with the natural characteristics of the three-dimensional solid of the geological body.


## 4.3. The Business Process of the Three-Dimensional Stratigraphic Space Hierarchical Structure Model

In order to form the triangulation layer, different geological layers must be managed, the points with the same attribute value are placed in the same layer, and they are divided into three exemptions. Due to the limited number of points with the same flexural value, interpolation and fitting of spatial discrete data are also needed. This article intends to use the kriging interpolation method to achieve. After all the strata are drawn, these strata need to be stitched together at once, and the stitching between strata also needs to be connected by triangles. After filling in the colors, the geological body model has been initially realized. Finally, there are related interfaces for the scalability of the model under the interactive function module. The business process of the 3D geological model is shown in Figure6.Figure 6
Business flow chart of the 3D hierarchical structure model of subway tunnel.


## 4.4. Interpolation Analysis of Three-Dimensional Stratigraphic Space

It is necessary to find the 10 nearest sample sites of the strata to be interpolated, obtain the actual measured value of each stratigraphic sample site closest tot, and repeat the steps until all stratigraphic sites are estimated by the expansion method. Table 1 shows the cross-validation results of the reduction and expansion methods using IDW spatial interpolation.Table 1
Cross-validation results of reduction and expansion using IDW spatial interpolation.
Interpolation methodInspection standardsReductionExpansion methodIDWMAE1.5361.562RMSIE2.16332.1655IDWMAE197.6204.6RMSIE275.6286.8It takes all the original points as the detection points, and their inclination angle value as the standard value. It compares the inclination angle values at these points with the trends generated by the above methods and calculates the inclination angle value with the occurrence point. Figure7 is the result of intersecting the ground plane of the subway tunnel.Figure 7
Intersection result of subway tunnel ground plane.After the above-mentioned interpolation, the triangulation network at each level and the triangulation network at partial fault levels have elevation information. This kind of situation can also occur when the formation is missing or the formation is pinched. Such an intersection can be corrected by adjusting the elevation of the formation. As shown in Figure7, after processing, in some areas, the elevation values of stratum 2 and 3 are equal to the elevation value of stratum 4. In essence, the formations 2 and 3 are missing, and after cross-processing, some formations can be treated as pinching.Table2 shows the relative error indicators of the inclination angle of different methods.Table 2
Various indicators of relative error of inclination angle of different methods.
MethodMaximumMinimumSumMeanSDTIN262.42.011193.3760.636572.048592Higher-order surface function58.334.88488.1526.171413.285446Tension spline function55.01561.04262255.164213.1916213.24285It can be seen that the relative error of the inclination angle of the tension spline function has the smallest average value and standard deviation compared with other methods. Therefore, we believe that the trend surface generated by the spatial interpolation method of the tension spline function is more accurate. When this method is used to interpolate the occurrence of the strata, the interpolation accuracy can be effectively improved.


## 4.5. Comparison and Analysis of Interpolation Algorithms

We will test the two interpolation methods introduced in this article and compare the interpolation results with the test values at each checkpoint. The accuracy of various interpolation methods is determined by the calculation error, the average absolute error, and the root mean square. For attack indicators, the interpolation accuracy of the IDW interpolation method is generally higher than that of the conventional kriging method, and the average absolute error and root mean square error are relatively small and superior. The result of the verification point test value is shown in Figure8.Figure 8
Verification point test value result.The degree of density deviation from the test value of the verification point is shown in Figure9.Figure 9
The degree of density deviation from the test value of the verification point.In the interpolation calculation, the inverse distance weighting method and the ordinary kriging method are two extremely widely used methods. The above two methods can simulate the spatial distribution characteristics of the erosion factor attributes in the buffer, but the accuracy of the interpolation data needs to be considered.The interpolation results vary with the density of the mesh. When the control point mesh is sparsely divided, the interpolation result has good smoothness, but the accuracy is low. Conversely, as the control mesh is finer, the interpolation result will be closer to the known discrete sampling points in its neighborhood, the interpolation accuracy is higher, and the original surface and the interpolation surface have good fusion. However, unrestricted fine control meshing also has certain drawbacks. First, as the meshing density increases, the amount of calculation will increase. On the other hand, after the mesh is subdivided to a certain extent, there will be few data points in each divided area, and the interpolation may cause local peaks or even serious distortion of the interpolation results. Figure10 is a comparison between the original data and the spatial self-oblique variance interpolation data.Figure 10
Original data and spatial self-oblique variance interpolation algorithm.The contrast space is self-oblique variance interpolation method, although the inverse distance weighting method has the advantages of low computational cost and small storage space. However, when the data distribution is extremely uneven, the interpolation results obtained by using this method often cannot meet the accuracy requirements in actual work. Therefore, the use of inverse distance-weighted interpolation is generally considered when data collection is dense and evenly distributed. In addition, it can be known from the incompleteness of the sampled data and the realization principle of the interpolation method. The interpolation data obtained by the inverse distance weighting method has the characteristics of being greater than the minimum value of the data and less than the maximum value of the data. Using the spatial self-oblique variance method for interpolation avoids the “bull eye” phenomenon and the “bulge” phenomenon. In the same slope area, the effect of using spatial self-oblique variance interpolation is smoother, while the surface obtained by inverse distance weighted interpolation appears uneven. Table3 shows the comparison results between the original data and the interpolated data.Table 3
Comparison of the original data and the data obtained by each interpolation.
IDOriginal coordinatesOriginal elevationOrdinary krigingImprove krigingABCCC15120−3032.55−3035.19−3038.4720522.64−3045.63−3045.22−3056.133015542.84−3051.17−3059.85−3088.29402821.26−3060.46−3071.63−3132.68503830.96−3170.27−3211.24−3255.60………………812762.4425882.1424.65418.42419.068131523.1726262.4477.92478.10479.338142260.832677.6524.51520.51532.448153744.222893.5543.16589.27609.53The realization principle of the spatial self-oblique variance interpolation algorithm is to use a polynomial fitting method to generate a smooth interpolation surface and reduce the error between the interpolation surface and the original discrete points through multiple fittings. Therefore, the interpolation result obtained by this method has higher accuracy. In places with less data, the general trend of the surface can be maintained, and in places with a large number, surfaces closer to scattered data can be generated. In terms of interpolation efficiency, spatial self-oblique variance interpolation is suitable for large-scale data interpolation. Inverse distance-weighted interpolation is generally used for medium- and small-scale data interpolation, and spatial self-oblique variance interpolation has obvious speed advantages. The statistical information of various interpolation transformations is shown in Table4.Table 4
Data transformation statistics.
StatisticsNo transformationLogarithmic transformationPower exponent transformationMean3974.68.19639077240Std. dev.1543.20.43156752000Kurtosis2.64612.87414.0166Skewness0.4316−0.524031.1702Among them, kurtosis and skewness can reflect whether the data obey normal distribution. Kurtosis is used to describe the index of the height of data distribution, skewness is used to describe the symmetry of data distribution, and the kurtosis and skewness of normal distribution should be equal to 0.


## 4.6. Comparison Results of Various Interpolation Algorithms

Kriging interpolation method is to find the optimal, linear, and fair geometric features and spatial structure after considering the sample point distribution, density, size, and estimated spatial distribution position of the sample points. The normal distribution and quality of the data are the key to the success of the kriging interpolation method. Figure11 is a comparative analysis of the performance of the elevation interpolation results obtained by different interpolation algorithms on the same stratum surface.Figure 11
Performance comparison of elevation interpolation.Table5 is a comparison of the results of the interpolation of the three-dimensional attribute data volume for each difference algorithm.Table 5
Comparison of the results of the difference algorithm in the interpolation of the three-dimensional attribute data volume.
IDWKrigingSmooth krigingImprove algorithmMaximum error137.2583.2139.3743.23Minimum error1.343.013.8227.87Average error96.6262.4529.1638.98Among them, the sliding kriging algorithm is an interpolation method that finds the global variable range and performs the calculation of the local kriging equations according to the known variable range. The above data shows that the overall error of the improved interpolation algorithm is smaller than the inverse distance weighted interpolation. Compared with the inverse distance weighted interpolation method, it also estimates the unknown points by calculating the weighted average of the known points. The kriging method is based on the spatial variogram, and the variogram is determined by the existing spatial sampling points. It makes full use of the spatial structure and randomness of data points, and the interpolation effect is more realistic. The inverse distance weighting method is an interpolation method that constructs an approximation function to fit the original point set. This method has a certain error, and the interpolation result is too smooth due to the influence of the grid density, and some details are lost. Compared with the inverse distance weighting method and the ordinary kriging method, the improved interpolation algorithm in this paper improves the performance by 23.4%. And the average error value of the algorithm is 38.98, which has been greatly improved compared with other algorithms.


## 5. Discussion

This paper studies the regularized grid data obtained by interpolation to better reflect the continuous distribution of variables in space when forming curves, surfaces, and dividing blocks. The thesis first studied inverse distance weighted interpolation and kriging interpolation and compared and analyzed the experimental results. Taking into account the spatial structure of the sampled data, it conducts in-depth research on the kriging interpolation method based on geostatistics and improves a sliding neighborhood kriging method based on variable range. Experiments prove that the improved method can perform large-scale volume interpolation of three-dimensional attribute data. In actual work, the changes of regionalization variables are complex and often manifested as anisotropy (that is, the material properties and the distribution characteristics in space change with the change of direction, showing different properties in different directions). Generally, when kriging interpolation is performed, the anisotropy problem is solved by fitting the variogram in each direction and finding the fit structure. This article implements kriging interpolation based on the assumption of isotropy. In the future, consideration of anisotropy needs to be added to the realization of the kriging interpolation method.


## 6. Conclusions

In the text, based on BIM and data fusion technology, the 3D stratum spatial interpolation algorithm for subway tunnels is studied. For the proposed method, a lot of research work has been done on spatial interpolation and the visualization of regular data fields, and a sliding neighborhood kriging method based on variable range is improved. The improved kriging algorithm can interpolate the three-dimensional attribute data volume and the interpolation effect is good. On the basis of regular kriging interpolation, a variable amplitude sliding kriging interpolation method is researched and implemented. This method uses the variable range obtained from all sampling points as the side length to divide the study area into a square grid and solve the kriging equations in the local cells obtained by the division. Experiments show that the above algorithm has the problem of slow interpolation speed or unable to complete the interpolation operation during the interpolation process of the three-dimensional attribute data volume. This article has made improvements to it, and the improvements include the realization of the local variable range of the variogram and the flooding method to fill the neighborhood point set. The speed of the improved algorithm is significantly increased, and it can be applied to the interpolation of three-dimensional attribute data volume. However, there are still shortcomings in the research. Due to the uncertainty of space interference, different research topics and different interference methods will lead to different results. Different survey data, different regions, different time and space scales, and optimal interpolation methods are also different. There is no optimal interpolation method under absolute conditions. Due to the uncertainty of interference, there is no suitable interpolation method for any spatial value, only the optimal interpolation method under certain conditions.



---

*Source: 3359402-2022-03-28.xml*