Identifying shale sweet spots by 3D microseismic imaging

Several studies suggest that shale sweet spots are likely associated with low Poisson’s ratio in the shale layer. Compared with conventional geophysical techniques with active seismic data, it is straightforward and cost effective to delineate the distribution of 3D Poisson’s ratio using microseismic data. In this study, we develop a method for determining microseismic event locations, 3D P-wave velocities, and the Poisson’s ratio model with data recorded from downhole monitoring arrays. The method combines the improved 3D traveltime tomography, which can invert P and S arrivals for 3D P-wave velocity and Poisson’s ratio structures simultaneously, and a 3D grid search approach for event locations in an iterative fashion. The traveltime tomography directly inverts the Poisson’s ratio strucure instead of calculating the Poisson’s ratio from P- and S-wave velocities that are inverted by conventional traveltime tomography separately. The synthetic results and analysis suggust that the proposed method is capable of directly recovering the true Poisson’s ratio model. We also apply the method to a field dataset, and the results show that it helps delineate the structure of the reservoir and identify shale sweet spot.


INTRODUCTION
Sweet spot is geologically the concurrence of several favorable geologic parameters such as fracturing, thermal history, gas content, reservoir thickness, or matrix rock properties, and economically areas that were most likely to be the developed parts of a continuous accumulation (Schenk, 2001). Norton and Maxwell (2013) suggest that a shale sweet spot where it may offer high gas production potential should correspond to low Poisson's ratio. With the increase of light hydrocarbon saturation, the compressional wave velocity decreases and shear wave velocity increase, thereby resulting into low velocity ratio (Vp/Vs ratio) or Poisson's ratio. They use the results of the Poisson's ratio distribution from surface seismic data interpretation, obtained by conventional seismic processing, to infer such sweet spots in shale.
Poisson's ratio or velocity ratio (Vp/Vs ratio) is a potential indicator for the reservoir lithology and fluid properties (Holbrook et al., 1988;Ikwuaor, 2006). In the field of active seismic exploration, many studies (i.e., Tatham, 1976;McCormack et al., 1985) adopt the Vp/Vs ratio for hydrocarbon identification and stratigraphic interpretations based on the reflection seismic data. Anderson and Lines (2008) review three methods estimating Vp/Vs ratio that can be used to extract rock-properties information from seismic data, including: AVO analysis, followed by post-stack inversion of the AVO attributes (Goodway, et al., 1997); prestack inversion of the Pwave gathers directly to impedances (Hampson et al., 2005); joint prestack inversion of the P-wave data and converted wave data directly to impedance (Hampson et al., 2005).
Although these techniques, utilizing the reflection seismic data, are recognized as Non-peer reviewed EarthArXiv preprint 4 effective tools to identify the sweet spot area in the reservoir, the high cost and long processing time limit the use of the techniques routinely (Martakis et al., 2003;Martakis et al., 2006).
In the field of passive seismology, 3D passive seismic traveltime tomography has become a well established technique to help understanding of seismotechonics and seismogenic processes over large areas (Thurber et al., 1995;Eberhart-Phillips and Michael, 1998;Chiaraba and Amato, 2003). Different aspects of the method have been reviewed in Thurber (1986), Kissling (1988), and Iyer and Hirahara (1993). Passive seismic or microseismic traveltime tomography has been also successfully applied in the hydrocarbon and geothermal exploration. Several studies introduce the role of microseismicity in the reservoir monitoring (e.g., Rutledge et al., 1998;Maxwell and Urbancic, 2001;Durham, 2003), and demonstrate its potential to delineate the reservoir for a relatively low cost compared to conventional 3D seismic surveys (Durham, 2003). Valoroso et al. (2008) apply 4D microseismic tomography to detect space-time dependency in response to fluid pressure. Zhang et al. (2009) show that microseismic traveltime tomography has a great potential for reservoir imaging and property estimation using induced seismicity of an oil field. Tselentis et al. (2011) present a case study in which high-resolution 3D Poisson's ratio variations are inferred by estimating the P-and S-wave velocities with passive seismic tomographic method. Barthwal et al. (2014) apply traveltime tomography to a microseismic dataset acquired during a mining operation and achieve a reasonable 3D velocity structure. Many studies (i.e., Chatterjee et al., 1985;Julian et al., 1996;Lees and Wu, 2000;Lees, 2002;Hauksson and Unruh, Non-peer reviewed EarthArXiv preprint 5 2007;Kummerow et al., 2012;Zang et al., 2014) also suggest that the tomographic measurement of Poisson's ratio or Vp/Vs ratio anomaly from earthquake or induced sesimicity is a promising technique for identifying geothermal resources.
Microsesimic activity induced by hydraulic fracturing is often monitored with geophone arrays in nearby monitoring wells. To infer shale sweet spots characterized by low Poisson's ratio, the effort involves the determination of microseismic event locations and 3D velocity models or a Poisson's ratio model with microseismic arrival time picks. It is well known that microseismic event locations are coupled with velocity models (Maxwell, 2009). To address the coupling problem, Thurber (1992) introduces joint inversion and altenating inversion with parameter separation, both methods can produce accurate models. Many studies (Thurber, 1986;Iyer and Hirahara, 1993;Thurber and Rabinowitz, 2000;Zhang et al., 2009;Concha et al., 2010;Jansky et al., 2010;Zhou et al., 2010;Tselentis et al., 2011) apply the alternating strategy, estimating a velocity model by fixing event parameters and then updating the velocity model based on the event location parameters, to solve the event-location problem. The approach does not require weight balancing tbetween event locations and model parameters. In this study, we take similar ideas but constrain 3D models with layer interfaces, using 3D Vp and Poisson's ratio grid parameters to represent the area, and allowing arbitrary Vp and Poisson's ratio variations within each layer. Directly inverting Poisson's ratio instead of calculating Poisson's ratio from Vp and Vs results allows applying contraints directly on Poisson's ratio and producing more reliable results (Thurber et al., 1993;Thurber and Eberhart-Phillips, 1999). The number of layers and layer thickness are Non-peer reviewed EarthArXiv preprint 6 fixed following well-log interpretation. A 3D grid search (Rodi and Toksӧz, 2001) is applied to search for the event locations in the velocity models (note that Vs is created by Vp and Poisson's ratio), followed by an improved 3D traveltime tomography to update Vp and Poisson's ratio models for updating event locations, and then the process is iterated until minimum misfits for both grid search and traveltime tomography are achieved. We shall describe the methodology in details and validate the method with two numerical and one field data examples.

