# Prediction of Favorable Carbonate Reservoirs under Extremely Thick Salts via Poststack Facies-Controlled and Prestack Zoeppritz Equation Inversions in the Santos Basin of Brazil

Authors: Sheng Zhang; Suoliang Chang; Handong Huang; Yinping Dong; Youyi Shen; Yaneng Luo; Baoheng Zhu
Journal: Geofluids
 (2020)

---

## Abstract

Subsalt carbonate rocks in Brazil’s deepwater region have broad prospects for oil and gas exploration and development. Due to the low-frequency bandwidth of the seismic data and the poor signal quality for this kind of reservoir target, there is a demand for accurate seismic prediction methods. We employ the facies-controlled inversion using a low-pass filter matrix to ensure the accuracy of the low frequency and to improve the robustness of the inversion results. We integrated the concept of adaptive regularization constraint of the Zoeppritz equation into the generalized linear inversion theory framework, which overcomes the shortcomings of the approximate equation. Making full use of the large angle prestack seismic information, Zoeppritz equation inversion improves the accuracy of the inversion results. The application of this method in carbonate reservoirs under extremely thick salts in the Santos Basin of Brazil indicates the feasibility and practicality of the proposed integrated prediction method.


---

## Body

## 1. Introduction

Oil and gas exploration in Brazil has experienced a significant increase in activity since the first major presalt discoveries were announced in 2007 [1–3]. One of the largest offshore basins in Brazil, the Santos Basin, is a passive continental margin basin with the second largest reserves of oil and gas in South America [4]. The total proven oil reserves were 67.1 billion barrels, while the natural gas reserves were 1.63×108 cubic meters, meaning that Santos Basin has great potential for future production increases (according to the International Energy Agency). In this area, the water depth is greater than 2000 m, and the target layer is covered with a thick salt caprock, which seriously restricts the exploration and development process [5, 6]. The exploration and development of oil and gas in deepwater areas are risky and costly, and the prediction of deepwater carbonate reservoirs, especially if they are located under extremely thick salts, is a great challenge due to the complex depositional system and the limited level of exploration [7]. Therefore, it is of great significance to adopt appropriate seismic inversion and reservoir prediction methods for carbonate formations in this area.Seismic acoustic impedance inversion has the capability of linking seismic data to geology and well log data. Due to the high signal-to-noise ratio of poststack data and its low computational cost, the poststack acoustic impedance inversion is widely used in the industry for characterizing the reservoir properties [8–11]. However, the poststack inversion result is highly dependent on the accuracy of the initial low-frequency model [12, 13]. Large error in the low-frequency model can deviate the inverted model far away from the true solution. Omitting the low frequency of the seismic data is often compensated by interpolating the well log data [14, 15]. In deepwater carbonate reservoirs under extremely thick salt settings with sparse wells, this low-frequency impedance model may be inaccurate. Moreover, the poststack inversion cannot accurately identify the fluid properties when the acoustic impedance of the deep reservoir overlaps.Compared to a poststack inversion, the advantage of a prestack inversion lies in the comprehensive consideration of the amplitude information of multiple offsets and the simultaneous acquisition of multiple inversion attributes [16–19]. Since Connolly [20] proposed elastic impedance, the amplitude variation with the offset method and elastic impedance inversion have been developed rapidly [21–23]. The conventional prestack seismic inversion mostly adopts the approximation equation of the exact Zoeppritz equation under the assumption of weak anisotropy and the small incident angles, with the increase in the incident angle, the error between them increases, but the large angle information is very important for the accurate inversion of the density parameters [24–26].To accurately predict the properties of carbonate reservoirs under extremely thick salts, we employ a low-pass filtering matrix, and a seismic facies constraint is applied to the seismic inversion to obtain a robust low-frequency trend. Compared with other regular constraints and petrophysical models, the seismic facies constraints are relatively insensitive to the weight coefficients, which helps to improve the inversion efficiency on the basis of guaranteeing validity and accuracy. At the same time, to reliably identify the properties of the fluid, in this paper, the exact Zoeppritz equation is introduced into the generalized linear inversion theory framework, and the iterative equation of the prestack inversion is derived. The adaptive regular constraint is added to improve the inversion accuracy of the elastic parameters, and the application of real data results proves the high accuracy and strong reliability of the proposed method, which can guide the exploration and development of similar oil and gas reservoirs.


## 2. Methodology

### 2.1. Poststack Facies-Controlled Seismic Inversion

For small angle incidence and weak impedance differences at the interface, the relationship between the reflection coefficient and the acoustic impedance can be expressed as [19, 27]
(1)rt=12dlnztdt,where rt and zt represent the reflection coefficient and acoustic impedance, respectively.Based on the small reflectivity assumption, equation (1) can be rewritten in the form of a matrix (Russell et al., 2006; [28])
(2)r=12Dlnz,where r=r1,r2,⋯,rN−1T, z=z1,z2,⋯,zNT, and N represent the reflectivity coefficient, acoustic impedance, and sampling number of a single seismic trace, respectively, and
(3)D=−11−11⋱⋱−11represents the difference matrix. The convolutional model of d can be written as follows:
(4)d=Wr+e=WDm+e,where m=1/2lnz represents the intermediate variable, W is the wavelet convolution matrix, and e is a noise item.Traditional poststack inversion is often carried out by minimizing the following objective function:(5)Om=dobs−WDm22,where dobs is the actual observed seismic data.Usually, the noise term obeys a Gaussian distribution. Considering the parameter correlation, the Cauchy trivariate distribution is adopted as a prior distribution.(6)Pm∣I=∏i=1N2Ψ−1/2π21+mTΦim2,where the matrix Ψ is a 3N×3N scale matrix, which has a role equivalent to the role of the correlation matrix for Gaussian statistics, N is the number of samples of seismic data, and Φi=DiTΨ−1Di.Based on the Bayesian framework, the posterior probability density function can be expressed as(7)Pm∣d,I=12πσnexp−d−WDmTd−WDm2σn2∏i=1N2Ψ−1/2π21+mTΦim2,where σn is the variance of the noise and T represents the transposition. Then, the objective function can be defined as
(8)Jm=d−WDmTd−WDm+2σn2∑i=1Nln1+mTΦim.Due to the frequency band-limited characteristic of the seismic data, directly solving for the objective function is strongly unstable [29–31]. In addition, incorrectly estimated low-frequency components may lead to erroneous inversion results [24, 32]. In the proposed study, the low-pass filtering process is expressed as a matrix, and the L2 norm is employed to minimize the difference between the low-frequency components of the inversion result and the logging data. The mathematical expression of the facies-constrained function (i.e., the low-frequency constraint) is
(9)ConstLowFre=m¯‐fm22=αIPSIPHIPm−lnLIPTSIPHIPm−lnLIP,where αIP is the constraint weight coefficients for the acoustic impedance, and SIP is the low-pass filter matrices of the acoustic impedance terms, which perform the low-pass filtering of the inversion parameters. HIP is the integral matrices of the acoustic impedance and can be expressed as
(10)HIP=10⋯011⋯0⋮⋮⋯⋮11⋯1N×N,(11)lnLIP=lnIP2IP1lnIP3IP1⋯lnIPN+1IP11×NT,where SIPHIPm represents the low-frequency components of the natural logarithms of the inverted results of the acoustic impedance. This constraint does not affect the middle and high-frequency components of the inversion result. Therefore, the seismic facies constraint is essentially different from the smooth regular and other constraints.Taking the construction of a low-pass filter matrix to the P-wave velocity as an example, its expression is(12)SIP=F‐1QF,where F and F‐1 are the positive and inverse Fourier transform matrices, respectively, and Q is the diagonal matrix composed of the Hanning window functions.Appending the facies-controlled constraint term to equation (9), we obtain a new objective function
(13)Jm=d−WDmTd−WDm+2σn2∑i=1Nln1+mTΦim+αIPSIPHIPm−lnLIPTSIPHIPm−lnLIP.Solving the partial derivative of the model parameterm in equation (13) and setting it equal to zero yields
(14)m=DTWTWD+2σn2Λ+αIPHIPTSIPTSIPHIP−1DTWTWD+αIPSIPHIPm−lnLIP,where Λ is an inconsistent weighted matrix expressed as
(15)Λ=∑i=1K2Φixy1+mTΦim,where ·xy represents the element of row x and column y of the matrix.


