-
Notifications
You must be signed in to change notification settings - Fork 103
/
shmultitapercse.html
115 lines (104 loc) · 8.68 KB
/
shmultitapercse.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
<!DOCTYPE HTML PUBLIC "-//W3C//DTD HTML 4.01//EN"
"http://www.w3.org/TR/html4/strict.dtd">
<html>
<head>
<title>SHTOOLS - Localized spectral analysis</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="../../localized.html" class="dir">Localized Spectral Analysis</a></p>
<h1 id="shmultitapercse">SHMultiTaperCSE</h1>
<p>Perform a localized multitaper cross-spectral analysis using spherical cap windows.</p>
<h1 id="usage">Usage</h1>
<p>call SHMultiTaperCSE (<code>mtse</code>, <code>sd</code>, <code>sh1</code>, <code>lmax1</code>, <code>sh2</code>, <code>lmax2</code>, <code>tapers</code>, <code>taper_order</code>, <code>lmaxt</code>, <code>k</code>, <code>alpha</code>, <code>lat</code>, <code>lon</code>, <code>taper_wt</code>, <code>norm</code>, <code>csphase</code>, <code>exitstatus</code>)</p>
<h1 id="parameters">Parameters</h1>
<dl>
<dt><code>mtse</code> : output, real*8, dimension (<code>lmax</code>-<code>lmaxt</code>+1)</dt>
<dd>The localized multitaper cross-power spectrum estimate. <code>lmax</code> is the smaller of <code>lmax1</code> and <code>lmax2</code>.
</dd>
<dt><code>sd</code> : output, real*8, dimension (<code>lmax</code>-<code>lmaxt</code>+1)</dt>
<dd>The standard error of the localized multitaper cross-power spectral estimates. <code>lmax</code> is the smaller of <code>lmax1</code> and <code>lmax2</code>.
</dd>
<dt><code>sh1</code> : input, real*8, dimension (2, <code>lmax1</code>+1, <code>lmax1</code>+1)</dt>
<dd>The spherical harmonic coefficients of the first function.
</dd>
<dt><code>lmax1</code> : input, integer</dt>
<dd>The spherical harmonic bandwidth of <code>sh1</code>.
</dd>
<dt><code>sh2</code> : input, real*8, dimension (2, <code>lmax2</code>+1, <code>lmax2</code>+1)</dt>
<dd>The spherical harmonic coefficients of the second function.
</dd>
<dt><code>lmax2</code> : input, integer</dt>
<dd>The spherical harmonic bandwidth of <code>sh2</code>.
</dd>
<dt><code>tapers</code> : input, real*8, dimension (<code>lmaxt</code>+1, <code>k</code>)</dt>
<dd>An array of the <code>k</code> windowing functions, arranged in columns, obtained from a call to <code>SHReturnTapers</code>. Each window has non-zero coefficients for a single angular order that is specified in the array <code>taper_order</code>.
</dd>
<dt><code>taper_order</code> : input, integer, dimension (<code>k</code>)</dt>
<dd>An array containing the angular orders of the spherical harmonic coefficients in each column of the array <code>tapers</code>.
</dd>
<dt><code>lmaxt</code> : input, integer</dt>
<dd>The spherical harmonic bandwidth of the windowing functions in the array <code>tapers</code>.
</dd>
<dt><code>k</code> : input, integer</dt>
<dd>The number of tapers to be utilized in performing the multitaper spectral analysis.
</dd>
<dt><code>alpha</code> : input, optional, real*8, dimension (3)</dt>
<dd>The Euler rotation angles used in rotating the windowing functions. <code>alpha(1)=0</code>, <code>alpha(2)=-(90-lat)*pi/180</code>, <code>alpha(3)=-lon*pi/180</code>. Either <code>alpha</code> or <code>lat</code> and <code>lon</code> can be specified, but not both. If none of these are specified, the spectral analysis will be centered at the north pole.
</dd>
<dt><code>lat</code> : input, optional, real*8</dt>
<dd>The latitude in degrees of the localized analysis. Either <code>alpha</code> or <code>lat</code> and <code>lon</code> can be specified but not both. If none of these are specified, the spectral analysis will be centered at the north pole.
</dd>
<dt><code>lon</code> : input, optional, real*8</dt>
<dd>The longitude in degrees of the localized analysis. Either <code>alpha</code> or <code>lat</code> and <code>lon</code> can be specified, but not both. If none of these are specified, the spectral analysis will be centered at the north pole.
</dd>
<dt><code>taper_wt</code> : input, optional, real*8, dimension (<code>k</code>)</dt>
<dd>The weights used in calculating the multitaper spectral estimates and standard error. Optimal values of the weights (for a known global power spectrum) can be obtained from the routine <code>SHMTVarOpt</code>.
</dd>
<dt><code>norm</code> : input, optional, integer, default = 1</dt>
<dd>1 (default) = 4-pi (geodesy) normalized harmonics; 2 = Schmidt semi-normalized harmonics; 3 = unnormalized harmonics; 4 = orthonormal harmonics.
</dd>
<dt><code>csphase</code> : input, optional, integer, default = 1</dt>
<dd>1 (default) = do not apply the Condon-Shortley phase factor to the associated Legendre functions; -1 = append the Condon-Shortley phase factor of (-1)^m to the associated Legendre functions.
</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>SHMultiTaperCSE</code> will perform a localized multitaper cross-spectral analysis of two input functions expressed in spherical harmonics, <code>SH1</code> and <code>SH2</code>. The maximum degree of the localized multitaper power spectrum estimate is <code>lmax-lmaxt</code>, where <code>lmax</code> is the smaller of <code>lmax1</code> and <code>lmax2</code>. The coefficients and angular orders of the windowing coefficients (<code>tapers</code> and <code>taper_order</code>) are obtained by a call to <code>SHReturnTapers</code>. If <code>lat</code> and <code>lon</code> or <code>alpha</code> is specified, then the symmetry axis of the localizing windows will be rotated to these coordinates. Otherwise, the localized spectral analysis will be centered over the north pole.</p>
<p>If the optional array <code>taper_wt</code> is specified, then these weights will be used in calculating a weighted average of the individual <code>k</code> tapered estimates (<code>mtse</code>) and the corresponding standard error of the estimates (<code>sd</code>). If not present, the weights will all be assumed to be equal. When <code>taper_wt</code> is not specified, the mutltitaper spectral estimate for a given degree will be calculated as the average obtained from the <code>k</code> individual tapered estimates. The standard error of the multitaper estimate at degree l is simply the population standard deviation, <code>S = sqrt(sum (Si - mtse)^2 / (k-1))</code>, divided by sqrt(<code>k</code>). See Wieczorek and Simons (2007) for the relevant expressions when weighted estimates are used.</p>
<p>The employed spherical harmonic normalization and Condon-Shortley phase convention can be set by the optional arguments <code>norm</code> and <code>csphase</code>; if not set, the default is to use geodesy 4-pi normalized harmonics that exclude the Condon-Shortley phase of (-1)^m.</p>
<h1 id="references">References</h1>
<p>Wieczorek, M. A. and F. J. Simons, Minimum-variance multitaper spectral estimation on the sphere, J. Fourier Anal. Appl., 13, doi:10.1007/s00041-006-6904-1, 665-692, 2007.</p>
<h1 id="see-also">See also</h1>
<p><a href="shmultitaperse.html">shmultitaperse</a>, <a href="shreturntapers.html">shreturntapers</a>, <a href="shreturntapersm.html">shreturntapersm</a>, <a href="shmtvaropt.html">shmtvaropt</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="../../localized.html" class="dir">Localized Spectral Analysis</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>