# Seismicity related to the hydraulic stimulation of GRT1, Rittershoffen, France

O. Lengline, M. Boubacar, J. Schmittbuhl

► To cite this version:

O. Lengline, M. Boubacar, J. Schmittbuhl. Seismicity related to the hydraulic stimulation of GRT1, Rittershoffen, France. Geophysical Journal International, 2017, 208 (3), pp.1704-1715. 10.1093/gji/ggw490. hal-01924555

HAL Id: hal-01924555

https://hal.science/hal-01924555

Submitted on 19 Oct 2021

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_2.png" /></figure>

# Seismicity related to the hydraulic stimulation of GRT1, Rittershoffen, France

O. Lengliné, M. Boubacar* and J. Schmittbuhl

IPGS/EOST, Université de Strasbourg/CNRS, 5 rue René Descartes, F-67084 Strasbourg, France. E-mail: lengline@unistra.fr

Accepted 2016 December 30. Received 2016 December 12; in original form 2016 July 29

## SUMMARY

The Rittershoffen deep geothermal reservoir, in Northeastern France, is well characterized and has been extensively studied by a multidisciplinary approach. A hydraulic stimulation for the development of the geothermal reservoir was performed in June 2013. This injection of fluid led to seismic activity which was closely monitored by a dedicated set of seismic stations. The seismic sequence started during the injection but showed an unusual long quiet period of 4 d after shut-in before the occurrence of a second swarm of events. Here we take the opportunity of this well-monitored activity to gain insight into the geomechanical factors favouring the development of induced earthquakes. We apply a template matching approach and a relative relocation procedure to obtain a precise estimate of the geometries of the activated structures. Our approach shows that the induced events during the injection took place on two parallel planar structures. It shows that details of the seismicity generally obtained from borehole seismic network are achievable from surface network when an appropriate analysis is performed. The development of this induced seismicity is in good agreement with the known stress field and failure criterion proposed for the reservoir. In particular, the orientation of the activated structure, associated focal mechanisms and the overpressure needed to initiate the seismic activity are all in line with the geomechanical model of the area. The swarm of delayed events, 4 d after shut-in, can be well explained by considering an aseismic slip on the imaged fault and the related static stress transfer. We therefore suggest that the ability to monitor local slow aseismic movements at depth, in conjunction with precise tracking of the seismicity, is of primary importance to understand induced earthquake activity.

**Key words:** Fracture and flow; Induced seismicity; Rheology and friction of fault zones.

## 1 INTRODUCTION

Understanding how the injection of fluids in the crust is linked to the occurrence of earthquakes is of primary importance for the assessment of induced seismicity. Tracking and understanding this seismicity is important as it poses a risk to the nearby population and the infrastructures (Giardini 2009) and also provides information about the reservoir mechanics and the fluid pathway (Zang et al. 2014). However, modelling the direct impact of fluid injection on induced events is not straightforward as it involves the coupling of various processes. Indeed the nucleation of earthquakes is affected both by the existing fractures and faults within the reservoir, the in situ stress field, the strength of the rock reservoir, the mechanical and chemical impact of the fluid flow, and the temperature variations, factors that all interact with one another (e.g. Blanpied et al. 1995). Most of these parameters are usually not known pre-

cisely nor monitored continuously at depth during the injection such that it is difficult to observe their influence on the development of the earthquake activity. The clear identification of the role of each of these parameters on the seismicity is also restricted by the limited resolution that can be achieved on earthquake detection and locations. Indeed, the structures activated by fluid injection are usually of very small scale compared to the uncertainties associated with the earthquake locations. It thus requires a dedicated seismological network and processing in order to record and locate precisely the associated seismic activity and highlight the geometry of the active structures.

Despite these difficulties some mechanisms have been tested and proved as having an influence on the development of the seismicity. The mechanical influence of the pore-pressure, P, on the triggering of earthquakes is generally evaluated from a Mohr–Coulomb failure criterion, where slip on an interface occurs when

$$ \tau - \mu_i(\sigma_n - P) - C_0 > 0 \quad (1) $$

<footnote>* Now at: INRS-Quebec, Institut National de la Recherche Scientifique, 490 rue de la Couronne, Québec, G1K9A9, Canada.</footnote>
<figure><img src="image_3.png" /><figcaption>**Figure 1.** Map of the area around the Rittershoffen geothermal site (from OpenStreetMap). The geothermal wells of the Soultz site are visible in red on the north west of the map. The principal stress directions deduced from Cornet *et al.* (2007) from imaging in these well are shown with thick black arrows. The names of the seismic stations are indicated above each instrument (triangles). The seismicity (orange and blue dots) defines a 1.5 km long cloud mostly located west of the injection well. We observe that events that occurred during the second swarm (blue dots) are mostly present in the northern end of this cloud.</figcaption></figure>

where $\sigma_n$ is the normal stress on the fault, $\tau$ the shear stress, $\mu_i$ the internal friction coefficient of the interface and $C_0$ is the cohesion of the material. Note that for pre-existing weak interfaces like joints or faults, the friction coefficient and the cohesion are typically smaller than for intact materials (Scholz 2002). The experiment performed in Rangely oil field, Colorado, USA, showed that the onset of the earthquake activity was linked to a threshold injection overpressure $P$ that makes the left hand term in eq. (1) positive (Raleigh *et al.* 1976). Cornet *et al.* (2007) also showed a nice example from Soultz-sous-Forêts where seismicity onset and location are well predicted by the simple description of eq. (1). It shows that knowing the state of stress and the strength of the activated faults allows good prediction of the onset and the fate of the seismicity. However, other examples of induced seismicity exist where the earthquake activity presents a more complex dynamics that the one described by the eq. (1) alone. Indeed, in some cases, the diffusion of the pore pressure, elastic stress redistribution caused by the earthquakes or the seismic movement within the reservoir has shown to influence the seismicity (Shapiro *et al.* 2002; Bourouis & Bernard 2007; Schoenball *et al.* 2012). It therefore appears that the induced earthquakes could be linked to various mechanisms.

Here we analyse the earthquake activity that develops during the hydraulic injection in the well GRT1 in Rittershoffen, France. This well is part of a deep geothermal project located at 6 km east of Soultz-sous-Forêts in Northern Alsace and developed by the ECOGI joint-venture. For this purpose, at the end of 2012, a first well (GRT1) was drilled to 2580 m depth through Triassic sediments and into the crystalline basement. In order to enhance the reservoir permeability, a hydraulic stimulation was performed in the GRT1 well in June 2013. The hydraulic stimulation in GRT1 lasted 2 d (2013 June 27 and 28) and was recorded by a dedicated seismic network (Maurer *et al.* 2015; Baujard *et al.* 2017). The seismic