### 2.2. Prestack Exact Zoeppritz Equation Seismic Inversion

Prestack inversion primarily utilized the Zoeppritz equation [33], based on the generalized linear inversion theory, assuming V=VP,VS,ρ are the parameters to be inverted, where VP, VS, and ρ are the P-wave velocity, S-wave velocity, and density terms, respectively. Then, the inversion target function was established with the least squares method
(16)fV=SV−D=∑i=0n−1SViΔ−Di2⟶min,where D is the actual angle gather of the seismic data, SVΔ=W∗RV represents the expected seismic model response, W is the wavelet Toeplitz matrix, and RV is the P-wave reflection coefficient variation with the angle calculated by the Zoeppritz equation.The Taylor expansion ofSVΔ=W∗RV at the initial model and the omitting of higher-order terms above the second order are as follows:
(17)SVΔ=W∗RV=SVi+∑k=0n−1ΔVk∂Si∂Vk.Substituting equation (17) into equation (16) and letting GV=∑k=0n−1∂Si/∂Vk, Δdi=DVi−SVi, the inversion objective functional transforms to
(18)fV=∑i=0n−1GVΔV−Δd2=∑i=0n−1GVΔV−ΔdTGVΔV−Δd⟶min.To obtain the minimum value of equation (18), which are the first-order partial derivatives of the parameter variations at both ends of equation (18), set it equal to zero, and these yield
(19)f′ΔV=2GTGΔV−2GTΔd=0,and the parameter iteration step-length is
(20)ΔV=GTG−1GTΔd,and the inversion parameters can be calculated iteratively according to the following equation:
(21)Vi=Vi−1+ΔV=Vi−1+GTG−1GTΔd.Because the first-order partial derivative matrix (Jacobi matrix) is a large sparse ill-conditioned matrix, the condition number ofGTG is larger, which makes the inversion unstable. To improve the stability of the iteration equation, based on the Tikhonov regularization idea, a small positive number (regularization factor) is added to the diagonal line of the ill-conditioned matrix to make the equation invertible; that is, equation (20) becomes
(22)Vi=Vi−1+GTG+λΙ−1GTΔd,and the regularization factor can be obtained by using the covariance matrix of the parameters according to Huang and Luo [34] and Zhang et al. [19]
(23)λ=σ2Cx,where σ2 is the adjusting factor, and σ2 is positively correlated with the noise level. Cx is the covariance matrix of the parameter perturbation.Then, the inverse iteration equation (22) can now be expressed as
(24)Vi=Vi−1+GTG+σ2Cx2I−1GTΔd,where GV=∑k=0n−1∂Si/∂Vk is derived from the above generalized linear inversion theory, and for S=W∗R, the wavelet matrix W is independent of the parameter V to be solved, which yields
(25)GV=W∗∂R1∂V1∂R1∂V2⋯∂R1∂VN∂R2∂V1∂R2∂V2⋯∂R2∂VN⋮⋮⋱⋮∂RN∂V1∂RN∂V2⋯∂RN∂VN=W∗∂R1∂VP1∂R1∂VS1∂R1∂ρ1∂R1∂VP2∂R1∂VS2∂R1∂ρ2000⋯00000∂R2∂VP2∂R2∂VS2∂R2∂ρ2∂R2∂VP3∂R2∂VS3∂R2∂ρ3⋯00⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯000⋯⋯0∂RN−1∂VPN−1∂RN−1∂VSN−1∂RN−1∂ρN−1∂RN−1∂VPN∂RN−1∂VSN∂RN−1∂ρN000⋯⋯⋯000∂RN∂VPN∂RN∂VSN∂RN∂ρN.The different ways of solving the Jacobi matrixG are the main differences of the different prestack inversion methods. The traditional prestack inversion method uses the Aki-Richards and other approximation equations to derive the partial derivatives of the reflection coefficients relative to the parameter V. The exact Zoeppritz equation inversion in this paper directly uses the partial derivatives of the reflection coefficients in the exact Zoeppritz equation to solve equation (24).The Zoeppritz equation can be written in matrix form(26)AR=B,where
(27)A=sinα1cosβ1−sinα2cosβ2cosα1−sinβ1cosα2sinβ2cos2β1−vs1vp1sin2β1−ρ2⋅vp2ρ1⋅vp1cos2β2−ρ2⋅vs2ρ1⋅vp1sin2β2sin2α1vp1vs1cos2β1ρ2⋅vs22⋅vp1ρ1⋅vs12⋅vp2sin2α2−ρ2⋅vp1⋅vs2ρ1⋅vs12cos2β2,B=−sinα1cosα1−cos2β1sin2α1.Equation (25) can be written as a display form
(28)R=A−1B,by calculating the derivative of R in equation (27) for V, which yields
(29)∂R∂Vpi=A−1−∂A∂VpiA−1B+∂B∂Vpi,∂R∂Vpi+1=A−1−∂A∂Vpi+1A−1B+∂B∂Vpi+1,∂R∂Vsi=A−1−∂A∂VsiA−1B+∂B∂Vsi,∂R∂Vsi+1=A−1−∂A∂Vsi+1A−1B+∂B∂Vsi+1,∂R∂ρi=A−1−∂A∂ρiA−1B+∂B∂ρi,∂R∂ρi+1=A−1−∂A∂ρi+1A−1B+∂B∂ρi+1.MatrixG can be obtained by calculating the partial derivatives of the parameters in matrices A and B, respectively. The matrix G obtained during each iteration is brought into equation (23) to solve the inversion result after the iteration, and the iteration is terminated after reaching a certain number of iterations.


## 2.1. Poststack Facies-Controlled Seismic Inversion

