A two-step method to apply Xu–Payne multi-porosity model to estimate pore type from seismic data for carbonate reservoirs

Carbonate reservoirs exhibit strong heterogeneity in the distribution of pore types that can be quantitatively characterized by applying Xu–Payne multi-porosity model. However, there are some prerequisites to this model the porosity and saturation need to be provided. In general, these application conditions are difficult to satisfy for seismic data. In order to overcome this problem, we present a two-step method to estimate the porosity and saturation and pore type of carbonate reservoirs from seismic data. In step one, the pore space of the carbonate reservoir is equivalent to a single-porosity system with an effective pore aspect ratio; then, a 3D rock-physics template (RPT) is established through the Gassmann’s equations and effective medium models; and then, the effective aspect ratio of pore, porosity and fluid saturation are simultaneously estimated from the seismic data based on 3D RPT. In step two, the pore space of the carbonate reservoir is equivalent to a triple-porosity system. Combined with the inverted porosity and saturation in the first step, the porosities of three pore types can be inverted from the seismic elastic properties. The application results indicate that our method can obtain accurate physical properties consistent with logging data and ensure the reliability of characterization of pore type.


Introduction
A significant portion of the world's oil and gas reserves resides in carbonate rocks (Alotaibi et al. 2010;Bagrintseva 2015). Unlike conventional siliciclastic rocks, carbonate rocks undergo complicated depositional and diagenetic processes. They have strongly heterogeneous pore system, such as vugs, moldic, intergranular, intragranular, intercrystal pores and micro-cracks. Complex pore structure has a very significant effect on elastic properties. Laboratory measurement results have demonstrated that the difference in pore type can lead to a great difference in P-wave velocity of up to 40% at certain same porosity (Anselmetti and Eberli 1993;Xu and Payne 2009). Such high heterogeneity in these porous rocks brings great uncertainties for carbonate reservoir prediction and makes oil and gas recovery difficult during development, which has been identified as the crucial parameter needed to resolve the scatter in seismic and petrophysical relationships (Anselmetti and Eberli 1999;Weger et al. 2009;Sun et al. 2015;Jin et al. 2017). Therefore, the pore-type estimation from seismic data is helpful to describe the heterogeneity of carbonate reservoirs, predict relatively high-quality reservoir, and further improve the accuracy of reservoir prediction.
Quantitative estimation of pore type from geophysical measurements is frequently attained by using various rock-physics models that relate rock elastic properties to geophysical observations (Kumar and Han 2005;Xu and Payne 2009;Zhao et al. 2013;Yuan et al. 2018). The idea of an effective pore aspect ratio has been used in rockphysics modeling and shear velocity prediction for complex reservoirs (Sams and Andrea 2001;Li et al. 2013). Sams and Andrea (2001) found that the results of numerical modeling with a spectrum of pore aspect ratios can be matched using the same porosity and a single aspect ratio in a clastic reservoir by numerical simulation. Li et al. Edited by Jie Hao and Xiu-Qiu Peng (2013) applied the differential effective medium (DEM) model to invert the effective pore aspect ratio from P-wave velocity and used it to predict shear wave velocity. Xu and Payne (2009) developed a new model to estimate the porosities of various pore types from well-logging data under the condition that the reservoir porosity and saturation are provided. These application conditions usually can be met for well-logging data because the porosity and saturation might be estimated from other logging data such as resistivity and neutron data. However, for seismic data there are not enough data to be applied. Some researchers attempt to apply Xu-Payne model (Xu and Payne 2009) to seismic data. Zhao et al. (2013) proposed to apply a neural network inversion method to only estimate the porosity and then estimate pore type in carbonate heavy oil reservoirs from seismic data. However, seismic elastic properties such as P-wave velocity are a function of porosity, saturation and pore type. That is, these reservoir parameters together have a comprehensive impact on the elastic parameters; therefore, one of them such as porosity cannot be recovered alone without concerning the others. Adelinet and Le Ravalec (2015) applied the effective medium theory to make an inference about the porosity and fracture from seismic wave impedance. Their method assumes that carbonate pore space contains two types: One is equant pore, and the other is crack, and these two porosity inclusions are filled with water, and then takes advantage of the available impedance to quantify the spatial variations in the equant porosity and fracture density. Therefore, there are some limitations in the current methods. For example, some methods assumed that pores are filled with water or heavy oil. Some methods simplified the pore system. At present, few scholars have proposed a solution to simultaneously invert the porosity, saturation and pore-structure parameters by means of the Xu-Payne model.
In this paper, we present a two-step method to estimate porosity and saturation as well as pore-type distributions from seismic data by using only seismic rock-physics models. Firstly, we summarize our rock-physics modeling method for complex carbonate reservoirs. Then, we design a workflow on how to give the equivalent way dealing with the pore space of carbonate reservoirs and estimate the reservoir parameters. By treating the pore space of carbonate reservoirs as a single-porosity system, we simultaneously estimate the porosity, fluid saturation and effective pore aspect ratio; after that, by treating the pore space of carbonate reservoirs as a triple-porosity system, we apply Xu-Payne model to estimate the porosities of three pore types from seismic data; furthermore, by numerical experiments on theoretical models we prove that the elastic properties of porous rock can be equivalently simulated by a singleporosity system with an effective pore aspect ratio. Finally, the results of estimating pore type of a dolomite reservoir from real data are given to demonstrate the effectiveness of the proposed method.

