Abstract:
A method for enhancing the vertical resolution of a nuclear well logging tool is set forth. A detector of nuclear radiation is moved along a well borehole to produce measurement data r(z) as a function of depth in the borehole. The data r(z) and a system response function g(z) are transformed into spatial frequency domain (f) counterparts of r(z) and g(z), the counterparts being R(f) and G(f). A window W(f) in the spatial frequency domain based on the function R(f) is defined. That smoothly eliminates the high frequency noise. A new profile function P(f) is obtained by dividing the filtered R(f) by the function G(f) to deconvolve R(f). The inverse Fourier transform of P(f) yields a vertically enhanced profile p(z) in the depth domain.

Description:
This application is a continuation of application Ser. No. 07/593,450, filed Oct. 5, 1990 now abandoned. 
    
    
     BACKGROUND OF THE DISCLOSURE 
     This invention relates to nuclear well logging and, more particularly, to methods for enhancing the vertical resolution of nuclear logging instruments. 
     Nuclear well logging tools usually employ detector systems such as gas filled tubes or scintillation detectors for counting either neutrons or gamma rays. These devices are event sensitive in the sense that an event such as passage of a gamma ray or neutron through the detector generates a discrete output or count of its occurrence. The events may occur at high frequency or low frequency depending on the source of the radiation being detected and other factors. But, in any case, the counts of the detector occur in a random, unpredictable time sequence and are summed, accumulated or integrated over a time or depth interval. Thus, the data for a particular time or depth interval is represented by a number of counts having a statistical uncertainty associated with it. 
     The nuclear logging tool outputs are typically sampled on a depth basis in most commercial well logging systems. Material intervening between the nuclear source and the detector of course alters the vertical response of the device. The system response function which is a characteristic of a particular logging tool geometry, collimation, source detector spacing, shielding, etc. also impacts the vertical response of the device. Generally, the vertical response has been thought of as being limited primarily by the source-detector spacing in the past. Movement of the logging tool through the borehole as measurements are made also blurs response to adjacent formations and further integrates or averages physical properties of the formation sought to be measured such as gamma ray activity, formation density and formation porosity. 
     BRIEF DESCRIPTION OF THE INVENTION 
     The techniques and methods of the present invention propose to enhance the vertical resolution of any common type of nuclear well logging tool or device by an appropriate processing of its output combining deconvolution and windowing performed in the spatial frequency domain. The Fast Fourier Transform (FFT) is used to bridge the device output between the depth (z) and the spatial frequency (f) domains. 
     More specifically, an assumed Gaussian system response function in conjunction with a modified Blackman window is utilized to enhance the signal and reduce statistical noise. A brief summary of the application of the technique may be described as follows. The FFT is used to transfer from the depth domain (z) the device output response and the Gaussian system response function into their spatial frequency domain (f) counterparts. An appropriately modified Blackman window is chosen in the spatial frequency domain based on the FFT of the original tool response as a function of depth (i.e., the log). A profile of the windowed convolved response function is then calculated in the frequency domain. Finally, an inverse FFT of the windowed convolved response function of the previous step is obtained to determine p(z), the vertically enhanced log profile back in the depth (z) domain. 
     The methods of the present invention may best be understood by reference to the following detailed description taken in conjunction with the appended drawings. It will be appreciated that the description and drawings may admit to other embodiments of the invention and are to be taken as illustrative only and not as limitative on the concepts of the invention. 
    
    
     BRIEF DESCRIPTION OF THE DRAWINGS 
     FIGS. 1A through H depict various frequency domain windows and the associated FFTs; 
     FIGS. 2A through H depict various frequency domain windows and the associated FFTs; 
     FIGS. 3A and B depict a Gaussian system function and a box-cusp function; 
     FIGS. 4A, B, C and D depict simulated data and deconvolved data with and without noise; 
     FIGS. 5A and B depict simulated data with and without noise and related deconvolved results; 
     FIGS. 6A, B, C and D depict compensated neutron porosity tool data in processed and unprocessed form; 
     FIGS. 7A, B, C and D depict natural gamma ray tool data in processed and unprocessed form; 
     FIGS. 8A, B, C and D depict compensated density logging tool data in processed and unprocessed form; 
     FIGS. 9A and B depict gamma ray spectral tool data in processed and unprocessed form; and 
     FIG. 10 is a schematic diagram depicting a typical well logging system. 
    
    
     DETAILED DESCRIPTION OF THE PREFERRED EMBODIMENT 
     Referring initially to FIG. 10 of the drawings, a typical well logging system layout is shown schematically. A well borehole 10 penetrates earth formations 11, 12, 13 and 14 and is filled with a borehole fluid 15. A nuclear well logging instrument 16 is suspended in the borehole 10 by a winch arrangement (not shown) via a well logging cable 17. The logging cable 17 is an armored, load bearing type having from one to several inner conductors for communicating electrical signals to and from the downhole tool 16. 
     The downhole tool contains a source S of nuclear radiation which may be typically neutron or gamma rays, a longitudinally spaced detector D, an electronic section E for control of the tool and a telemetry section TLM for sending signals to surface equipment and for receiving control signals from the surface equipment. As the tool 16 is moved vertically through the borehole 10, characteristics of the earth formations 11, 12 13 and 14 may be measured and sent via the cable 17 to the telemetry receiver 18 at the surface. 
     The telemetry receiver 18 decodes and formats the data appropriately for input to a computer 19 wherein processing according to concepts of the present invention may be performed by appropriate programs utilizing the methods of the present invention. The computer 19 will then output processed results to a recorder 20 or a tape unit 21 for recording as a function of borehole 10 depth. Depth information is supplied to the computer 19 and the recorder 20 from a sheave or measuring wheel depicted schematically at 22. Movement of the cable 17 across the wheel 22 generates depth pulses or information in a manner known in the art. 
     The log measurement process blurs the formation and averages the properties over some finite depth. The blurring effect mainly depends on the system response function. Deconvolution is a technique used to reduce the blurring effect of the detector, thereby improving the vertical resolution of the log measurement. A simple and applicable Gaussian function is used as the system response function to be deconvolved with the log. Because the deconvolution generally amplifies the high-frequency noise, it must be truncated. However, directly truncating the high-frequency noise introduces some ripples on the log. These ripples are reduced by a smooth truncation procedure called &#34;windowing&#34;, but at the cost of losing some signal. Window characteristics in the depth and the spatial frequency domain show how to filter the noise or ripples and preserve the maximum vertical resolution. 
     The vertical resolution of a logging tool depends on the system response function, which is defined as the depth domain response to a delta function input (an idealized, infinitely thin layer of formation). If the system response function can be well defined, it may be possible to sharpen or enhance the log measurement by deconvolution. Unfortunately, nuclear data are inherently noisy and the exact system response function for nuclear tools is usually unknown. 
     In the case of γ-ray logging, it has been shown by J. G. Conaway and P. G. Killeen, in &#34;Quantitative Determination from Gamma Ray Logs by Application of Digital Time Series Analysis,&#34; Geophysics, v. 43, pp 1204-1221, that a very accurate representation of the intrinisic system function is the convolution of a cusp and a box function which is related to the detector spacing. 
     On the other hand, a Gaussian system function with a variance being a function of the average migration length of the formation was found to be a reasonably accurate system function for neutron porosity tools, see Michael P. Smith, &#34;Enhanced Vertical Resolution Processing of Dual Spaced Neutron and Density Tools Using Standard Shop Calibrations and Borehole Compensation Procedures,&#34; SPWLA 1990. In the former case, the spectra of the Gaussian and the cusp-box function differ only in the high frequency range which we intend to truncate in our procedure. Consequently, a Gaussian function is a reasonable approximation to the system function for nuclear tools. 
     In conjunction with this system function, we have chosen a windowing procedure that suppresses the noise while maintaining maximum vertical resolution. This leads to a Gaussian function, ##EQU1## which was used as the system response function to be deconvolved with the data and the modified Blackman window was used as a filter for all the data processing. The exponential coefficient α in the Gaussian function of Equation 1 is equal to 1.66/FWHM, where FWHM is the full width of the Gaussian at its half maximum. The value of FWHM can be estimated from the tool response. For example, FWHM is about half of the transition length (10% to 90% point) at a boundary. 
     Two sets of simulated data with various thickness beds were generated to verify the effect of the deconvolution and windowing technique. This technique is applied to the laboratory induced γ ray spectrum and the log data from neutron porosity tool, density tool, and natural gamma ray measurement to demonstrate its capability. 
     The Fast Fourier transform (FFT) makes deconvolution and windowing performed in the spatial frequency domain possible and more efficient. A signal, the response of the measuring instrument, is represented mathematically as a function of an independent variable (such as time or depth). In well logging the independent variable is usually depth; therefore, log data are referred to as depth domain signals. Through the FFT, the signal in the depth domain can be converted to the spatial frequency domain. 
     The responses of the logging measurement in the depth domain, denoted as r(z), can be described as convolving the formation profile, p(z), with the intrinsic system function, s(z) where z is the depth. Mathematically; ##EQU2## where * represents a convolution. 
     For linear invariant systems, the previous convolution product in the depth domain can be replaced by a simple multiplication product in the spatial frequency domain. Although it is well known that the system functions for most nuclear logging tools vary with the formation properties, it will be assumed that the error made in maintaining invariance is small compared to the statistical variations of the signal. Therefore, the relationship in Equation (2) can be described in the spatial frequency domain as: 
     
         R(f)=P(f) . S(f)                                           (3) 
    
     Capital letters (e.g., Equation 3) indicate the corresponding spatial frequency domain function. The spatial frequency f is defined in cycles per unit depth. 
     The symbols, r(z) and R(f), p(z) and P(f), s(z) and S(f) are called Fourier transform pairs, where the transform process connecting a pair is denoted as 
     r(z) R(f); 
     p(z) P(f); and 
     s(z) S(f). 
     Fourier transform pairs satisfy the following relationships: ##EQU3## 
     Theoretically, the formation profile in the spatial frequency domain can be obtained by simply dividing the spatial frequency response of the tool by the transfer system function (the Fourier transform of the system function): ##EQU4## The actual formation profile in the depth domain is obtained by taking the inverse transform of the above, which is deconvolution. One problem may occur when the system response function is not well defined or noise exists in the log. There is an optimum filter to be used to deconvolve the signal in a noisy environment. The Weiner filter, W(f), is given by Equation 7: ##EQU5## where the numerator is the complex conjugate and N(f) is the spectra of the noise. 
     It is apparent that, whenever the noise power (|N(f)| 2 ) is large compared with the input power |P(f)| 2 , the Weiner filter vanishes. On the other hand, if the input power at some frequencies is much greater than the corresponding noise power, then, W approaches the reciprocal of S(f) which agrees with Equation (6). Since the input power and noise spectra are usually unknown in a logging environment, the implementation of the Weiner filter is impractical. 
     The filter used in the methods of the present invention is surprisingly useful and easily understood. Furthermore, the results from an application of this filter are less noisy and possesses much better definition than the original data. If G(f) is the Gaussian transfer function to be deconvolved with the log data, then the deconvolved results in the frequency domain P&#39;(f) can be expressed as: ##EQU6## The actual profile can be completely recovered (P&#39;(f)=P(f)), only if the Gaussian function is identical to the intrinsic system function and no noise exists in the log. The noises will be amplified for the high spatial frequency f at which G(f) approaches zero. A window in the spatial frequency domain should be used to reduce these high-frequency noises by setting an upper frequency limit and to compensate for differences between the intrinsic system function and the Gaussian function. Then, Equation (6) can be revised as: ##EQU7## where, W(f) is called the windowing transfer function. The upper limit frequency is called cutoff frequency. W(f)=0, when the frequency f is equal to or larger than cutoff frequency. W(f)/G(f) is sometimes called the deconvolution filter. The procedures for this technique are summarized as follows: 
     (1) Use the FFT to transfer measurement response, r(z), and Gaussian system function, g(z), to the spatial frequency domain as R(f) and G(f); 
     (2) Choose the proper window, W(f), in spatial frequency domain; 
     (3) Calculate P(f)=W(f)R(f)/G(f) for discrete frequencies; and 
     (4) Take the inverse FFT of P(f) to obtain p(z). 
     As mentioned previously, windowing is used to suppress high-frequency noises. A step function truncation of the Fourier series will lead to the well-known Gibbs phenomenon, (i.e., a fixed percentage overshoot and ripple around a discontinuity). In the theory of Fourier series, the Gibbs phenomenon can be moderated through the use of a less abrupt truncation of the Fourier series. Since the multiplication of Fourier coefficients by a window corresponds to convolving the signal response with the Fourier transform of the window, a design objective is to find a finite window whose Fourier transform has relatively small side lobes and narrow main lobe. 
     FIGS. 1A-H and 2A-H show several different kinds of windows in the spatial frequency domain on the left and their Fourier transform on the right of the views. From the relation of Fourier transform pairs, Equation (5) can be written for window pairs as: ##EQU8## i.e., the peak height of the main lobe on the right plot will be equal to the area under the window curve on the left. The rectangular window at the top left corresponds to the step function truncation of the Fourier series. The Fourier transform of the rectangular window, shown at the top right, has a narrow main lobe, but has side lobes which contain a large part of the total energy and decay very slowly. 
     The other three windows are modified Blackman windows, which are represented as: ##EQU9## where n i  is the frequency sample number at the start of the transition zone, and n c  is the frequency sample number at the cutoff frequency. The fourth one is generally called the Blackman window, which is the special case of the modified Blackman window with n i  =0. 
     In FIGS. 1A-H, the four windows have the same cutoff frequency (n c  =100) with different starting points for the transition zone. The area under each window curve in the spatial frequency domain is decreased from the top to the bottom. The peak height of the main lobe in the transform is proportional to the area under the window curve. The larger the peak height of the main lobe is, the narrower the main lobe width will be. The peak height of the side lobes depends on the degree of sharpness in the spatial frequency domain transition zone. The less sharp the transition is, the smaller the side lobes will be. 
     The windows with different cutoffs shown in FIGS. 2A-H have the same area under the curves. The peak height of the main lobe in their transforms is the same (i.e., the vertical resolution in the depth domain is the same for these four windows). The side lobes diminish for the (unmodified) Blackman window. The Blackman window provides the best performance in terms of obtaining the maximum peak height of the main lobe and least peak height of the side lobes, but only if its cutoff frequency is 2.5 times that of the rectangular window. Unfortunately, in real logs, because of the resolution limitation of the detector and noise in the signal, it is difficult to extend the cutoff and frequency of the Blackman window to get maximum resolution without increasing the noise. A satisfactory compromise between the resolution and ripples (noise) is obtained by using the modified Blackman window with suitable adjustments of n i  and n c . 
     The dotted line in FIGS. 4A-D shows an arbitrary intrinsic formation profile in counts per second (CPS) versus depth. The solid curve in FIG. 4A is the noise-free simulated data which was obtained by convolving the formation profile sampled at 0.1 feet with the Gaussian system function (FWHM=1.0 feet) shown in FIG. 3A. The solid curve in FIG. 4B is the deconvolution result of the noise-free data. The same Gaussian function with more samples was used as the deconvolver. The result shows the vertical resolution of all the beds is enhanced. The center bed readings reach the true value, except the 0.5&#39; bed. Notice that there are no overshoots. 
     FIG. 4C illustrates the case where random noise, based on counts per sample assuming a logging speed of 30&#39;/min, was added to the noise-free simulated data. FIG. 4D indicates that the vertical resolution is enhanced and the noise is reduced. 
     Another simulated log profile is denoted via FIG. 5A by the dotted line. The dashed curve in this figure is the noise-free data, which was obtained by convolving the profile sampled at 0.25&#39; with the box-plus-cusp functions (dashed line in FIG. 3B). The deconvolver is the Gaussian function with a FWHM of 4 samples (solid line in FIG. 3B). The deconvolved results of the noise-free data (solid curve) shows a very good enhancement of the vertical resolution for all the beds. The essential reason for using the box-plus-cusp system function, which is a better approximation to the system function, was to indicate the robustness of the Gaussian deconvolution procedure. 
     FIG. 5B shows the same formation with noisy data and the deconvolved results. Except for a slight overshoot on the 2&#39; bed at a depth of 30&#39;, which is caused by noise amplification, the vertical resolution of all the beds was improved. 
     The results of the simulated data show that the deconvolution and windowing technique can be successfully used to enhance the vertical resolution and reduce the noise if the system function is well defined. 
     FIGS. 6A-D show the raw data from a compensated neutron tool near and far detectors at 0.1&#39; sample rate and the deconvolution results with a modified Blackman window. A Gaussian system function with a FWHM of 10 and 7 samples was used for far and near data, respectively. The deconvolved logs are sharper and cleaner than the raw data logs. The thin beds at depths 1852&#39;, 1869&#39; and 1879&#39; show the discernible enhancement. 
     FIGS. 7A-D shows two sets of raw data within the same depth region from the natural gamma ray tool at a 0.1 feet sample rate and the deconvolved results with a modified Blackman window. A Gaussian system function with a FWHM of 7 samples was used. Notice, in the deconvolved results, the vertical resolution is enhanced, the noise is reduced and the thin beds become sharper. 
     FIGS. 8A-D show the raw data of the near and the far detectors from a compensated density log (CDL) tool at a 0.1 feet sample rate and the deconvolved results. Sodium iodide detectors are used for the CDL tool. A Gaussian system function with a FWHM of 6 and 5 samples was used for far (FCDL) and near (NCDL) data, respectively. The log is cleaner and possesses some noticeable thin bed enhancements at the depths 1924.8&#39;, 1926.2&#39;, 1926.9&#39;, and 1927.6&#39;. 
     In the top of FIG. 9 is the raw data of a laboratory induced γ-ray spectrum obtained in a 20% water-filled sandstone formation with an 8&#34; borehole. The deconvolved spectrum using a modified Blackman window is shown in the bottom of the same figure. The FWHM value of the Gaussian system function is 5 sample units. Notice, that the deconvolved result provides a better resolved set of data with reduced statistical variations. 
     The results of the simulated, laboratory and nuclear log data show that deconvolution and windowing technique can enhance the vertical resolution of the thin beds, sharpen the boundary of the beds (thicker than 2&#39;) without degrading the center bed reading, and reduce the noise. A simple Gaussian function is appropriate as the system response function for the nuclear tools described here. Moreover, the bandwidth of the modified Blackman window is readily obtained from the Fourier transform of the recorded log data. 
     The foregoing description may make other changes and modifications of the techniques of the present invention apparent to those of skill in the art. The aim of the appended claims is to cover all such changes and modifications as fall within the true spirit and scope of the invention.