For small angle incidence and weak impedance differences at the interface, the relationship between the reflection coefficient and the acoustic impedance can be expressed as [19, 27]
(1)rt=12dlnztdt,where rt and zt represent the reflection coefficient and acoustic impedance, respectively.Based on the small reflectivity assumption, equation (1) can be rewritten in the form of a matrix (Russell et al., 2006; [28])
(2)r=12Dlnz,where r=r1,r2,⋯,rN−1T, z=z1,z2,⋯,zNT, and N represent the reflectivity coefficient, acoustic impedance, and sampling number of a single seismic trace, respectively, and
(3)D=−11−11⋱⋱−11represents the difference matrix. The convolutional model of d can be written as follows:
(4)d=Wr+e=WDm+e,where m=1/2lnz represents the intermediate variable, W is the wavelet convolution matrix, and e is a noise item.Traditional poststack inversion is often carried out by minimizing the following objective function:(5)Om=dobs−WDm22,where dobs is the actual observed seismic data.Usually, the noise term obeys a Gaussian distribution. Considering the parameter correlation, the Cauchy trivariate distribution is adopted as a prior distribution.(6)Pm∣I=∏i=1N2Ψ−1/2π21+mTΦim2,where the matrix Ψ is a 3N×3N scale matrix, which has a role equivalent to the role of the correlation matrix for Gaussian statistics, N is the number of samples of seismic data, and Φi=DiTΨ−1Di.Based on the Bayesian framework, the posterior probability density function can be expressed as(7)Pm∣d,I=12πσnexp−d−WDmTd−WDm2σn2∏i=1N2Ψ−1/2π21+mTΦim2,where σn is the variance of the noise and T represents the transposition. Then, the objective function can be defined as
(8)Jm=d−WDmTd−WDm+2σn2∑i=1Nln1+mTΦim.Due to the frequency band-limited characteristic of the seismic data, directly solving for the objective function is strongly unstable [29–31]. In addition, incorrectly estimated low-frequency components may lead to erroneous inversion results [24, 32]. In the proposed study, the low-pass filtering process is expressed as a matrix, and the L2 norm is employed to minimize the difference between the low-frequency components of the inversion result and the logging data. The mathematical expression of the facies-constrained function (i.e., the low-frequency constraint) is
(9)ConstLowFre=m¯‐fm22=αIPSIPHIPm−lnLIPTSIPHIPm−lnLIP,where αIP is the constraint weight coefficients for the acoustic impedance, and SIP is the low-pass filter matrices of the acoustic impedance terms, which perform the low-pass filtering of the inversion parameters. HIP is the integral matrices of the acoustic impedance and can be expressed as
(10)HIP=10⋯011⋯0⋮⋮⋯⋮11⋯1N×N,(11)lnLIP=lnIP2IP1lnIP3IP1⋯lnIPN+1IP11×NT,where SIPHIPm represents the low-frequency components of the natural logarithms of the inverted results of the acoustic impedance. This constraint does not affect the middle and high-frequency components of the inversion result. Therefore, the seismic facies constraint is essentially different from the smooth regular and other constraints.Taking the construction of a low-pass filter matrix to the P-wave velocity as an example, its expression is(12)SIP=F‐1QF,where F and F‐1 are the positive and inverse Fourier transform matrices, respectively, and Q is the diagonal matrix composed of the Hanning window functions.Appending the facies-controlled constraint term to equation (9), we obtain a new objective function
(13)Jm=d−WDmTd−WDm+2σn2∑i=1Nln1+mTΦim+αIPSIPHIPm−lnLIPTSIPHIPm−lnLIP.Solving the partial derivative of the model parameterm in equation (13) and setting it equal to zero yields
(14)m=DTWTWD+2σn2Λ+αIPHIPTSIPTSIPHIP−1DTWTWD+αIPSIPHIPm−lnLIP,where Λ is an inconsistent weighted matrix expressed as
(15)Λ=∑i=1K2Φixy1+mTΦim,where ·xy represents the element of row x and column y of the matrix.


## 2.2. Prestack Exact Zoeppritz Equation Seismic Inversion

Prestack inversion primarily utilized the Zoeppritz equation [33], based on the generalized linear inversion theory, assuming V=VP,VS,ρ are the parameters to be inverted, where VP, VS, and ρ are the P-wave velocity, S-wave velocity, and density terms, respectively. Then, the inversion target function was established with the least squares method
(16)fV=SV−D=∑i=0n−1SViΔ−Di2⟶min,where D is the actual angle gather of the seismic data, SVΔ=W∗RV represents the expected seismic model response, W is the wavelet Toeplitz matrix, and RV is the P-wave reflection coefficient variation with the angle calculated by the Zoeppritz equation.The Taylor expansion ofSVΔ=W∗RV at the initial model and the omitting of higher-order terms above the second order are as follows:
(17)SVΔ=W∗RV=SVi+∑k=0n−1ΔVk∂Si∂Vk.Substituting equation (17) into equation (16) and letting GV=∑k=0n−1∂Si/∂Vk, Δdi=DVi−SVi, the inversion objective functional transforms to
(18)fV=∑i=0n−1GVΔV−Δd2=∑i=0n−1GVΔV−ΔdTGVΔV−Δd⟶min.To obtain the minimum value of equation (18), which are the first-order partial derivatives of the parameter variations at both ends of equation (18), set it equal to zero, and these yield
(19)f′ΔV=2GTGΔV−2GTΔd=0,and the parameter iteration step-length is
(20)ΔV=GTG−1GTΔd,and the inversion parameters can be calculated iteratively according to the following equation:
(21)Vi=Vi−1+ΔV=Vi−1+GTG−1GTΔd.Because the first-order partial derivative matrix (Jacobi matrix) is a large sparse ill-conditioned matrix, the condition number ofGTG is larger, which makes the inversion unstable. To improve the stability of the iteration equation, based on the Tikhonov regularization idea, a small positive number (regularization factor) is added to the diagonal line of the ill-conditioned matrix to make the equation invertible; that is, equation (20) becomes
(22)Vi=Vi−1+GTG+λΙ−1GTΔd,and the regularization factor can be obtained by using the covariance matrix of the parameters according to Huang and Luo [34] and Zhang et al. [19]
(23)λ=σ2Cx,where σ2 is the adjusting factor, and σ2 is positively correlated with the noise level. Cx is the covariance matrix of the parameter perturbation.Then, the inverse iteration equation (22) can now be expressed as
(24)Vi=Vi−1+GTG+σ2Cx2I−1GTΔd,where GV=∑k=0n−1∂Si/∂Vk is derived from the above generalized linear inversion theory, and for S=W∗R, the wavelet matrix W is independent of the parameter V to be solved, which yields
(25)GV=W∗∂R1∂V1∂R1∂V2⋯∂R1∂VN∂R2∂V1∂R2∂V2⋯∂R2∂VN⋮⋮⋱⋮∂RN∂V1∂RN∂V2⋯∂RN∂VN=W∗∂R1∂VP1∂R1∂VS1∂R1∂ρ1∂R1∂VP2∂R1∂VS2∂R1∂ρ2000⋯00000∂R2∂VP2∂R2∂VS2∂R2∂ρ2∂R2∂VP3∂R2∂VS3∂R2∂ρ3⋯00⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯⋯000⋯⋯0∂RN−1∂VPN−1∂RN−1∂VSN−1∂RN−1∂ρN−1∂RN−1∂VPN∂RN−1∂VSN∂RN−1∂ρN000⋯⋯⋯000∂RN∂VPN∂RN∂VSN∂RN∂ρN.The different ways of solving the Jacobi matrixG are the main differences of the different prestack inversion methods. The traditional prestack inversion method uses the Aki-Richards and other approximation equations to derive the partial derivatives of the reflection coefficients relative to the parameter V. The exact Zoeppritz equation inversion in this paper directly uses the partial derivatives of the reflection coefficients in the exact Zoeppritz equation to solve equation (24).The Zoeppritz equation can be written in matrix form(26)AR=B,where
(27)A=sinα1cosβ1−sinα2cosβ2cosα1−sinβ1cosα2sinβ2cos2β1−vs1vp1sin2β1−ρ2⋅vp2ρ1⋅vp1cos2β2−ρ2⋅vs2ρ1⋅vp1sin2β2sin2α1vp1vs1cos2β1ρ2⋅vs22⋅vp1ρ1⋅vs12⋅vp2sin2α2−ρ2⋅vp1⋅vs2ρ1⋅vs12cos2β2,B=−sinα1cosα1−cos2β1sin2α1.Equation (25) can be written as a display form
(28)R=A−1B,by calculating the derivative of R in equation (27) for V, which yields
(29)∂R∂Vpi=A−1−∂A∂VpiA−1B+∂B∂Vpi,∂R∂Vpi+1=A−1−∂A∂Vpi+1A−1B+∂B∂Vpi+1,∂R∂Vsi=A−1−∂A∂VsiA−1B+∂B∂Vsi,∂R∂Vsi+1=A−1−∂A∂Vsi+1A−1B+∂B∂Vsi+1,∂R∂ρi=A−1−∂A∂ρiA−1B+∂B∂ρi,∂R∂ρi+1=A−1−∂A∂ρi+1A−1B+∂B∂ρi+1.MatrixG can be obtained by calculating the partial derivatives of the parameters in matrices A and B, respectively. The matrix G obtained during each iteration is brought into equation (23) to solve the inversion result after the iteration, and the iteration is terminated after reaching a certain number of iterations.


