TY - JOUR AU - Yu,, Gang AB - Abstract To achieve the economic and effective development of shale gas, it is not only necessary to determine favourable geological areas but also to pay more attention to favourable engineering areas. Geologically favourable areas are measured on their total organic carbon (TOC) content, thickness of high-quality shale and fractures, and favourable engineering areas are assessed on brittleness and fractures, but in most cases fractures occur in both geological and engineering favourable areas. The Zhao104 well is located in the Sichuan Basin in China, in the Longmaxi (LMX) Formation, and has an abundant shale gas content based on drilling and well log interpretation. It is necessary to clarify the distribution of shale gas favourable areas in this area, which plays a vital role in determining horizontal wells. In this paper, petrophysics and three-dimensional seismic attributes are used to predict the favourable shale reservoir area. Based on the petro-physical model, the most sensitive elastic parameters of the TOC content and brittleness are determined, a relationship model of the sensitive elastic parameters with TOC content and brittleness is created from logs and pre-stack inversion is used to characterize the shale brittleness, TOC content and thickness of high-quality shale. Different scale seismic fracture attributes are used to predict the strength and direction of fractures. The field data show that combined petrogeophysics and 3D seismic attributes can effectively predict a favourable shale reservoir area, reducing the need for exploration and development risks. petro-physical modelling, pre-stack inversion, P-wave anisotropy, 3D seismic attributes, gas reservoir 1. Introduction Shale pertains to a reservoir with an ultralow porosity and permeability; its occurrence mode, accumulation pattern of natural gas and development model are remarkably different from those of conventional oil and gas reservoirs. How to identify favourable shale reservoir areas has been a very important task in the exploration and development of shale gas reservoirs in North America (Sena et al.2011; Jaiswal et al.2014; Heath-Clarke et al.2016; Convers et al.2017). Rickman et al. (2008) defined the rock brittleness index with Young's modulus E and Poisson's ratio σ, and calculated the brittleness index. It has been shown by their research that the brittleness will be higher if the Young's modulus is larger and Poisson's ratio is lower. Other scholars have also defined rock brittleness from mineral components (quartz, calcite, clay content etc.) and combined elastic parameters (Young's modulus × density (Eρ), Lame constant × density (λρ), shear modulus × density (μρ) etc.), and different brittleness prediction methods have been compared, thus showing that the different methods have their own advantages and should be selected according to the geological conditions of the actual field (Rybacki et al.2016; Liu et al.2017, 2019b; Kivi et al.2017; Verma et al.2016). Sena et al. (2011) have proposed one method that integrates pre-stack seismic inversion and logging to identify favourable areas, while geo-mechanical properties and stress fields can be obtained from seismic data. Tavella et al. (2013) have provided a relationship between the elastic parameters and brittleness of a well, and we can extrapolate this relationship to predict TOC and brittleness distribution characteristics from a well to a three-dimensional area. A favourable area identification workflow for shale gas has been established in North America (Chopra et al.2012; Carpenter 2018). Pore pressure research is also beginning to pay attention to shale gas reservoirs (Banik et al.2014; Dasgupta et al.2016). In comparison with foreign countries, shale gas exploration and development in China have achieved certain positive effects in recent years. Several shale gas demonstration zones have been established, and the North American shale gas favourable area workflow has been applied to the Sichuan Basin in China achieving suitable results (Wang et al.2015,; Yu et al.2017; Zeng et al.2017, 2018; Chen et al.2018; Zhao et al.2018; Liu et al.2019a). At the same time, Chinese scholars have made certain improvements according to the previous North American favourable area workflow because shale gas reservoir geological conditions are very different between America and China. Yu et al. (2014) have paid more attention to shale gas rock physics modelling, including logging curve correlation and mineral disturbance analysis, which can directly affect the accuracy of favourable area prediction. Zhang et al. (2014) have introduced a new TOC evaluation technique, which has considered the influence of hydrodynamic conditions and sedimentation environments on TOC contents. This method has improved the accuracy of TOC prediction compared with a conventional TOC calculation method proposed by Passey et al. (1990). In this work, petro-physical modelling and three-dimensional seismic data are used to predict the shale reservoir fracture development, TOC, brittleness and high-quality shale thickness. First, the most sensitive elastic parameters of the TOC content and brittleness are determined based on a petro-physical model, which relies on the elastic parameter product of the Eρ. Then, we create a relationship model of the sensitive elastic parameters and the TOC content and brittleness from logs. Second, the relationship model is applied to the entire area by pre-stack elastic inversion and the Eρ attribute is used to characterize the shale brittleness, TOC content and thickness of high-quality shale. The seismic attributes and P-wave anisotropy are used to predict the strength and direction of fractures based on post-stack and pre-stack seismic data, respectively. Finally, we integrate the results of shale reservoir fracture development, TOC content, brittleness and high-quality shale thickness, and determine the favourable areas. 2. Geological overview The study area is the Longmaxi (LMX) Formation shale reservoir in the well area of Zhao104, Sichuan Basin, China. The well area of Zhao104 is located at the junction of the Sichuan Basin and the Yunnan-Guizhou Plateau, and the tectonic stress field and geological structure are complex. To determine the potential of the shale gas resources and obtain shale gas evaluation parameters in this area, parameter well Zhao104 was deployed in October 2010. Figure 1 shows a comprehensive evaluation map of the Zhao104 well. The gas content was suitable from the LMX Formation to the Wufeng (WF) Formation, and some cores were taken for later study including petro-physical analysis. The thickness of the high shale gas segment was approximately 40 m. In addition, the TOC content was relatively high, and the highest TOC content located at the bottom of the LMX Formation was nearly 5%. The gas logging data indicated that this interval was the best shale gas segment, which created a good start for shale gas exploration and development in this area. In 2011, 100 km2 of 3D seismic acquisitions were implemented in the Zhao104 area, and 3D seismic data were used to predict the favourable shale areas. The purpose of this study was to clarify the distribution of shale gas favourable areas in the area, thereby providing geophysical evidence for the next exploration and development stages. Figure 1. Open in new tabDownload slide Comprehensive logging evaluation map of the Zhao104 well. Figure 1. Open in new tabDownload slide Comprehensive logging evaluation map of the Zhao104 well. 3. Analysis of the exploration difficulties The target layer is the Silurian LMX Formation. The study area has undergone multi-phase intense orogenic tectonic movement since the Indosinian period, and the area is in a steep mountain region; the fault development and shale gas exploration and development levels are very low, so the research area of favourable areas is confronted with certain challenges. On the one hand, how to accurately and effectively identify fractures is a challenge due to the multi-phase tectonic movement and complex geological conditions, which have caused fractures to develop. On the other hand, the depth of the target layer varies greatly, from 500 to 5000 m, and the high-quality shale layer is approximately 30 m. In addition, the exploration and development levels of the study area are low. There is only one parameter well, Zhao104, and one horizontal well, YSH1–1, in the early stage of exploration and development. There is little information, which increases the multiple solutions to the favourable area prediction results. Before 3D seismic exploration, vertical well Zhao104 and horizontal well YSH1–1 were deployed in the area, which were used to obtain evaluation parameters of the study area. Hydraulic fracturing was implemented for horizontal well YSH1–1 during shale gas development, but the shale gas production was not very high. To accelerate the progress of shale gas exploration and development in this area and to reduce the risk of horizontal well failure, 3D seismic exploration was carried out in the area; the acquisition observation system used a 30 m × 30 m bin, and the aspect ratio in the target layer reached 0.8, which was aimed at meeting the theory of horizontal transverse isotropy (HTI) fracture prediction; the latter would help to solve the difficulty of micro-fracture identification. Figure 2 shows the surface conditions of the study area, and figure 2a shows the elevation of the study area. The red and yellow characters represent high altitudes, the green and grey characters represent low altitudes and the elevation difference of the whole area is very large. Figure 2b shows the detailed elevation difference of one line in the study area. The maximum elevation difference was nearly 700 m. Figure 3 shows the logging response of high-quality shale thickness characteristics for well Zhao104. The thickness of the entire LMX Formation was approximately 200 m, the lithology was mainly shale and the underlying layer was the WF Formation. The lithology was mainly shale with a high organic carbon content, which was between 35 and 40 m; the layer below the WF Formation was the Baota (BT) Formation and its lithology was limestone. High-quality shale was characterized by three highs and three lows; namely, high acoustic time difference, high gamma ray, high TOC content, low density, low neutron density and low resistivity. Figure 2. Open in new tabDownload slide Surface conditions: (a) the elevation of the study area and (b) the elevation of one inline in the study area. Figure 2. Open in new tabDownload slide Surface conditions: (a) the elevation of the study area and (b) the elevation of one inline in the study area. Figure 3. Open in new tabDownload slide Logging responses of the high-quality shale thickness characteristics for the Zhao104 well. Figure 3. Open in new tabDownload slide Logging responses of the high-quality shale thickness characteristics for the Zhao104 well. 4. Method and technical strategy There is only one well in this area; we have done a lot of work to ensure the reliability of subsequent predictions, which is different from the conventional methods. We get more logging parameters for this well as far as possible, and carry out logging curve corrections to ensure the reliability of logging input information and to improve the reliability of prediction. In the single well prediction, we first use element capture logging to get the mineral component, and then we use the actual TOC, quartz core data to constrain or calibrate the TOC or quartz calculation, which can increase the reliability of that brittleness or TOC we get from a single well. We use quartz contents to represent shale brittleness. Elemental capture logging and core data is used to constrain TOC and quartz calculation after the logging correction. Compared with the conventional TOC and brittleness calculation method, TOC calculation is mainly based on Passey's et al. (1990) method or its modification, the TOC method key is the reliability of acoustic time logging, resistivity logging. The brittleness calculation is mainly based on Young's modulus and Poisson's ratio (Rickman et al.2008). For the method of brittleness calculation, its key is to reliability of the P-wave velocity, S-wave velocity and density curve. Because the brittleness calculation is from the Young's modulus and Poisson's ratio, and the Young's modulus and Poisson's ratio are obtained by P-wave velocity, S-wave velocity and density curve. If the initial p-wave velocity, s-wave velocity and density curve are not corrected, the accuracy of the subsequent brittleness prediction is unreliable. If curve correction is not done, reasonable results cannot be guaranteed. Besides, conventional methods cannot directly calculate TOC and brittleness mineral content, they only use a general logging curve to calculation TOC and brittleness because many wells do not carry out element capture logging and core data. From a single well prediction, logging curves and rock physics modelling are some necessary strategies to ensure the reliability of results. There are many factors affecting the quality of logging curves, such as wellbore collapse, logging time differences and missing well curves. At present, intersection analysis of the logging curve is an effective solution to the wellbore collapse problem and missing log curves. The wellbore collapse correction is mainly to analyse the relationship between the measured value and the theoretical value, and then optimize the theoretical model to correct the outliers into the theoretical model. In this paper, the purpose of curve correction is to eliminate the impact of wellbore collapse; intersection analysis is used to carry out logging curve correction analysis in three steps: Step 1: Establish different theoretical sandstone analysis templates. Create four kinds of models, (i) 100% quartz by the Raymer Model, (ii) 80% quartz plus 20% clay based on the Unconsolidated Model, (iii) 100% clay based on the Reuss Model and (iv) 100% clay based on the Unconsolidated Model. Step 2: Making an intersection analysis between the original density curve and original P-wave velocity. We can select the target depth P-velocity and density curves to analyse. Step 3: Identify the abnormal curve segment and correct it. The actual logging sampling point that falls down within the theoretical line is the corrective sampling point, and the outside of the theoretical line is the abnormal point. The theoretical model that is most suitable for the research area is preferred, and the abnormal point is corrected based on the preferred theoretical model. There are two types of fracture prediction methods; one is based on post-stack seismic attributes, i.e. coherence curvatures etc.; and the other one is based on pre-stack P-wave azimuthal anisotropy, i.e. amplitude variation with azimuth, frequency variation with azimuth etc. The post-stack seismic data prediction method is effective for large-scale fracture (greater than quarter wavelength) prediction; however, post-stack seismic data is not able to obtain small-scale fracture information (less than or equal to quarter wavelength) due to seismic resolution limitations, which means it cannot be imaged on less than a quarter wavelength scale. Therefore, for small-scale fracture prediction, P-wave azimuthal anisotropy is used to identify small-scale fractures based on HTI theory, which assumes that the anisotropy in the rock is caused by a set of directional vertical fractures; the fast P-wave is parallel to fractures, the slow P-wave is perpendicular to fractures and $$\begin{equation}V(\beta ) = {V_0} + \alpha \cdot \cos 2\beta, \end{equation}$$ (1) where |$V(\beta )$| is the P-wave velocity in different azimuth, |${V_0}$| is the average value of the azimuth velocity, |$\alpha $| is the velocity difference between the azimuth extreme velocity and the azimuth average velocity and |$\beta $| is angle between the shot azimuth and fractures strike. The elliptic velocity equation (1) is determined by three parameters, |${V_0}$|⁠, |$\alpha $| and |$\beta $|⁠. It is only necessary to know the velocity values of the three azimuths, then we can solve the equation and get the parameters. The long axis of the velocity azimuth ellipse represents the direction of the fractures, while the ellipse oblations represent the magnitude of fractures, described as follows $$\begin{equation}\eta = ({\sigma _1} + {\sigma _2})/({\sigma _1} - {\sigma _2}),\end{equation}$$ (2) where |$\eta $| indicates the fracture development magnitude, and |${\sigma _1}$| and |${\sigma _2}$| are the long axis and short axis of the ellipse, respectively. With regard to the brittleness and TOC, we can carry out elastic impedance inversion and construct sensitive elastic parameters based on petrophysical modelling. This article applies the equation of elastic impedance deduced on the basis of Connolly's work (1999): $$\begin{equation}EI\,(\theta ) = V_p^{(1 + {{\sin }^2}\theta )} \cdot V_s^{( - 8K{{\sin }^2}\theta )} \cdot {\rho ^{(1 - 4K{{\sin }^2}\theta )}},\end{equation}$$ (3) where |$EI(\theta )$| is the elastic impedance that will remarkably vary as the incidence angle changes; |${V_p}$| is the P-wave velocity; |${V_s}$| is the S-wave velocity, |$\rho $| is the medium density and |$\theta $| is the incident angle |$K = {( {{V_p}/{V_s}} )^2}$|⁠. It is unfavourable to compare elastic impedance values at different angles, so the elastic impedance is normalized to obtain the following: $$\begin{equation}EI\,(\theta ) = {\alpha _0}\,{\rho _0}\left[ {{{\left( {\frac{{{V_p}}}{{{\alpha _0}}}} \right)}^a} \cdot {{\left( {\frac{{{V_s}}}{{{\beta _0}}}} \right)}^b} \cdot {{\left( {\frac{\rho }{{{\rho _0}}}} \right)}^c}} \right],\end{equation}$$ (4) where $$\begin{eqnarray} \left\{\begin{array}{@{}l@{}} a = (1 + {{\sin }^2}\theta )\\ b = - 8K{{\sin }^2}\theta \quad\quad\quad,\\ c = (1 - 4K{{\sin }^2}\theta ) \end{array}\right. \end{eqnarray}$$ (5) and |${\alpha _0}$|⁠, |${\beta _0}$| and |${\rho _0}$| have been introduced in this equation as reference values, representing the P-wave velocity, S-wave velocity and density from well logging, which normalizes the impedance of the elastic waves to the acoustic impedance scale. Then, we can obtain Eρ by $$\begin{equation}E\rho = {I_s}^2\frac{{3{I_p}^2 - 4{I_s}^2}}{{{I_p}^2 - {I_s}^2}},\end{equation}$$ (6) where, |${I_p}$| and |${I_s}$| are P-impedance and S-impedance, respectively. Eρ can improve the results credibility. From a single well to a plane prediction, the conventional main method for obtaining TOC and brittleness is a fitting relationship between the TOC (or brittleness) and elastic parameters from a single well. The elastic parameters can be obtained by pre-stack inversion, so this kind of method can predict TOC and brittleness. If there is more well data, the fit relationship between the TOC (or brittleness) is better and more reliable. In this study area, as there is only one well in this area, the sample number is very small, so the conventional fit method does not get reliable results. Furthermore, if we use the fit method, we need to get elastic parameters (P-impedance, S-impedance, density, Young's Module etc.) by pre-stack inversion. If the fit relationship includes a density item (e.g. density, passion ration, Young's module), it will have an iterative error due to inversion density leading to poor reliability. Compared with conventional methods, the proposed method in this paper can directly characterize TOC and brittleness by using Eρ obtained from only P-impedance and S-impedance, then it can avoid iterative errors (in equation 3). So, we use the Eρ intersection analysis method to predict brittleness and TOC. In this study area, due to the difficulties of fracture development, small high-quality shale thickness and low accuracy of favourable area prediction, this paper proposes the following technical workflow solution. In the southern marine shale area, we used post-stack coherence and curvature attributes to identify large-scale fractures, and we used the pre-stack P-wave anisotropy to identify small-scale fractures. In this area, we combined two types of different-scale method to evaluate the fracture characteristics in this area, which could reduce the uncertainty of a single fracture prediction method. For high-quality shales and low accuracy favourable areas, we propose carrying out mineral disturbance analysis and petro-physical modelling to determine the sensitive elastic parameters of the TOC content and brittleness. We create a relationship model between the sensitive elastic parameters, TOC content and brittleness from logs. The relationship model is applied to the entire area by pre-stack inversion, and the sensitive elastic parameter results can be used to characterize the shale brittleness, TOC and thickness of high-quality shale, which could reduce iteration errors caused by conventional fitting due to the lack of wells, thus improving the accuracy of favourable shale gas area prediction. The main technical workflow of this paper is shown in figure 4. The main steps are as follows: Step 1: Logging curves and rock physics modelling: well logging curves are analysed and corrected, which is the basis of this study. Core data and logging curves are used to carry out mineral component inversion, and a shear wave prediction model is established for the study area. On this basis, different mineral models are established by changing the contents of quartz, TOC and clay, which are used to analyse the variation in elastic parameters corresponding to the different models and determine the sensitive elastic parameters of the TOC content and brittleness in the study area. Step 2: Fracture prediction: conventional post-stack attributes is used to identify large-scale fractures, and the P-wave anisotropy is used to predict small-scale fractures from the pre-stack gathers (Ruger 1998; Marfurt et al.1999; Gray & Head 2000). Different-scale fracture information is integrated to finely characterize the fracture characteristics of the shale reservoirs. Step 3: Pre-stack elastic parameter inversion: pre-stack elastic inversion is performed by well logging information that is constrained based on rock physics modelling, and the reservoir P-impedance, S-impedance and density values are obtained (Connolly 1999; Whitcombe 2002). The sensitive elastic parameters of the TOC content and brittleness are established in Step 1 by combining the P-impedance, S-impedance and density values, which is the foundation of the TOC content and brittleness predictions. Step 4: TOC/brittleness predictions: based on Step 1, intersection analyses are conducted between the TOC content or brittleness and sensitive elastic parameters to determine the sensitive elastic parameter ranges corresponding to high TOC or brittleness values. Then, the high TOC volume can be distinguished according to the sensitive elastic parameter value ranges. The brittleness prediction can also be implemented in this way. Step 5: Prediction of the high-quality shale thickness: using the high TOC volume determined in Step 4, we can make horizon tracks along the top and bottom of the high TOC volume, the TOC time thickness can be obtained and then the thickness of high-quality shale can be realized by time-depth conversion. Figure 4. Open in new tabDownload slide Shale gas favourable area prediction workflow. Figure 4. Open in new tabDownload slide Shale gas favourable area prediction workflow. Finally, we integrate the TOC content, brittleness, high-quality shale thickness and fracture characteristics from Steps 1 to 5 to comprehensively analyse the favourable shale areas, which are often characterized by a high TOC levels, high brittleness values, large thicknesses of high-quality shale and developed fractures. 5. Shale gas favourable area identification Before the favourable area prediction, a petro-physical analysis of the Zhao104 well was carried out in the study area, which included logging curve correction, shear wave prediction, mineral component calculation and mineral disturbance analysis. Logging curve correction was the basis of rock physics modelling, which ensured the reliability of the petrophysical analysis (Yu et al.2014). Figure 5 shows the intersection analysis of the P-wave velocities and density; green and red lines are computed for 100% quartz and 80% quartz plus 20% clay cases using the Raymer equation (Makar & Kamel 2011) and the Unconsolidated Model (Xu & White 1996), respectively. Blue and brown lines represent the result of Reuss equation and the Unconsolidated Model, respectively, calculated for the 100% clay case. Black points represent original data, while the coloured data is the final. Black points (figure 5) that fall beyond the model lines show poor measurement due to borehole conditions. One of the blue ovals shows the Wufeng formation (organic rich interval). The red theoretical model line that goes through the middle of the points representing Wufeng FR shows that this interval can be predicted using this model. Figure 5. Open in new tabDownload slide Logging curve intersection analysis of the P-wave velocities and density. Figure 5. Open in new tabDownload slide Logging curve intersection analysis of the P-wave velocities and density. Figure 6 compares the curves before and after correction. Projecting the sample points beyond the theoretical values in figure 5 onto the lithology section in figure 6, it was found that the points beyond the theoretical value range were located in the logging section corresponding to relatively large calliper changes (marked by the dotted rectangle in figure 6). Therefore, it was confirmed that the calliper changes led to the unreliability of the P-wave and S-wave well logs. After the calliper was corrected, the reliability of the basic data was guaranteed. Figure 6. Open in new tabDownload slide Curve comparison of the Zhao well logging curves before and after correction (red is the original logging curve and black is the corrected curve). Figure 6. Open in new tabDownload slide Curve comparison of the Zhao well logging curves before and after correction (red is the original logging curve and black is the corrected curve). S-wave prediction is very important for favourable shale gas areas, but we have carried out cross-dipole acoustic logging in this area; thus, we do not need to discuss the shear wave prediction. Mineral component calculation and disturbance analysis are key to shale gas petrophysical modelling (Yu et al.2014). The premise of mineral disturbance analysis is that the mineral components can be obtained. In the area of well 104, conventional well logging with the density, P-wave velocity, porosity and resistivity for the whole well was not only carried out, but also element logging and dipole transverse wave logging of the target layer was conducted. The clay content, TOC content, quartz content and calcite content of well 104 can be calculated through element capture logging. Figure 7 shows the prediction of the mineral content in well 104. In the mineral component prediction, a core was used for calibration (shown by the red dot in figure 7) to ensure the reliability of the predicted mineral composition. On this basis, by changing the contents of the TOC and brittle mineral components in the shale reservoirs, the influence of the mineral component content change on the reservoir elastic parameters was studied, and which elastic factors were the most sensitive to the TOC content and brittleness was determined. In shale reservoirs, the content of quartz is often much higher than those of other brittle minerals, so we chose quartz to represent the brittle minerals and used the quartz content to characterize brittleness. Two types of model were created at the target layer, including six models (excluding the original formation), as follows: Figure 7. Open in new tabDownload slide Mineral component calculations for Zhao104. Figure 7. Open in new tabDownload slide Mineral component calculations for Zhao104. Four TOC models: Model 1: The original formation, which has no changes to the mineral content. Model 2: The TOC content is decreased by 3%, and the clay content is increased by 3%. Model 3: The TOC content is increased by 6%, and the clay content is decreased by 6%. Model 4: The TOC content is increased by 9%, and the clay content is decreased by 9%. Four quartz models: Model 1: The original formation, which has no changes to the mineral content. Model 2: The quartz content is decreased by 5%, and the clay content is increased by 5%. Model 3: The quartz content is increased by 5%, and the clay content is decreased by 5%. Model 4: The quartz content is increased by 10%, and the clay content is decreased by 10%. Every model corresponds to a kind of elastic parameter; for the different models, the elastic parameters were different. The differences between the elastic parameters of different models were compared, and the difference factor was defined as the sensitivity (named the sensitivity factor, SF), which was the ratio of the variation value to the original value; the larger the ratio was, the more sensitive the elastic parameter was to the TOC content or brittleness. Therefore, we could determine the most sensitive elastic factor of the TOC content and brittleness by the different mineral model disturbance analyses. Figures 8 and 9 show the TOC and quartz perturbation analyses, respectively. The results showed that when the TOC or brittle mineral contents were changed, the corresponding elastic parameter changes were different. Table 1 shows the difference factor statistics of the mineral component disturbance analyses for TOC and quartz. The statistical analysis results showed that for the TOC content, λρ was the most sensitive to TOC and Eρ came second. For quartz, μρ was the most sensitive to the brittleness and Eρ came second. Thus, Eρ was very sensitive to both the TOC content and brittleness and Eρ was considered to be the most sensitive elastic factor to the TOC content and brittleness. In subsequent research, the relationship among the TOC content, brittleness and Eρ was established. The most sensitive elastic factor, that is the Eρ volume, was constructed by pre-stack inversion from the well to the 3D area, which was a basic factor for TOC and brittleness plane prediction. Figure 8. Open in new tabDownload slide TOC perturbation analysis. Figure 8. Open in new tabDownload slide TOC perturbation analysis. Figure 9. Open in new tabDownload slide Quartz perturbation analysis. Figure 9. Open in new tabDownload slide Quartz perturbation analysis. Table 1. Mineral component disturbance analysis difference factor statistics (Notes: SF = (variation value/original value), unit = %). Sensitivity Factor σ Eρ μρ λρ TOC 3.22 9.02 8.47 13.84 Quartz 2.06 4.00 4.91 0.23 Sensitivity Factor σ Eρ μρ λρ TOC 3.22 9.02 8.47 13.84 Quartz 2.06 4.00 4.91 0.23 Open in new tab Table 1. Mineral component disturbance analysis difference factor statistics (Notes: SF = (variation value/original value), unit = %). Sensitivity Factor σ Eρ μρ λρ TOC 3.22 9.02 8.47 13.84 Quartz 2.06 4.00 4.91 0.23 Sensitivity Factor σ Eρ μρ λρ TOC 3.22 9.02 8.47 13.84 Quartz 2.06 4.00 4.91 0.23 Open in new tab Figure 10 shows the seismic amplitude profile across the Zhao104 well. The bottom of the LMX Formation is characterized by a strong amplitude seismic response with continuous events, which means that the deposition environment is stable. Figure 11 shows the characteristics of the fracture distribution at the bottom of the LMX Formation from post-stack seismic attributes. Figure 11 parts a and b are the coherence attributes and minimum curvature attributes, respectively. Coherence and curvature attributes can describe the large-scale fractures of the target layer and the directions of the fractures are consistent with the principal structures, which are comprehensive responses for fracture systems. The post-stack coherence and curvature attributes indicate that there are many sets of large-scale fractures at the bottom of the LMX Formation, the northwest and east-west directions are the principal fracture directions; the northwest direction is auxiliary and the east-west directional fractures divide the study area into two parts, which is a result of several structures squeezing each other. The 104 and YSH1–1 wells are located in the mutual compression zone with several sets of structures, coherence attributes and curvature attributes, thus indicating that the fracture systems around wells 104 and YSH1–1 are undeveloped. Figure 12 is a coherence attribute section across well 104, where black indicates that the fractures are relatively developed, and white indicates that the fractures are not developed. It can be seen that the fractures are not developed at the bottom of the LMX Formation. Figure 13 shows the amplitude slices from four different azimuth data along the bottom of the LMX Formation; the amplitudes in figure 13a–d are 0–37°, 37–90°, 90–135° and 135–180°, respectively. There is a large amplitude difference between the different azimuths. This phenomenon is the basis for the P-wave anisotropy fracture prediction. In addition, the aspect ratio of the target layer is close to 0.8 in this area, which means that it is feasible to use P-wave anisotropy to predict the fractures. Figure 14 shows the strength and direction of the fracture obtained by P-wave anisotropy based on a pre-stack gather. The length and direction of the line represent the strength and direction of the fracture, respectively. The colour also reflects the strength of the fracture and the coherence attribute is the base map in figure 14. The P-wave anisotropic fracture prediction results show that the fracture system near Zhao104 and the horizontal well YSH1–1 is not developed, and the fracture direction of the Zhao104 well is mainly east-west. Compared with the post-stack coherence and curvature attributes, the principal fracture system distribution characteristics based on the pre-stack data agree with the coherence and curvature attributes. The difference is that the P-wave anisotropy method can reflect more details in addition to the large-scale fracture system, as shown by the red circle in figure 14. The red circle area indicates that the fractures are relatively developed, but in the coherence and curvature attribute slices the fracture systems are not very developed. Clearly, the P-wave anisotropy can characterize smaller-scale fractures. These differences are very important in later hydraulic fracturing. According to the prediction results of the pre-stack and post-stack fracture prediction methods, horizontal wells YSH1–1 and Zhao104 are located in the squeezed zone with several structures. In this area, these wells are subjected to multiple sets of tectonic stresses, but the fracture system near the wells is not developed. Figure 10. Open in new tabDownload slide Seismic amplitude profiles across the 104 well. Figure 10. Open in new tabDownload slide Seismic amplitude profiles across the 104 well. Figure 11. Open in new tabDownload slide (a) Fracture predictions of the coherence and (b) minimum curvature attributes and in the LMX Formation. Figure 11. Open in new tabDownload slide (a) Fracture predictions of the coherence and (b) minimum curvature attributes and in the LMX Formation. Figure 12. Open in new tabDownload slide Coherence attribute profile across the 104 well. Figure 12. Open in new tabDownload slide Coherence attribute profile across the 104 well. Figure 13. Open in new tabDownload slide Amplitude slices from different azimuths: (a) 0–37°, (b) 37–90°, (c) 90–135° and (d) 135–180°. Figure 13. Open in new tabDownload slide Amplitude slices from different azimuths: (a) 0–37°, (b) 37–90°, (c) 90–135° and (d) 135–180°. Figure 14. Open in new tabDownload slide Fracture strength and direction of the LMX Formation based on the P-wave anisotropy. Figure 14. Open in new tabDownload slide Fracture strength and direction of the LMX Formation based on the P-wave anisotropy. Figure 15 shows the TOC intersection analysis and horizontal well section display. Figure 15a shows the three-parameter intersection analysis of Eρ, Poisson's ratio and TOC. According to the intersection analysis, the high TOC shale is clearly separated from the surrounding rock, which lays the foundation for pre-stack inversion identification of high-quality shale. We can determine the Eρ range, where the TOC content is higher than 4% (figure 15a; marked with a red rectangle). Then, the value range of the red rectangle is projected onto the Eρ profile (figure 15b) and the high TOC (≥4%) distribution can be obtained. After determining the range of TOC values higher than 4%, the 3D geological body carving technique is used to identify the high TOC volumes. The top and bottom of the high TOC volume are tracked and we can obtain the TOC time thickness. Then, a time-depth conversion can be used to achieve the thickness of high-quality shale. Figure 16 shows the plane distribution of the target layer TOC and high-quality shale thickness. Figure 16 parts a and b show the high TOC (≥4%) plane prediction by the previously mentioned intersection analysis and thickness, respectively. From the results of TOC distribution and high-quality thickness, the area of TOC ≥ 4% in the study area is relatively small, and the plane distribution is also not continuous. Furthermore, we analysed the brittleness distribution characteristics of the study area. Figure 17 is the Fusion display of the high-quality shale and seismic amplitude profiles. Obviously, TOC on the left of the well is more enriched than right. Figure 18 shows the brittleness intersection analysis and plane characteristics of the study area. Figure 18a shows the three-parameter analysis of Eρ, Poisson's ratio and quartz content. Similar to the TOC intersection analysis described above, when the quartz content is higher than 60%, the corresponding Eρ value can be obtained, which indicates the high brittleness zone as shown by the rectangle marked in figure 18a. After we obtain the Eρ value range corresponding to the high brittleness obtained by the intersection analysis, the range is projected to the Eρ attribute profile and the brittleness body can be carved similar to that of the TOC mentioned earlier. Figure 18b shows a greater than 60% brittleness plane distribution. It is generally believed that a brittleness zone value greater than 60% is relatively high, which is suitable for hydraulic fracturing (Rickman et al.2008). We have found that there is only a small proportion with a brittleness greater than 60% in the study area, and the brittleness distribution is not continuous. Wells 104 and YSH1–1 are located in a relatively low brittleness area. Although horizontal well YSH1–1 is located in a favourable TOC area, the brittleness is relatively low. As the brittleness of the horizontal well zone is relatively low and fractures have not developed above, the two factors may result in the shale reservoirs not forming effective network communication during hydraulic fracturing, which may be the main reason why shale gas production is unsatisfactory for YSH1–1. Figure 19 shows the Eρ profile across YSH1–1. The P-wave velocity is significantly higher, which corresponds to the high quartz content. Additionally, the P-wave velocity is significantly lower in the high TOC section. The high quartz content region and high TOC region have significant differences in Eρ attribute profiles, and a high quartz content (high brittleness) corresponds to a high (red-yellow in figure 19) Eρ attribute value, whereas a high TOC content corresponds to a low (blue in figure 19) Eρ attribute value. The Eρ profile seismic response is consistent with the P-wave velocity log characteristic, which indicates that using the Eρ attribute to describe the TOC content and brittleness is reasonable. Figure 15. Open in new tabDownload slide TOC intersection analysis and horizontal well section display: (a) three-parameter intersection analysis of Eρ, Poisson's ratio and TO, and (b) Eρ profile for horizontal well YSH1–1. Figure 15. Open in new tabDownload slide TOC intersection analysis and horizontal well section display: (a) three-parameter intersection analysis of Eρ, Poisson's ratio and TO, and (b) Eρ profile for horizontal well YSH1–1. Figure 16. Open in new tabDownload slide (a) Plane distributions of the TOC content and (b) high-quality shale thickness. Figure 16. Open in new tabDownload slide (a) Plane distributions of the TOC content and (b) high-quality shale thickness. Figure 17. Open in new tabDownload slide Fusion display of the high-quality shale and seismic amplitude profiles. Figure 17. Open in new tabDownload slide Fusion display of the high-quality shale and seismic amplitude profiles. Figure 18. Open in new tabDownload slide (a) Brittleness intersection analysis and (b) plane characteristics of the Zhao104 well area. Figure 18. Open in new tabDownload slide (a) Brittleness intersection analysis and (b) plane characteristics of the Zhao104 well area. Figure 19. Open in new tabDownload slide Eρ profile response versus P-wave velocity log characteristic. Figure 19. Open in new tabDownload slide Eρ profile response versus P-wave velocity log characteristic. The shale gas favourable area mainly considers several key parameters: high TOC content, large thickness of high-quality shale, relatively developed fractures and relatively high brittleness. Additionally, conventional factors, such as burial depth and structure, must also be considered. Figure 20 shows the comprehensive analysis of the 104 area, parts a–d are high TOC plane distributions, high-quality shale thickness distributions, high brittleness distributions and micro-fracture distributions, respectively. We have selected relatively favourable areas for every evaluation parameter, and the red dotted circle is the relatively favourable area corresponding to every evaluation parameter (figure 20). Finally, all evaluation parameters are integrated, and we attempt to find overlapping areas of these evaluation parameters; two favourable areas, named A and B (figure 20a), are preferably selected in the study area. Overall, the areas that meet several evaluation parameters simultaneously (e.g. high TOC content, large high-quality shale thickness, relatively developed fractures and relatively high brittleness), are very few. Wells 104 and YSH1–1 are at the edge of the favourable areas, where the TOC content and high-quality shale thickness are relatively suitable, but the brittleness is relatively low and fractures are not developed, which may be the reasons that the horizontal well production is not very good. In the study, there is also a zone, named C in figure 20a, which has a high TOC content, large high-quality shale thickness, relatively developed fractures and high brittleness, as shown in the red ellipse in figure 20a, specifically in the area where the fractures are significantly better than those in A and B. However, the C area is very close to a large fault, which is not favourable for the preservation of a shale gas reservoir. Furthermore, C is located at the boundary of the study area, and the reliability of the favourable area prediction results is poor; thus, area C is not the preferred favourable area zone. Figure 20. Open in new tabDownload slide Comprehensive analysis of the favourable area prediction in the Zhao104 well area: (a) high TOC plane distribution, (b) high-quality shale thickness distribution, (c) high brittleness distribution and (d) fracture distribution. Figure 20. Open in new tabDownload slide Comprehensive analysis of the favourable area prediction in the Zhao104 well area: (a) high TOC plane distribution, (b) high-quality shale thickness distribution, (c) high brittleness distribution and (d) fracture distribution. In the subsequent exploration and deployment stages, because the favourable area is relatively small based on the results, the structure of the area is complex and, considering the economic feasibility of shale gas development, commercial horizontal well group deployment has not been carried out in this area. The 3D seismic favourable area prediction results have directly determined the subsequent horizontal well location plans, which has greatly reduced the exploration risk. 6. Discussion and conclusions Petro-physical modelling forms the basis for favourable shale gas area predictions, and mineral disturbance analysis can help us determine the most sensitive elastic factors in the TOC content and brittleness. The logging curve correction is the basis for rock physics modelling. Based on post-stack seismic data, coherence and curvature attributes can be used to characterize the large-scale fracture distribution in shale reservoirs. The P-wave anisotropy can be used to finely describe small-scale fracture characteristics, and the comprehensive use of these two methods can better characterize the fracture characteristics of shale reservoirs. The Eρ attribute is constructed by seismic pre-stack elastic inversion. Eρ, Poisson's ratio and brittleness or TOC content have been analysed to determine the Eρ range, which corresponds to high brittleness or high TOC levels. The Eρ attribute is used to indicate the characteristics of brittleness or TOC distributions. TOC and high-quality shale thickness predictions show that the area where the TOC content is higher than 4% is relatively small and not continuous in this area. Additionally, the area where the brittleness is greater than 60% is also very small and not continuous, and the favourable areas are relatively small in the entire study area. Well 104 and horizontal well YSH1–1 are located at the edge of the favourable area; the TOC content and high-quality shale thickness are relatively suitable, but the brittleness is relatively low and the fractures are not developed, which may be the reasons that the shale gas production has not been very good for the horizontal well. In addition, well 104 and horizontal well YSH1–1 are not located in the best TOC area and the largest thickness of high-quality shale, specifically for well 104, which is located in a relatively poor TOC area. As shown in figure 17, we find that the TOC content is very high to the left of the 104 well. If well 104 can be moved to the left area, it may be more beneficial to well deployment. The thicknesses of high-quality shale are approximately 10–40 m for well 104 and approximately 10–20 m for horizontal well YSH1–1, which is not located in the area with the largest thickness of high-quality shale. If the two well positions can be moved to the left or down to the dotted area, this may be better. But whether the left part of the 104 well is more beneficial to well deployment is not only determined by TOC and thickness. Also, even if we have taken many measures to make our results more convincing, there is only one well in the study area, which has brought us more solutions. We also need to consider the other geological parameters and get more well information to continue analysis. It is worth remembering that the elastic sensitive parameters of the TOC content and brittleness are different for different shale gas areas. Although the Eρ attribute has achieved suitable results in this area, the attribute may not be applicable in other shale gas areas, so petro-physical modelling must be carried out according to the actual fields. Acknowledgements The authors acknowledge the support of the National Natural Science Foundation of China (Grant No. 41874168), the National Science and Technology Major Project of China (Grant No. 2017ZX05018004) and the Sichuan Science and Technology Program for Distinguished Young Scholars (Grant No. 2019JDJQ0053). The authors confirm that any identifiable participants in this study have given their consent for publication. Conflict of interest statement. None declared. References Banik N. , Koesoemadinata A. , Wager C. , Inyang C. , Agarwal V. , Priezzhev I. , 2014 . Predrill prediction of subsalt pore pressure from seismic impedance . The Leading Edge , 33 , 400 – 412 . Google Scholar Crossref Search ADS WorldCat Carpenter C. , 2018 . Unconventional-reservoir characterization with azimuthal seismic diffraction imaging . Journal of Petroleum Technology , 70 , 92 – 94 . Google Scholar Crossref Search ADS WorldCat Chopra S. , Sharma R.K. , Keay J. , Marfurt K.J. , 2012 . Shale gas reservoir characterization workflow , In: The 82nd Society of Exploration Geophysicists Annual Meeting , Las Vegas, USA , 4–9 November . Google Preview WorldCat COPAC Chen S. , Zhao W.Z. , Zeng Q.C. , Yang Q. , He P. , Gai S.H. , Deng Y. , 2018 . Quantitative prediction of total organic carbon content in shale-gas reservoirs using seismic data: a case study from the Lower Silurian Longmaxi Formation in the Chang Ning gas field of the Sichuan Basin, China . Interpretation , 6 , 153 – 168 . Google Scholar Crossref Search ADS WorldCat Connolly P. , 1999 . Elastic impedance . The Leading Edge , 18 , 438 – 452 . Google Scholar Crossref Search ADS WorldCat Convers C. , Hanitzsch C. , Curia D. , Davis T. , Tura A. , 2017 . Elastic parameter estimation for the identification of sweet spots, Vaca Muerta Formation, Neuquén Basin, Argentina . The Leading Edge , 36 , 948a1 – 948a10 . Google Scholar Crossref Search ADS WorldCat Dasgupta S. , Chatterjee R. , Mohanty S.P. , 2016 . Prediction of pore pressure and fracture pressure in Cauvery and Krishna-Godavari basins, India . Marine and Petroleum Geology , 78 , 493 – 506 . Google Scholar Crossref Search ADS WorldCat Gray D. , Head K. , 2000 . Fracture detection in Manderson field: a 3D AVAZ case history . The Leading Edge , 19 , 1214 – 1221 . Google Scholar Crossref Search ADS WorldCat Heath-Clarke M. , Taylor K. , David H. , Anthony F. , Fred H. , Matthew H. , Andy M. , 2016 . The characterization of unconventional reservoirs in the Bowland sequence using onshore 3D seismic data, Cleveland Basin, UK . First Break , 34 , 45 – 52 . WorldCat Jaiswal P. , Varacchi B. , Ebrahimi P. , Dvorkin J. , Puckette J. , 2014 . Can seismic velocities predict favourable areas in the Woodford Shale? A case study from McNeff 2–28 Well, Grady County, Oklahoma . Journal of Applied Geophysics , 104 , 26 – 34 . Google Scholar Crossref Search ADS WorldCat Kivi I.R. , Zare-Reisabadi M. , Saemi M. , Zamani Z. , 2017 . An intelligent approach to brittleness index estimation in gas shale reservoirs: a case study from a western Iranian basin . Journal of Natural Gas Science and Engineering , 44 , 177 – 190 . Google Scholar Crossref Search ADS WorldCat Liu X.W. , Guo Z.Q , Liu C. , Liu Y.W. , 2017 . Anisotropy rock physics model for the Longmaxi shale gas reservoir, Sichuan Basin, China . Applied Geophysics , 14 , 886 – 898 . WorldCat Liu J. , He Z.L. , Liu X.W. , Huo Z.Z. , Guo P. , 2019a . Using frequency-dependent AVO inversion to predict the ‘sweet spots’ of shale gas reservoirs . Marine and Petroleum Geology , 102 , 283 – 291 . Google Scholar Crossref Search ADS WorldCat Liu W. , He Z.H. , Cao J.X. , Zhang J.J. , Xu G. , Wan X.P. , Yu G. , 2019b . Integrated geology sweet spots and micro-seismic monitor to optimize reservoir stimulation-a case for shale gas, China , In: The International Petroleum Technology Conference , Beijing, China , 26–28 March, IPTC 19197 Google Preview WorldCat COPAC Makar K.H. , Kamel M.H. , 2011 . An approach for minimizing errors in computing effective porosity in reservoir of shaly nature in view of Wyllie-Raymer-Raiga relationship . Journal of Petroleum Science and Engineering , 77 , 386 – 392 . Google Scholar Crossref Search ADS WorldCat Marfurt K.J. , Sudhaker V. , Gersztenkorn A. , Crawford K.D. , Nissen S.E. , 1999 . Coherence calculations in the presence of structural dip . Geophysics , 64 , 104 – 111 . Google Scholar Crossref Search ADS WorldCat Passey Q.R. , Creaney S. , Kulla J.B. , Moretti F.J. , Stroud J.D. , 1990 . A practical model for organic richness from porosity and resistivity logs . AAPG Bulletin , 74 , 1777 – 1794 . WorldCat Rickman R. , Mullen M. , Petre E. , Grieser W.V. , Kundert D. , 2008 . A practical use of shale petro physics for stimulation design optimization: all shale plays are not clones of the Barnett shale , In: The 2008 SPE Annual Technical Conference and Exhibition , Denver, Colorado, USA , 21–24 September , Society of Petroleum Engineer , SPE 115258 . Google Preview WorldCat COPAC Ruger A. , 1998 . Variation of P-wave reflectivity with offset and azimuth in anisotropic media . Geophysics , 63 , 935 – 947 . Google Scholar Crossref Search ADS WorldCat Rybacki E. , Meier T. , Dresen G. , 2016 . What controls the mechanical properties of shale rocks? – Part II: brittleness . Journal of Petroleum Science and Engineering , 144 , 39 – 58 . Google Scholar Crossref Search ADS WorldCat Sena A. , Castillo G. , Chesser K. , Voisey S. , Estrada J. , Carcuz J. , Carmona E. , Hodgkins P. , 2011 . Seismic reservoir characterization in resource shale plays: stress analysis and sweet spot area discrimination . The Leading Edge , 30 , 758 – 764 . Google Scholar Crossref Search ADS WorldCat Tavella J. , Moirano J. , Zambrano J.C. , Zanca G. , García M.V. , Soto H.S , Gutierrez C.G. , 2013 . Integrated characterization of unconventional Upper Jurassic Reservoir in Northern Mexico , In: The 75th EAGE Conference & Exhibition incorporating SPE EUROPEC , London, UK , 10–13 June . Google Preview WorldCat COPAC Verma S. , Zhao T. , Marfurt K.J. , Devegowda D. , 2016 . Estimation of total organic carbon and brittleness volume . Interpretation , 4 , 373 – 385 . Google Scholar Crossref Search ADS WorldCat Wang K.N. , Sun Z.D. , Dong N. , 2015 . Prestack inversion based on anisotropic Markov random field-maximum posterior probability inversion and its application to identify shale gas sweet spots . Applied Geophysics , 12 , 533 – 544 . Google Scholar Crossref Search ADS WorldCat Whitcombe D.N. , 2002 . Elastic impedance normalization . Geophysics , 67 , 60 – 62 . Google Scholar Crossref Search ADS WorldCat Xu S. , White R.E. , 1996 . A physical model for shear wave velocity prediction . Geophysical Prospecting , 44 , 687 – 717 . Google Scholar Crossref Search ADS WorldCat Yu G. , Zhang Y.S. , Paola N. , Azizov I. , 2014 . Rock physics diagnostics and modeling for shale gas formation characterization in China , In: Beijing 2014 International Geophysical Conference & Exposition , Beijing, China , 21–24 April . Google Preview WorldCat COPAC Yu G. , Zhang Y.S. , Wang X.M. , Liang X. , Liu W. , Guo R. , Strecker U. , Smith M. , 2017 . Shale gas reservoir characterization and sweet spot prediction in China . First Break , 35 , 59 – 63 . WorldCat Zhao D.F. , Guo Y.H. , Yin S. , Ren C.Y. , Wang Y.J. , 2018 . Prediction of geo-mechanical sweet spots in a tight gas sandstone reservoir: A case study of lower Permian strata in the southern Qinshui Basin, China . Interpretation , 7 , 207 – 219 . Google Scholar Crossref Search ADS WorldCat Zeng H.L. , Wang W. , Liang Q.S. , 2017 . Seismic expression of delta to deep-lake transition and its control on lithology, total organic content, brittleness, and shale-gas sweet spots in Triassic Yanchang Formation, southern Ordos Basin, China . Interpretation , 5 , SF1 – SF14 . Google Scholar Crossref Search ADS WorldCat Zeng Q.C. et al. . 2018 . Quantitative prediction of shale gas sweet spots based on seismic data in Lower Silurian Longmaxi Formation, Weiyuan area, Sichuan Basin, SW China . Petroleum Exploration and Development , 45 , 422 – 430 . Google Scholar Crossref Search ADS WorldCat Zhang Y.S. , Yu G. , Liang X. , 2014 . The Marine Shale Gas-rich Zones Geophysical Prediction and Evaluation Techniques in South China , In: The 76th EAGE Conference and Exhibition , Amsterdam , 16 – 19 , June . Google Preview WorldCat COPAC © The Author(s) 2019. Published by Oxford University Press on behalf of the Sinopec Geophysical Research Institute. This is an Open Access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted reuse, distribution, and reproduction in any medium, provided the original work is properly cited. TI - Combined petrophysics and 3D seismic attributes to predict shale reservoirs favourable areas JF - Journal of Geophysics and Engineering DO - 10.1093/jge/gxz060 DA - 2019-10-01 UR - https://www.deepdyve.com/lp/oxford-university-press/combined-petrophysics-and-3d-seismic-attributes-to-predict-shale-bayFAyxA09 SP - 974 VL - 16 IS - 5 DP - DeepDyve ER -