-
Notifications
You must be signed in to change notification settings - Fork 103
/
makegeoidgrid.html
116 lines (105 loc) · 8.26 KB
/
makegeoidgrid.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
<!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="../../../index.html">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="../../f95-routines.html" class="dir">Fortran 95</a> > <a href="../../gravmag.html" class="dir">Gravity and Magnetics</a></p>
<h1 id="makegeoidgrid">MakeGeoidGrid</h1>
<p>Create a global map of the geoid.</p>
<h1 id="usage">Usage</h1>
<p>call MakeGeoidGrid (<code>geoid</code>, <code>cilm</code>, <code>lmax</code>, <code>r0</code>, <code>gm</code>, <code>potref</code>, <code>omega</code>, <code>r</code>, <code>gridtype</code>, <code>order</code>, <code>nlat</code>, <code>nlong</code>, <code>interval</code>, <code>lmax_calc</code>, <code>a</code>, <code>f</code>, <code>exitstatus</code>)</p>
<h1 id="parameters">Parameters</h1>
<dl>
<dt><code>geoid</code> : output, real*8, dimension(<code>nlat</code>, <code>nlong</code>)</dt>
<dd>A global grid of the height to the potential <code>potref</code> above a sphere of radius <code>r</code> (or above a flattened ellipsoid if both <code>a</code> and <code>f</code> are specified). The number of latitude and longitude points depends upon <code>gridtype</code>: (1) <code>lmax+1</code> by <code>2lmax+1</code>, (2) <code>2lmax+2</code> by <code>2lmax+2</code>, (3) <code>2lmax+2</code> by <code>4lmax+4</code>, or (4) <code>180/interval+1</code> by <code>360/interval+1</code>.
</dd>
<dt><code>cilm</code> : input, real*8, dimension (2, <code>lmax</code>+1, <code>lmax</code>+1)</dt>
<dd>The real spherical harmonic coefficients (geodesy normalized) of the gravitational potential referenced to a spherical interface of radius <code>r0pot</code>.
</dd>
<dt><code>lmax</code> : input, integer</dt>
<dd>The maximum spherical harmonic degree of the gravitational-potential coefficients. For <code>gridtype</code>s 1, 2 and 3, this determines the number of latitudinal and longitudinal samples.
</dd>
<dt><code>r0</code> : input, real*8</dt>
<dd>The reference radius of the spherical harmonic coefficients.
</dd>
<dt><code>gm</code> : input, real*8</dt>
<dd>The product of the gravitational constant and mass of the planet.
</dd>
<dt><code>potref</code> : input, real*8</dt>
<dd>The value of the potential on the chosen geoid, in SI units.
</dd>
<dt><code>omega</code> : input, real*8</dt>
<dd>The angular rotation rate of the planet.
</dd>
<dt><code>r</code> : input, real*8</dt>
<dd>The radius of the reference sphere that the Taylor expansion of the potential is performed on. If <code>a</code> and <code>f</code> are not specified, the geoid height will be referenced to this spherical interface.
</dd>
<dt><code>gridtype</code> : input, integer</dt>
<dd>The output grid is (1) a Gauss-Legendre quadrature grid whose grid nodes are determined by <code>lmax</code>, (2) an equally sampled <code>n</code> by <code>n</code> grid used with the Driscoll and Healy (1994) sampling theorem, (3) ar a similar <code>n</code> by 2<code>n</code> grid that is oversampled in longitude, or (4) a 2D Cartesian grid with latitudinal and longitudinal spacing given by <code>interval</code>.
</dd>
<dt><code>order</code> : input, integer</dt>
<dd>The order of the Taylor series expansion of the potential about the reference radius <code>r</code>. This can be either 1, 2, or 3.
</dd>
<dt><code>nlat</code> : output, integer</dt>
<dd>The number of latitudinal samples.
</dd>
<dt><code>nlong</code> : output, integer</dt>
<dd>The number of longitudinal samples.
</dd>
<dt><code>interval</code>: optional, input, real*8</dt>
<dd>The latitudinal and longitudinal spacing of the output grid when <code>gridtype</code> is 4.
</dd>
<dt><code>lmax_calc</code> : optional, input, integer</dt>
<dd>The maximum degree used in evaluating the spherical harmonic coefficients. This must be less than or equal to <code>lmax</code>.
</dd>
<dt><code>a</code> : optional, input, real*8, default = <code>r0</code></dt>
<dd>The semi-major axis of the flattened ellipsoid that the output grid <code>geoid</code> is referenced to. The optional parameter <code>f</code> must also be specified.
</dd>
<dt><code>f</code> : optional, input, real*8, default = 0</dt>
<dd>The flattening <code>(R_equator-R_pole)/R_equator</code> of the reference ellipsoid. The optional parameter <code>a</code> (i.e., <code>R_equator</code>) must be specified.
</dd>
<dt><code>exitstatus</code> : output, optional, integer</dt>
<dd>If present, instead of executing a STOP when an error is encountered, the variable exitstatus will be returned describing the error. 0 = No errors; 1 = Improper dimensions of input array; 2 = Improper bounds for input variable; 3 = Error allocating memory; 4 = File IO error.
</dd>
</dl>
<h1 id="description">Description</h1>
<p><code>MakeGeoidGrid</code> will create a global map of the geoid, accurate to either first, second, or third order, using the method described in Wieczorek (2007; equation 19-20). The algorithm expands the potential in a Taylor series on a spherical interface of radius <code>r</code>, and computes the height above this interface to the potential <code>potref</code> exactly from the linear, quadratic, or cubic equation at each grid point. If the optional parameters <code>a</code> and <code>f</code> are specified, the geoid height will be referenced to a flattened ellipsoid with semi-major axis <code>a</code> and flattening <code>f</code>. The pseudo-rotational potential is explicitly accounted for by specifying the angular rotation rate <code>omega</code> of the planet.</p>
<p>It should be noted that this geoid calculation is only strictly exact when the radius <code>r</code> lies above the maximum radius of the planet. Furthermore, the geoid is only strictly valid when it lies above the surface of the planet (it is necessary to know the density structure of the planet when calculating the potential below the surface).</p>
<p>The geoid can be computed on one of four different grids: (1) a Gauss-Legendre quadrature grid (see <code>MakeGridGLQ</code>), (2) A <code>n</code> by <code>n</code> equally sampled grid (see <code>MakeGridDH</code>), (3) an <code>n</code> by 2<code>n</code> equally spaced grid (see <code>MakeGridDH</code>), or (4) A 2D Cartesian grid (see <code>MakeGrid2D</code>). This routine uses geodesy 4-pi normalized spherical harmonics that exclude the Condon-Shortley phase. This can not be modified.</p>
<h1 id="references">References</h1>
<p>Driscoll, J.R. and D.M. Healy, Computing Fourier transforms and convolutions on the 2-sphere, Adv. Appl. Math., 15, 202-250, 1994.</p>
<p>Wieczorek, M. A. Gravity and topography of the terrestrial planets, Treatise on Geophysics, 10, 165-206, 2007.</p>
<h1 id="see-also">See also</h1>
<p><a href="makegrid2d.html">makegrid2d</a>, <a href="makegridglq.html">makegridglq</a>, <a href="makegriddh.html">makegriddh</a>, <a href="makegravgriddh.html">makegravgriddh</a>, <a href="makegravgradgriddh.html">makegravgradgriddh</a></p>
<p class="dir">
> <a href="../../../index.html" class="dir">Home</a> > <a href="../../documentation.html" class="dir">Documentation</a> > <a href="../../f95-routines.html" class="dir">Fortran 95</a> > <a href="../../gravmag.html" class="dir">Gravity and Magnetics</a></p>
<table class="footer2">
<tbody>
<tr>
<td class="c1"><a href="https://lagrange.oca.eu/?lang=en">Laboratoire Lagrange</a></td>
<td class="c2"><a href="https://www.oca.eu/?lang=en">Observatoire de la Côte d'Azur</a></td>
<td class="c3">© 2017 <a href="https://github.com/SHTOOLS">SHTOOLS</a></td>
</tr>
</tbody>
</table>
</div>
</body>
</html>