# -----------------------------------------------------------------------------
# Bjorn-Gustaf J. Brooks, January 19, 2009
# -----------------------------------------------------------------------------
The "routines" directory contains the Octave/Matlab routine "pcm_locate.m" 
that is the basis for the paper: "Applying Wavelet and Fourier Transform 
Analysis to Large Geophysical Datasets". The code is simple enough so that 
it can easily be translated into another language. 

There are various other example IDL and Octave/Matlab scripts that I have 
adapted from other sources in order to provide practical examples of 
geophysical data analysis. There are however many other packages for 
spectral analysis. For example, GNU Octave provides an extra routine, 
specgram, (http://octave.sourceforge.net/doc/f/specgram.html) that can be 
used to produce nice spectrograms.

# -----------------------------------------------------------------------------
# LOCATE GRID CELLS THAT SHOW ANOMALOUS CHANGES IN PERIODICITY ----------------
"pcm_locate.m" is the basis of the main paper, which can be used to search 
through millions of data points of climate output. Care should be taken 
with regard to the frequency (f) over which you search each time series 
because braoder frequency ranges are much more computationally expensive.

# -----------------------------------------------------------------------------
# GENERATE FOURIER TRANSFORM POWER SPECTRUM -----------------------------------
Generate DFTPS the same as those shown in Figure 1 of the paper using 
the Octave/Matlab script contained in "routines/" directory "dft_h_all.m". 
This will plot the DFTPS of the first column from the data file 
"pcm_B05_12_hpt_0001_monthly.dat" and the DFTPW of the random walk and 
random block Monte Carlo simulations. Specify the variable "numtest" a 
suitable number of iterations for your MC trials. 10,000 trials will 
take about a half hour.

# -----------------------------------------------------------------------------
# GENERATE WAVELET SPECTROGRAM ------------------------------------------------
Generate a spectrogram the same as those shown in Figure 3 of the paper 
using the IDL scripts in the "routines/" directory "wv.pro" and 
"wavelet.pro", which are modified from Torrence and Compo's routines 
(available at: http://paos.colorado.edu/research/wavelets/). 
"wv.pro" will produce a spectrogram from data file 
"pcm_B05_12_hpt_0001_monthly.dat", by making a call to wavelet.pro, 
which does the convolution.

# -----------------------------------------------------------------------------
# COMPUTE MONTE CARLO SIGNIFICANCES -------------------------------------------
An example Octave/Matlab script, called "mc_sig_check.m", can be found in 
the "routines/" directory "mc_sig_check.m" is adapted from the MC methods 
published in Cornette's (2007 CISE) paper. "mc_sig_check.m" will calculate 
the significances of the soil moisture, precipitation, and 
temperature data from "pcm_B05_12_hpt_0851_monthly.dat" and send them to 
standard out. 10,000 trials should take nearly an hour for all 3. The 
routine will only search in the frequencies where the 12 month cycle 
should occur, and will count the number of times a peak generated by MC 
permutations exceeds that of the original time series. The significance 
is reported as the the fraction MC simulations that exceed the peak in 
the original series divided by the total number of trials. To 
change the frequency range that is tested, modify:
"[a b]=max(p(4501:6000));"

# -----------------------------------------------------------------------------

Enjoy!