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

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

Construction of tomographic functionals 303

and taking into account that conditions ϕ

out

+∞

= 0 , ∂/∂t

−∞

= 0, we

obtain

: p



out

|ϕ



∂

∂t

out

⊗

∂

∂t

. (10.20)

Thus, formulas (10.15), (10.18) and (10.20) determine the interaction operators S

and S

for the ﬁelds ϕ

and ϕϕ

out

. Physical interpretation of the tomographic

functionals is that their integral kernels are spatial functions of the inﬂuence of varia-

tions of the ﬁelds of the required medium parameters on the sampling experimental

values of the wave ﬁeld of the sounding signal. The support of the tomographic

functional is determined by the region of interaction of the ﬁelds

and ϕϕ

out

. This

region is certainly bounded on account of physical conditions of the tomographic

experiment (a ﬁnite temporal interval of measurements and a ﬁnite velocity of the

ﬁeld ϕ

out

propagation) .

The norm of the tomographic functional is determined by the amplitude of the

inﬂuence function related to the interaction operator. In Fig. 10.2, the diagrams

of the scattering of the plane wave by the elementary perturbations of the shear

module µ and the mass density ρ are represented (the perturbation of the Lame

parameter λ has an isotropic diagram). The diagrams characterize the interaction

of the ﬁelds

and ϕ

out

, if their sources are placed far from the scattering object

(near-ﬁeld region). The ﬁeld is generated by a “source” with the time dependence

determined by the apparatus function of the seismic recording channel. Since the

frequency band of the apparatus function is limited, the interaction of the ﬁelds ϕ

and

out

, which is always a convolution in the time domain (and a product in the

frequency domain), involves only the overlapping spectral components of the ﬁelds

ϕϕ

and ϕ

out

The analysis of tomographic functionals provides the possibility of mathemat-

ically designing a tomographic experiment by monitoring of the parameters in-

volved in tomographic functional (the localization of the source and receivers of

the ﬁelds ϕ

and ϕ

out

in space, the shape of the sounding signal ϕϕ

, and the

choice of the time interval for recording the ﬁeld ϕϕ

out

). In this case, the character

of the interaction is determined only by the inﬂuence of the medium parameters

upon the physics of propagation of the sounding signal. For example, if the ﬁeld

is pure rotational (the transversal wave in an isotropic reference medium), then

the tomographic functional p

becomes zero. Accordingly, such a tomographic ex-

periment obviously provides no information on the ﬁeld of the parameter λ. A

physically informative tomographic experiment is characterized not only by the

norms of tomographic functionals, but also by the domain of the overlapping of

their support.

304 STATISTICAL METHODS OF GEOPHYSICAL DATA PROCESSING

Fig. 10.2 The diagrams of the scattering of shear and longitudinal wave by the elementary per-

turbations. To the right of the spatial diagrams their sections in the plane (x, z) are represented;

and a

are respectively the normal to the wave front and the unit amplitude vector of the

incident q-wave (q = p, s).

Construction of tomographic functionals 305

Consider the case where the reference medium may be treated as homogeneous.

Let the variations of the elastic parameters in a homogeneous isotropic medium be

such that |∇∇λ|, |

∇

∇µ| and δρ are small enough. We consider two cases where the

ﬁeld ϕ

is either a part of the pure longitudinal wave or the pure transversal wave.

In the ﬁrst case, we introduce the function

) = ρ

− 1

where, in addition, ρ

δ(v

) = δ(λ + 2µ). Writing the action of the perturbation

operator upon the ﬁeld ϕϕ

≡ ϕϕ

, where ϕϕ

∇

∇ × ϕ

≡ 0, we obtain

δLϕ

= δλ∇∇∇ ·

+ 2δµ∆ϕ

= (δλ + 2δµ)∆ϕ

= δ(λ + 2µ)∆ϕ

. (10.21)

Making an observation that at the space points where the sources are absent the

equation



∂

∂t

− (λ

+ 2µ

)∆



= 0 ,

holds true, we replace ∆

with (ρ

/λ

+ 2µ

) ∂

