Examensarbete vid Institutionen för geovetenskaper
Degree Project at the Department of Earth Sciences
ISSN 1650-6553 Nr 485
Analysis of Seismic Data Acquired in the Hverahlíð Geothermal Area
Agnieszka Stoch
INSTITUTIONEN FÖR GEOVETENSKAPER
D E P A R T M E N T O F E A R T H S C I E N C E S
Examensarbete vid Institutionen för geovetenskaper
Degree Project at the Department of Earth Sciences
ISSN 1650-6553 Nr 485
Analysis of Seismic Data Acquired in the Hverahlíð Geothermal Area
Agnieszka Stoch
The work for this thesis was carried in cooperation with Reykjavik Energy.
ISSN 1650-6553
Copyright © Agnieszka Stoch
Published at Department of Earth Sciences, Uppsala University (www.geo.uu.se), Uppsala, 2020
Abstract
Analysis of Seismic Data Acquired in the Hverahlíð Geothermal Area Agnieszka Stoch
Volcanic rifting environments, such as in Iceland, are challenging for conventional seismic reflection methods using active surface seismic sources. This study demonstrates the potential of a novel technique, called Virtual Reflection Seismic Profiling (VRSP) for imaging reflections in geothermal regions, like Hverahlíð, where a dense seismic array recorded a number of local microearthquakes for cross-correlation. Uppsala University, in collaboration with Reykjavik Energy, recorded seismicity in Hverahlíð using both seismometers and geophones. Acquired data were processed using the VRSP method, which applies seismic interferometry only to selected events, in this thesis local microearthquakes. Cross-correlation of the signal from a microearthquake recorded at one of the stations, which acts as a virtual source, with a ghost reflection recorded by the remaining receivers, produces a virtual shot gather. Stacking each station’s result, for all available events, and following a conventional multichannel processing sequence resulted in two stacked seismic images. Potential reflections observed in the obtained sections could be associated with major feed zones identified in the area by the borehole measurements. Eight dynamite explosions were processed with a conventional seismic reflection method, as a complementary source to the microearthquakes. In the produced stacked seismic section two potential reflections could be observed. Results from both passive and active datasets were 3D visualised to verify whether the reflections correspond to each other between sections.
Two horizons were traced throughout all three stacked sections. One more interface appeared on two images obtained from processing the passive data. This study shows promising results for using natural sources to image the subsurface in a challenging environment.
Keywords: Hverahlíð geothermal area, reflection seismic, virtual reflection seismic profiling, virtual shot gather, cross-correlation, autocorrelation
Degree Project E in Geophysics, 1GE029, 30 credits Supervisor: Christopher Juhlin
Department of Earth Sciences, Uppsala University, Villavägen 16, SE-752 36 Uppsala (www.geo.uu.se) ISSN 1650-6553, Examensarbete vid Institutionen för geovetenskaper, No. 485, 2020
The whole document is available at www.diva-portal.org
Popular science summary
Analysis of Seismic Data Acquired in the Hverahlíð Geothermal Area Agnieszka Stoch
Iceland, located on the Mid-Atlantic Ridge, characterised by rift zones with active faulting, and above a mantle plume, is a challenging environment for imaging the subsurface with the conventional seismic reflection method. Artificially generated seismic waves at the surface have difficulty in penetrating the lithology, because they are reflected at shallow depths due to horizontal lava flows. Strong, irregular scattering effects cause additional difficulty. This study demonstrates the potential of a novel technique, called Virtual Reflection Seismic Profiling (VRSP) for imaging reflections in Hverahlíð geothermal area, SW Iceland. This method uses selected events with clear signals to image the subsurface, in this project it were local microearthquakes. Upgoing waves created at a subsurface source, such as a microearthquake, are more advantageous to use than artificially generated seismic waves at the surface, since they are distributed at various depths. VRSP has also the financial advantage of using natural sources or ambient noise instead of expensive artificial surface sources.
Uppsala University, in collaboration with Reykjavik Energy, deployed a dense seismic array to continuously record seismicity at Hverahlíð geothermal reservoir, in the southern part of Hengill central volcano, SW Iceland. Additionally, eight dynamite explosions were made as a complementary source to the microearthquakes. In this thesis, virtual shot gathers were produced by cross-correlating the signal from a microearthquake recorded at one of the stations, which acted as a virtual source, with a ghost reflection recorded by the remaining receivers. For all available events, results from each station were stacked and after following a conventional multichannel processing sequence two seismic images were produced with potential reflections. The dynamite explosions were processed using a conventional seismic reflection method and in the produced image two reflections can be observed. Results from both were 3D visualised to verify whether the reflections correspond to each other between sections. Two reflections were traced throughout all three sections. One more interface appeared on two images obtained from processing the microearthquake data. Potential reflections observed in the obtained sections could be associated with intrusions, magmatic bodies, or magma related fluids within a fissure swarm. This study shows promising results for using natural sources to image the subsurface in a challenging environment, such as in Iceland.
Keywords: Hverahlíð geothermal area, reflection seismic, virtual reflection seismic profiling, virtual shot gather, cross-correlation, autocorrelation
Degree Project E in Geophysics, 1GE029, 30 credits Supervisor: Christopher Juhlin
Department of Earth Sciences, Uppsala University, Villavägen 16, SE-752 36 Uppsala (www.geo.uu.se) ISSN 1650-6553, Examensarbete vid Institutionen för geovetenskaper, No. 485, 2020
The whole document is available at www.diva-portal.org
Table of Contents
1. Introduction ... 1
2. Motivation and objectives ... 2
3. Background ... 3
3.1 Study area: Hverahlíð, Iceland ... 3
3.2 Previous research ... 4
4. Data acquisition ... 5
5. Methodology ... 7
5.1 Virtual Reflection Seismic Profiling ... 7
5.2 Used software ... 8
5.3 Pre-processing ... 9
5.3.1 Dynamite shots ... 9
5.3.2 Microearthquakes ... 9
5.4 Processing the dynamite shots... 11
5.4.1 Survey geometry and elevation statics ... 11
5.4.2 Data filtering ... 12
5.4.3 NMO correction and stack ... 14
5.5 Processing the microearthquakes ... 14
5.5.1 Elevation statics and data filtering ... 15
5.5.2 Autocorrelation ... 16
5.5.3 Cross-correlation ... 16
5.5.4 Survey geometry ... 17
5.5.5 NMO correction and stack ... 17
6. Results ... 19
6.1 Dynamite shots ... 19
6.2 Microearthquakes ... 20
6.2.1 Autocorrelation ... 20
6.2.2 Cross-correlation ... 21
6.3 Comparison ... 25
7. Discussion ... 25
8. Conclusions ... 27
9. Acknowledgments ... 28
10. References ... 29
1
1. Introduction
Iceland’s unique location on the Mid-Atlantic Ridge, characterised by rift zones with active faulting, and above a mantle plume, creates a challenging environment for imaging with the conventional seismic reflection method due to strong, irregular scattering effects. Semi-horizontal lava flows at relatively shallow depths cause additional difficulty, because of their high impedance contrasts. These factors make surface seismic sources less suitable for exploration of geothermal reservoirs in the Icelandic volcanic rifting environment (Gudmundsson, 2019). Recent improvements in seismic analysis, as well as technological advances in instrumentation make it possible to extract information about the subsurface with a coherently sampled wavefield. To do so, Virtual Reflection Seismic Profiling is applied. This technique, based on seismic interferometry, instead of using the measured ambient noise field, includes only selected events with clear signals. For this project local microearthquakes were processed to produce seismic sections. Artificially generated seismic waves at the surface have difficulty in penetrating the lithology, because they are reflected at shallow depths due to semi-horizontal lava flows. Upgoing waves created at a subsurface source, such as a microearthquake, might be more advantageous to use since they are distributed at various depths. Furthermore, events of magnitude M=
0, 1, 2 correspond to confined explosions of 15, 250 and 5000 kg of dynamite, respectively (Chael, 2009). Deploying dense arrays for an extended period of time allows to coherently sample the wavefield in the study area.
A pilot study was conducted at Hverahlíð geothermal reservoir, in the southern part of Hengill central volcano, SW Iceland. Uppsala University, in collaboration with Reykjavik Energy, deployed a dense seismic array consisting of both seismometers and geophones to continuously record seismicity.
Obtained data were used for noise studies and cross-correlation with stacking. Additionally, eight dynamite explosions were made as a complementary source to the microearthquakes.
This thesis focuses on the passive data recorded for one month and the seismic data obtained from the detonations. Analysis of the passive data included extracting 2-minute-long traces, corresponding to the time window of the microearthquakes located by the COSEISMIC network, and creating Virtual Shot Gathers using the Virtual Reflection Seismic Profiling method. Stacked images were compared to artificial source results processed using a conventional seismic reflection method.
2
2. Motivation and objectives
The aim of this project was to test how the Virtual Reflection Seismic Profiling method works in the challenging environment of the Hverahlíð area and assessing the viability of this approach for further exploration of the geothermal reservoir.
The aim of this thesis was to analyse the seismic data acquired in the Hverahlíð geothermal area.
The main objectives were:
• analysis of dynamite shot gathers for identification of phases, application of corrections and stacking,
• using passive data to generate Virtual Shot Gathers by autocorrelation and cross-correlation of selected microearthquakes,
• applying corrections to Virtual Shot Gathers and stacking for each station.
3
3. Background
3.1 Study area: Hverahlíð, Iceland
Iceland is located above a mantle plume on the Mid-Atlantic Ridge, which separates the American and Eurasian plates. The ridge separating the island is characterised by rift zones with active faulting and consists of two segments: in the south - the Reykjanes Ridge, and in the north – the Kolbeinsey Ridge.
Volcanism on the plate boundary in Iceland extends from southwest to the north (Denk et al., 2011).
Hverahlíð is located in the Hengill region, southwest Iceland, which is one of the most extensive geothermal areas in Iceland (Hardarson et al., 2009). This currently active zone is located at the triple junction of the South Iceland Seismic Zone (SISZ), the Reykjanes Peninsula (RP) and the Western Volcanic Zone (WVZ) as shown in Figure 1. This system is comprised of fissure swarms, volcanoes and grabens, into which lavas and hyaloclastites gather, creating highlands in central Hengill (Franzson et al., 2010). The Hengill region is dominated by NE-SW faults intersected by easterly striking features.
Three fissure eruptions that occurred 9000, 5000 and 2000 years ago are considered postglacial volcanism, and affected the permeability of the area (Harðarson, 2014). The two younger active fissure swarms, striking NE-SW, are the main outflow channel for the geothermal exploration at the Nesjavellir field and at Hellisheidi (Harðarson, 2014). Groundwater is heated up by intrusions in the roots of the volcano, resulting in an upflow caused by buoyancy forces. Those upflow zones, situated underneath the Hengill region, explain the geothermal activity in the area (Nielsson & Franzson, 2010). The volcano-tectonic activity in the Hengill area is related to rifting, along with induced seismicity caused by drilling and fluid injection into the geothermal fields (Li et al., 2019).
The Hverahlid system is situated in the southern part of the Hengill region, outside of the highlands, but within the fissure swarm. It is an example of an ice-confined lava margin that formed during eruptions underneath ice, close to the edge of a glacier (Sæmundsson, 2010). West of Mt. Skálafell lies a deep crater at the peak of the lava shield that formed during the last glacial period. The highest part of the main margin is located northwest, close to the crater. The shield lava flowed down to the east, creating a lowland where crops of hyaloclastite breccia can be seen within lava scoured by ice (Sæmundsson, 2010). This high-temperature system is comprised of predominantly two rock types:
hyaloclastites created in subglacial eruptions, producing highlands, and lava series formed during interglacial periods, accumulated in the lowlands (Helgadóttir et al., 2010). Basaltic hyaloclastites form while magma cools and piles up at the base of the glacier and have relatively low permeability despite their high porosity (Helgadóttir et al., 2010).
4
Figure 1. Geological map of Iceland. Volcanic zones: Reykjanes Ridge (RR), Reykjanes Peninsula (RP), Western Volcanic Zone (WVZ), Mid-Iceland Volcanic Zone (MVZ), Northern Volcanic Zone (NVZ), Eastern Volcanic Zone (EVZ), Vestmanna Islands (VI), South Iceland Seismic Zone (SISZ), Tjörnes Fracture Zone (TFZ). Red dots mark high-temperature sectors, orange circle - The Hengill volcanic system (from Hardarson et al., 2010).
3.2 Previous research
Extensive studies have been conducted in the Hengill region to outline the geothermal anomaly. Several geophysical methods have been used in the area, including aeromagnetic surveys (Bjornsson & Hersir, 1981), gravimetry (Thorbergsson et al., 1984, Árnason et al., 1987), seismic refraction (Pálmason, 1971) and passive seismic (Foulger, 1984). Most information about the geothermal reservoir was obtained through resistivity studies (Franzson et al. 2010). In 1980, Hersir used the Schlumberger configuration and in 1981, dipole-dipole, together with Björnsson (Hersir, 1980; Bjornsson & Hersir, 1981).
Throughout the years numerous transient electromagnetic soundings were made, as well as magnetotelluric surveys (Árnason et al., 2010).
Geophysical exploration, geological mapping and geochemical sampling indicated that the low resistivity structure covers approximately 110 km2. The Hengill geothermal field has been utilized in the Nesjavellir and Hellisheidi areas to produce hot water and generate electricity for several years, Hverahlid and Bitra are being developed currently. Nevertheless, continuous studies show that there are several geothermal fields that could be exploited in the future (Hersir et al., 2009).
5
4. Data acquisition
Measurements in the Hverahlid geothermal area were made with 13 short period seismometers (8 Lennartz 5s sensors and 5 Lennartz 1s sensors, all accompanied by Reftek digitizers), 2 three-component broadband Reftek seismometers, 2 three-component broadband Guralp seismometers and 30 vertical 5 Hz SmartSolo geophones with autonomous in-built recording. Geometry of the deployed seismic array is presented in Figure 2. Acquisition of data started on the 2nd of July 2019 and continued until 19th of September 2019. This project focuses on measurements made between the 7th of August and 2nd of September. Spacing between stations operating during this time varied from 45m to 160m. Data acquired by the Guralp station SUDA had to be excluded from the analysis due to GPS clock failure.
Figure 2. Survey area showing the location of the three-component intermediate broadband and bandwidth seismometers (orange triangles), and vertical 5 Hz SmartSolo geophones (blue triangles) in the Hverahlid geothermal area. Red lines mark an apparent fault, green circles indicate boreholes.
Two types of sources were used for data processing: dynamite explosions and microearthquakes.
Artificial source detonations took place on the 28th of August in 8 locations presented in Figure 3A.
Each shot hole was 76 mm wide and cased for the uppermost 3 m. Seven shot holes were 10 m deep and one was 8 m deep. Five of the shots were made using 5 kg of dynamite, two with 2.5 kg and one with 3
6
kg. The weight was reduced in three wells due to close location to the deep production wells. During dynamite explosions, the 40 stations shown in Figure 3A were operating.
A search through the catalogue of automatically located events from the SIL seismic system of the Icelandic Meteorological Office (IMO) for the project’s timeframe listed 7 microearthquakes within 2 km radius from the centre of the array. In order to obtain better results, it was necessary to use more events. The COSEISMIC network on Hellisheidi operated by ISOR provided locations of 39 events within the same radius. The earthquakes extracted from the SIL stations were also recorded with the COSEISMIC seismometers, but with slightly different locations. Since the COSEISMIC network is denser and operates more stations in the study area, events provided by ISOR were analysed (Figure 3B). Magnitudes of the recorded events varied from 0.05 to 2.2.
Figure 3. (A) Geometry of the seismic array recording during the dynamite explosions (yellow crosses) on the 28th of August 2019. (B) Map of the receivers recording between the 2nd of July and the 2nd of September 2019 with plotted locations of the microearthquakes (red dots) used for data processing. SmartSolo geophones are represented by blue triangles and seismometers by orange triangles.
7
5. Methodology
5.1 Virtual Reflection Seismic Profiling
Seismic signals can be produced by a range of sources: from impact devices such as a sledgehammer or weight drop, through explosives, airguns to vibrators. Data are collected by geophones that convert seismic energy to electrical voltage. The seismic reflection method is based on processing the back- scattered energy from sub-surface discontinuities. The signal generated on the surface, at known time, radiates through the earth. It is reflected back to the receiver when a contrast in the acoustic impedance appears across the interface.
Seismic interferometry, introduced by (Claerbout 1968), is a technique that allows imaging the subsurface using the ambient noise field. This method, also known as daylight imaging, assumes sources to be randomly distributed underneath the reflectors. Mathematical investigations conducted by Wapenaar et al., (2002) and Wapenaar (2004) broaden the use of this concept, by showing that the Green’s function of any heterogeneous medium can be recovered from the cross-correlation of wavefields recorded at two receivers, were one of them acts as a virtual source position.
This prominent and powerful tool can be especially used in regions with complex terrain, or where the monitoring is conducted permanently, as for reservoir production (Holzner et al., 2005) and in CO2
sequestration surveillance (Xu et al., 2012). Several projects have been based on this method. Yu and Schuster (2006) performed analysis where the source of seismic energy was a drill bit. Draganov et al., (2009) used it to retrieve reflection profiles from surface-wave noise. Seismic interferometry can also be applied at a larger scale. It was used by Nakata et al., (2014) to image the structure of the lithosphere in Tibet and by Gorbatov et al., (2013) to retrieve Moho reflections. This project follows the application of seismic interferometry to local microearthquakes and is referred to as Virtual Reflection Seismic Profiling, based on the study conducted by Kim et al., (2017) for the magma reflection imaging in Krafla, Iceland.
In order to retrieve the Green’s function by cross-correlation, a few assumptions need to be made.
Principally, the ambient noise sources must be uncorrelated and randomly distributed, and the measurements should be made over a long enough period to record a sufficient number of events (Cheraghi et al., 2015). Equation 1 is used to retrieve the seismic impulse response.
{G
p,q(x
A,x
B,t)+ G
p,q(x
A,x
B,-t)}
*a(t) ≈
∑ 𝑣𝑁 𝑝𝑖𝑖
(x
A,-t)
*𝑣𝑞𝑖(x
B,t)
(1)The left-hand side of the equation represents the convolution of the noise source time function autocorrelation a(t) with the Green’s function. Gp,q(xA,xB,t) represents the casual part and Gp,q(xA,xB,-t) the acasual (time-reversed) part of the Green’s tensors for the receivers positions xA and xB. The right-
8
hand side of the equation represents the sum over N available microearthquake records. It is carried out on the cross-correlation of the i-th particle velocity component recorded in the p and q directions.
Virtual shot gathers are created by cross-correlating the signal recorded at a station, which acts as a virtual source, with remaining receivers. This redatums ghost reflections to surface source locations.
Figure 4 presents a conceptualization of this method. Receiver A and B, located at the surface, measure a wavefield created by a microearthquake underneath. The signal recorded at station B is the ghost reflection, which is energy reflected from the surface travelling down to an arbitrary reflector and back to the receiver. Correlating the two signals from A and B produces a seismic impulse response at B as if it came from the virtual source at A. Virtual shot gathers are created by selecting one station at a time and cross-correlating a recorded trace with the remaining receivers’ response. This produces correlation panels at all receiver positions for each microearthquake recording available. The autocorrelation of the noise source time function is compensated by deconvolving results with a wavelet extracted from the autocorrelation trace or by normalising traces so that the source energy is equally weighted (Draganov et al., 2009; Cheraghi et al., 2015).
Figure 4. Conceptualization of the Virtual Reflection Seismic Profiling method. Cross-correlation of the signal from a microearthquake recorded at the station A with a ghost reflection recorded by the station B, produces a seismic impulse response at B as if it came from the virtual source at A, a virtual shot gather (modified from Kim et al., 2017).
5.2 Used software
List of programmes used in this thesis:
• GLOBE ClaritasTM Seismic Processing Software
• OpendTect – Open Seismic Interpretation Software
• QGIS 3.4 Madeira
• SAC v101.6a software distributed by IRIS DMC
• the Seismic Waveform Tool (SWFT), a code designed by Lawrence Livermore National Laboratory
• SeisWare Coordinate Converter
9
5.3 Pre-processing
Data for this project had to be examined and prepared for further processing. Station elevations were obtained from the LMÍ Hæðarlíkan 2016 DEM, available through the National Land Survey of Iceland- Landmælingar Íslands portal (https://www.lmi.is). All data were converted using SeisWare, from ISN93/ Lambert 1993 coordinate system into UTM zone 27N to match the GLOBE ClaritasTM software requirements.
5.3.1 Dynamite shots
Thirty geophones and ten seismometers recorded dynamite detonations on the 28th of August (Figure 3A). The SEGY data from the SmartSolo instruments consist of 2-minute-long traces with 30001 samples in each record and a sampling rate of 4 ms.
SAC formatted data obtained from seismometers were 24 hours long with a sampling rate of 5 ms.
In order to create one seismic dataset containing data from both types of instruments, SAC files had to be pre-processed to fit the geophone data. At first, 2-minute-long traces corresponding to the time of the 8 dynamite shots were extracted using the Seismic Waveform Tool (SWFT), a code designed by Lawrence Livermore National Laboratory. Given a list of events the Seismic Data Processor module extracted selected traces for each seismometer. Output files were converted to SEGY format and resampled to 4 ms using GLOBE ClaritasTM Seismic Processing Software. Merging the traces with the SmartSolo data showed a time shift in the extracted events. To align the seismic data correctly, dynamite explosions had to be extracted once more. The same process was repeated using the SAC v101.6a software distributed by IRIS DMC. The produced traces had to be shifted down 18 s to account for the difference in UTC and GPS time. Ten seconds, containing each shot, were extracted for further processing. There was an evident discrepancy on all shots between data recorded by geophones and nearby seismometers, so traces were shifted to match P-wave arrivals on the closest SmartSolo geophone. Further analysis showed that traces from seismometers had reversed polarity compared to the SmartSolo stations, so all of them were flipped before processing.
5.3.2 Microearthquakes
The SEGY files from the SmartSolo geophones contained one day of passive data with 720 two-minute long traces. In order to extract data corresponding to the time of the microearthquakes, traces first had to be appropriately sorted. Afterwards, adequate hour, minute, and second values were added into the headers.
Seismometers deployed for the survey had gaps in the recordings, so available data content was determined by the time of the microearthquakes. Analysis showed that the number of files available for further processing varied from ten up to fifteen seismometers depending on the event (Figure 3B). Data
10
recorded by seismometers were extracted, converted, resampled and had their polarity reversed in the same manner as for the dynamite shots. Two-minute traces were cut down to five seconds, to contain only recorded microearthquakes. An example of a microearthquake recording is presented in Figure 5.
Figure 5. Microearthquake recorded by the deployed seismic array, on the 24th of August 2019 at 5:54:28.
11
5.4 Processing the dynamite shots
Artificial source recordings were processed using a conventional seismic reflection procedure. An outline of the processing flow is presented in Figure 6.
Figure 6. Outline of the dynamite shots processing sequence.
5.4.1 Survey geometry and elevation statics
Creating an image of the sub-surface requires a geometry database, which will relate SHOTID and CHANNEL numbers, stored in the recorded trace headers, to surface locations and elevations. Point coordinates of each station and dynamite shot were added into the headers. That information allowed the calculation of offsets between source and receivers. Then, a geometry database was created.
Following that, traces were grouped into common depth point (CDP) gathers in two variants: 2D with 32 bins and 3D with 8 bins.
To account for near surface velocity variations and changes in the elevation of sources and receivers, refraction static corrections had to be calculated. This approach vertically shifts traces to a flat datum.
It required picking the earliest seismic arrivals on all traces for each shot. Having that information along with the geometry database, enabled use of the Refstat application, part of the GLOBE ClaritasTM software. Refstat, based on the geometry information, produced an initial model that uses ray tracing of estimated theoretical times of first arrivals. Results were compared with the picked first breaks. After that, the difference between the two is minimised, by running the inversion iteratively. The output file consists of surface-consistent residual time delays.
Pre-processing Geometry
database Static correction Bandpass filter
SPEW - spectral
whitening AGC (200 ms) FDFILT AGC (500 ms)
SMUTE CDP sort NMO correction Stack
FXDECON AGC (300 ms)
12
To achieve the best result, it was necessary to compare the influence of different static corrections on the raw data. Using trial parameters, modifying the velocities of the initial model, and running inversions for 2D and 3D geometry produced numerous static shift files. The most satisfactory static correction (Figure 7B) was obtained through 6 iterations with initial velocities of 900m/s and 3500 m/s in the starting model and using 3D geometry.
Figure 7. Example of a raw shot gather (A), and post application of the static correction (B).
5.4.2 Data filtering
Prior to velocity analysis data had to be cleared of noise to enhance potential reflections and emphasise the body waves. Based on the analysis of the amplitude spectrum graph (Figure 8), dynamite shots were filtered using a bandpass filter (1-6-80-120 Hz) to attenuate the contribution of the surface waves, which have high amplitudes and low frequency. Following that, a spectral whitening function was applied and automatic gain control (AGC) with a 200 ms window. To make the overall results clearer and remove the over-whitening, a zero-phase filter with a quasi-trapezoidal amplitude spectrum was used (3-6-20- 30 Hz). Then AGC was repeated using a 500 ms window. Traces were shifted 700 ms, to position the first arrivals approximately at 200 ms two-way traveltime. Afterwards, to prepare the shot gathers for stacking, the first arrival energy was surgically muted. Four out of eight processed shots are shown in Figure 9 with red arrows indicating potential reflections.
13
Figure 8. Amplitude spectrum graph for each recorded dynamite shot.
Figure 9. Four examples of dynamite shot gathers with potential reflections marked by red arrows.
14 5.4.3 NMO correction and stack
Distance between the receiver and the source influences the arrival time of a reflection. Increasing offset creates a time delay. Removing this effect is done by applying the normal move-out correction (NMO) to the CDP gathers. Using too high or too low velocity will not flatten the reflection. Having that in mind, several constant velocities ranging from 4000 m/s to 7000 m/s at 100 m/s intervals were tested.
A file with time and velocity pairs was constructed after analysing the effect of each constant velocity correction on the stacked data (Figure 10). Another NMO correction was created in the Claritas Velocity Analysis application. Using all produced corrections, dynamite shot data were then stacked and compared. A post-stack FX deconvolution filter was applied to attenuate random noise and improve event continuity. This complex Wiener deconvolution is used on the seismic section transformed to FX space. Filtering for each frequency is performed in the X direction. Afterwards, the section is converted back to TX space. Additionally, AGC with a 300 ms window was applied to enhance the reflections.
Figure 10. Velocity model for the dynamite shots created through analysis of the influence ofthe normal move- out correction with constant velocities ranging from 4000 m/s to 7000 m/s at 100 m/s intervals.
5.5 Processing the microearthquakes
Local microearthquakes recorded by the seismic array with magnitudes ranging from 0.05 to 2.2 were analysed as energy sources. Only the vertical component data were considered in order to minimise the S wave influence. Virtual Reflection Seismic Profiling was used to produce stacked seismic images of the subsurface. Figure 11 presents the simplified processing sequence.
15
Figure 11. Simplified processing sequence used for the passive data.
5.5.1 Elevation statics and data filtering
Static corrections created in the Refstat application, based on first breaks picked on microearthquakes and created for this data set’s geometry, did not bring satisfactory results. Since 40 stations used for the dynamite shots overlapped with the passive data processing, the static corrections used previously for the active source were again applied with slight modification. Time shifts for the missing 5 stations were approximated to match neighbouring receivers. The adjusted file was then used for the static correction.
To reduce contributions of the surface waves and simultaneously emphasise body waves, microearthquake recordings were bandpass filtered between 6 and 80 Hz. Next, a spectral whitening function was applied to the whole data set. Following that, DC bias was removed. Traces then were normalised, to ensure that their amplitudes have the same root-mean-square amplitude level for each microearthquake. The normalising process compensated for the effect of the autocorrelation of the noise source time function a(t) in equation 1. Prior to autocorrelation and cross-correlation AGC with a 200 ms window was applied.
16 5.5.2 Autocorrelation
Three-second long autocorrelation sections were calculated for all receivers, within each microearthquake record. Then, each stations’ results were stacked. To emphasize coherent arrivals, a zero-phase filter with a quasi-trapezoidal amplitude spectrum (3-6-20-30 Hz) was applied, together with AGC using a 200 ms window.
5.5.3 Cross-correlation
Ghost reflections were redatumed to surface source locations by cross-correlating the signal recorded at a selected station, representing a virtual source, with the remaining receivers’ response. This procedure required use of the VIBCORR module, where an auxiliary trace of the station acting as a virtual source was cross-correlated with all other available receivers within each microearthquake record. Repeating that procedure for one station at the time, for each event, produced 3-second-long correlation panels.
Afterwards, all available microearthquakes were stacked for each station, giving 45 virtual shot gathers.
To check results at negative time cross-correlation was repeated with the pilot trace shifted 500 ms up.
Figure 12 presents a virtual shot gather using receiver 107 as the virtual source.
Figure 12. Virtual shot gather created with receiver 107 as a virtual source. Red arrows indicate a potential reflection.
17 5.5.4 Survey geometry
The passive data processing required creating a geometry database for the microearthquakes.
Coordinates of each virtual source receiver were added into the headers from a text file. Following that, traces were grouped into 32 CDP bins along the longest continuous part of the receiver line as presented in Figure 13A (line 1). Additionally, a geometry database with 26 CDP bins (line 2) was created perpendicular to the other line (Figure 13B).
Figure 13. (A) Positions of the 32 CDP bins (line 1) created along the longest continuous part of the receiver line.
(B) Positions of the 26 CDP bins (line 2) created perpendicular to the line 1.
5.5.5 NMO correction and stack
Prior to velocity analysis for the passive data, a zero-phase filter with a quasi-trapezoidal amplitude spectrum (3-6-20-30 Hz) was applied, to clear virtual shot gathers of noise. Then, AGC was used with a 500 ms window.
Several NMO files with constant velocities ranging from 4000 m/s to 7000 m/s at 200 m/s intervals were tested. Arbitrary time-velocity pairs, results from the Claritas Velocity Analysis application and automatically created files in VELSPEC and NMOPICK modules were also tried. The influence of different constant normal moveout velocities on an exemplary virtual shot gather are presented in Figure 14 together with a result from the Claritas Velocity Analysis application. Results were compared after sorting the data into common depth point gathers and stacking. To improve event continuity and attenuate random noise a post-stack FX deconvolution filter was applied.
Additionally, this procedure was applied to microearthquakes located south from the dynamite shots to check if the results improved. The obtained stacked seismic section exhibited poorer resolution and no enhancement of the reflections, so it was not considered for further analysis.
18
Figure 14. Influence of different NMO corrections with constant velocities on a virtual shot gather with receiver 105 as a virtual source: (A) 4000 m/s, (B) 4800 m/s, (C) 5600 m/s, (D) 6400 m/s, (E) 7000 m/s and the result from the Claritas Velocity Analysis application (F). The green arrow indicates the direct wave, the red arrow points to the potential reflection.
19
6. Results
6.1 Dynamite shots
A seismic array, consisting of both geophones and seismometers, recorded eight dynamite explosions as energy sources. Data were analysed and processed using a conventional seismic reflection procedure to obtain a subsurface image. Due to uncertainties in the velocity model, results were not converted from time to depth.
Figure 15 shows a stacked seismic section produced using the NMO correction velocity presented in Figure 10. There are visible parts of a reflection, labelled C1, between 400 ms and 600 ms two-way traveltime, dipping from high to low number CDPs. A second marked reflection, C3, appears at approximately 700 ms two-way traveltime. The discontinuity in the C1 and C3 reflections, between 114 and 120 CDP, could be explained by the apparent fault marked in Figure 2.
Figure 15. Stacked seismic section obtained from processing the dynamite shot data. The C1 reflection is marked with green arrows and C3 with yellow arrows.
20
6.2 Microearthquakes
6.2.1 Autocorrelation
Autocorrelations of each trace within 39 microearthquake records were calculated and stacked for each receiver to produce the seismic section presented in Figure 16. The generated image imitates “an unmigrated stacked, single-fold, zero source-receiver offset survey” (Kim et al., 2017). If the microearthquakes were located underneath each station or were dispersed stochastically, the autocorrelation results could be used as a base for identifying virtual reflections.
A few scattered, visible coherent arrivals are marked with yellow arrows in the produced autocorrelation panel (Figure 16). Since there was a limited number of available events, that were not appropriately distributed, the visible coherency does not provide enough information to positively identify virtual reflections on this seismic section.
Figure 16. Seismic section produced by stacking each receiver’s autocorrelation results for all microearthquake recordings. The yellow arrows indicate coherent arrivals.
21 6.2.2 Cross-correlation
The same microearthquake recordings, used for the autocorrelation procedure, were used to produce Virtual Shot Gathers with the Virtual Reflection Seismic Profiling method. Initially, data were stacked along the longest continuous part of the receiver line (Figure 13A, line 1). On the processed seismic section (Figure 17), a coherent event labelled C1 can be observed between 300 ms to 400 ms, dipping from high to low number CDPs. To verify whether C1 is a potential reflection or a processing artifact, a perpendicular CDP line was created (Figure 13B, line 2). Data were processed and stacked with different trial velocities in the same manner as the previous line. In the obtained stacked image (Figure 18), a corresponding event is visible in the same time interval, dipping from low to high number CDPs.
This would indicate that C1 could be a reflection. To confirm that this interface is consistent on both seismic sections, the crossing lines were 3D visualised using Opendtect software (Figure 19). This seismic interpretation software allowed for tracking of each reflection in the subsurface image and extrapolating it to the adjacent line. The C1 horizon was traced on both CDP profiles with a green line.
In Figure 17, another visible reflection, C2, at 500 ms is shown. At the perpendicular line in Figure 18, a corresponding horizon, potentially C2, is also visible between CDP 108 and 114. On both sections in the 3D view (Figure 19), C2 is traced with a blue line. Possible parts of the C3 interface, identified in the dynamite shot section (Figure 15), are marked with a yellow line.
22
Figure 17. Stacked seismic section along the longest continuous part of the receiver line (line 1) obtained through the Virtual Reflection Seismic Profiling method. Reflection C1 is marked with green arrows, C2 with blue arrows and C3 with yellow arrows.
23
Figure 18. Stacked seismic section, line 2, perpendicular to line 1, obtained through the Virtual Reflection Seismic Profiling method. Reflection C1 is marked with a green arrow, the visible parts of the C2 and C3 reflections with a blue and yellow arrow, respectively.
24
Figure 19. 3D visualization of the cross-correlation results. Line 1 was created along the longest continuous part of the receiver’s line, line 2 is perpendicular to it. The green horizon traces the C1 reflection, the blue horizon marks the C2 reflection and the yellow horizon traces the C3 reflection.
25
6.3 Comparison
The unmigrated stacked section obtained from processing the recorded dynamite explosions was compared with the cross-correlation results. Stacked CDP lines 1 and 2 are visualised in Figure 20, together with the dynamite result (line 3). Reflection C1 (green line) was tracked on all stacked images, as well as parts of C3 (yellow line). The C2 interface is only visible on lines 1 and 2.
7. Discussion
Acquiring detailed subsurface images with the reflection seismic method is difficult in the challenging Icelandic environment. High impedance contrast of the semi-horizontal lava flow at shallow depths, irregular scattering effects and high seismic attenuation considerably affect the penetration of seismic waves (Gudmundsson, 2019).
Despite poor resolution and an uncertain velocity model, some potential reflections can be observed in the obtained subsurface images of the Hverahlíð geothermal area (Figures 19 and 20). These potential interfaces could be associated with a magmatic body or magma related fluids including entrapped brines, steam, CO2 and SO2 (Kim et al., 2017). Observed reflectivity could also be related to intrusions or to the transition from hyaloclastites and Quaternary lava flows to Tertiary lava flows and dikes (Wagner et al., 2019). To determine whether the identified interfaces are real or are a processing artifact, dynamite shot gathers (Figure 9), as well as virtual shot gathers (Figure 12) were examined. Parts of the reflections were identified in both active and passive datasets reducing the likelihood that the sub-horizontal events are artifacts.
Kim et al., (2017) showed the benefit of correlating obtained stacked seismic sections with a well- logging measurement to confirm that a visible reflection was from a known magmatic body. Based on that a synthetic VRSP section could be computed.
The seismic array deployed for this study did not directly overlap with any boreholes, so the obtained sections could not be correlated with known lithology. However, measurements made in four boreholes drilled north from the seismic array (Figure 2), indicated the presence of three feed zones at depth ranges:
700-900 m, 1000-1300 m and 1800-2000 m (Okbatsion, 2010). For line 1 the estimated depth to reflection C1 is 900 m, C2 is 1.2 km and C3 is 1.8 km, assuming a velocity of 5000 m/s. Line 2 shows C1 dipping north-west and intersecting line 3 at 1.2 km depth. Analysis of the active and passive data indicates that the observed reflectivity could correspond to the existing feed zones. A more detailed study with a 3D recording set-up, deployed for a longer time span to record more microearthquakes underneath the seismic array could identify and correctly position the reflections in order to rule out the possibility of the results being a processing artifact. Borehole measurements with sonic logs and/or borehole seismic could provide necessary information for building an appropriate velocity model.
Additionally, applying seismic interferometry with multidimensional deconvolution to coda waves, if implemented, could provide an additional result to validate the reflections.
26
Figure 20. 3D visualization of the cross-correlation results (line 1 and 2) compared with the image obtained from dynamite shots (line 3). The green horizon traces the C1 reflection, the blue horizon marks the C2 reflection and the yellow horizon marks the C3 reflection.
27
Although VRSP is a promising method, there are a few issues that must be considered. If the distribution of microearthquakes is not uniform, in regard to the azimuth, the cross-correlation can produce artifacts. Other complications can arise from polarity changes in the focal mechanism, converted phases of the seismic waves and S wave contamination affecting the cross-correlation function.
8. Conclusions
This study aimed to image the subsurface beneath the Hverahlíð geothermal area and test how the Virtual Reflection Seismic Profiling method works in this challenging environment. A dense seismic array recorded dynamite explosions and local microearthquakes and the obtained datasets were processed using a conventional seismic reflection method and Virtual Reflection Seismic Profiling, respectively.
The processing workflow used in this study produced a stacked image from the dynamite shot data with a visible reflection, C1 at approximately 500 ms two-way traveltime, and C3, at 700 ms two-way traveltime (Figure 15). Selected local microearthquakes recorded by the deployed array and processed using the VRSP method resulted in two images of the subsurface also with visible reflections, C1 and C2, and parts of C3 (Figure 19). The identified interfaces may be associated with three major feed zones, related to intrusions, magmatic bodies, or magma related fluids within a fissure swarm. A 3D recording set-up with a dense seismic array intersecting a borehole and deployed for an extended period of time could be used to further confirm the presence of the reflections and verify this interpretation.
VRSP has the financial advantage of using natural sources or ambient noise instead of artificial surface sources, typically used in conventional seismic reflection surveys. However, the active source data served two important purposes, (1) it provided the static corrections that improved the coherency of the reflections in the passive processing and (2) it provided an additional section with potential reflections that could be compared to the VRSP processing.
28
9. Acknowledgments
First and foremost, I would like to express my sincere appreciation to my supervisor, Christopher Juhlin, for his guidance throughout this project. His support, advice and feedback were invaluable for this research to be completed. I am grateful to have been given the opportunity by Chris to join the field works in Orebro and Kiruna, where I gained hands on experience using different geophysical methods.
I would like to express my sincere gratitude to Reykjavik Energy who made this research possible.
My sincere thanks to Ólafur Gudmundsson for help with the data pre-processing. His insight and knowledge of the subject guided me throughout this project work. Special thanks to Alireza Malehmir for being a great teacher throughout courses, field works and who encouraged me to be a part of the organising committee of IGSC 2019. Thanks to him I joined the SEG Uppsala Student Chapter, through which I had an opportunity to attend SEG Annual Meeting in San Antonio, Texas. I am also grateful to Alex Hobé for valuable suggestions and help with SIL seismic system of the Icelandic Meteorological Office. I would also like to thank ISOR Icelandic GeoSurvey for providing me the data from the COSEISMIC network.
Heartfelt thanks to my “office-wife”, Salóme Jórunn Bernharðsdóttir, especially for help with converting the seismic data, but also for keeping me motivated throughout the process of writing, and for the fun, board-game evenings we shared. I also wish to thank Jack Hewitt for all the helpful advice and proof-reading my work. Finally, I would like to give my deepest appreciation to my parents, my sister, and my dear friends, who provided their unfailing support and relentless encouragement throughout the years. Thank you.
29
10. References
Árnason, K., Haraldsson, G.I., Johnsen, G.V., Þorbergsson, G., Hersir, G.P., Sæmundsson, K., Georgsson, L.S., Rögnvaldsson, S.Th., Snorrason, S.P. (1987). Nesjavellir-Ölkelduháls, surface exploration in 1986. Orkustofnun report in Icelandic, OS-87018/JHD-02, pp. 112.
Árnason, K., Eysteinsson, H. & Hersir, G.P. (2010). Joint 1D inversion of TEM and MT data and 3D inversion of MT data in the Hengill area, SW Iceland. Geothermics, vol. 39 (1), pp. 13–34.
Bjornsson, A. & Hersir, G.P. (1981). Geophysical reconnaissance study of the Hengill high-temperature geothermal area, SW-Iceland. Trans. - Geotherm. Resour. Counc.; (United States), vol. 5. Available at:
https://www.osti.gov/biblio/6741513 [2020-05-12]
Chael, E.P. (2009). Local magnitudes of small contained explosions. SANDIA REPORT 2009-7941.
DOI: https://doi.org/10.2172/973338
Cheraghi, S., Craven, J.A. & Bellefleur, G. (2015). Feasibility of virtual source reflection seismology using interferometry for mineral exploration: A test study in the Lalor Lake volcanogenic massive sulphide mining area, Manitoba, Canada: Interferometry for mineral exploration. Geophysical Prospecting, vol. 63 (4), pp. 833–848.
Claerbout, J.F. (1968). Synthesis of a layered medium from its acoustic transmission response.
Geophysics, vol. 33 (2), pp. 264–269.
Denk, T., Grímsson, F., Zetter, R. & Símonarson, L.A. (2011). Introduction to the Nature and Geology of Iceland. Late Cainozoic Floras of Iceland. Dordrecht: Springer Netherlands, pp. 1–29.
Draganov, D., Campman, X., Thorbecke, J., Verdel, A. & Wapenaar, K. (2009). Reflection images from ambient seismic noise. Geophysics, vol. 74 (5), pp. A63–A67.
Foulger, G.R. (1984). Seismological studies at the Hengill geothermal area, SW Iceland. PhD thesis, University of Durham, pp. 313.
Franzson, H., Gunnlaugsson, E., Árnason, K., Steingrímsson, B. & Harðarson, B.S. (2010). The Hengill Geothermal System, Conceptual Model and Thermal Evolution. Proceedings World
Geothermal Congress, Bali, Indonesia, 25-29 April 2010.
Gorbatov, A., Saygin, E. & Kennett, B.L.N. (2013). Crustal properties from seismic station autocorrelograms. Geophysical Journal International, vol. 192 (2), pp. 861–870.
Gudmundsson, O. (2019). Unpublished proposal for a pilot microseismic reflection survey at Hverahlíð.
Harðarson, B.S. (2014) Geothermal exploration of the hengill high-temperature field. Paper presented at Short Course IX on Exploration for Geothermal Resources, Kenya, pp. 8.
Hardarson, B.S., Einarsson, G.M., Franzson, H. & Gunnlaugsson, E. (2009). Volcano-Tectonic- Geothermal Interaction at the Hengill Triple Junction, SW Iceland. Transactions - Geothermal Resources Council, vol. 33, pp. 42-47.
Helgadóttir, H.M., Sn, S.Ó., Harðarson, B.S., Einarsson, G.M. & Franzson, H. (2010). Geology and Hydrothermal Alteration in the Reservoir of the Hellisheiði High Temperature System, SW-Iceland.
Proceedings World Geothermal Congress, Bali, Indonesia, 25-29 April 2010.
30
Hersir, G.P., Árnason, K. & Steingrímsson, B. (2009). Exploration and development of the Hengill geothermal field. Paper presented at Short Course on Exploration for Geothermal Resources, El Salvador. pp. 12.
Hersir, G.P. (1980). Electric and electromagnetic measurements across the Mid-Atlantic ridge in southwest Iceland, with special reference to the high temperature area of Hengill. Mag scient thesis, University of Aarhus, Denmark. pp.165.
Holzner, R., Eschle, P., Zürcher, H. & Dangel, S. (2005). Applying microtremor analysis to identify hydrocarbon reservoirs. First Break, vol. 23(5). pp. 41-46.
Kim, D., Brown, L.D., Árnason, K., Águstsson, K. & Blanck, H. (2017). Magma reflection imaging in Krafla, Iceland, using microearthquake sources: Magma Imaging with Microearthquakes. Journal of Geophysical Research: Solid Earth, vol. 122 (7), pp. 5228–5242.
Li, K.L., Abril, C., Gudmundsson, O. & Gudmundsson, G.B. (2019). Seismicity of the Hengill area, SW Iceland: Details revealed by catalog relocation and collapsing. Journal of Volcanology and Geothermal Research, vol. 376, pp. 15–26.
Nakata, N., Snieder, R. & Behm, M. (2014). Body-wave interferometry using regional earthquakes with multidimensional deconvolution after wavefield decomposition at free surface. Geophysical Journal International, vol. 199 (2), pp. 1125–1137.
Nielsson, S. & Franzson, H. (2010). Geology and Hydrothermal Alteration of the Hverahlid HT-System, SW-Iceland. Proceedings World Geothermal Congress, Bali, Indonesia, 25-29 April 2010.
Okbatsion, M.A. (2010). Analysis of temperature and pressure characteristics of the Hverahlíd geothermal field in the Hengill geothermal system, SW-Iceland. Geothermal Training Programme, Report 6, Reykjavík, Iceland.
Pálmason, G. (1971). Crustal structure of Iceland from explosion seismology. Societas Scientiarum Islandica, Rit, vol. 40, pp. 187.
Sæmundsson, K. (2010). 28. Hverahlíð - Glacier dammed lava flow. Iceland Geosurvey. Available:
https://en.isor.is/28-hverahlid-glacier-dammed-lava-flow [2020-05-12]
Wagner F, Tryggvason A, Gudmundsson Ó. (2019). Double-difference relocation and local-earthquake tomography of the Vp, Vs, and Vp/Vs structure in Hengill, SW Iceland. Manuscript in preparation.
Wapenaar, K. (2004). Retrieving the Elastodynamic Green’s Function of an Arbitrary Inhomogeneous Medium by Cross Correlation. Physical Review Letters, vol. 93 (25), pp. 254301.
Wapenaar, K., Draganov, D., Thorbecke, J. & Fokkema, J. (2002). Theory of Acoustic Daylight Imaging Revisited. Proceedings of 2002 SEG Annual Meeting, January 1, 2002. Society of Exploration Geophysicists.
Xu, Z., Juhlin, C., Gudmundsson, O., Zhang, F., Yang, C., Kashubin, A. & Lüth, S. (2012).
Reconstruction of subsurface structure from ambient seismic noise: an example from Ketzin, Germany:
Subsurface structure from seismic noise. Geophysical Journal International, vol. 189 (2), pp. 1085–
1102.
Yu, J. & Schuster, G.T. (2006). Crosscorrelogram migration of inverse vertical seismic profile data.
Geophysics, vol. 71 (1), pp. S1–S11.
31
Thorbergsson, G., Magnússon, I.Th., Gunnarsson, Á., Johnsen, G.V., Björnsson, A. (1984). Geodetic and gravity surveys in the Hengill area 1982 and 1983. Orkustofnun report in Icelandic, OS- 84003/VOD-03 B, pp. 58.
Examensarbete vid Institutionen för geovetenskaper ISSN 1650-6553