Xu-Payne Model
Xu-Payne model (Xu and Payne 2009) has been widely used in rock-physics modeling and inversion on carbonate rocks (Zhao et al. 2013;Neto et al. 2014). It partitions the pore space of carbonate rock into clay and non-clay pores. The non-clay pores are divided into three types of elliptical pores with different aspect ratios. These three representative pore types are defined as follows: (a) Interparticle and intercrystal pores are defined as reference pores; they have moderate pore aspect ratios which are between 0.12-0.15; (b) microfractures and micro-cracks are defined as compliant cracks; they have lower aspect ratios which are usually below 0.02; (c) moldic pores and vugs are defined as stiff pores; they have higher aspect ratios which are usually over 0.7. The porosity of clay pores that are regarded as micropores with bound water is the product of the total porosity and the shale volume. Xu and Payne (2009) proposed a rock-physics modeling workflow for their model. There are two main aspects. One is to estimate the elastic moduli of dry rock frame using the effective medium model. The other is to calculate the velocities of low-frequency, water-saturated rock using Gassmann's formulae (Mavko et al. 2009). When clay micropores exist, they are added to the matrix using the differential effective medium as a part of "solid" (Xu and Payne 2009).

Dry frame moduli for porous rock
The elastic properties of carbonate rocks have been simulated using the effective medium models, such as Kuster-Toksöz theory and DEM theory, with different concentrations of pores and specified pore shapes or pore types (Kumar and Han 2005;Sayers 2008;Xu and Payne 2009;Zhao et al. 2013;Li and Zhang 2018). The classical DEM process and the Kuster-Toksöz theory applied in the original Xu-Payne model are used to calculate effective elastic moduli for dry rocks (Kuster and Toksöz 1974;Xu and Payne 2009). Compared with the Kuster-Toksöz theory, the DEM theory simulates the effective elastic properties of the composite more accurately and efficiently even if the volume percentages of the inclusion are large (Norris 1985;Berryman and Berge 1996;Li and Zhang 2014). At present, the DEM theory has been extended from the single-porosity medium to the multiple-porosity medium without the problem of the added order of the inclusions (Li and Zhang 2014). The ordinary differential equations for dry porous rocks can be written as (Li and Zhang, 2014): where y is porosity, K * (y) and G * (y) are, respectively, the effective bulk and shear moduli of the composite, i is the volume percentage of the ith inclusion in total porosity , P * i and Q * i are called as the polarization factors of the effective bulk and shear moduli K * and G * , respectively; they depend on the aspect ratios of the inclusions (Mavko et al. 2009). N is the number of the inclusions. When N = 1 , they are the same as the classical DEM theory for dry rocks (Mavko et al. 2009;Li and Zhang 2014). When N = 3 , they become the triple-porosity model which is suitable to predict the effective elastic properties for Xu-Payne model.
The effective elastic properties for dry rocks might be obtained not only from the numerical solutions of differential Eqs. 1 and 2 by integrating them from porosity y = 0 with values K * (0) = K m and G * (0) = G m to the desired highest value y = , but also from their analytical approximation formulae (Li and Zhang 2014).

