background image

GEOLOGICA CARPATHICA, JUNE 2005, 56, 3, 285–294

Seismic reflection investigations for gas storage in aquifers

(Mura Depression, NE Slovenia)


1, 2


Environmental Agency of the Republic of Slovenia, Dunajska 47, SI-1000 Ljubljana, Slovenia;


University of Ljubljana, Faculty of Natural Sciences and Engineering, Askerèeva 12, SI-1000 Ljubljana, Slovenia

(Manuscript received January 7, 2004; accepted in revised form June 16, 2004)

Abstract: Two antiform structures in the Mura Depression were investigated as the most promising in Slovenia for the

construction  of  an  underground  gas  storage  facility  in  aquifers.  Seventeen  reflection  profiles  with  a  total  length  of

157 km were recorded, and three boreholes were drilled at their locations. Structural models based on interpretation of

the reflection seismic data, were constructed for the two main horizons (the pre-Tertiary basement and the Badenian/

Sarmatian boundary). Evaluation of different seismic velocity data was carried out to establish lateral velocity variations

in order to perform correct time-to-depth conversion. The porous rock in the Pecarovci structure is a 70 m thick layer of

dolomite, occurring at a depth of 1900 m, whereas layers of marl, several hundred meters thick, represent the imperme-

able cap rock. Due to faults, the Dankovci structure, at a depth of 1200 m where the reservoir rocks consist of thin layers

of conglomerate and sandstone, was found to be less reliable. 1D synthetic seismograms were used to correlate the

geological and seismic data at the borehole locations. Raytracing modelling was applied to confirm lateral continuity of

some horizons and to improve the structural interpretation at the locations of faults, which are critical factors for the

storage of gas.

Key words:  Pannonian  Basin,  Mura  Depression,  raytracing  modelling,  synthetic  seismograms,  gas  storage,  seismic



Slovenia  imports  almost  all  the  gas  it  needs.  The  supply

through  the  pipelines  is  fairly  constant  throughout  the  year,

but consumption is subject to seasonal changes. For this rea-

son, gas must be stored in the summer months in order to per-

mit higher consumption during the winter. For economic and

safety  reasons,  the  storage  of  natural  gas  is  reasonable  only

when underground (Dussaud 1989). There are four main types

of underground storage: salt caverns, abandoned mines, aqui-

fers or depleted oil or gas fields, and hard rock caverns (Gaus-

sens 1986).

The geological structure of Slovenia only permits storage in

aquifers (Fig. 1). Geological investigations for such facilities

have been ongoing for more than 15 years. The goal has been

to find an appropriate antiform structure, at a depth of between

500 and 2000 m, composed of porous (reservoir) rock with an

impermeable covering layer (Tek 1989; Flanigan 1995). In the

first stage, 13 different locations (Sadnikar 1993) were investi-

gated and two of them (Pecarovci and Dankovci) in the Mura

Depression were selected for further exploration. Geophysical

methods, especially reflection seismics, had an important role

in the evaluation of possible suitable locations.

The  Pecarovci  and  Dankovci  structures  are  located  in  the

Mura  Depression  (Gosar  2005),  on  the  slope  of  the  Murska

Sobota massif, which dips towards the Radgona depression. In

this area, the depth to the pre-Tertiary basement (Pt horizon) is

between  1800  and  2000 m.  The  possible  collector  rocks  for

gas  storage  are  the  Mesozoic  carbonates  in  the  pre-Tertiary

basement  and  the  thin  layers  of  porous  conglomerates  and

sandstones above the discordant boundary between the Bade-

nian  and  Sarmatian  layers  (KB  horizon)  inside  the  Tertiary

sediments (Turk 1993).

Seismic reflection investigations

A dense net of reflection seismic profiles was recorded at

the  Pecarovci-Dankovci  location.  For  the  construction  of  a

structural  model  in  an  area  measuring  8×8 km  (Fig. 2),  we

used 17 profiles. The total length of these profiles is 157 km

and their length inside the modelled area is 94 km. Geofizika

Zagreb performed the field data acquisition. Most of the pro-

files were recorded with line explosive sources (Geoflex and

Fig. 1. Schematic drawing of gas storage facility in an aquifer (af-

ter Gaussens 1986).

background image

286                                                                                                 GOSAR

Primacord),  end-offset  layout  and  a  Texas  Instruments  DSF

IV seismograph. The distance between geophone groups was

40 m, and CMP coverage was 24- or 30-fold. Three profiles

(Pec-1v-89,  Dan-1v-89  and  Pec-Dan-1v-89)  were  recorded

with Vibrosise source, split offset layout, a Texas Instruments

DSF V seismograph, 30 m geophones group spacing and 24-

fold coverage. Linear geophone arrays composed of 24 SM-4/

UB  10  Hz  geophones  were  used.  Standard  data  processing

Fig. 2. Location map of the seismic profiles and boreholes at the Pecarovci–Dankovci location.

was performed by INA-Naftaplin in Zagreb and by Western

Geophysical in London.

In  the  region  of  the  Pecarovci  and  Dankovci  structures,

three deep boreholes have been drilled to date (Fig. 2). Two of

them (Dan-1 and Pec-1) have reached the pre-Tertiary bedrock

at a depth of between 1800 and 1950 m, while Dan-3 has ter-

minated below the KB horizon at a depth of 1400 m (Sadnikar

1993). In the older Dan-1 borehole drilled for oil exploration,

background image


only basic electrical well-logging was performed. In the Dan-3

and Pec-1 boreholes, a complete set of well-logging measure-

ments was performed, including sonic velocity measurements.

In the Dan-3 and Pec-1 boreholes, down-hole seismic velocity

measurements  were  also  performed.  In  the  Mesozoic  dolo-

mites, overlying metamorphic complex thermal water with the

temperature of 102 °C was drilled in the Pec-1 borehole at a

depth of 1850 m (Sadnikar 1993). This borehole (Fig. 3) has

drilled  the  Mura,  Lendava  and  Murska  Sobota  Formations

(Gosar 2005). The main horizons interpreted in all seismic re-

flection profiles in the region are the KB horizon, which corre-

sponds to the discordant boundary between the Badenian and

Sarmatian sediments, and the Pt horizon, which corresponds

to the pre-Tertiary carbonate or metamorphic bedrock.

1D synthetic seismograms

Seismic  modelling  in  one  dimension  (1D)  is  a  commonly

applied tool for the correlation of seismic (time) data with geo-

logical or well-logging (depth) data (Sheriff & Geldart 1995).

This is important because, for the seismic data, lower vertical

resolution limited by the wavelength of the signal is character-

istic, but good lateral coverage along the profile is also appar-

ent. On the other hand, data from boreholes has good vertical

resolution, but is limited in lateral extent (Neidell 1981). In the

case of 1D modelling, the geological structure near the bore-

hole is approximated by horizontal layers. If the wavelength

Fig. 3.  Lithologic  column,  interval  velocity,  density,  acoustic  impedance,  reflection  coefficient  and  two  way  time  diagrams  for  the

Pec-1 borehole.

of the signal is small compared to the distance between adja-

cent interfaces, then good separated reflections are obtained.

But,  in  cases  where  the  layers  are  thin  with  respect  to  the

wavelength, the seismic trace is a result of the interference of

signals reflected from several interfaces. The theoretical verti-

cal resolution of seismic data is 1/4 of the signal wavelength,

whereas in practice it is not greater than 3/8 of the wavelength

(Widess 1973).

At  the  Pecarovci–Dankovci  location,  very  thin  layers  and

the problem of interference close to the KB horizon were en-

countered. 1D modelling was therefore applied in order to im-

prove the interpretation of the KB horizon, and to correlate the

seismic and borehole data. For this modelling, synthetic seis-

mograms were constructed for boreholes Pec-1 and Dan-3 us-

ing the Vista software package (Sis 1990). The reflectivity se-

ries was computed from the sonic and density log data. The

sonic log was corrected on the basis of the down-hole mea-

surements. In the lowest part, close to the Pt horizon, the re-

sults  of  laboratory  measurements  of  velocity  on  cores  were

also used. In the Pec-1 borehole, 49 layers of different acous-

tic impedance were distinguished (Fig. 3). The highest reflec-

tion  coefficients  corresponded  to  the  top  and  bottom  of  the

thin conglomerate layers above the KB horizon. We compared

the synthetic seismograms with two seismic profiles measured

close to the borehole location. Comparing the synthetic traces

with  the  Pec-Dan-1v-89  profile  (Fig. 4b),  which  is  trending

SW–NE, good correlation can be observed at 0.5 s, near the

KB  horizon  (0.9 s),  and  at  the  Pt  horizon  (1.35 s),  but  poor

background image

288                                                                                                 GOSAR

correlation between 0.7 and 0.9 s. If the same synthetic traces

are compared with the Pec-1v-89 profile, which is trending S–

N (Fig. 4a), there is also a good correlation at both main hori-

zons and at 0.5 s, but two or three good reflections between

0.7 and 0.85 s not visible on the previous profile can also be

observed. It can be concluded that, in this area, a high degree

of velocity anisotropy is encountered. 1D modelling was also

used to test what influence the frequency content of the input

signal has on the vertical resolution of the seismic data to sup-

port advanced processing of the data (Yilmaz 1987).

Raytracing modelling

With  2D  raytracing  modelling  stacked  seismic  sections,

field shot records and CMP (common midpoint) gathers were

simulated  (Fagin  1991).  2D  modelling  was  performed  using

the  Sierra  Quik  package  (Sierra  1990)  on  the  models  con-

structed with the Sierra Mimic program. Quik programs use

the asymptotic raytracing theory methods to find the path of

seismic energy between the source and receivers (May & Cov-

Fig. 4. Comparison of the synthetic seismogram for the Pec-1 borehole with the profiles Pec-1v-89 (a) and Pec-Dan-1v-89 (b).

ey 1981). At each intersection of the ray with a horizon, the

program computes the time and the reflection coefficient. In

each layer, the program uses straight rays even if the velocity

varies. This approximation is good enough for most models.

At interfaces, the rays are refracted according to Snell’s law.

In layers with variable velocity, the direction of the ray is cal-

culated in following manner. When the ray enters the layer,

the local velocity is used to determine direction. The ray then

continues straight to the next interface, where a new local ve-

locity is used to compute the direction in the lower layer. The

result  of  a  simulation  is  a  spike  section  of  reflection  coeffi-

cients. The convolution of a spike section with the input wave-

let results in a synthetic seismogram. Using this method, fully

complex  seismic  amplitudes  can  be  obtained  in  geological

structures of almost arbitrary configuration.

By  2D  modelling  of  stacked  seismic  sections,  an  attempt

was made to confirm the continuity of some reflections related

to the KB horizon and the interpretation of faults at the Pt ho-

rizon (Gosar 1995). The results of the simulation are presented

for the characteristic profile, Pec-Dan-1v-89 (Fig. 5). The ba-

sis for interpretation and construction of the time model was

background image


Fig. 5. The unmigrated (a), and migrated (b) seismic profile Pec-Dan-1v-89.

background image

290                                                                                                 GOSAR

the migrated seismic section (Figs. 5b and 6). When evaluat-

ing  the  structural  interpretation  of  the  Pt  horizon,  it  was

proved that better results are obtained with modelling of the

unmigrated  seismic  section  (Fig. 5a),  where  faults  are  more

Fig. 6. Line drawing interpretation of the seismic profile Pec-Dan-1v-89.

Fig. 7. Normal incidence raytracing for the profile Pec-Dan-1v-89.

evident because of diffractions. The input model consists of

nine layers of different velocity (Figs. 6 and 7). The two thin

layers of higher velocity represent the conglomerate sequenc-

es above the KB horizon. Between the Pecarovci and Dankov-

background image


Fig. 8. Synthetic seismograms for normal-incident raytracing shown in Fig. 7 (a) and for the corrected model (b).

ci structures there is a fault at which one sequence terminates

and the second is displaced. Above the Pt horizon, there is a

50 m  thick  layer  of  breccia.  Normal  incidence  raytracing,

which simulates the unmigrated seismic section, is shown on

the depth model of this profile in Fig. 7 and the corresponding

synthetic seismogram in Fig. 8a. By comparing the synthetic

and  the  original  seismic  section  (Fig. 5a),  it  was  concluded

that the structure on the NE side of the Pecarovci antiform is

more complex. To prove this, a new model of the Pt horizon

was constructed at this location with two normal faults on the

NE side of the Pecarovci antiform instead of only one. This

model was simplified to only one interface because it was rec-

ognized that the upper layers do not affect the rays significant-

ly.  The  synthetic  seismogram  for  this  simulation  (Fig. 8b)

showed  better  correspondence  with  the  seismic  section

(Fig. 5a). Therefore, we concluded that the Pecarovci structure

has, on its NE side, at least two normal faults. It is possible

also that the structure is even more complex.

Structural models

For structural modelling based on seismic reflection data, a

square area, 8×8 km in size, was selected. Seventeen profiles

with a total length of 94 km, and data from three boreholes,

were used (Fig. 2). The aim of the structural modelling was to

construct time and depth maps of the two most important hori-

zons,  that  is  the  KB  —  the  top  boundary  of  the  Badenian

background image

292                                                                                                 GOSAR

rocks, and the Pt — pre-Tertiary basement. Computer assisted

contouring  of  interpreted  surfaces,  taking  into  account  fault

traces, was applied.

First, a detailed analysis of the available velocity data was

carried  out  to  enable  correct  time-to-depth  conversion.  Four

types of velocity data were used: a — velocity analysis from

seismic data; b — down-hole measurements in boreholes; c —

sonic logs; and, d — laboratory measurements on cores. The

velocity function was based on the sonic log data, which were

corrected using the down-hole measurements, and fitted to the

datum plane of the seismic profiles (150 m above sea level).

Laboratory  measurements  on  cores  were  used  in  the  deeper

part of the Pec-1 borehole, where no other data was available.

Reflection  velocity  analysis  data  was  used  mainly  to  deter-

mine  lateral  velocity  changes.  The  velocity  isolines  for  the

characteristic profile Pec-Dan-1v-89 measured in the SW–NE

direction  showed  no  significant  lateral  velocity  variations  in

the upper part of the section, until a two-way time 1.0 s was

reached. However, there was a slight increase in velocity in

the NE direction, between 1.0 and 1.5 s. Established velocity

variation was applied for time-to-depth conversion.

The  structural  model  of  two  main  horizons  (KB  and  Pt),

showing  two-way  traveltime  isochrons  of  reflected  waves,

was first constructed (Figs. 9 and 10). A structural time map of

the KB horizon (Fig. 9) shows a closed antiform structure at

Dankovci  confined  by  the  980 ms  (milliseconds)  isochron,

which has a peak at 950 ms. Five faults cut the structure, but

they do not indicate significant slips. The reservoir rocks in

this  horizon  are  thin  layers  of  conglomerate  and  sandstone.

Small  quantities  of  oil  and  gas  were  found  in  these  layers,

proving the tightness of the cap rocks. On the other hand, be-

cause of the thin layers (a couple of meters thick), a fault could

easily separate two parts of the layer and reduce the volume of

the reservoir. At the Pecarovci location, there is no closed an-

tiform structure in the KB horizon.

On the Pt horizon (Fig. 10), there are closed antiform struc-

tures at both locations. The reservoir rock is a layer of porous

dolomite,  approximately  70 m  thick,  overlying  metamorphic

rocks  in  the  basement.  At  Dankovci,  it  is  confined  by  the

1460 ms isochron, while the top of the structure is at 1350 ms.

The area of the closed part is 5.42 km


. At Pecarovci, the top

of the structure is at 1350 ms and the antiform is confined by

the 1400 ms isochron. The highest point of the opening is on

the SW side. The area of the closed part is 1.576 km


. The pre-

liminary interpretation with a single fault on the NW side of

the anticline is presented in the model in Fig. 10. According to

raytracing modelling, the new interpretation with two normal

faults is shown in Fig. 11, together with a prognostic geologi-

cal profile. A minor reverse fault is also seen on the SW side

of the structure.

On  the  basis  of  the  time-to-depth  conversion,  structural

depth  maps  were  constructed  for  both  horizons.  The  depth

from  the  datum  plane  to  the  KB  horizon  is  from  1100 m  to

1200 m,  and  the  depth  to  the  Pt  horizon  is  from  1900 m  to

2000 m (Fig. 11b)

Of the three antiform structures, the Pecarovci (Pt) structure

was selected for further investigations as the most promising

for the construction of the gas storage facility. The structure at

Dankovci (Pt) is too big for the desired storage volume, while

the structure at Dankovci (KB) was found to be less reliable

because of faults (Sadnikar 1993).


On the basis of the constructed models, the available storage

volume in the Pecarovci (Pt) structure (Fig. 11) was estimated.

The calculations were performed using the Evasit Program for

the evaluation of porous gas storage facilities (Gaz de France

1990).  The  input  data  comprised  the  areas  of  closed  isoch-

Fig. 9. Time structural map of the KB horizon,  equidistance = 20 ms,

x — digitized points.

Fig. 10. Time structural map of the Pt horizon, equidistance = 20 ms,

x — digitized points.

background image


