Commit 7192e225 authored by Liam Healy's avatar Liam Healy
Browse files

FFT Unpack

Function #'unpack will switch between the various unpacking generic
functions, depending on the specified unpack type and whether the
vector length is a power of 2.  New file unpack.lisp has all the
unpacking functions.
parent 39c888a4
;; Functions for fast fourier transforms on real data.
;; Sumant Oemrawsingh, Sat Oct 31 2009 - 20:12
;; Time-stamp: <2009-11-01 23:04:46EST fft-half-complex.lisp>
;; Time-stamp: <2009-11-03 23:02:17EST fft-half-complex.lisp>
;; /usr/include/gsl/gsl_fft_halfcomplex.h
;; /usr/include/gsl/gsl_fft_halfcomplex_float.h
......@@ -33,22 +33,6 @@
:documentation
"Inverse FFT for a half-complex radix-2 vector")
(defmfun fft-half-complex-radix2-unpack
((h-complex vector)
&key (stride 1) (n (expt 2 (floor (log (size h-complex) 2))))
(complex-output (eltcase single-float (make-marray '(complex single-float) :dimensions n)
t (make-marray '(complex double-float) :dimensions n))))
("gsl_fft_halfcomplex" :type "_radix2_unpack")
(((c-pointer h-complex) :pointer) ((c-pointer complex-output) :pointer) (stride sizet) (n sizet))
:definition :generic
:element-types :float
:inputs (h-complex complex-output)
:outputs (complex-output)
:return (complex-output)
:documentation
"Convert an array of half-complex coefficients as returned by
real-fft-radix2-transform, into an ordinary complex array.")
;; Mixed Radix general-N functions
(defmobject fft-half-complex-wavetable-double-float
......@@ -131,20 +115,3 @@
:documentation
"Forward FFT for a half-complex vector")
(defmfun fft-half-complex-unpack
((h-complex vector)
&key (stride 1) (n (size h-complex))
(complex-output (eltcase single-float (make-marray '(complex single-float) :dimensions n)
t (make-marray '(complex double-float) :dimensions n))))
("gsl_fft_halfcomplex" :type "_unpack")
(((c-pointer h-complex) :pointer) ((c-pointer complex-output) :pointer) (stride sizet) (n sizet))
:definition :generic
:element-types :float
:inputs (h-complex complex-output)
:outputs (complex-output)
:return (complex-output)
:documentation
"This function converts an array of half-complex
coefficients as returned by fft-real-transform, into an ordinary complex
array. It fills in the complex array using the symmetry z_k = z_{n-k}^* to
reconstruct the redundant elements.")
;; Functions for fast fourier transforms on real data.
;; Sumant Oemrawsingh, Sun Oct 25 2009 - 16:35
;; Time-stamp: <2009-11-01 22:41:14EST fft-real.lisp>
;; Time-stamp: <2009-11-03 23:02:17EST fft-real.lisp>
;; /usr/include/gsl/gsl_fft_real.h
;; /usr/include/gsl/gsl_fft_real_float.h
......@@ -65,21 +65,3 @@
:documentation
"This function allocates a workspace for a real float transform of length
n.")
(defmfun fft-real-unpack
((vector vector)
&key (stride 1) (n (size vector))
(complex-output
(eltcase single-float (make-marray '(complex single-float) :dimensions n)
t (make-marray '(complex double-float) :dimensions n))))
("gsl_fft_real" :type "_unpack")
(((c-pointer vector) :pointer) ((c-pointer complex-output) :pointer) (stride sizet) (n sizet))
:definition :generic
:element-types :float
:inputs (vector complex-output)
:outputs (complex-output)
:return (complex-output)
:documentation
"This function converts a single real array into an equivalent complex
array (with imaginary part set to zero), suitable for fft-complex
routines.")
;; Unpack functions for FFT vectors.
;; Sumant Oemrawsingh, Sun Oct 25 2009 - 16:35
;; Time-stamp: <2009-11-03 23:20:35EST unpack.lisp>
(in-package :gsl)
;; /usr/include/gsl/gsl_fft_real.h
;; /usr/include/gsl/gsl_fft_real_float.h
;; /usr/include/gsl/gsl_fft_complex.h
;; /usr/include/gsl/gsl_fft_complex_float.h
(defmfun fft-real-unpack
((vector vector)
&key (stride 1)
(output
(eltcase single-float
(make-marray '(complex single-float) :dimensions (size vector))
t
(make-marray '(complex double-float) :dimensions (size vector)))))
("gsl_fft_real" :type "_unpack")
(((c-pointer vector) :pointer)
((c-pointer output) :pointer)
(stride sizet) ((size vector) sizet))
:definition :generic
:element-types :float
:inputs (vector output)
:outputs (output)
:return (output)
:export nil
:index unpack
:documentation
"This function converts a single real array into an equivalent complex
array (with imaginary part set to zero), suitable for fft-complex
routines.")
(defmfun fft-half-complex-radix2-unpack
((vector vector)
&key (stride 1)
(output
(eltcase single-float
(make-marray '(complex single-float) :dimensions (size vector))
t
(make-marray '(complex double-float) :dimensions (size vector)))))
("gsl_fft_halfcomplex" :type "_radix2_unpack")
(((c-pointer vector) :pointer)
((c-pointer output) :pointer)
(stride sizet) ((size vector) sizet))
:definition :generic
:element-types :float
:inputs (vector output)
:outputs (output)
:return (output)
:export nil
:index unpack
:documentation
"Convert an array of half-complex coefficients as returned by
real-fft-radix2-transform, into an ordinary complex array.")
(defmfun fft-half-complex-unpack
((vector vector)
&key (stride 1)
(output
(eltcase single-float
(make-marray '(complex single-float) :dimensions (size vector))
t
(make-marray '(complex double-float) :dimensions (size vector)))))
("gsl_fft_halfcomplex" :type "_unpack")
(((c-pointer vector) :pointer)
((c-pointer output) :pointer)
(stride sizet) ((size vector) sizet))
:definition :generic
:element-types :float
:inputs (vector output)
:outputs (output)
:return (output)
:export nil
:index unpack
:documentation
"This function converts an array of half-complex coefficients as
returned by fft-real-transform, into an ordinary complex array. It
fills in the complex array using the symmetry z_k = z_{n-k}^* to
reconstruct the redundant elements.")
(export 'unpack)
(defun unpack
(vector &rest args &key (stride 1) (unpack-type 'real) &allow-other-keys)
(let ((pass-on-args (copy-list args)))
(remf pass-on-args :unpack-type)
(if (eq unpack-type 'complex)
(if (power-of-2-p (floor (size vector) stride))
(apply 'fft-half-complex-radix2-unpack pass-on-args)
(apply 'fft-half-complex-unpack pass-on-args))
(apply 'fft-real-unpack pass-on-args))))
;; Definition of GSLL system
;; Liam Healy
;; Time-stamp: <2009-11-01 23:01:16EST gsll.asd>
;; Time-stamp: <2009-11-03 23:23:15EST gsll.asd>
(when (asdf:find-system :fsbv nil)
(pushnew :fsbv *features*))
......@@ -140,6 +140,7 @@
(:file "fft-real")
(:file "fft-half-complex" :depends-on ("fft-real"))
(:file "forward")
(:file "unpack")
(:file "fft-interface" :depends-on ("fft-complex" "fft-half-complex" "fft-real"))
(:file "fft-interface-example")))
(:module random
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment