Geology Reference
In-Depth Information
C Initialize DFT
. DFT(I)=(1.D0,0.D0)
DO 16 J=1,K
AJ=DFLOAT(J)
ARG=-OMEGA*FREQ(I)
ARG=ARG*AJ*DT
DFT(I)=DFT(I)+G(J)*DCMPLX(DCOS(ARG),DSIN(ARG))
16 CONTINUE
15 CONTINUE
C Find and write out MEM spectrum.
DO 17 I=1,301
SME(I)=DT*PEP/(DFT(I)*DCONJG(DFT(I)))
WRITE(2,18)FREQ(I),SME(I)
18
FORMAT(1X,D15.5,D26.15)
17
CONTINUE
END
The spectral density estimated by the maximum entropy method (MEM), for
a prediction error filter of length 500, is shown in Figure 4.7. Both the Chandler
resonance and the annual motion are well resolved and there is the suggestion
of a weaker retrograde annual motion. Compared with the spectra produced by
singular value decomposition, shown in Figure 4.5, and that found for the pole path
interpolated by natural splines found by conventional discrete Fourier transform
methods, shown in Figure 4.6, the maximum entropy spectrum does indeed show
much improved resolution.
A more detailed plot of the Chandler and annual features of the maximum en-
tropy spectral density is shown in Figure 4.8. The annual motion peaks at 0.992 cpy
with corresponding period about 3 days longer than a year. Interestingly, the much
weaker retrograde annual motion shown in Figure 4.7 peaks at a period a few days
shorter than a year.
10 5
10 4
10 3
10 2
10
1
10 -1
10 -2
-1.5
- 1.0
- 0.5 - 0.0 0.5 1.0
Frequency in cycles per year (cpy)
1.5
Figure 4.7 The maximum entropy spectral density found by a prediction error
filter of length 500.
Search WWH ::




Custom Search