Inverse Vening Meinesz formula and de¯ection-geoid formula:
applications to the predictions of gravity and geoid
over the South China Sea
C. Hwang
Department of Civil Engineering, National Chiao Tung University, 1001 Ta Hsueh Road, Hsinchu 30050, Taiwan, ROC Fax: +886 3 5716257; e-mail: hwang@geodesy.cv.nctu.edu.tw
Received: 7 April 1997 = Accepted: 7 January 1998
Abstract. Using the spherical harmonic representations of the earth's disturbing potential and its functionals, we derive the inverse Vening Meinesz formula, which converts de¯ection of the vertical to gravity anomaly using the gradient of the H function. The de¯ection-geoid formula is also derived that converts de¯ection to geoidal undulation using the gradient of the C function. The two formulae are implemented by the 1D FFT and the 2D FFT methods. The innermost zone eect is derived. The inverse Vening Meinesz formula is em-ployed to compute gravity anomalies and geoidal undulations over the South China Sea using de¯ections from Seasat, Geosat, ERS-1 and TOPEX//POSEIDON satellite altimetry. The 1D FFT yields the best result of 9.9-mgal rms dierence with the shipborne gravity anomalies. Using the simulated de¯ections from EGM96, the de¯ection-geoid formula yields a 4-cm rms dierence with the EGM96-generated geoid. The predicted gravity anomalies and geoidal undulations can be used to study the tectonic structure and the ocean circulations of the South China Sea.
Key words De¯ection-geoid formula Gravity anomaly Inverse Vening Meinesz formula Satellite altimetry South China Sea
1 Introduction
Since the publication of the Vening Meinesz formula (Vening Meinesz 1928), little attention has been paid to its inverse formula, which converts de¯ections of the vertical to gravity anomalies. This is mainly because measurements of de¯ections are not widely available. With the advent of satellite altimetry, however, de¯ections of the vertical become available in the oceans and the inverse Vening Meinesz formula can be useful if one wishes to compute
marine gravity from satellite altimetry. Marine de¯ections of the vertical can be derived from altimeter-measured geoidal undulations (if the sea-surface topography is properly removed) and the use of de¯ection as data type can reduce many systematic errors in satellite altimetry (Hwang 1997). Indeed, a frequency-domain version of the inverse Vening Meinesz formula exists in the literature, e.g., Haxby et al. (1983), Hwang and Parsons (1996), Sandwell and Smith (1997). A space-domain version has also been derived in, e.g., Molodenskii et al. (1962, Eq. III. 2.11). This paper attempts to derive the inverse Vening Meinesz formula for all cases using a spectral representa-tion approach. The de¯ecrepresenta-tion-geoid formula, which con-verts de¯ections of the vertical to geoidal undulations, will be also derived using the same approach. Practical methods for implementing the two formulae will be presented. As an example, the two formulae will be employed to compute the gravity anomalies and the geoidal undulations over the South China Sea using the de¯ections of the vertical from Seasat, Geosat, ERS-1 and TOPEX/POSEIDON satellite altimetry.
2 Fundamentals
First we brie¯y review some of the basic equations in physical geodesy necessary for the derivations of the inverse Vening Meinesz formula and the de¯ection-geoid formula. The earth's disturbing potential T can be expanded into a series of spherical harmonics as T r; /; k GMr X1 n2 R r nXn m0 X1 a0 Ca nmYnma /; k 1
where r; /; k are the spherical coordinates (geocentric distance, geocentric latitude and longitude), R is the earth's mean radius, Ca
nmare the harmonic coecients and Ynma are
the fully normalized spherical harmonics de®ned as Ya
nm /; k RSnm /; k Pnm sin / cos mk; a 0 nm /; k Pnm sin / sin mk; a 1
with Pnm sin / being the fully normalized associated
Legendre function (Heiskanen and Moritz 1967). On the sphere of radius R, geoidal undulation can be expressed as N /; k T c0jrR R X1 n2 Xn m0 X1 a0 Ca nmYnma /; k 3
and gravity anomaly as Dg /; k ÿoTor ÿ2rT jrR c0 X1 n2 n ÿ 1 Xn m0 X1 a0 Ca nmYnma /; k 4 where c0GMR2 5
is the mean gravity. Spherical harmonic representations of other functionals of the earth's disturbing potential can be found in Rummel and van Gelderen (1995). Geoidal undulation and gravity anomaly are two scalar functions derived from the earth's disturbing potential. De¯ection of the vertical, however, is a vector function and can be expressed as
z /; k n g ÿ1 RrN /; k ÿX1 n2 Xn m0 X1 a0 Ca nmrYnma /; k 6 where n and g are the north-south component and west-east component of the de¯ection vector, respectively, and r is the gradient operator on the sphere de®ned as:
r o o/ o cos /ok 7 Thus the basis functions of the de¯ection vector are rYa
nm, rather than Ynma .
Finally we recall a variant of Green's formula (Meissl 1971, p. 12) ZZ rrf rgdr ÿ ZZ rf D gdr ÿZZ rgD f dr 8
where f and g are two arbitrary functions de®ned on the unit sphere, and D is the Laplace surface operator
(Courant and Hilbert 1953) de®ned as Dcos /1 o/o cos /o/o cos12/ o2 ok2 9 For the surface spherical harmonics, we have
DYa
nm /; k n n 1Ynma /; k 0 10
Using Eqs. (8) and (10) and the orthogonality relation-ship of fully normalized spherical harmonics (Heiskanen and Moritz 1967, p. 31), we obtain the result
ZZ rrY a nmrYsrbdr 4pn n1; if n s and m r and a b 0; if n6 s or m 6 r or a 6 b 11
3 The inverse Vening Meinesz formula
The key to deriving the inverse Vening Meinesz formula is to look for a suitable kernel function for converting de¯ection of the vertical to gravity anomaly. Based on Meissl's (1971) approach, we introduce the kernel function H de®ned as H wpq X1 n2 2n 1 n ÿ 1 n n 1 Pn cos wpq 12 where p and q are two points on the unit sphere with a spherical separation of wpq, so that (see Fig. 1)
cos wpqsin /psin /q cos /pcos /qcos kqÿ kp 13
The Legendre polynomial Pn cos wpq in Eq. (12) can be
further decomposed into the series: Pn cos wpq2n 11 Xn m0 X1 a0 Ya nm /p; kpYnma /q; kq 14
which is termed decomposition formula by Heiskanen and Moritz (1967), see also Hobson (1965). Integrating the scalar products of rH and rN over the unit sphere and using Eq. (11), we getZZ
rrqH wpq rqN qdrq R ZZ r X1 n2 n ÿ 1 n n 1 Xn m0 X1 a0 Ya nm prqYnma q " # X1 n2 Xn m0 X1 a0 Ca nmrqYnma q " # drq RX1 n2 n ÿ 1 n n 1 Xn m0 X1 a0 Ca nmYnma p ZZ r rqY a nm q rqYnma q drq 4pRX1 n2 n ÿ 1Xn m0 X1 a0 Ca nmYnma p 15
Fig. 1. Spherical distance wpqbetween points p and q, and components
Comparing Eqs. (4) and (15), we have Dg p c0
4pR ZZ
rrqH wpq rqN qdrq 16
As shown in the Appendix, the closed form of H is H wpq 1 sinwpq 2 log sin3 w2pq 1 sinwpq 2 ! 17 Furthermore, rqH wpq dwdH pqrqwpq dH dwpq owpq o/q owpq cos /qokq ! 18 With Eq. (17) the derivative of H with respect to wpq is
H0 dH dwpq ÿ coswpq 2 2 sin2 wpq 2 cos wpq 2 3 2 sinw2pq 2 sinwpq 2 1 sinw2pq 19 Dierentiating Eq. (13) with respect to /q and kq, we
have (cf. Heiskanen and Moritz 1967, p. 113) ÿ sin wpqowo/pq
q cos /qsin /pÿ sin /qcos /pcos kqÿ kp
ÿ sin wpqowokpq
q ÿ cos /pcos /qsin kqÿ kp
20 Referring to the spherical triangle in Fig. 1, the following relationships hold:
sin wpqcos aqp cos /qsin /pÿ sin /qcos /pcos kqÿ kp
sin wpqsin aqp ÿ cos /psin kqÿ kp 21
Thus owpq
o/q ÿ cos aqp;
owpq
okq ÿ cos /qsin aqp 22
Inserting Eqs. (18) and (22) into Eq. (16) we ®nally get the inverse Vening Meinesz formula
Dg p ÿ4pc0 ZZ rH 0 owpq o/q owpq cos /okq ! nqgqdrq c0 4p ZZ rH 0 n qcos aqp gqsin aqp ÿ drq 4pc0 ZZ rH 0e qpdrq 23
where eqp is the de¯ection component at point q in the
direction of the azimuth aqp (Heiskanen and Moritz
1967, p. 187), or simply the longitudinal de¯ection component. The meaning of the inverse Vening Meinesz formula is: assuming that everywhere on the unit sphere the north-south and west-east de¯ection components are known, the gravity anomaly at any given point can
be obtained by integrating the products of H0 and the
longitudinal de¯ection components over the unit sphere.
Figure 2 shows the function H0, which changes
rap-idly as w approaches zero. The ®rst zero crossing of H0
occurs at w 43. When w is small, we have the
as-ymptotic representation: H0 w ÿ 2
w2 24
It is noted that the asymptotic representation of H0 is
equal to that of dS=dw (the derivative of Stokes' function) and agrees with the asymptotic representation of the kernel function in Eq. (III.2.111) of Molodenskii et al. (1962). Furthermore, H w 2
w as w approaches
zero, so the asymptotic representations of the function H and Stokes' function are identical.
4 The de¯ection-geoid formula
Next we derive a formula for converting de¯ection of the vertical to geoidal undulation. The derivation is almost the same as that for the inverse Vening Meinesz formula. First, we introduce the kernel function C C wpq X1 n2 2n 1 n n 1Pn cos wpq 25 Then, ZZ rrqC wpq rqN qdrq R ZZ r X1 n2 1 n n 1 Xn m0 X1 a0 Ya nm prqYnma q " # X1 n2 Xn m0 X1 a0 Ca nmrqYnma q " # drq
RX1 n2 1 n n 1 Xn m0 X1 a0 Ca nmYnma p ZZ rrqY a nm q rqYnma qdrq 4pRX1 n2 Xn m0 X1 a0 Ca nmYnma p 26
Thus, the geoidal undulation at point p can be obtained by integrating the scalar products of rH and rN over the unit sphere:
N p 4p1 ZZ
rrqH wpq rqN qdrq 27
Using Eqs. (22) and (27) we obtain the de¯ection-geoid formula N p 4pR ZZ r dC dwpqeqpdrq 28
According to the Appendix, the closed form of C is C wpq ÿ2 log sinw2pqÿ32cos wpqÿ 1 29
and C0 dC dwpq ÿ cot wpq 2 3 2sin wpq 30
which agrees with the result of Molodenskii et al. (1962, Eq. III.2.8) if the summation of the series in Eq. (25) starts from n 1. Figure 2 also shows the funciton C0.
The ®rst zero crossing of C0 occurs at w 70:5. The
asymptotic representation of C0 when w is small is:
C0 w ÿ2
w 31
5 Computations by 1D FFT: rigorous implementations We propose two computational schemes for the inverse Vening Meinesz formula and the de¯ection-geoid for-mula when the north-south and west-east components of de¯ections are given on a regular grid. The ®rst scheme is based on the one-dimensional fast Fourier transform (1D FFT) method, which, given regularly gridded data, can rigorously implement a surface integral such as Eq. (23) or (28). In such a scheme, gravity anomalies or geoidal undulations at the same parallel are computed simultaneously by FFT as (cf. Haagmans et al. 1993) Dg/p kp N/p kp ( ) c0 R D/Dk 4p X/n /q/1 Xkn kqk1 H0 Dkqp C0 Dk qp
ncoscos aqp gcossin aqp
c0 R D/Dk 4p Fÿ11 ( X/n /q/1 F1 H0 Dkqp cos aqp F1 C0 Dkqp cos aqp F1 ncos F1 H0 Dkqp sin aqp F1 C0 Dkqp sin aqpF1 gcos ) where ncos n cos /; gcos g cos /; Dkqp kqÿ kp; D/
and Dk are grid intervals in the directions of latitude and longitude, and F1is the 1D FFT. Since all quantities
are real-valued, we can compute the Fourier transforms of two real-valued arrays simultaneously to save computer time (Hwang 1993). Speci®cally, let h k and g k, k 0; . . . ; n ÿ 1, be the two real-valued arrays to be Fourier transformed. We ®rst form the complex array y k as
y k h k i g k; k 0; . . . ; n ÿ 1 33 where i pÿ1. Let Y k be the Fourier transform of y k, we have H 0 Re Y 0; and H k 1 2 Re Y k Y n ÿ k 1 2i Im Y k Y n ÿ k for k 1; . . . ; n ÿ 1 G 0 Im Y 0; and G k 1 2Re Y k Y n ÿ k ÿ1 2i Im Y k Y n ÿ k for k 1; . . . ; n ÿ 1 34
where H k and G k are the Fourier transforms of h k and g k; respectively, and Re(.) and Im(.) are the real and the imaginary parts of a complex number. In practice, the complex array holding ncosand gcos, and the
complex array holding H0 Dk
qp cos aqp and
H0 Dk
qp sin aqp (or C0 Dkqp cos aqp and C0 Dkqp
sin aqp are Fourier transformed. Taking advantage of
the gridded data, azimuth and spherical distance can be calculated as
tan aqp ÿ cos /psin Dkqp
ÿ sin /qÿ /p 2 sin /qcos /psin2 Dk2qp
35 sin2 wqp 2 sin2 D/qp 2 sin2 Dkqp 2 cos /qcos /p 36 where D/qp /qÿ /p.
6 Computations by 2D FFT: planar approximations The second computational scheme is based on the planar approximations of the two formulae, and hence the two-dimensional fast Fourier transform (2D FFT). In a local rectangular, x ÿ y coordinate system, the surface element and the spherical distance can be approximated as R2dr
q dxqdyq and wqp qRqp, with
qqp being the planar distance. With the asymptotic
representation in Eq. (24), the inverse Vening Meinesz formula becomes
Dg xp; yp 4pc0 ZZ D ÿ2R2 q3 qp nq qqpcos aqp gq qqpsin aqp dxRqdy2 q ÿc0 2p ZZ D ypÿ yq xpÿ xq2 ypÿ yq2 3 2nqdxqdyq ( ZZ D xpÿ xq xpÿ xq2 ypÿ yq2 3 2gqdxqdyq ) ÿc2p0 y x2 y232 ! n x x2 y232 ! g ( ) 37 where * is the convolution operator and D is the data domain. Schwarz et al. (1990) show that
F xy x2 y2ÿ3=2 ÿ2pi uv u2 v2ÿ3=2 38
where F is the 2D FT, and u and v are the spatial frequencies. Thus the relationship between gravity anomaly and de¯ection components in the frequency domain is
DG u; v ic0
u2 v2
p vX u; v uE u; v 39
where DG; X and E are the Fourier transforms of Dg, n, g, respectively. Equation (39) can also be found in, for example, Hwang and Parsons (1996), Sandwell and Smith (1997), Haxby et al. (1983), who derived this formula using dierent approaches.
The planar approximation of the de¯ection-geoid formula reads N xp; yp 4pR ZZ D ÿ2R q2 qp nq qqpcos aqp gq qqpsin aqp dxRqdy2 q ÿ12p ZZ D ypÿ yq xpÿ xq2 ypÿ yq2nqdxqdyq ( ZZ D xpÿ xq xpÿ xq2 ypÿ yq2gqdxqdyq ) ÿ12p x2 yy 2 n x2 yx 2 g 40 The Fourier transforms of x
x2y2and x2yy 2, which do not
exist in the literature, are now derived. We begin with the de®nite integral found in Gradshteyn and Ryzhik (1994, p. 782): Z 1 0 log 1 a2 r2 p ÿ log rJ0 brrdr 1 b2 1 ÿ eÿab; a > 0; b > 0 41
Let a in Eq. (41) approach in®nity so that 1 pa2 r2 K remains a constant for any r, and
eÿab 0. Then, Z 1 0 log K r J0 brrdr b12 42
This means that the Hankel transform of log K
r is 2pq12, with q2 u2 v2. Furthermore, x y 1 x2 y2 xy 1 r2 ÿ o ox o oy ! log Kr 43 By the dierentiation theorem of Fourier transform (see, e.g., Mesko 1984), we have
F xy 1 x2 y2 ÿi uv 1 u2 v2 44
With Eqs. (40) and (44), we obtain the relationship between geoidal undulation and de¯ection components in the frequency domain
N u; v 2p u2i v2vX u; v uE u; v 45 which can also be found in Olgiati et al. (1995, Eq. 6). Using Eq. (39) and the geoid-gravity spectral relation-ship [see, e.g., Schwarz et al. (1990)], one can also derive Eq. (45).
7 The innermost zone eects
At zero spherical distance the kernel function H0and C0
become singular and the azimuth is unde®ned. Thus we must account for the innermost zone eect (Heiskanen and Moritz 1967). First we consider such an eect in the inverse Vening Meinesz formula. The de¯ection compo-nents at the neighbourhood of point p can be expanded into the Taylor series (see Fig. 3)
Fig. 3. Local rectangular coordinates for the innermost zone eect in a cap of radius s0
nq np xnx yny12ÿx2nxx 2xynxy y2nyy gq gp xgx ygy12ÿx2gxx 2xygxy y2gyy 46 where nxonox; ny onoy; nxxo 2n ox2; nyyo 2n oy2 and nxy o 2n oxoy.
Retaining only the linear terms in Eq. (46) and assuming that the innermost zone is circular, with Eq. (24) we have Dgi4pc0 Z 2p a0 Z s0 s0 ÿ2 s2 np s sin a nx s cos a ny
cos a p gp s sin a gx s cos a gy
sin a psdsda c0 2p ny Z 2p a0cos 2adaZ s0 s0ds gx Z 2p a0sin 2adaZ s0 s0ds s02c0 ny gx 47
Using a similar derivation, the innermost zone eect in the case of the de¯ection-geoid formula is
Nis 2 0
4 ny gx 48
Thus the innermost zone eects for gravity anomaly and geoidal undulation depend on the gradients of the de¯ection components. For discrete data ny and gx can
be obtained by numerically dierentiating n and g along the y and x directions, respectively. If the planar grid intervals are Dx and Dy, the radius of the innermost zone may be approximated by s0 DxDy p r 49
8 Applications: gravity and geoid over the South China Sea from satellite altimetry
As an application of the inverse Vening Meinesz formula, we computed marine gravity anomalies over the South China Sea (de®ned domain: 5 latitude
25, 105 longitude 125 using the 1D FFT and
2D FFT methods and the remove-restore procedure. The EGM96 geopotential model (Lemoine et al. 1997) to degree 360 was used as the reference ®eld. The altimeter data used are from Seasat, Geosat/ERM, Geosat/GM, ERS-1/35-day, ERS-1/GM and TOPEX/ POSEIDON. The sea-surface topography of Levitus (1982) is subtracted from the altimeter sea surface heights before generating the de¯ections of the vertical. At the centre of the South China Sea the average altimeter data density is 1560 points in 1 1, and at
the continental borders the densities drop sharply. We used the method of least-squares collocation and the covariance functions of de¯ections derived by Hwang
and Parsons (1995) to grid the de¯ections of dierent azimuths into the north-south and west-east compo-nents at a 20 20interval. We used a 1-border to avoid
bad results at the edges. Further, 100% zero paddings were applied to data arrays and kernel arrays to avoid edge eects in convolutions by FFT. Table 1 shows the comparisons between the shipborne gravity anomalies and the gravity anomalies derived from the 1D FFT, the 2D FFT, and Sandwell and Smith's (1997) methods. The shipborne gravity anomalies were provided by the National Geophysical Data Center (NGDC). A total of 180297 shipborne gravity anomalies were used for the comparisons. The 1D FFT produces a slightly better result than the 2D FFT. Table 1 also shows the CPU time ratio between a given method and the 2D FFT method on a Sun Sparc 20 machine. The 1D FFT requires more than doubled computer time than the 2D FFT. Employing the innermost zone eect improves the result. We also did tests in other areas, and found that the innermost zone eect always improves the result. Furthermore, the result from the 1D FFT is about 30% better than the recently published altimeter-derived gravity anomalies from Sandwell and Smith (1997). Figure 4 shows the predicted gravity anomalies over the South China Sea (1D FFT plus innermost zone eect). In Fig. 4, the outline of the basin of the South China Sea is clearly visible. A median valley-like feature running from the southwest to the northeast is probably the spreading centre of the South China Sea, which has been identi®ed by, e.g., Briais et al. (1993), using geomagnetic data. The predicted gravity anomalies can be further used to interpret the tectonic structure of the South China Sea.
Because generally there are no measured geoidal undulations at sea, we used the simulation approach of Tziavos (1996) to evaluate the performances of the 1D and 2D FFT methods that implement the de¯ection-geoid formula. First, over the South China Sea we generated gridded north-south and west-east compo-nents of de¯ections, and geoidal undulations using harmonic coecients from EGM96 (Lemoine et al. 1997) at a 7:50 7:50 interval. Only harmonic
coe-cients between degrees 181 and 360 were used, so that the remove-restore procedure need not be used. Be-cause going from de¯ection to geoid is a smoothing process, we used a 5 border to avoid edge eects. The
north-south and west-east de¯ection components were then used to compute geoidal undulations. The statis-Table 1. RMS dierences between the shipborne and the predicted gravity anomalies over the South China Sea and the CPU time ratios (aIE: innermost zone eect)
Method RMS
dierence (mgal)
CPU time ratio
1D FFT/IEa 9.90 2.8
1D FFT/no IE 10.06 2.4
2D FFT 10.11 1.0
tics of the dierences between the EGM96-generated undulations, considered as the ``ground truth'', and the computed undulations from various methods are given in Table 2. From Table 2, the 1D FFT again gives a better result than the 2D FFT. Also, the use of the innermost zone eect signi®cantly reduces the dier-ence between the `true' and the predicted geoids. Having done these tests, we then used the de¯ection data as used in predicting the gravity anomalies to compute the geoidal undulations over the South China Sea. The result is shown in Fig. 5. The predicted geoid can be used to study the ocean circulations over the South China Sea, which recently have received con-siderable attention from the oceanographers in South-east Asia.
9 Conclusion
In this paper we showed the detailed derivations of the inverse Vening Meinesz formula and the de¯ection-geoid formula using the spherical harmonic representa-tions of the functionals of the earth's disturbing potential, and for each we presented the 1D FFT and 2D FFT methods for computations. In all cases the 1D FFT yields better results than 2D FFT, but the former needs nearly doubled computer times. Over the South China Sea, by the inverse Vening Meinesz formula and the 1D FFT we derived a set of gravity anomalies better than that from Sandwell and Smith (1997) when comparing with the shipborne gravity anomalies. Using the simulated de¯ections from EGM96 over the South Fig. 4. Grey-shaded relief map of the pre-dicted gravity anomalies over the South China Sea, with illumination from the north-west
Table 2. Statistics of the EGM96 geoid and the dierences between the predicted and the EGM96 geoids over the South China Sea (unit: m;
aIE: innermost zone eect)
Case mean min. max. std. dev. RMS
EGM96 deg 181 to 360 0.000 )2.998 3.812 0.542 0.542
EGM96 - 1D FFT/IEa 0.038 )0.010 0.100 0.014 0.041
EGM96 - 1D FFT/no IE 0.038 )0.055 0.169 0.021 0.043
China Sea, the de¯ection-geoid formula yields a 4-cm accuracy. The predicted gravity anomalies and geoidal undulations over the South China Sea are freely available to all scientists. Interested readers please send e-mail to hwang@geodesy.cv.nctu.edu.tw.
Acknowledgements. This research is funded by the National Science Council of Republic of China, under contract NSC86- 2611-M-009-001-OS, with the project title: Surveying gravity and bathy-metry over the South China Sea using satellite altibathy-metry. The au-thor is grateful to Roger Haagamans and an anonymous reviewer for their important suggestions.
Appendix: derivations of closed forms of H and C functions
Consider the generating function of Legendre poly-nomials (Hobson 1965) f k 1 1 ÿ 2k cos w k2 p X1 n0 knP n cos w A1
The series in Eq. (A1) is absolutely and uniformly convergent if k < 1 (Hotine 1969, p. 310). The case k 1 will be of conditional convergence except at w 0. In the following derivations, we exclude the point w 0 in all results. Setting k 1 in Eq. (A1), we have
X w X1 n0 Pn cos w 1 2 ÿ 2 cos w p 1 2 sinw2 A2
Using Eq. (A1) and the fact that P0 cos w 1, we have
f k ÿ 1 k X1 n1 knÿ1P n cos w A3
Integrating Eq. (A3) with respect to k between the limits k, 0, and using the result in Gradshteyn and Ryzhik (1994, p. 101), we have X1 n1 Z k 0 k nÿ1dkP n cos w X1 n1 1 nknPn cos w Z k 0 f k ÿ 1 k dk A4
ÿ log 2p1 ÿ 2k cos w k2 2 ÿ 2k cos w log 4
Setting K=1, we get Y w X1
n1
1
nPn cos w ÿ log sinw2 1 sinw2
A5 Furthermore, integrating Eq. (A1) with respect to k between the limits k, 0 and using the result in Gradshteyn and Ryzhik (1994, p. 99), we have
Fig. 5. The predicted geoid over the South China Sea, contour interval is 0.5 m
X1 n0 Z k 0 k ndkP n cos w X1 n0 1 n 1kn1Pn cos w Z k 0 f kdk log 2p1 ÿ 2k cos w k2 2k ÿ 2 cos w
ÿ log 2 ÿ 2 cos w A6 Setting k 1, we have Z w X1 n0 1 n 1Pn cos w log 1 sinw 2 sinw 2 ! A7 It is easy to see that H w and C w are linear combinations of the three basic in®nite series, X w; Y w and Z w, namely, H w X1 n2 2n 1 n ÿ 1 n n 1 Pn cos w X1 n2 2 ÿ1 nÿ 2 n 1 Pn cos w 2 X w ÿ P0ÿ P1 ÿ Y w ÿ P1 ÿ 2 Z w ÿ P0ÿ12P1 1 sinw2 log sin3 w 2 1 sinw2 ! A8 and C w X1 n2 2n 1 n n 1Pn cos w X1 n2 1 n 1 n 1 Pn cos w Y w ÿ P1 Z w ÿ P0ÿ12P1 ÿ2 log sinw 2ÿ32cos w ÿ 1 A9 where P0 1; P1 cos w. References
Briaiis A, Patrait P, Tapponnier P (1993) Updated interpretation of magnetic anomalies and sea ¯oor spreading stages in the South China Sea: implications for the Tertiary tectonics of Southeast Asia. J Geophys Res 98: 6299±6328
Courant R, Hilbert D (1953) Methods of mathematical physics, vol I. Interscience, New York
Gradshteyn IS, Ryzhik IM (1994) Table of integrals, series, and products (5th edn). Academic Press, New York
Haagmans R, de Min E, van Gelderen M (1993) Fast evaluation of convolution integrals on the sphere using 1D FFT, and a comparison with existing methods for Stokes integral. Manuscr Geod 18: 227±241
Haxby WF, Karner GD, Labreque JL, Weissel JK (1983) Digital images of combined oceanic and continental data sets and their used in tectonic studies. EOS Trans Am Geophys Un 64: 995± 1004
Heiskanen WA, Moritz H (1967) Physical geodesy. Freeman WH, San Francisco
Hobson EW (1965) The theory of spherical and ellipsoidal har-monics (2nd edn). Chelsea, New York
Hotine M (1969) Mathematical geodesy. ESSA Monogr no 2, US Dept Commerce, Washington DC
Hwang C (1993) Fast algorithm for the formation of normal equations in a least-squares spherical harmonic analysis by FFT. Manuscr Geod 18: 46±52
Hwang C (1997) Analysis of some systematic errors aecting al-timeter-derived geoid gradient with applications to geoid de-termination over Taiwan. J Geod 71: 113±130
Hwang C, Parsons B (1995) Gravity anomalies derived from Sea-sat, GeoSea-sat, ERS-1 and TOPEX/POSEIDON altimetry and ship gravity: a case-study over the Reykjanes Ridge. Geophys J Int 122: 511±568
Hwang C, Parsons B (1996) An optimal procedure for deriving marine gravity from multi-satellite altimetry. Geophys J Int 125: 705±718
Lemoine FG, et al. (1997) The development of the NASA GSFC and NIMA joint geopotential model. In: Segawa J, Fwimoto H (eds) Proc Int Symp Gravity, geoid, and marine geodesy Levitus S (1982) Climatological atlas of the world ocean. NOAA
professional paper 13, US Dept Commerce. Rockville, Md Meissl P (1971) A study of covariance functions related to the
earth's disturbing potential. Rep no 151, Dept Geod Surv, Ohio State University, Columbus
Mesko AM (1984) Digital ®ltering: applications in geophysical exploration for oil. Akademiai Kiado, Budapest
Molodenskii MS, Eremeev VF, Yurkina MI (1962) Methods for study of the external gravitational ®eld and ®gure of the earth. Works of Central Research Institute of Geodesy, Aerial Pho-tography and CarPho-tography, no 131, Moscow
Olgiati A, Balmino G, Sarrailh M, Green CM (1995) Gravity anomalies from satellite altimetry: comparison between com-putation via geoid heights and via de¯ections of the vertical. Bull Geod 69: 252±260
Rummel R, van Gelderen M (1995) Meissl scheme ± spectral characteristics of physical geodesy. Manuscr Geod 20: 379±385 Sandwell D, Smith WHF (1997) Marine gravity anomaly from Geosat and ERS-1 satellite altimetry. J Geophys Res 102: 10039±10054
Schwarz KP, Sideris MG, Forsberg R (1990) The use of FFT techniques in physical geodesy. Geophys J Int 100: 485±514 Tziavos IN (1996) Comparisons of spectral techniques for geoid
computations over large regions. Manuscr Geod 70: 357±373 Vening Meinesz FA (1928) A formula expressing the de¯ection of
the plumb-line in the gravity anomalies and some formulae for gravity ®eld and the gravity potential outside the geoid. Proc Koninkl Ned Akad Wettenschap. vol 31