METHOD
The Poisson's ratio or velocity ratio is an important indicator for the reservoir lithology and fluid properties (Ikwuaor, 2006). And P-wave velocity helps characterize the reservoir structure (Tselentis et al., 2011). Also, conventional technique taking the ratio of P-and S-wave velocities cannot ensure that the area is well resolved in which the Pand S-wave raypaths are not identical (Eberhart-Phillips, 1990). In this section, we therefore attempt to describe an alternating method to estimate the 3D P-wave velocity as well as Poisson's ratio and correpsonding event locations using the P-and S-wave arrival times recorded by the receivers in borehole. Analogous to the alternating strategy, we separate the whole inversion system into two parts: one is to solve for the event locations and origin time using a 3D grid search approach, the other one is to invert for the 3D P-wave velocity and Poisson's ratio models simultaneously in a modified traveltime tomography. Hence, the modified traveltime tomography and the 3D grid search are iteratively employed in the alternating way until the arrival time Non-peer reviewed EarthArXiv preprint 7 misfit residual becomes negligible. In the following, we shall explicitly introduce the employed approaches in the alternating method.
We apply the 3D grid search method to optimize the event location parameters before updating the P-wave velocity and Poisson's ratio, since the grid search approach does not require a initial guess of the locations and is able to avoid trapping into a local minimum. In microseismic data processing, the observed arrival time ij t from a microseismic event i to a receiver j in well is formulated as,   Popovici, 1999) to calculate the synthetic P-and S-wave traveltimes. Note that in the S-wave raytracing, the utilized S-wave velocity model is converted from the initial or updated P-wave velocity and Poisson's ratio models, which are obtained in the aftermentioned traveltime tomography.
During the searching process, the azimuth and incident angles are also utilized to constrain the event locations. After determining the event parameters using 3D grid search, we utilize the searched event locations to invert for the P-wave velocity and Poisson's ratio models with the computed traveltimes of P-and S-wave. Based on a priori information, the linear approximation of the traveltime residuals can be expressed to the perturbation of P-wave velocity or Poisson's ratio model parameters, r is the P-or S-wave arrival time residuls, obs ij T is the observed P-or S-wave arrival time data achieved in eqaution (4), cal ij T represents the calculated arrival times based on the a priori information,   is the origin time perturbation, x  refers to the event location perturbation, m  is the perturbation to the P-wave velocity or Poisson's ratio. For P-and S-wave, the equation (5) is therefore expressed as follows, (7) can be transformed step by step as follows, By combining the equation (6) and (8), we invert P-and S-wave traveltimes for Pwave velocity and Poisson's ratio models at the same time. From the perspective of the inverse problem, we can define the objective function of the traveltime tomography as, To solve the above equation and obtain the model parameter updates, we utilize the Gauss-Newton method to linearize the equation (9) and apply the conjugate gradient method to iteratively invert for Vp and Poisson's ratio updates. Hence, the traveltime tomography for Vp and Poisson's ratio models can be expressed as, or Poisson's ratio ( m  ), which corresponds to the partial derivatives in equation (6) and (8) (10). We shall test and show the validity of our inversion in the synthetic test.
We describe the approches of 3D grid search and traveltime tomography. During the performance of the alternating method, both approches will be iteratively employed to update the corresponding event locations as well as Vp and Poisson's ratio models.
To clearly understand the alternating method, we summerize the entire procedure as follows:

SYNTHETIC EXAMPLES
We conduct two synthetic cases to test the efficiency of our proposed method. In the first case, two observation wells along with a red cloud with 726 events are distributed in Figure 1. In each well there are 20 receivers with the spacing of 20 m. With the assumption of flat layer interfaces, we design two layered grid models (Vp and Vs) with an abnomal volume embeded into each model. Then we derive the true Poisson's ratio model based on the Vp and Vs models. As shown in Figure 2(a-c) and 3(a-c), the anomalies are supposed as the zone of hydraulic fracturing in a stage. During the fracturing treatment, we assume that the red event cloud (see Figure 1) is triggered by fractures with shear or tensile failure, which may enhance the porosity and permeability decrease the Poisson's ratio significantly. The stimulated region or anomaly is usually veiwed as the shale sweet spot. We obtain the initial models by perturbing the true Vp and Vs layered models by 5% and 10%, respectively. The initial models are assumed from the sonic logs or calibrated models. The initial Poisson's ratio model is further obtained based on initial Vp and Vs models.
After velocity models and survey geometry are set up, we apply the 3D wavefront raytracing to calculate the traveltime data with the true event locations in the true models. The synthetic data are viewed as the observed or picked arrival time data in the field situation. We invert the Vp and Poisson's ratio and utilize the grid search method to refine the event location during each loop, and the procedure could be terminate until the event locations are unchanged. Figure 1 shows that the event locations obtained by the grid search method are well recovered in comparison with the true locations. The maximum location error is about 30 m. We recognize that the error of the event location in Y axis is smaller than the location in X axis due to the constraint of two receiver array.
Besides the event locations, the final results, including Vp and Poisson's ratio models, are inverted with the layer constraint of the initial models. Figure   During the whole procedure, we calculate the data misfit for every loop, which includes one 3D grid search for event location and one traveltime tomography for Vp and Poisson's ratio models simultaneously. Figure 5a shows the convergence curve of the data misfit over four loops. We only show misifit of four loops because additional loops could not contribute to the data misfit decrease. Additionally, in each traveltime tomography ten iterations are performed to obtain stable solutions. In the last loop, the traveltime misfit of each iteration is shown in Figure 5b. In the traveltime tomgoraphy of each loop, we input the initial models, instead of the Vp and Poisson's ratio solutions inverted by the previous loop, into the traveltime tomography, because it may avoid the the bias in the previous models during the traveltime tomography.
Under some circumstances, the stage of hydraulic fracturing is away from the monitoring wells, which may result in poor raypath coverage. To investigate whether our method works on the situation or not, we design another case in which the anomaly is set up farther away from the observation wells than before ( Figure 6). The anomaly is still assumed as the sweet spot. In Figure 6, the true red cloud is also regarded as the event cloud generated by some fractures with shear or tensile failure in the area of hydraulic fracturing. Similar to the previous case, we invert for the location and Vp as well as Poisson's ratio results with our proposed procedure. Figure 6 shows the comparison of true and searched event locations on different views, which reveals that the event locations are reasonably recovered. The maximum error is about 80 m in the location results. Figure 7 shows the comparsion of different maps of true, initial, and inversion Vp models. The comparison demonstrates that there is a reasonable recovery in the abnormal area. Compared with true Poisson's ratio model (see Figure 8a-c), the inverted Poisson's ratio model (Figure 8g-i) presents that the abnormal area (marked by black rectangle) has a reasonable recovery. The convergence curves of the data misfits both in the whole procedure and in the traveltime tomography of last loop are displayed in Figure 9. Hence, our method also works in this situation, although it is not better than the previous case.

