基于IGS跟踪站的大面积矿区GNSS变形监测-有色金属在线

简介概要

基于IGS跟踪站的大面积矿区GNSS变形监测

来源期刊：中国有色金属学报(英文版)2014年第2期

论文作者：卞和方张书毕张秋昭郑南山

文章页码：514 - 519

关键词：矿区沉降；卫星导航系统(GNSS)；抗差卡尔曼滤波(RKF)；IGS连续跟踪站

Key words：mining subsidence; global navigation satellite system (GNSS); Robust Kalman Filtering (RKF); IGS station

摘要：为了精确地对大面积矿区进行沉降监测，以周围IGS跟踪站为参考点，建立高精度卫星导航系统(GNSS)变形监测网。鉴于IGS跟踪站的非线性运动，结合速度场信息及周解坐标，给出抗差卡尔曼滤波模型(RKF)，并用于确定IGS跟踪站不同历元对应的坐标基准。通过皖北矿区沉降监测实例对该方法进行验证。结果表明：RKF模型确定的坐标基准优于预测模型及IGS分析中心提供的周解坐标；当观测时间大于4 h时，对应的监测精度可以达到毫米级。该方法可以高效、准确地监测大面积矿区沉降。

Abstract: In order to monitor large-area mining subsidence accurately, a high-precision global navigation satellite system (GNSS) monitoring network was established based on the nearby international GNSS service (IGS) stations taken as reference points. Given the non-linear motions of IGS stations, the robust Kalman filtering (RKF) model was presented to determine the datum of multi-period monitoring network considering the velocity and weekly solution of IGS stations. The theory proposed was applied to monitoring mining subsidence in northern Anhui coal mine in China. According to the case study, the RKF model to establish monitoring datum is better than the prediction method and the weekly solution from IGS analysis centers (ACs), and the corresponding precision of deformation can reach up to millimeter level with 4 h observation. The research provides an efficient and accurate approach for monitoring large-area mining subsidence.

Trans. Nonferrous Met. Soc. China 24(2014) 514-519

Monitoring large-area mining subsidence by GNSS based on IGS stations

He-fang BIAN¹, Shu-bi ZHANG^{1, 2}, Qiu-zhao ZHANG¹, Nan-shan ZHENG^{1, 2}

1. School of Environment and Spatial Information, China University of Mining and Technology, Xuzhou 221116, China;

2. Key Laboratory of Land Environment and Disaster Monitoring of SBSM, China University of Mining and Technology, Xuzhou 221116, China

Received 17 November 2012; accepted 3 November 2013

Key words: mining subsidence; global navigation satellite system (GNSS); Robust Kalman Filtering (RKF); IGS station

1 Introduction

Global navigation satellite system (GNSS) techniques, especially the global positioning system (GPS), are playing an increasing role in deformation monitoring, such as crustal deformation, bridge- deformation monitoring and mining subsidence [1-6]. In general, reference point, which is usually established in stable geological condition and easy to long preserve, is essential for deformation monitoring [5,6]. For a large area of mining, it is difficult to establish stable reference points within the mine deformation area.

As we all know, the international GNSS service (IGS) stations are characterized by stable geologic condition, easy to obtain observation data and continuous observation [7-10]. Obviously, the IGS stations can be used as reference points of known coordinates to precisely determine the datum of large-range deformation monitoring with long observation duration. However, the long observation duration limited the wide application of IGS stations in practical engineering. In this work, the feasibility of monitoring large-area mining subsidence is investigated by using IGS stations as reference points. Three methods associated with the establishment of deformation datum for multi-period monitoring network are studied. In the end, an efficient approach for monitoring large-area mining subsidence is demonstrated. The technology was successfully applied in northern Anhui coal mine in China.

2 Feasibility of IGS stations as reference stations for mining subsidence monitoring

Nowadays, IGS stations are quite scarce on global scale. For example, there are only 9 permanent stations in China mainland. Therefore, double difference (DD) model is adopted to process the baseline of the monitoring network combined GPS and GLONASS. With a constellation of 31 operational GPS satellites and 24 GLONASS satellites, compared with only GPS observation, the number of total observations is increased by a factor of 1.8 [11].

2.1 Model of combined GPS/GLONASS

For GPS double difference model, the observation equations are [9,12]

(1)

where is a double difference operator with respect to receivers kl and satellites ij; L₁ and L₂ are phase observations for both frequencies respectively; ρ is the distance between satellite and receiver; λ₁, λ₂ are the wavelengths of L₁ and L₂, respectively; I is the ionospheric refraction; T is the tropospheric refraction; n is the integer ambiguity; ε is the measurement noise and the errors which cannot be modeled.