## 3. Synthetic Data Examples

To verify the feasibility and applicability of the proposed integrated inversion method for estimating the carbonate reservoirs under extremely thick salts, we use synthetic seismic data to test the inversion strategy. The logging and seismic data used are from the Lula oil field, Santos Basin, Brazil. Figure1 shows the inverted results of the different inversion methods; Figure 1(a) shows the well logging curves of the study area, where the black solid line is the actual well curve, and the dotted line is the low-frequency components; Figures 1(b) and 1(c) represent the smooth regular constrained inversion with values of 103 and 107, respectively. The black solid line is for the well curve, and the black dotted line is for its low-frequency components; the red solid line is for the inversion results, and the blue dotted line is for its low-frequency components. We find that under the condition of small weight, although the details of the inversion result are valuable, the low-frequency background deviation is more severe, resulting in the distortion of the inversion result (Figure 1(b)). Figure 1(c) shows the inversion result under the large weight. The low-frequency component of the inversion result can be well matched with the actual well, and the low frequency of the actual well is restored. However, the smoothing effect of the inversion result is obvious, and the inversion resolution is difficult to guarantee. Figures 1(d) and 1(e) represent the facies-controlled constrained inversion result with facies weight values of 103 and 106, respectively. The black solid line is for the well curve, and the black dotted line is for its low-frequency components; the red solid line is for the inversion results, and the blue dotted line is for its low-frequency components. The inversion results show that the inversion effect is minimally affected by the different facies-controlled weighting factors and the inversion resolution and fidelity are high, which can be accurate to the measured data.Figure 1
(a) The well logging curve of the study area (black solid line for the well curve, dotted line for its low-frequency components). (b) and (c) represent the smooth regular constrained inversion with a value of 10-3 and 10-7, respectively (black solid line for the well curve, black dotted line for its low-frequency components; red solid line for the inversion results, blue dotted line for its low-frequency components). (d) and (e) represent the facies-controlled constrained inversion result with facies weight values of 103 and 106, respectively (black solid line for the well curve, black dotted line for its low-frequency components; red solid line for the inversion results, blue dotted line for its low-frequency components).Figure2 shows the inverted results of the Aki approximate equation inversion method and the exact Zoeppritz equation inversion method. We can clearly observe that the exact Zoeppritz equation inversion method has higher resolution results, especially the density inversion result, which is better than the Aki approximated equation inversion.Figure 2
Synthetic test of the different prestack inversion methods. The black solid line represents the real data; the green solid line represents the initial model; the blue line represents the Aki approximate equation inversion result; the red line is the exact Zoeppritz equation prestack inversion result.The test results show that the proposed integrated method has a higher resolution and correlation coefficient than the traditional inversion method, which can prove the superiority of the proposed method.


## 4. Real Data Application

The Lula oil field is located on the San Paulo plateau in the Santos Basin in southeastern Brazil, as shown in Figure3(a). The basin is adjacent to Campos Basin in the north and bounded by the Charcot Sea seamount in the south, which covers an area of 32.7×104km2. The study area is Block S of the Lula oil field. The average water depth of the oilfield is 2000 m-2400 m, and the oil-bearing area is 802 km2. As shown in Figure 3(b), the basic structure of the oil field is a fault anticline structure of basement uplift, and it is a set of carbonate rocks under the cover of a thick salt rock.Figure 3
(a) Location map of the Lula oil field, Santos Basin, Brazil. (b) Schematic diagram of the geological environment in the Lula oil field.
(a)
(b)As shown in Figure4(a), the structural characteristics of the study area show a pattern of high east-west flanks and low troughs in the middle. There are two NE-SW uplift areas in the study area, the middle of which is an intraplatform depression. The east side is a steep slope, and the west side is a gentle slope. The drilling data reveal that the structural uplift area is the main oil and gas accumulation area. In the study area, there are two large northeast-southwest faults and more than 30 other faults (Figure 4(b)).Figure 4
(a) The paleogeography of the target area and (b) fault development profile in the study area.
(a)
(b)The Santos Basin, as a typical passive continental margin basin in the Atlantic Ocean, has three stages of tectonic evolution that have an impact on the formation of oil and gas resources, and the three stages include the rift stage, the transition stage, and the drift stage [35]. A set of thick salt rocks in the middle divides the basin into two sequences, the subsalt sequence and suprasalt sequence. As shown in Figure 5, the subsalt strata in the study area can be further divided into the depression stage and a rift stage, and the depression stage can be divided into an upper depression sequence, middle depression sequence, and a lower depression sequence, corresponding to the BVE100, BVE200, and BVE300 layers, respectively. The rift stage can be divided into an upper rift sequence, middle rift sequence, and lower rift sequence, corresponding to the ITP, PIC, and CAM layers, respectively.Figure 5
Map of the presalt stratigraphy in the Santos Basin.At the top of the BVE group, there is a layer of salt rock with a thickness of more than 2000 meters and a layer of high-velocity gypsum that is 10-20 meters thick in the middle. A large number of stromatolite have been developed in the tectonic uplift zone, and the grainstone and spherulite have been developed in the lower slope zone. These lithological properties are indicative of a good kind of reservoir; in the BVE300 layer, the ITP group also developed muddy limestone, layered limestone, and wackestone. Their physical properties are indicative of a poorer type of reservoir, and these are the second and third types of reservoirs in the BVE group. On the top layer of the BVE layer, a set of microscopic reef facies is deposited. The reef beach is built on the hard reef base of the eastern uplift and is distributed in strips. The gentle slope zone has relatively weak water dynamics and mainly develops granular limestone and layered limestone.From the well-tie seismic cross section in the study area, the top of the reef has mostly strong reflection, with obvious hilly characteristics, and the interior has chaotic reflection or weak reflection. The top of the organic beach is relatively flat, with good continuity, disordered internal reflection, and weak energy. The reflection characteristics of the eastern uplift area are vague and difficult to distinguish. It can be defined as the reef-beach complex, as shown in Figure6.Figure 6
The well-tie seismic cross section in the study area.In the study area, the reef-shoal facies and intershoal low facies developed in the subsalt carbonate rocks, which have strong reservoir heterogeneity and great lateral variation. Among them, the reef-shoal facies mainly developed stromatolite, granular limestone, and granulated limestone, with good reservoir properties andporosity>10%. The muddy rocks and ash-bearing limestones are mainly developed in the intershoal facies, and the physical properties of the reservoir are relatively poor, with porosity≤6%. According to the observations of the core and outcrop, the reservoirs are divided into three types. The reservoirs of type I are reef beaches, the reservoirs of type II are beach margins, and the reservoirs of type III are intershoal. Based on the above classification criteria, the 14 wells in the study area are divided into the three types of reservoirs. These reservoirs can be further divided into three types according to their porosity, where type I is a microbial reef shoal reservoir (porosity>10%), type II is an intershoal lows and intraplatform shallow reservoir (porosity 5%~10%), and type III is a semideep lake-deep lake marl reservoir (porosity<5%). The type I reservoirs are characterized by high porosity, low P/S wave velocity ratio, low Poisson ratio, and low Lamé coefficient for multiplication density. The type III reservoirs have opposite characteristics than type I, and type II reservoirs have characteristics that are centered between type I and type III, as shown in Figure 7.Figure 7
Three reservoir classification types, porosity andVP/VS and Poisson’s ratio intersection diagram.The study found that the lithology types of the subsalt carbonate rocks in the Lula oil field are rich and diverse, and they have both of the lithological characteristics of large carbonate rocks and unique reservoir characteristics. In particular, the unique stromatolite reservoirs in the area formed the organic reef beach sedimentary system of the BVE group.Through the analysis of logging curves in the study area, we realize the logging responses of the different lithologies. By using the velocity parameters and natural gamma intersection analysis, the reservoirs and nonreservoirs can be identified well, and the high-quality reservoirs and poor reservoirs can be further identified. The results of the intersection analysis of the rock physical parameters show that the P-wave velocity (acoustic impedance) can effectively distinguish the limestone and mudstone, gypsum, etc., which proves that the velocity is sensitive to the reservoir, but the lithology of the limestone and dolomite is superimposed (Figure8). From the analysis of the key petrophysical intersection of the logging data, it can be seen that using the density-P/S wave velocity ratio and density-Poisson ratio to build a rock physics template of the reservoir can better distinguish the three types of reservoirs, as shown in Figure 9. In the different target layers, the average density of type I reservoirs ranges from 2.36 to 2.62 g/cm3, the average P/S wave velocity ratio ranges from 1.4 to 1.62, and the average Poisson ratio ranges from 0.08 to 0.18. The average density of type II reservoirs ranges from 2.55 to 2.72 g/cm3, the average P/S wave velocity ratio ranges from 1.5 to 1.8, and the average Poisson ratio ranges from 0.12 to 0.29. The average density of type III reservoirs ranges from 2.65 to 2.8 g/cm3, the average P/S wave velocity ratio ranges from 1.62 to 1.87, and the average Poisson ratio ranges from 0.19 to 0.33. These parameters can provide the relevant parameter constraints for subsequent multiparameter interactive reservoir thickness predictions.Figure 8
Different lithological velocities and natural gamma intersections in the study area.Figure 9
The different parameter intersection map of the study area.


