A 3D Image Reconstruction Model for Long Tunnel Geological Estimation

Long tunnels often collapse during the construction period. To ensure personnel safety, the geological characteristics must be predicted before tunnel face excavation. In this study, the ground-penetrating radar (GPR) technique is introduced to obtain information regarding the tunnel excavation face at a certain interval. 0e amplitude of the radar echo signal is expressed as a function of the position and travel time. A B-scan strategy is selected for the GPR to obtain tunnel information. A frequencydomain (w-k) focusing algorithm, namely, a synthetic aperture radar focusing algorithm, is applied to focus scattered radar signals to obtain focused images. A low-pass filter is designed to remove noises from the original signals.0e contours of target objects are extracted from the background information using the edge detection technique. Space coordinate values of the objects are converted to polar coordinates using the Hough transform algorithm for 3D modeling. Visual C++ and AutoCAD are combined to develop a 3D CADmodel to help managers in controlling the construction process. 0e system creates 3D visualization model images and evaluates the geological characteristics of the tunnel excavation faces.0e Taigu Tunnel located in the Shanxi Province of China is taken as a case study. A procedure for the geological analysis of this tunnel is introduced in detail, and a 3D image model is built. 0e results show that the 3D model can help predict rock compositions and locate potential hazards. Moreover, it has better accuracy than conventional models and can be applied to similar transportation construction projects.


