Arch. Min. Sci., Vol. 60 (2015), No 2, p. 645656 Electronic version (in color) of this paper is available: http://mining.archives.pl DOI 10.1515/amsc-2015-0042 SYLWIA TOMECKA-SUCHO*, HENRYK MARCAK* UYCIE ATRYBUTÓW GPR DO WYZNACZANIA REJONÓW ZAGROONYCH POJAWIENIEM SI PUSTEK POEKSPLOATACYJNYCH Sinkholes which occur in regions of old mine workings increase the risk to building and transport safety. Geophysical surveys, particularly with the use of ground penetrating radar (GPR), can help to locate underground voids which migrate towards the surface before they transform into sinkholes. The mining region in Upper Silesia, Poland was selected to test the method. The test was carried out on the profile at which sinkhole appeared few months after measurements. It can be assumed that the development of deformations in the ground was preceded by hydraulic and geomechanical processes, which directly caused this event. To identify the cause of the sinkhole formation exactly in this place in which it is located we carried out interpretation of GPR measurements through the calculation of GPR signals attributes such as instantaneous phase, instantaneous amplitude envelope, envelope derivative, envelope second derivative. The difference between two similar recorded data can be interpreted as a result of existence of hydraulic channels. On reflection, it appears that GPR signals attributes can be an important tool not only in the location of a cavity voids, but also can help in understanding the mechanisms of formation of the sinkholes. Keywords: GPR method, sinkholes, GPR attributes, geophysics in mines Osiadajce plytkie zroby po wyeksploatowanym wglu mog by przyczyn rozwoju procesów zapadliskowych na powierzchni ziemi. Takie zapadliska pojawiaj si nagle i stanowi due zagroenie dla ludzi, zwierzt i obiektów budowlanych. Do lokalizacji ,,wdrujcych pustek" w kierunku powierzchni ziemi mona wykorzysta metody geofizyczne, w szczególnoci metod georadarow. Przeprowadzono badania testowe w rejonie Sierszy na Górnym lsku, gdzie równie powstaj zapadliska górnicze. Prowadzono badania na profilu pomiarowym w czasie bezporednio poprzedzajcym pojawienie si zapadliska. Mona przyj, e obserwacje georadarowe byly prowadzone w czasie, kiedy rozwijal si geomechaniczny i hydrauliczny proces przygotowania deformacji zapadliskowej. W interpretacji materialu pomiarowego zastosowano metod atrybutów sygnalów georadarowych takich jak faza chwilowa, pierwsza * AGH UNIVERSITY OF SCIENCE AND TECHNOLOGY, FACULTY OF GEOLOGY, GEOPHYSICS AND ENVIRONMENTAL PROTECTION, DEPARTMENT OF GEOPHYSICS, AL. A. MICKIEWICZA 30, 30-059 KRAKOW, POLAND. E-MAIL: tomecka@ agh.edu.pl, email@example.com pochodna i druga pochodna sygnalu analitycznego. Anomalne wartoci atrybutów wskazuj na rozwój kanalów hydraulicznych, które byly bezporedni przyczyn powstania zapadliska. Z pracy wynika, e badanie rozkladu atrybutów analitycznego sygnalu georadarowego moe pomóc w zrozumieniu procesu tworzenia zapadliska i identyfikowa istnienie takiego procesu. Slowa kluczowe: metoda georadarowa, pustki wdrujce, atrybuty georadarowe 1. Introduction A sinkhole is a natural hole in the Earth surface caused by karst processes or collapse of old mine workings close to the surface. Generally, such holes form in the ground in a short period of time, and can have catastrophic consequences when they appear under buildings or transport tracks. There are many examples demonstrating the efficiency of using geophysical methods for forensic evaluation of sinkhole potential (Dobecki et al., 2006). In particular, the georadar method has been very useful in locating the developing destruction leading to sinkholes (Al-fares et al., 2002; Beres et al., 2001; Chamberlain, 2000; Cokun, 2012; Jeng et al., 2012; Kofman et al., 2006; McMechan et al., 1998). In Poland shallow copper and zinc mining has been carried out in the northern part of the Upper Silesia Mining Basin for hundreds of years. It is a cause of sinkholes in this region. The subsidence of rocks overlying old mine workings has formed upward-moving loose volumes surrounded by fractured zones (Marcak, 1999). Hydraulic channels are formed, water removes soil and rock particles from the upper part of the overburden to the bottom of the hydraulic system, causing gradual and slow migration of the loose volume towards the ground surface and thus creating the risk of sinkholes. Like in the results obtained in other parts of the world, these near-surface loose volumes can be detected using geophysical methods. In Poland different geophysical methods were used for locating moving voids in mining areas, such as microgravity (Fajklewicz, 1985), geoelectrical methods (Zakolski, 1974) and ground penetrating radar (GPR) (Marcak et al., 2006, 2008, 2010). Changes in physical properties, such as bulk density, electric resistivity or dielectric permittivity, associated with underground voids and their surroundings can be detected using geophysical methods in order to locate the centre of the void. If the location of the centre of the void is known, the risk of subsidence can easily be reduced with a relatively small amount of work and at a low cost by filling the void with rock materials. Previous experience gained from geophysical surveys confirms that GPR can be used to locate and to determine the parameters associated with underground voids migrating towards the ground surface (Batayneh et al., 2002; Beres et al., 2001; Chamberlain, 2000; Marcak et al., 2006, 2008, 2010). It encourages further research on the application of advanced data processing focused on extracting specific information in identifying affected zones and voids. This article demonstrates that GPR data collected in an area underlain by a developing underground void can be interpreted using GPR attributes. These attributes correspond with seismic attributes and they indeed facilitate the extraction of useful information from survey data. 2. GPR method The GPR method uses high-frequency electromagnetic waves (10-2000 MHz) in order to image subsurface structures based on changes in electric parameters such as resistivity and dielectric permittivity. These parameters can be directly related to geological, physical and chemical properties of soils. GPR is a method which uses antennas to send a pulsating signal underground and records reflected and diffracted waves in time. The record from the receiver is known as a radiogram (Carcione et al., 2000; Daniels, 1996). The most important elements of GPR are the antennas as well as the transmitter and receiver which control them. The frequency spectrum, of the impulse sent by the transmitting antenna as well as location of the transmitting and receiving antennas significantly influence the depth of penetration and resolution of data. The transmitter sends an impulse into the subsurface, while the receiver records the echo reflected and diffracted from geological interfaces and anomalous underground features. The length of a single impulse ranges from 0.5 to 100 ns which results in short impulses with a wide amplitude spectrum. Thus the resolution of the method is high. The GPR impulse has a relatively high peak power (up to several hundred watts). 3. Dielectric permittivity of materials surrounding a migrating underground void Underground voids and affected zones surrounding them are usually filled with water. This results in alteration of the dielectric permittivity within the void and its surrounding zone. Table 1 summarises the electromagnetic parameters, i.e. dielectric permittivity (r), electric conductivity () and velocity of electromagnetic signal (v) as well as attenuation coefficient (), associated with the near-surface strata in the area of Siersza. TABLE 1 Summary of electromagnetic parameters of geological materials associated with voids and zones surrounding it (Annan, 1999, 2001; Carcione et al., 2000) Relative dielectric permittivity r (-) Electric conductivity (mS/m) Velocity of electromagnetic signal v (m/ns) Attenuation coefficient of strata (dB/m) Material Air Fresh water Mineralised water (average value related to anion and cation concentration in water) Sand, sandstone (Dry-wet) Loam (Dry-wet) 1 80 80 330 540 0.0 0.010.5 1000 0.011.0 21000 0.30 0.03 0.03 0.170.05 0.130.05 0.0 0.010.1 1000 0.00.3 1100 Empirical research (Saarenketo, 1998) indicates that the electric properties of soil are strongly influenced by the water content. The dielectric permittivity of soil can be calculated on the basis of the CRIM model (Carcione, 2007) expressed by the following relationship: ( ) = vk (1) where vk is content by volume of three components: air (a), water (w) and rock skeleton (s). There are the following relationships: va = (1 S ) vw = S vs = 1 where is the porosity, S is the water saturation coefficient, and k () is the dielectric permittivity of individual rock components. Changes in dielectric permittivity cause changes in the propagation of radar waves. Reflection of radar waves happens particularly on interfaces between dry and water-saturated zones, while the reflection coefficient R (when position of transmitter and receiver is the same) can be determined based on the formula (Carcione et al., 2006): R= (2) saturated dry saturated + dry (3) where dry is the dielectric permittivity of dry soil, and saturated is the dielectric permittivity of water-saturated soil. 4. GPR traces The transmitting antenna sends a signal underground. In the subsurface it is reflected from interfaces associated with changes in dielectric permittivity. Waves are scattered within diffracting strata and form scattered waves. Reflected and scattered waves which travel in the direction opposite to the direction of the propagated GPR signal are recorded using the receiving antenna. During recording the signal is changed into a series of samples representing one column in radar records. Thus, a radiogram represents a matrix of samples in which rows show amplitude values of the electric field recorded for the same travel time across traces, while columns show recorded traces (a column is a record of the amplitude in time for a single location of the receiving and transmitting antennas, while rows are records of amplitudes for the same time on different traces). These matrices will be transformed in this article in order to extract information included in recorded data. GPR traces can be treated as the result of the convolution of spikes (corresponding to interface locations in time scale and proportional to their reflection coefficients) with the signal function (wavelets) (Fig. 1). The wavelet depends on the wave source and on the properties of the subsurface through which the signal travels (Sheriff et al., 1995). Properties of the signal function can be highlighted by calculating attributes. Fig. 1. Schematic explanation of the nature of GPR traces. GPR trace = series of Dirac functions* convoluted with amplitude function (trace = spike*wavelet), * convolution 5. Attributes An attribute can be any property of wave records which is useful in recognizing the geological structure. A large number of attributes are used in interpretation of seismic data (Chopra et al., 2007). The proper special class of GPR attributes has been introduced based on analytical signals (Taner et al., 1979). The instantaneous class of GPR attributes, like seismic attributes, can be calculated based on complex signals. Complex trace F(t) is defined as: F t f t i f t ~ ~ (4) where f(t) is the real part of the analytical signal trace, and its imaginary part f (t) is calculated using the Hilbert transform. f t f t (5) This way, a series of attributes can be calculated based on the analytical trace: amplitude envelope ~ 2 A t f 2 t f 2 t 1 (6) The GPR amplitude envelope is a positive value and is a measure of the total signal energy. In the convolution model the amplitude of the spikes depends on the signal energy, and according to formula (3) is directly dependent on the dielectric permittivity contrast in the reflective boundary. Consequently, this attribute in our examples depends on the difference in water satura- tion between two layers producing signal reflection at their interface. Changes in this parameter indicate changes in the distribution of dielectric permittivity and distribution of reflection coefficient which is related to dielectric permittivity. The following attributes were calculating: instantaneous phase ~ f t t arctan f t (7) This attribute can be used to correlate the phases of the signal between traces in GPR records, thus to eliminate natural amplitude decrease with depth. Phase continuity is one of the main parameters used in the interpretation of GPR records. instantaneous frequency d = const dt (8) Instantaneous frequency is a parameter which indicates changes in mean signal frequency and shows changes in the lithology of the investigated strata. envelope derivative R t = dA dt (9) The envelope derivative is an attribute which indicates changes in physical properties of the subsurface. It is also sensitive to absorption of energy of the wave field. It indicates discontinuities in the subsurface. second derivative of envelope D t = d2A dt 2 (10) This attribute indicates locations of changes in physical properties in the investigated subsurface. It is also a good indicator of lateral changes in electrical properties of strata. This parameter seems to be useful in differentiating the water saturation in the investigated areas. The entire GPR record is a conglomeration of different information that can be broken down into more basic information with the use of the attributes. In particular, the amplitude envelope attribute is directly related to the distribution of dielectric permittivity contrast along the reflection border, while the instantaneous phase attribute correlation demonstrates the continuity along geological borders and local discontinuities. Also, the two attributes of envelope derivative and second derivative of envelope help in evaluating the horizontal changes in the GPR anomalies. In summary it can be stated that the attributes extract specific information from the record field. Not all of the attributes provided information in a specific investigation example in Siersza, discussed in the paper. However, some of them were very useful. 6. GPR survey in the Siersza area The Siersza area has been intensively mined since the 18th century. First, zinc and copper ores were extracted, and later mining focused on deeper coal deposits. Dewatering of operating mines and flooding of abandoned old mine workings led to the creation of depressions and sinkholes. GPR data were collected from an area where shortly after the survey a sinkhole appeared at the surface. The survey was carried out using bistatic radar RAMAC GPR with a 250 MHz antenna. The results of the GPR measurements are shown in Fig. 2. The depth scale is constructed for the wave velocity which in the subsurface is equal to 10 cm/ns. The section where the surface wave is recorded (depth of between 0.0 and 0.6 m) has been excluded from the analysis along with data recorded below 100 ns, which due to radar wave attenuation is not interpretable. Fig. 2. Radiogram from the area where between x = 21 m and x = 28 m a sinkhole appeared six months later. A is the area where the sinkhole appeared, B the anomalous area without the sinkhole The survey data were used to identify unambiguously the saturated zone. The observed anomalies were divided into two parts marked A and B. Initially they were similar, but six months after the initial measurements a sinkhole appeared in region marked as A. Below 60 meters on the echogram shown in Fig. 2 several lines appeared, which can be results of diffraction. The migration of data removes most of these effects with the exception of one, which is considered in the paper as the real effect of the geological structure. 7. Attributes calculated based on data from the Siersza survey line The distribution of the following attributes was calculated: · Instantaneous phase (Fig. 3), · Amplitude envelope (Fig. 4), · Envelope derivative (Fig. 5), · Envelope second derivative (Fig. 6). Water flowing through this channel weakened the overlying zone to such an extent that it collapsed. The channel appearance is the condition for construction the sinkhole. The question is if the channel location can be recognized from the GPR data. Fig. 3. Distribution of the instantaneous phase on the Siersza profile. The line was conducted along the correlated signals of the attribute Fig. 4. Distribution of the amplitude envelope on the Siersza profile Fig. 5. Distribution of the envelope derivative on the Siersza profile In Fig. 3, which shows the distribution of the instantaneous phase, a probably hydraulic channel is visible, which joins area A with the underlying decreased-density zone. The distribution of all three attributes: the amplitude envelope (Fig. 4), the envelope derivative (Fig. 5), and the envelope second derivative (Fig. 6), indicate the presence of both anomalies, A and B, similarly as shown by initial data (Fig. 2). A comparison of Figs 5 and 6 reveals that anomaly A is smooth and deeper, while anomaly B is jagged. One may notice the linear part of the left edge on anomaly A, more clearly visible on pictures of derivatives. The envelope derivative (Fig. 5) clearly contours the triangle of anomaly A where the sinkhole was created. However, anomaly B has a condensed character only in the shallow zone, while deeper there are only individual strips. Fig. 6. Distribution of the envelope second derivative on the Siersza profile 8. Attribute processing The attribute distribution can be further transformed by normalising its values. Because the attributes have their physical interpretation while normalised attributes are informative parameters, when interpreting normalised values one should ensure the processing of the transformed data began with relatively simple procedures. Normalised attributes were subtracted one from another and only positive values were included in the final result (those matrix elements were shown which were strengthened with the new attribute). Results of such transformations were averaged on the net 7 × 7 elements and are shown for recorded data below 40 meters depth. Figure 7 presents such differences between the normalised distribution of the envelope second derivative and envelope first derivative. The results of calculations show that the marked left part of anomaly A indicates the stronger linear part of anomaly. This is a strongest linear element in this picture and it can representing the formulating sliding surface. As mentioned earlier, from the geochemical point of view the formation of a sinkhole must be preceded by washout of material from the floor of the future sink-hole. The anomalies shown on pictures (Fig. 7 and Fig. 8) create a system of potential channels, through which water can flow from the earth surface to subterranean aquifers. Our interpretation of results seems to be confirmed by the fact of appearance of a sinkhole in this center of anomaly A within six month. Concentration of water could obviously be expected in the area of anomaly B, but it is much more scattered, forming random local appearances of huge anomaly values. In contrast, in the region of left edge of anomaly A the concentrations re of linear character. This indicates that higher water saturation and more intense water migration should be expected in this part of anomaly A. It also indicates the possibility of the existence of a second channel near the one recognised based on the distribution of the instantaneous phase, and it could have contributed to the creation of the sinkhole. This channel was also recognised based on the difference between the envelope second derivative and the survey data (Fig. 8). It shows that both of the channels could have been connected and they controlled the distribution of stress in the studied part of the subsurface. Weakened strata in the left part of anomaly A, where material was transported downward through the channels, weakened the area of anomaly A shown in Figs 8 and 3. As a result the sinkhole was created here. Thus, the interpretation of attributes can be used not only for the location of loose volumes but also the recognition of the mechanisms which finally led to the development of the sinkhole. Fig. 7. Positive values of the difference between normalised values of the envelope second derivative and the envelope derivative Fig. 8. Positive values of the difference between normalised values of the second derivative of amplitude, attribute and the survey results shown in Fig. 2 9. Summary The emergence of the sink-hole is preceded by the changes in the local distribution of the soil strength, which is mostly the result of leaching of particles in zones of increased water circulation from the surface to underground workings. If the area of loosening is relatively large it is followed by a sinkhole. The predicting of the places in which it may arise sinkhole is comes down to location the deepest most saturated areas in the researched region. In the Siersza region two anomalies (A) and (B) are similar but the sinkhole appeared only in the area of the first anomaly. Interpretation of GPR attributes allow distinguish between the two anomalies showing that in the area of the first the saturated area is deeper than in the area of the second anomaly and the position of the probable strongest water flow stream can be located. The introduction of attributes to the interpretation of GPR data helps to present information included in survey data in a much more complete way than that based only on radiograms. As was shown using data from the Siersza area, the analysis of attribute distribution was utilised in order to recognise the mechanisms of sinkhole development in the mining area of Siersza. Calculating attributes from GPR data also provides an opportunity to use modern computer methods (Tadeusiewicz, 2011; Tadeusiewicz et al., 2014; Szymczyk et al., 2014), such as image recognition or neural networks, in modelling the subsurface investigated with the GPR method. Acknowledgements The project was funded by the National Research Centre based on Agreement No. UMO-2011/01/B/ST7/06178, awarded based on decision No. DEC-2011/01/B/ST7/06178.
Archives of Mining Sciences – de Gruyter
Published: Jun 1, 2015