### 4.1. Poststack Facies-Controlled Inversion to Distinguish the Lithology

The target layer is buried deeply with a large number of blank reflections in the middle. The salt rocks divide the seismic profile into two sets of carbonate sedimentary systems, i.e., upper and lower salt systems. Due to the absorption of salt rock, a large amount of reflected wave energy is absorbed, and the high-frequency signal disappears. The seismic signal quality of the carbonate reservoirs under the salt is very poor, which brings great difficulties and challenges to the inversion and reservoir prediction.We intercepted a section of a two-dimensional section from the target layer for inversion testing. Well A is selected as the logging constraint. Four sequence boundaries are divided in the target interval. Figure10 shows the inversion result of well 1, where the vertical resolution is high and the lithological change trend from the high to the low to the high acoustic impedance is obvious. Because of the strong heterogeneity of the carbonate rocks, the transverse variation trend of the carbonate rocks can also be clearly seen from the section. In the log interpretation, the red label is interpreted as an oil reservoir, which corresponds well to the high acoustic impedance in the inversion profile. Generally, the distribution of the carbonate reservoirs that cannot be recognized by seismic profiles can be clearly depicted using facies-controlled inversion, and the rich geological signals contained in the weak seismic signals can be clearly observed.Figure 10
Poststack acoustic impedance inversion profile of well 1.Comparing the inversion results with the log interpretation, the lithology characteristics are in good agreement with the logging interpretation. As shown in Figure11, the high-acoustic impedance gypsum and low-acoustic impedance salt rock are clearly distinguished. The high-acoustic impedance limestone reservoirs and low-acoustic impedance limestones have distinct characteristics in their inversion profiles, and their inversions highlight the different lithological characteristics, which facilitates the prediction and interpretation of the subsequent reservoir distributions. The organic reef is one of the favorable reservoirs in the Lula oil field. The seismic profiles cannot identify the hidden organic reef characteristics, and only some weak internal reflection can be observed. The inversion section clearly depicts the overlapping characteristics of the organic reef, which makes the hidden signals clearly observed.Figure 11
The acoustic impedance inversion profiles of the different lithologies.Figure12 shows the inverted results of the organic reef facies and the bioclastic limestone of the shallow facies in the study area. It is difficult to observe the lithological changes in the seismic profiles, but lithological differences such as reef and limestone can be clearly identified in the inversion profiles. The lithology of the target stratum is complex and diverse. The limestone, dolomite, gypsum, salt rock, and mudstone are all well distributed. The lithological changes can be seen from the inversion profiles. The characteristics of the organic beach and biolimestone reservoirs clearly clarify the lithological characteristics of the reservoirs and can identify the favorable areas. It is necessary to establish a lithological interpretation model for reservoir prediction based on the characteristics of the seismic amplitude and inversion.Figure 12
The acoustic impedance inversion profile of the organic reef and bioclastic limestone reservoir.


### 4.2. Prestack Exact Zoeppritz Equation Inversion for Identifying the Favorable Reservoirs

According to the petrophysical analysis, it is difficult to identify the type of favorable reservoirs only by the poststack inversion. Because the density parameters are insensitive to the seismic data with small and medium incidence angles, it is possible to obtain the reservoir density more accurately by the prestack exact Zoeppritz equation inversion using seismic data with large incidence angles [28, 36]. Then, we used the proposed method to identify the type of favorable reservoirs. Figure 13 shows the prestack exact Zoeppritz equation inversion results of the P-wave velocity, S-wave velocity, density, and the ratio of P- and S-wave of well 4. On the whole, it is considered that the upper section of BVE100 has a set of high-quality reservoirs with good continuity. The limestone in the BVE200 interval has a weak hydrocarbon potential. The reservoir of BVE300 is well distributed, but its continuity is poor. The lateral variation of the limestone in BVE300 is large, and its distribution is discontinuous.Figure 13
The prestack exact Zoeppritz equation inversion results of the (a) P-wave velocity, (b) S-wave velocity, and (c) density of well 4.In terms of the inversion effect, the prestack exact Zoeppritz equation inversion method is more convincing than the poststack single parameter inversion. The inverted three-parameter data volume can obtain more favorable reservoir-sensitive parameters through multiattribute interactive calculations. Vertically, the oil reservoirs are divided into three layers. Well W4 has a thin reservoir thickness but good oil- and gas-bearing capacity. In Figure13, the red arrow indicates the location of the reservoir, where No. 1 and No. 3 are type I and II reservoirs, respectively, which are the oil layers, and No. 2 is a type III reservoir, which is a dry layer; the lower density between the BVE200 and BVE300 is the salt rock. By comparing the drilling test data of the wells with the inversion results, the drilling reservoir is consistent with the oil- and gas-bearing interval in the inversion results, which has proven that the proposed inversion methodology is reliable and accurate. Based on the above elastic parameters and combining with other more petrophysical-sensitive parameters to distinguish between type I and type II reservoirs.Figure14 is a comparison of the P-S wave velocity ratio and Poisson ratio inversion with the logging reservoir classification. In the log reservoir classification interpretation, red and yellow represent class I and II reservoirs, respectively, and green represents class III reservoirs. The comparison shows that the inversion results are in good agreement with the logging reservoir classification results, which indicates that the reservoir can be comprehensively divided by using multiple parameters combined with the petrophysical statistical results.Figure 14
Comparison of the P-S wave velocity ratio and Poisson ratio inversion with the logging reservoir classification.