Rock-physics inversion in carbonate reservoirs
When we use Xu-Payne model to estimate the pore type of carbonate reservoirs, we must firstly define the reference line (Xu and Payne 2009) which is calculated by applying rockphysics model and the known porosity and saturation under the assumption that the whole pore spaces of porous rock contain only reference pores. It means that we have to provide the porosity and saturation for Xu-Payne model when it is applied to seismic data. Since the complex pore system of the carbonate reservoir makes the relationship between the velocity and porosity widely scattered (Rafavich et al. 1984;Eberli et al. 2003;Sayers 2008;Baechle et al. 2009), it would lead to a large error in the estimation of porosity and saturation if the effect of pore structure was not taken into account in the adopted rock-physics model (Li and Zhang 2018). Li and Zhang (2018) proposed an inversion scheme for porosity and saturation to solve this problem based on three dimensions rock-physics template (3D RPT). Here, we draw on their research ideas and present a two-step workflow that only applies rock-physics model to estimate pore-type distribution of carbonate reservoirs from seismic data ( Fig. 1). In step one, the pore space of carbonate reservoirs is treated as a single-porosity system with the effective aspect ratio of pore; then a 3D RPT on P-and S-wave impedances and density with the effective aspect ratio of pore is established through the Gassmann's equations and effective medium models, porosity and fluid saturation; then, the effective aspect ratio of pore, porosity and fluid saturation are simultaneously estimated from the seismic data based on 3D rockphysics template (Li and Zhang 2018). In step two, the pore space of carbonate reservoirs is treated as a triple-porosity system just as the Xu-Payne model defines; then combined with the inverted porosity and saturation obtained in the first step, three types of pores which include reference pores, stiff pores and cracks can be inverted by the means of Xu-Payne method. The seismic elastic data in Fig. 1 generally include P-and S-wave velocities and density, which are inverted from pre-stack gather data based on the amplitude variation with offset (AVO) theory.

The effectiveness on porous rock taken as a single-porosity system
Assuming that the weighted summation of the polarization factors P * i and Q * i of inclusion i with the pore aspect ratio i in Eqs. 1 and 2 can be approximated by the average polarization factors P * I and Q * I with an effective pore aspect ratio I , that is, I , the differential Eqs. 1 and 2 become:  Fig. 1 Diagram workflow of a two-step method to estimate porosity and saturation as well as pore-type distribution from seismic data for porous rock 1 3 The above equations indicate that the elastic properties for porous rocks can be simulated by DEM equations for a single-porosity system with an effective pore aspect ratio I . We give a comparison of simulation results between Eqs. 1-4. Here, we imitate the Xu-Payne model and take dolomite as the host medium, let K m = 89 GPa and G m = 37 GPa. We assume that there are three pore types in the rock, which are the stiff pores with the pore aspect ratio of 0.9, the reference pores with the pore aspect ratio of 0.15 and the micro-cracks with the pore aspect ratio of 0.005, respectively. We also assume that three different materials exist within the pore spaces, which are oil, water and dry void, respectively. Then, the fourth-order Runge-Kutta scheme can be used to get numerical solutions of the differential Eqs. 1-4. The step size here used in numerical integration is set to Δy = 0.1 . In Fig. 2, the scattered points stand for the simulated elastic moduli using Eqs. 1 and 2 for different pore volumes. The red points represent the simulated values for the volume contents of the stiff pores, the reference pores and micro-cracks, which are 89.9%, 10% and 0.1%, respectively, indicating that the stiff pores in the pore space of porous rock are dominant, while the blue points represent the simulated values for the volume contents, which are 10%, 89.9% and 0.1%, respectively, indicating that the reference pores in porous rock are dominant; The pink points represent the simulated values for the volume contents, which are 45%, 50% and 5%, respectively, indicating that micro-cracks in porous rock pores are dominant. The curves represent (c) Triple porosity system for (10%, 89.9%, 0.1%) Triple porosity system for (45%, 50%, 5%) Single porosity system for α = 0.420 Single porosity system for α = 0.155 Single porosity system for α = 0.070 Triple porosity system for (89.9%, 10%, 0.1%) Triple porosity system for (10%, 89.9%, 0.1%) Triple porosity system for (45%, 50%, 5%) Single porosity system for α = 0.420 Single porosity system for α = 0.155 Single porosity system for α = 0.070 Triple porosity system for (89.9%, 10%, 0.1%) Triple porosity system for (10%, 89.9%, 0.1%) Triple porosity system for (45%, 50%, 5%) Single porosity system for α = 0.420 Single porosity system for α = 0.155 Single porosity system for α = 0.070 Triple porosity system for (89.9%, 10%, 0.1%) Triple porosity system for (10%, 89.9%, 0.1%) Triple porosity system for (45%, 50%, 5%) Single porosity system for α = 0.420 Single porosity system for α = 0.155 Single porosity system for α = 0.070 the simulated elastic moduli using Eqs. 3 and 4 for different effective pore aspect ratios. The red curve stands for the simulated values for the effective pore aspect ratio eff = 0.42 ; the blue curve stands for the simulated values for eff = 0.155 ; the pink curve stands for the simulated values for eff = 0.07 . Note that these curves for both dry rocks and the water-and oil-saturated rocks match well with the scattered points. It is fully illustrated that the elastic properties of multiple-porosity medium can be simulated by an equivalent single-porosity medium with an effective aspect ratio.

