We have seen a number of examples of Lisp-Stat's built-in functions, now let us see how to define our own. A common task in geophysics is to compute the power spectrum or auto-correlation of a time series. Let's see how we might do this with Lisp-Stat. First let's make up some data and look at the real part of its spectrum.
>(def x (normal-rand 256)) X (plot-lines (iseq (length x)) (realpart (fft x)))
We won't bother showing this plot here, but you can convince yourself quickly that the spectrum is flat except for the oscillations about the (zero) mean.
Now, (help 'fft)
tells us that Lisp-Stat's fast
Fourier transform routine returns the un-normalized
Fourier transform of its first argument, and the inverse
fft if its optional second argument is true.
Thus the power spectrum and the autocorrelation are given by
> (def power (realpart (* (fft x) (conjugate (fft x))))) POWER > (def autocor (realpart (fft power t))) AUTOCOR > ; which we can normalize using the built-in max function (def normalized (/ autocor (max autocor))) NORMALIZED
Now this is fine as far as it goes, but if we're going to be computing the power and autocorrelation all the time, we might as well define our own lisp functions to do this. Let's start by making a function to compute the power spectrum.
> (defun power (z)
(realpart (* (fft z) (conjugate (fft z))))
)
POWER
To invoke this we just say (power x)
where x is an array.
The syntax for defun
is
(defun <name>
<parameters> <body>)
where name is the symbol used as the function name, parameters is a list of symbols making up the arguments, and body is one or more expressions making up the body of the function.
>(defun autocor (z) (realpart (fft (power z) t)) ) AUTOCOR
With this function we could define, calculate and plot the autocorrelation of a random time series by:
> (def x (normal-rand 256)) X > (plot-lines (iseq (length x)) (autocor x))
But we can go one step further in making our autocorrelation function useful, by defining an optional argument that lets us decide whether we want the autocorrelation normalized or not.
> (defun autocor (z &key normalize)
(let ((auto (realpart (fft (power z) t))))
(if normalize (/ auto (max auto)) auto)
)
)
AUTOCOR
The function let (along with let*) is used
in Lisp to define local variables. Thus auto
is used defined only inside this function. The
&key means that the string ``normalize'' defines
an optional keyword; see [Tie90] for more
details.
Now we can compute the un-normalized or normalized autocorrelation with the same Lisp function, to wit:
> (def ax (autocor x)) AX > (def nax (autocor x :normalize t)) NAX > (plot-lines (iseq (length ax)) ax) > (plot-lines (iseq (length nax)) nax)
with the results as shown in Figure 1.6.
Figure 1.6: Ordinary and normalized autocorrelation of an
uncorrelated Gaussian process.