### 4.3. The Integrated Prediction of Favorable Reservoirs

The lithology of the target layer is complex and diverse, with limestone, dolomite, gypsum, salt rock, and mudstone well distributed. To clarify the lithology characteristics of the reservoir and find the favorable areas, it is necessary to establish a lithological combination model for reservoir prediction. The prestack and poststack inversion multiparameter joint interpretation mode can reduce the number of solutions of the prediction and can provide reference and criteria for fine reservoir characterization in the study area.Figure15 is the interpretation model of the lithologic assemblage. The reef reservoirs show domal seismic reflection, and the poststack inversion results show high-speed characteristics. Therefore, the morphological characteristics of the reef reservoirs can be depicted, but the hydrocarbon-bearing characteristics of the reef reservoirs cannot be predicted. It is necessary to further integrate with the prestack inversion results and interactively interpret the results with the multiple parameters.Figure 15
Prediction models of the different lithologies and lithological combinations.In terms of the density inversion results, we found that the reef reservoirs of the different layers have medium and low density characteristics and have medium and low P-S wave velocity ratios, while the Poisson ratios and Lamé coefficient have multiplying densities. Combining these factors, the fluid-bearing characteristics of the reef reservoirs can be further estimated.The biological beach is mainly composed of bioclastic limestone, laminated limestone, and a small amount of micrite limestone. The poststack inversion results show that the P-wave velocity of the biobank reservoirs is relatively high, ranging from 5600 to 6000 m/s, and the Poisson ratio of the prestack inversion is more than 0.26. Compared with the reef reservoirs, the physical properties of such reservoirs are generally poor.In the poststack seismic inversion section, the interaction layer of limestone and marl is mostly strip-shaped on the inversion section, and the velocity is characterized by high-speed interaction with medium-low velocity. Meanwhile, it can be seen that the algal stromatolite and the globular limestone are block and narrow strip features. In the prestack inversion profiles, the algal stromatolites and spherical limestones show obvious characteristics of low P-S wave velocity ratio, and the interbedded limestone and marl show interaction characteristics of high P-S wave velocity ratio and low P-S wave velocity ratio. The P/S wave velocity ratio of this limestone reservoir is usually less than 1.8, which is the main reservoir combination model in the study area.The poststack inversion results show that the velocity of the nonreservoir and poor reservoir is relatively high, ranging from 5500 m/s to 6000 m/s. In the prestack inversion profile, the corresponding P-S wave velocity is relatively high and has poor physical properties.Based on the prestack and poststack inversion results and the multiparameter constrained reservoir thickness extraction method, we mapped the planar thickness distribution of the favorable reservoirs in BVE100, BVE200, and BVE300. The results are shown in Figure16. The favorable reservoirs of type I+II in BVE100 are concentrated in two uplift areas. The thickest reservoir is developed in well 9-LL-7 area with a thickness of 70 m, followed by well 8-LL-1D and 9-RJS-686 areas with maximum thicknesses of 50 m and 25 m, respectively, and the average thickness is 20 m (Figure 16(a)). A type III reservoir in the BVE100 member is developed locally with a small thickness. Reservoirs are concentrated in the well 9-RJS-686D area. Reservoir development in the 9-LL-12D well area is relatively small (Figure 16(b)). The type III reservoir in the BVE200 member has a small development range and thin thickness. The average thickness of the reservoir is approximately 20 m, the thickest is 30 m and the thinnest is approximately 5 m (Figure 16(c)). The development of favorable type I+II reservoirs in BVE300 is between layer BVE100 and layer BVE200. The reservoirs are concentrated in the uplift area with strip distribution and tend to pinch out in the groove and slope area with thin reservoir thickness. The maximum thickness of well 1-RJS-628A is approximately 30 m and that of well 9-LL-12D and well 8-LL-5 is approximately 15 m, with an average thickness of approximately 10 m (Figure 16(d)).Figure 16
Favorable reservoir thickness plane distribution prediction.
(a)
(b)
(c)
(d)The overall structure pattern of the study area is “two uplifts and one shallow” along the NE direction, and the northwest is a gentle slope monoclinic structure, as shown in Figure17(a). There are two NE-SW uplift areas with an interplatform depression in the middle. The east side of well 9-LL-7 is a steep slope, while the west side of well 9-RLS-686D is a gentle slope. Drilling data reveal that the structural uplift area is the main oil and gas accumulation area.Figure 17
(a) Contour map and trap distribution of the top interface structure of the BVE100 layer, (b) reef-shoal facies reservoir and structural overlap map, and (c) sedimentary facies map.
(a)
(b)
(c)There are five structural traps in the main target layer of BVE100. In the trap, wells of 7-LL-28D, 8-26D, 3-674, 8-5, 9-665, 9-7, and 9-12D all encounter oil layers, and seven wells have been drilled outside the trap. Except for the two wells of 8-LL-54D and 8-62D, all other five wells met the reservoir. The analysis shows that the structural pattern of the “two uplifts and one shallow” makes a large amount of oil and gas accumulate in the high part of the uplift. In addition to the oil and gas drilling in the structural traps, the oil reservoirs are also drilled in the high-quality reservoir development areas outside the traps. The reservoir is controlled by multiple factors, such as structure, lithology, and physical properties.The sedimentary facies map is compiled based on the seismic facies, paleo-tectonics, and well logging data, as shown in Figure17(b), and on the reef-shoal facies reservoir and structure superimposed map, as shown in Figure 17(c). It can be seen that two sets of reef-bank complexes are developed in the Lula oil field, in which the reef nucleus is the main body of the algae stratified reef, the algal stromatolite microbial limestone, and the interbedded layer of the granulite and limestone; the microbial limestone is mainly developed before the reef. The velocity distribution is relatively uniform; the post reef is mainly composed of granular limestone, as the velocity is high, and the stratification is better; the hydrodynamics of the shoal are weak, mainly distributing high-speed compact limestone and layered limestone.To further analyze the main controlling factors of the oil and gas migration in the study area, seismic fracture detection is carried out based on edge detection and ant body attributes. Figure18 shows the fracture detection results of layer BVE 100 and ITP. We find that the study area is mainly distributed with two large east-north (red line) and west-south direction faults and more than 30 faults of other sizes and scales. The fracture development area is near the two main faults, and the fault direction is mainly to the northeast.Figure 18
Seismic fracture detection results of the BVE100 and ITP layers.
(a)
(b)Figure19 compares the plane and section of the seismic fracture detection results. Fractures are widely developed on the eastern side of well 7-LL-53D and on the western side of well 7-LL-28D. There are two NE-SW fracture systems on both sides of the uplift area: the distribution of faults with different sizes, ground barriers, and anticlines. During the transition period, the sedimentation is stable, and the fault and fractures are less developed. During the rift period, the faults are widely developed, and the tectonic movement is intense. It can be seen that faults are widely developed in the Lula oil field during the rift period, and the hydrocarbon accumulation is controlled by the NE-trending Horst formed during this period. A set of source rocks in the early rift period migrated along the fault and eventually (red arrow) aggregated in the depression period, indicating that faults and fractures are the main channels for hydrocarbon migration.Figure 19
Comparison of the plane and section of the seismic fracture detection results.Three potential areas are predicted based on the reservoir, structure, and oil testing. As shown in Figure20, potential area 1 is located east of well 9-RJS-665, potential area 2 is located east of well 1-RJS-628A, and potential area 3 is located southwest of 7-LL-28D. The three areas are all located in the structural mountain ridge zone, where reservoirs are well developed, and industrial oil flow is obtained in the low part of the wells, with obvious exploration potential.Figure 20
The results of the integrated prediction of the potential area distribution.


