Research on Tunnel Construction Monitoring Method Based on 3D Laser Scanning Technology

: A tunnel is a symmetric structure working under uncertain load. In the traditional method of monitoring and measurement of tunnel deformation, only a small number of symmetric measuring point displacements are engaged to represent the overall deformation of the tunnel, and it is difﬁcult to fully reﬂect the displacement changes of the tunnel during the construction process. In this paper, a deformation monitoring method is proposed for the tunnel construction process based on 3D laser scanning technology. Initially, 3D laser scanning was used to obtain scattered clouds of irregular structures considering the New Austrian Tunneling Method (NATM). Afterward, the modiﬁed B-spline interpolation and greedy triangulation were used to ﬁt the surfaces. Moreover, the normal vector matrix was innovatively applied to express the deformation of the tunnel ﬁtting surface, which solved the problem of scattered point clouds in 3D laser scanning. A normal vector was obtained by the intersection of the normal of one ﬁtted surface with the other. Subsequently, the maximum entropy method was introduced to form the probability density function of the normal vector, and the 1% and 5% probability eigenvalues were used to analyze the overall deformation trend of the tunnel. Finally, the eigenvalue of 5% probability, which was less affected by construction uncertainties, was selected for analysis. This way, the analysis and prediction method for overall tunnel deformation was established.