/∂t

in (10.21), which yields

δL

= δν

∂

∂t

The corresponding tomographic functional and the interaction operator are of the

form

: p

= hϕ

|ϕ

out

i = −

∂

∂t

⊗

∂

∂t

ϕϕ

out

. (10.22)

In the second case, we introduce the function

) = ρ



− 1



Then

ρδ(v

) = δµ .

The action of the perturbation operator upon the ﬁeld ϕ

≡ ϕϕ

, where

∇

∇·ϕ

≡

0, is expressed as

δLϕ

= δµ∆ϕ

By using the equation



∂

∂t

− µ∆



= 0

for the solenoidal ﬁeld ϕϕ

, we obtain

δL

= δν

∂

∂t

306 STATISTICAL METHODS OF GEOPHYSICAL DATA PROCESSING

The tomographic functional for ν

is of the form

: p

= hϕϕ

|ϕ

out

i = −

∂

∂t

ϕϕ

⊗

∂

∂t

out

. (10.23)

Despite the fact that the operators S

and S

given by (10.22) and (10.23), respec-

tively, and determining the interaction of the ﬁelds

and ϕ

out

, are identical, the

supports of the corresponding tomographic functionals are diﬀerent, which stems

from the fact that the velocities of propagation of the compressional (v

) and shear

) waves are diﬀerent. To illustrate the conﬁguration of the support of the to-

mographic functional, we consider the plane incoming wave. Fig. 10.3(a) shows

the characteristic cones corresponding to the ﬁeld of the point source ϕϕ

out

, that

is generated by the receiver at the moment t

and propagates in the reverse time.

The outer conical surface corresponds to the propagation with the compressional

velocity, and the inner one corresponds to the propagation with the shear velocity.

The ﬁeld is concentrated between these conical surfaces in the three-dimensional

space R

(for simplicity, we assume that the “source” is simulated by the δ-function

in time, and, therefore, in Fig. 10.3(a) there is no action of the convolution over

time). The normal to the front of the plane wave is oriented in the direction op-

posite to that of the vector e

. In Fig. 10.3(a), the ﬁeld ϕ

represents the plane

compressional wave with ﬁnite signal length. The support of the tomographic func-

tional (Fig. 10.3(a)) is limited in the space by two surfaces. The outer one is the

paraboloid of revolution that is the space projection of the section of the outer char-

acteristic cone (v

) by the plane front of the incident P -wave. The inner surface

is the ellipsoid of revolution that is the space projection of the section of the inner

characteristic cone (v

) by the plane front of the incident P -wave. The symmetry

axes of both surfaces limiting the support of the tomographic functional coincide

with the vector e

In Fig. 10.3(c), the incident ﬁeld ϕ

is the plane S-wave. In this case, the outer

conical surface limiting the support of the functional (Fig. 10.3(d)) is a hyperboloid

of revolution, whereas the inner one is a paraboloid of revolution. The symmetry

axis is e

as in the previous case.

10.3.3 The transport equation of the stationary sounding signal

The transport equation arises in studying the spatial structure of the atmosphere of

the Earth and other planets, in studying plasma objects in a wide range of physical

experiments where the measurement data are quadratic functions of wave ﬁelds,

and in describing passages of particle ﬂows through a medium.

Let the ﬁeld of the sounding signal be

ϕ(x, n) , n ∈ Ω = {n : n ∈ R

, |n| = 1}.

The transport equations for the reference (L

) and perturbed (L

) media in the

Construction of tomographic functionals 307

Fig. 10.3 The characteristic cones, corresponding to the ﬁeld of the point source

out

. The ﬁeld

, with ﬁnite signal duration, is a plane compressional wave (a), (b) and a plane shear wave (c),

(d); (b), (d) are the supports of the tomographic functional.

case of isotropic scattering are as follows:

ϕ = (n, ∇∇ϕ) + α

(x)ϕ − σ

(x)

Ω

ϕ(x, n)dn, (10.24)

ϕ = (n,

∇

∇ϕ) + α(x)ϕ − σ(x)

