Multidimensional recursive filter preconditioning in geophysical estimation problems

Multidimensional examples

Our first multidimensional example is the SeaBeam dataset, a result of water bottom measurements from a single day of acquisition. SeaBeam is an apparatus for measuring water depth both directly under a ship and somewhat off to the sides of the ship's track. The dataset has been used at the Stanford Exploration Project for benchmarking different strategies of data interpolation. The left plot in Figure 12 shows the original data. The right plot shows the result of (unpreconditioned) missing data interpolation with the Laplacian filter after 200 iterations. The result is unsatisfactory, because the Laplacian filter does not absorb the spatial frequency distribution of the input dataset. We judge the quality of an interpolation scheme by its ability to hide the footprints of the acquisition geometry in the final result. The ship track from the original acquisition pattern is clearly visible in the Laplacian result, which is an indication of a poor interpolation method.

seabdat
Figure 12. On the left, the SeaBeam data: the depth of the ocean under ship tracks; on the right, an interpolation with the Laplacian filter.

We can obtain a significantly better image (Figure 13) by replacing the Laplacian filter with a two-dimensional prediction-error filter estimated from the input data. The result in the left plot of Figure 13 was obtained after 200 conjugate-gradient iterations. If we stop after 20 iterations, the output (the right plot in Figure 13) shows only a small deviation from the input data. Large areas of the image remain unfilled. At each iteration, the interpolation process progresses only to the length of the filter.

seabold
Figure 13. SeaBeam interpolation with the prediction-error filter. The left plot was taken after 200 conjugate-gradient iterations; the right, after 20 iterations.

Inverting the PEF convolution with the help of the helix transform, we can now apply the inverse filtering operator to precondition the interpolation problem. As expected, the result after 200 iterations (the left plot in Figure 14) is similar to the result of the corresponding unpreconditioned (model-space) interpolation. However, the output after just 20 iterations (the right plot in Figure 14) is already fairly close to the solution.

seabnew
Figure 14. SeaBeam interpolation with the inverse prediction-error filter. The left plot was taken after 200 conjugate-gradient iterations; the right, after 20 iterations.

For a more practical test, we chose the North Sea seismic reflection dataset, previously used for testing azimuth moveout and common-azimuth migration (Biondi, 1997; Biondi et al., 1998). Figure 15 shows the highly irregular midpoint geometry for a selected in-line and cross-line offset bin in the data. The data irregularity is also evident in the bin fold map, shown in Figure 16. The goal of data regularization is to create a regular data cube at the specified bins from the irregular input data, which have been preprocessed by normal moveout without stacking.

cmp-win Figure 15. Midpoint distribution for a 50 by 50 m offset bin in the 3-D North Sea dataset.

fold-win Figure 16. Map of the fold distribution for the 3-D data test.

The data cube after normalized binning is shown in Figure 17. Binning works reasonably well in the areas of large fold but fails to fill the zero fold gaps and has an overall limited accuracy.

bin-win
Figure 17. 3-D data after normalized binning.

For efficiency, we perform regularization on individual time slices. Figure 18 shows the result of regularization using bi-linear interpolation and smoothing preconditioning (data-space regularization) with the minimum-phase Laplacian filter (Fomel et al., 2002). The empty bins are filled in a consistent manner but the data quality is distorted because simple smoothing fails to characterize the complicated data structure. Instead of continuous events, we see smoothed blobs in the time slices. The events in the in-line and cross-line sections are also not clearly pronounced.

smo2-win
Figure 18. 3-D data regularized with bi-linear interpolation and smoothing preconditioning.

We can use the smoothing regularization result to estimate the local dips in the data, design invertible local plane-wave destruction filters (Fomel, 2001), and repeat the regularization process. Inverse interpolation with plane-wave data-space regularization is shown in Figure 19. The result is noticeably improved: the continuous reflection events become clearly visible in the time slices. Despite the irregularities in the input data, the regularization result preserves both flat reflection events and steeply-dipping diffractions. Preserving diffractions is important for correct imaging of sharp edges in the subsurface structure (Biondi et al., 1998).

For simplicity, we assumed only a single local dip component in the data. This assumption degrades the result in the areas of multiple conflicting dips, such as the intersections of plane reflections and hyperbolic diffractions in Figure 19. One could improve the image by considering multiple local dips. Fomel (2002) describes an alternative offset-continuation approach, which uses a physical connection between neighboring offsets instead of assuming local continuity in the midpoint domain.

int4-win
Figure 19. 3-D data regularized with cubic B-spline interpolation and local plane-wave preconditioning.

The 3-D results of this paper were obtained with an efficient 2-D regularization in time slices. This approach is computationally attractive because of its easy parallelization: different slices can be interpolated independently and in parallel. Figure 20 shows the interpolation result for four selected time slices. Local plane waves, barely identifiable after binning (left plots in Figure 20), appear clear and continuous in the interpolation result (right plots in Figure 20). The time slices are assembled together to form the 3-D cube shown in Figure 19.

winslice
Figure 20. Selected time slices of the 3-D dataset. Left: after binning. Right: after plane-wave data regularization. The data regularization program identifies and continues local plane waves in the data.

Multidimensional recursive filter preconditioning in geophysical estimation problems