During subway operation, various factors will cause long-term land subsidence, such as the vibration subsidence of foundation soil caused by train vibration load, incomplete consolidation deformation of foundation soil during tunnel construction, dense buildings and structures in the vicinity of the tunnel, and changes in water level in the stratum where the tunnel is located. The monitoring of long-term land subsidence during subway operation in high-density urban areas differs from that in low-density urban construction areas. The former is the gathering point of the entire urban population. There are many complex buildings around the project, busy road traffic, high pedestrian flow, and less vegetation cover. Several existing items require monitoring. However, monitoring distance is long, and providing early warning is difficult. This study uses the 2.8 km operation line between Wulin Square station and Ding’an Road station of Hangzhou Subway Line 1 as an example to propose the integrated method of DInSAR-GPS-GIS technology and the key algorithm for long-term land subsidence deformation. Then, it selects multiscene image data to analyze long-term land subsidence of high-density urban areas during subway operation. Results show that long-term land subsidence caused by the operation of Wulin Square station to Ding’an Road station of Hangzhou Subway Line 1 is small, with maximum subsidence of 30.64 mm, and minimum subsidence of 11.45 mm, and average subsidence ranging from 19.27 to 21.33 mm. And FLAC3D software was used to verify the monitoring situation, using the geological conditions of the soil in the study area and the tunnel profile to simulate the settlement under vehicle load, and the simulation results tended to be consistent with the monitoring situation.
Subways constructed in soft soil areas, such as Shanghai and Hangzhou in China, are currently in operation. By the end of 2020, 45 cities in mainland China have opened subway lines with a total length of 7978.19 km [
Long-term land subsidence monitoring in high-density urban areas during subway operation differs from those in low-density urban construction areas, urban suburbs, and other areas. A high-density urban area is a commercial, financial and trade, entertainment, and cultural center of a city. It is also the gathering point of the urban population. Many existing items require monitoring, including the main structure of a foundation pit, and buildings and structures along with subways, pipelines, and bridges [
Given the limitation in the period, most current domestic and foreign studies on long-term land subsidence during subway operation have focused on body subsidence and dynamic response caused by train vibration and the prediction of long-term land subsidence. In terms of body subsidence and dynamic response, Leclerc et al. [
In recent years, interferometric synthetic aperture radar (InSAR) technology has been increasingly applied to the topographic survey. This technology exhibits the following characteristics: noncontact survey, unnecessary monitoring of control network setting, sub-millimeter accuracy, high efficiency, low cost, wide-coverage, no weather influence, and high spatial resolution. Therefore, InSAR has also been applied to the monitoring of land subsidence during subway operations. Fadhillah et al. [
In this study, the 2.8 km operating line from Wulin Square station to Ding’an Road station of Hangzhou Metro Line 1 is used as an example to analyze the long-term subsidence in high-density urban areas during metro operation by InSAR-GPS-GIS technology.
With a total length of 2.8 km, the region from Wulin Square station to Ding’an Road station of Hangzhou Subway Line 1 is selected as the research area, including Wulin Square station, Fengqi Road station, Longxiangqiao station, and Ding’an Road station. The study area is the largest commercial street in Hangzhou. Yan’an Road connects the three largest business circles of Hangzhou (Wulin, Hubin, and Wushan Business Circles) in series, and the population of the study area is large. Many important buildings, such as Hangzhou Tower, Yintai Wulin Department Store, and Hangzhou Department Store, are distributed on both sides of the road. A total of 20 scenes Sentinel-1A satellite data are downloaded for processing in the research area. For InSAR techniques, the atmospheric delay benefit is one of the most significant error sources and cannot be completely eliminated, so various researchers find ways to reduce the errors. This study uses GPS technology for precise coordinates and multi-scene image data for error reduction.
The DEM data used in this study are obtained from Shuttle Radar Topography Mission (SRTM) 90 m data. The converted DEM data in standard format are shown in
POD Precise Orbit Ephemerides data are adopted to remove systematic errors caused by orbital errors. Sentinel-1 satellite POD are imaged into time format file S1A_OPER_AUX_POEORB_OPOD_20170628T121621_V20170607T225942_20170609T005942. The empirical orthogonal function (EOF) indicates that the imaging satellite is S1A, the release date of POD is June 28, 2017, and the data imaging date is June 08, 2017.
Combined with the research area, the key coordinates of the region from Wulin Square station to Ding’an Road station are located, and the coordinates of the key points are obtained, as shown in
Serial number | GPS point | Longitude | Latitude | Note |
---|---|---|---|---|
1 | Point 1 | 120.170928 | 30.278203 | Wulin Square station |
2 | Point 2 | 120.170353 | 30.268285 | Fengqi Road station |
3 | Point 3 | 120.170496 | 30.260612 | Longxiangqiao station |
4 | Point 4 | 120.171287 | 30.254123 | Jianhengli |
5 | Point 5 | 120.174305 | 30.251814 | Ding’an Road station |
6 | Point 6 | 120.16101 | 30.278702 | Tianmu Mountain Road and Huancheng West Road Intersection |
7 | Point 7 | 120.163022 | 30.266351 | Shengtangguan Pavilion |
8 | Point 8 | 120.168915 | 30.256868 | West Lake Lotus Area |
9 | Point 9 | 120.178617 | 30.270406 | North Zhonghe Road station |
10 | Point 10 | 120.181132 | 30.25431 | Qinghua Hotel |
The ArcMap component of ArcGIS Desktop is used for data processing to import experimental data and display them inline chart form through ArcMap software to illustrate the law of long-term land subsidence changing with time.
The DInSAR-GPS-GIS integrated method is primarily based on the multi-view integration platform, as shown in
The flowchart of DInSAR deformation processing is shown in
The obtained satellite data of Sentinel-1A for two scenes are used as an example, namely 20160601 and 20161128. 20160601 satellite data are selected as master images, whereas 20161128 satellite data are selected as slave images. The intensity map of the entire scene image in the research area is presented in
The theoretical height and theoretical deformation accuracy of baseline estimation in the research area are shown in
The baseline estimation results in the research area are as follows:
Normal Baseline (m) = 59.401, Critical Baseline min-max (m) = [−6452.217]−[6452.217], Range Shift (pixels) = 2.734, Azimuth Shift (pixels) = −8.031, Slant Range Distance (m) = 879019.561, Absolute Time Baseline (Days) = 180, Doppler Centroid diff. (Hz) = −31.079, Critical min-max (Hz) = [−486.486]−[486.486], 2π Ambiguity height (InSAR) (m) = 261.718, 2π Ambiguity displacement (DInSAR) (m) = 0.028, 1 Pixel Shift Ambiguity height (Stereo Radargrammetry) (m) = 21984.307, 1 Pixel Shift Ambiguity displacement (Amplitude Tracking) (m) = 2.330, Master Incidence Angle = 39.623, Absolute Incidence Angle difference = 0.004.
The baseline estimation results indicate that the spatial baseline of the two-scene master-slave satellite data is 59.401 m, which is considerably less than the critical baseline of 6452.217 m. It satisfies the accuracy condition. The time baseline is 180 days. The surface deformation represented by a phase change period during DInSAR processing is 0.028 m.
Image alignment is the first step that needs to be performed in InSAR data processing. This process requires that the same pair of pixel points of the two images correspond to the same pixel cell on the ground according to the coordinate mapping relationship between the corresponding points of the two images, and the complex values corresponding to the corresponding image elements in the two images are multiplied by the complex values, and the result of the two-phase conjugate multiplication is the interferogram. The selected polarization mode is VV, Range Looks is 5, Azimuth Looks is 1, Grid Size for Suggested Looks is 20, and Coregistration with DEM is True. The Global parameter sets Make TIFF True. The obtained interferogram after flattening in the research area is shown in
Three filtering methods, namely, adaptive, boxcar, and Goldstein, are available. Considering that the adaptive method is suitable for high-resolution data, the boxcar method uses the frequency of local interference fringes to optimize the filter and maintain small interference fringes as much as possible. Despite the variability of the filter in the Goldstein method, it improves the definition of interference fringes in the research area and reduces the decorrelation noise caused by spatial or temporal baselines. Therefore, the Goldstein method is selected to process the filtering of the research area, and the interferogram INTERF_out _fint (as shown in
Phase change takes 2π as its period. As long as the phase change exceeds 2π, the phase will start and the cycle will begin again. Therefore, phase unwrapping involves unwrapping the phase after flattening and filtering, such that it corresponds to the linearly changing terrain information, and solving the 2π ambiguity problem. Three types of unwrapping methods, namely, region growing, minimum cost flow (MCF), and Delaunay MCF, are available. The region growing method does not set an extremely high coherence threshold to leave adequate free growth space. The phase mutation part exists as an unwrapped island on the unwrapped image. This method reduces the errors caused by phase mutation. The MCF method can achieve better results than the region growing method when unwrapping is difficult due to a large area of low coherence or other factors that restrict growth. This method uses a square grid, considers all the pixels in the image, and masks the pixels whose coherence is less than the threshold. The Delaunay MCF method is different from the MCF method because it does not consider all the pixels in an image but only those whose coherence is greater than the threshold value. Moreover, it uses Droney triangular grids instead of square grids. Therefore, the MCF method is selected as the unwrapping method, 1.0 is selected for the unwrapping decomposition level, and 0.2 is selected for the unwrapping coherence threshold. The phase unwrapping results in INTERF_out_upha are shown in
GCPs should be selected away from the deformation area, phase jump area with unwrapping errors, and residual terrain phase. When dealing with surface deformation caused by DInSAR, a GCP can be selected at a position far from the deformation area in the absence of orbital errors. If orbital errors occur, then multiple GCPs should be considered as stable reference points from which the error phase can be calculated and removed by the program. The selection of GCPs is shown in
Refinement and reflattening involve the calculation of orbit refinement and phase shift, the elimination of possible slope phase, and the correction of satellite orbit and phase shift. Polynomial refinement is selected as the refinement method. The algorithm is robust and can be used even when the baseline is small. The Refinement Res Phase Poly Degree selects 3, which is the number of polynomials used to estimate the phase slope in the deflattening process. If the number of control point inputs is small, then the number of times will automatically decrease. The number 3 indicates that the phase slope with a constant phase shift in the range and azimuth directions will be corrected.
After the refinement and reflattening process, the optimized results are as follows:
ESTIMATE A RESIDUAL RAMP
Points selected by the user = 40, Valid points found = 40, Extra constrains = 2, Polynomial Degree choose = 3, Polynomial Type: = k0 + k1 * rg + k2 * az, Polynomial Coefficients (radians): k0 = 24.7173808990, k1 = 0.0006048432, k2 = −0.0046972519, Root Mean Square error (m) = 63.2093312825, Mean difference after Remove Residual refinement (rad) = −0.1101105943, Standard Deviation after Remove Residual refinement (rad) = 1.5099697025.
The interferogram INTERF_out_reflat_fint after reflattening is shown in
Phase-to-displacement conversion and geocoding combine the absolutely calibrated and unwrapped phases with the composite phase and convert them into deformation data and geocode them into a mapping coordinate system. By default, the system obtains deformation in the LOS direction. The product coherence threshold is 0.2, and the phase whose coherence is greater than this value is changed to the deformed value. Vertical displacement selects False, and the deformation in the vertical direction is not calculated. Slope displacement selects False, and the deformation in the slope direction is not calculated. Displacement custom direction selects False. Azimuth angle selects 0. The inclination angle selects 0. X dimension (m) selects 20. Y dimension (m) selects 20. Dummy removal selects True to mask the areas outside the image.
The deformation in the LOS direction is shown in
20160601 Satellite Data Analysis
The deformation in the LOS direction _slc_out_disp is shown in
The comparison between satellite DEM and deformation area in the research area is presented in
According to the location of Hangzhou Subway Line 1 and the surrounding location of the soil geological conditions, as shown in
Number | Name of material | Thickness | Young’s modulus | Poisson’s ratio | Weight density | Angle of internal friction | Cohesion |
---|---|---|---|---|---|---|---|
1 | Fill | 3 | 20 | 0.3 | 19.9 | 12 | 25.3 |
2 | Powdery clay | 6 | 17 | 0.3 | 19.3 | 10 | 27.2 |
3 | Silty clay | 7.5 | 13.5 | 0.3 | 18.8 | 8 | 17.5 |
4 | Powdery clay | 5 | 15 | 0.3 | 18.9 | 11 | 26 |
5 | Chalky soil with chalky sand | 38.5 | 40 | 0.3 | 19.6 | 30 | 11.6 |
6 | Tunnel lining | / | 35000 | 0.2 | 28 | / | / |
7 | Tunnel bed | / | 20000 | 0.2 | 21 | / | / |
The train dynamic load is the key to subway environmental vibration response analysis. Considering the limitation of the test conditions of the subway train dynamic load real measurement method and the limitation of the empirical formula method, the train one-track model is established using the US five-level track upset spectrum to obtain the train dynamic load; through a series of steps such as modeling, grid calculation cell division, parameter setting, constraint boundary conditions, establishing force field, calculating yield criterion establishment, etc., the vertical stress, maximum principal stress and displacement clouds of the tunnel-soil 3D finite element model are obtained and the results are analyzed, respectively.
Analysis of vertical stress and maximum principal stress results
The train will be subject to vertical force and horizontal force at the same time during operation, and the horizontal force usually dominates during high-speed travel.
From
The main stress is positively correlated with the tunnel depth and soil depth, and the maximum main stress is usually at the bottom, but due to the excavation of the tunnel, the stress concentration phenomenon is generated on the surrounding soil, resulting in a rapid increase in stress, and with the completion of excavation, the stress concentration phenomenon gradually disappears, and the stress near the tunnel only increases to a certain extent compared to the distal end, and the overall trend of stress distribution becomes larger with the increase of depth, and the maximum main stress also appears at The maximum principal stresses also appear at the bottom of the soil, and the results presented in the model calculation clouds agree with the basic theory to a large extent. The maximum principal stress diagram is shown in Displacement results analysis
The train vibration load will cause deformation, and the ground surface will collapse under long-term train vibration load, uneven settlement, and other geohazard problems.
From
Using the 2.8 km operation line between Wulin Square station and Ding’an Road station of Hangzhou Subway Line 1 as an example, this study proposes the integrated method of DInSAR-GPS-GIS technology, including Sentinel-1A satellite data, satellite POD, DEM data, GPS control points, and ArcGIS processing methods. This study presents the key algorithm for long-term land subsidence deformation during subway operation, including DInSAR deformation processing flow, baseline estimation, master-slave image file polarization mode selection, interferogram generation, filtering and coherence calculation, phase unwrapping, control point selection, orbit refinement and reflattening, phase transformation deformation, and geocoding. Multiscene image data, such as 20160601, 20161128, 20170608, and 20170912 satellite data, are selected to analyze long-term subsidence in high-density urban areas during subway operation. The long-term subsidence caused by the operation of Hangzhou Metro Line 1 from Wulin Square station to Ding’an Road station was found to be smaller in the study section compared to the long-term maximum subsidence of 240 mm in Nanjing Metro Line 10, China. The maximum subsidence is 30.64 mm, the minimum subsidence is 11.45 mm, and the average subsidence ranges from 19.27 to 21.33 mm. Based on the basic data of the study area, FLAC3D software was used to simulate the settlement of the study area subjected to vehicle loading, and it was found that the simulation results tended to be consistent with the monitoring by DInSAR-GPS-GIS technology, thus proving the reasonableness and accuracy of the algorithm.