Efficient inversion of near singular geophysical signals

Details Efficient inversion of near singular geophysical signals Efficient inversion of near singular geophysical signals

2000-01-14

2002-11-05

Lefkowitz, Edward (Department: 2862)

Data processing: measuring, calibrating, or testing

Measurement system in a specific environment

Earth science

C703S008000

Reexamination Certificate

active

06477467

ABSTRACT:

BACKGROUND OF THE INVENTION
1. Field of the Invention
The present invention relates generally to seismic signal processing and, more particularly, to a method for efficient inversion of near singular geophysical signals.
2. Description of the Related Art
A common data processing activity in exploration geophysics is inverting geophysical signals. This can arise in Wiener filtering in its various forms, for example, deconvolution, cross equalization, adaptive filtering and echo cancellation. Such signals are characterized by a convolutional equation system and inverting the signals are all variations on the convolutional model least squares problem. These problems are discussed in Treitel, S. and Lines, L. R., “Linear Inverse Theory and Deconvolution” Geophysics, vol. 47, p p. 1153-1159, 1982 and in Lines, L. R. and Treitel, S., “Tutorial: A Review of Least-Squares Inversion and its Applications to Geophysical Problems”, Geophysical Prospecting, Vol. 32, pp.159-186, 1984. For geophysical signals, the corresponding systems are often near singular (ill-conditioned) or rank-deficient. Classic methods for inverting such systems require stabilization methods such as adding white random noise, using a &ggr; damping factor, spectral shaping or constraint, or other so-called regularization methods. For the convolutional model, this is usually done in the normal equation form using Levinson recursion on the Toeplitz autocorrelation matrix, as will be explained below.
In this class of problems, one is given an input signal i(t), an output signal o(t) and a finite impulse response filter u(t) having the relationship:
i
(
t
)=i {circle around (x)}u(
t
)=
o
(
t
)  (1)
Here the symbol {circle around (x)} indicates convolution and the vectors i, u and o are each functions of time t. It is typically desired to find the particular solution filter or wavelet û(t) which which gives the minimum error e(t) in the least squares sense. The solution to the least squares minimization problem is found by first defining the error vector by
e
(
t
)=
i
(
t
)
{circle around (x)}ûl (
t
)−
o
(
t
)  (2)
and then minimizing the mean squared scalar error &egr; defined by the inner product of the error vector e from Eq. (2) with itself
&egr;=
e
H
·e.
  (3)
Here the superscript
H
indicates the Hermetian transpose operator or complex conjugate, and the dot · indicates the vector inner product.
For a filter vector u(t) of length n and an output vector o(t) of length m in Eq. (1), the input vector i(t) will be of length m+n−1. Convolution of the input vector i with filter vector u may be represented by multiplication of u by an m×n matrix T, where T is the Toeplitz matrix defined by
T
(
t,&tgr;
)=
i
(
n−&tgr;+t−
1)  (4)
for t=0, . . . m+n−2 and &tgr;=0, . . . n−1.
FIG. 1
shows the matrix T. Eq. (4) allows Eq. (1) to be rewritten in matrix form as
T·u=o
  (5)
Eq. (5) is known as the direct form equations of the least squares problem. In geophysical signal processing, the number of data observations, represented by rows m, is usually larger than the number of physical parameters, represented by columns n. Thus the direct form equations, Eq. (5), are typically a system of m overdetermined equations.
The least squares solution û may be found from the normal form equations, an exactly determined Toeplitz n×n system
(
T
H
·T
)·
û=T
H
·o
  (6)
Eq.(6) is Hermetian Toeplitz if i(j)=0 for j=0 . . . , n−2 and j=m, . . . , m+n−2. This is shown in
FIG. 2
, which indicates a column circular matrix.
FIG. 1
then represents the transient-free case, while
FIG. 2
represents the more typical fully-windowed case.
Eq. (6) is typically simply written as
A·û=c
  (7)