activity related to the GRT1 hydraulic stimulation was processed in real-time and gave rise to a first seismicity catalogue composed of a total of 212 events, from 2013 June 27 to July 4 (Maurer *et al.* 2015). The catalogue reveals that the seismicity stopped shortly after injection but started again after four completely quiet days on July 2, in the form of an intense seismic swarm that lasted less than one day. The occurrence of this second swarm after a quiet period is relatively rare and poses the question of its mechanical origin. In order to understand how this second swarm developed several days after the injection was stopped, we apply a dedicated set of tools to recover and locate the earthquakes that occurred during this sequence as precisely as possible. We see that, given the current state of stress of the reservoir, the strength of the rock mass forming the reservoir and the orientation of the main activated structures, we can propose a logical geomechanical model for the development of the observed induced and triggered seismicity.

## 2 DATA

The seismic network around the Rittershoffen geothermal field at the time of the injection in June 2013 was composed of 18 surface stations in total (Maurer *et al.* 2015; Fig. 1). The recording sampling frequency of these stations is distributed as follows: 9 at 150 Hz, 4 at 100 Hz and 5 at 300 Hz. All stations were continuously recording the seismic signal. Most of the sites have three components sensors and, depending on the site, these are broad-band or short-period sensors. The material differences are the result of the various equipments installed in the area over time: the Soultz permanent network, the Rittershoffen permanent network and the Rittershoffen temporary network. These networks provide a dense station coverage with 12 stations located less than 5 km away from the injection point. The
<table><caption>**Table 1.** 1-D velocity model of the Rittershoffen area used for the initial location of the events.</caption><thead><tr><th>Depth (m)</th><th>V<sub>P</sub> (km s<sup>-1</sup>)</th></tr></thead><tbody><tr><td>0</td><td>1.8</td></tr><tr><td>500</td><td>2.4</td></tr><tr><td>600</td><td>2.5</td></tr><tr><td>850</td><td>2.9</td></tr><tr><td>1050</td><td>3.3</td></tr><tr><td>1450</td><td>3.5</td></tr><tr><td>1650</td><td>4.4</td></tr><tr><td>2250</td><td>4.8</td></tr><tr><td>2350</td><td>4.9</td></tr><tr><td>2550</td><td>5.6</td></tr><tr><td>15 000</td><td>5.9</td></tr></tbody></table>

region at the south of the injection point, however, is not as well uniformly covered by seismic sensors because of permit issues that limited their installation.

We first process the seismic signals by performing a Short Term Average / Long Term Average (STA/LTA) approach (Earle & Shearer 1994) to detect possible seismic events over the whole time period. All signals are first filtered between 10 and 40 Hz, the STA is computed over 1 s and the LTA over 5 s of signal. We set a detection criterion when at least three stations reach the STA/LTA threshold fixed at 2. This leads to 946 detections. We then review manually all detections to discard false detections and only keep clear seismic events, resulting in a total of 686 earthquakes. We finally pick visible *P*- and *S*-wave arrival times for these events and obtain on average eight *P*-wave picks per event. The location of these events, based on the picking arrival times, is then obtained by the software HYPOINVERSE (Klein 2002) using a 1-D gradient velocity model of the area (Table 1). This velocity model is derived from a sonic log that was performed in GRT1 well (Maurer et al. 2015). The Vp/Vs ratio varies between 1.76 and 1.77. A final set of 674 seismic events could be located with this approach. Earthquakes that are discarded at the location stage are those with too few picks to be located.

## 3 TEMPLATE MATCHING DETECTION

The processing described above, STA/LTA detection combined with the location based on absolute arrival time picks, shows that almost all of the seismicity occurs less than 1 km from the injection point (Fig. 1). We also observe that events of the second swarm show an offset compared to the events that occurred during the injection. However, the lack of resolution in the location of the detected earthquakes prevents any further detailed analysis at this stage. The low resolution is partly related to the difficulty of picking the arrival times of the noisy seismic signals and the supposedly small scale of the activated structures.

The seismicity we are investigating took place in a densely populated area (cf. Fig. 1) with numerous noise sources masking the earthquake signals (Lehujeur et al. 2015). Not only can this cause imprecise picking arrival times but this could also lead to a reduced number of events in the final catalogue. Indeed, a low signal-to-noise ratio at a sufficient number of stations, or events occurring close in time, will not result in a positive detection by the STA/LTA algorithm. It is also possible that some events that were detected by the STA/LTA procedure are actually not located confidently because of too few or inconsistent picks resulting in the removal of

these events during latter processing stages. All these factors result in numerous missed events. In order to improve the seismic catalogue and to recover events that may have been missed during the previous procedure, we employ a template matching approach. The template matching approach has proved to be successful in recovering previously missed events in various contexts (e.g. Peng & Zhao 2009; Helmstetter et al. 2015; Lengliné et al. 2016). It has also been successfully employed to track induced seismicity (Skoumal et al. 2014; Huang & Beroza 2015).

The template matching approach uses a known seismic signal in order to detect newer events similar to the tested one. This technique has the ability to detect events even when the signal to noise ratio at an individual station is lower than 1. It however requires the knowledge of some template events to compare the continuous seismic signal with. The choice of these events can affect the final set of newly detected events as we will detect mainly similar events to the tested ones. It implies that the newly detected events are not located too far from the template event, otherwise waveforms will not be similar and stacking the correlation signals at different stations after correcting from travel times between stations will not result in a coherent stack. One could consider all the detected events as possible template events. This approach will be quite unnecessarily time consuming as many detected events are similar such that it will akin to processing the same data several times thus obtaining the same results. Furthermore, because some identified events have a low signal-to-noise ratio, it is possible that at some stations the correlation detector will mainly correlate with the noise source and not the earthquake signal. In this case it will produce multiple false detections. To overcome these difficulties, we consider here a different approach. We first group all the 674 located events into clusters of similar waveforms. We compute the correlation matrix between all events around the *P*-wave on 128 samples long window (i.e. 1.28 s) filtered between 8 and 25 Hz. The frequency range 8–25 Hz has been chosen based on the analysis of the spectrogram of the continuous data. This frequency range is best for filtering out most of the noise present in the data while keeping the earthquake signature and thus maximizing the signal to noise ratio. This range typically corresponds to the range where the corner frequency of the events is expected to fall within. All signals are first re-sampled at 100 Hz for all stations. If the mean normalized correlation coefficient at least at three stations is higher than 0.8, we associate the two events into the same cluster. We tested that reducing the correlation coefficient threshold for linking events in a same group do not lead to an improvement of the final number of detected events. We finally keep 13 clusters that all have at least two events. For each group we stack the waveforms of all events in the group in order to create a synthetic waveform representing the average waveform of the events of the group. Events are first aligned on the *P*-wave arrival before stacking is performed at all possible stations. We thus create synthetic template signals for the 13 groups of events. This approach is somewhat similar to the subspace detector, but here we consider only the first vector of the singular value decomposition that represent the stack of the events (Harris & Paik 2006; Maceira et al. 2010; Barrett & Beroza 2014).