Porosity, saturation and effective pore aspect ratio inversion based on a single-porosity system
The simultaneous inversion of porosity and saturation as well as effective pore aspect ratio can be carried out by using the 3D rock-physics templates (Li and Zhang 2018). Li and Zhang (2018) give concrete steps to realize the inversion: First determine the distribution ranges of the porosity, saturation and pore aspect ratio from logging data in the study area, and generate the 3D mesh in the reservoir parameters space, which is made up of the porosity, saturation and pore aspect ratio; ensure that all the grid nodes are assigned to a specific saturation, porosity and pore aspect ratio; then, employ the DEM Eqs. 3-4 and fluid distribution model for each grid node to compute its bulk and shear moduli of dry rock and fluid moduli under the given pore aspect ratio, porosity and saturation, and then calculate the P-and S-wave velocities and density of saturated rocks by means of the Gassmann's equations, until all the nodes in the grid have been calculated; after that, build a 3D RPT in the elastic parameters space such as Ip-Is-ρ or λρ-μρ-ρ; superpose the elastic parameters measured by logging tools onto the built 3D RPT to verify the validity of the 3D RPT; finally, project these elastic parameters inverted from seismic gather onto the newly built 3D RPT, and use the least square method to search the closest grid point for the inverted elastic data; choose the pore aspect ratio and porosity as well as saturation corresponding to the closest mesh point as the inversion results of reservoir results.

Pore-type inversion based on a multi-porosity system
Once the porosity and saturation are determined, the volume percentages of three pore types (reference pores, stiff pores and cracks) can be estimated from the seismicinverted elastic properties by following the inversion scheme of Xu and Payne (Xu and Payne 2009;Zhao et al. 2013).
The key inversion steps include: Firstly, define the aspect ratios of the three pore types that are reference pores, stiff pores and cracks in the carbonate reservoir. Then, assume that the whole pore spaces of carbonate rock contain reference pores alone and apply the DEM analytical models for single-porosity system and the Gassmann's equations to compute P-wave velocity V reference p at given porosity, saturation and aspect ratio. Next, determine whether the pore type of reservoir belongs to crack pore type or cavity pore type by comparing the measured P-wave velocity Vp with the estimated P-wave velocity V reference p . Finally, apply nonlinear global optimization algorithm to find the best estimate for the pore volume percentages of both cracks and reference pores or both cavities and reference pores by minimizing the error between theoretical predictions and geophysical measurements based on the dual-porosity system.

