Three-dimensional seismic characterization and imaging of the Soda Lake geothermal field
TThree-dimensional seismic characterization andimaging of the Soda Lake geothermal field
Kai Gao ∗ ,1 , Lianjie Huang , and Trenton Cladouhos Los Alamos National Laboratory, Geophysics Group, MS D446, Los Alamos,NM 87545, USA Cyrq Energy, Inc., 4010 Stone Way North, Suite 400, Seattle, WA 98103, USA
Abstract
Accurate characterization of subsurface geophysical properties and detection of the faultsystem are essential for geothermal energy exploration and production. The Soda Lake geother-mal field is in western Nevada with a complex fault system. Previous seismic character-ization only produced a low-resolution, smooth velocity model along with a simple, con-ceptual fault model. Using optimized correlation-based full-waveform inversion, wavefield-separation-based reverse-time migration, and automatic fault detection techniques, we present3D seismic characterization for the Soda Lake geothermal field using 3D surface seismic dataacquired with Vibroseis sources. We obtain 3D high-resolution velocity, density, and acous-tic impedance models, 3D seismic images with different grid spacings, and a high-resolutionfault system. Consistency check between the constructed faults and currently active injectionand production geothermal wells verifies that our seismic inversion and imaging results anddetected faults are reliable. These results can provide valuable information for optimizing wellplacement and geothermal energy production at the Soda Lake geothermal field.
Keywords:
Full-waveform inversion, fault detection, reverse-time migration, Soda Lake geother-mal field, surface seismic data ∗ Corresponding Author: [email protected] (K. Gao); [email protected] (L. Huang) a r X i v : . [ phy s i c s . g e o - ph ] A ug Introduction
Geothermal energy is an increasingly important component in the renewable energy sector ofthe United States. In the continental U.S., western states and Hawaii host most of the geothermalenergy resources. The Soda Lake geothermal field is owned and operated by Cyrq Energy, Inc.,and is located eight miles north-west of Fallon, Nevada, one of the most geothermal-rich states.Figure 1 shows the location of the Soda Lake geothermal field. The Soda Lake geothermal plantis a binary geothermal electric-generating facility. Soda Lake I operated from 1987 through 2018,and Soda Lake II operated from 1991 through 2019. Figure 1 shows the location of the Soda Lakegeothermal field in the U.S. continent and its location in western Nevada.In this study, we aim to provide a high-resolution, reliable, 3D seismic characterization for theSoda Lake geothermal field using 3D surface seismic data acquired in 2009. Previous studies con-structed a low-resolution velocity model for the Soda Lake geothermal field, and built a conceptualfault model based on geophysical imaging and geological analysis. Our goal through this studyis to reveal the complexity of subsurface geophysical properties in this region with our recentlydeveloped seismic-waveform inversion algorithm, and to construct a high-resolution fault systemusing our advanced seismic migration imaging and fault construction algorithms in a deterministicand systematic manners.The 3D surface seismic data acquired at the Soda Lake geothermal field contain a total of8,321 compressional-to-compressional (PP)-component common-shot gathers. Figure 2 shows thedistributions of sources (red dots) and receivers (blue triangles) in this 3D seismic survey. Eachcommon-shot gather covers a surface area of up to approximately × km . The entire surveyhas a clockwise azimuth angle of 33.97 ◦ , measuring from the positive Y direction (i.e., the Northdirection). There exist some “holes” in the source and receiver distributions where there are nodata acquired. The source lines are along the NE-SW direction with a line interval of 235.5 m.The receiver lines are along the NW-SE direction with a line interval of 167.5 m. The inline sourceinterval is 33.5 m, and the inline receiver interval is 67 m. The surface area covered by all thesources and receivers is approximately × km .To obtain high-resolution subsurface medium parameter models for the Soda Lake geothermalfield, we perform full-waveform inversion (FWI) of the acquired 3D surface seismic data. FWI is2 nonlinear inversion method for estimating subsurface medium properties by minimizing the dif-ference between observed seismic waveforms and synthetic waveforms (Tarantola, 1984; Virieuxand Operto, 2009). First developed in 1980s, FWI is becoming a standard seismic inversion toolin both theoretical research and industry, particularly for highly complex geology. Suffering fromhigh nonlinearity and the cycle-skipping issue, conventional FWI usually fails to converge to cor-rect results for field seismic data, even for complex synthetic models. We employ our recentlydeveloped FWI algorithm that uses a different waveform match criterion in our algorithm to mit-igate these difficulties as much as possible, even though not completely. Accompanied with ouradvanced, parallel full-waveform inversion codes that run on Los Alamos National Laboratory’ssuper-computing platform, we are able to perform a 3D full-waveform inversion for the Soda Lakegeothermal field to obtain a set of high-resolution subsurface medium parameter models up to2.5 km in depth.A subsurface structural image is one of the most important products from seismic characteriza-tion, which provides essential information for geothermal well placement and production. Previousstudies have produced a subsurface image using conventional ray-based Kirchhoff migration that isincapable of handling complex structures. To improve subsurface imaging, we apply reverse-timemigration to the 3D surface seismic data using our FWI-inverted 3D velocity model. Reverse-time migration (RTM) is an advanced seismic imaging technique for imaging complex structures(McMechan, 1983; Chang and McMechan, 1987). Using the PP seismic data separated from the3D multi-component seismic data, we perform a wavefield-separation-based low-artifact RTM (Feiet al., 2015) to obtain a high-resolution subsurface structural image for the Soda Lake geothermalfield. We subsequently perform a fault-enhancing processing on the 3D RTM image and an auto-matic fault detection to delineate faults from the 3D image volume.Through this comprehensive seismic characterization, we achieve a set of high-resolution andreliable 3D subsurface medium parameter models, a 3D structural image, and a 3D fault system forthe Soda Lake geothermal field. Our inversion and imaging results reveal a complex fault systemunravelling many faults that are not in the previous fault model. We find that current geothermalproduction wells either penetrate through the faults we detect from the image volume, or are fairlyclose to these faults. The consistency between the faults and currently active injection/productionwells validates the accuracy and reliability of our inversion and imaging results and detected faults.3ur high-resolution 3D subsurface medium property models and 3D fault system image can pro-vide valuable information for optimizing well placement and geothermal energy production at theSoda Lake geothermal field.Our paper is organized as follows: In the Methodology section, we briefly describe our seismicinversion and imaging methods applied to the 3D surface seismic data from the Soda Lake geother-mal field. In the Results section, we present and analyze our seismic inversion and imaging results.Particularly, we perform a check on accuracy and reliability of our imaging and fault detectionresults by plotting currently active geothermal injection and production wells in the 3D space. Wesummarize our findings in the Conclusions section. The seismic data we use in this study are pre-processed by Geokinetics, an industrial data pro-cessing company. Their data processing procedure includes random noise attenuation, groundrollattenuation, surface-consistent deconvolution and amplitude correction, spherical expanding am-plitude compensation, and optionally automatic gain control or time-variant spectral whitening,etc. Considering the complexity of 3D seismic data processing for such a large dataset, we per-form no additional data processing except frequency-domain filtering and offset-based data selec-tion based on our needs for our 3D seismic inversion and imaging algorithms briefly described inthis section.
Seismic processing by Geokinetics produced a smooth P-wave velocity model based on mi-gration velocity analysis (MVA) and basalt body building. Figures 3-5 display the smooth P-wavevelocity model at three different slicing locations. The MVA also reveals a high-velocity basaltbody at the center of the Soda Lake geothermal field, at the depth range of approximately 0.5 to1 km. All the other regions of the initial velocity model is very smooth and contain almost no anyhigh-wavenumber model perturbations that indicate either faults or sedimentary reflectors.We employ FWI to derive high-resolution subsurface medium property models for the SodaLake geothermal field. FWI is a nonlinear inversion approach to estimating medium properties4sing both the amplitude and phase information of the full seismic wavefield (Tarantola, 1984;Plessix, 2006; Virieux and Operto, 2009). In its simplest form, FWI is a L -norm nonlinear opti-mization problem: ψ ( m ) = min m || u − d || , (1)where m is the model parameter, d is the observed seismic waveform, and u is the syntheticseismic waveform. FWI is applicable in either acoustic or elastic media. In our seismic charac-terization, we use its acoustic form and invert for both the P-wave velocity V p and the density ρ ,i.e., m = ( V p , ρ ) . Inverting the density model in addition to the P-wave velocity model facilitatesthe FWI to achieve better amplitude match. Without including density inversion, some of the re-flections caused by acoustic impedance contrasts would be attributed entirely to velocity contrastsof the model, leading to “over-update” of the P-wave velocity model. In this case, visually, the P-wave velocity would contain fairly obvious high-wavenumber perturbations that are geologicallyless plausible.FWI in the form of equation (1) is highly nonlinear and difficult to converge because of cycle-skipping issue, particularly for noisy field seismic data. To alleviate these difficulties, we adopt anoptimized correlation misfit function in our seismic inversion (Choi and Alkhalifah, 2016): ψ ( m ) = min m ||C τ ( u , d ) − C τ ( d , d ) || , (2)where C is the cross-correlation operation, and τ is the time lag of the cross-correlation. Our studiesshow that this correlation-based FWI misfit function usually leads to much better convergence andreliable results compared with conventional FWI with equation (1).In FWI, the inversion gradients associated with the medium parameters are computed using theadjoint-state method (Plessix, 2006; Virieux and Operto, 2009). Our FWI is based on the followingfirst-order form of the acoustic-wave equation: ∂p∂t + K ∇ · v = f, (3) ρ ∂ v ∂t + ∇ p = 0 , (4)where K = ρV p is the bulk modulus of the medium, and f is the source term. We invert for V p and ρ simultaneously, therefore we compute the gradients associated V p and ρ based on the wave5quation system as ∇ V p ψ = − (cid:88) N s ,N r (cid:90) T max ρV p ∂p∂t p † dt, (5) ∇ ρ ψ = (cid:88) N s ,N r (cid:90) T max ∂ v ∂t v † dt, (6)where p = p ( x , t ) is the source pressure wavefield, p † is the adjoint-state pressure wavefield, v and v † are the source and adjoint particle velocity wavefields, respectively, N s and N r are the numbersof sources and receivers, respectively, and T max is the maximum propagation time of the wavefield.The adjoint-state wavefields p † and v † are the solutions to the adjoint-state wave equation withthe adjoint source being the external source term, solved in a reverse-time manner. The adjointsource is computed based on equation (2). One can refer Choi and Alkhalifah (2016) for thedetailed expression of the adjoint source term. High-resolution subsurface structural images can reveal faults that are crucial for optimizinggeothermal well placement and production. Conventional imaging such as ray-based Kirchhoffmigration can produce satisfactory images for simple structures, but usually fails to provide clearand reliable images when the subsurface geological structures are complex, particularly for faultinggeothermal fields.Reverse-time migration (RTM) is the industrial state-of-the-art imaging technique for imagingcomplex structures (McMechan, 1983; Chang and McMechan, 1987). RTM uses full wavefield toform subsurface images. In acoustic media, the PP image is formed by the zero time-lag cross-correlation between the source pressure wavefield and the receiver pressure wavefield. To properlyattenuate low-wavenumber artifacts caused by high medium contrasts, as at the Soda Lake geother-mal field containing a high-velocity basalt body, we employ the wavefield-separation-based RTMimaging condition (Fei et al., 2015): I pp = (cid:88) N s ,N r (cid:90) T max [ p s p r − H z ( p s ) H z ( p r ) − p s H z ( q r ) − H z ( p s ) q r ] dt, (7)where N s and N r are the numbers of source and receivers, respectively, T max is maximum wave-field propagation time, and p s = p s ( x , t ) and p r = p r ( x , t ) are the source and receiver wavefields,6espectively, with x being the spatial location and t being the time. The wavefield q r = q r ( x , t ) = g ( x ) ∗ H t ( d ( t )) is an auxiliary reverse-time-propagated wavefield with the temporal Hilbert trans-formed seismic data as the source term, g = g ( x ) is Green’s function, and d = d ( t ) is therecorded data at the receiver locations. H z denotes the Hilbert transform in the depth direction.The auxiliary dataset H t ( d ( t )) is computed prior to RTM imaging. Directly imaging faults with steep dips is a challenging task using seismic data (e.g., Tan andHuang, 2014), which requires fault-generated seismic scattering wavefields well preserved in data.For seismic characterization at the Soda Lake geothermal field, we alternatively use post-imagingautomatic fault detection method to delineate faults from the 3D seismic image volume producedusing the aforementioned RTM algorithm.We delineate the subsurface fault system at the Soda Lake geothermal field using the optimalfault surface voting method (Wu and Fomel, 2018). The method first automatically picks a set ofsparse seed points from an initial input fault attribute image, and then uses them to construct theoptimal surface patch based on global maximum attribute values. The method creates the finalfault attribute map, such as the fault likelihood map, from the smoothed attribute maps based oncollected accumulation scores. Finally, the method forms fault surfaces based on the computedfault strikes, dips and probabilities. One can refer to Wu and Fomel (2018) for algorithmic details.We then employ the detected faults to further enhance the RTM PP image with a fault-preserving,nonlinear anisotropic diffusion filtering method (Fehmers and H¨ocker, 2003; Wu and Guo, 2018).We process the PP image produced using RTM with the following nonlinear diffusion-type partialdifferential equation: ∂I∂t = ∇ ( ε D ∇ I ) , (8)where I = I ( x ) is the structural image, ε = ε ( x ) is the spatial coherence information, e.g.,the detected faults, and D = D ( x ) is the spatially heterogeneous anisotropic diffusion tensorconstructed from the structural tensor of the image I ( x ) . Equation (8) can effectively suppressrandom noises, improve lateral continuity, while enhance faults of the image. One can refer toFehmers and H¨ocker (2003) and Wu and Guo (2018) for algorithmic details.7 Results
The initial MVA P-wave velocity model built by Geokinetics as shown in Figures 3-5 satisfiessome kinematic properties of the seismic data. To accommodate our FWI of the 3D surface seismicdata from the Soda Lake geothermal field, we resample the initial velocity model with a 10 mregular grid spacing in three spatial directions. This resampling results in a regular-grid initialP-wave velocity model of 603 sample points in both the X- and Y-directions, and 255 samplingpoints in the depth direction. We build an initial density model using Gardner’s rule (Gardner et al.,1974) as ρ = 310 × V . p with a unit of kg/m .Three-dimensional FWI is a computationally intensive inversion problem. Considering thatthe survey contains over 8,300 common-shot gathers, and the data have a relatively high signal-to-noise ratio after processing, we limit the dominant frequency of inversion to 10 Hz, and we displaythe inversion results at the 13th iteration.Figures 6-8 display our FWI-inverted P-wave velocity model at three different slicing locationsand view angles. We find apparent layer-structured model perturbations in both the shallow anddeep regions compared with the initial smooth velocity model shown in Figures 3-5. Specifically,at the two slicing positions shown in Figures 7 and 8, we observe some faulting discontinuitiescutting through the layered structures across the entire model, even at positions that are not beneaththe high-contrasted basalt body.These layer-structured model perturbations and inferred faults are even more clear on our FWI-inverted density model shown in Figures 9-11. For example, in Figure 10, at inline receiver posi-tions spanning from approximately 1.5 km to 5 km, we observe some clear faulting discontinuitiesthat break the layers. Some of these faults are not beneath the basalt body, indicating that theyare not artifacts caused by the high medium property contrast of the basalt body. We find similarfaulting structures in Figure 11.The absolute values of the inverted density model might not be accurate out of two possi-ble reasons. First, the seismic data input for the inversion is not truly acoustic, but is PP dataseparated from the acquired elastic multi-component data. Therefore, their amplitude might notbe completely matched with the synthetic acoustic data. Second, we employ a cross-correlation8ype misfit function expressed in equation (2). This misfit function facilitates better convergencewhen using field seismic data, yet does not require a strict absolute amplitude match between thesynthetic and the observed waveforms. Nevertheless, this inverted “pseudo-density” model, asdescribed in our Methodology, prevents the inversion producing an over-updated velocity model.In this sense, we consider the inverted density model still has reasonable relative contrasts andaccurate structures, although it may not be accurate in terms of absolute values.With both the FWI-updated P-wave velocity and density models, we compute their correspond-ing acoustic impedance model shown in Figure 12, where we find clear discontinuities cuttingthrough the layers. These features can be reasonably interpreted as faults in this region. In addi-tion, we find clear layered structures in the acoustic impedance model, which generate reflectionsin the observed seismic data. These structures are obviously missing in the initial model shown inFigures 3-5. The spatial variations of the inverted models indicate the complexity of the subsurfacestructures at the Soda Lake geothermal field.Figure 13 depicts the convergence of the relative data misfit of our FWI over a total of 13iterations, showing an obvious misfit reduction even with the complexity of both the seismic dataand the geological structures in this area. Additional updates may further reduce the data misfit.In Figure 14, we compare among the observed data (Figure 14a), the synthetic data in the initialmodel (Figure 14b) and the synthetic data in the FWI-inverted model (Figure 14c), for a randomlyselected common-shot gather in the survey. The synthetic data in the initial model are obviouslymismatched with the observed data. Particularly, the reflection signals after approximately 1 s inthe observed data are completely missing in the synthetic data. By contrast, the match between thesynthetic data in the FWI-inverted model and the observed data is clearly improved. Seismic wave-forms in Figure 14c before and after 1 s closely resemble those in the observed data in Figure 14a.Other common-shot gathers have also a similar data match improvement after our FWI. We subsequently perform 3D reverse-time migration of the 3D surface seismic data from theSoda Lake geothermal field using the FWI-inverted velocity model. The center frequency of thesource wavelet used in RTM is also 10 Hz, the same as that used in our FWI.9igures 15-17 show the structural image up to 2.5 km in depth of the Soda Lake geothermalfield at three different slicing positions. We also plot the detected faults along with the images, withcolors of the faults representing fault probability. The horizontal slices on the top-left panels ofFigures 15-17 show that major faults in this area are along the diagonal direction of the geometry.Interpreted based on the map shown in Figure 1, these major faults are approximately along thenorth-south direction, with a small to moderate azimuth angle. On the vertical slices, we find thelayers in this region are well imaged, even blow the high-contrasted basalt body. Faults in thisregion have steeply dipping angles, as indicated by the fault images in Figures 15-17. Some of thethese faults penetrate through the basalt body, indicating that the faulting in this region occurredafter the basalt body was formed in the geological history.Figure 18 shows the detected faults from the structural image shown in Figure 15 at threedifferent view angles. Consistent with those shown in Figures 15-17, the colors of the fault surfacesrepresent the fault probability. Although the fault probability is not full for every spatial point onthe fault surfaces, most fault surfaces have moderate to high probability, indicating the detectedfaults are reliable.There are several currently active injection and production wells at the Soda Lake geothermalfield. To validate the accuracy and reliability of our imaging and fault detection results, we vi-sualize the wells, the structural images, and the detected faults in their true 3D spatial positions,as displayed in Figure 19. Figures 19a and b show several image slices superimposed with thedetected faults (in white-blue colors). We place the injection wells (green-colored tubes) and pro-duction wells (red-colored tubes) in the 3D space. We find that all the injection and productionwells either run through the detected faults or are fairly close to one or two faults. For instance, inFigure 19a, a curved production well in the center of the image penetrates exactly a major fault onthe image. The production well on the right of the figure also penetrates a location where two faultsintersect. Figure 19b shows that that the green-colored injection well on the right of the figure islocated adjacent to a major fault, and part of the injection well overlies the fault path.We show the constructed fault surfaces in Figures 19c-f, along with injection and productionwells plotted in the same 3D scene. Clearly, all the production wells either penetrate one or twofaults, or are very close to a fault surface. Some of the production wells penetrate through a samefault twice. The consistency between the currently active injection and production wells and our10onstructed faults from our 3D image volume demonstrates that our imaging and fault detectionresults are very close to the realistic geology in this region. The detected fault surfaces consist ofa complex fault system in this region.Our preceding results in this paper reveal complex subsurface structures and medium propertyvariations up to 2.5 km at the Soda Lake geothermal field, with a spatial grid interval of 10 malong lateral and vertical directions. We aim to reveal the complex fault system with a higherspatial resolution using a grid spacing of 6.7 m in the horizontal directions and 2.5 m in the depthdirection, up to 1 km in depth for this area. We use a source wavelet with a center frequency of30 Hz for RTM in this fine grid, enabling us to resolve fine layers for this geothermal field, and toobtain a high-resolution fault system.Figures 20-22 show the image volume with a horizontal grid spacing of 6.7 m and a vertical gridspacing of 2.5 m, superimposed by the detected faults from this image volume. We observe thatthe detected faults have similar strikes with those shown in Figures 15-17. These high-resolutionimages and detected faults further verify the geological complexity of the near-surface region upto 1 km in depth.We further construct fault surfaces using the high-resolution image volume shown in Fig-ures 20-22. Figure 23 show the constructed fault surfaces at three different view angles. Thesefault surfaces have relatively higher overall fault probability compared with those associated withthe 10-m-grid-spacing image shown in Figure 18. Some of these faults are not properly detectedfrom the 10-m-grid-spacing image. We also find interlacing fault surfaces in Figure 23, which fur-ther verify that the fault system at the Soda Lake geothermal field is fairly complex, and requireshigh-resolution images to reveal those faults.Similar to the examination on the fault-well consistency for the 10-m-grid-spacing image, weshow in Figure 24 that all currently active injection and production wells either run through or arevery close to our detected faults. This fault-well consistency, along with that associated with the10-m-grid-spacing image volume, further verifies the accuracy and reliability of our subsurfacefault imaging at the Soda Lake geothermal field.11
Conclusions
We have conducted a 3D, high-resolution seismic characterization for the Soda Lake geother-mal field using full-waveform inversion and reverse-time migration of 3D surface seismic data. Wehave obtained updated velocity, density, and acoustic impedance models using 3D full-waveforminversion of the pre-processed 3D PP seismic data. We have also performed high-resolution,wavefield-separation-based reverse-time migration to obtain high-resolution 3D structural images,including a 10-m-grid-spacing image up to 2.5 km in depth, and a high-resolution image with avertical grid spacing of 2.5 m up to 1 km. We have detected faults from these image volumesand constructed corresponding fault surfaces, revealing the complex fault system at the Soda Lakegeothermal field. A careful check on the consistency between the constructed fault surfaces andcurrent active injection and production wells validate that our seismic inversion and imaging re-sults and detect faults are accurate and reliable. These results can provide valuable information foroptimizing well placement and geothermal energy production. Future work aims at using multi-component elastic seismic data to conduct isotropic and anisotropic elastic full-waveform inversionto reveal anisotropic characteristics of the Soda Lake geothermal field.
This work was supported by the U.S. Department of Energy (DOE) Geothermal TechnologiesOffice through the Los Alamos National Laboratory (LANL). LANL is operated by Triad NationalSecurity, LLC, for the U.S. DOE National Nuclear Security Administration (NNSA) under Con-tract No. 89233218CNA000001. This research used resources provided by the LANL InstitutionalComputing Program supported by the U.S. DOE NNSA under Contract No. 89233218CNA000001.We thank Xinming Wu and Benxin Chi for helpful discussions, and David Li for his review of themanuscript. 12 eferences
Chang, W.-F., and G. A. McMechan, 1987, Elastic reverse-time migration: Geophysics, , no.10, 1365–1375, doi: 10.1190/1.1442249.Choi, Y., and T. Alkhalifah, 2016, An optimized correlation-based full waveform inversion:Conference Proceedings of 78th EAGE Conference and Exhibition, , no. 1, 1–5, doi:10.3997/2214-4609.201600642.Fehmers, G. C., and C. F. W. H ¨ocker, 2003, Fast structural interpretation with structure-orientedfiltering: Geophysics, , no. 4, 1286–1293, doi: 10.1190/1.1598121.Fei, T. W., Y. Luo, J. Yang, H. Liu, and F. Qin, 2015, Removing false images in reverse time mi-gration: The concept of de-primary: Geophysics, , no. 6, S237–S244, doi: 10.1190/geo2015-0289.1.Gardner, G. H. F., L. W. Gardner, and A. R. Gregory, 1974, Formation velocity and density: The di-agnostic basics for stratigraphic traps: Geophysics, , no. 6, 770–780, doi: 10.1190/1.1440465.McMechan, G. A., 1983, Migration by extrapolation of time-dependent boundary values: Geo-physical Prospecting, , no. 3, 413–420, doi: 10.1111/j.1365-2478.1983.tb01060.x.Plessix, R.-E., 2006, A review of the adjoint-state method for computing the gradient of a func-tional with geophysical applications: Geophysical Journal International, , no. 2, 495–503,doi: 10.1111/j.1365-246X.2006.02978.x.Tan, S., and L. Huang, 2014, Least-squares reverse-time migration with a wavefield-separationimaging condition and updated source wavefields: Geophysics, , no. 5, S195–S205, doi:10.1190/geo2014-0020.1.Tarantola, A., 1984, Inversion of seismic reflection data in the acoustic approximation: Geo-physics, , no. 8, 1259–1266, doi: 10.1190/1.1441754.Virieux, J., and S. Operto, 2009, An overview of full-waveform inversion in exploration geo-physics: Geophysics, , no. 6, WCC1–WCC26, doi: 10.1190/1.3238367.Wu, X., and S. Fomel, 2018, Automatic fault interpretation with optimal surface voting: Geo-physics, , no. 5, O67–O82, doi: 10.1190/geo2018-0115.1.Wu, X., and Z. Guo, 2018, Detecting faults and channels while enhancing seismic structural andstratigraphic features: Interpretation, , no. 1, T155–T166, doi: 10.1190/INT-2017-0174.1.13 ist of Figures × km . . . . . 173 Slices and a 3D view of the initial P-wave velocity model produced by Geokineticsat slicing position 1. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 184 Slices and a 3D view of the initial P-wave velocity model produced by Geokineticsat slicing position 2. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 195 Slices and a 3D view of the initial P-wave velocity model produced by Geokineticsat slicing position 3. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 206 Slices and a 3D view of the updated P-wave velocity model produced using ourFWI at slicing position 1 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 217 Slices and a 3D view of the updated P-wave velocity model produced using ourFWI at slicing position 2. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 228 Slices and a 3D view of the updated P-wave velocity model produced using ourFWI at slicing position 3. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 239 Slices and a 3D view of the updated density model produced using our FWI atslicing position 1. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2410 Slices and a 3D view of the updated density model produced using our FWI atslicing position 2. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2511 Slices and a 3D view of the updated density model produced using our FWI atslicing position 3. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2612 Slices and a 3D view of the updated acoustic impedance model produced using ourFWI. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2713 Data misfit convergence over a total of 13 iterations in our FWI. . . . . . . . . . . 2814 Panels (a)-(c) show the observed data, the synthetic data in the initial model shownin Figure 3, and the synthetic data in the FWI-updated model shown in Figure 6,respectively. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2915 Slices and a 3D view of the structural image of the Soda Lake geothermal field upto 2.5 km in depth produced with our fault-enhancing RTM algorithm, superim-posed with the detected faults at slicing position 1. . . . . . . . . . . . . . . . . . 3016 Slices and a 3D view of the structural image of the Soda Lake geothermal field upto 2.5 km in depth produced with our fault-enhancing RTM algorithm, superim-posed with the detected faults at slicing position 2. . . . . . . . . . . . . . . . . . 3117 Slices and a 3D view of the structural image of the Soda Lake geothermal field upto 2.5 km in depth produced with our fault-enhancing RTM algorithm, superim-posed with the detected faults at slicing position 3. . . . . . . . . . . . . . . . . . 3218 Panels (a)-(c) show the constructed fault surfaces from the RTM structural imageshown in Figure 15 at three different view angles. Colors of the fault surfacesrepresent the fault probability and are consistent with those in Figures 15-17. . . . 33149 Panels (a)-(b) show the detected faults for the RTM structural image shown inFigure 15 at two different view angles, along with active injection (green) andproduction (red) wells at the Soda Lake geothermal field. Panels (c)-(f) show theconstructed fault surfaces with currently active injection and production wells inthis geothermal field. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3420 Slices and a 3D view of the structural image up to 1 km in depth produced withour fault-enhancing RTM algorithm, superimposed by the detected faults at slicingposition 1. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3521 Slices and a 3D view of the structural image up to 1 km in depth produced withour fault-enhancing RTM algorithm, superimposed by the detected faults at slicingposition 2. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3622 Slices and a 3D view of the structural image up to 1 km in depth produced withour fault-enhancing RTM algorithm, superimposed by the detected faults at slicingposition 3. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3723 Panels (a)-(c) show the detected faults from the structural image volume shown inFigure 20 at three different view angles. Colors of the fault surfaces represent faultprobability and are consistent with those in Figures 20-22. . . . . . . . . . . . . . 3824 Panels (a)-(b) show the constructed fault surfaces from the structural image shownin Figure 20 at two different view angles, along with active injection (green) andproduction (red) wells in this area. Panels (c)-(f) depict the constructed fault sur-faces with currently active injection and production wells at the Soda Lake geother-mal field. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3915igure 1: Location of the Soda Lake geothermal field in the western Nevada.16igure 2: Source and receiver distributions of the 3D surface seismic survey at the Soda Lakegeothermal field. The survey contains 8,321 common-shot gathers. Each common-shot gather hasa coverage area of up to approximately × km .17igure 3: Slices and a 3D view of the initial P-wave velocity model produced by Geokinetics atslicing position 1. 18igure 4: Slices and a 3D view of the initial P-wave velocity model produced by Geokinetics atslicing position 2. 19igure 5: Slices and a 3D view of the initial P-wave velocity model produced by Geokinetics atslicing position 3. 20igure 6: Slices and a 3D view of the updated P-wave velocity model produced using our FWI atslicing position 1 . 21igure 7: Slices and a 3D view of the updated P-wave velocity model produced using our FWI atslicing position 2. 22igure 8: Slices and a 3D view of the updated P-wave velocity model produced using our FWI atslicing position 3. 23igure 9: Slices and a 3D view of the updated density model produced using our FWI at slicingposition 1. 24igure 10: Slices and a 3D view of the updated density model produced using our FWI at slicingposition 2. 25igure 11: Slices and a 3D view of the updated density model produced using our FWI at slicingposition 3. 26igure 12: Slices and a 3D view of the updated acoustic impedance model produced using ourFWI. 27igure 13: Data misfit convergence over a total of 13 iterations in our FWI.28 a)(b)(c) Figure 14: Panels (a)-(c) show the observed data, the synthetic data in the initial model shown inFigure 3, and the synthetic data in the FWI-updated model shown in Figure 6, respectively.29igure 15: Slices and a 3D view of the structural image of the Soda Lake geothermal field upto 2.5 km in depth produced with our fault-enhancing RTM algorithm, superimposed with thedetected faults at slicing position 1. 30igure 16: Slices and a 3D view of the structural image of the Soda Lake geothermal field upto 2.5 km in depth produced with our fault-enhancing RTM algorithm, superimposed with thedetected faults at slicing position 2. 31igure 17: Slices and a 3D view of the structural image of the Soda Lake geothermal field upto 2.5 km in depth produced with our fault-enhancing RTM algorithm, superimposed with thedetected faults at slicing position 3. 32 a)(b)(c)
Figure 18: Panels (a)-(c) show the constructed fault surfaces from the RTM structural image shownin Figure 15 at three different view angles. Colors of the fault surfaces represent the fault proba-bility and are consistent with those in Figures 15-17.33 a) (b)(c) (d)(e) (f)
Figure 19: Panels (a)-(b) show the detected faults for the RTM structural image shown in Figure 15at two different view angles, along with active injection (green) and production (red) wells at theSoda Lake geothermal field. Panels (c)-(f) show the constructed fault surfaces with currently activeinjection and production wells in this geothermal field.34igure 20: Slices and a 3D view of the structural image up to 1 km in depth produced with ourfault-enhancing RTM algorithm, superimposed by the detected faults at slicing position 1.35igure 21: Slices and a 3D view of the structural image up to 1 km in depth produced with ourfault-enhancing RTM algorithm, superimposed by the detected faults at slicing position 2.36igure 22: Slices and a 3D view of the structural image up to 1 km in depth produced with ourfault-enhancing RTM algorithm, superimposed by the detected faults at slicing position 3.37 a)(b)(c)
Figure 23: Panels (a)-(c) show the detected faults from the structural image volume shown inFigure 20 at three different view angles. Colors of the fault surfaces represent fault probability andare consistent with those in Figures 20-22. 38 a) (b)(c) (d)(e) (f)a) (b)(c) (d)(e) (f)