Commit 976eac7b by Sumant Oemrawsingh

Added some extra functions for use with FFTs

```fft-frequency-vector can be used to make a vector that contains the actual
sample frequencies (so you don't have to guess or write your own).

fft-shift and fft-inverse-shift shift the zero frequency to the center (for
easy plotting), and the inverse (for proper inverse FFT).

These functions are comparable to fftfreq, fftshift and ifftshift in numpy (if
anyone cares...).```
parent 7a78c863
 ... ... @@ -407,6 +407,19 @@ GSL doesn't have them.
• invert-matrix finds the inverse of a matrix and uses GSL's LU decomposition functions.
• IEEE floating point number analysis.
• fft-frequency-vector returns a vector where the sample frequencies are contained. If you perform an FFT on a vector of a given size and :sample-size, this vector will contain the sample frequencies in order. If the :shifted keyword is T, then the frequencies are ordered in ascending order.
• fft-shift returns a copy of a vector where the zero frequency has been shifted to the center; the frequency components will be sorted according to their frequency, in ascending order. Optionally, a :stride can be provided.
• fft-inverse-shift performs the inverse action of fft-shift; the zero and positive frequency components are shifted to the beginning, so that the resulting vector is suitable for an inverse FFT. Optionally, a :stride can be provided.
• Status

... ...
 ;; Extras for FFT. ;; Sumant Oemrawsingh, Wed Nov 25 2009 - 22:03 (in-package :gsll) (defun fft-frequency-step (size &optional (sample-spacing 1)) "Return the frequency step for the FFT of a vector with the given size and assuming the given sample spacing." (/ (* size sample-spacing))) (defun fft-highest-frequency (size &optional (sample-spacing 1)) "Return the highest frequency for the FFT of a vector with the given size and assuming the given sample spacing." (if (evenp size) (values (/ sample-spacing 2) (/ size 2)) (- (/ sample-spacing 2) (/ (fft-frequency-step size sample-spacing) 2)))) (defun fft-frequency-split (size) "Return the index where the highest frequency component is located in a vector of the given size, after applying an FFT. The second value that is returned, is the most negative frequency divided by the sample frequency step." (if (evenp size) (let ((n (/ size 2))) (values n (- 1 n))) (let ((n (/ (- size 1) 2))) (values n (- n))))) (export 'fft-frequency-vector) (defun fft-frequency-vector (element-type size &key (sample-spacing 1) shifted) "Make and return a vector that contains the sample frequencies of an FFT that has been applied to a vector with the given size and :sample-spacing. If :shifted is T, then the vector will contain the sample frequencies after FFT-SHIFT has been applied to the result of an FFT." (let* ((df (fft-frequency-step size sample-spacing)) pos-f neg-f) (multiple-value-bind (max-f-index min-f-index) (fft-frequency-split size) (setf pos-f (loop for i from 0 to max-f-index collect (coerce (* i df) element-type)) neg-f (loop for i from min-f-index below 0 collect (coerce (* i df) element-type))) (make-marray element-type :initial-contents (if shifted (nconc neg-f pos-f) (nconc pos-f neg-f)))))) (export 'fft-shift) (defun fft-shift (vector &key (stride 1)) "Return a copy of a vector that is the result of an FFT, with the zero and positive frequencies shifted to the center and end, so that the data is suitable for e.g. plotting." (let* ((size (size vector)) (n (floor size stride)) (shifted (copy vector)) (split (fft-frequency-split n)) (zero-pos (if (evenp n) (- split 1) split))) ;; Positive frequencies (loop for i from zero-pos to (- n 1) do (setf (maref shifted (* i stride)) (maref vector (* (- i zero-pos) stride)))) ;; Negative frequencies (loop for i from (+ split 1) to (- n 1) do (setf (maref shifted (* (- i split 1) stride)) (maref vector (* i stride)))) shifted)) (export 'fft-inverse-shift) (defun fft-inverse-shift (vector &key (stride 1)) "Return a copy of a vector where the zero and positive frequency components are shifted to the beginning, i.e. ordered so that it is suitable for an inverse FFT." (let* ((size (size vector)) (n (floor size stride)) (inv-shifted (copy vector)) (split (fft-frequency-split n)) (zero-pos (if (evenp n) (- split 1) split))) ;; Positive frequencies (loop for i from zero-pos to (- n 1) do (setf (maref inv-shifted (* (- i zero-pos) stride)) (maref vector (* i stride)))) ;; Negative frequencies (loop for i from (+ split 1) to (- n 1) do (setf (maref inv-shifted (* i stride)) (maref vector (* (- i split 1) stride)))) inv-shifted))
 ... ... @@ -143,6 +143,7 @@ (:file "select-direction") (:file "unpack") (:file "discrete") (:file "extras") (:file "example"))) (:module random :depends-on (init data) ... ...
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!