Introduction
China has the largest and widest network of traffic tunnels in the world. By the end of 2018, the number of road tunnels alone is expected to reach 17,738/17.236 million linear meters, of which 1,058 will be extralong tunnels per 4.7065 million linear meters and 4,315 will be long tunnels per 7.421 million linear meters [1,2]. e construction of such tunnels is associated with significant construction difficulties and risks, including high geostresses and many geological problems such as high geotemperature, high hydraulic pressure, special adverse strata, high-pressure water-rich karsts, large deformation of soft rock, and large deformation of the Archaean metamorphic hard rock [3,4]. Tunnel collapse, roof caving, water inrush, and/or other accidents increase the construction difficulty and construction cost; the safety of personnel is another major concern [5][6][7]. Visualizing a geological structure through 3D modeling is considered an effective method for predicting rock composition and determining potential hazards.
Geological structure prediction is an essential routine work in the excavation of tunnels, providing important prior information that can help ensure a safe, economical, and efficient tunneling [8][9][10][11]. Several methodologies have been developed for the geological prediction of tunnel rocks, such as surface investigation, geological analysis, geophysical prospecting, and rock mapping. ese efficient and unbiased methods for hazardous inspection can be grouped into two categories: destructive and nondestructive methods [12]. As early as the 1950s, geological structures were being analyzed by directly drilling the tunnel faces [13]. is destructive method can help access the rock samples ahead of the tunnels and the surrounding rocks. However, it is timeconsuming and costly. Nondestructive methods for geological prediction during tunnel construction saw their emergence in Europe in the 1980s and America and China in the 1990s. Unlike destructive methods, nondestructive methods are designed to provide information about the geological properties without deteriorating the material microstructure, making them suitable for the geological forecasting of tunnel faces. Nowadays, geological prospecting has been developed, ranging from destructive methods such as drilling and pilot tunnel excavation to nondestructive methods such as mechanical, electromagnetic, and electrical methods [14]. Each approach has its advantages and disadvantages. Zhang et al. recommended the seismic reflection method, ground-penetrating radar (GPR), infrared water detection, a transient electromagnetic method, and advanced probe borehole drilling to predict the location, size, and distribution of karst features ahead of tunnel faces [15]. Yuan et al. conducted a comprehensive investigation through geophysical exploration, deep hole drilling, and ground stress testing to determine the geological and hydrogeological conditions of a tunnel [16]. e mechanical methods include negative apparent velocity, seismic, ultrasonic, and acoustic emission methods. Wang et al. detected hidden structures in a tunnel surrounded by rocks using the pulse reflection method [17]. Jiang et al. applied a joint seismic travel time and waveform inversion method to image a buried tunnel with concrete walls and a void space inside [18]. Zhang et al. proposed an improved shallow soil-rock interface seismic reflection method for a geotechnical field investigation [19]. However, this seismic method is not ideal as it is difficult to accurately determine the seismic velocity of the rocks around the tunnel walls. ese methods may yield inaccurate geological inferences, which will hamper the classification of the wall-rock characteristics. Considering the drawbacks of the seismic method for hydrogeological surveys, Liu et al. used the seismic and resistivity methods to predict tunnel collapse and water inrush due to fracture zones and water-bearing structures [20]. Pazzi et al. conducted a comprehensive geophysical survey using resistivity tomography, microgravity method, and single-station seismic noise measurement to explore the cause of a massive sinkhole called Tiankeng [21]. Computer techniques have been successfully implemented in the construction industry for decision-making and monitoring.
Electrical methods include electromagnetic (EM) and electrical resistivity (ER) investigations. e inspection requires special measuring devices and complex investigations of the information contained in the measured signals. e receiver records the signal strength at certain points on the near-surface to depths greater than 300 ft. e ratio of the voltage to the current is the resistance. When the measurement is made over a homogeneous surface, the apparent resistivity is equal to the true resistivity of the ground. Resistivity imaging is a technique used to detect and characterize water bodies in open fractures or faults based on the resistivity contrast between the water bodies and surrounding materials [22]. Ismail established a 3D filtered tunnel model of the weathering grade, quality grade, and heat capacity of rock mass based on surface survey data and electrical resistivity tomography technology [23]. Recently, the GPR, which is an electromagnetic method, has been considered the most useful technique to search for metals and minerals at relatively shallow depths (approximately 15-30 m). A GPR comprises a radar antenna that transmits and receives wide-frequency electromagnetic energy impulses at frequencies between 25 MHz and 1 GHz. e different dielectric constants of the medium will produce electromagnetic wave reflection and the amplitude of the received electromagnetic waves. e frequency characteristics, subsurface structure, and stratigraphic characteristics can be identified using geophysical technology with high observational accuracy. is method has been successfully applied to map subsurface geological structures and groundwater contaminants. In 1988, Daniels et al. published a tutorial paper entitled "Introduction to subsurface radar," giving a good overview of the GPR technology. e GPR is widely used in traffic infrastructure and tunnel detection because it can accurately and continuously detect the contours of road surfaces [24]. Li et al. used GPR to measure tunnel lining thickness [25]. Li et al. comprehensively applied GPR and Geo-3d technology to survey the geometric characteristics of karst caves and their spatial relationship with tunnels [26]. Alani et al. and Zhang et al. studied the relationship between different detection frequencies of GPR and imaging resolution and detection depth [27,28]. Wei et al. found that integrating the data obtained by multiple GPRs and seismic sensors can help recover the 3D position of possible events in front of the mining face [29].
is study is focused on controlling geological hazards in advance during the evacuation period using nondestructive methods.
e GPR technology is used to capture information regarding the rock structures of excavation faces. To extract the features of the tunnel faces, radar reflection is transformed into digital images using synthetic aperture radar (SAR) algorithm. e rock conditions can be easily interpreted and conceived in a 3D environment. Typically, two-dimensional maps are applied to show the locations and directions of the weak planes, and the excavation face conditions are determined by manual visual investigation for most tunnels. It is a costly and time-consuming process and even extremely dangerous. In this regard, Zhao et al. conducted a 3D simulation on the complex failure evolution of a brittle fault zone and quantitatively evaluated the potential failure mode [30]. Liu et al. established an overall high-precision 3D model by comprehensively considering the topography, folds, and major faults and revealed the zoning phenomenon of the in situ stress along the line [31]. Vanneschi et al. established a 3D model of a geological structure by combining 3D laser scanning technology and photographic image technology [32]. Xiong et al. proposed a 3D multiscale geological modeling method for the risk assessment in tunnel engineering by coupling the surrounding rock geological body model and a tunnel model [33]. With the development of information management systems, the use of 3D CAD models has become increasingly common in construction projects. A 3D CAD model can help visualize a project in 3D, and tunnel-related characteristics can be recognized using a computer system, thus enabling the efficient control of the construction process. Tonini et al. proposed a geological modeling method based on field survey data and reconstruction of a geological section perpendicular to a geological structure (fold) axis (geological section) using a CAD software system [34]. um proposed a method based on a digital geographic information system (GIS) to conduct a geological survey and realized the automation of geological and geomechanical mapping in the process of tunnel excavation [35]. e rest of this paper is organized as follows. First, the tunnel face detection technique GPR is briefly introduced. Subsequently, the B-scan strategy for tunnels is reviewed. In the third part, we discuss the SAR algorithm and Hough transform for 3D image model construction. Computerized software for 3D tunnel modeling is illustrated in the fourth part. In the discussion part, the Taigu Tunnel project undertaken by China is taken as a case study to demonstrate the efficiency of the software. Finally, the conclusions are summarized.

