/
pycilmplusdh.html
97 lines (86 loc) · 6.68 KB
/
pycilmplusdh.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
<!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="../../images/favicon.ico">
</head>
<body>
<div class="main">
<p class="centeredimage"><img src="../../images/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="https://github.com/SHTOOLS/SHTOOLS/releases">DOWNLOAD</a></td>
<td class="selected"><a href="../../documentation.html">DOCUMENTATION</a></td>
<td><a href="../../faq.html">FAQ</a> </td>
</tr>
</tbody>
</table>
<p class="dir">
> <a href="../../../index.html" class="dir">Home</a> > <a href="../../documentation.html" class="dir">Documentation</a> > <a href="../../python-routines.html" class="dir">Python</a> > <a href="../../pygravmag.html" class="dir">Gravity and Magnetics</a></p>
<h1 id="cilmplusdh">CilmPlusDH</h1>
<p>Calculate the gravitational potential exterior to relief referenced to a spherical interface using the finite-amplitude algorithm of Wieczorek and Phillips (1998).</p>
<h1 id="usage">Usage</h1>
<p><code>cilm</code>, <code>d</code> = pyshtools.CilmPlusDH (<code>gridin</code>, <code>nmax</code>, <code>mass</code>, <code>rho</code>, [<code>lmax</code>, <code>n</code>, <code>dref</code>, <code>sampling</code>])</p>
<h1 id="returns">Returns</h1>
<dl>
<dt><code>cilm</code> : float, dimension (2, <code>lmax</code>+1, <code>lmax</code>+1)</dt>
<dd>The real spherical harmonic coefficients (geodesy normalized) of the gravitational potential corresponding to constant density relief referenced to a spherical interface of radius <code>d</code>.
</dd>
<dt><code>d</code> : float</dt>
<dd>The mean radius of the relief in meters.
</dd>
</dl>
<h1 id="parameters">Parameters</h1>
<dl>
<dt><code>gridin</code> : float, dimension (<code>nin</code>, <code>sampling</code>*<code>nin</code>)</dt>
<dd>The radii of the interface evaluated on a grid, determined by a call to <code>MakeGridDH</code>.
</dd>
<dt><code>nmax</code> : integer</dt>
<dd>The maximum order used in the Taylor-series expansion used in calculating the potential coefficients.
</dd>
<dt><code>mass</code> : float</dt>
<dd>The mass of the planet in kg.
</dd>
<dt><code>rho</code> : float</dt>
<dd>The density contrast of the relief in kg/m^3.
</dd>
<dt><code>lmax</code> : optional, integer, default = <code>n/2-1</code></dt>
<dd>The maximum spherical harmonic degree of the output spherical harmonic coefficients. <code>lmax</code> must be less than or equal to <code>n/2-1</code>.
</dd>
<dt><code>n</code> : optional, integer, default = <code>nin</code></dt>
<dd>The number of samples in latitude when using Driscoll-Healy grids. For a function bandlimited to <code>lmax</code>, <code>n=2(lmax+1)</code>.
</dd>
<dt><code>dref</code> : optional, float</dt>
<dd>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 <code>gridin</code>.
</dd>
<dt><code>sampling</code> : optional, integer, default determined by dimensions of <code>gridin</code></dt>
<dd>If 1 the output grids are equally sampled (<code>n</code> by <code>n</code>). If 2, the grids are equally spaced (<code>n</code> by 2<code>n</code>).
</dd>
</dl>
<h1 id="description">Description</h1>
<p><code>CilmPlus</code> will calculate the spherical harmonic coefficients of the gravitational potential exterior 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 input relief <code>gridin</code> must correspond to absolute radii. The parameter <code>nmax</code> 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 <code>gridin</code>. However, if the optional parameter <code>dref</code> is specified, this will be used instead as the reference radius.</p>
<p>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 <code>gridin</code> or <code>dref</code>) raised to the nth power, i.e., <code>(gridin-d)\*\*n</code>. As such, if the input function is bandlimited to degree <code>L</code>, the resulting function will be bandlimited to degree <code>L*nmax</code>. This subroutine implicitly assumes that the <code>gridin</code> has an effective spherical harmonic bandwidth greater or equal to this value. (The effective bandwidth is equal to <code>n/2-1</code>) If this is not the case, 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 <code>L</code>, though this should be verified for each application. Thus, if the input function is considered to be bandlimited to degree <code>L</code>, the function should be evaluated on a grid corresponding to a maximum degree of about <code>3*</code>L.</p>
<p>This routine uses geodesy 4-pi normalized spherical harmonics that exclude the Condon-Shortley phase.</p>
<h1 id="references">References</h1>
<p>Wieczorek, M. A. and R. J. Phillips, Potential anomalies on a sphere: applications to the thickness of the lunar crust, J. Geophys. Res., 103, 1715-1724, 1998.</p>
<h1 id="see-also">See also</h1>
<p><a href="pycilmplusrhohdh.html">cilmplusrhohdh</a>, <a href="pycilmminusdh.html">cilmminusdh</a>, <a href="pycilmminusrhohdh.html">cilmminusrhohdh</a>, <a href="pymakegriddh.html">makegriddh</a></p>
<p class="dir">
> <a href="../../../index.html" class="dir">Home</a> > <a href="../../documentation.html" class="dir">Documentation</a> > <a href="../../python-routines.html" class="dir">Python</a> > <a href="../../pygravmag.html" class="dir">Gravity and Magnetics</a></p>
<table class="footer2" summary = "SHTOOLS; Fortran and Python spherical harmonic transform software package">
<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">© 2016 <a href="https://github.com/SHTOOLS">SHTOOLS</a></td>
</tr>
</tbody>
</table>
</div>
</body>
</html>