Background In the last two annual reports on CO2 storage monitoring, we discussed the feasibility... more Background In the last two annual reports on CO2 storage monitoring, we discussed the feasibility of using different geophysical options for subsurface CO2 monitoring, and described an innovative approach for near-continuous subsurface monitoring. We also discussed special challenges and requirements with continuous and leak detection of CO2. We chose seismic methods for subsurface monitoring primarily because they are broadly applicable in a variety of geological settings and sensitive to changes in CO2 saturation. To further develop seismic methods, for this particular application, we have proposed adaptive sparse data acquisition to address the problem of cost-effective time-lapse monitoring. We focused our work on leak detection and have proposed data gathering and processing concepts for dynamic imaging and model-based inversion. Developments reported in previous reports included a suite of seismic modeling tools to simulate monitoring and a small-scale field example on seismic...
Abstract We present preliminary results from an on-going geophysical investigation of the former ... more Abstract We present preliminary results from an on-going geophysical investigation of the former DOE Pinellas site in Pinellas County, Florida, a site with confirmed NAPL contamination. The goal of this work was to demonstrate the combined use of high-resolution crosswell seismic and crosswell radar data for characterization of NAPL distribution in near-surface environments. Borehole conductivity, gamma, and cone penetrometry measurements provided secondary information to constrain lithology. ...
We describe a novel framework for PDE (partial-differential-equation)-constrained full-waveform i... more We describe a novel framework for PDE (partial-differential-equation)-constrained full-waveform inversion (FWI) that estimates parameters of subsurface flow processes, such as rock permeability and porosity, using time-lapse observed data. The forward modeling couples flow physics, rock physics, and wave physics models. For the inverse modeling, we handle the back-propagation of gradients by an intelligent automatic differentiation strategy that offers three levels of user control: (1) At the wave physics level, we adopt the discrete adjoint method in order to use our existing high-performance FWI code; (2) at the rock physics level, we use built-in operators from the $\texttt{TensorFlow}$ backend; (3) at the flow physics level, we code customized PDE operators for the potential and nonlinear saturation equations, which highly resemble neural network layers. These three levels of gradient computation strike a good balance between computational efficiency and programming efficiency, and when chained together constitute the coupled inverse system. We use numerical experiments to demonstrate: (1) the three-level coupled inverse problem is superior in terms of accuracy to a traditional strategy; (2) it is able to simultaneously invert for parameters in empirical relationships such as the rock physics models; (3) and, the inverted model can be used for reservoir performance prediction and reservoir management/optimization purposes.
We use the first arrival traveltime to correct the phase distortion in a nonlinear wave equation ... more We use the first arrival traveltime to correct the phase distortion in a nonlinear wave equation inversion scheme.This improves the precision of tomographic reconstruction of a velocity structure with large variations and helps solve the ill-posed problem of wave equation inversion.When the variation of the velocity distribution is large, general non-linear wave equation inversions are very ill-posed and for such strong nonlinear we can not obtain a correct inversion.One of main reasons is that the calculated and observed phase of the wavefield differs greatly if the initial model is far from the true model.This leads to highly mismatched phase between the calculated and the observed wave field.This is so-called “Cycle Skipping” problem in the full waveform inversion.The phase mismatch is even more pronounced if a high operating frequency is employed in order to increase resolution.To address this problem, we use the first arrival to “demodulate” the wave field in the frequency domain with a goal of restoring the phase of wave field.Then we minimize an objective function consisting of so called “demodulated” wave field to solve wave equation inversion problem.In this way, we find that the inversion is much improved, and when the velocity perturbation in a complicated model reaches 35%, we can still obtain a good inversion.A computer simulation shows that our method is very robust for acoustical wave inversion with good reconstruction precision.
Q-compensated imaging combines attenuation (Q-1) estimation with attenuation-compensated reverse ... more Q-compensated imaging combines attenuation (Q-1) estimation with attenuation-compensated reverse time migration (Q-RTM). We review the approach (theory and algorithm), and present results from a field application of Q-compensated imaging using borehole seismic data. The goal of Q-compensated imaging is to restore pre-stack amplitude loss, improve signal bandwidth, and improve pre-stack balance of amplitude and phase in an effort to improve stack coherency and subsequently image resolution. To this end, we developed and tested a consistent approach for Q estimation and forward-time and reverse-time simulation for waves in constant Q media. This approach incorporates robust and practical algorithms for Q-compensated imaging. The approach is tested on a crosswell field dataset from west Texas, USA where improved spatial resolution and improved spatial continuity can be observed in the Q-compensated image.
Sparse coding can be applied to train an overcomplete dictionary on time-lapse seismic data or im... more Sparse coding can be applied to train an overcomplete dictionary on time-lapse seismic data or images. The learned dictionary generally consists of sparse representations of one or more images. We then use such sparse representations, along with L 1-regularization techniques, to predict missing values in seismic images by solving an inverse problem. The practical outcome of the proposed methodology can be a significant reduction in field operational costs by requiring only sparse instead of dense surveys, and by integrating in the seismic images the information captured by the learned dictionary from previous time-lapse and baseline images. A synthetic example is presented to test the method.
In signal processing we often perform time shifts. Considering the duality of time domain and fre... more In signal processing we often perform time shifts. Considering the duality of time domain and frequency domain representations, we may also want to shift the spectrum of a signal. The frequency shift transform is such a process. It is equivalent to complex modulation or demodulation of a signal (Oppenheim, 1987). The FST applied to a spectrum S(f) is represented by S(f+f& where fOcan be positive or negative. Here, we should assume that S(f) is bandlimited so that we have free space along the frequency axis to shift the spectrum (see Figures la and lb). time domain signals with the same bandwidth corresponding to S(f) and S(f+fo) exhibit the same envelope, but different oscillations exist under the envelope. It is well known that the peak to sidelobe level of a correlated sweep depends on the bandwidth of the signal. In practical terms, the required bandwidth is often stated in octaves, e.g., 2-3 octaves are needed for a “good” correlation. It is not so well known that the envelope of the correlated sweep depends only on the absolute bandwidth Af and not the number of octaves of bandwidth. Therefore, the envelope of a sweep from 1000 Hz to 1500 Hz (less than one octave) has the same envelop as a sweep from 5 Hz to 505 Hz (over six octaves). The fact that the high frequency sweep is narrowband whereas the low frequency signal is wideband can be either an advantage or a disadvantage depending on the application. One direct application of FST is to extract the envelope of a signal. In this case FST acts as demodulation. Another possible usage is to solve some aliasing problems. In this case we reduce the dominant frequency of a signal by shifting the spectrum to lower frequencies. Then we can use a larger interval when sampling this signal. This feature is specially useful for spatial aliasing because it is expensive or difficult to increase the spatial sampling rate. The FST is a linear time-variant transform, Therefore, when applying FST to signals with different timedelays, we find that these signals after FST exhibit phase distortion. Special efforts must be _ made to remove this distortion. It should be pointed out that the complex FST is reversible: when necessary we can take a reversed FST that restores the signal to the original style.
An approach for quasi-continuous, geophysical time-lapse monitoring with sparse seismic data is p... more An approach for quasi-continuous, geophysical time-lapse monitoring with sparse seismic data is proposed. This approach takes advantage of the small changes in the seismic property of a geological reservoir that are expected to occur in a small time interval. The goal of this approach is to obtain high temporal and spatial resolution in reconstructed, time-lapse geophysical images using comparable resources that would have provided high spatial but low temporal resolution images with conventional approaches. This is done by acquiring spatially sparse data at small time intervals. In this case, a spatially sparse dataset refers to that dataset which is a small fraction (as little as 5%) of what would be acquired to reconstruct a high spatial resolution tomographic image of the subsurface. The high spatial resolution obtained by the proposed approach occurs because unrecorded data are predicted from future and past data. With high temporal and spatial resolution, early detection of important reservoir changes is more likely to occur.
A stochastic approach to seismic inversion using the ensemble Kalman filter (EnKF) is proposed. S... more A stochastic approach to seismic inversion using the ensemble Kalman filter (EnKF) is proposed. Seismic depth and time image data are used as the input for EnKF stochastic seismic inversion. The sonic log is used to estimate source wavelet and create initial models for the inversion, which provides an efficient integration of sonic log data and seismic data. We use both travel time and waveform data for the inversion and obtain the absolute seismic velocity instead of the relative impedance. EnKF can continuously update the model using time-lapse data. A synthetic example is used to demonstrate the possible application to seismic monitoring.
We present the results of five experiments investigating the dielectric properties of granular ma... more We present the results of five experiments investigating the dielectric properties of granular materials partially saturated with trichloroethylene (TCE), a common dense non-aqueous contaminant. Previous research has investigated the radar signatures of similar solvents in controlled field experiments but no core-scale measurements have verified the appropriate petrophysical model. Broadband dielectric measurements were performed using a time domain reflectometry (TDR) system coupled to a solvent-compatible coaxial transmission line. Two synthetic samples and three natural aquifer samples were fully saturated with water and then subjected to an axial TCE injection until breakthrough was observed. The resulting dielectric measurements show good agreement with the empirical complex refractive index model (CRIM) allowing a reasonable prediction of the radar reflectivities and transmission velocities expected in field surveys targeting pools of similar non-aqueous contaminants.
Subsurface imaging with common‐source cross‐hole data can be achieved using prestack reverse‐time... more Subsurface imaging with common‐source cross‐hole data can be achieved using prestack reverse‐time migration. The algorithm consists of extrapolation of the recorded wave field, application of the excitation‐time imaging condition, and postprocessing of the resulting image with a low‐pass wavenumber filter. The wavenumber filter removes the artifact associated with the direct arrival; this artifact is not separable from the scattered data before migration because, in the cross‐hole geometry, they significantly overlap in time, space, and wavenumber. Migration of synthetic data produces the best possible results, but images produced by migration of scale‐model data are not greatly inferior. Apparently, acceptable images can be obtained from a surprisingly few sources, if these sources are located sufficiently far apart to give independent information and the recording aperture is sufficiently wide.
IEEE Transactions on Geoscience and Remote Sensing, Jul 1, 1987
We present both theory and numerical simulations of diffraction tomography for arrays of line sou... more We present both theory and numerical simulations of diffraction tomography for arrays of line sources. The approach taken is applicable to a wide range of imaging problems where discrete sources and receivers are located near the object to be imaged rather than in its far field. As such, these new results tie into the vast body of research on inverse scattering, which to date is based largely on planewave sources. Our derivation implicitly includes a method for synthesizing plane waves; this method leads to inversion formulas based on the generalized projection-slice theorem. Although related techniques for synthesizing plane waves from discrete arrays are known, it is helpful to have available a theory that incorporates the synthesis directly into the scattering theory. The formulation is presented for propagating fields satisfying the Born and Rytov approximations in weakly inhomogeneous media and provides a convenient means for treating both forward and inverse scattering problems. Through a numerical example, we illustrate two important features of diffraction tomography inversion, i.e., 1) the effects of limited view and 2) results of probing with different signal frequencies. The example utilizes data generated by an exact forward modeling technique thus providing strong evidence supporting the usefulness of weak scattering approximations for inversion problems.
A crosswell dataset, collected at Conoco' s Borehole Test Facility in Oklahoma, was first process... more A crosswell dataset, collected at Conoco' s Borehole Test Facility in Oklahoma, was first processed using a tomography algorithm based on an isotropic velocity model. The resulting 2-D velocity tomogram was found to have huge artifacts believed to be caused by nonuniform ray coverage and possibly anisotropy. A comparison of the zero offset crosswell velocities (Horizontal Path Log) with the log-derived sonic velocities (Vertical Path Log) indicated nearly 25% P-wave anisotropy in some formations. Despite the artifacts, a good interpretation of the 2-D isotropic tomogram was made with the help of synthetic data. The data were then inverted using an algorithm that incorporates a model of elliptical transverse isotropy. This anisotropy model produced an unambiguous image for two components of velocity that simultaneously match the sonic log and the crosswell data. Both inversions support a final interpretation of 1-D vertical stratification of layers, some of which exhibit significant P-wave anisotropy. This interpretation is consistent with the sonic logs, crosswell data, and available geological information.
Background In the last two annual reports on CO2 storage monitoring, we discussed the feasibility... more Background In the last two annual reports on CO2 storage monitoring, we discussed the feasibility of using different geophysical options for subsurface CO2 monitoring, and described an innovative approach for near-continuous subsurface monitoring. We also discussed special challenges and requirements with continuous and leak detection of CO2. We chose seismic methods for subsurface monitoring primarily because they are broadly applicable in a variety of geological settings and sensitive to changes in CO2 saturation. To further develop seismic methods, for this particular application, we have proposed adaptive sparse data acquisition to address the problem of cost-effective time-lapse monitoring. We focused our work on leak detection and have proposed data gathering and processing concepts for dynamic imaging and model-based inversion. Developments reported in previous reports included a suite of seismic modeling tools to simulate monitoring and a small-scale field example on seismic...
Abstract We present preliminary results from an on-going geophysical investigation of the former ... more Abstract We present preliminary results from an on-going geophysical investigation of the former DOE Pinellas site in Pinellas County, Florida, a site with confirmed NAPL contamination. The goal of this work was to demonstrate the combined use of high-resolution crosswell seismic and crosswell radar data for characterization of NAPL distribution in near-surface environments. Borehole conductivity, gamma, and cone penetrometry measurements provided secondary information to constrain lithology. ...
We describe a novel framework for PDE (partial-differential-equation)-constrained full-waveform i... more We describe a novel framework for PDE (partial-differential-equation)-constrained full-waveform inversion (FWI) that estimates parameters of subsurface flow processes, such as rock permeability and porosity, using time-lapse observed data. The forward modeling couples flow physics, rock physics, and wave physics models. For the inverse modeling, we handle the back-propagation of gradients by an intelligent automatic differentiation strategy that offers three levels of user control: (1) At the wave physics level, we adopt the discrete adjoint method in order to use our existing high-performance FWI code; (2) at the rock physics level, we use built-in operators from the $\texttt{TensorFlow}$ backend; (3) at the flow physics level, we code customized PDE operators for the potential and nonlinear saturation equations, which highly resemble neural network layers. These three levels of gradient computation strike a good balance between computational efficiency and programming efficiency, and when chained together constitute the coupled inverse system. We use numerical experiments to demonstrate: (1) the three-level coupled inverse problem is superior in terms of accuracy to a traditional strategy; (2) it is able to simultaneously invert for parameters in empirical relationships such as the rock physics models; (3) and, the inverted model can be used for reservoir performance prediction and reservoir management/optimization purposes.
We use the first arrival traveltime to correct the phase distortion in a nonlinear wave equation ... more We use the first arrival traveltime to correct the phase distortion in a nonlinear wave equation inversion scheme.This improves the precision of tomographic reconstruction of a velocity structure with large variations and helps solve the ill-posed problem of wave equation inversion.When the variation of the velocity distribution is large, general non-linear wave equation inversions are very ill-posed and for such strong nonlinear we can not obtain a correct inversion.One of main reasons is that the calculated and observed phase of the wavefield differs greatly if the initial model is far from the true model.This leads to highly mismatched phase between the calculated and the observed wave field.This is so-called “Cycle Skipping” problem in the full waveform inversion.The phase mismatch is even more pronounced if a high operating frequency is employed in order to increase resolution.To address this problem, we use the first arrival to “demodulate” the wave field in the frequency domain with a goal of restoring the phase of wave field.Then we minimize an objective function consisting of so called “demodulated” wave field to solve wave equation inversion problem.In this way, we find that the inversion is much improved, and when the velocity perturbation in a complicated model reaches 35%, we can still obtain a good inversion.A computer simulation shows that our method is very robust for acoustical wave inversion with good reconstruction precision.
Q-compensated imaging combines attenuation (Q-1) estimation with attenuation-compensated reverse ... more Q-compensated imaging combines attenuation (Q-1) estimation with attenuation-compensated reverse time migration (Q-RTM). We review the approach (theory and algorithm), and present results from a field application of Q-compensated imaging using borehole seismic data. The goal of Q-compensated imaging is to restore pre-stack amplitude loss, improve signal bandwidth, and improve pre-stack balance of amplitude and phase in an effort to improve stack coherency and subsequently image resolution. To this end, we developed and tested a consistent approach for Q estimation and forward-time and reverse-time simulation for waves in constant Q media. This approach incorporates robust and practical algorithms for Q-compensated imaging. The approach is tested on a crosswell field dataset from west Texas, USA where improved spatial resolution and improved spatial continuity can be observed in the Q-compensated image.
Sparse coding can be applied to train an overcomplete dictionary on time-lapse seismic data or im... more Sparse coding can be applied to train an overcomplete dictionary on time-lapse seismic data or images. The learned dictionary generally consists of sparse representations of one or more images. We then use such sparse representations, along with L 1-regularization techniques, to predict missing values in seismic images by solving an inverse problem. The practical outcome of the proposed methodology can be a significant reduction in field operational costs by requiring only sparse instead of dense surveys, and by integrating in the seismic images the information captured by the learned dictionary from previous time-lapse and baseline images. A synthetic example is presented to test the method.
In signal processing we often perform time shifts. Considering the duality of time domain and fre... more In signal processing we often perform time shifts. Considering the duality of time domain and frequency domain representations, we may also want to shift the spectrum of a signal. The frequency shift transform is such a process. It is equivalent to complex modulation or demodulation of a signal (Oppenheim, 1987). The FST applied to a spectrum S(f) is represented by S(f+f& where fOcan be positive or negative. Here, we should assume that S(f) is bandlimited so that we have free space along the frequency axis to shift the spectrum (see Figures la and lb). time domain signals with the same bandwidth corresponding to S(f) and S(f+fo) exhibit the same envelope, but different oscillations exist under the envelope. It is well known that the peak to sidelobe level of a correlated sweep depends on the bandwidth of the signal. In practical terms, the required bandwidth is often stated in octaves, e.g., 2-3 octaves are needed for a “good” correlation. It is not so well known that the envelope of the correlated sweep depends only on the absolute bandwidth Af and not the number of octaves of bandwidth. Therefore, the envelope of a sweep from 1000 Hz to 1500 Hz (less than one octave) has the same envelop as a sweep from 5 Hz to 505 Hz (over six octaves). The fact that the high frequency sweep is narrowband whereas the low frequency signal is wideband can be either an advantage or a disadvantage depending on the application. One direct application of FST is to extract the envelope of a signal. In this case FST acts as demodulation. Another possible usage is to solve some aliasing problems. In this case we reduce the dominant frequency of a signal by shifting the spectrum to lower frequencies. Then we can use a larger interval when sampling this signal. This feature is specially useful for spatial aliasing because it is expensive or difficult to increase the spatial sampling rate. The FST is a linear time-variant transform, Therefore, when applying FST to signals with different timedelays, we find that these signals after FST exhibit phase distortion. Special efforts must be _ made to remove this distortion. It should be pointed out that the complex FST is reversible: when necessary we can take a reversed FST that restores the signal to the original style.
An approach for quasi-continuous, geophysical time-lapse monitoring with sparse seismic data is p... more An approach for quasi-continuous, geophysical time-lapse monitoring with sparse seismic data is proposed. This approach takes advantage of the small changes in the seismic property of a geological reservoir that are expected to occur in a small time interval. The goal of this approach is to obtain high temporal and spatial resolution in reconstructed, time-lapse geophysical images using comparable resources that would have provided high spatial but low temporal resolution images with conventional approaches. This is done by acquiring spatially sparse data at small time intervals. In this case, a spatially sparse dataset refers to that dataset which is a small fraction (as little as 5%) of what would be acquired to reconstruct a high spatial resolution tomographic image of the subsurface. The high spatial resolution obtained by the proposed approach occurs because unrecorded data are predicted from future and past data. With high temporal and spatial resolution, early detection of important reservoir changes is more likely to occur.
A stochastic approach to seismic inversion using the ensemble Kalman filter (EnKF) is proposed. S... more A stochastic approach to seismic inversion using the ensemble Kalman filter (EnKF) is proposed. Seismic depth and time image data are used as the input for EnKF stochastic seismic inversion. The sonic log is used to estimate source wavelet and create initial models for the inversion, which provides an efficient integration of sonic log data and seismic data. We use both travel time and waveform data for the inversion and obtain the absolute seismic velocity instead of the relative impedance. EnKF can continuously update the model using time-lapse data. A synthetic example is used to demonstrate the possible application to seismic monitoring.
We present the results of five experiments investigating the dielectric properties of granular ma... more We present the results of five experiments investigating the dielectric properties of granular materials partially saturated with trichloroethylene (TCE), a common dense non-aqueous contaminant. Previous research has investigated the radar signatures of similar solvents in controlled field experiments but no core-scale measurements have verified the appropriate petrophysical model. Broadband dielectric measurements were performed using a time domain reflectometry (TDR) system coupled to a solvent-compatible coaxial transmission line. Two synthetic samples and three natural aquifer samples were fully saturated with water and then subjected to an axial TCE injection until breakthrough was observed. The resulting dielectric measurements show good agreement with the empirical complex refractive index model (CRIM) allowing a reasonable prediction of the radar reflectivities and transmission velocities expected in field surveys targeting pools of similar non-aqueous contaminants.
Subsurface imaging with common‐source cross‐hole data can be achieved using prestack reverse‐time... more Subsurface imaging with common‐source cross‐hole data can be achieved using prestack reverse‐time migration. The algorithm consists of extrapolation of the recorded wave field, application of the excitation‐time imaging condition, and postprocessing of the resulting image with a low‐pass wavenumber filter. The wavenumber filter removes the artifact associated with the direct arrival; this artifact is not separable from the scattered data before migration because, in the cross‐hole geometry, they significantly overlap in time, space, and wavenumber. Migration of synthetic data produces the best possible results, but images produced by migration of scale‐model data are not greatly inferior. Apparently, acceptable images can be obtained from a surprisingly few sources, if these sources are located sufficiently far apart to give independent information and the recording aperture is sufficiently wide.
IEEE Transactions on Geoscience and Remote Sensing, Jul 1, 1987
We present both theory and numerical simulations of diffraction tomography for arrays of line sou... more We present both theory and numerical simulations of diffraction tomography for arrays of line sources. The approach taken is applicable to a wide range of imaging problems where discrete sources and receivers are located near the object to be imaged rather than in its far field. As such, these new results tie into the vast body of research on inverse scattering, which to date is based largely on planewave sources. Our derivation implicitly includes a method for synthesizing plane waves; this method leads to inversion formulas based on the generalized projection-slice theorem. Although related techniques for synthesizing plane waves from discrete arrays are known, it is helpful to have available a theory that incorporates the synthesis directly into the scattering theory. The formulation is presented for propagating fields satisfying the Born and Rytov approximations in weakly inhomogeneous media and provides a convenient means for treating both forward and inverse scattering problems. Through a numerical example, we illustrate two important features of diffraction tomography inversion, i.e., 1) the effects of limited view and 2) results of probing with different signal frequencies. The example utilizes data generated by an exact forward modeling technique thus providing strong evidence supporting the usefulness of weak scattering approximations for inversion problems.
A crosswell dataset, collected at Conoco' s Borehole Test Facility in Oklahoma, was first process... more A crosswell dataset, collected at Conoco' s Borehole Test Facility in Oklahoma, was first processed using a tomography algorithm based on an isotropic velocity model. The resulting 2-D velocity tomogram was found to have huge artifacts believed to be caused by nonuniform ray coverage and possibly anisotropy. A comparison of the zero offset crosswell velocities (Horizontal Path Log) with the log-derived sonic velocities (Vertical Path Log) indicated nearly 25% P-wave anisotropy in some formations. Despite the artifacts, a good interpretation of the 2-D isotropic tomogram was made with the help of synthetic data. The data were then inverted using an algorithm that incorporates a model of elliptical transverse isotropy. This anisotropy model produced an unambiguous image for two components of velocity that simultaneously match the sonic log and the crosswell data. Both inversions support a final interpretation of 1-D vertical stratification of layers, some of which exhibit significant P-wave anisotropy. This interpretation is consistent with the sonic logs, crosswell data, and available geological information.
Uploads
Papers by Jerry Harris