where A=T
H
·T is the n×n autocorrelation or Gram matrix and c=T
H
·o is the n×1 cross-correlation vector. The typical method for solving Eq. (7) has been Levinson recursion, which will solve an n×n Toeplitz system proportional to n
2
multiply-add steps, that is, has order O(n
2
) computational complexity. The computations of matrix A and vector c will then add another (m+1)−n steps to these calculations.
The convolutional model least squares problems encountered in geophysical signal processing are often ill-conditioned. A measure of the condition of a matrix T is the condition number cond(T). The condition number of a matrix T is defined as the ratio of the largest to the smallest singular values &sgr;(i) in the singular value decomposition of T or the square root of the ratio of the largest to the smallest eigenvalues &lgr;(i) in the eigenvalue decomposition of the corresponding Gram matrix A. Thus
cond(
T
H
·T
)=cond(
T
)
2
=cond(
A
)=(&sgr;
max
/&sgr;
min
)
2
=&lgr;
masx
/&lgr;
min
  (8)
where
&sgr;
max
=maximum singular value of T,
&sgr;
min
=minimum singular value of T,
&lgr;
max
=maximum eigenvalue of A, and
&lgr;
min
=minimum eigenvalue of A.
Thus the normal form equations given in Eqs. (6) and (7) have a condition number which is the square of the condition number of the direct form equations of Eq. (5).
Matrices with larger condition numbers are more difficult to accurately invert. This suggests that it would be advantageous to solve the direct form of the least squares problem in Eq. (5) without resorting to conversion to the less well-conditioned normal form of Eqs. (6) and (7). This is usually done by using any number of matrix decomposition techniques, all of which are O(m·n
2
) in computational complexity. The advantage to the normal form equation approach for Toeplitz systems is reduced computations at the expense of numerical stability.
The matrix A is defined as singular if &lgr;
min
=0, and as nearly singular or ill-conditioned if |&lgr;
min
|<&dgr;, for an appropriately selected small &dgr;. The eigenvalues of Gram matrices are nonnegative. For a practical lower limit, one may choose &dgr;=n·&lgr;
max
·eps where eps is the limit of relative floating point accuracy on the computer used. For IEEE standard 754 double precision floating point arithmetic, eps=2
−52
≈10
−16
The maximum conditioning number defining ill-conditioning, c
max
, may be estimated by substituting |&lgr;
min
|<&dgr; into Eq. (8), yielding c
max
=1/(n·eps). Then a nearly singular or ill-conditioned system is defined by cond(A)>c
max
.
For a well-conditioned matrix A, a least squares solution û to Eq. (7) is typically found with a matrix decomposition technique such as QR or LU matrix decomposition. The ill-conditioned matrix case, where cond(A)>c
max
, may be solved using singular value decomposition. The rank-deficient case may also be solved by the singular value decomposition approach. Another approach for an ill-conditioned matrix A is to use a small damping factor satisfying &ggr;<<|&lgr;
max
| and &dgr;<&ggr;<|&lgr;
min
|. This leads to solving the new system
(
A+&ggr;·I
)·
û=c
  (9)
where I is the identity matrix. The singular value decomposition is the stablest approach, but also the slowest approach.
In geophysical applications, the spectral content of the input signal i(t) is often characterized by regions of low amplitude caused by natural shaping, filtering, and noise or reverberation. Such spectral gaps or notches make spectral inversion very difficult without resorting to stabilization methods such as adding white random noise, using a &ggr; damping factor as in Eq.(9), spectral shaping or constraint, or other so-called regularization methods. The solutions, wavelets, obtained by traditional methods are often ill-behaved and as well may not give useful results. Whereas using methods such as singular value decomposition lead to excellent results, they are generally too computationally ex

Affiliated with

Also associated with

Rating

Efficient inversion of near singular geophysical signals does not yet have a rating. At this time, there are no reviews or comments for this patent.

If you have personal experience with Efficient inversion of near singular geophysical signals, we encourage you to share that experience with our LandOfFree.com community. Your opinion is very important and Efficient inversion of near singular geophysical signals will most certainly appreciate the feedback.

Rate now

Comments { 0 }

Profile ID: LFUS-PAI-O-2919808