Ω

ϕ(x, n)dn (10.25)

where α is the absorption coeﬃcient, β is scattering coeﬃcient, Θ

||α

(x), σ

(x)||, ν(δθ) = δθ). The interaction operator has two components: S

, S

], and the tomographic functional is represented as





= [









]. The

308 STATISTICAL METHODS OF GEOPHYSICAL DATA PROCESSING

interaction operators S

and S

are easily found from the following representation

of the perturbation operator for equations (10.24) and (10.25):

δL = δL

+ δL

δL

: δL

ϕ = δαϕ , (10.26)

δL

: δL

ϕ = δσ

Ω

ϕ(x, n

) dn

. (10.27)

The tomographic functionals are written in the form



|δα



Ω

out

δL

dn ,



|δσ



Ω

out

δL

dn .

Substituting δL

in the form (10.26) and δL

in the form (10.27), we obtain



(x)





out

|ϕ



Ω



out

|ϕ



Ω

, (10.28)







out

|ϕ



Ω



out

Ω

(x, n

) dn



Ω

. (10.29)

Examining formulas (10.28) and (10.29), we see that, in the isotropic case, the

action of the tomogarphic functional corresponding to the absorption coeﬃcient

α, reduces to the integration of the product of the ﬁeld ϕ

and ϕ

out

over the

full solid angle, whereas the action of the tomographic functional of the scattering

coeﬃcient is determined by the product of the integrals of the ﬁelds ϕ

and ϕ

out

over the full solid angle. We note that the ﬁeld ϕ

, occurring in the expressions for

the tomographic functional is the solution of the transport equation (10.24) in the

reference medium. The ﬁeld ϕ

out

satisﬁes the conjugate equation

−(n,

∇

∇ϕ) + α

(x) − σ

(x)

Ω

ϕ(x, n

) dn

= h ,

where h is the ﬁeld generated by the source.

10.3.4 The diﬀusion equation

In the geoelectrics it is common to consider the quasistationary ﬁelds E and H. In

the case of the absence of charges, these ﬁelds satisfy the homogeneous diﬀusion

equation

∆ϕ − σ

∂

∂t

ϕ = 0 ,

Construction of tomographic functionals 309

where

ϕ coincides with either E or H. The evolutionary operators in the reference

) and perturbed (L

) media are as follows

= ∆ − σ

(x)

∂

∂t

= ∆ − σ (x)

∂

∂t

where σ(x) = σ

(x) + δσ(x). The perturbation operator δL

can be represented as

δL

= δσ∂/∂t. Then the tomographic functional takes the form

= ϕ

⊗

∂

∂t

out

. (10.30)

We note that ϕϕ

and ϕ

out

occurring in (10.30) satisfy the equation

ϕϕ

out

= 0 , L

∗

= h ,

where

∗

= ∆ + σ

(x)

∂

∂t

10.4 Ray Tomographic Functional in the Dynamic and Kinematic

Interpretation of the Remote Sounding Data

The interpretation of the data of remote sounding can be carried out by two ways.

As the basic data for the tomography experiment can be used the phase character-

istic of the sounding signal (in this case we deal with the ray tomography) or the

full ﬁeld of the sounding signal (it corresponds to the diﬀraction tomography). In

the latter case the wave ﬁelds structure carries an information about space domain

which is limited by kinematic reason, instead of that in the ray tomography, the

information about medium properties is located in the vicinity of the bent ray. The

kinematic restriction in the diﬀraction tomography appears because of a ﬁnite time

function of the source, if before the switching-on of the source, the medium has been

in rest. It is possible to develop the tomographic interpretation using the station-

ary wave ﬁelds. It must take into account both kinematic and dynamic restrictions,

moreover, the interference principle of the wave ﬁeld shaping allows us to construct

its spectral amplitude as a superposition of the ﬁelds having the ray structure.

Let us show forming of wave amplitude spectrum on the example of acoustic

ﬁeld excited by a point source located inside the uniform layer with thickness h and

satisﬁed the homogeneous boundary conditions on the upper and lower bounds:

(∇

+ k

)ϕ = −δ(z − z

)

δ(ρ)

2πρ

, (10.31)



z=0

= 0 , ϕ



z=h

= 0 .

Let us consider the solution (10.31) without the consideration of the boundary

conditions, i.e. ϕ

(x) = e

ikl

/(4πl) (l — the ray segment), which joins the source

310 STATISTICAL METHODS OF GEOPHYSICAL DATA PROCESSING

point (0,z

) and the observation point x = (ρ, z)). The rays are reﬂected by the

upper bound z = 0 and the reﬂection is determined by the boundary condition

ϕ|

z=0

= 0. In this case the reﬂection coeﬃcient is equal −1. The reﬂected ray can

be represented as the ray generated by some virtual source with a mirror position

relatively to the plane z = 0 in the point z = −z

. Then the reﬂected wave can

be described by the formula e

ikl

/(4πl), where l is the ray segment between virtual

source and observation point. In the case of lower bound, the reﬂection coeﬃcient is

equal +1, and the virtual source is located in the point 2h −z

. Multiple reﬂection

of the waves can be described using a set of rays, when each ray has an own virtual

source located in the point 2h

+ z

, where n = 0, ±1, ±2, . . . . The wave ﬁeld is

the sum of the incident wave and multiple reﬂected waves. Each reﬂected wave is

the spherical wave, generated by appropriate virtual source (mirror-image method):

ϕ(x) =

4π

∞

n=−∞

(−1)



ikl

↑

−

ikl

↓



, (10.32)

where

↑

= [ρ

+ ((z − z

) + 2hn)

]

1/2

is the ray, which goes from the source toward the upper bound;

↓

= [ρ

+ ((z + z

) + 2hn)

]

1/2

is the ray, which goes from the source toward the upper bound.

Let us show as simplest ray expansion (10.32) is connected with more compli-

cated approximation of the wave ﬁeld vertically inhomogeneous medium. Taking

into account the axial symmetry of the problem, we write the wave equation for the

Hankel image of the ﬁeld ˜ϕ(ka, z):

˜ϕ(ka, z) = 2π

∞

(kaρ)ϕ(ρ, z)ρdρ ,

where J

(akρ) is the Bessel function of zero order;

2π

∞

(kaρ)

∂

∂ρ



∂ϕ

∂ρ



dρ +

∂

ϕ(ka, z)

∂z

+ k

(z) ˜ϕ(ka, z) = δ(z − z

) ; (10.33)

˜ϕ(ka, 0) = 0 ; ∂ ˜ϕ/∂z − (ka, h) = 0 ;

n(z) is the refraction index. Taking into account that

2π

∞

(kaρ)

∂

∂ρ



∂ϕ

∂ρ



dρ = −k

˜ϕ(ka, z)

and we obtain another form for the equation (10.33):

∂

˜ϕ(ka, z)

∂z

+ k

(z) − a

] ˜ϕ(ka, z) = −δ(z − z

) (10.34)

Construction of tomographic functionals 311

with the same boundary condition. In the case of the solution of the equation

(10.34), which satisﬁes the ﬁrst boundary condition, we introduce the notation

˜ϕ

(ka, z

). And for the solution which satisﬁes the second boundary condition we

introduce the notation ˜ϕ

(ka, z

), then we can write general solution in the form

˜ϕ(ka, z) = ˜ϕ

(ka, z

) ˜ϕ

(ka, z

)/W (ka) , (10.35)

where

W =



˜ϕ

∂ ˜ϕ

∂z

∂ ˜ϕ

∂z



is the Wronskian, the symbols z

, z

mean, that z > z

and z < z

correspondingly.

We represent ˜ϕ

(ka, z

) and ˜ϕ

(ka, z

) using up-going and down-going wave

ﬁelds

˜ϕ

(ka, z

) = U(ka, z

) + r

(ka, z

) D(ka, z

) ,

˜ϕ

(ka, z

) = D(ka, z

) + r

(ka, z

) U (ka, z

) .

Where U is the solution of the homogeneous equation for up-going waves in the