Introduction
Tunnels are symmetric structures working under uncertain loads. The concept of the New Austrian Tunneling Method (NATM) was proposed in the 1960s, and it is still considered to be a basic method for tunnel construction in weak and broken surrounding rock areas. Its core principle is to apply a thin layer of flexible shotcrete and bolt support close to the surrounding rock near the excavation surface to control the deformation and stress release of the surrounding rock. To implement the dynamic design of the tunnel, NATM requires vigilant tunnel monitoring and measurement, which are crucial tools for controlling deformation.
Monitoring and measurement of the stability and deformation of the surrounding rock, surface, and supporting structure are critical to collecting necessary data during tunnel construction. On-site measurements of the surrounding rock can be used to monitor the structural response, evaluate the structural safety in real-time, maintain structure safety and symmetry in nature [1], determine the timing of secondary lining, substantiate the rationality of the selected support method, and guide the design and construction of tunnels [2]. As shown in Figure 1, traditional monitoring and measurement use indices such as vault subsidence and convergence deformation [3] to characterize the deformation of the surrounding rock of the tunnel. Currently, monitoring and measurement are usually Three-dimensional laser scanning technology uses the principle of laser ranging to make a contactless measurement and can directly obtain massive 3D point clouds with irregular spatial distribution that can be used in reverse engineering [6]. Additionally, it has a certain degree of penetrability and is less affected by the environment. It has played an increasingly important role in global change, smart cities, resource surveys, environmental monitoring, basic mapping, and other pertinent fields [7][8][9].
Three-dimensional laser scanning technology has been gradually employed in tunnel construction owing to its fast detection speed, large quantity, high level of automation, strong environmental adaptability, and good adaptability to tunnel deformation measurement. Further, it allows for a comparison between the point cloud or the fitted point cloud surface with the design surface and helps in detecting the over-excavation and under-excavation of the tunnel [10]. Three-dimensional laser scanning technology has been widely used in the construction and operation of shield tunnels. The overall procedure of tunnel deformation using this technology includes the projection of the point cloud on the tunnel section, fitting the 3D laser scanner into the tunnel outline, calculating the vertical and horizontal diameters, and judging the tunnel deformation through the change in ellipticity [11]. The geometric information of the tunnel is contained in the detected point cloud during tunnel construction as per NATM, and this technology is also being adopted gradually for the detection of over-and under-excavation of tunnels. The tunnel point cloud is analyzed and studied using direct graphic comparison [12]. However, 3D laser scanning technology is rarely used in tunnel monitoring and measurement schemes. Therefore, this paper focuses on solving two major obstacles, i.e., the comparison of randomly distributed point clouds and the denoising of abnormal points. Additionally, considering that 3D laser scanning technology is an important means of reverse engineering, Three-dimensional laser scanning technology uses the principle of laser ranging to make a contactless measurement and can directly obtain massive 3D point clouds with irregular spatial distribution that can be used in reverse engineering [6]. Additionally, it has a certain degree of penetrability and is less affected by the environment. It has played an increasingly important role in global change, smart cities, resource surveys, environmental monitoring, basic mapping, and other pertinent fields [7][8][9].
Three-dimensional laser scanning technology has been gradually employed in tunnel construction owing to its fast detection speed, large quantity, high level of automation, strong environmental adaptability, and good adaptability to tunnel deformation measurement. Further, it allows for a comparison between the point cloud or the fitted point cloud surface with the design surface and helps in detecting the over-excavation and under-excavation of the tunnel [10]. Three-dimensional laser scanning technology has been widely used in the construction and operation of shield tunnels. The overall procedure of tunnel deformation using this technology includes the projection of the point cloud on the tunnel section, fitting the 3D laser scanner into the tunnel outline, calculating the vertical and horizontal diameters, and judging the tunnel deformation through the change in ellipticity [11]. The geometric information of the tunnel is contained in the detected point cloud during tunnel construction as per NATM, and this technology is also being adopted gradually for the detection of over-and under-excavation of tunnels. The tunnel point cloud is analyzed and studied using direct graphic comparison [12]. However, 3D laser scanning technology is rarely used in tunnel monitoring and measurement schemes. Therefore, this paper focuses on solving two major obstacles, i.e., the comparison of randomly distributed point clouds and the denoising of abnormal points. Additionally, considering that 3D laser scanning technology is an important means of reverse engineering, this paper employs the cloud platform to collect data and self-compiled software for analysis and processing.
In order to solve the analysis problem of the scattered point cloud, B-spline interpolation and a greedy triangulation were used to fit the tunnel surfaces. Afterward, the surfaces of the tunnel were compared before and after the measurement. This paper innovatively used the method of forming a normal vector matrix for comparison and introduced the Symmetry 2022, 14, 2065 3 of 17 maximum entropy method to create the probability density function of the normal vector and extract the deformation eigenvalues. The deformation eigenvalues were then used to create a curve that could assess and forecast tunnel deformation. The flowchart is shown in Figure 2. this paper employs the cloud platform to collect data and self-compiled software for analysis and processing.
In order to solve the analysis problem of the scattered point cloud, B-spline interpolation and a greedy triangulation were used to fit the tunnel surfaces. Afterward, the surfaces of the tunnel were compared before and after the measurement. This paper innovatively used the method of forming a normal vector matrix for comparison and introduced the maximum entropy method to create the probability density function of the normal vector and extract the deformation eigenvalues. The deformation eigenvalues were then used to create a curve that could assess and forecast tunnel deformation. The flowchart is shown in Figure 2.

Fitting Surface
Surface fitting is the basis for conducting scattered point cloud analysis, whereas scattered point clouds cannot be compared. After the curved surface is fitted, the comparison between the curved surfaces can be realized to understand the deformation of the tunnel.
Surface fitting of a point cloud is the process of fitting a discrete point cloud onto a surface for general analysis. The accuracy of the surface fitting has direct implications for the calculation outcomes [13]. The degree of fitting optimization can be effectively improved by selecting the right fitting method to obtain the reference surface for flatness calculation. The choice of an appropriate fitting model is crucial for surface fitting because

Fitting Surface
Surface fitting is the basis for conducting scattered point cloud analysis, whereas scattered point clouds cannot be compared. After the curved surface is fitted, the comparison between the curved surfaces can be realized to understand the deformation of the tunnel.
Surface fitting of a point cloud is the process of fitting a discrete point cloud onto a surface for general analysis. The accuracy of the surface fitting has direct implications for the calculation outcomes [13]. The degree of fitting optimization can be effectively improved by selecting the right fitting method to obtain the reference surface for flatness calculation. The choice of an appropriate fitting model is crucial for surface fitting because it has a direct impact on the precision of parameter solution and surface fitting in the subsequent data processing.

