/
cilmplus.html
220 lines (157 loc) · 10.1 KB
/
cilmplus.html
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
<!DOCTYPE HTML PUBLIC "-//W3C//DTD HTML 4.01//EN"
"http://www.w3.org/TR/html4/strict.dtd">
<html>
<head>
<title>SHTOOLS - Gravity and magnetics routines</title>
<meta http-equiv="Content-Type" content="text/html; charset=UTF-8">
<link rel="stylesheet" type="text/css" href="../CSS/sh.css">
<link rel="icon" type="image/vnd.microsoft.icon" href="../Figures/favicon.ico">
</head>
<body>
<div class="main">
<p class="centeredimage"><img src="../Figures/logo.jpg" width=894 height=135 alt="SHTOOLS --- Tools for working with spherical harmonics"></p>
<table class="menu">
<tbody>
<tr>
<td><a href="http://shtools.ipgp.fr/">HOME</a></td>
<td><a href="http://sourceforge.net/projects/shtools/">DOWNLOAD</a></td>
<td class="selected"><a href="../documentation.html">DOCUMENTATION</a></td>
<td><a href="https://sourceforge.net/p/shtools/discussion/">FORUM</a></td>
<td><a href="../faq.html">FAQ</a> </td>
</tr>
</tbody>
</table>
<p class="dir">
> <a href="../../SHTOOLS.html" class="dir">Home</a> > <a href="../documentation.html" class="dir">Documentation</a> > <a href="../gravmag.html" class="dir">Gravity and Magnetics</a></p>
<PRE>
<!-- Manpage converted by man2html 3.0.1 -->
<B>CILMPLUS(1)</B> SHTOOLS 2.9 <B>CILMPLUS(1)</B>
</PRE>
<H2 class="man">CilmPlus</H2 class="man"><PRE>
CilmPlus - Calculate the gravitational potential of relief referenced to a spherical
interface using the finite amplitude algorithm of Wieczorek and Phillips (1998).
</PRE>
<H2 class="man">SYNOPSIS</H2 class="man"><PRE>
SUBROUTINE CilmPlus ( CILM, GRIDIN, LMAX, NMAX, MASS, D, RHO, GRIDTYPE, W, ZERO, PLX, N, DREF
)
REAL*8 CILM(2, LMAX+1, LMAX+1), GRIDIN(LMAX+1, 2*LMAX+1), MASS, D, RHO,
REAL*8, OPTIONAL W(LMAX+1), ZERO(LMAX+1), PLX(LMAX+1, (LMAX+1)*(LMAX+2)/2), DREF
INTEGER LMAX, NMAX, GRIDTYPE
INTEGER, OPTIONAL N
</PRE>
<H2 class="man">DESCRIPTION</H2 class="man"><PRE>
<B>CilmPlus</B> will calculate the spherical harmonic coefficients of the gravitational potential
that correspond to constant density relief referenced to a spherical interface. This is
equation 10 of Wieczorek and Phillips (1998), where the potential is strictly valid only when
the coefficients are evaluated at a radius greater than the maximum radius of the relief. The
relief is input as a grid, whose type is specified by GRIDTYPE (1 for Gauss-Legendre
quadrature grids, 2 for NxN Driscoll and Healy sampled grids, and 3 for Nx2N Driscoll and
Healy sampled grids). The input relief GRIDIN must correspond to absolute radii. The
parameter NMAX is the order of the Taylor series used in the algorithm to approximate the
potential coefficients. By default, the relief and spherical harmonic coefficients will be
referenced to the mean radius of GRIDIN. However, if the optional parameter DREF is
specified, this will be used instead as the reference radius.
It is important to understand that as an intermediate step, this routine calculates the
spherical harmonic coefficients of the relief (referenced to the mean radius of GRIDIN or
DREF) raised to the nth power, i.e., (GRIDIN-D)**n. As such, if the input function is
bandlimited to degree L, the resulting function will be bandlimited to degree L*NMAX. This
subroutine implicitly assumes that the GRIDIN has an effective spherical harmonic bandwidth
greater or equal to L*NMAX. (The effective bandwidth is equal to LMAX for GRIDIN=1, and is
specified by N/2-1 for GRIDIN=2 or 3.) If this is not the case, then aliasing will occur. In
practice, for accurate results, it is found that the effective bandwidth needs only to be
about three times the size of L, though this should be verified for each application. Thus,
if the input function is considered to be bandlimited to degree L, the function should be
evaluated on a grid corresponding to a maximum degree of about 3*L. Aliasing effects can be
minimized by using Driscoll and Healy Nx2N grids.
If the input grid is evaluated on the Gauss-Legendre points, it is necessary to specify the
optional parameters W and ZERO, or W and PLX, which are calculated by a call to <B>PreCompute</B>.
In contast, if Driscoll-Healy grids are used (NxN or Nx2N), it is necessary to specify the
optional parameter N.
If memory is not an issue, the algorithm can be speeded up considerably when using Gauss-
Lengendre grids by inputing the optional array PLX (along with W) of precomputed associated
Legendre functions on the Gauss-Legendre nodes. Both of these variables are computed by a
call to <B>PreCompute</B>.
</PRE>
<H2 class="man">ARGUMENTS</H2 class="man"><PRE>
CILM (output) REAL*8, DIMENSION (2, LMAX+1, LMAX+1)
The real spherical harmonic coefficients (geodesy normalized) of the gravitational
potential corresponding to constant density relief referenced to a spherical
interface of radius D.
GRIDIN (input) REAL*8, DIMENSION (LMAX+1, 2*LMAX+1) for GRIDTYPE 1, DIMENSION (N, N) for
GRIDTYPE 2, DIMENSION (N, 2*N) for GRIDTYPE 3
The radii of the interface evaluated on a grid corresponding to a function of
maximum spherical harmonic degree LMAX. This is calculated by a call to either
<B>MakeGridGLQ</B> or <B>MakeGridDH</B>.
LMAX (input) INTEGER
The maximum spherical harmonic degree of the output spherical harmonic
coefficients. This degree also determines the dimension of the input relief GRIDGLQ
for GRIDTYPE 1. (As a general rule, this should be about twice the spherical
harmonic bandwidth of the input function.) For Driscoll-Healy grids, LMAX must be
less than or equal to N/2 - 1.
NMAX (input) INTEGER
The maximum order used in the Taylor-series expansion used in calculating the
potential coefficients. As a rule, this should be about 4.
MASS (input) REAL*8
The mass of the planet in kg.
D (output) REAL*8
The mean radius of the relief in meters.
RHO (input) REAL*8
The density contrast of the relief in kg/m^3.
GRIDTYPE (input) INTEGER
1 = Gauss-Legendre grids, calculated using <B>PreCompute</B> and <B>MakeGridGLQ</B>. 2 = Equally
sampled Driscoll-Healy grids, NxN, calculated using <B>MakeGridDH</B>. 3 = Equally spaced
Driscoll-Healy grids, Nx2N, calculated using <B>MakeGridDH</B>.
W (input) REAL*8, OPTIONAL DIMENSION (LMAX+1)
The weights used in the Gauss-Legendre quadrature, which are required for GRIDTYPE
= 1. These are calculated from a call to <B>PreCompute</B> (or alternatively, <B>PreGLQ</B>).
ZERO (input) REAL*8, OPTIONAL, DIMENSION (LMAX+1)
The nodes used in the Gauss-Legendre quadrature over latitude for GRIDTYPE 1,
calculated by a call to <B>PreCompute</B>. One of PLX or ZERO must be present when
GRIDTYPE = 1, but not both.
PLX (input) REAL*8, OPTIONAL, DIMENSION (LMAX+1, (LMAX+1)*(LMAX+2)/2)
An array of the associated Legendre functions calculated at the nodes used in the
Gauss-Legendre quadrature for GRIDTYPE 1. These are determined from a call to
<B>PreCompute</B>. One of PLX or ZERO must be present when GRIDTYPE = 1, but not both.
N (input) INTEGER, OPTIONAL
The number of samples in latitude when using Driscoll-Healy grids. For a function
bandlimited to LMAX, N = 2 (LMAX + 1). This is required for GRIDTYPE = 2 and 3.
DREF (input) REAL*8, OPTIONAL
The reference radius to be used when calculating both the relief and spherical
harmonic coefficients. If this is not specified, this parameter will be set equal
to the mean radius of GRIDIN.
</PRE>
<H2 class="man">NOTES</H2 class="man"><PRE>
This routine uses geodesy 4-pi normalized spherical harmonics that exclude the Condon-
Shortley phase; This can not be modified.
This routine requires the fast Fourier transform library <B>FFTW</B>, which is available at
<http://www.fftw.org>.
</PRE>
<H2 class="man">SEE ALSO</H2 class="man"><PRE>
<B>cilmplusrhoh(1)</B>, <B>shexpandglq(1)</B>, <B>makegridglq(1)</B>, <B>precompute(1)</B>, <B>preglq(1)</B>, <B>glqgridcoord(1)</B>,
<B>makegriddh(1)</B>
<http://shtools.ipgp.fr/>
</PRE>
<H2 class="man">REFERENCES</H2 class="man"><PRE>
Wieczorek, M. A. and R. J. Phillips, Potential anomalies on a sphere: applications to the
thickness of the lunar crust, <B>J.</B> <B>Geophys.</B> <B>Res.</B>, 103, 1715-1724, 1998.
</PRE>
<H2 class="man">COPYRIGHT AND LICENSE</H2 class="man"><PRE>
Copyright 2012 by Mark Wieczorek <wieczor@ipgp.fr>.
This is free software; you can distribute and modify it under the terms of the revised BSD
license.
SHTOOLS 2.9 2012-11-06 <B>CILMPLUS(1)</B>
</PRE>
<p class="dir">
> <a href="../../SHTOOLS.html" class="dir">Home</a> > <a href="../documentation.html" class="dir">Documentation</a> > <a href="../gravmag.html" class="dir">Gravity and Magnetics</a></p>
<table class="footer2" summary = "Mark Wieczorek">
<tbody>
<tr>
<td class="c1"><a href="http://www.ipgp.fr/">Institut de Physique du Globe de Paris</a></td>
<td class="c2"><a href="http://www.sorbonne-paris-cite.fr/index.php/en">University of Sorbonne Paris Cité</a></td>
<td class="c3">© 2014 <a href="http://www.ipgp.fr/~wieczor">Mark Wieczorek</a></td>
</tr>
</tbody>
</table>
</div>
</body>
</html>