vicinity of z

; D is down-going wave. At that the normalization of U and D is

chosen as follows

W [U, D] = −2ik ,

U(ka, 0)

D(ka, 0)

, r

∂

D(ka, h)

∂

U(ka, h)

are the reﬂection coeﬃcients at upper and lower bounds.

Taking into account that

W ( ˜ϕ

, ˜ϕ

) = (1 − r

)W (U, D) = 2ik(1 − r

)

and using the Hankel image ˜ϕ(ka, z) (10.35) we can write the ﬁnal representation

of the wave ﬁeld

ϕ(ρ, z) =

4π

∞

(kaρ)[U (ka, z

) + r

D (ka, z

)]

× [D (ka, z

) + r

U (ka, z

)](1 − r

)

−1

a da .

Factoring of the binomial (1 −r

)

−1

, we obtain the multiple scattering represen-

tation (by analogy with the ray representation)

ϕ(ρ, z) =

∞

n=0

4π

∞

(ka, ρ)[(r

)

U (ka, z

)D (ka, z

)

+ (r

)

n+1

U (ka, z

) D (ka, z

) + r

n+1

D(ka, z

)

× D(ka, z

) + r

n+1

U(ka, z

)U(ka, z

) ]a da. (10.36)

Here the ﬁrst addend describes the wave which come out from the source in upward

direction, then after n reﬂections for each bound it propagates as down-going wave;

312 STATISTICAL METHODS OF GEOPHYSICAL DATA PROCESSING

the second addend describes similar wave with (m + 1) reﬂections for each bound;

the third and the fourth addends are interpreted as the wave, which goes out from

the source in deﬁnite direction, then after n reﬂections from one bound and n + 1

reﬂections from other one it propagates in the opposite direction. In the consid-

ered cases as reﬂection we understand the refraction phenomenon also, appropriate

turning point is described by the equality n

(z) = a

(see Eq. (10.34)).

The wave ﬁeld (10.36) is generated by the point source, located in the point

ρ = 0 and z = z

, i.e. the solution is the Green function. The ray interpretation

of the Green function, considered above, is based on the ray expansion of the wave

ﬁeld:

G(x, x

) =

(x, x

) exp[iωτ

(x, x

)] ,

where g

is the amplitude of the wave ﬁeld for j-th ray; τ

is the eikonal.

10.5 Construction of Incident and Reversed Wave Fields

in Layered Reference Medium

The tomographic functional (10.10) includes the wave ﬁelds ϕ

and ϕ

out

, which are

the solutions of the direct problems. In the case of the arbitrary reference medium

such solutions can be connected with the great mathematical and computing diﬃ-

culty. As usual, the solution of the direct problem can be simpliﬁed if the reference

medium has some type of the symmetry. The most symmetric is, obviously, uniform

inﬁnite medium allows an analytic representation of the wave ﬁelds ϕ

and ϕ

out

The next class of the medium is a class of piecewise-layered media. Such media have

a wide spread occurrence owing to shell structure of the Earth and near-Earth space.

Under the local description the spherical symmetry layers can be approximated by

plane parallel layers. Let us describe the algorithms for calculation of the ﬁelds ϕ

and ϕ

out

in the case of the layered medium (Brekhovskikh, 1960; Petrashen et al.,

1985; Kennett, 1983).

The problem of propagation of elastic, acoustic and electromagnetic waves with

zero-initial conditions can be reduced to the solution of the diﬀerential matrix equa-

tion of the ﬁrst order (Ursin, 1983; Ursin and Berteussen, 1986):

∂

∂x

˜ϕ = β

L ˜ϕ = β



˜ϕ



. (10.37)

So, if the equation Lϕ = 0 describes propagation of elastic waves in an isotropic

medium with the parameters ρ, λ, µ, then the Fourier transform by lateral coordi-

nates (x

, x

) and the temporal frequency ω:

ϕ(ω, p

, p

, x

) =

∞

−∞

dt dx

ϕ(t, x

, x

)

× exp[iωt − ip

− ip

] . (10.38)