GPR Theory
e precision and reliability of conventional detection methods were improved in this study. Since tunnel rocks have varied shape, size, and depth and are filled with a variety of materials, the selection of the geophysical technique for the investigation depends on many variables.
Underground excavation is associated with many problems related to safety. Since complete records of tunnel construction are unavailable, common radar images are a valuable source to interpret geological information [36]. In this study, the GPR was selected to inspect the different layers of a tunnel; this method uses radio waves to probe the natural geological material. In a typical GPR measurement process, a transmitter and a receiver are deployed in a fixed geometry and are moved over the surface to detect any reflections from subsurface features, as shown in Figure 1.
e GPR antenna moves horizontally across the surface of the ground. e transmitting antenna radiates short high-frequency electromagnetic pulses into the ground, where they are refracted, diffracted, and reflected as they encounter variations in the dielectric permittivity and electric conductivity. e most common operation mode of the GPR for detecting tunnel face geological structures is the reflection mode, whereby the traces of returned waves are collected either continuously or in stations along a line, thus creating a time cross section or a profile image of the subsurface. e velocity of the waves varies with the physical and chemical properties of the material through which they travel.
Signals transmit electromagnetic pulses from the surface antennas into the ground. e time elapsed between sending the pulses and receiving them back at the surface is measured [37]. When the travel times of the energy pulses are measured and their velocity through the ground is known, the distance can be accurately determined. e radar travel time is calculated using the following equation: where z is the distance to the reflectors, in m; t is the two-way travel time, in ns; x is the distance between the antenna and receiver; and v is the velocity of the electromagnetic waves in the medium, which is calculated using the following equation: where c is the propagation velocity of the electromagnetic waves in vacuum (0.29979 m/ns); ε r is the relative permittivity of the subsurface material; and μ r is the relative permeability of the magnetic medium. e electromagnetic wave velocity is equal to the product of the wavelength λ and the frequency f as expressed in the following equation: e transmission frequency f in a GPR survey ranges from 25 to 100 MHz.
When the electromagnetic wave propagation changes due to the relative permittivity under geological conditions, the electromagnetic waves will produce reflection and transmission, and their energy allocation is closely related to the abnormal changes in the electromagnetic wave reflection coefficient. e reflection coefficient is calculated using the following equation: where y is the reflection coefficient of the electromagnetic waves at the interface and ε 1 and ε 2 are the relative dielectric constants of the first and second layers, respectively.

GPR Scanning Strategy
e radar feedback signal is processed to form an amplitude image of the tunnel faces ahead via different scanning strategies. e algorithms can be divided into A-scan, B-scan, and C-scan depending on the scanning geometries and directions. A-scan screening involves recording two parameters: amplitude and pulse transit time. e transit time of the radar echoes is used to fix the distance of the Journal of Advanced Transportation objects and tells the observer of certain interfaces along a line by one dimension. However, 1D image information is not useful for tunnel face evaluation. Many researchers prefer at least a 2D image than a simple one-dimensional image because it can be used to take an image of a cross section through the object.
Hence, we employed the B-scan GPR algorithm to show weak planes. e scattered data were recorded from several spatial positions on a straight scan by moving antennas on the tunnel surfaces and looking downward. As shown in Figure 2, the antenna moves along the x direction, and a series of pulse-echo measurements are acquired for all (x, y) positions over a plane, where z is the space depth. e individual reflected waves received are digitized into a reflection trace, and when many traces are stacked next to each other, a 2D vertical profile can be produced along the transect. e set of A-scans is gathered to form a 2D dataset b(x, t), called the B-Scan, as shown in Figure 3.
In the B-scan, the vertical axis represents the screening time and the horizontal axis represents the radar scan traces [38]. During this process, a transmission and reception phenomenon occurs at each observation position, and then, the recorded data are combined to produce a 2D space-time GPR image.

GPR Focusing Algorithms
Because of the beam width of the transmitting and receiving antennas, the reflections on a structure spread over a broad region in the recorded data. Many focusing techniques have been developed to correct defocusing images.
e SAR algorithm is applied to original GPR B-scan images to obtain true information about the location, shape, and size of targets. Furthermore, the SAR focusing algorithm is efficient and does not introduce additional losses to the digital signal acquisition system. Since basic B-scan signals in the time domain cannot handle waves of varying speed, the frequency-domain w-k method based on the SAR is introduced in this study.
In this study, a 2D Fourier transform is applied to interpolate the signal in the Fourier space. As shown in Figure 3, the antenna moves along the X axis. e magnitude and phase of the EM scatterings are recorded at each discrete point. Suppose the target point P is located at (xi and zi), the recorded data Es(x, ω) in the frequency domain can be written as follows [39][40][41]: where ω � 2πf is taken for convenience of nomenclature, ρ is the complex amplitude of the scattered electric field of the point target, and v is the velocity of the electromagnetic waves in the medium, as expressed in equation (2). By applying 1D Fourier transform to Es(x, ω) from space x to spatial frequency k x , we can obtain the w-k domain representation of the scattered data using the following equation: Two-dimensional Es ′ (k x , k z ) data can be obtained by interpolating Es(k x , ω) onto a (k x , k z ) mesh. ereafter, the final focused GPR image can be obtained by transforming the data from the (k x , k z ) domain to the (x, z) domain via the following equation: where e(x, z) is the focused GPR image in the spatial-spatial domain after the w-k domain SAR algorithm.

