Arch. Min. Sci., Vol. 59 (2014), No 2, p. 307321 Electronic version (in color) of this paper is available: http://mining.archives.pl DOI 10.2478/amsc-2014-0022 JAN BLACHOWSKI*1, ADAM CHRZANOWSKI**, ANNA SZOSTAK-CHRZANOWSKI*,** APPLICATION OF GIS METHODS IN ASSESSING EFFECTS OF MINING ACTIVITY ON SURFACE INFRASTRUCTURE ZASTOSOWANIE METOD GIS W OCENIE WPLYWU DZIALALNOCI GÓRNICZEJ NA INFRASTRUKTUR NA POWIERZCHNI Tilt (T), curvature (K) and horizontal strain () in ground subsidence troughs are the basic deformation parameters, which are used in the assessment of mining effects on surface infrastructure. The parameters can be determined from mathematical functions describing the continuous displacement field. The latter can be obtained by the least squares fitting of selected displacement functions to results of three-dimensional monitoring of horizontal and vertical displacements at discrete points. A methodology based on spatial data modelling in Geographic Information Systems (GIS) facilitates the above process as demonstrated on the example of a mining area in Canada. Polish guidelines regarding classification of mining risk categories based on the values of these parameters have been used in the example. Keywords: mining, deformation parameters, mining damage, GIS, map algebra Nachylenie (T) , krzywizna (K) oraz odksztalcenie poziome () w niecce obnieniowej to podstawowe wskaniki uywane do oceny wplywu dzialalnoci górniczej na infrastruktur na powierzchni. Wskaniki te mona wyznaczy poprzez wpasowanie odpowiedniej funkcji matematycznej opisujcej pole przemieszcze. Pole to mona otrzyma poprzez wpasowanie wybranych funkcji przemieszcze do wyników pomiarów (3D) poziomych i pionowych przesuni punktów sieci obserwacyjnej. Metodyka oparta na modelowaniu danych przestrzennych w systemach informacji geograficznej (GIS) wspomaga ten proces co zademonstrowano na przykladzie terenu górniczego w Kanadzie. Wykorzystano polskie wytyczne w zakresie kategorii terenów górniczych oparte na wartociach tych wskaników. Slowa kluczowe: górnictwo, wskaniki deformacji, szkody górnicze, GIS, algebra mapy * ** FACULTY OF GEOENGINEERING, MINING AND GEOLOGY, WROCLAW UNIVERSITY OF TECHNOLOGY, NA GROBLI 15, 50-421 WROCLAW, POLAND CANADIAN CENTRE FOR GEODETIC ENGINEERING, UNIVERSITY OF NEW BRUNSWICK, FREDERICTON, NB, E3B 5A3 CANADA CORRESPONDING AUTHOR: email@example.com 1. Introduction Effects of ground subsidence on the surface infrastructure in mining areas are usually evaluated by using values of tilt (T), curvature (K) or radius of curvature R = 1/K, and horizontal strain () parameters in the subsidence trough. The values of these parameters serve as criteria for evaluating potential damages to the objects located on the surface. If the functions describing horizontal and vertical displacement fields are known then the parameters can be calculated from the relationships (1): T x K x w x x , T y K y w y y 2 w y y 2 v y y (1) 2 w x x 2 u x x where u, v and w are the functions describing the components of displacement fields in the x, y and z directions of a 3-D coordinate system. Depending on the available information on geological, mining and geo-mechanical conditions, as well as on the results of deformation monitoring, the functions describing the continuous displacement field may be determined by one of the following methods: · Deterministic modeling of ground deformation (Szostak-Chrzanowski & Chrzanowski, 2008; Szostak-Chrzanowski et al., 2011); · Empirical modeling using, for instance, Knothe's theory of ground subsidence (Knothe, 1984; Malinowska & Hejmanowski, 2010; Sroka, 2009); · Three-dimensional (3-D) monitoring surveys using least squares fitting of a selected displacement function (e.g., higher order polynomials or spline functions) to the results at discrete points (Chen, 1983; Chrzanowski et al., 1986). Deterministic and empirical modeling are applicable at the design stage and during the mining operation. Monitoring surveys give horizontal and vertical displacements at discrete points during the mining operation. They may serve as a verification of the predicted deformations and as a source of information on the impact of mining activity on the surface infrastructure. The recent advances in spatial analysis methods, particularly spatial data modelling techniques employed in Geographic Information Systems (GIS) may help to determine the continuous displacement field and the resulting deformation parameters: T, K, and . The latter information is then used to classify the mining area into categories of potential risk to the surface structures. This paper discusses the proposed methodology based on the GIS technology, followed by a case study of evaluating effects of mining activity on the surface infrastructure in a mining area in Canada. 2. Methodology Based on GIS Technology The methodology of determining deformation parameters, used to assess mining influence on surface infrastructure is based on the analytical functions of GIS including interpolation, spatial data modelling (cartographic modelling) and map algebra. 2.1. Cartographic Modelling Cartographic modeling is a method of spatial data processing that generates new spatial data (i.e. map layers) based on input spatial data and operations performed on them (Heywood et al., 2005). The methodology uses the concept of Map Algebra (Tomlin, 2008) or Mapematics (Berry, 1987), where spatial data regarded as spatial variables and spatial operations equivalent to algebraic operations are used to develop a process model in GIS. Operation of these functions (e.g. extraction, overlay, generalization, reclassification, interpolation, etc.) can be controlled by additional parameters such as: range, mask, cell size of output data. New spatial data (variables) can be obtained in the result of a single operation or a sequence of operations that constitute a cartographic model. Cartographic modelling can be performed on raster type data where real objects, phenomena or processes are represented as a matrix of regular orthogonal cells (pixels) or vector type data representing these objects as points, lines or polygons defined by nodes with set coordinates. The raster type is usually used to represent continuous data (e.g. land use, terrain elevation or precipitation) and can store integer or floating point values. The vector type better represents discrete objects such as borders, technical infrastructure networks and other things. The application of the proposed method for deformation analysis is performed in two stages: (1) determination of continuous displacement fields from geodetic measurements of vertical and horizontal displacements at discrete points for a given period of time using interpolation techniques of GIS. (2) calculation of deformation parameters T, K and using map algebra for the displacement surfaces obtained in stage one. The results are stored in a geographical database as raster surfaces representing values of these parameters for individual cells or as isolines representing equal values of a given parameter. These data are used to prepare maps of mining terrain risk categories using raster reclassification techniques and to assess the potential of damage to surface infrastructure. 2.2. Determination of Continuous Horizontal and Vertical Displacement Fields Continuous fields of displacements can be obtained from measured displacements at discrete points by using selected spatial interpolation functions. These interpolation methods can be global, i.e. applying a single mathematical function to all sample points, or local, i.e. applying a single mathematical function repeatedly to subsets of sample points. The interpolating functions can produce either surfaces passing exactly through the observed points (i.e. be exact) or surfaces that adjust these values to fit a trend in the data (i.e. be approximate). In addition, some interpolation methods generate smooth surfaces between observed points (i.e. are gradual) or produce stepped surfaces (i.e. are abrupt). Thus, not all are suitable for estimation of every variable or phenomenon. Commercial and public domain GIS software offer a range of different interpolation techniques. The Spatial and Geostatistical analytical extensions of the ArcGIS v.10 software (ESRI, 2010a, 2010b) used in this study, have the following interpolation methods: Natural Neighbour, Inverse Distance Weighted, Trend, Spline, Radial Basis Function and Kriging. The last being a geostatistical interpolator. The Natural Neighbour is a method that uses the closest subset of input samples to a query point and applies weights to them based on proportionate areas to interpolate a value (Sibson, 1981). It is a local interpolator that will not estimate values above or below the observed quantitates (ESRI, 2010a). The Inverse Distance Weighted (IDW) method estimates values using a linearly weighted combination of a set of sample points. The weight is a function of inverse distance (ESRI, 2010a). It is a local and exact interpolator that does predict values above or below the maximum and minimum input values, as in the case of the Natural Neighbour, which may generate abnormal peaks and troughs in the surface. Therefore, both methods are unsuitable for our study. The Trend technique fits a surface defined by a mathematical function (a polynomial) through all the data points so that the difference between the estimated value at a sample point and its original value is minimized. It is a global and approximate method. One must remember that the higher the order of the polynomial the more complicated is the generated surface, i.e. first order polynomial function fits an inclined plane through the data points but a third order polynomial will fit a surface curved in two dimensions and so on. The Spline is a similar method to the trend but it works on a subset of local points. The spline function estimates values fitting a least curvature surface that passes, unlike in the trend method, directly through the input points (Franke, 1982). Therefore it is a local and exact method of interpolation. The Radial Basis Function is a special form of the spline function (ESRI, 2010b). It is a local and exact interpolator. Thus, the three methods: Trend, Spline, and Radial base Function can estimate values above the maximum or below the minimum measured displacement values (Kowalczyk et al., 2010). The Kriging function represents a group of geostatistical methods, where distance or direction reflect spatial correlations between measurement points (autocorrelation) that can be used to explain variation in the surface. It has the form of linear regression that is used to determine optimum value of the estimator (Olea, 1999). Kriging is a multistep process that comprises several components, such as: examining the data, calculating the empirical semi-variogram, fitting a model to the empirical values, generating the matrices of kriging equations, and solving them to obtain a predicted value and the error (uncertainty) associated with it for each location in the output surface. Unlike in, for example, the IDW method, the weights of the surrounding measured values used to derive a prediction for an un-sampled location are based not only on the distance between the measured points and the prediction location but also on the overall spatial arrangement of the measured points. To use the spatial arrangement in the weights, the spatial autocorrelation must be quantified (ESRI, 2010b). Assessment of particular interpolation functions for estimation of vertical and horizontal displacement fields is based on the interpolation mean error and root mean square error values. Preliminary tests on various interpolation methods have indicated that the smallest errors of partial vertical and horizontal displacement fields estimations for optimized interpolation param- eters, that is the ones producing the smallest differences between measured and predicted values determined with cross-validation technique (Hofierka et al., 2007), have been obtained for the spline, radial basis function and the kriging methods. These three interpolation techniques have been selected for further study. The displacement field surfaces calculated with these methods have been used as input data for determination of the deformation parameters in the case study described in section 3 of this presentation. 2.3. Determination of Deformation Parameters Using GIS The particular deformation parameters, for the resulting displacement fields, are calculated using the surface analysis and the map algebra GIS functions. The Tilt parameter surfaces are calculated applying the Slope spatial analysis function of ArcGIS software. This function calculates the first derivative of the subsidence surface. The value of each cell of the output raster representing slope is calculated based on the value and size of the input raster cell and the neighbouring eight cells on the moving window basis. This allows to determine the tilt from the following equation (2): dz dz T dx dy where: dz -- is the rate of change in the x direction, dx dz -- is the rate of change in the y direction. dy The Horizontal Strain parameter values are calculated applying the Slope function of ArcGIS (2). Because of the nonlinear characteristic of the horizontal displacements, first component displacement fields in the x and y directions respectively are interpolated using the functions described in part 2.2. and then combined to form the complete horizontal displacement surface using map algebra expression (3) and (4). These horizontal displacement field surfaces form the input data to calculate the horizontal strain parameter. S X 2 Y2 Y X (3) (4) (2) a tan 1 where: S -- is the displacement, a -- is the direction. The Curvature parameter values are calculated applying the Curvature spatial analysis function of ArcGIS to interpolated vertical displacement field surfaces treated as input spatial variables. This function calculates the second derivative of the function describing the subsidence surface. The value of each cell of the output raster representing curvature field is calculated based on the value and size of the input raster cell and the neighbouring eight cells on the cell by cell basis (ESRI, 2010a). For each location, a fourth order polynomial is formed (5): Z Ax 2 y 2 Bx 2 y Cxy 2 Dx 2 Ey 2 Fxy Gx Hy I where: Z -- is the height of the input surface, x, y -- are the coordinates of particular location. The coefficients are calculated from the 3 ×3 cell surface and curvature (K) is calculated from the equation (ESRI, 2010a): K = 2 × (D + E) × 100 where: Z1 Z 2 2 Z5 D L2 Z3 Z 4 2 Z5 E L2 (6) (5) (7) (8) Z1 to Z4 -- are the values of cells adjacent to the cell Z5, for which curvature is calculated, L -- is the dimension of cells. To obtain the radius (R) of curvature parameter, map algebra statement is used to calculate the reciprocal of curvature values. 2.4. Classification of Mining Areas into Categories of Risk for Surface Infrastructure Once the values of T, K, and are determined, they are used in a classification of mining areas into risk categories for the surface infrastructure. Here, use of Polish guidelines (GIG, 2000) is recommended. According to Polish regulations, mining areas are classified into six categories (Table 1) of a potential risk for the surface infrastructure as a function of the magnitude of the deformation parameters. TABLE 1 Risk classification of mining areas (GIG, 2000) Mining ground category Tilt T [mm/m] Radius of curvature K [km] Horizontal strain [mm/m] 0 I II III IV V T 0.5 0.5 < T 2.5 2.5 < T 5 5 < T 10 10 < T 15 15 < T 40 |R| 20 |R| <40 12 |R| <20 6 |R| <12 4 |R| <6 |R| < 4 || 0.3 0.3 < || 1.5 1.5 < || 3 3 < || 6 6 < || 9 9 < || 3. Case Study The described methodology has been used to assess a potential influence of withdrawal of large deposits of evaporates and underlying natural gas on the surface infrastructure, mainly small residential houses and farm buildings, in a mining area in Canada. 3.1. Monitoring Data Over the last 20 years, the mine has maintained a 3-D monitoring system consisting of geodetic leveling and GPS positioning. Configuration of the monitoring network has been changing annually due to destruction of some points by road and pipeline construction works and due to addition of monitoring points according to the progress of mining. Fig. 1 shows the outline of mining workings and distribution of the monitoring points on the surface in 2010. Fig. 1. Geodetic monitoring network (2010) and outline of underground workings (after Chrzanowski et al., 2012) Over the period of 15 years (1995-2010), the maximum vertical and horizontal displacements reached 0.75 m and 0.45 m respectively. Average standard deviation of horizontal and vertical displacements has been estimated at ±7 mm. In order to evaluate the potential effects of the observed ground deformation on the surface buildings, the analysis has been performed in 3 stages: (1) Determination of continuous vertical and horizontal displacement fields; (2) Determination of deformation parameters; and (3) Classification of the mining area into risk categories. The changeable configuration of the monitoring network limited the number of common points for displacement calculations between measurement epochs. For example, for the whole 1995-2010 period, the number of common points for vertical displacements has been limited to 46 and for horizontal displacements to 13 points only. In order to increase the number of common points and increase the accuracy of displacement field calculations, the 1995-2010 period has been divided into three periods. To determine total vertical displacement field the total period has been divided into three segments of: 1995-2000, 2000-2005 and 2005-2010 with 78, 78 and 56 common points respectively. Similarly, to determine total horizontal displacement field, three periods: 1995-2003, 2003-2006 and 2006-2010 have been used with 22, 28 and 23 joint points respectively. 3.2. Determination of Continuous Displacement Fields In order to get a better insight into the application of the cartographic methodology, three interpolation functions, namely spline, radial basis, and kriging have been used throughout the analysis. The interpolated surfaces representing vertical displacement fields in three consecutive periods have been saved in a raster format with cell size 10 by 10 m. These surfaces have been then added together using map algebra to obtain total vertical displacement fields in the 1995-2010 period for the three interpolation methods. The output maps provide subsidence values in the analysed space with horizontal accuracy corresponding to the pixel size (10 ×10 m). Fig. 2 shows vertical displacement fields for the 15 years of mining obtained from the three interpolation methods. Fig. 2. Comparison of vertical displacement fields in the 1995-2010 period determined with three interpolation methods Analogously, but taking into consideration direction of the displacements as described in part 2.3, partial horizontal displacement fields have been determined using the three interpolation functions. Surfaces representing horizontal displacement fields in three consecutive periods have been added together using map algebra to obtain total horizontal displacement fields. 3.3. Determination of Deformation Parameters In the second stage of the analysis, deformation parameters T, K and have been calculated and surfaces representing each parameter values have been generated using the interpolated horizontal and vertical displacement field surfaces. The maximum derived values of the tilt parameter during the 1995 to 2010 period range from 1.5 mm/m for the subsidence surface interpolated with the ordinary kriging function, 1.8 mm/m for the spline function and 2.5 mm/m for the radial basis function. The contours representing tilt parameter values and location of underground workings are shown in Fig. 3. The derived maximum absolute values of the radius of curvature during the 1995 to 2010 period range from 38.6 km for the subsidence surface interpolated with the ordinary kriging function, 30.4 km for the spline function and 23.4 km for the radial basis function. The curvature parameter maps are shown in Fig. 4. Fig. 3. Maps of the tilt parameter determined for surfaces obtained with three interpolation techniques Fig. 4. Maps of the radius of curvature parameter determined for surfaces obtained with three interpolation techniques The horizontal strain parameter values have been calculated applying the Slope and Map Algebra of ArcGIS functions to the three interpolated horizontal displacement field surfaces treated as input spatial variables. Directions of the horizontal displacement have been taken into account using non-linear interpolation technique as described in (Williams, 1999). The maximum derived absolute values of the horizontal strain parameter during the 1995 to 2010 period range from 0.7 mm/m for the horizontal displacement field surface interpolated with the spline function, 0.8 mm/m for the ordinary kriging function and 1.0 mm/m for the radial basis function. The horizontal strain parameter maps are shown in Fig. 5. 3.4. Accuracy Analysis The accuracy of determination of vertical and horizontal displacement fields has been affected by the limited number of input points and their uneven distribution. Calculation of total displacements as a sum of partial displacements for the three consecutive periods has increased the number of input points and improved their distribution. Partial displacement fields have been generated for optimized interpolation parameters producing the smallest differences between measured and estimated vertical displacement values determined with the cross-validation tech- Fig. 5. Maps of the horizontal strain parameter determined for surfaces obtained with three interpolation techniques nique (Hofierka et al., 2007). Fitting of the partial displacement fields with the three interpolation methods has been assessed based on mean error and root mean squared error statistics. In the presented example, the accuracy of the total displacement maps and deformation parameter maps may be affected by propagation of errors of partial displacement maps. Thus, the total displacement error St is equal to (9): 2 St S12 Sn (9) where Si are the partial displacement errors. With the aim of obtaining quantitative information on the accuracy of calculations, a set of test points with known values of displacements in the 1995-2010 period has been used. In the case of vertical displacement fields, 46 points with displacement values measured between 1995 and 2010 have been used as test points to determine the difference between the displacement values estimated with the three interpolation methods and the measured displacements. The results have been presented in Table 2. It shows the number of test points falling in consecutive 5% accuracy intervals. For the spline and radial basis interpolation function methods, 80% of the estimated values are within 5% of the measured values. For the kriging interpolation methods, 65% of the compared values are within 5%. The maximum differences between the estimated and measured values for single points have been 26% for the radial basis function, 29% for the spline function and 40% for the kriging function. TABLE 2 Analysis of displacement field interpolation accuracy Error of displacement determination Spline Number of test points RBF Kriging < 5% 5.1-10% 10.1-15% > 15% The differences between extreme values of deformation parameters obtained for surfaces determined by three interpolation methods are up to 40% for the tilt and curvature parameters between the radial basis function and the kriging methods. The differences between the spline and kriging methods reach up to 20%. 3.5. Assessment of the Deformation Effects on Surface Buildings Table 3 summarizes maximum values of tilt (T), radius of curvature (R), and horizontal strain () in the area of ground subsidence of the investigated mine. TABLE 3 Maximum values of deformation parameters calculated with cartographic modeling Interpolation Parameter Spline Radial basis function Kriging Tilt [mm/m] Radius of curvature [km] Horizontal strain [mm/m] Due to a lack of Canadian specifications regarding the resistance of buildings to ground deformation in mining areas, the aforementioned Polish guidelines have been used in the assessment of the potential mining effects. According to the guidelines summarized in Table 1, the obtained maximum values of deformation parameters correspond to risk category I. Figure 6 shows as an example, the classification map of the investigated area obtained from the evaluation of the tilt parameters. Fig. 7 shows a typical type of residential buildings in the investigated area. The experience gained in mining areas in Poland indicates that small (max. two stories) residential buildings, should not show any effects of ground deformations of category II or even III. Thus, the ground deformation in the investigated area can be considered as being harmless to the houses on the surface. More detailed analysis of the resistance of various structures to the ground deformation could be performed by using point scoring method (Kwiatek, 2007) in Fig. 6. Map of risk categories using determined ground tilts Fig. 7. Typical residential building in the investigated area the evaluation of individual buildings according to their dimensions, type of foundations etc. or assessing the correlation between the rate of building wear and the parameters of surface deformations (Firek & Wodynski, 2011) 4. Conclusions The method based on geographic information systems (GIS) in the determination of deformation parameters from scattered monitoring results has been successfully applied in a mining area in Canada. Three interpolation techniques: spline, radial basis function and kriging have been compared and evaluated. The maximum differences between the estimated and the measured values of displacements at test points have reached 30%. Each of these local and exact interpolation techniques has different characteristics. The kriging method produces the mildest surface of the three and being a geostatistical technique, allows for a more complex analysis of interpolation errors. The radial basis function method may generate more abrupt changes near the input points (ESRI, 2010b). In the authors opinion, based on the assessment of the accuracy, in the presented exercise, the spline interpolation method has produced the most consistent results. The problem of a small number of monitoring points in the case study in the analysed period of 15 years has been overcome by performing the analysis over three shorter periods and summing up the results. This has almost doubled the number of available measurements. According to the Polish guidelines, adopted for the purpose of this study, the obtained values of deformation parameters in the tested mining area, should not produce any damage to the structural integrity of the investigated properties on the surface. Acknowledgments This work has been partly financed from the National Science Centre Project UMO-2012/07/B/ ST10/04297 titled "Development of a numerical method of mining ground deformation modelling in complex geological and mining conditions" realised in the Institute of Mining Engineering at the Wroclaw University of Technology (Poland). References Berry J.K., 1987. Fundamental operations in computer-assisted map analysis. International Journal of Geographic Information Systems. International Journal of Geographical Information Systems, Vol. 1, No 2, p. 119-136. Chen Y.Q., 1983. Analysis of Deformation Surveys A Generalised Method. Ph.D. dissertation. Department of Geodesy and Geomatics Engineering. Technical Report, No 94, University of New Brunswick, Canada. Chrzanowski A., Chen Y.Q., Romero P., Secord J.M., 1986. Integration of Geodetic and Geotechnical Deformation Surveys in the Geosciences, Tectonophysics, 130, p. 369-383. Chrzanowski A., Blachowski J., Szostak-Chrzanowski A., 2012. Ocena wplywu eksploatacji górniczej na infrastruktur powierzchniow w kompleksowych warunkach geologiczno-górniczych z zastosowaniem metod GIS. [In:] Kowalski A. (Ed.), Bezpieczestwo i ochrona obiektów budowlanych na terenach górniczych. GIG, Katowice, p. 34-45; ESRI, 2010a. ArcGIS User Manual. Professionals Library Spatial Analyst, ESRI Press, USA. ESRI, 2010b. ArcGIS User Manual. Professionals Library Geostatistical Analyst, ESRI Press, USA. Firek K., Wodyski A., 2011. Qualitative and quantitative assessment of mining impacts influenceon traditional development in the mining areas. Arch. Min. Sci., Vol. 56, No 2, p. 179-188 Franke R., 1982. Smooth Interpolation of Scattered Data by Local Thin Plate Splines. Computer and Mathematics with Applications, Vol. 8, No 4, p. 273-281. GIG, 2000. Zasady oceny moliwoci prowadzenia podziemnej eksploatacji górniczej z uwagi na ochron obiektów budowlanych. Glówny Instytut Górnictwa, Instrukcja No 12. Heywood I., Cornelius S., Carver S., 2006. An Introduction to Geographic Information System, 3rd Edition, Pearson, Prentice Hall, Essex, England. Hofierka J., Cebecauer T., Súri M., 2007. Optimisation of Interpolation Parameters Using Cross-validation. Peckham R.J., Jordan G. (Ed.), Digital Terrain Modelling, Lecture Notes in Geoinformation and Cartography, Springer, 67-82. Knothe S., 1984. Prognozowanie wplywów eksploatacji górniczej, Wyd. lsk, Katowice. Kowalczyk K., Rapiski J., Mróz M., 2010. Analysis of Vertical Movements Modelling Through Various Interpolation Techniques. Acta Geodyn. Geomater., Vol. 7, No 4 (160), p. 399-409. Kwiatek J., 2007. Obiekty budowlane na terenach górniczych. Glówny Instytut Górnictwa. Katowice. Malinowska A., Hejmanowski R., 2010. Building damage risk assessment on mining terrains in Poland with GIS application. International Journal of Rock mechanics & Mining Sciences, Vol. 47, p. 238-245. Olea R.A., 1999. Geostatistics for Engineers and Earth Scientists. Kluwer Academic Publishers, Boston, USA. Sibson R., 1981. A Brief Description of Natural Neighbor Interpolation, In: Interpolating Multivariate Data. New York: John Wiley & Sons, p. 21-36. Sroka A., 2009. Influence of time-based parameters of longwall panel exploitation on objects inside rock mass and on the surface. Arch. Min. Sci., Vol. 54, No 4, p. 819-826. Szostak-Chrzanowski A., Chrzanowski A., Hiroji A., 2011. Ground Subsidence Modeling based on Stress Redistribution in Rock Mass, Journal GIG, Glówny Instytut Górnictwa, Katowice 2011, No 2/1/2011, p. 524-536. Szostak-Chrzanowski A., Chrzanowski A., 2008. Canadian Contributions to Monitoring and Physical Interpretation of Ground Subsidence. Markscheidewesen. Jahrgang 115 (2008), No 3, October, p. 14-20. Tomlin C.D., 2008. Cartographic Modelling. Kemp, K. K (Ed), Encyclopedia of Geographic Information Science, Sage. Williams R., 1999. Non-Linear Surface Interpolations: Which Way Is The Wind Blowing. Proceedings of the 1999 ESRI International User Conference, USA.
Archives of Mining Sciences – de Gruyter
Published: Jun 1, 2014