Similarly, the GLONASS observation equations can be expressed as [9]

(2)

where is the double-difference ambiguity with respect to receivers kl and satellites ij; is the single-difference ambiguity with respect to receivers kl and satellites j; λⁱ is the carrier wavelength of satellite i; Δλ^ij is the carrier wavelength difference between satellites j and i where λⁱ=λ^j+Δλ^ij.

2.2 Feasibility verification

Data from three Chinese stations (BJFS, SHAO and WUHN) of the IGS network are selected for the reference points of monitoring network. The locations of the three IGS stations are given in Fig. 1. The distances between monitoring area to the reference points are between 400 and 700 km. The observations of TY01 monitoring point in mining area are available for 8 h. The distances between TY01 point to the three IGS stations (BJFS, SHAO and WUHN) are 681 km, 478 km and 417 km, respectively.

The root mean square (RMS) of TY01 station for different time-lengths (8 sessions) in double difference model is shown in Fig. 2, while the three IGS station (BJFS, SHAO and WUHN) are taken as fixed stations. In order to validate the true estimated accuracy of TY01 station for each session, the coordinates bias between different sessions estimated results and the true value (weekly solution) is illustrated in Fig. 2. It is obvious that 4 h GPS observation can achieve 5 mm positioning accuracy from Fig. 2. Therefore, the 4 h observation is sufficient for millimeter-level GPS deformation monitoring in large-area mine.

Fig. 1 IGS stations distribution for reference points of mining deformation monitoring

Fig. 2 RMS of TY01 station for different time-lengths in double difference model (a) and coordinates bias (b)

3 Robust Kalman filtering determining coordinates of reference points

3.1 Analysis of IGS station time-series

The international terrestrial reference frame (ITRF), as a realization of the international terrestrial reference system (ITRS), is presented by a set of station positions and linear velocities [12]. In order to test whether the coordinate variation of IGS stations is a constant velocity motion, the weekly solutions of SHAO station form 1512 to 1628 GPS weeks (December 31,2008 to March 23, 2011) are used to explain. From Fig. 3, it can be seen that from the weekly solutions, the coordinate variation is not completely linear, especially in y- and z-direction. Therefore, it is difficult to extrapolate the coordinates of IGS station accurately according to the velocity model and local geodetic datum.

Fig. 3 Coordinate changes of SHAO station

3.2 Robust Kalman filter model (RKF)

Considering the movement of IGS station as a constant speed c₀, the linear dynamic system of station coordinate can be given by [13,14]

X(t)=X(t-1)+c₀t+w_t (3)

Assuming sampling interval Δt=1 (week), t_k=k·Δt=k, and denote X(t_k)=X_k, the discrete-time system for Eq. (3) can be expressed as

X_k=X_k_-1+c₀+w_k (4)

Equation (3) is called state equation. And w_k is the dynamic noise, E(w_k)=0, k∈{1, 2, …, n}. X_k denotes the 3×1 vector of unknown random parameters for station coordinate at the epoch k.

Measurement equation:

L_k=X_k+e_k (5)

where L_k is the measured value, i.e., weekly solutions coordinates; e_k is the measurement noise. Combining Eqs. (4) and (5), Kalman filter equation can be expressed as [5,13]

(6)

with

(7)

where X(k/k-1) is the predicted state estimate; E_k is the prediction residual; J_k is the Kalman gain; Q_k is the dynamic noise covariance; R_k is the measurement noise covariance.

Because IGS station coordinate from GNSS observations may be affected by the observation environment, e.g., equipment changes, data problems, earth quakes and other geophysical processes, the Robust Kalman filtering (RKF) model is built to find outliers and discontinuities through equivalent Kalman gain [5,13],

(8)

where is the equivalent Kalman gain; k₀ and k₁ are the robust parameters, k₀=2.5-3.5, k₁=3.5-4.5.

The initial values of X(k/k) and D_X(k/k) can be determined by solving the first-period GNSS monitoring network, and the dynamic noise can be obtained through linear fitting the time series of IGS stations. In addition, the three-dimensional coordinate accuracy of the weekly solutions is 2-3 mm [11,15,16]. According to the international GPS service products, the velocity of the three IGS stations is shown in Table 1.

Table 1 Velocity of BJFS, SHAO and WUHN stations in ITRF2008 (mm/y)