Signal Processing.
Synthetic aperture radar images contain speckle noise that degrades the image quality, and the noise can significantly impact the performance of the radar in judging the target [42]. Although weak planes can be detected and located by the SAR algorithm, the exact material depth in the layers should be determined using the signal process techniques. In this study, a low-pass filter was used to detect images and suppress noise. All the noise signals are considered high-frequency noise and have a smoothing effect for images. It is better to perform contour detection automatically after the low-pass filter process. Edges define the boundaries between the regions in an image, which helps with segmentation, and the object contour position of each geological layer is calculated using an edge detection recognition algorithm. To evaluate the number of fringes quickly and exactly, an extreme point tracking method is proposed in this study. A 3 × 3 window is applied to the images with a certain threshold. e threshold is recalculated as the mean of nine intensity values of the pixels in the window. If a pixel has an intensity value greater than this threshold, it is set to 1; otherwise, it is set to 0. e original image is converted to a binary pattern with a 3 × 3 window. Sixteen possible edgelike patterns arise in 3 × 3 windows, as shown in Figure 4.
If the binary window obtained above matches any of the sixteen masks, the center pixel of the window is set to be an edge pixel. e previous steps are then repeated for each pixel in the image as the center pixel of the window, and all edges can be determined. To remove false edges due to noise, the variance for each 3 × 3 window is calculated and then compared with a global threshold based on the level of noise in the image. If the value is greater than the threshold, it is retained as an edge. If it is not greater than the threshold, it is removed. Figure 5 shows the schematic diagram of the algorithm.

ree-Dimensional CAD Model.
With the development of computer technology, 3D CAD models are being applied to visualize processes during the construction phase, instead of 2D images. ree-dimensional tunnel modeling not only involves geometrical data reconstruction but also texture extraction. e Hough transform technique is used to calculate the tunnel position parameters in space, which can be used to isolate the features of a particular shape within an image. To simplify the 3D visualizing process, all the geological layers are considered regular shapes. Hence, the classical Hough transform is applied to 3D tunnel modeling since it is concerned with the identification of characteristics in the image.
We take the simplest linear transform as an example; a straight line can be expressed by y � ax + b in image space and can be graphically plotted for each pair of image points (x, y). In the Hough transform, a straight line is considered a function of the slope a and the intercept parameter b. Hence, the straight line can be represented as a point (a, b). For computational reasons, it is better to use a different pair of parameters, namely, r and θ, for the lines in the Hough transform, as shown in Figure 6. e parameter r represents the distance between the line and the origin, and θ is the angle of the vector from the origin to this closest point. Using this parameterization, we can express the equation of a line as follows: where r ∈ R, θ ∈ [0, 2π], r ≥ 0, and for straight lines in two dimensions, we can rewrite the equation as r(θ) � x 0 cos θ + y 0 sin θ. Hence, it is possible to associate with each line of a pair.

Journal of Advanced Transportation
To obtain a 3D model of a tunnel, the Hough transform concept is extended to 3D space. Suppose any plane in the X-Y-Z space corresponds to a ρ-θ-φ space and considering a particular point A (x a, y a , z a ) in the X-Y-Z space, planes passing through this point can be expressed as in the following equation: ρ � x a sin θ cos ϕ + y a sin θ sin ϕ + z a cos θ, where (ρ, θ, φ) is a vector from the origin to the nearest point on the plane. In the study, all the contour points of the SAR images in the ρ-θ-φ space were converted to the spatial domain using the Hough transform technique. AutoCAD and Visual C++ software package were used to model the 3D tunnel image. Tunnel shaped data are obtained by on-site measurement at each cross-section along the X axis. e appropriate values of the length are determined based on the evacuation process. ree-dimensional models are generated from the cross-section of the tunnels obtained using the software. e software maps the entire model of the images and extracts the tunnel texture. Figure 7 shows the procedure for establishing a 3D tunnel model.