Bicubic Spline Interpolation Surface Construction
The cubic B-spline interpolation method is the process of obtaining a curve function group by mathematically solving the three bending moment equations through a smooth curve having a series of shape value points.
In the 2D plane curve interpolation, the cubic spline function on nodes x 0 , x 1 , ..., x n is defined as S(x) ∈ C2[a, b]; it is a cubic polynomial in each small interval [x j , x j+1 ], where a = x 0 < x 1 < · · · < x n = b defines subjected node. If the function value Y j = f (X j )(j = 0, 1, · · · , n) is given on the node x j , and the cubic spline function is established as S(x j ) = y j (j = 0, 1, · · · , n), then S(x) is referred to as the cubic spline interpolation function.
For m + n + 1 scattered points in the space, the surface is solved by the point cloud composed of scattered points; Q h (h = 0, 1, 2, . . . , m + n) represents the set of these scattered points. Thus the n-th spline curve can be defined as follows [14]: where Y k,n (t) is the basic function of the n-th degree spline function, and t represents the node parameter, i.e., t ∈ [0, 1]. In spline interpolation, when n = 3 and k = 0, 1, 2, 3, the polygon can be represented by 4 points, i.e., Q 0 , Q 1 , Q 2 , Q 3 . The cubic spline is a cubic polynomial defined as follows [15]: Thus, the expression of the cubic spline is: (6) where i = 0, 1, 2, . . . , m.
Following the principle of cubic spline interpolation in 2D space, a method is introduced for 3D space to obtain the interpolation surface of the bicubic spline. For this purpose, considering the n × m-th degree parametric surface as the calculation surface and (n + 1) × (m + 1) scattered points in the space, the multi-degree polynomial composed of i and j can be used for the solution, which can be solved in the form of a matrix. The space grid is formed by connecting two adjacent points in the point set with line segments, which is used as the control network and can be defined as follows.
where Q ij is the control vertex of the surface, u and v represent parameters in two directions, i.e., 0 ≤ u ≤ 1 and 0 ≤ v ≤ 1; Y in (u) and Y jm (v) represent the spline-based functions in u and v directions on the plane coordinates. The characteristic mesh of the bicubic B-spline surface is formed when m = n = 3, and the equation of the bicubic B-spline surface can be given as follows: Symmetry 2022, 14, 2065

Greedy Projection Triangulation Algorithm
The Delaunay triangulation method is tedious and difficult to converge. The greedy algorithm does not consider the overall situation, and it only seeks the local optimal solution [16][17][18]. Therefore, in the point cloud fitting, the greedy algorithm is considered a local fitting method. On the other hand, the number of point clouds is extremely large; therefore, the greedy algorithm eliminates the need to exhaust all possibilities to find the optimal solution and requires extensive computing time. For each point cloud subset, the problem is divided into smaller subsets, and the best fitting results from those subsets are used, since the greedy algorithm is capable of obtaining the local optimal result. Given that the correlation between the point clouds is not very high, the greedy algorithm is feasible for fitting the 3D laser scanning point clouds.
Greedy projection triangulation involves projecting the point cloud from the point cloud subset into a plane along the normal direction of the point cloud data and then triangulating to determine the topological relationship of each point using the Delaunay triangulation method [19,20]. By selecting a sample triangle as the initial surface, the boundary of the surface is continuously expanded and extended until all the points that are geometrically and topologically correct are connected. This way, a complete triangular mesh surface is formed. The specific steps are as follows [21]: The KD tree structure is used to determine the nearest points; this step is called nearest neighbor search.

2.
The plane and the surface are made roughly tangential to project the point cloud's normal information to the plane in blocks.

