Troyan V., Kiselev Y. Statistical Methods of Geophysical Data Processing

Подождите немного. Документ загружается.

Methods of transforms and data analysis 373

The subscripts in ν and ω can be omitted, because after proceeding to limit it does

not has any value. The sum

f(t) =

∞

n=−∞

exp{inω

transforms to the integral

f(t) =

∞

−∞





∞

−∞

f(t) exp{iωt}dt





exp{iωt}dω/2π. (12.5)

The inner integral from the right hand side of the equality (12.5) is called the Fourier

transform of the function f (t):

f(ω) =

∞

−∞

f(t) exp{−iωt}dt, (12.6)

and outer integral from the right hand side of the equality (12.5) is called the inverse

Fourier transform of the function

f(ω):

f(t) =

2π

∞

−∞

f(ω) exp{iωt}dω. (12.7)

There are other kinds of the Fourier and inverse Fourier transforms, for example,

with multiplier

1/2π in both Fourier transform and inverse Fourier transform.

The suﬃcient condition for the existence of the Fourier transform of the function

f(t) consists in a boundedness of the integral

∞

−∞

|f(t)|dt.

If the function f(t) is real one, then inverse Fourier transform can be written

down as

f(t) =

2π

∞

−∞

R(ω) cos ωtdω −

2π

∞

−∞

X(ω) sin ωtdω, (12.8)

where

R(ω) =

∞

−∞

f(y) cos ωydy (12.9)

is cosine transform of f (t),

−X(ω) =

∞

−∞

f(y) sin ωydy (12.10)

374 STATISTICAL METHODS OF GEOPHYSICAL DATA PROCESSING

is sine transform of f(t). In addition the following formula is valid

f(ω) = R(ω) + iX(ω) = A(ω) exp{iγ(ω)},

where

A(ω) =



(ω) + X

(ω)



1/2

(12.11)

is a amplitude spectrum,

γ(ω) = arctg [X(ω)/R(ω)] (12.12)

is a phase spectrum.

In the case of multivariable function f(t

, t

, . . . , t

), if the following integral

exists

∞

−∞

∞

−∞

. . .

∞

−∞

|f(τ

, τ

, . . . , τ

)|dτ

dτ

. . . dτ

then the direct

C(ω

, . . . , ω

) =

(2π)

1/2

∞

−∞

∞

−∞

. . .

∞

−∞

f(τ

, τ

, . . . , τ

)

× exp





−i

j=1





dτ

. . . dτ

and inverse

f(t

, . . . , t

) =

(2π)

1/2

∞

−∞

∞

−∞

. . .

∞

−∞

C(ω

, ω

, . . . , ω

)

× exp





−i

j=1





dω

. . . dω

Fourier transforms are determined.

12.2 Laplace Transform

The Laplace transform is closely connected with the Fourier transform. If we do

not use generalized functions, then many functions, for example, sin at and cos at

functions, do not have the Fourier image, because the direct Fourier integral is

divergent one. But if we multiply appropriate function by exp (−σ|t|), where σ is

a real positive value, which satisﬁes the condition lim

t→±∞

{exp(−σ|t|)f (t)} = 0, then

the function {exp(−σ|t|)f(t)} has the Fourier image. Such transform is called the

Methods of transforms and data analysis 375

Laplace transform of f(t). If f (t) = 0 at t < 0, then we obtain the one-side Laplace

transform

F (s) =

∞

f(t) exp(−st)dt, (12.13)

where s = σ + iω, in addition the following condition lim

t→±∞

{exp(−s|t|)f(t)} = 0

must be satisﬁed. The inverse Laplace transform, can be written as

f(t) =

2πi

σ+i∞

σ−i∞

F (s) exp(st)ds, (12.14)

where the path of integration is a straight line, which is parallel to the imaginary

axis and placed in the right half-plane of the complex plane.

In many cases Laplace transform is more preferable than Fourier transform, for

example, in the cases when the function F (s) from the transform of (12.14) has the

pole particularities in the imaginary axis.

12.3 Z-Transform

Z-transform is a special kind of the transform, which is useful for calculations with

the functions represented by their discrete sequences.

Let us consider the sequence of numbers b

, which is a discrete representation

of the continuous function f(t) at given time points t

= k∆ (k = . . . , −1, 0, 1, . . . ),

where ∆ is a sampling interval. Let us write the discrete image f

(t) of the function

f(t) as

(t) =

∞

k=−∞

δ(t − k∆), (12.15)

where δ(t) is delta function. The Fourier transform of f

(t) can be written as the

following sum

B(ω) =

∞

−∞

(t) exp(−iωt)dt =

∞

k=−∞

exp(−ikω∆). (12.16)

Introducing the notation Z = exp(−iω∆) we can write B(ω) (12.16) as polynomial

B(Z) =

∞

k=−∞

, (12.17)

where B(Z) is Z-transform of f

(t)

Once Z is determined as Z = exp(iω∆), then a polynomial with reversed sign of Z is obtained

in comparison with the sign of the polynomials which are considered above.

376 STATISTICAL METHODS OF GEOPHYSICAL DATA PROCESSING

12.4 Radon Transform for Seismogram Processing

The direct Radon transform reads as

s(p, τ ) =

∞

−∞

u(x, px + τ)dx

(Fig. 12.1). The inverse Radon transform is determined as

Fig. 12.1 The direct and inverse Radon transform. The straight line in the space-time plane (a)

and the point in the plane angle coeﬃcient – reﬂection time (b). The point in the plane angle

coeﬃcient – reﬂection time (c) and the line in space-time plane (d).

u(x, t) = −

∞

−∞

dτ

s(p, τ )dp = −P HD[s(p, τ)],

u(x, t) = −

∞

−∞



−

πτ

∗

∂

∂τ

s(p, τ )



dp.

Let us write the two-dimensional Fourier transform of the seismogram (Fig. 12.2):

u(k

, ω) =

∞

−∞

u(x, t) exp[i(ωt − k

x)]dxdt.

Using the formula k

= pω, we can write the spectrum u(k

, ω) in the form

u(pω, ω) =

∞

−∞

u(x, t) exp[i(ωt − pωx)]dxdt,

Methods of transforms and data analysis 377

Fig. 12.2 The line in the space-time plane (a) and in the wave number – frequency plane (b).

Fig. 12.3 The window function for the integration in k-domain.

or, by introducing variable t = τ + px, we obtain the following expression:

u(pω, ω) =

∞

−∞

iωτ

∞

−∞

u(x, τ + px)dx

dτ

∞

−∞

iωτ

s(p, τ )dτ = F R[u(x, t)],

where s(p, τ) =

−∞

u(x, τ + px)dx is the Radon transform of the seismogram

u. Thus, the two-dimensional Fourier transform can be represented as the Radon

transform and the one-dimensional Fourier transform.

The ﬁltration in the domain wavenumber–frequency (k − ω) reduces to the in-

tegration with the limits k

min

and k

max

(Fig. 12.3):

u(x) =

2π

max

min

u(k) exp(ikx)dk,

where k = |ω|p, or

u(x) =

|ω|

2π

max

min

u(p) exp(iωpx)dp (dk = |ω|dp).

The complete ﬁltration of the seismogram in (ω, k) plane can be written down as

378 STATISTICAL METHODS OF GEOPHYSICAL DATA PROCESSING

following (Fig. 12.4)

u(x, t) =

4π

max

min

max

min

u(ω, p)|ω|exp[−iω(t − px)]dpdω

2π

max

min



2π

max

min

exp[−iω(t − px)] × F R[u(x, t)]|ω|dω



2π

max

min

R[u(x, t)] ∗ (F

−1

[|ω|])dp.

Fig. 12.4 The domain of integration in (ω, k) plane.

12.5 Gilbert Transform and Analytical Signal

Let the function S(t) is a bounded function with the norm L

, i.e.

S(t) ∈ L

(−∞, ∞),

∞

−∞

|S(t)|

dt < ∞, p ≥ 1.

In this case it is possible to deﬁne the direct

Q(t) = −

V.P.

∞

−∞

S(τ)

t − τ

dτ = −S(t) ∗

πt

and inverse Gilbert transform

S(t) =

V.P.

∞

−∞

Q(τ)

t − τ

dτ,

where V.P. denotes the Cauchy principal value of the integral at the point t = τ.

Let us calculate the direct Fourier transform of the Gilbert transform Q(t):

Q(ω) =

∞

−∞

Q(t) exp(iωt)dt = −

∞

−∞

S(τ)

∞

−∞

exp(iωt)

τ − t

dtdτ.

Methods of transforms and data analysis 379

After change of variable τ −t = u we obtain the simple expression for the spectrum:

Q(ω) = −

∞

−∞

S(τ) exp(iωτ )

∞

−∞

exp(iωu)

dudτ

= −S(ω)

∞

−∞

exp(iωu)

= −S(ω)

∞

sin ωu

du = −iS(ω) sign ω,

because (Fig. 12.5)

∞

sin ωu

du =

sign ω.

We should note, that |Q(ω)| = |S(ω)|, arg Q(ω) = arg S(ω) ± π/2.

Fig. 12.5 The graphic representation of the function sign ω.

The function Q(ω) is a quadrature ﬁlter, which changes a phase of the signal by

±π/2, but does not change its amplitude.

Let Z(t) is a complex (temporal) function or an analytical signal:

Z(t) = S(t) + iQ(t), or Z(t) = a(t) exp(iϕ(t)),

where

S(t) = ReZ(t) = a(t) cos ϕ(t),

Q(t) = ImZ(t) = a(t) sin ϕ(t).

The amplitude and phase function can be represented as following

a(t) = (S

(t) + Q

(t))

1/2

and ϕ(t) = arctan(Q(t)/S(t)).

The amplitude and phase function are used for the interpretation of seismic data.

380 STATISTICAL METHODS OF GEOPHYSICAL DATA PROCESSING

12.6 Cepstral Transform

The convolution can be realized as a multiplication of the spectra of the initial

signals, which are obtained with the help of the direct Fourier transform (from

the time domain to the spectrum domain). The cepstral transform (from the time

domain to the cepstral domain) allows to substitute the convolution of the signals

for the summation of their cepstral images. Besides, in some cases a separation

of the signals is better to implement in the cepstral domain than in the spectral

domain, that can leads to more eﬀective ﬁltration.

The inverse Fourier transform of the logarithm of the frequency spectrum yields

a cepstrum

f(ζ) = 1/2π

∞

−∞

ln[

f(ω)] exp(iωζ)dω. (12.18)

The transform from the time domain to the cepstral domain can be implemented

per four steps:

f(ω) =

∞

−∞

f(t) exp(−iωt)dt,

f(ω) = |

f(ω)|exp[iϕ(ω)],

3) ln