## 4.1. Poststack Facies-Controlled Inversion to Distinguish the Lithology

The target layer is buried deeply with a large number of blank reflections in the middle. The salt rocks divide the seismic profile into two sets of carbonate sedimentary systems, i.e., upper and lower salt systems. Due to the absorption of salt rock, a large amount of reflected wave energy is absorbed, and the high-frequency signal disappears. The seismic signal quality of the carbonate reservoirs under the salt is very poor, which brings great difficulties and challenges to the inversion and reservoir prediction.We intercepted a section of a two-dimensional section from the target layer for inversion testing. Well A is selected as the logging constraint. Four sequence boundaries are divided in the target interval. Figure10 shows the inversion result of well 1, where the vertical resolution is high and the lithological change trend from the high to the low to the high acoustic impedance is obvious. Because of the strong heterogeneity of the carbonate rocks, the transverse variation trend of the carbonate rocks can also be clearly seen from the section. In the log interpretation, the red label is interpreted as an oil reservoir, which corresponds well to the high acoustic impedance in the inversion profile. Generally, the distribution of the carbonate reservoirs that cannot be recognized by seismic profiles can be clearly depicted using facies-controlled inversion, and the rich geological signals contained in the weak seismic signals can be clearly observed.Figure 10
Poststack acoustic impedance inversion profile of well 1.Comparing the inversion results with the log interpretation, the lithology characteristics are in good agreement with the logging interpretation. As shown in Figure11, the high-acoustic impedance gypsum and low-acoustic impedance salt rock are clearly distinguished. The high-acoustic impedance limestone reservoirs and low-acoustic impedance limestones have distinct characteristics in their inversion profiles, and their inversions highlight the different lithological characteristics, which facilitates the prediction and interpretation of the subsequent reservoir distributions. The organic reef is one of the favorable reservoirs in the Lula oil field. The seismic profiles cannot identify the hidden organic reef characteristics, and only some weak internal reflection can be observed. The inversion section clearly depicts the overlapping characteristics of the organic reef, which makes the hidden signals clearly observed.Figure 11
The acoustic impedance inversion profiles of the different lithologies.Figure12 shows the inverted results of the organic reef facies and the bioclastic limestone of the shallow facies in the study area. It is difficult to observe the lithological changes in the seismic profiles, but lithological differences such as reef and limestone can be clearly identified in the inversion profiles. The lithology of the target stratum is complex and diverse. The limestone, dolomite, gypsum, salt rock, and mudstone are all well distributed. The lithological changes can be seen from the inversion profiles. The characteristics of the organic beach and biolimestone reservoirs clearly clarify the lithological characteristics of the reservoirs and can identify the favorable areas. It is necessary to establish a lithological interpretation model for reservoir prediction based on the characteristics of the seismic amplitude and inversion.Figure 12
The acoustic impedance inversion profile of the organic reef and bioclastic limestone reservoir.


## 4.2. Prestack Exact Zoeppritz Equation Inversion for Identifying the Favorable Reservoirs

According to the petrophysical analysis, it is difficult to identify the type of favorable reservoirs only by the poststack inversion. Because the density parameters are insensitive to the seismic data with small and medium incidence angles, it is possible to obtain the reservoir density more accurately by the prestack exact Zoeppritz equation inversion using seismic data with large incidence angles [28, 36]. Then, we used the proposed method to identify the type of favorable reservoirs. Figure 13 shows the prestack exact Zoeppritz equation inversion results of the P-wave velocity, S-wave velocity, density, and the ratio of P- and S-wave of well 4. On the whole, it is considered that the upper section of BVE100 has a set of high-quality reservoirs with good continuity. The limestone in the BVE200 interval has a weak hydrocarbon potential. The reservoir of BVE300 is well distributed, but its continuity is poor. The lateral variation of the limestone in BVE300 is large, and its distribution is discontinuous.Figure 13
The prestack exact Zoeppritz equation inversion results of the (a) P-wave velocity, (b) S-wave velocity, and (c) density of well 4.In terms of the inversion effect, the prestack exact Zoeppritz equation inversion method is more convincing than the poststack single parameter inversion. The inverted three-parameter data volume can obtain more favorable reservoir-sensitive parameters through multiattribute interactive calculations. Vertically, the oil reservoirs are divided into three layers. Well W4 has a thin reservoir thickness but good oil- and gas-bearing capacity. In Figure13, the red arrow indicates the location of the reservoir, where No. 1 and No. 3 are type I and II reservoirs, respectively, which are the oil layers, and No. 2 is a type III reservoir, which is a dry layer; the lower density between the BVE200 and BVE300 is the salt rock. By comparing the drilling test data of the wells with the inversion results, the drilling reservoir is consistent with the oil- and gas-bearing interval in the inversion results, which has proven that the proposed inversion methodology is reliable and accurate. Based on the above elastic parameters and combining with other more petrophysical-sensitive parameters to distinguish between type I and type II reservoirs.Figure14 is a comparison of the P-S wave velocity ratio and Poisson ratio inversion with the logging reservoir classification. In the log reservoir classification interpretation, red and yellow represent class I and II reservoirs, respectively, and green represents class III reservoirs. The comparison shows that the inversion results are in good agreement with the logging reservoir classification results, which indicates that the reservoir can be comprehensively divided by using multiple parameters combined with the petrophysical statistical results.Figure 14
Comparison of the P-S wave velocity ratio and Poisson ratio inversion with the logging reservoir classification.


## 4.3. The Integrated Prediction of Favorable Reservoirs

