-
Notifications
You must be signed in to change notification settings - Fork 103
/
shbiasadmitcorr.doc
88 lines (79 loc) · 4.2 KB
/
shbiasadmitcorr.doc
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
Calculate the expected multitaper admittance and correlation spectra associated
with the input global cross-power spectra of two functions.
Usage
-----
admit, corr = SHAdmitCorr (sgt, sgg, stt, tapers, [lmax, lwin, k, mtdef,
taper_wt])
Returns
-------
admit : float, dimension (lmax-lwin+1)
The biased admittance spectrum obtained using the localized (cross-)power
spectra of Sgt and Stt.
corr : float, dimension (lmax-lwin+1)
The biased correlation spectrum obtained using the localized (cross-)power
spectra of Sgt, Stt, and Sgg.
Parameters
----------
sgt : float, dimension (lmaxin1+1)
The global cross-power spectrum of the functions G and T.
sgg : float, dimension (lmaxin2+1)
The global power spectrum of the function G.
stt : float, dimension (lmaxin3+1)
The global power spectrum of the function T.
tapers : float, dimension (lwinin+1, kin)
The spherical harmonic coefficients of the spherical cap localizing windows.
Each column corresponds to the non-zero coefficients of a single angular
order. Since all that is important is the power spectrum of each window, the
exact angular order is not important. These are generated by a call to
SHReturnTapers or SHReturnTapersM.
lmax : optional, integer, default = min(lmaxin1, lmaxin2, lmaxin3)
The maximum spherical harmonic degree to use of the input spectra sgt, sgg,
stt.
lwin : optional, integer, default = lwinin
The spherical harmonic bandwidth of the localizing windows.
k : optional, integer, default = kin
The number of localizing windows to use. Only the first k columns of tapers
will be employed, which corresponds to the k best-concentrated localizing
windows.
mtdef : optional, integer, default = 1
1 (default): Calculate the multitaper spectral estimates Sgt, Sgg and Stt
first, and then use these to calculate the admittance and correlation
functions. 2: Calculate admittance and correlation spectra using each
individual taper, and then average these to obtain the multitaper admittance
and correlation functions.
taper_wt : optional, float, dimension (k), default = -1
The weights to apply to each individual windowed spectral estimate. The
weights must sum to unity and are obtained from SHMTVarOpt. The default
specifies that taper weights are not used.
Description
-----------
Given the global cross-power spectra Sgt, Sgg and Stt of two functions G and T,
SHBiasAdmitCorr will calculate the expected multitaper admittance and
correlation spectra associated with the two global functions. This routine
expects as input a matrix containing the spherical harmonic coefficients of the
localizing windows, which can be generated by a call to SHReturnTapers or
SHReturnTapersM. Only the k best-concentrated localization windows will be
employed when calculating the biased cross-power spectra. The maximum calculated
degree of the output biased admittance and correlation spectra corresponds to
lmax-lwin, as it is assumed that the input cross-power spectra beyond lmax are
unknown, and not zero.
Two manners of calculating the localized admittance and correlation spectra are
possible according to the value of the optional parameter mtdef. In case 1, the
multitaper cross-power spectra of Sgt, Sgg, and Stt are first calculated, and
from these, the admittance and correlation spectra are formed. In case 2, the
biased admittance and correlation spectra are calculated for each individual
taper, and these are then averaged to obtain the biased multitaper admittance
and correlation spectra.
The default is to apply equal weights to each individual windowed estimate of
the spectrum, but this can be modified when mtdef is 1 by specifying the weights
in the optional argument taper_wt. The weights must sum to unity and can be
calculated by SHMTVarOpt.
References
----------
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.
Simons, F. J., F. A. Dahlen and M. A. Wieczorek, Spatiospectral concentration on
the sphere, SIAM Review, 48, 504-536, doi:10.1137/S0036144504445765, 2006.
Wieczorek, M. A. and F. J. Simons, Localized spectral analysis on the sphere,
Geophys. J. Int., 162, 655-675.