AVO study on a known geothermal reservoir located in the fractured carbonate formations of the pre-Cenozoic basement, Northwest Hungary

Amplitude Versus Offset (AVO) analysis has been successfully applied in the hydrocarbon industry for more than three decades. This effective tool of the joint seismic and well-log data processing enables to predict and analyze fluid saturated porous geological formations. AVO methodology is based on the anomalous behavior of the pre-stack reflected amplitudes observed from fluid bearing rocks. However, the potential of AVO methodology is still unexploited in geothermal exploration, although the lithology and rock physical properties are very similar. In this study, we summarize the theoretical backgrounds and calculate synthetic AVO responses of a known geothermal reservoir located in the fractured carbonates of the Mesozoic basement. We demonstrate that the AVO response of a deep geothermal reservoir can be quite different from the amplitude response observed from a hydrocarbon bearing clastic formation. AVO attributes of the investigated geothermal reservoir are presented, and the potential of its detection by seismic amplitude data are discussed.


Introduction
was the first who demonstrated that Amplitude Versus Offset (AVO) responses from both the top and base of a porous gas sand show anomalous (increasing) trends unlike the AVO responses observed from other lithological boundaries (showing normal decreasing trends). AVO analysis was introduced in the hydrocarbon exploration 1 3 in the late 80's and after that the number of productive wells increased significantly not only in the Gulf of Mexico (Allen and Peddy 1993) but all around the world. Later, several authors proved that utilizing the pre-stack amplitude variations, AVO analysis works also for fractured carbonate formations to indicate hydrocarbon bearing porous zones (Harvey 1993;Lynch et al. 1997;Li et al. 2003). Parallel with the above studies, the question came up whether AVO can be applied for other lithological investigations beyond the hydrocarbon exploration, for example to estimate rock properties in the depth rage of the middle crust (Pratt et al. 1993;Simon 1998). Hajnal 2000 andTakács et al. 2021 went even deeper and demonstrated the successful utilization of AVO processing and analysis focused on the vicinity of the crust-mantle boundary beneath south-east Hungary.
After Goodway (2001) discussed detailed relations between the Lamé parameters (λ and μ) and the lithology, and Russell (2006) pointed out the importance of elastic impedance inversion (EII) for lithological interpretation, the above question was finally answered. The elastic parameters (incompressibility and rigidity) are much more characteristic for lithological interpretation than the P-wave velocity or the S-wave velocity in itself. The new seismic tool (EII), based on the joint inversion of the pre-stack reflected P-and S-wave amplitude data, was added to the toolbox of successful hydrocarbon exploration.
In the last decade, only a few authors studied the benefits of AVO analysis for geothermal exploration, however the number of publications is growing in this topic. For example Aleardi and Mazzotti (2012) investigated a geothermal reservoir located in the fractured intrusive basement of the Larderello-Travale geothermal field, Italy. They concluded that AVO responses from the massive rocks of the basement can be different from the well-known amplitude responses (Ostrander 1984) taking in account for the hydrocarbon exploration. The authors declared that "we need to derive a new AVO attribute which may highlight fracture locations in this peculiar rock type". Russell (2020) also underlined the potential of AVO analysis and inversion in the geothermal exploration and concluded that adaptation of the methodology might be very promising for geothermal reservoir development. Recently Allo et al. (2021) published their study on the utilization of neural networks in the characterization of a carbonate geothermal reservoir located in the Paris Basin, France.

AVO theory
AVO procedure applies true amplitude processed pre-stack seismic data analyzing the reflected P-wave amplitudes as a function of the source-receiver offset (or incident angle). The offset-angle conversion of the CDP offset gathers is performed utilizing the velocity field obtained from seismic velocity analysis or the available sonic logs (Hampson-Russell Software Services Ltd. 2004). In this way, the rock physical parameters Pand S-wave velocity (V P and V S ) and density (ρ) above and below the seismic interfaces can be obtained utilizing a resistant and robust inversion algorithm. An advantage of the AVO method is that no direct S-wave recording is necessary to get S-wave velocity information from the seismic gathers. On the other hand, several other useful elastic parameters can be derived beside the above mentioned three parameters obtained from the AVO procedure. Those derived elastic parameters are the velocity ratio (V P /V S ), the Poissons's ratio (σ), the incompressibility (λ), and the rigidity (μ).