The lithology of the target layer is complex and diverse, with limestone, dolomite, gypsum, salt rock, and mudstone well distributed. To clarify the lithology characteristics of the reservoir and find the favorable areas, it is necessary to establish a lithological combination model for reservoir prediction. The prestack and poststack inversion multiparameter joint interpretation mode can reduce the number of solutions of the prediction and can provide reference and criteria for fine reservoir characterization in the study area.Figure15 is the interpretation model of the lithologic assemblage. The reef reservoirs show domal seismic reflection, and the poststack inversion results show high-speed characteristics. Therefore, the morphological characteristics of the reef reservoirs can be depicted, but the hydrocarbon-bearing characteristics of the reef reservoirs cannot be predicted. It is necessary to further integrate with the prestack inversion results and interactively interpret the results with the multiple parameters.Figure 15
Prediction models of the different lithologies and lithological combinations.In terms of the density inversion results, we found that the reef reservoirs of the different layers have medium and low density characteristics and have medium and low P-S wave velocity ratios, while the Poisson ratios and Lamé coefficient have multiplying densities. Combining these factors, the fluid-bearing characteristics of the reef reservoirs can be further estimated.The biological beach is mainly composed of bioclastic limestone, laminated limestone, and a small amount of micrite limestone. The poststack inversion results show that the P-wave velocity of the biobank reservoirs is relatively high, ranging from 5600 to 6000 m/s, and the Poisson ratio of the prestack inversion is more than 0.26. Compared with the reef reservoirs, the physical properties of such reservoirs are generally poor.In the poststack seismic inversion section, the interaction layer of limestone and marl is mostly strip-shaped on the inversion section, and the velocity is characterized by high-speed interaction with medium-low velocity. Meanwhile, it can be seen that the algal stromatolite and the globular limestone are block and narrow strip features. In the prestack inversion profiles, the algal stromatolites and spherical limestones show obvious characteristics of low P-S wave velocity ratio, and the interbedded limestone and marl show interaction characteristics of high P-S wave velocity ratio and low P-S wave velocity ratio. The P/S wave velocity ratio of this limestone reservoir is usually less than 1.8, which is the main reservoir combination model in the study area.The poststack inversion results show that the velocity of the nonreservoir and poor reservoir is relatively high, ranging from 5500 m/s to 6000 m/s. In the prestack inversion profile, the corresponding P-S wave velocity is relatively high and has poor physical properties.Based on the prestack and poststack inversion results and the multiparameter constrained reservoir thickness extraction method, we mapped the planar thickness distribution of the favorable reservoirs in BVE100, BVE200, and BVE300. The results are shown in Figure16. The favorable reservoirs of type I+II in BVE100 are concentrated in two uplift areas. The thickest reservoir is developed in well 9-LL-7 area with a thickness of 70 m, followed by well 8-LL-1D and 9-RJS-686 areas with maximum thicknesses of 50 m and 25 m, respectively, and the average thickness is 20 m (Figure 16(a)). A type III reservoir in the BVE100 member is developed locally with a small thickness. Reservoirs are concentrated in the well 9-RJS-686D area. Reservoir development in the 9-LL-12D well area is relatively small (Figure 16(b)). The type III reservoir in the BVE200 member has a small development range and thin thickness. The average thickness of the reservoir is approximately 20 m, the thickest is 30 m and the thinnest is approximately 5 m (Figure 16(c)). The development of favorable type I+II reservoirs in BVE300 is between layer BVE100 and layer BVE200. The reservoirs are concentrated in the uplift area with strip distribution and tend to pinch out in the groove and slope area with thin reservoir thickness. The maximum thickness of well 1-RJS-628A is approximately 30 m and that of well 9-LL-12D and well 8-LL-5 is approximately 15 m, with an average thickness of approximately 10 m (Figure 16(d)).Figure 16
Favorable reservoir thickness plane distribution prediction.
(a)
(b)
(c)
(d)The overall structure pattern of the study area is “two uplifts and one shallow” along the NE direction, and the northwest is a gentle slope monoclinic structure, as shown in Figure17(a). There are two NE-SW uplift areas with an interplatform depression in the middle. The east side of well 9-LL-7 is a steep slope, while the west side of well 9-RLS-686D is a gentle slope. Drilling data reveal that the structural uplift area is the main oil and gas accumulation area.Figure 17
(a) Contour map and trap distribution of the top interface structure of the BVE100 layer, (b) reef-shoal facies reservoir and structural overlap map, and (c) sedimentary facies map.
(a)
(b)
(c)There are five structural traps in the main target layer of BVE100. In the trap, wells of 7-LL-28D, 8-26D, 3-674, 8-5, 9-665, 9-7, and 9-12D all encounter oil layers, and seven wells have been drilled outside the trap. Except for the two wells of 8-LL-54D and 8-62D, all other five wells met the reservoir. The analysis shows that the structural pattern of the “two uplifts and one shallow” makes a large amount of oil and gas accumulate in the high part of the uplift. In addition to the oil and gas drilling in the structural traps, the oil reservoirs are also drilled in the high-quality reservoir development areas outside the traps. The reservoir is controlled by multiple factors, such as structure, lithology, and physical properties.The sedimentary facies map is compiled based on the seismic facies, paleo-tectonics, and well logging data, as shown in Figure17(b), and on the reef-shoal facies reservoir and structure superimposed map, as shown in Figure 17(c). It can be seen that two sets of reef-bank complexes are developed in the Lula oil field, in which the reef nucleus is the main body of the algae stratified reef, the algal stromatolite microbial limestone, and the interbedded layer of the granulite and limestone; the microbial limestone is mainly developed before the reef. The velocity distribution is relatively uniform; the post reef is mainly composed of granular limestone, as the velocity is high, and the stratification is better; the hydrodynamics of the shoal are weak, mainly distributing high-speed compact limestone and layered limestone.To further analyze the main controlling factors of the oil and gas migration in the study area, seismic fracture detection is carried out based on edge detection and ant body attributes. Figure18 shows the fracture detection results of layer BVE 100 and ITP. We find that the study area is mainly distributed with two large east-north (red line) and west-south direction faults and more than 30 faults of other sizes and scales. The fracture development area is near the two main faults, and the fault direction is mainly to the northeast.Figure 18
Seismic fracture detection results of the BVE100 and ITP layers.
(a)
(b)Figure19 compares the plane and section of the seismic fracture detection results. Fractures are widely developed on the eastern side of well 7-LL-53D and on the western side of well 7-LL-28D. There are two NE-SW fracture systems on both sides of the uplift area: the distribution of faults with different sizes, ground barriers, and anticlines. During the transition period, the sedimentation is stable, and the fault and fractures are less developed. During the rift period, the faults are widely developed, and the tectonic movement is intense. It can be seen that faults are widely developed in the Lula oil field during the rift period, and the hydrocarbon accumulation is controlled by the NE-trending Horst formed during this period. A set of source rocks in the early rift period migrated along the fault and eventually (red arrow) aggregated in the depression period, indicating that faults and fractures are the main channels for hydrocarbon migration.Figure 19
Comparison of the plane and section of the seismic fracture detection results.Three potential areas are predicted based on the reservoir, structure, and oil testing. As shown in Figure20, potential area 1 is located east of well 9-RJS-665, potential area 2 is located east of well 1-RJS-628A, and potential area 3 is located southwest of 7-LL-28D. The three areas are all located in the structural mountain ridge zone, where reservoirs are well developed, and industrial oil flow is obtained in the low part of the wells, with obvious exploration potential.Figure 20
The results of the integrated prediction of the potential area distribution.


## 5. Conclusions

We have proposed a strategy for the integrated prediction of carbonate reservoirs under extremely thick salts via poststack facies-controlled and prestack exact Zoeppritz equation inversion in the Santos Basin of Brazil. To reduce the errors of the poststack inversion caused by the low-frequency loss of seismic data due to the extremely thick salts, we use a low-pass filter matrix to establish the facies-controlled constraints to ensure the accuracy of the low-frequency information of the poststack inversion results and improve the robustness of the inversion. Meanwhile, we make the inversion results conform to the geological characteristics. To make full use of the seismic data with a large incident angle and solve the problem that poststack inversion cannot distinguish the reservoir types, we introduce the exact Zoeppritz prestack inversion method with an adaptive regularization to reduce the error of the density inversion to ensure the reliability of predicting favorable reservoirs.Combining the application of prestack and poststack inversions based on the characteristics of the favorable Lula carbonate reservoirs can not only characterize the shape of reefs and other reservoirs but the distribution of reservoirs also be depicted accurately, which provides important support for exploration. The quantitative reservoir prediction results show that the BVE100 reservoir in the Lula area is concentrated in two uplift areas and western gentle slope areas, with large reservoir thickness, fast lateral change, and obvious exploration potential. Future research should combine the drilling data and the oil test and exploration data to the calibrated prediction results to reduce the risk of overlooking carbonate reservoirs under extremely thick salts in the exploration of the Santos Basin. The integrated methods enable us to predict the distribution of the favorable carbonate reservoirs under extremely thick salts, which can be used as an important guide for the further exploration of carbonate reservoirs under extremely thick salts.



---

*Source: 6205185-2020-01-23.xml*