Fig. 11. Planned underground gas storage facility in the Pecarovci (Pt) structure: a — time structural map, b — depth structural map,

c — prognostic geological profile A–A’.

rones and corresponding seismic velocities in the top rock. As

the working volume, fifty percent of the total volume was tak-

en. Using this data, the working volume in the Pecarovci (Pt)

structure  was  estimated  on  275–300  million  m


(n).  This  is

above the desired minimum volume of 200 million m


(n) se-

lected  for  the  project  of  gas  storage  in  Slovenia  (Sadnikar


Among  all  evaluated  locations  in  Slovenia,  the  antiform

structure in the pre-Tertiary basement at Pecarovci was proved

to  be  the  most  promising  for  the  construction  of  an  under-

ground gas storage facility. Its structure is defined by seven

seismic profiles and one borehole. To prove the structural in-

terpretation of seismic data and to test the hydrogeological pa-

rameters of the reservoir layer and the impermeable cap rock,

another four boreholes are planned. The main disadvantage of

this structure is the great depth of the storage layer, which re-

quires compressors of higher power, and higher costs during


Acknowledgments:  The  author  is  grateful  to  M.  Bielik,  F.

Hubatka and J. Šefara for constructive reviews of the article.

Part of this work was conducted while the author was receiv-

background image

294                                                                                                 GOSAR

ing a grant from the University of Trieste. The assistance and

many  valuable  suggestions  of  R.  Nicolich  are  gratefully  ac-

knowledged. The author is indebted also to D. Ravnik, J. Sad-

nikar and S. Kranjc for their help and constructive remarks.


Dussaud M. 1989: Review of World Wide Storage Projects-France.

In: Tek M.R. (Ed.): Underground Storage of Natural Gas. Klu-

wer Academic Publishers, Dordrecht, 23–29.

Fagin S.W. 1991: Seismic modelling of geological structures. Soci-

ety of Exploration Geophysicists, Tulsa, 1–288.

Flanigan  O.  1995:  Underground  gas  storage  facilities,  design  and

implementation. Gulf Publishing, Huston, 1–179.

Gaussens P. 1986: Stockages souterrains de gas, Titre XIII. In: Man-

uel pour le transport et la distribution du gaz. Association tech-

nique de l’industrie du gaz en France, Paris, 1–333.

Gaz  de  France  1990:  Evasit,  PC  Version,  User’s  Manual.  Gaz  de

France, Paris, 1–25.

Gosar  A.  1995:  Modelling  of  seismic  reflection  data  for  under-

ground  gas  storage  in  the  Pecarovci  and  Dankovci  structures-

Mura depression. Geologija 37–38, 483–549 (in Slovenian).

Gosar A. 2005: Geophysical and structural characteristics of the pre-

Tertiary basement of the Mura Depression (SW Pannonian Ba-

sin, NE Slovenia). Geol. Carpathica 56, 2, 103–112.

May B.T. & Covey J.D. 1981: An inverse ray method for computing

geological structures from seismic reflections: zero offset case.

Geophysics 46, 4, 268–87.

Neidell N.S. 1981: Stratigraphic modelling and interpretation: geo-

physical principles and techniques. Education Course Note Se-

ries 13. AAPG, Tulsa, 1–141.

Sadnikar J. 1993: Investigations for underground gas storage in Slo-

venia. Rud.-Metalur. Zbor. 49, 1–2, 149–167 (in Slovenian).

Sheriff  R.E.  &  Geldart  L.P.  1995:  Exploration  seismology.  Cam-

bridge University Press, Cambridge, 1–592.

Sierra 1990: Quik raytracing, Release 3.2, User Manual. Sierra Geo-

physics, Seattle, 1–497.

SIS 1990: Vista 6.5, User Manual. Seismic Image Software, Calgary,


Tek  M.R.  1989:  Underground  storage  of  natural  gas.  Kluwer  Aca-

demic Publishers, Dordrecht, 1–458.

Turk  V.  1993:  Reinterpretation  of  chronostratigraphic  and  lithos-

tratigraphic correlations in the Mura depression. Rud.-Metalur.

Zbor. 49, 1–2, 149–167 (in Slovenian).

Widess  M.B.  1973:  How  thin  is  a  thin  bed?  Geophysics  38,  6,


Yilmaz  O.  1987:  Seismic  data  processing.  Society  of  Exploration

Geophysicists, Tulsa, 1–526.