Figure 4 demonstrates the three methods used for determining coordinate data of IGS station (form 1512 to 1628 GPS weeks), and y coordinate of SHAO station as an example. For the prediction method, the coordinate datum is extrapolated from the velocity and time span, and compared with the other two methods. the maximum discrepancy may be more than 10 mm.

3.3 Uniform data

The first-period GNSS monitoring network provides an initial coordinate datum which can be used to compare with the other period monitoring. Therefore, the latter period GNSS monitoring results should be unified into the initial coordinate datum using BURSA model [16,17]. After the latter period monitoring results are unified into the initial coordinate data, deformation of each point in the monitoring network can be achieved. A detailed technical scheme for monitoring large-area mining subsidence is shown in Fig. 5.

Fig. 4 RKF model determining y coordinate of SHAO station (Residual=RKF model-weekly solution)

Fig. 5 Technical scheme for GNSS large-area mine subsidence monitoring

4 Experimental

The technique proposed was applied to monitoring mining subsidence for Huaibei Mining Group Co., Ltd., which is located in northern Anhui province in China. The company has over 20 production mines and more than 50 years of coal mining history. When monitoring the whole mine area subsidence, it is difficult to establish stable points as monitoring reference points within the large deformation area. By establishing a GNSS monitoring network combined with the nearby IGS stations, the RKF model is introduced to determine the monitoring datum for large-area deformation.

4.1 Initial coordinate data establishment

In order to monitor the mining subsidence accurately, a research project was carried out to monitor the mine deformation by GNSS from March 4, 2009. The three IGS stations (BJFS, SHAO and WUHN) are taken as reference points, and the location of monitoring points with the three IGS stations is shown in Fig. 1. The detailed GNSS monitoring network is shown in Fig. 6, with coverage of 8000 km², and the averaged distance between monitoring points is about 35 km.

Fig. 6 Location of GNSS monitoring network points

As can be seen in Fig. 6, the GNSS monitoring network consists of 19 points: 13 mine monitoring points, 5 mine control points and one bedrock mark, which is built on the stable bedrock. Considering the stable bedrock underground 70 m depth, the bedrock mark is relatively stable in the deformation area. The first period GNSS observation was carried out at March 4, 2009, and the observation time of each monitoring point was more than 4 h, with an average of 5 h. The Bernese 5.0 software was used to process the observation data by double difference model, and the positioning precision achieves millimeter level with N, E directions 2-3 mm, and U direction less than 5 mm.

4.2 Determination of reference datum

According to the concept of ITRF, the realization includes station position and velocity, and the IGS station coordinates are changing for different epochs. Consequently, the coordinates of three reference stations (BJFS, SHAO and WUHN) of the first period observation are different from those of the second period observation which was carried out at November 10, 2010. When solving the second-period GNSS monitoring network, it needs to determine the coordinate datum accurately from the current epoch. Therefore, the three methods are used in the example. One typical method is to extrapolate station coordinates from reference epoch (the first observation time: March 4, 2009) to current epoch by applying ITRF velocity field [11], which is called prediction method. The other two methods are weekly solution and RKF model respectively. For the weekly solution method, the reference epoch is GPS week 1521, and the current epoch is GPS week 1609. The coordinate changes of the three reference stations are shown in Table 2 for the three methods.

Table 2 Coordinate changes of three reference stations (mm)

From Table 2, the coordinate datum of the second-period GNSS monitoring network can be calculated by

x_2nd=x_1st+Δx (9)

where x_1st is the initial coordinate datum and Δx is the coordinate change. Then, fixing the three reference stations, all the points coordinates of the GNSS monitoring network refer to the current epoch (November 10, 2010) can be obtained based on high precision baselines.

4.3 Result and analysis

After the second-period monitoring results are unified into the initial coordinate data, the deformation of each monitoring network can be achieved. The geodetic height of LQYK point and the height variation for different methods are presented in Table 3.

Table 3 Comparison of geodetic height of LQYK point for different methods

Due to the fact that the LQYK point is bedrock mark, the height variation should be zero theoretically. The result in Table 3 of height variation determined by RKF model shows good agreement with the theoretic value: the height variation is -6.8 mm. Compared with -18.6 mm height variation of the prediction method, the weekly solution method has a greater accuracy. From these results, it can be conclude a that RKF model is the best method for determining the coordinate data, while the prediction method is the worst.