Rock physical basics
Any Amplitude Versus Offset (AVO) response depends on the P-wave (V P ) and S-wave (V S ) velocities, and the density (ρ) of the porous reservoir rock-which involves the matrix material, the pore space, and the fluid fills in the pores. A common way to look the ratio of V P and V S is given by the Poisson's ratio (σ): Poisson's ratio depends on the solid-rock material, the volume and shape of the pore space, and the fluid saturation (Hampson-Russell Software Services Ltd. 2004;Chopra and Castagna 2014). Since the rock material is dependent on the mineral composition, it influences the Poisson's ratio (Fig. 1). It is also important to note that in case of two-phase fluid saturated rocks, the Poisson's ratio is decreasing with increasing gas content (Fig. 2). The presence of only a very few amount of gas in the pores will significantly decrease the V P / V S ratio and consequently the Poisson's ratio (σ). The reason is that P-wave velocity drops very suddenly even in the presence of a small amount of gas content, but S-wave velocity does not change significantly because it does not travel in fluids, it travels only in the rock matrix. Gassmann (1951) and Biot (1956) established independently their developments for the wave propagation in fluid saturated rocks (Fig. 3). They derived expressions for the saturated bulk (K) and shear (μ) moduli and applied them in the regular equations of wave velocities (V P and V S ):  (Chopra and Castagna 2014). Note that mineral composition and porosity of the rocks have much more effect on the Vp/Vs and Poisson's ratio parameters than the P-T conditions 1 3

Biot-Gassmann equations
where is the density In the Biot-Gassmann equations, the shear modulus does not change with varying fluid saturation at a constant porosity (μ sat = μ dry ) and the bulk modulus is defined as follows: where sat distinguishes the saturated rock, dry the dry rock, m the rock matrix, fl the fluid, and ϕ is the porosity. A rearranged version of the above equation provides a more instinctive formula by Mavko et al. (1998):

3
The bulk modulus of the solid rock matrix (K m ) is usually taken from published data measured on drill core samples. The typical value for sandstones is 40 GPa and it is 60 GPa for limestones. The bulk modulus for the fluid fill can be calculated by the following relation obtained from the volume average equation: K fl is the bulk modulus of the fluid, K w is the bulk modulus of the water, K hc is the bulk modulus of the hydrocarbon, and S w is the water saturation. The Biot-Gassmann theory summarized above provides a good understanding of the elastic parameters in fluid saturated geological formations.

Zoeppritz equations and their approximations
Theoretical studies of the reflection and transmission coefficients were provided by Zoeppritz (1919) to determine the amplitude variations of the reflected and transmitted P-and S-waves with the incident angle at a single seismic interface for the entire range of incident angles (between 0° and 90°). In case of an incident P-wave, the equation is written as follows: Seismic waves generated at a single interface in case of incident P-wave (modified after Chopra and Castagna 2014) where R and T are the reflection and transmission coefficients, V P , V S , and ρ are the P-and S-wave velocities and the density, while θ and ϕ are the angles of incidence and transmission for all generated P-and S-waves (Fig. 4).
In practice, several simplified approaches of the above matrix equation are applied to calculate theoretical AVO responses and to perform AVO inversion of the reflected P-wave data. The most popular approaches are the Aki-Richards and Shuey's equations (Richards and Frasier 1976;Shuey 1985). Those approximations work properly in the range of incident angles that usually can be observed by a regular surface seismic survey.
Based on the Aki-Richards theory (in the range of 0°-40°), the reflection coefficient (R) can be calculated by a three-term function of the square of sinus and tangent incident angle (θ) at a single seismic interface: where the constants can be expressed as follows while θ is the incident angle, α is the P-wave velocity, β is the S-wave velocity, and ρ is the density of the lithological formations. According to Shuey's approach (in the range of 0°-30°), the third term C (curvature) can be neglected, and his approximation results in a simple linear function with the two very basic AVO attributes A intercept and B gradient (see also Fig. 5): It should be noted that the 'intercept' is related with the zero-offset reflection coefficient and the 'gradient' is very sensitive for the porosity and fluid content.