f(ω) = ln |

f(ω)| + iϕ(ω),

f(ζ) = 1/2π

∞

−∞

[ln |

f(ω)| + iϕ(ω)] exp(iωζ)dω.

(12.19)

To return in the time domain we should implement the inversion of the transforms

(12.19):

1) f (ω) =

∞

−∞

f(ζ) exp(−iωζ)dζ,

f(ω) = exp[f(ω)],

3) f (t) = 1/2π

∞

−∞

f(ω)|exp(iωt)dω.

(12.20)

The variable ζ is called the cepstral frequency. By analogy with a ﬁltering in

the time and frequency domains, the ﬁltering in the cepstral domain is called the

cepstral ﬁltering.

Let’s note that during passing to the cepstral domain an important step consists

in the phase determination ϕ(ω). It is due to the periodic property of the tangent,

which has the period nπ. Therefore the phase value ϕ(ω) can be determined ac-

curate within nπ. If the phase ϕ(ω) is assumed a continuous value, then at their

calculation we should use, for example, a criterion of the smallness of the absolute

value of the phase function.

Methods of transforms and data analysis 381

12.7 Bispectral Analysis

In Sec. 1.10.5 it was introduced the correlation function for a real random function.

Let us consider an autocovariance function of the second order

b(τ

, τ

) for the real