FIELD EXAMPLE
We further apply our method to a field dataset. To gain enough raypath coverage and ensure the quality of the picked data, we select 156 effective events recorded by both receiver arrays in stage 5. Figure 10 exhibits the recordings of 2 receiver arrays of one event, and the arrivals of P and S waves can be clearly picked up. We construct the initial Vp and Vs layer models from the sonic logs and convert them to 3D Vp and Vs grid models. Based on the grid models, we calculate the Poisson's ratio model. The map views of 3D initial Vp and Poisson's ratio models are displayed in Figure 12(a-c) and 13(a-c), respectively. With the proposed procedure, we obtain the event location and Vp as well as Poisson's ratio model solutions until the event locations remain unchanged. Figure 11a shows the 3D distribution of our event locations (blue dots) and event location results (green dots) provided by a service company. The comparison demonstrates that the average location difference is about 30 m. We also project both event location results on different views (Figure 11b-d). We observe that our event locations remain a similar distribution as green ones, even our event locations are more convergent than the green ones. Seven red stars in Figure 11 represent the perforations in the fifth stage. From the X-Y plan view, the event locations are agreed to perforations, which is attributed to the constraint of two observation wells. From the X-Z or Y-Z plan view, the event location on the depth seems consistent with that of events from the service company based on the calibrated layer models.
The different views of Vp and Poisson's ratio inversion results are shown in Figure   12(d-f) and 13(d-f), respectively. From various plan views, the relative low Poisson's ratio region (marked by black rectangle) is observed. Around the depth of 1.9 Km, the Poisson's ratio is relatively the lowest in the region marked by black rectangle. Based on the Poisson's ratio distribution, we could delineate the area of the sweet spot and estimate the stimulated volume of the hydrocarbon. We also provide the data misfit curve in Figure 14 to show the convergence of our method. We execute many loops in the whole procedure, and observe that the data misfit is almost unchanged after four loops. The data misifts of four loops are presented in Figure 14a. We also display the misfit cuvers of P-and S-wave in the traveltime tomography of the last loop.