The template matching approach is performed by correlating the template signals with the continuous signal at four stations (E15, BETS, RITT and KUHL) representing 10 channels in total (for the first three stations, we consider the signals on the three components while at station KUHL only the vertical component is available) at each time step (0.01 s). Template window signals are 2.56 s long and start 0.5 s before the *P*-wave pick. Correlations are performed after filtering the template signal and the windowed continuous signal
<figure><img src="image_4.png" /><figcaption>**Figure 2.** Cumulative number of earthquakes (dashed black line) and hourly rate of events (dashed red line) detected by the STA/LTA procedure and validated during the location. The two injection periods are in light grey and dark grey. We observe that no earthquake occurred during the second injection and that a burst of earthquake took place on the 2nd of July, after 4 d of quiescence. The plain lines refer to the number of events obtained after detection by the template matching approach.</figcaption></figure>

between 8–25 Hz. We use the overlap-add algorithm and the FFTW library for a fast implementation of the correlation computation (Harris & Paik 2006; Frigo & Johnson 1998).

The selected stations and channels are the four closest stations from the injection point with the best signal-to-noise ratio and are distributed with a good azimuthal coverage relative to the injection point. After the correlation coefficient is obtained at all individual channels for a given template event, these vectors are time-shifted according to the travel time differences between the stations. This step ensures that the resulting stack of the correlation coefficient vectors from these different stations is only coherent if the detected signal originates from a source proximal to the template event. In order not to be too restrictive using this procedure and possibly excluding new events not exactly co-located with the template event, we apply a maximum filtering over a duration of 0.1 s before stacking the correlation signals at the different stations. The maximum filtering simply transforms the correlation signal by replacing every value by the maximum in its 0.1 s range (Van Herk 1992; Gil & Werman 1993). It thus permits the detection of somewhat similar waveforms to the template but not located exactly in the same area. Stacking is performed over the 10 channels for each time step. We consider all parts of this averaged correlation signal where the correlation coefficient is higher than 0.4 as a possible detection. We visually check that extracted waveforms do correspond to earthquake signals. Reducing this correlation coefficient threshold produce some new events that could not be ambiguously distinguished from noise. If multiple templates detect the same event we only consider one detection associated with the template for which the correlation coefficient is the highest over the considered time-window. We run this template matching approach for the seven days covering the injection period and the second swarm (from the 2013 June 27 until 2013 July 3). This leaves us with 1395 events from the original set of 13 earthquake templates. We notably recover 93 per cent of the original catalogue’s events. For all of these new detections we extract the seismic waveform on all stations around the P-wave arrival for 2.56 s. Our template detection approach con-

firms that only one earthquake occurred in the injection area during the nearly 4 d following the shut-in. Seismicity resumed on 2013 July 2 with a strong burst with 225 events in this improved version of the catalogue (see Fig. 2).

## 4 RELATIVE RELOCATIONS

Because of the high similarity between the detected events (Fig. 3), it is possible to obtain a precise relative relocation of these events based on double-difference algorithm. We compute travel time delays between all the 1395 detected events. Delays are computed on 128 sample long windows filtered between 8 and 25 Hz and centered on the *P*-wave arrivals for vertical components, and on *S*-wave arrivals for horizontal components. We keep all traveltime delays associated with a correlation coefficient higher than 0.6 for the relocation process. In total we observe a much higher signal-to-noise ratio and similarity of the signal on the horizontal components. This is attested by the number of delays finally retained: 7 255 735 from *S* waves (both *E* and *N* components) and 2 306 002 from *P*-waves. The relocation is performed from the Software HYPODD (Waldhauser & Ellsworth 2000) where we set the initial position of all the events to be located on the open hole borehole trajectory at 2300 m depth. We are able to relocate confidently 1393 out of the 1395 detected events.

The relocated data set shows with finer details that the seismicity took place on two geometrically independent faults in the reservoir (Fig. 4 and Supporting Information Fig. S1). We acknowledge that the overall geometries of the apparent structures we are resolving with the relocation can be slightly affected by the velocity model as this will modify the computed ray geometries (Kinnaert et al. 2016). There is a clear offset between these two quasi-parallel structures. The earthquakes taking place on the structure in the south are well fitted by a plane of azimuth N025°E and dipping 74° towards west (Fig. 5). The second structure cannot be as well fitted by a single plane. We find that the azimuth varies between N005°E and N025°E.
<figure><img src="image_5.png" /><figcaption>**Figure 3.** Example of the waveforms of all 1395 detected events at two different stations. All waveforms are first filtered between 8 and 25 Hz and normalized by their maximum amplitude. Top: waveforms aligned on the S-wave arrival (around 1.5 s) recorded at station E15 on the north component. Bottom: waveforms aligned on the P-wave arrival (around 0.5 s) at the station KUHL on the vertical component. The vertical black line in both figures refers to the last event that occurred during the injection time period. We can notice the distinct waveforms for the events that took place during the second burst.</figcaption></figure>

The sparsity of the events on this second structure prevents us to obtain a reliable estimation of its dip (Fig. 4). Events that occurred during the injection (induced events) are mostly located on the fault to the south and the events of the second swarm on the fault to the north. The open-hole section of the borehole extends between 1920 and 2560 m but we observe that most earthquakes occur between 2200 and 2400 m depth. There exists a slight shift of nearly 20 m between the location of the first earthquakes and the borehole. Such an absolute shift is however not constrained from the relative relocation procedure and it is likely that the imaged structure actually intersects the borehole. Based on borehole imaging of well GRT1 we constrain the first earthquake that occurred during the injection to coincide with the well trajectory at 2368 m depth where