3.
The procedure is completed by projecting the points to the plane, triangulating the points on the plane, reflecting the topological relationship into the space, and concluding the continuous cycle of surface reconstruction.

1.
The point cloud and surface are compatible with the greedy triangulation method, which can accurately express the tunnel's local surface. When the density of the point cloud is high, the effect of the abrupt change in the normal angle can be overlooked because of the small change in the normal angle between adjacent triangles. A major surface fitting method for 3D laser scanning monitoring and measurement can be used in future research, given the high acceptance of this method in the engineering field.

2.
The B-spline interpolation method can restore the tunnel girdle area with minimal curvature by fitting the tunnel face with a polynomial in the whole area. At the same time, it maintains high goodness of fit. However, the fitting surface of the vault (area with large curvature) has a distorted position, and the partition fitting method is used to improve the fitting degree. The normal angle is continuously changing; the fitting surface is smooth and continuous at high orders. However, the fitting surface does not tend to overlap with the point cloud and has a tendency to disregard the local variations of the tunnel. Therefore, engineers and technicians have a low level of acceptance for this method; thus, it can be used as an auxiliary fitting method in future studies and compared to the greedy triangulation method for validation.

Deformation Hypothesis
A full-fitting surface is created after tunnel fitting. However, since point clouds cannot be tracked by 3D laser scanners, the authors were unable to determine the direction in which each point in the tunnel was moved. The initial support surface of the tunnel is also quite complicated, which seriously affects the monitoring and measurement of the tunnel. The tunnel comprises the surrounding rock and initial support to form an arched structure. The arched structure is a compression bending member or an axial compression member. Compressive instability failure features of the tunnel arch ring under pressure are seen in the case of the original support of the tunnel being destroyed. Therefore, it is assumed that the deformation direction of the initial support surface of the tunnel is in the normal direction of the fitted surface. According to the traditional monitoring and measurement method, the monitoring of vault settlement and peripheral convergence also aligns with the assumption that the tunnel is deformed along the normal direction. In order to determine the process of tunnel deformation, make predictions of deformation and realize early warnings of tunnel arch bridge instability, this article compares and analyzes the two-point cloud fitting surfaces before and after using the normal vector.

Observation Point Arrangement of the Fitting Surface
Cloud data of 3D laser scanning are discrete and randomly distributed and require a continuous function model [22]. Therefore, the point cloud data of the initial support of the tunnel detected each time are fitted to obtain the fitting surface function of the point cloud, i.e., S = f (x, y, z). In order to simplify the calculation, the point cloud data are selected at preset time intervals, and the first, second and nth data fitting surfaces are set as According to the same circumferential angle and longitudinal distance, the points are arranged on the curved surface of the previous tunnel fitting as the monitoring points for monitoring and measurement. This ensures the uniform arrangement of the tunnel surface points on the plane composed of the tunnel central axis and the tunnel vault axis. The observation points needed for monitoring and measuring are uniformly distributed and are assigned to matrix A in the point set, as follows.

Calculation of Initial Support Deformation of Tunnel
Starting from the observation point matrix B and assuming that the initial support of the tunnel only creates displacement in the normal direction, the normal of the first data fitting surface, S 1 = f 1 (x, y, z), is determined. Based on the partial derivative of the binary function, the partial derivative at the point (x 0 , y 0 ) can be given as follows: Further, the first data fitting surface S 1 = f 1 (x, y, z) is an implicit function, i.e., F(x, y, z), which implies that ∂z ∂x ∂z ∂y Thus, the normal equation for the observation point a i,j (x 0 , y 0 ) on the first data fitting surface S 1 = f 1 (x, y, z) can be given as follows: The observation point is based on the normal f i,j of the first data fitting surface; thus, the normal matrix, F, is the absolute value matrix of the normal vector starting from the point a i,j as follows: Moreover, matrix B is defined as the set of intersections between the normal and the two curved S n = f n (x, y, z) surfaces in the following equation; the distance from A to B can be roughly calculated as the deformation of the initial support of the tunnel.
As shown in Figure 3, the tunnel deformation length is the distance between the normal line and the intersection points of S 1 = f 1 (x, y, z) and S n = f n (x, y, z) surfaces. The displacement of S n is set to be positive when S n is inside S 1 , and negative when S n is outside S 1 . The normal vector matrix D is formed by summarizing the observation points of the displacement.
As shown in Figure 3, the tunnel deformation length is the distance between the normal line and the intersection points of 1 1 , , S f x y z  ( ) and n n , , S f x y z  ( ) surfaces. The displacement of Sn is set to be positive when Sn is inside S1, and negative when Sn is outside S1. The normal vector matrix D is formed by summarizing the observation points of the displacement.