Survey of the study area
Our workflow has been applied in porosity, saturation and pore-type prediction of Cambrian Longwangmiao formation dolomite reservoirs in Anyue gas field in southwestern China's Sichuan Basin. Both matrix porosity and permeability of this dolomite reservoir are low. The porosity and permeability of the core are 2-15% and 0.001-10 mD, respectively. The reservoir porosity has strong heterogeneity. A large number of cavities and micro-cracks are developed in the reservoir. Pore types in the study area are classified into three types: intergranular and intragranular pores, dissolved cavities and fractures (Zhou et al. 2015;Li and Zhang 2018). Figure 3 Figure 4 shows the analysis of 3D rock-physics template on well X12. The parameters and display modes used in the rock-physics modeling are same as what Li and Zhang (2018) did. They took the bulk and shear moduli of rock matrix as K m = 89 GPa, G m = 37 GPa, respectively, and the density of rock matrix is taken as ρ m = 2.8 g/cm 3 , and the bulk modulus of water is taken as K w = 2.2 GPa. The blue surface represents the trend of P-and S-wave impedances and density with the pore aspect ratio and porosity for water saturation S w = 100%, and the brown surface represents the trend of P-and S-wave impedances and density with the pore aspect ratio and porosity for water saturation S w = 0%. The yellow curve plotted on each surface denotes the one whose porosities are constant, and the light blue curve denotes the one whose pore aspect ratios are constant. Light green surface indicates the trend of P-and S-wave impedances and density with porosity and saturation for the effective pore aspect ratio α = 0.1. The blue surface for fully watersaturated rock is always situated above the brown for fully gas-saturated rock. This may due to the fact that the density of the water-bearing rock is always greater than that of the gas-bearing rock. The points on the cross-plot are color coded with water saturation in Fig. 4a and porosity in Fig. 4b, respectively (Li and Zhang 2018). Note that most of the samples whose water saturations that interpreted from logging data are below 20% situate near the brown surface, and are drawn in red colors. Only very few samples whose water saturations that inferred from well-logging data are more than 80% situate near the blue surface, and are drawn in blue colors. The reservoirs with high gas saturation have high porosities. More 3D RPT analysis results for other two wells which come from the same the study area are found in the article of Li and Zhang (2018). These results show that the predictions from the 3D RPTs are fully coincident with the log interpretation and real drilling results. Hence, they can be effectively used to characterize the P-and S-wave impedances and density with porosity and water saturation and pore aspect ratio (Li and Zhang 2018). Figure 5 shows the inversion results of the pore aspect ratio, porosity and saturation from logging data based on the 3D RPT in Fig. 4. In Fig. 5, the black curves indicate the measured P-and S-wave velocities, density, and the inferred porosity and saturation from logging data, respectively. The red curves represent, respectively, the inverted effective pore aspect ratio, porosity and water saturation. We also display the lithology derived from the field by the petrophysicist in the figure. It shows that most of the reservoir is dolomite, with only a small amount of clay whose volume content is almost constant and a very small amount of quartz minerals at its top and bottom. Observe that the porosity of well X12 is high, the maximum porosity is up to 15%, the gasbearing property is good, and the water saturation is generally less than 20%. Both the inverted porosities and water saturation based on the 3D RPT are in well coincidence with the inferred ones from logging data, and also coincide with those from the actual drilling. The inverted effective pore aspect ratios for high porosity reservoir are around 0.2. Analysis of 3D rock-physics template on logging data for dolomite reservoir in gas well X12 with data points is coded with water saturation (a) and porosity (b), respectively Figure 6 shows the pore-type prediction results for gas well X12 following the solution of Xu and Payne (Xu and Payne 2009). The blue colors represent the inferred stiff pores, the green colors represent the inferred reference pores, and the red colors represent the inferred cracks. It can be seen that there are remarkable variations in both pore types and the corresponding effective pore aspect ratios α eff , which are shown in Fig. 4 Fig. 6 Pore-type inversion results for well X12. Area of red, green and blue colors represent volume concentrations of crack-like, reference and stiff pores in the pore space, respectively reference porosity and stiff porosity indicate that there are high porosities at a depth of 4622.1 m. There are several cracks in the bottom thin slice, and the inverted crack porosity indicates that there are high porosities at a depth of 4666.7 m. Figure 7a-c shows the prediction results of the porosity, saturation and effective pore aspect ratio from the inverted seismic elastic data through wells X12, X203 and X8 based on the 3D RPT for a single-porosity system. In addition to well X12 mentioned above, well X8 is also an industrial gas well, while well X203 is a water-bearing well whose reservoir has low porosity (Li and Zhang 2018). Figure 7a-c shows the inferred porosity, saturation and pore aspect ratio profiles, respectively. It can be found that both industrial gas wells X8 and X12 situate at locally high points and their porosities and gas saturations are relatively high, while water well X203 situates at the secondary high point and its porosities and gas saturations are relatively low. It is obvious that the pore aspect ratios of industrial gas well X8 are large, which indicates the pores are relatively harder, while the pore aspect ratios of well X203 are lower, which means the pores are relatively softer. The pore aspect ratios of well X12 are smaller than those of well X8. Figure 8a-c demonstrates the pore-type predictions from the inverted elastic properties combining with the inverted porosity and saturation displayed in Fig. 7a, b. One can observe that the dolomite reservoirs are extremely strong and heterogeneous. The inferred reference porosity distribution as shown in Fig. 8a indicates that the interparticle pores are the dominant pore type in this dolomite reservoir. The inferred stiff porosity distribution as shown in Fig. 8c indicates that the vuggy and moldic pores are mainly developed around wells X8 and X12. The porosities of reference pores and stiff pores of industrial gas well X12 are relatively larger than those of well X8, whereas the inferred compliant porosity distribution as shown in Fig. 8b indicates that cracks are not well developed in the whole cross-well profile, only developed around well X203.