According to the height variation in Table 3, the maximum discrepancy among the three methods is -11.8 mm which can be considered a systematic bias of each monitoring point. Therefore, the prediction method which considers only the velocity field of IGS station may result in systematic error of the coordinate data. Furthermore, RKF model, integrating the velocity field and weekly solution of IGS station, can be used to determine coordinate data accurately. In the end, the height variation of LQYK point by RKF model shows the external precision of GNSS monitoring network and also reaches millimeter level.

5 Conclusions

1) GNSS monitoring network can be used to monitor large-area mining subsidence with millimeter level precision based on the nearby IGS stations when the observation time is more than 4 h.

2) The RKF model to determine coordinate data of monitoring network may be better than the prediction method and the weekly solution, according to our case study.

3) A detailed technical scheme for GNSS monitoring the large-area mine subsidence, which is characterized by no reference points in the mine deformation area, has been presented and successfully applied in northern Anhui coal mine, showing significant improvement in efficiency.

References

[1] ZHANG X H, Andersen O B. Surface ice flow velocity and tide retrieval of the Amery ice shelf using precise point positioning [J]. J Geod, 2006, 80:171-176.

[2] LU Wei-cai, CHENG Shi-guang, YANG Hai-sheng, LIU Da-ping. Application of GPS technology to build a mine-subsidence observation station [J]. Journal of China University Mining and Technology, 2008, 18(3): 377-380.

[3] GAO Jing-xiang, LIU Chao, WANG Jian, LI Zeng-ke, MENG Xiang-chao. A new method for mining deformation monitoring with GPS-RTK [J]. The Chinese Journal of Nonferrous Metals, 2011, 21(6): 659-664. (in Chinese)

[4] ZHANG Qin, HUANG Guang-wen, Ding Xiao-guang, Wang Li, Zhao Chao-ying. Research and realization of high-precision GPS datum, considering plate movement, stability and system errors [J]. Chinese Journal of Geophysics, 2009, 52(12): 3158-3165. (in Chinese)

[5] YU Xue-xiang, LU Wei-cai. Robust Kalman filtering model and its application in GPS monitoring networks [J]. Acta Geodaitica et Cartographica Sinica, 2001, 30(1): 27-30. (in Chinese)

[6] WanG Jian, PENG Xiang-guo, XU Chang-hui. Coal mining GPS subsidence monitoring technology and its application [J]. Mining Science and Technology, 2011, 21(4):463-467.

[7] Altamimi Z, Collilieux X, Metivier L. ITRF2008: an improved solution of the international terrestrial reference frame [J]. Journal of Geodesy, 2011, 83: 457-473.

[8] YAO Yi-bin. On the non-linear motion of IGS station [J]. Geo-spatial Information Science, 2007, 10(4): 240-244.

[9] Dach R, Hugentobler U, Fridez P, Meindl M. Bernese GPS Software Version 5.0 Tutorial [M]. Bern: University of Bern, 2007: 38-122.

[10] Bei Jin-zhong, LI Yu-lin, Zhang Peng, Jiang Zhi-hao. Accuracy of GPS fiducial station coordinate and velocity field and their time series [J]. Geomatics and Information Science of Wuhan University, 2004, 29(9): 763-767. (in Chinese)

[11] Collilieux X, Metivier L, Altamimi Z, Dam TV, Ray J. Quality assessment of GPS reprocessed terrestrial reference frame [J]. GPS Solutions, 2011, 15(3): 219-231.

[12] Ray J, Dong D, Altamimi Z. IGS reference frames: status and future improvements [J]. GPS Solutions, 2004, 8(4): 251-266.

[13] Yang Yuan-xi. Adaptive navigation and kinematic positioning [M]. Beijing: Surveying and Mapping Press, 2006: 51-94. (in Chinese)

[14] Chai Hong-zhou, Cui Yue, Ming Feng. The determination of Chinese mainland crystal movement model using least-squares collocation [J]. Acta Geodaitica et Cartographica Sinica, 2009, 38(1): 61-65. (in Chinese)

[15] Collilieux X, Dam T V, Ray J, Coulot D, Me'tivier L, Altamimi Z. Strategies to mitigate aliasing of loading signals while estimating GPS frame parameters [J]. Journal of Geodesy, 2012, 86: 1-14.

[16] Ferland R, Piraszewski M. The IGS-combined station coordinate, earth rotation parameters and apparent geocenter [J]. J Geod, 2009, 83: 385-392.

[17] Collilieux X, Altamimi Z, Coulot D, Dam TV, Ray J. Impact of loading effects on determination of the international terrestrial reference frame [J]. Advances in Space Research, 2010, 45: 144-154.