Discussion
is paper discusses a three-dimensional modeling system of tunnel construction, which can establish the three-dimensional model of tunnel and then ensure the safety of tunnel construction. e following case study demonstrates the benefits of the 3D modeling system for tunnel construction projects. e 3D modeling system was developed using Autodesk AutoCAD 2010 and Visual C++ tool. e system is used to analyze the rock geological layers along the drilling direction and visualize the drilling process. Potential hazards are identified in advance to ensure safety during construction. e project is the construction of the Shanxi Long Tunnel, a section of the Taigu Tunnel, which is the secondlongest tunnel in China. It is located in the east side of the Lu-Liang mountainous area in Shanxi Province and passes through coal seams, goaves, and other adverse geological environments, as shown in Figure 8. e tunnel is characterized by a split pattern design. e left way is 13654 m long, and the right way is 13570 m. e construction sections start from K0 +9 00 and end at K7 + 550 m. e tunnel altitude is between 887.92 and 1385.58 m.  e angle between the mainline and the horizontal projection is 5.9°. e project was to be completed in December 2014. Owing to the absence of a reference coordinate system for the long tunnel, a hand measuring wheel was made to run along the tunnel to paint marks on the wall at 15 m intervals to increase the positioning accuracy.
Conventional geological forecasting involves tunnel rock monitoring, hydrogeologic monitoring, geological mapping, and rock category identification. In this study, a SIR-3000 series GPR manufactured in the US was selected to detect the geological layer, and its frequency is 100 MHz. A continuous reading of the radar detection signal was completed manually at the marked points. It took approximately 15-20 min to capture the images for one round. e radar images were processed in the control center at the entrance of the tunnel. Figure 9(a) shows the cross-sectional view of five testing points. Figure 9(b) shows the radar antenna. Figure 9(c) shows the marking on the wall for the testing points. Figure 9(d) shows the recorded original signals. e two-dimensional datasets obtained using the GPR in the time-domain were collected, including the original echo signals from the tunnel surfaces. However, it was difficult to extract the target from the background signals shown in Figure 10. e SAR was applied to focus the original beams. e radar speed and echo depth were analyzed simultaneously. Figure 10(a) shows the SAR image. e straight lines indicate the radar beams. A low-pass filter was used for noise reduction. Figure 10(b) shows the image after denoising. An edge detection technique was used to determine the contours of each geological layer. Figure 10(c) shows the edge extraction process. It can be seen from the Figure 10(c) that there are many noise signals, which need to be processed. e lines in the images indicate noise generated due to the coupling of the transmitting and receiving antenna, and an adaptive method is used to eliminate this noise, as shown in Figure 10(d).
OPGL package and C++ combined with AutoCAD were used to develop the computer software to visualize the process of tunnel excavation. e design data were stored in the database in advance and compared with the process in practice. e points in the 2D spaces were projected into a 3D space by Hough transform. Figure 10 illustrates the historical excavation process and construction process. e left part is the 3D tunnel model of the geological structure. e red color indicates the rock layers, and the purple color indicates the tunnel lining. e green color indicates the excavation face. e cross-sectional view is displayed in the right part of the software interface. e right and left holes can be displayed based on user requirement, as shown in Figure 11.
It can be seen from Figure 11 that this method can more accurately estimate the geological characteristics of the tunnel, which is beneficial to ensure safety during the actual construction of a long tunnel.

Conclusions
e hazard risk during the construction period of a tunnel should be identified in advance. In this study, the GPR detection technique was applied to obtain the geotechnical   Journal of Advanced Transportation features of a tunnel. e original radar signals were collected at five points. After the analysis, the following conclusions were obtained: (1) e B-scan strategy was used to obtain 2D images from the cross-sectional view when the antenna moves along the x direction. e SAR algorithm was applied to the original signals to focus all the beams at certain points. (2) Since the basic B-scan signals in the time domain cannot handle waves of varying speed, a 2D Fourier transform was applied to interpolate the signal into the Fourier space based on the SAR algorithm. (3) e signal processing consists of two phases: lowpass filtering and edge detection. In this study, a low-pass filter was used to detect images and suppress the noise. All the noise signals were considered high-frequency components, with a smoothing effect for images. e contour position of each geological layer was calculated using edge detection recognition. (4) A 3D modeling software was developed using Visual C++ for a better analysis of the tunnel construction process. e 3D model produces animations that serve as accurate representations to predict the construction process, thus facilitating the control of the construction process and identification of geological risks. (5) e Taigu Tunnel project was taken as an example to demonstrate the efficiency of the developed software. Different rock layers could be visualized using the 3D model, including the tunnel cross-sectional view. e 3D computerized model can also be applied to other construction projects.

Data Availability
e data used to support the findings of this study are available from the corresponding author upon request.

Conflicts of Interest
e authors declare that there are no conflicts of interest regarding the publication of this paper.