Examples for normal and anomalous AVO responses
The lithological boundaries, in most cases, generate decreasing (normal) AVO response but some of them produce increasing (anomalous) trend. In Fig. 6, a normal amplitude response is presented that usually can be observed from a boundary located between solid geological formations (with no porosity). In Fig. 7, an anomalous AVO response is demonstrated that usually can be detected from a boundary located between solid and fluid saturated porous rocks. At this point, it is important to define the terms of 'decreasing' and 'increasing' AVO trends. Decreasing and increasing amplitude trends in AVO analysis are meant regarding the absolute values of the reflected amplitudes. The amplitude variation with offset is obviously decreasing in Fig. 6 because positive amplitudes are decreasing with the increasing offset. However, the amplitude variation presented in Fig. 7 is increasing because the negative amplitudes (regarding their absolute values) are increasing with the increasing offset (all reflection coefficient is negative).
The P-wave reflection coefficient modelling was carried out with the utilization of the exact solution of the Zoeppritz's matrix equations (Young and Braile 1976). For the polarity, the North American convention was applied, in other words a positive reflection coefficient represents a peak in the seismic traces. It is also important to note that in Fig. 7, the zero-offset reflection coefficient is a negative value thus the AVO response is considered as an increasing trend because the absolute value of the amplitudes (or the reflected energy) is increasing.

AVO modelling based on well-log data available in the study area
Before starting a reliable AVO analysis, an important step is modelling the probable AVO responses of the study area (Masri and Takács 2020). We calculated theoretical AVO trends to study whether fluid saturated porous rocks of the area of investigation can cause any AVO anomaly. Our main task was to study the possible AVO responses from the top and base of the hot water bearing fractured zone located inside the Triassic carbonate basement. The lithological model and its rock physical parameters obtained from the well log data is presented in Table 1 where the productive (hot water bearing) fractured zone with a thickness of 54 m is highlighted in the green colored row (Table 2).
Resistivity, P-and S-wave velocity, and density logs from the most informative deviated hole after an 11-point median filtering are displayed in Fig. 8. We note that only the resistivity and P-wave logs were measured in the well, the S-wave velocity and density logs were calculated based on the Castagna et al. (1985) and the Geertsma and Smit (1961) equations. The correlation between the measured resistivity and P-wave logs is excellent and both the calculated logs show the expected decreased values inside the fractured zone between the green markers of FR and T1 (top and base of the fractured zone). Decreasing values below the green T1 marker are explained with the presence of Triassic calcareous marl which was drilled down to the end of hole (see Table 1).
The blocked rock physical parameters (V P , V S , σ, and ρ) of the model, which was utilized in the subsequent AVO modelling, are demonstrated in Fig. 9. It must be noted that a Poisson's ratio (σ) value of 0.13 was set up for the fractured zone (between the markers FR and T1); considering high formation pressure in the water saturated carbonates that   True amplitude recovery (TAR) See Table 3 for the desired aims and effects Trace editing (KILL, MUTE) Trace editing was carried out formerly when the data was processed for conventional imaging Surface consistent amplitude balancing See Table 3 for the desired aims and effects Air blast attenuation Attenuated velocity range: 200-700 m/s (in 50 m/s steps), filter to enhance air blast: 10-15-25-30 Hz Surface consistent deconvolution See Table 3 for the desired aims and effects Bandpass filtering 12-16-60-65 Hz Refraction statics Refraction Statics were calculated formerly when the data was processed for conventional imaging Velocity analysis Velocity analysis was carried out formerly when the data was processed for conventional imaging Residual statics Time varying trim statics using a sliding window with 400 ms length Kirchhoff PSTM migration See Table 3 for the desired aims and effects can make the rock matrix weak and decrease the S-wave velocity (V S ) in this way. The calculated zero offset reflectivity and the Zoeppritz synthetics are also presented in this figure. The synthetic traces show the variations of theoretical reflected amplitudes with the incident angle between 0° and 40°. An anomalous effect in the wave field caused by the hot water saturated fractured zone is clearly seen inside the dark green ellipse. A more detailed insight is provided in Fig. 10 showing all theoretical AVO responses of the model. We concluded that only the light green AVO curve from the base of the fluid saturated fractured zone (T1) shows anomalous behavior compared with the other curves. It has negative intercept, positive gradient, and very importantly an 'angle of crossover' less than 30°. This feature (low 'angle of crossover') distinguishes this AVO response from the other ones. Thus, a combination of the intercept, gradient, and the 'angle of crossover' as a new AVO attribute might be a helpful indicator (beside the conventional AVO attributes) for mapping the hot water bearing deep fractured zone in this study area. We intend to calculate this new attribute soon, see preliminary conclusion of Aleardi and Mazzotti (2012) mentioned in the introduction part of this study.