In our synthetic examples, we show the feasibility of mapping the sweet spot or low
Poisson's ratio distribution using two observation wells, however, the event location and velocity solutions are obviously better in the first case than in the other case due to better ray coverage. As the distance becomes larger between events and receivers, it will be more difficult to recover the event locations and structural models. We also test our method in the case of one monitoring well, the results reveal that the event location and velocity could not be resolved effectively, due to the trade-off between the event locations and structural models. We futher test our method with three observation wells in the form of triangle distribution, the event location and velocity results are recovered very well. Hence, with the borehole monitoring, the poor raypath coverage resulted from the one well may not recover the event locations and velocity models. However, when using at least two observation wells with suitable positions relative to the area of located events, our method is promising for reconstructing Vp and Poisson's ratio models.
Because the 3D grid search is used to locate the events, we do not need any a priori event location. The initial velocity models can be constructed directly from the sonic logs. Upon the high-quality observed data and initial models are prepared well, our method is able to derive a reliable Poisson's ratio model and build the sweet spot map.
In a complicated geological situation, the P and S ray may differ significantly so that based on the Vp and Vs models the Poisson's ratio result may be biased (Nicholson and Simpson, 1985). Under this circumstance, conventional techniques, such as taking the ratio of Vp and Vs (Nakamura et al., 2003) and inverting S-P times ratio (Chiarabba and Amato, 2003;Zhang et al., 2009), are inaccurate to estimate Poisson's ratio or Vp/Vs ratio models. However, in this study the raypaths of P-and S-wave are not assumed to be identical, our method is therefore not contrained by the poorly resolved portions of the Vp model. In the synthetic test, we have shown that our traveltime tomography can produce more accurate P-wave velocity and Poisson's ratio solutions than the traditional technique taking the ratio of Vp and Vs. It is noted that, since Pand S-wave raypaths are similar for the short distance between events and receivers in our real cases, our solutions do not show an obvious improvement compared to that of traditional technique.
We recognize that the induced event cluster is clearly associated with low Poisson's ratio marked by black rectangle in Figure 13. The low Poisson's ratio anomaly could be interpreted due to changes of the rock properties of the shale when the fractures are triggered during high-pressure fluid injection. The fracturing treatment may improve the porosity and permeability of the target zone, meanwhile, the hydrocarbon may fill in the fractures and then the Poisson's ratio would be decreased.
Hence, the low Poisson's ratio may be adopted as a proxy for the hydrocarbon or sweet spot prediction. In the field example, due to limitation of effective data, we only employ our method to delineate the sweet spot distribution in the area of stage 5. We believe that our method is also helpful to build sweet spot map for other stages or mult-stage treatments. The produced sweet spot map can be utilized to estimate the stimulated reservoir volume during the hydraulic fracturing.

CONCLUSIONS
In this study, we propose an alternating method, which combines the 3D traveltime tomography with the 3D grid search, to recover 3D Vp as well as Poisson's ratio models and event locations in 3D. Our synthetic tests show that the proposed method is feasible to reconstruct 3D Vp and Poisson's ratio models directly in the event area. Our synthetic results indicate that the low Poisson's ratio area, assumed as the sweet spot, achieves a good recovery in the Poisson's ratio model. In the field example, we obtain the reasonable location results and a promising Poisson's ratio model. The location results could be used to characterize the fractures. The low Poisson's ratio anomaly may also be used to estimate the stimulated reservoir volume. Therefore, with sufficient ray coverage, our alternating method may obtain the distribution of the Poisson's ratio, which should be helpful to identify the sweet spot in shale reservoirs.     The circle represents the data misfit of each loop, and the triangle and square are the traveltime misfits over each iteration of P and S waves, respectively.    loop. The circle represents the data misfit of each loop, and the triangle and square are the traveltime misfits over each iteration of P and S waves, respectively.

Figure 10
The different recordings of the same microseismic event from (a) eightreceiver array and (b) twelve-receiver array, respectively. The P-and S-wave arrivals display a high signal-to-noise ratio and are easily picked (red and blue respectively) on all 20 receivers.        The circle represents the data misfit of each loop, and the triangle and square are the traveltime misfits over each iteration of P and S waves, respectively.