Maximum Entropy Analysis of Tunnel Deformation Maximum Entropy Analysis
Entropy S(x) is used to quantitatively describe the uncertainty or information content of random events or variables [23,24]. The principle of maximum entropy [9] refers to that among all the many probability density functions that satisfy the given constraints; the probability density function with the largest information entropy is considered to be the best (i.e., has the smallest deviation) [25][26][27].
where p(x i ) is the probability against the random variable value x i , and N is the number of samples. When x is a continuous random variable, entropy can be defined by the following Equation: where f (x) is the probability density function of the random variable distribution.
The following conditions are satisfied considering the normal vector of each observation point of the tunnel ∆ as a continuous random variable and the probability density function f (∆) as a function of it.
where n i is the i-th order origin moment of the statistical sample of ∆, which can be determined by the calculation of the statistical sample.
In order to realize the entropy of the ∆(i.e., Y(∆)) for the maximum value under the conditions of Equations (23) and (24), the following equation can be specially derived, which is referred to as Lagrange's Equations (26) and (27): where β 0 , β 1 , β 2 , · · · , β N is the Lagrange multiplier series. When

∂F(∆)
∂ f (∆) = 0, the entropy Y(∆) reaches the maximum value, which implies that the following condition is satisfied: Thus, the probability density function of the normal vector ∆ of the tunnel observation point represented by the maximum entropy theory can be obtained as follows: The undetermined coefficients β 0 , β 1 , β 2 , · · · , β N can be solved by the i+1-element i-th degree equation system using Equations (23) and (24). The probability density function of the normal vector of the tunnel observation point is shown in Figure 4. To improve the convergence of the calculation, this paper adopts the fourth-order moment of the maximum entropy method (i.e., N = 4) [28][29][30][31][32]. This can yield the peak value and position of the This paper uses the ±1% probability feature deformation value in the probability density function of the initial support normal vector of the tunnel as the characteristic deformation value of the tunnel monitoring and measurement analysis. In addition, it takes the ±5% probability feature deformation values for the auxiliary analysis [33][34][35]. Employing the maximum entropy approach's fourth moment may fully extract sample information generated using this method. Thus, Equations (23) and (24) can yield the following results:

Example Verification
In this paper, two sets of point cloud data measured in a railway tunnel were fitted, and the normal vector matrix of the fitted surface was obtained by the method discussed in Section 3. This paper used the maximum entropy method to analyze and process the normal vector probability density function. The fitting methods involved greedy triangulation and B-spline interpolation approaches.