a major structure was observed (see Section 7.1). The earthquake activity during the injection progressively migrates upward and to the south along the identified fault plane (see Fig. 6). It thus mostly corresponds to an asymmetric migration relative to the injection point. At the end of the injection, the seismicity illuminates a planar 300 m long structure with a vertical extent of 200 m. The average along strike migration speed of the seismicity is of the order of 12 m hr⁻¹. We observe that some events at the end of the injection, before shut-in, are located in the northern structure mainly activated 4 d later (see Fig. 4 around latitude 48.8985°N). Earthquakes that occurred during the second swarm appear to occur on a plane with quite a similar orientation as those activated during the injection. We also observe that there exists an approximately 100 m shift between
<figure><img src="image_6.png" /><figcaption>**Figure 4.** Map of the 1393 detected and relocated earthquakes (dots). The orange dots refer to earthquakes that took place during the injection time period and the blue ones are those that took place during the second swarm, 4 d later. The star indicates the first earthquake that occurred during the injection. The purple line shows the trajectory of the injection well GRT1 with the blue termination indicating the open-hole section of the well. The grey rectangle marks the best fitting plane to the seismicity during the injection. The polarities of the two events of the injections with the clearest P-wave arrivals are reported. We also show the two best fitting mechanisms.</figcaption></figure>

this plane and the one inferred from the location of earthquakes
during injection (see Figs 4 and 6). It is possible that the two faults
imaged by the seismicity are part of *en-echelon* system as commonly
observed for faults in the Rhine graben (Schumacher 2002). We also
observe that earthquakes during both crises are all located around
the same depth interval and that almost very few events are located
below 2400 m. The second burst of events that occurred 4 d after
shut-in started at 48.900°N, at the middle of the visible structure
where we observe a change of direction of the event locations. The
earthquakes during this second burst appear to migrate quickly from
the initiation point. This migration occurs in both directions from

<figure><img src="image_7.png" /><figcaption>**Figure 5.** Cross-section of the detected and relocated earthquakes that occurred during the injection (orange dots). The cross-section line is displayed in Fig. 4. The star marks the location of the first earthquake of the sequence which is positioned at the location where the main permeable structure intersects the well (Vidal *et al.* 2016). Positive distance indicates east. The purple line indicates the trajectory of the injection well GRT1. In all the displayed depth range of this figure, the well is an open hole. The grey line is the orientation of the best fitting plane of the seismic cloud.</figcaption></figure>

the first earthquake of this second burst. The depth range of the
events in this second burst is similar to the events of the first swarm.

5 FOCAL MECHANISM

Obtaining a focal mechanism based on first motion polarities for
the detected events is not an easy task as most of these events show
a low signal-to-noise ratio. Furthermore, the station coverage is not
optimal for deducing the focal mechanisms of the recorded earth-
quakes owing to the absence of stations in the south. However, it
is still possible to estimate the mechanisms of a few large events
with a higher signal amplitude. We focus on the two events that
occurred during the injection with the largest signal to noise ratio at
most stations and with the most clear first motions. At each station
we visually determine the polarity of the signals. We compute the
ray parameters for the two events based on their relocation. We ob-
serve that the inferred orientation of the fault plane from the spatial
distribution of the seismicity is in agreement with the identified
<figure><img src="image_8.png" /><figcaption>**Figure 6.** Projected positions of all relocated events (circles) along a line with an azimuth of N25° corresponding to the azimuth of the best fitting plan of the injection-induced events. The star indicates the first earthquake that took place during the injection. The fading of the orange and blue colours represents the time since the first event of each swarm. Note the nonlinear scale for the time of the events during the second swarm related to the fact that numerous earthquakes during this swarm occurred in less than 1 hr.</figcaption></figure>

polarities (Fig. 4). It suggests that these events took place on the aforementioned fault plane and with a pure strike-slip mechanism, although the polarities can still be as well fitted by adding a slight normal component to the slip vector.

## 6 RELATIVE MOMENT

We estimate the relative moments between all earthquakes of the sequence to gain insight into any possible moment variations during the investigated period. For all relocated events we employ a Singular Value Decomposition to compute earthquake’s relative moment (Rubinstein & Ellsworth 2010). The coefficients of the first basis vector for each event represent the amplitude of the coherent part between each tested waveform and can then be taken as representative of a relative moment. We employ this approach to compute the relative moments at four stations with a good azimuthal coverage around the injection well (E17, E15, BETS and RITT). The waveforms are not filtered and the singular value decomposition is performed on both horizontal components at each station. The final estimate of the relative moment is taken for each event as the median value obtained for all traces. All relative moments are normalized to the event with the lowest moment, here fixed arbitrarily to 1. We observe that the moment distribution for the events that occurred during the injection (first swarm) and during the second swarm are clearly distinct (Fig. 7). Earthquakes that took place during the injection have lower moments. Their distribution can be described above the completeness moment by a high B value where B is defined as

$$N(M_0) = a M_0^{-B} \quad (2)$$

and $N(M_0)$ is the number of earthquakes with a relative moment higher than $M_0$ and a is a pre-factor. We observe that the distribution of moments for the events that occurred during the second swarm is well compatible with a typical $B = 2/3$ value as observed for earthquakes worldwide (Scholz 2002). It suggests that the events of the second swarm have different characteristics compared to the one induced during the injection (first swarm). The cumulative seismic

<figure><img src="image_9.png" /><figcaption>**Figure 7.** Distribution of relative moments for the earthquakes during the injection (blue circles) and for the ones of the second swarm (orange circles). All moments are normalized by the moment of the smallest event. A power law decay with an exponent $B = -2/3$ is drawn as a black dashed line for reference.</figcaption></figure>

moment released by the event during the second swarm is 1.5 times larger than the seismic moment released by the events during the injection.

## 7 DISCUSSION

### 7.1 Borehole fault zone and seismicity fault plane

Well logging was performed in the well GRT1 both before and after the hydraulic injection. It reveals that a major structure associated
with the largest temperature anomaly and shows that a significant permeability enhancement is located at 2368 m depth. Acoustic imaging allows the determination of the orientation of this fault zone which strikes N175°E, and dips 65°W (Vidal et al. 2016). Geological logs show that this major fault zone intersects the borehole in the granitic basement close to the bottom of the upper altered zone of the granite. This major permeable structure (its estimated width is 24 cm) is then supposedly the main fluid pathway for the injected water. We thus hypothesize that the main structure oriented N25°E and dipping 70°W as evidenced by the earthquake locations could correspond to this main permeable structure imaged in the well bore. This is also in agreement with the first motion polarities. The slight variation of the azimuth of the inferred fault plane between seismicity and acoustic imaging could result from the uncertainty associated with the strike determination in the acoustic image but also from the different scale of the two types of measurement. Indeed, the orientation given by fitting the seismic cloud provides a large-scale global orientation while the borehole image only gives a local determination of the geometry of the fault structure. It confirms that the identified earthquakes do not occur as an isotropic 3-D cloud but on a pre-identified planar fault.

## 7.2 Regional stress field and induced seismicity onset