random function ξ(t):

b(τ

, τ

) =

∞

−∞

ξ(t)ξ(t + τ

)ξ(t + τ

)dt. (12.21)

After Fourier transform of b(τ

, τ

), we obtain the two-dimensional spectrum

B(ω

, ω

) =

∞

−∞

∞

−∞

b(τ

, τ

) exp{−i(ω

+ ω

)}dτ

dτ

(12.22)

of the autocorrelation function, at that the next inequality is valid

B(ω

, ω

) = f(ω

)f(ω

∗

(ω

+ ω

where f(ω) is the Fourier transform of ξ(t), symbol ∗ denotes a complex conjugation.

Let us note, that the function b(τ

, τ

) from the formula (12.21) is a symmetric one:

b(τ

, τ

) = b(τ

, τ

) = b(−τ

, τ

− τ

)

= b(τ

− τ

, −τ

) = b(−τ

, τ

, −τ

) = b(τ

− τ

, −τ

from that follows the next equalities for the spectral function B(ω

, ω

B(ω

, ω

) = B(ω

, ω

) = B(ω

, −ω

− ω

)

= B(−ω

− ω

, −ω

) = B(ω

, −ω

− ω

) = B(−ω

− ω

, −ω

Since we consider here a real process ξ(t), the next equality is valid also

B(ω

, ω

) = B

∗

(−ω

, −ω

Let us consider the relation between the phase function of a linear ﬁlter and

the phase characteristic of the bispectral function of the output signal. Let us the

process x

is represented as the discrete convolution of the ﬁlter w

with the “input”

process g

∞

s=0

t−s

If the process g

is the stationary and uncorrelated one, then the second order

autocovariance of the process x

can be written as

b(τ

, τ

) = γ

∞

s=0

s+τ

, (12.23)

Here we do not consider the living conditions of the autocovariance function of the second order.

In action, the random function have a limited length, therefore the existence of similar function

does not raise doubts.

382 STATISTICAL METHODS OF GEOPHYSICAL DATA PROCESSING

where γ

= M[g

], w

= 0, s < 0. We represent spectrum B

(ω

, ω

) of the

autocovariance function b(τ

, τ

) from the relation (12.23) as multiplication of the

ﬁlter W (ω) spectral functions:

(ω

, ω

) = γ

∞

=−∞

∞

=−∞

∞

s=−∞

s+τ

× exp[−i(ω

+ ω

)] = γ

W (ω

)W (ω

)W (−(ω

+ ω

))

= γ

W (ω

)W (ω

∗

(ω

+ ω

where

W (ω) =

∞

s=−∞

exp(−iωs),







−π ≤ ω

≤ π,

−π ≤ ω

≤ π,

−π ≤ ω

+ ω

≤ π.

We note, that the spectrum of the second order autocovariance function of x

accu-

rate within multiplier coincides with spectrum of the second order autocovariance

function of the ﬁlter w

(ω

, ω

) = γ

(ω

, ω

If we denote the phase of the spectra of the second order covariance function of the

ﬁlter and the ﬁlter phase function as ψ(ω

, ω

) and ϕ(ω) respectively

(ω

, ω

) = G

(ω

, ω

) exp(iψ(ω

, ω

)),

W (ω) = H(ω) exp(iϕ(ω),

then we can write the equality

ψ(ω

, ω

) = ϕ(ω

) + ϕ(ω

) − ϕ(ω

+ ω

). (12.24)

So, the calculation of the spectrum phase ψ(ω

, ω

) of the second order autocovari-

ance function of x

, it is possible to ﬁnd a spectrum phase of the ﬁlter ϕ(ω).

12.8 Kalman Filtering

The recurrence procedure, which has been considered in Sec. 6.11, can be assumed

as a basis of the construction of the algorithm of the dynamic (Kalman) ﬁltration.

A shape of the seismic signals can be changed during the propagation. Such changes

can be taken into account using the following model:

+ ε

, ε

∈ N (0, σ

= F

ii−1

i−1

+ Gδ

i−1

, δ

∈ N(0, σ

where

is an amplitude vector for i-th sample:







(1)

(M)







i−1







(1)

i−1

(M)

i−1







, G =