3D Laser Scanning Inspection and Deterministic Analysis of Traditional Monitoring Quantity
Considering that the spline interpolation method surface is relatively continuous and stable, and the influence of accidental deformation such as debris interference and blasting damage in the tunnel can be excluded [36,37], the deformation direction is compared using the traditional monitoring and measurement method and the B-spline interpolation method. The monitoring and measuring parts of two traditional techniques in the laser scanning interval, No. 01 and No. 02, were chosen for this study. Table 1 shows differences in the deformation direction of the vault area during the same period using the traditional monitoring method. This phenomenon is also established in the 3D laser scanning detection results. Figure 5 presents the 3D laser scanning inspection, from which it can be interpreted that the deformations of the vault and the arch foot are more complicated, and the positive and negative deformations appear alternately, which is also substantiated by the deformation data of the 14th-15th days in Table  1. The comparison demonstrates that both approaches reflect the same tunnel deformation. On the other hand, the traditional monitoring method employs section This paper uses the ±1% probability feature deformation value in the probability density function of the initial support normal vector of the tunnel as the characteristic deformation value of the tunnel monitoring and measurement analysis. In addition, it takes the ±5% probability feature deformation values for the auxiliary analysis [33][34][35]. Employing the maximum entropy approach's fourth moment may fully extract sample information generated using this method. Thus, Equations (23) and (24) can yield the following results: where ∆ 1 and ∆ 2 can be used as the characteristic deformation value of tunnel monitoring and measurement. Considering that the deformation values are mostly concentrated in the peak probability density function, the deformation values ∆ 1max and ∆ 2max can be used as reference indicators for tunnel monitoring and measurement.

Example Verification
In this paper, two sets of point cloud data measured in a railway tunnel were fitted, and the normal vector matrix of the fitted surface was obtained by the method discussed in Section 3. This paper used the maximum entropy method to analyze and process the normal vector probability density function. The fitting methods involved greedy triangulation and B-spline interpolation approaches.