The state of stress in the investigated area has been extensively discussed in a number of studies mainly from measurements performed in the wells of the Soultz-sous-Forêts geothermal site (e.g. Klee & Rummel 1993). As this geothermal zone is located less than 7 km away from the current site and as no major active geological structure is located between the two sites, we will assume that they both are under the influence of the same regional stress field, that is, as proposed in Cornet et al. (2007)

$$S_H \approx S_v \quad (3)$$

$$S_h = 0.54S_v \quad (4)$$

$$S_v = 33.8 + 0.0255(z - 1377) \quad (5)$$

$$P_p = 0.9 + 0.0098z, \quad (6)$$

where z is the depth (in metres) and stress are given in MPa. $S_H$ is the maximum horizontal stress, $S_h$ the minimum horizontal stress and $S_v$ is the vertical stress. The direction of $S_H$ has been found to vary with depth in the area. As the analysed earthquakes occurred over a limited depth range around 2300 m depth, we can consider that the orientation of $S_H$ is constant with an azimuth of N170°E ±10° (Cornet et al. 2007). The analysis of the orientation of the breakouts in the well GRT1 confirms that the direction of the maximum horizontal stress is nearly N-S and gives a direction of $S_H$ N180°E (Vidal et al. 2016).

The stress tensor of eqs (3)–(6) implies that it would require an overpressure in the well of 7.7 MPa for the fluid pressure to be higher than $S_h$ and thus create hydraulic fractures. The overpressure in the well did not exceed 3.1 MPa during the injection which ruled out such a scenario. Furthermore, the seismicity during this stimulation is not aligned in the direction of $S_H$ (perpendicular to $S_h$) which also confirms that the recorded events do not result from tensile hydraulic fractures. As we show below, the increase of overpressure first promotes the sliding on a pre-existing interface and at much lower stress that the one needed to propagate a mode I fracture.

<figure><img src="image_10.png" /><figcaption>**Figure 8.** Mohr circle for the initial state of stress at 2360 m for the Rittershoffen geothermal reservoir (black circle) and for the perturbed state of stress for a 2.65 MPa overpressure in the well (blue circle). The two lines refer to the two proposed failure criteria: (i) for the intact fault slipping aseismically (green line), and (ii) for the seismic patches (orange line).</figcaption></figure>

The stress field of eqs (3)–(6) promotes the occurrence of strike-slip, normal faulting. For a typical friction coefficient $\mu = 0.85$, the optimally oriented plane for rupture under a strike-slip environment has an azimuth of N25°E. It therefore suggests that the imaged strike-slip fault corresponds to a nearly optimally oriented structure with a slip vector in agreement with the regional stress field. We can represent the Mohr diagram corresponding to the stress field at the depth of the fault plane where it intersects the borehole (2368 m; Fig. 8). If the image log of the permeable structure at 2368 m actually corresponds to the same structure where seismicity is observed, as hypothesized, it then suggests that this structure was existing prior to the injection. We first consider the fault to be at the friction equilibrium defined by $\tau = \mu\sigma_n$, where $\mu = 0.85$ is the typical friction coefficient at low normal stress (Byerlee 1978). We observe that even in the absence of perturbation caused by fluid injection, the fault is nearly critical. Only a very small amount of stress perturbation will promote slip on this interface. It suggests that the fault starts to slip at the very beginning of the injection, for a very small overpressure, as predicted by Byerlee’s law. In the absence of seismicity it is suggested here that most of the slip occurred aseismically during this period.

We notice that the earthquake activity does not start immediately at the time of the injection. We record the first earthquake at 15:25 UTC on 2013 June 27. This time corresponds to a change of the injection flow rate which increases at that time from 27 to 40 l s⁻¹. This flow rate increase is linked to an increase of the downhole overpressure (measured at 1920 m deep) which reached 2.5 MPa at that time (Baujard et al. 2017). The fate of the earthquake activity at the time of the injection occurs around 06:30 UTC on the 28 of June which corresponds to the time where the flow rate is decreasing and the measured downhole overpressure goes below 2.5 MPa. Cornet et al. (2007) proposed that the failure condition for the nucleation of earthquakes in Soultz is well modelled by a Coulomb criterion (eq. 1) considering an internal friction coefficient $\mu_i = 0.96$ and a small cohesion of $C_0 = 1$ MPa. Such a value of the internal friction coefficient close to 1.0 is typical for a wide variety of rocks (Carmichael 1982). The orientation of the optimally planes for failure are 23.0° from the direction of $S_H$ considering this internal friction coefficient. We compute that the fluid overpressure needs
<figure><img src="image_11.png" /><figcaption>**Figure 9.** Distance from the injection point as a function of time for all earthquakes that occur during the injection (orange circles). The injection point is taken on the borehole path at the depth of the first earthquake. The green curve shows the maximum pressure reached during each of the injection step. The purple curve indicates the radius on the disk based on eq. (9) while the dark curve is the front of the diffusing pressure for a diffusivity of $5 \times 10^{-2} \text{ m}^2 \text{ s}^{-1}$.</figcaption></figure>

to increase to 2.65 MPa in order for the Mohr envelope to reach such a failure criterion (Fig. 8). The development of the induced seismicity during the injection very well obeys this failure criterion. When the fluid pressure reaches the 2.5–2.65 MPa level the patches of the fault representing seismic asperities (with a higher cohesion and an internal friction coefficient) are activated and the seismic activity starts. The same process applies at the end of the injection when seismicity stops as soon as the Coulomb failure criterion is not met and the fault stops to slip when the fluid overpressure goes back to 0.

## 7.3 Seismicity migration during the injection

We propose two main hypotheses to explain the migration of the seismic activity that occurs during the injection. In a first model, we propose that the water pressure resulting from the injection is propagating in an open fault. It leads to a stable crack propagation that is related to the volume of the injected water. If we suppose that this injected volume, $V$ is expanding over a disk of radius $r$ and thickness $h$, then

$$V(t) = \pi h r^2(t) \qquad (7)$$

where $h$ is the thickness of the fault where the fluid is propagating and is supposed constant. The injected volume can be computed from the different flow rates $Q(t)$ used in the various steps of the injection procedure,