Discussion
It is important to note that the proposed method to infer the porosity and saturation as well as volume percentage or porosity of each pore type in pore system from seismic   Fig. 7 Sections of the predicted porosity (a), water saturation (b) and effective pore aspect ratio (c) from seismic data for dolomite reservoir data is achieved by employing the Gassmann's equations and the differential effective medium theory of multiple-porosity system. It is generally known that for a porous system, Gassmann's equations have the following basic assumptions (Mavko et al. 2009;Adam et al. 2006): (a) The rock is macroscopically homogeneous and isotropic porous medium; (b) the pore space is connected; (c) a frictionless fluid is filled in the pores; (d) the rock-fluid system is closed; (e) when the rock is excited by a wave, the relative motion between the solid rock and the fluid is negligibly small compared to the motion of the whole saturated rock itself; (f) the pore fluid does not interact with the solid in a way that would soften or harden the frame. Therefore, when using Gassmann's equations, the above assumptions must be taken into account. The assumption (a) on macroscopic homogeneity and isotropy of rock ensures that the wavelength is long compared to the grain and pore sizes. Most rocks can generally satisfy this assumption for seismic frequency band . On the other hand, the DEM theory assumes that the pores in porous medium disconnected each other and fluids cannot move through pores; this approach simulates very high-frequency saturated rock behavior suitable for ultrasonic laboratory condition (Mavko et al. 2009). However, some researchers have testified that there are almost no differences between elastic properties of rocks compared to their measurements at high and low frequencies in laboratory experiments (Mavko et al. 2009;Adam et al. 2006). For low-frequency seismic condition, a practical solution is to calculate first the dry rock bulk and shear moduli by applying effective medium model, then use the Gassmann's equations to calculate the bulk and shear moduli of saturated rocks (Mavko et al. 2009;Kumar and Han 2005;Xu and Payne 2009;Zhao et al. 2013;Adelinet and Ravalec 2015). Carbonate rocks have strong elastic frames and very complex pore systems compared to siliciclastic rocks. Because a significant number of pores in carbonate rocks is not well connected, carbonate rocks not always meet Gassmann's basic assumption on pore connectivity. The applicability of Gassmann's theory in carbonate rocks has been targetedly investigated, and an inspiring conclusion from those researches is that Gassmann's equations are appropriate to simulate the saturated bulk and shear moduli of carbonate rocks (Sharma et al. 2006;Sharma et al. 2011;Grechka 2009;Gegenhuber and Pupos 2015;Yan and Han 2016). By using numerical tests, Grechka (2009) . 8 Sections of pore-type distribution inverted from seismic data at the target area: a reference porosity, b micro-cracks porosity and c stiff porosity 1 3 disconnected porosity is unlikely to cause significant errors in infill substitution when the pore aspect ratios are greater than about 0.2. If the pores are aligned or the bulk and shear moduli of the substituted infills are similar, the limits of practical applicability of Gassmann's predictions extend to smaller aspect ratios. Therefore, theoretically, Gassmann's assumption (c) on pore connectivity is not a necessary condition (Grechka 2009;Yan and Han 2016). Our method to invert for pore types of reservoir from seismic data is based on P-and S-wave impedances and density inverted from pre-stack seismic data. Therefore, the quality of pore-type estimation depends strongly on the accuracy of pre-stack seismic inversion which is a relatively mature technique and has been widely applied in seismic reservoir prediction (Fatti et al 1994;Buland and More 2003;Hampson et al. 2005;Wang et al. 2019). It is generally believed that P-and S-wave data can be unambiguously resolved by pre-stack seismic inversion, but density, which is crucial for hydrocarbon detection, cannot (Tarantola 1986;Debski and Tarantola 1995). This is mainly because the amplitude is not sensitive enough to density, the offset limited data tend to result in unstable results of the three-parameter AVO inversion, and small noise may lead to a huge inversion error (Tarantola 1986;Debski and Tarantola 1995). In order to ensure the reliability of density inversion, a large number of studies have been carried out on the correlation between P-wave impedance and both S-wave impedance and density. Hampson et al. 2005 proposed a new method to invert P-and S-wave velocities and density from pre-stack seismic data simultaneously by assuming the linear relationship between the logarithm of P-wave impedance and both S-wave impedance and density. Their method is effective in both model data and actual seismic data Swan (2007) proposed a background trend correction method based on the correlation coefficient between P-and S-wave impedance. The P-wave reflectivity obtained by AVO inversion is used to correct the S-wave reflectivity obtained together, so that the correlation between them can match with the correlation of theoretical models or log data. Following their ideas, we correlate the inverted S-wave impedance and density by using the P-wave impedance obtained simultaneously to make sure that the correlation coefficients between the inverted P-wave impedance with both S-wave impedance and density are matched with the corresponding correlation coefficients in well log data. Figure 9 shows the result of pre-stack threeterm simultaneous inversion. The P-and S-wave impedances  . 9 Sections of the inverted P-wave impedance (a), S-wave impedance (b) and density (c) which are used to infer the reservoir parameters and densities of the two gas-bearing wells in the target interval (between the horizons E1g and E1l) have relatively low values, while the corresponding values of the water wells are relatively high. Figure 10 shows the cross-plots between P-wave impedance and both S-wave impedance and density for logging data of the three well in target strata and seismic inversion data adjacent to the borehole. The correlation coefficients among the inverted seismic data are close to those among the well-logging data. Our inversion method ensures the accuracy of the inverted reservoir parameters through correlation coefficient matching correction.