True amplitude pre-processing of the 3D seismic data
AVO analysis of real seismic data requires true amplitude pre-processing on the recorded seismic traces (Mazzotti and Mirri 1991). We must recover the relative amplitudes of the reflections which are distorted by the wave propagation and the data acquisition (for example energy attenuation and variations of the source and receiver coupling). On the Fig. 8 Logs from a productive deviated geothermal well (modified after Masri and Takács 2020). MD: measured depth [m]. Blue and red logs: measured resistivity and P-wave velocity; yellow logs: calculated S-wave velocity and density; green FR and T1 markers: top and base of the fractured zone inside the Triassic basement other hand, we must avoid any data processing steps that can change the trace-by-trace relative amplitudes (for example Automatic Gain Control or whole Trace Equalization). In that way, we can get seismic gathers with reliable relative amplitude variations (AVO responses) that are related with only the changes of lithology and fluid content. After a careful seismic data pre-processing, the next procedure is to create AVO attributes from the pre-processed CDP (Common Depth Point) gathers. AVO attributes are very helpful for lithological interpretation, the two most popular ones are the 'product' and the 'scaled Poisson's ratio change' derived from the basic AVO attributes 'intercept' and 'gradient' defined in Fig. 5 and discussed in Chapter 6.
After considering the above requirements, we applied the following processing sequence on the 3D seismic data by the utilization of ProMAX software (Landmark Graphics Corporation 2001;Yilmaz 2001) to get appropriate CDP gathers for the subsequent AVO analysis: The most beneficial steps of the above true amplitude pre-processing sequence from the viewpoint of AVO analysis are summarized in Table 3.

Seismic well tie
Well log data is measured in the depth domain and its conversion to the seismic two-waytime (TWT) will not give a good fit to the seismic time section for the first try. The reason is that depth to time conversion of any well log is performed usually by the acoustic logs and the velocity calculated from them differs from the seismic velocity. One of the reasons 1 3 for the difference is the different frequency used for acoustic logging (~ 10,000-30,000 Hz) and for the surface seismic observation (~ 10-100 Hz). On the other hand, the drilling weakens the rocks therefore acoustic logging can provide a bit lower velocity compared with the seismic velocity in some intervals. On the third hand, acoustic logs usually start below the topographic surface (focused for the target zone) thus they do not provide velocity information for the upper part of the geological model. The missing velocity information is estimated by an extrapolation of the measured acoustic data up to the surface.
Because of the above reasons, seismic well tie is necessary to get proper time converted well logs to insert them in the surface seismic section. This procedure means the shifting, stretching, and squeezing the logs based on the similarity between the synthetic seismic traces obtained from the acoustic and density logs (acoustic impedance) and the real seismic data (Fig. 11). After finalizing an accurate seismic well tie, the correlation coefficient between the synthetic and real seismic traces was 0.807, which can be considered as an appropriately high value. The correlation coefficient was determined by a built process of the Hampson-Russell software calculating cross-correlation function between the time series (synthetic and real traces). The closer the cross-correlation value is to 1, the more closely the compared traces are identical.  Table 1 and are displayed with different colors. Note the light green anomalous trend from the bottom of fracture zone (T1) with a low value of 'angle of crossover' (less than 30°) In Fig. 12, the acoustic impedance log before and after the seismic well tie is inserted in the closest seismic crossline. The effect of the procedure can be clearly seen in the figure. Before the seismic well tie, the log does not fit to the seismic section at all but after the correlation the log shows a very good fit to the seismic data. The correlation (seismic well tie) was performed using the interactive panel of the AVO software by shifting, stretching, and squeezing the logs looking for the best correlation between the synthetic (blue) and measured (black) traces shown in Fig. 11. 1 3

