**2.1 Solution for the product of Moho depth and density contrast** ð Þ *D***Δ***ρ*

In H. Moritz' original publication [17] the problem is to determine the MD ð Þ *D* such that the compensation attraction (*AC*) fully compensates the Bouguer gravity anomaly. Here we employ this condition in the last part of Eq. (1), which can be written (cf. [7, 21])

$$
\delta \mathbf{g}^I = \delta \mathbf{g}^B + A\_C = \mathbf{0},\tag{2}
$$

where *δg<sup>B</sup>* is the Bouguer gravity disturbance (i.e., the free-air gravity disturbance after removal of the topographic attraction).

The VMM technique uses both gravimetric and seismic data in a least squares combination to determine the MD (*D*) and/or MDC ð Þ Δ*ρ* . The method assumes that the crust is in isostatic balance, implying that the isostatic gravity anomaly Δ*g<sup>I</sup>* � � and disturbance *δg<sup>I</sup>* � � vanish at each point on the Earth's surface as in Eq. (1) above. Note that the compensation attraction is a function of both MD and MDC. Approximating the Earth's surface by a sphere of radius *R,* one obtains after several manipulations of Eq. (2) the following equation in *D* for a constant Δ*ρ*:

$$RG\Delta\rho\left\|\left[K(\boldsymbol{\nu},s)d\sigma=f,\right.\tag{3}$$

where *G* is the gravitational constant, *K*ð Þ *ψ*, *s* is an integral kernel function with arguments *ψ* ¼geocentric angle between integration and computation points and *<sup>s</sup>* <sup>¼</sup> <sup>1</sup> � *<sup>D</sup>=R*, and *<sup>f</sup>* ¼ � *<sup>δ</sup>g<sup>b</sup>* <sup>þ</sup> *AC*<sup>0</sup> � �*=G*. Here *AC*<sup>0</sup> is zero-degree harmonic of the compensation attraction (which does not affect the Moho undulation). Eq. (3) is a

masses of the crust and its compensation along each vertical is assumed to be

*Geodetic Sciences - Theory, Applications and Recent Developments*

height is compensated by variations in the depth of the crust. That is, the mass excess of topography is compensated by the mass deficit of mountain roots in the upper mantle. In ocean areas anti-roots of mantle material compensates for the light

The P/H model assumes a constant depth of compensation of the solid Earth topography (including negative topography over oceans), while the density of the

Below we will present the least-squares theory for determining a combined VMM-seismic model for both MD and MDC. The theory is finally applied in a new

**2.1 Solution for the product of Moho depth and density contrast** ð Þ *D***Δ***ρ*

In H. Moritz' original publication [17] the problem is to determine the MD ð Þ *D* such that the compensation attraction (*AC*) fully compensates the Bouguer gravity anomaly. Here we employ this condition in the last part of Eq. (1), which can be

where *δg<sup>B</sup>* is the Bouguer gravity disturbance (i.e., the free-air gravity distur-

The VMM technique uses both gravimetric and seismic data in a least squares combination to determine the MD (*D*) and/or MDC ð Þ Δ*ρ* . The method assumes that the crust is in isostatic balance, implying that the isostatic gravity anomaly Δ*g<sup>I</sup>* � � and disturbance *δg<sup>I</sup>* � � vanish at each point on the Earth's surface as in Eq. (1) above. Note that the compensation attraction is a function of both MD and MDC. Approx-

imating the Earth's surface by a sphere of radius *R,* one obtains after several manipulations of Eq. (2) the following equation in *D* for a constant Δ*ρ*:

*σ*

where *G* is the gravitational constant, *K*ð Þ *ψ*, *s* is an integral kernel function with arguments *ψ* ¼geocentric angle between integration and computation points and

compensation attraction (which does not affect the Moho undulation). Eq. (3) is a

� �*=G*. Here *AC*<sup>0</sup> is zero-degree harmonic of the

*RG*Δ*ρ* ðð

bance after removal of the topographic attraction).

*<sup>s</sup>* <sup>¼</sup> <sup>1</sup> � *<sup>D</sup>=R*, and *<sup>f</sup>* ¼ � *<sup>δ</sup>g<sup>b</sup>* <sup>þ</sup> *AC*<sup>0</sup>

**26**

*<sup>δ</sup>g<sup>I</sup>* <sup>¼</sup> *<sup>δ</sup>g<sup>B</sup>* <sup>þ</sup> *AC* <sup>¼</sup> 0, (2)

*K*ð Þ *ψ*, *s dσ* ¼ *f*, (3)

Due to the elasticity of the Earth's crust these local models are not very realistic. Hence, [14] modified the A/H model by introducing a model with a regional compensation in which mass loads and unloads are balanced by a gentle bending or flexure of the crust over a regional area. [17] generalized the VM model from a regional to a global compensation with a spherical sea level approximation. [7] and, finally, [6] generalized the VMM model to allow for variations both in crustal density and depth. In this way the VMM can be seen as a generalization of both the A/H and P/H models with global isostatic compensations by variations of both

The A/H model assumes a constant crustal density, and variations in topographic

constant from place to place.

topography varies with topographic height.

mountain root and crustal density.

mass of the ocean.

global model.

**2. The VMM theory**

written (cf. [7, 21])

non-linear Fredholm integral equation of the first kind, which has the following first- and second-order solutions:

$$D\_1 = \frac{1}{4\pi\Delta p} \sum\_{n=0}^{\infty} \frac{2n+1}{n+1} \sum\_{m=-n}^{n} f\_{nm} Y\_{nm} \tag{4}$$

where *Ynm* is a fully-normalized spherical harmonic, *f nm* is the corresponding coefficient given by the Bouguer gravity disturbance *f*, and

$$(D\_2)\_P = (D\_1)\_P + \frac{\left(D\_1^2\right)\_P}{R} - \frac{1}{32R\pi} \iiint\_{\sigma} \left[\frac{\left(D\_1^2\right)\_Q - \left(D\_1^2\right)\_P}{\sin^3 \psi\_{PQ}}\right] d\sigma\_Q. \tag{5}$$

Here subscripts *P* and *Q* denote computation and integration points, respectively, *f nm* is the spherical harmonic coefficient of *f*. Note that the integral contributes significantly only locally around the computation point. The formula can be improved by a few steps of iteration:

$$D\_{P}^{k+1} = D\_{P}^{k} + \frac{\left(D\_{P}^{k}\right)^{2}}{R} - \frac{1}{32R\pi} \iint\limits\_{\sigma} \left[\frac{\left(D\_{Q}^{k}\right)^{2} - \left(D\_{P}^{k}\right)^{2}}{\sin^{3}\psi\_{PQ}}\right] d\sigma\_{Q}; k = 0, 1, 2, \dots,\tag{6}$$

where *D*<sup>0</sup> *<sup>P</sup>* ¼ *D*<sup>1</sup> at point *P* determined by Eq. (4).

As the isostatic balance of the crust is hardly valid for crustal blocks of diameter smaller than, say, 100 km ([32], p.195), the upper limit of the series in Eq. (4) of should not exceed *n*<sup>2</sup> =180. Also, as we shall see later, the low-degree harmonics in *D*1, say, below *n*<sup>1</sup> =10, are not contributing to the isostatic balance but are due to mass anomalies in the Earth's interior below the crust.

The integrals in Eqs. (5) and (6) are local, as the integrand quickly vanishes with distance away from the computation point. Hence a flat earth approximation may be relevant (See [6]).

If the MDC varies laterally, the following 2nd-order approximation of Eq. (3) can be found in the spectral domain (cf. [6]) when introducing the notation *χ* ¼ *D*Δ*ρ*

$$f\_{nm} = 4\pi \frac{n+1}{2n+1} \left[ \chi\_{nm} + \frac{n+2}{2R} (\chi D)\_{nm} \right] \tag{7}$$

and, after summing up, one obtains:

$$\chi = \sum\_{n=n\_1}^{n\_2} \left[ \frac{2n+1}{n+1} \frac{f\_n}{4\pi} - \frac{n+2}{2R} (\chi D)\_n \right],\tag{8}$$

where *f <sup>n</sup>* and ð Þ *χD <sup>n</sup>* are the Laplace harmonics

$$
\begin{pmatrix} f\_n \\ (\chi D)\_n \end{pmatrix} = \sum\_{m=-n}^{n} \binom{f\_{nm}}{(\chi D)\_{nm}} Y\_{nm}.\tag{9}
$$

Using the approximation

$$
\chi D \approx \chi^2 / \Delta \rho\_P,\tag{10}
$$

one obtains from (8) the iterative formula

$$\chi\_{\mathbb{P}}^{k+1} = \sum\_{n=n\_1}^{m\_1} \left[ \frac{2n+1}{n+1} \frac{f\_n}{4\pi} - \sum\_{m=-n}^{n} \frac{n+2}{2R\Delta\rho\_{\mathbb{P}}} \left(\chi^2\right)\_n^k \right]; \mathbb{k} = 0, 1, 2, \dots \tag{11}$$

where *χ*<sup>0</sup> *<sup>P</sup>* is the first-order solution:

$$\chi^0 = \frac{1}{4\pi} \sum\_{n=n\_1}^{n\_2} \frac{2n+1}{n+1} f\_n. \tag{12}$$

Where *dD* and *d*Δ*ρ* are the (unknown) corrections to the initial values *DP* and

*d*Δ*ρ*

Assuming that the observation errors are random with expectation zero and covariance matrix **Q**, the weighted least squares solution of this system becomes:

**A** � ��**<sup>1</sup>**

From this result, the adjusted MD and MDC for point *P* are obtained by:

As the first equation ð Þ *l*<sup>1</sup> is a linearization, it could make sense to iterate the adjustment procedure by replacing the previous initial values *DP* and Δ*ρ<sup>P</sup>* in Eq. (16) by their adjusted values *D*^,Δ^*ρ* and repeat the above computation procedure

First, the result of the least-squares procedure depends on the quality and weighting of the gravity and seismic observations. The weights should be selected as proportional to the inverse standard errors (STEs) of the observations squared. The STEs of seismic data is, hopefully, provided along with the data files. For the gravity data we derive the global mean STE in Section 3.1. In Section 3.2 we propagate the data errors to error estimates in the VMM least-squares results of Moho constituents. Finally, in Section 3.3 a method for validating the modeled Moho undulations is presented.

Assuming that there are no systematic errors and disregarding 2nd –order terms in Eq. (8), one obtains the error in *χ* by simple error propagation from Eq. (12):

2*n* þ 1

*<sup>n</sup>* <sup>þ</sup> <sup>1</sup> *df <sup>n</sup>*, (23)

*dcn* vuut , (24)

ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi

2*n* þ 1 *n* þ 1 � �<sup>2</sup>

X*n*2 *n*1

X*n*2 *n*¼*n*<sup>1</sup>

9 = ; <sup>¼</sup> <sup>1</sup> 4*π*

where *df <sup>n</sup>* is the error in *f <sup>n</sup>*. Then it follows that the global Root Mean Square

**3.1 The uncertainty in the gravimetric-isostatic observation equation**

*εχ* <sup>¼</sup> <sup>1</sup> 4*π*

� � and **<sup>L</sup>** <sup>¼</sup>

**A**<sup>T</sup>**Q**�**<sup>1</sup>**

*<sup>D</sup>*^ <sup>¼</sup> *DP* <sup>þ</sup> *dD*^ and <sup>Δ</sup>^*<sup>ρ</sup>* <sup>¼</sup> <sup>Δ</sup>*ρ<sup>P</sup>* <sup>þ</sup> *<sup>d</sup>*Δ^*ρ:* (22)

*<sup>P</sup>* ¼ 0, and *ε<sup>i</sup>* are the errors of the observations. In

**AX** ¼ **L** � **ε**, (19)

*l*1 *l*2 *l*3 3 7

**L***:* (21)

<sup>5</sup>*:* (20)

2 6 4

Δ*ρP*, *l*<sup>1</sup> ¼ *χ* � *χ<sup>P</sup>* � *I*

**A** ¼

2 6 4

until sufficient convergence.

**3. Uncertainty estimations**

Error (RMSE) of *χ* becomes

**29**

*RMSE*ð Þ¼ *χ*

1 <sup>4</sup>*<sup>π</sup> <sup>E</sup>* ðð

8 < :

*σ ε*2 *χ* � �*d<sup>σ</sup>*

where

*k*

*DOI: http://dx.doi.org/10.5772/intechopen.97449*

*<sup>P</sup>*, where *I*

matrix form the adjustment system can be written

3 7 0

*On Moho Determination by the Vening Meinesz-Moritz Technique*

<sup>5</sup>, **<sup>X</sup>** <sup>¼</sup> *dD*

**<sup>X</sup>**^ <sup>¼</sup> **<sup>A</sup>**<sup>T</sup>**Q**�**<sup>1</sup>**

Alternatively, we may present Eq. (11) by the iterative formula:

$$
\chi\_P^{k+1} = \chi\_P^0 + I\_P^k, \; k = 0, 1, 2, \dots, \tag{13}
$$

where

$$I\_P^k = \frac{1}{R\Delta\rho\_P} \iiint\_{\sigma} \frac{\left(\chi\_P^k\right)^2 - \left(\chi\_Q^k\right)^2}{\sin^3\psi} d\sigma\_Q. \tag{14}$$

Again, this integral is very local, which suggests the use of a flat-Earth approximation. Also, assuming that ð Þ *n*<sup>2</sup> þ 2 *D*0*=*ð Þ 2*R* <1, Eq. (7) leads to the approximate solution:

$$\chi\_P \approx \frac{1}{4\pi} \sum\_{n=n\_1}^{n\_2} \frac{2n+1}{n+1} \frac{f\_n}{1+(n+2)D\_0/(2R)}.\tag{15}$$

Note that the solution *χ<sup>P</sup>* is the product of the MD and MDC. If one of the parameters is known, the other can be determined by the equation. Hence, gravity data alone cannot be used to distinguish between the two Moho constituents. Hence, additional information, e.g., from seismic and/or geological data, is needed to separate the two. However, as we shall see later, usually such data is not taken for granted in the VMM technique, but the gravity data used in Eq. (8) is typically applied to improve a priori Moho constituents in a least-squares procedure.

The solution (8) can be derived from Eq. (1), and from the inversion of a 3-D Newton integral. See Appendix A.