Conclusions
In this paper, we proposed a new two-step method to estimate the porosity and saturation and quantitatively characterize spatial distribution of pore type from seismic data on carbonate reservoirs by using only seismic rock-physics models. In step one, the pore space of carbonate reservoirs was equivalent to a single-porosity system with an effective pore aspect ratio; then, a 3D RPT on P-and S-wave impedances and density with the effective aspect ratio of pore, porosity and fluid saturation was established through the  Fig. 10 Cross-plots between P-wave impedance and S-wave impedance for log data (a), seismic inversion data (b) adjacent to the borehole, and between P-wave impedance and density for log data (c), seismic inversion data (d) adjacent to the borehole Gassmann's equations and effective medium models; and then, the effective aspect ratio of pore, porosity and fluid saturation were simultaneously estimated from the seismic data based on 3D RPT. In step two, the pore space of carbonate reservoirs was equivalent to a triple-porosity system just as what the Xu-Payne model defines. Combined with the inverted porosity and saturation in the first step, the porosities of three pore types (reference pores, stiff pores and cracks) were inverted from the seismic elastic properties following the Xu-Payne's scheme.
The effectiveness on porous rock taken as a single-porosity system was testified by an effective equation derived from the differential effective medium theory and numerical experiments on theoretical models. The predicted elastic moduli for both dry rocks and the water-and oil-saturated rocks by the differential effective medium (DEM) model of multiple-porosity rock are well matched with the simulated results by the DEM model of single-porosity rock.
The applications to well-logging data and seismic data both indicate that the proposed inversion scheme can obtain more accurate physical properties and ensure the reliability of characterization of pore-type distribution. Meanwhile, the estimated pore type shows that the dolomite reservoir has strongly heterogeneous pore space in the target study area.