3D Laser Scanning Inspection and Deterministic Analysis of Traditional Monitoring Quantity
Considering that the spline interpolation method surface is relatively continuous and stable, and the influence of accidental deformation such as debris interference and blasting damage in the tunnel can be excluded [36,37], the deformation direction is compared using the traditional monitoring and measurement method and the B-spline interpolation method. The monitoring and measuring parts of two traditional techniques in the laser scanning interval, No. 01 and No. 02, were chosen for this study. Table 1 shows differences in the deformation direction of the vault area during the same period using the traditional monitoring method. This phenomenon is also established in the 3D laser scanning detection results. Figure 5 presents the 3D laser scanning inspection, from which it can be interpreted that the deformations of the vault and the arch foot are more complicated, and the positive and negative deformations appear alternately, which is also substantiated by the deformation data of the 14th-15th days in Table 1. The comparison demonstrates that both approaches reflect the same tunnel deformation. On the other hand, the traditional monitoring method employs section deformation judgment, which is incomplete. Three dimensional laser scanning can thoroughly exhibit the deformation, which is more suitable for correct judgment. deformation judgment, which is incomplete. Three dimensional laser scanning can thoroughly exhibit the deformation, which is more suitable for correct judgment.  Construction records indicate that local blasting was carried out at the under-excavated position on the left side of the tunnel on the 15th-17th days. It can be seen from Table 1 that the overall displacement of the tunnel to the left occurred due to the local blasting. In Figure 6, the yellow area of the tunnel indicates the deformation area to the outside of the tunnel, and the red area indicates the deformation to the inside of the tunnel. During this period, there was a noticeable displacement to the left, and it can be recognized qualitatively that the displacement direction detected by traditional monitoring was fundamentally consistent with the deformation direction acquired by 3D laser scanning detection analysis. Construction records indicate that local blasting was carried out at the under-excavated position on the left side of the tunnel on the 15th-17th days. It can be seen from Table 1 that the overall displacement of the tunnel to the left occurred due to the local blasting. In Figure 6, the yellow area of the tunnel indicates the deformation area to the outside of the tunnel, and the red area indicates the deformation to the inside of the tunnel. During this period, there was a noticeable displacement to the left, and it can be recognized qualitatively that the displacement direction detected by traditional monitoring was fundamentally consistent with the deformation direction acquired by 3D laser scanning detection analysis.  Table 2 compares the characteristic deformation values of the normal vector probability density function of 1% and 5% probability for the greedy triangulation and B-spline interpolation methods. It can be observed that the positive deformation (deformation to the outside of the tunnel) obtained by the B-spline interpolation method was significantly larger on the 15th-17th days and significantly smaller on the 20th−21st days. However, when the data were sufficiently concentrated, the computation did not converge on the 22nd-23rd days while utilizing the maximum entropy approach for data analysis. Moreover, observations could not be made by 3D laser scanning on the 16th day and traditional monitoring on the 17th-20th days due to construction interference in the tunnel [38]. Figures 7 and 8 compare the trend of the characteristic deformation value of the 1% probability feature deformation value of the 3D laser scanning with the deformation Figure 6. Schematic diagram of the deformation direction of the surface fitted by the B-spline interpolation method on the 15th-17th day (Note: the deformation direction of the golden part is outward, and the deformation direction of the red part is inward). Table 2 compares the characteristic deformation values of the normal vector probability density function of 1% and 5% probability for the greedy triangulation and B-spline interpolation methods. It can be observed that the positive deformation (deformation to the outside of the tunnel) obtained by the B-spline interpolation method was significantly larger on the 15th-17th days and significantly smaller on the 20th−21st days. However, when the data were sufficiently concentrated, the computation did not converge on the 22nd-23rd days while utilizing the maximum entropy approach for data analysis. Moreover, observations could not be made by 3D laser scanning on the 16th day and traditional monitoring on the 17th-20th days due to construction interference in the tunnel [38]. Figures 7 and 8 compare the trend of the characteristic deformation value of the 1% probability feature deformation value of the 3D laser scanning with the deformation value measured by traditional monitoring. The 3D laser scanning point cloud adapted well to the greedy triangulation and B-spline interpolation methods for surface fitting.   It can be seen from Figures 7 and 8 that the overall observation trend of the traditional monitoring method and 3D laser scanning were consistent with each other. However, 3D laser scanning yielded more conservative results than the traditional monitoring method.    It can be seen from Figures 7 and 8 that the overall observation trend of the traditional monitoring method and 3D laser scanning were consistent with each other. However, 3D laser scanning yielded more conservative results than the traditional monitoring method. It can be seen from Figures 7 and 8 that the overall observation trend of the traditional monitoring method and 3D laser scanning were consistent with each other. However, 3D laser scanning yielded more conservative results than the traditional monitoring method. Given the presence of steel arches in the tunnel's first support, the surrounding rock deformation was not necessarily coordinated with the deformation of the initial support. In addition, it can be observed that the analysis results of the surface fitted by the greedy triangulation method and the surface fitted by the B-spline interpolation method were consistent with marginal deviations at some points.

Exceedance Probability Characteristic Deformation Value Analysis
In this paper, the gradient change of probability feature deformation value, K n , was used to validate the results obtained by the greedy triangulation and the B-spline interpolation methods, as follows.
where S n is the probability feature deformation value of the specific transcendence probability on an nth day. Table 3 shows that the trend of K n of the 1% probability feature deformation value is approximately in line with that of the 5% probability feature deformation value. The result of the B-spline interpolation method showed significant inconsistency, and the greedy triangulation method's judgment results were more stable and suited for assessing the initial support deformation in engineering.  10 show a symmetrical connection between positive and negative deformations. Further, the 5% probability feature deformation value curve was compatible with the 1% probability feature deformation value curve, whose fluctuation of deformation value was substantially more steady. Similarly, the trends of observations obtained from the greedy triangulation and B-spline interpolation methods were comparable. The Bspline interpolation underwent a mutation for a few days; however, after the mutation, the deformation became consistent with the greedy triangulation method. Compare dates Greedy Trigonometry 1% probability feature deformation value B-spline interpolation 1% probability feature deformation value Greedy Trigonometry 5% probability feature deformation value B-spline interpolation 5% probability feature deformation value Compare dates exceeding probability rate of change(%) Greedy Trigonometry 1% probability feature deformation value B-spline interpolation 1% probability feature deformation value Greedy Trigonometry 5% probability feature deformation value B-spline interpolation 5% probability feature deformation value Figure 10. Comparison of negative feature deformation between greedy triangulation and B-spline interpolation.