AVO attributes obtained from pre-stack seismic amplitudes
Mathematical background and calculation of the basic AVO attributes (A-intercept and B-gradient), as well as their rock physical meanings were discussed in Chapter 2.3 where we presented how to get those two attributes from the pre-stack CDP gathers (Fig. 5). The basic AVO attributes (A and B) are rarely used in themselves because other useful prestack parameters with more lithological meanings can be derived from them. Currently, we believe that 'Poisson's ratio change' is one of the most helpful fluid indicators to detect fluid saturated porous geological formations. With some assumptions, it can be easily derived as A + B (Chopra and Castagna 2014). One of the 'Poisson's ratio change' AVO attribute cross-section selected from the 3D data cube in the vicinity of a productive well is presented in Fig. 13. The most important conclusion based on the above figure is that Poisson's ratio is decreasing at the top of the fractured zone (marked with FR), and it is increasing at the base of the known reservoir (marked with T1). In other words, the green colors at around FR mean negative change in the Poisson's ratio while the blue colors at around T1 mean positive change. This feature of the Poisson's ratio distribution is characteristic for the fluid saturated geological formations. The decreased acoustic impedance (red log) between FR and T1 also indicate the fluid saturated fractured zone that can be mapped along the presented AVO attribute section and the entire 3D volume.
The lithological and structural interpretation has not been completed yet because we intend to continue our study with the simultaneous inversion, as well as with the LMR inversion of the pre-stack 3D seismic data. Those procedures will provide additional rock physical models to investigate several more parameters (P-and S-wave acoustic impedances, density, Vp/Vs ratio, and the Lamé parameters) to complete the final interpretation of the study area.
All productive holes in the area of investigation were drilled before the AVO analysis, and the deviated well shown in Figs. 12 and 13 was planned remarkably well utilizing only the 3D migrated volume and some of its conventional seismic attributes (reflection strength and variance). However, based on the 'Poisson's ratio change' AVO attribute, Fig. 12 Acoustic impedance log (a) before and (b) after seismic well tie inserted in the seismic stack (time domain). Blue trace: well path, blue markers: tops of geological formations, red log: acoustic impedance, blue horizon: top of the presumed reservoir zone (with negative reflection coefficient) further important conclusions can be drawn. On the one hand, the well hit the edge of the Triassic fractured zone and evidently that is the most disturbed zone providing the highest porosity and/or permeability along the displayed section. On the other hand, very likely, the well passed through another porous zone in the Miocene sediments just beneath the M2 marker at about 1670 ms (see the decreasing values both in the Poisson's ratio and the acoustic impedance at that location). They targeted the deeper Triassic formations for a higher temperature geothermal reservoir.
As the last step of our present study, we visualized the calculated 'Poisson's ratio change' values on the surface of the presumed dolomitic fractured zone (Fig. 14). Three wells drilled prior the AVO analysis successfully hit the geothermal reservoir, they are also marked in the map. It can be clearly seen that the bottom of all three wells is located in the light blue area of the map showing negative values of the 'Poisson's ratio change'. Note that the horizontal extension of the map is only 2 × 2 km because we wanted to cut the running time of the pre-stack migration process before testing the AVO analysis. However, this map proves that 'Poisson's ratio change' AVO attribute is a useful tool to design wells not for only the hydrocarbon exploration but also for geothermal purposes. Based on these results, we believe that AVO analysis adopted from the  hydrocarbon exploration can also be a successful tool in the geothermal exploration to assign prospective drilling targets.

Conclusions
The modelling results presented in this study are characteristic for only the present study area because of its specific geological-lithological fabric. However, the methodology of AVO analysis providing Poisson's ratio change vertical sections and maps can be extended to other study areas of the geothermal exploration.
From the viewpoint of further investigations, several useful conclusions can be drawn from our present investigations. By modelling theoretical AVO responses, we demonstrated that the amplitude response of the studied geothermal reservoir located in the fractured carbonates (light green curve in Fig. 10) is different from the usually negative increasing response of the top of a gas bearing clastic sandstone (Fig. 7). Our result confirms Aleardi and Mazzotti's (2012) previous statement, even though they studied fractured geothermal reservoirs in the intrusive basement.