$$V(t) = \int_{0}^{t} Q(t') dt'. \qquad (8)$$

Combining eqs (7) and (8) leads to the expression of the radius of expanding injected volume

$$r(t) = \sqrt{\frac{\int_{0}^{t} Q(t') dt'}{\pi h}} \qquad (9)$$

Applying eq. (9) to the actual injection history leads to results presented in Fig. 9. The best fit of the model to the seismicity data is obtained when setting the thickness of the fault (the boundary of

the injected volume) to 12 mm. We show that during the first part of the injection, when we observe the migration of the seismicity, the radius $r$ defining the limit of the propagating crack in our model is expanding almost linearly with time.

Our second hypothesis to explain the migration of the seismicity is to consider that the pressure pulse is diffusing in a permeable medium. In this case, and hypothesizing that the applied perturbation to the medium can be approximated as a step function, Shapiro et al. (1997) showed that the pressure pulse is propagating away from the injection point as $r(t) = \sqrt{4\pi Dt}$, where $D$ is the hydraulic diffusivity of the surrounding medium. Applying this model to the seismicity recorded during the injection, we found that the best hydraulic diffusivity that fits the data is $7 \times 10^{-2} \text{ m}^2 \text{ s}^{-1}$. This value is of the same order of magnitude to the value deduced for the nearby Soultz reservoir during the 1993 injection ($5 \times 10^{-2} \text{ m}^2 \text{ s}^{-1}$; Shapiro et al. 2002). We observe that the first hypothesis, despite its simplicity, explains the data at least as well as the second one. It suggests that the development of the seismicity on the fault plane can be well explained knowing the failure criterion of the fault, the actual state of stress and the injection history.

## 7.4 The triggering of the second swarm—Coulomb stress transfer from continuous slow aseismic slip

It is difficult to explain the occurrence of the earthquakes of the second swarm only from the diffusion of the injected fluid. Indeed, in this case we would rather expect a continuous seismic activity over the entire period extending from the shut-in to the time of this second swarm. Such a prolonged activity after the end of injection is very often encountered in geothermal exploitation (Zang et al. 2014). Furthermore, as the over-pressure in the well reaches almost 0 MPa at the end of the 28th June it is difficult to propagate such a pressure front with no over-pressure at the source. This scenario is confirmed by the analysis of the moment distribution which points to a very fast decay of the moment distribution of earthquakes during the injection and to a much lower decay compatible with normal tectonic earthquakes for the events of the second swarm. It is tempting to link this change of the $B$ value to the fluid overpressure as observed in various reservoirs (e.g. Bachmann et al. 2012). It then indicates that the events of the second swarm are certainly related to a low fluid pressure. These observations point out that another process is responsible of the observed activity. We test the hypothesis that a large scale aseismic slip on the fault plane imaged by the seismicity could have been responsible for the triggering of the activity during this second swarm. Aseismic slip in geothermal reservoirs have been reported from various observational evidences: shift of the borehole, velocity variations and repeating earthquakes (Cornet et al. 1997; Bourouis & Bernard 2007; Calò et al. 2011). The effect of fluid injection on natural interfaces has been found to promote aseismic sliding (Guglielmi et al. 2015; Wei et al. 2015; De Barros et al. 2016). Such a slow slip within the reservoir would be coherent with low amplitude events observed during the injection compared to the post-injection events (e.g. Lengliné et al. 2014)).

We suppose that the whole structure imaged by the seismicity during the injection defines a single plane that slips slowly during the injection. It thus defines a plane of surface $A = 200 \text{ m} \times 300 \text{ m}$, and we suppose that the movement of this aseismic slip is in agreement with the one reported by the focal mechanism of earthquake taking place on this plane (pure left lateral strike slip faulting). From the dimension of this plane, an empirical relation for earthquakes gives a typical displacement of the order of 1 cm (Wells
& Coppersmith 1994). We assume that such an empirical relation also applies to slow slip although the stress drop involved in such events might be lower than those of earthquakes (e.g. Brodsky & Mori 2007). We note that it is difficult to constrain the actual size of the hypothesized aseismic slip with no broadband sensors at depth and that the approximate dimension of the plane given by the extension of the seismicity gives only a lower bound estimation. We thus look at the impact of uniform 1 cm of pure left-lateral slip on the defined fault in the variation of Coulomb stress changes on the nearby fault hosting the activity of the second swarm. As this second fault has a variable strike, we compute the Coulomb stress changes for two possible orientations of the receiver fault that spans the possible azimuth range (N005°E–N025°E) of the fault as observed from the seismicity. Coulomb stress changes are computed using the Coulomb 3 software (Lin & Stein 2004; Toda et al. 2005).

In each case we consider that the receiving fault is vertical but we check that assuming a dip of 74°, similar to the slipping fault, has no significant impact on the calculated stress changes. All stress changes are computed at a depth of 2300 m, corresponding to the depth of the earthquakes of the second swarm. We observe that the second swarm is mostly located in a region of stress increase. This stress increase is the most important when computed on structures with an almost N–S orientation. Such an orientation is similar to the strike of the northern branch of the fault that ruptured during the second swarm.

Although modest, these stress changes on the second fault can be responsible for the increase in seismicity observed 4 d later. Similar amplitude of stress increases are also found in the south of the slipping plane, but no earthquakes were recorded in this area. It is possible that no fault exists at that location or that faults are badly oriented for rupture such that the stress changes do not lead to an increase in seismicity. Several explanations can be invoked to explain the delay of the seismicity of the second swarm. First, the aseismic plane has been progressively slipping since the injection started and continued to slip after shut-in such that it reaches a 1 cm cumulative slip only after several days and then triggers the earthquake activity. However this explanation is not very plausible because it will be difficult to promote slip on this fault plane after shut-in as the fault pressure goes quickly to zero and thus moves the interface away from failure. A second hypothesis is that all the aseismic slip actually took place during the injection. We note that the stress increase on the second fault is only modest (Fig. 10). If this second fault is initially in the same state of stress as the one related to the injection, then it shows that the stress change of ~2 bars is not large enough to promote instantaneous failure on this plane. However, we can consider that this stress step actually causes a clock advance of the nucleating asperities on the fault plane that is slowly loaded (Gomberg et al. 1998; Perfettini et al. 2003). This stress step makes the frictional unstable patches accelerate instantaneously, such that they get closer to failure but still need more time before they actually reach instability. Then as the first event on this second plane occurs we can suppose that most of this part of the fault is now close to failure such the extra stress change caused by the elastic stress redistribution of the event itself triggers this short lived burst of activity during the second swarm.

# 8 IMPLICATION FOR POST-INJECTION SEISMICITY

Post-injection seismicity occurs frequently in many cases of induced seismicity (e.g. Häring et al. 2008). In most of these instances, the

<figure><img src="image_12.png" /><figcaption>**Figure 10.** Coulomb stress changes computed for a 1 cm homogeneous left-lateral slip on the structure outlined in white and delimited by the extent of the induced seismicity cloud. Orange and blue circles refer to earthquakes that occur during the injection and during the second swarm respectively. The two figures refer to two possible orientations of the receiver faults as indicated by the black line in the upper right corner.</figcaption></figure>

seismic activity following the shut-in is continuous and does not stop and resume at latter times as observed during the GRT1 stimulation. It is then difficult to discriminate if the post-injection seismicity results from the only effect of the pore-pressure diffusion or if aseismic movements within the reservoir could also participate to the triggering of earthquakes. The unique GRT1 induced seismicity sequence suggests that the post-injection events are well explained as the result of such slow slip events. It is therefore appealing to propose that aseismic movements are ubiquitous features of stimulated reservoir and can be responsible for the observed post-injection seismicity.

# 9 CONCLUSIONS

Using dedicated tools (template matching, relative relocation), we were able to obtain a precise image of the fault structures activated during the injection. The induced seismicity related to fluid injection appears to agree with a geomechanical model but it requires detailed knowledge of the investigated area: the existing fracture and fault network, the regional stress field amplitudes and orientations, the failure criterion for the existing fractures and faults, and the pore pressure at depth during the injection. It also requires a detailed analysis of the seismicity with precise location. Such information allows an estimation of when and where the seismicity would appear in the geothermal reservoir. As aseismic movements appear to have a significant role in the deformation of the reservoir and can redistribute the stress locally it is also important to monitor the displacement associated with such movements at depth. Unfortunately, geodetic measurements at the surface might not be sensitive
enough to resolve the small displacement on these structures. Some downhole broadband instruments or dedicated geodetic monitoring might be needed in the future to capture these slow movements.

## ACKNOWLEDGEMENTS

We deeply thank F. Cornet, M. Lehujeur, J. Vergne, M. Grunberg, A. Maggi, J. Vidal, V. Maurer, C. Baujard, A. Genter, N. Cuenot, E. Gaucher, X. Kinnaert, T. Kohl, M.-A. Rico and the Seismology group of Labex G-EAU-THERMIE-PROFONDE for fruitful discussions and M. Heap, L. Lau and M. Dolatshahi for grammatical assistance. We thank two anonymous reviewers for their comments. We also thank ECOGI for funding part of the monitoring network and for sharing some geoscientific data. This work was conducted in the framework of the Labex G-EAU-THERMIE-PROFONDE, which is co-funded by the French government under the program ‘Investissements d’Avenir’, a cross-border collaboration between industrial companies, ECOGI, and ES-Géothermie, academics (KIT, EOST) and the EEIG ‘Heat Mining’ of Soultz-sous-Forêts. We also thank the Geophysical Instrument Pool Potdam from the GFZ German Research Centre for Geosciences (Germany) for providing a set of temporary seismological units.

## REFERENCES

Bachmann, C.E., Wiemer, S., Goertz-Allmann, B. & Woessner, J., 2012. Influence of pore-pressure on the event-size distribution of induced earthquakes, *Geophys. Res. Lett.*, **39**, L09302, doi:10.1029/2012GL051480.

Barrett, S.A. & Beroza, G.C., 2014. An empirical approach to subspace detection, *Seismol. Res. Lett.*, **85**(3), 594–600.

Baujard, C., Genter, A., Dalmais, E., Maurer, V., Hehn, R., Rosillette, R., Vidal, J. & Schmittbuhl, J., 2017. Hydrothermal characterization of wells grt-1 and grt-2 in Rittershoffen, France: implications on the understanding of natural flow systems in the rhine graben, *Geothermics*, **65**, 255–268.

Blanpied, M.L., Lockner, D.A. & Byerlee, J.D., 1995. Frictional slip of granite at hydrothermal conditions, *J. geophys. Res.*, **100**(B7), 13 045–13 064.

Bourouis, S. & Bernard, P., 2007. Evidence for coupled seismic and aseismic fault slip during water injection in the geothermal site of Soultz (France), and implications for seismogenic transients, *Geophys. J. Int.*, **169**(2), 723–732.

Brodsky, E.E. & Mori, J., 2007. Creep events slip less than ordinary earthquakes, *Geophys. Res. Lett.*, **34**, L16309, doi:10.1029/2007GL030917.

Byerlee, J., 1978. Friction of rocks, *Pure appl. Geophys.*, **116**(4–5), 615–626.

Calò, M., Dorbath, C., Cornet, F. & Cuenot, N., 2011. Large-scale aseismic motion identified through 4-D P-wave tomography, *Geophys. J. Int.*, **186**(3), 1295–1314.

Carmichael, R.S., 1982. *Handbook of Physical Properties of Rocks*, CRC Press. Inc.

Cornet, F., Helm, J., Poitrenaud, H. & Etchecopar, A., 1997. Seismic and aseismic slips induced by large-scale fluid injections, in *Seismicity Associated with Mines, Reservoirs and Fluid Injections*, pp. 563–583, ed. Talebi, S., Springer.

Cornet, F., Bérard, T. & Bourouis, S., 2007. How close to failure is a granite rock mass at a 5km depth?, *Int. J. Rock Mech. Min. Sci.*, **44**(1), 47–66.

De Barros, L. et al., 2016. Fault structure, stress or pressure control of the seismicity in shale? Insights from a controlled experiment of fluid-induced fault reactivation, *J. geophys. Res.*, **121**, 4506–4522.

Earle, P.S. & Shearer, P.M., 1994. Characterization of global seismograms using an automatic-picking algorithm, *Bull. seism. Soc. Am.*, **84**(2), 366–376.

Frigo, M. & Johnson, S.G., 1998. FFTW: An adaptive software architecture for the FFT, in *Proceedings of the 1998 IEEE International Conference on Acoustics, Speech and Signal Processing*, vol. 3, pp. 1381–1384, IEEE.

Giardini, D., 2009. Geothermal quake risks must be faced, *Nature*, **462**(7275), 848–849.

Gil, J. & Werman, M., 1993. Computing 2-d min, median, and max filters, *IEEE Trans. Pattern Anal. Mach. Intell.*, **15**(5), 504–507.

Gomberg, J., Beeler, N., Blanpied, M. & Bodin, P., 1998. Earthquake triggering by transient and static deformations, *J. geophys. Res.*, **103**(B10), 24 411–24 426.

Guglielmi, Y., Cappa, F., Avouac, J.-P., Henry, P. & Elsworth, D., 2015. Seismicity triggered by fluid injection-induced aseismic slip, *Science*, **348**(6240), 1224–1226.

Häring, M.O., Schanz, U., Ladner, F. & Dyer, B.C., 2008. Characterisation of the Basel 1 enhanced geothermal system, *Geothermics*, **37**(5), 469–495.

Harris, D.B. & Paik, T., 2006. *Subspace Detectors: Efficient Implementation*, United States Department of Energy.

Helmstetter, A., Nicolas, B., Comon, P. & Gay, M., 2015. Basal icequakes recorded beneath an alpine glacier (glacier d'Argentière, Mont Blanc, France): evidence for stick-slip motion?, *J. geophys. Res.*, **120**(3), 379–401.

Huang, Y. & Beroza, G.C., 2015. Temporal variation in the magnitude-frequency distribution during the Guy-Greenbrier earthquake sequence, *Geophys. Res. Lett.*, **42**(16), 6639–6646.

Kinnaert, X., Gaucher, E., Achauer, U. & Kohl, T., 2016. Modelling earthquake location errors at a reservoir scale: a case study in the upper rhine graben, *Geophys. J. Int.*, **206**(2), 861–879.

Klee, G. & Rummel, F., 1993. Hydrofrac stress data for the european HDR research project test site soultz-sous-forêts, *Int. J. Rock Mech. Min. Sci. Geomech. Abstr.*, **30**(7), 973–976.

Klein, F.W., 2002. User's guide to HYPOINVERSE-2000, a Fortran program to solve for earthquake locations and magnitudes, Tech. Rep., US Geological Survey.

Lehujeur, M., Vergne, J., Schmittbuhl, J. & Maggi, A., 2015. Characterization of ambient seismic noise near a deep geothermal reservoir and implications for interferometric methods: a case study in northern Alsace, France, *Geotherm. Energy*, **3**(1), 1–17.

Lengliné, O., Lamourette, L., Vivin, L., Cuenot, N. & Schmittbuhl, J., 2014. Fluid-induced earthquakes with variable stress drop, *J. geophys. Res.*, **119**(12), 8900–8913.

Lengliné, O., Duputel, Z. & Ferrazzini, V., 2016. Uncovering the hidden signature of a magmatic recharge at Piton de la Fournaise volcano using small earthquakes, *Geophys. Res. Lett.*, **43**(9), 4255–4262.

Lin, J. & Stein, R.S., 2004. Stress triggering in thrust and subduction earthquakes and stress interaction between the southern San Andreas and nearby thrust and strike-slip faults, *J. geophys. Res.*, **109**, B02303, doi:10.1029/2003JB002607.

Maceira, M., Rowe, C., Beroza, G. & Anderson, D., 2010. Identification of low-frequency earthquakes in non-volcanic tremor using the subspace detector method, *Geophys. Res. Lett.*, **37**, L06303, doi:10.1029/2009GL041876.

Maurer, V., Cuenot, N., Gaucher, E., Grunberg, M., Vergne, J., Wodling, H., Lehujeur, M. & Schmittbuhl, J., 2015. Seismic monitoring of the Rittershoffen EGS project (Alsace, France), in *Proceedings of World Geothermal Congress*, Melbourne, Australia.

Peng, Z. & Zhao, P., 2009. Migration of early aftershocks following the 2004 Parkfield earthquake, *Nat. Geosci.*, **2**(12), 877–881.

Perfettini, H., Schmittbuhl, J. & Cochard, A., 2003. Shear and normal load perturbations on a two-dimensional continuous fault: 1. Static triggering, *J. geophys. Res.*, **108**(B9), 2408, doi:10.1029/2002JB001804.

Raleigh, C., Healy, J. & Bredehoeft, J., 1976. An experiment in earthquake control at Rangely, Colorado, *Science*, **191**, 1230–1237.

Rubinstein, J.L. & Ellsworth, W.L., 2010. Precise estimation of repeating earthquake moment: Example from Parkfield, California, *Bull. seism. Soc. Am.*, **100**(5A), 1952–1961.

Schoenball, M., Baujard, C., Kohl, T. & Dorbath, L., 2012. The role of triggering by static stress transfer during geothermal reservoir stimulation, *J. geophys. Res.*, **117**, B09307, doi:10.1029/2012JB009304.

Scholz, C.H., 2002. *The Mechanics of Earthquakes and Faulting*, Cambridge Univ. Press.
Schumacher, M.E., 2002. Upper Rhine Graben: role of preexisting structures during rift evolution, *Tectonics*, **21**, 6-1-6-17.

Shapiro, S.A., Huenges, E. & Borm, G., 1997. Estimating the crust permeability from fluid-injection-induced seismic emission at the KTB site, *Geophys. J. Int.*, **131**(2), F15–F18.

Shapiro, S.A., Rothert, E., Rath, V. & Rindschwentner, J., 2002. Characterization of fluid transport properties of reservoirs using induced microseismicity, *Geophysics*, **67**(1), 212–220.

Skoumal, R.J., Brudzinski, M.R., Currie, B.S. & Levy, J., 2014. Optimizing multi-station earthquake template matching through re-examination of the Youngstown, Ohio, sequence, *Earth planet. Sci. Lett.*, **405**, 274–280.

Toda, S., Stein, R.S., Richards-Dinger, K. & Bozkurt, S.B., 2005. Forecasting the evolution of seismicity in southern California: animations built on earthquake stress transfer, *J. geophys. Res.*, **110**, B05S16, doi:10.1029/2004JB003415.

Van Herk, M., 1992. A fast algorithm for local minimum and maximum filters on rectangular and octagonal kernels, *Pattern Recognit. Lett.*, **13**(7), 517–521.

Vidal, J., Genter, A. & Schmittbuhl, J., 2016. Pre-and post-stimulation characterization of geothermal well GRT-1, Rittershoffen, France: insights from acoustic image logs of hard fractured rock, *Geophys. J. Int.*, **206**(2), 845–860.

Waldhauser, F. & Ellsworth, W.L., 2000. A double-difference earthquake location algorithm: method and application to the northern Hayward fault, California, *Bull. seism. Soc. Am.*, **90**(6), 1353–1368.

Wei, S. et al., 2015. The 2012 Brawley swarm triggered by injection-induced aseismic slip, *Earth planet. Sci. Lett.*, **422**, 115–125.

Wells, D.L. & Coppersmith, K.J., 1994. New empirical relationships among magnitude, rupture length, rupture width, rupture area, and surface displacement, *Bull. seism. Soc. Am.*, **84**(4), 974–1002.

Zang, A., Oye, V., Jousset, P., Deichmann, N., Gritto, R., McGarr, A., Majer, E. & Bruhn, D., 2014. Analysis of induced seismicity in geothermal reservoirs–an overview, *Geothermics*, **52**, 6–21.

## SUPPORTING INFORMATION

Supplementary data are available at *GJIRAS* online.

**Figure S1.** Locations of the earthquakes in the original catalogue (blue circles) and the relocated events (red circles). The relocated events have been shifted to the centre of the seismic cloud formed by the original location of the earthquakes.

Please note: Oxford University Press is not responsible for the content or functionality of any supporting materials supplied by the authors. Any queries (other than missing material) should be directed to the corresponding author for the paper.