Discussion
Considering that the 3D laser scanning point cloud is in the scattered form, a surface fitting is required. Compared to the B-spline interpolation method, the results of the greedy triangulation method coincided well with the fitting surface, showing better Compare dates Greedy Trigonometry 1% probability feature deformation value B-spline interpolation 1% probability feature deformation value Greedy Trigonometry 5% probability feature deformation value B-spline interpolation 5% probability feature deformation value Compare dates exceeding probability rate of change(%) Greedy Trigonometry 1% probability feature deformation value B-spline interpolation 1% probability feature deformation value Greedy Trigonometry 5% probability feature deformation value B-spline interpolation 5% probability feature deformation value Figure 10. Comparison of negative feature deformation between greedy triangulation and B-spline interpolation.

Discussion
Considering that the 3D laser scanning point cloud is in the scattered form, a surface fitting is required. Compared to the B-spline interpolation method, the results of the greedy triangulation method coincided well with the fitting surface, showing better

Discussion
Considering that the 3D laser scanning point cloud is in the scattered form, a surface fitting is required. Compared to the B-spline interpolation method, the results of the greedy triangulation method coincided well with the fitting surface, showing better goodness of fit.
Further, when the density of the point cloud was high enough, the curvature of the edges of the triangles in the fitted surface was reduced. When the abrupt change in curvature is minimal, it is possible to assume that the tunnel displacement is in the normal direction of the surface.
Moreover, the normal vector was found to be perpendicular to the first fit. The displacement of the tunnel at this location (which was a vector itself) was the length of the line segment at the point of intersection of the normal vector and two fitted surfaces. A normal matrix can be formed by uniformly arranging normals on a surface. The normal vector matrix accurately expressed the deformation of each position of the tunnel. It is reasonable to assume that the tunnel surface deforms in the normal direction mostly. However, this assumption may not stand correct when the tunnel flatness is seriously insufficient; thus, it is necessary to ensure that the tunnel is flat enough before monitoring.
In this paper, the maximum entropy method was employed to form the probability density function of the normal vector. The authors reflected the deformation trend of the tunnel through the changing trend of the eigenvalues of the normal vector probability density function. This method can not only eliminate the abnormal points caused by the insufficient flatness of the tunnel surface but also fully extract reliable information about the tunnel deformation. The entire deformation of the tunnel may be accurately identified in space. The research also discovered that when the tunnel deformation is small, the normal vector is small, and the probability density function calculation of the maximum entropy method does not converge. If this continues, we may say that the tunnel deformation has stabilized. It was observed in the current study that the 1% probability feature deformation value and the 5% probability feature deformation value could make an accurate judgment. It may be deduced that when the exceedance probability is too small, it is more likely to be disturbed by construction factors. Furthermore, the 5% probability feature deformation value may more accurately depict the deformation of each section of the tunnel.
In future studies, different methods can be used to form the probability density function, and the change in the normal vector matrix can be identified through AI algorithms to judge the deformation trend of the evaluated tunnel.

Conclusions
Using the greedy triangulation method to fit the tunnel point cloud, the tunnel fitting surface completely coincided with the fitting surface. The normal vector matrix was created by making the normal of the fitted surface consistent and intersecting with another fitted surface, assuming that the tunnel surface was deformed towards the normal. In this paper, the normal vector matrix was used to characterize the tunnel deformation for the first time, and the probability density function of the normal quantity was formed by the reliability theory of the maximum entropy method. It is proposed to use the probability characteristic deformation value of the normal vector to analyze the trend of tunnel deformation and realize tunnel deformation prediction.

Data Availability Statement:
The data used to support the findings of this study are included in the article.