diff --git a/fast-fourier-transforms/backward.lisp b/fast-fourier-transforms/backward.lisp index e9866b4ff15e6ca96271d569ead820b31ff38b80..8b62af1d75d85d97a9707c518a8c994f67fa8166 100644 --- a/fast-fourier-transforms/backward.lisp +++ b/fast-fourier-transforms/backward.lisp @@ -126,16 +126,18 @@ (export 'backward-fourier-transform) (defun backward-fourier-transform (vector &rest args - &key decimation-in-frequency (stride 1) &allow-other-keys) + &key decimation-in-frequency (stride 1) non-radix-2 &allow-other-keys) "Perform a backward fast Fourier transform on the given vector. If the length of the vector is not a power of 2, and the user has a suitable wavetable and/or workspace, these can be supplied as - keyword arguments. If the vector is real, it is assumed to be - in half-complex form." + keyword arguments. If the vector is real, it is assumed to be in + half-complex form. If the length of the vector is a power of 2, use + of a non-radix-2 transform can be forced." (let ((pass-on-args (copy-list args))) (remf pass-on-args :half-complex) (remf pass-on-args :decimation-in-frequency) - (if (power-of-2-p (floor (size vector) stride)) + (remf pass-on-args :non-radix-2) + (if (and (not non-radix-2) (power-of-2-p (floor (size vector) stride))) (if (subtypep (element-type vector) 'real) (apply 'backward-fourier-transform-halfcomplex-radix2 vector pass-on-args) diff --git a/fast-fourier-transforms/example.lisp b/fast-fourier-transforms/example.lisp index 175bc67de117def6ac96b7710a08b9b4bc577b4b..a05d5d3bfccc723724e07c69572e8ea4adbd1342 100644 --- a/fast-fourier-transforms/example.lisp +++ b/fast-fourier-transforms/example.lisp @@ -146,32 +146,33 @@ (setf (grid:gref vector i) (coerce (/ (grid:gref vector i) length) element-type))) vector)) -(defun test-real-fft-noise (vector &key (stride 1)) +(defun test-real-fft-noise (vector &key (stride 1) non-radix-2) "Test forward and inverse FFT for a real vector, and return both results in unpacked form." (let* ((forward (forward-fourier-transform (realpart-vector vector :stride stride :init-offset 0) - :stride stride)) + :stride stride :non-radix-2 non-radix-2)) (inverse - (forward-fourier-transform (copy forward) :half-complex t :stride stride))) + (forward-fourier-transform (copy forward) :half-complex t :stride stride + :non-radix-2 non-radix-2))) (values (unpack forward :unpack-type 'complex :stride stride) (unpack (vector/length inverse :stride stride) :stride stride)))) -(defun test-complex-fft-noise (vector &key (stride 1)) +(defun test-complex-fft-noise (vector &key (stride 1) non-radix-2) "Test forward, inverse and backward FFT for a complex vector and return all three results." (let ((forward (forward-fourier-transform (copy-with-stride vector :stride stride :init-offset 2000) - :stride stride))) + :stride stride :non-radix-2 non-radix-2))) (values forward (inverse-fourier-transform (copy-with-stride forward :stride stride :init-offset 0) - :stride stride) + :stride stride :non-radix-2 non-radix-2) ;; in backward-fourier-transform, we could use copy instead of ;; ;; copy-with-stride (backward-fourier-transform (copy-with-stride forward :init-offset 2000) - :stride stride)))) + :stride stride :non-radix-2 non-radix-2)))) -(defun test-fft-noise (element-type size &key (stride 1)) +(defun test-fft-noise (element-type size &key (stride 1) non-radix-2) "A test of real forward and complex forward, reverse, and inverse FFT on random noise. Returns the result of the DFT forward Fourier transformation and the forward FFT, which should be the same, and the original vector and the @@ -188,7 +189,8 @@ :init-offset 3000)))) (if (subtypep element-type 'complex) (multiple-value-bind (forward inverse backward) - (test-complex-fft-noise random-vector :stride stride) + (test-complex-fft-noise random-vector :stride stride + :non-radix-2 non-radix-2) (values dft-random-vector ; DFT forward result for reference forward ; FFT forward result; should check @@ -196,7 +198,8 @@ inverse ; The inverse FFT applied to the forward result (vector/length backward :stride stride))) (multiple-value-bind (forward inverse) - (test-real-fft-noise random-vector :stride stride) + (test-real-fft-noise random-vector :stride stride + :non-radix-2 non-radix-2) (values dft-random-vector forward random-vector inverse))))) ;; (test-fft-noise '(complex double-float) 2 :stride 2) diff --git a/fast-fourier-transforms/forward.lisp b/fast-fourier-transforms/forward.lisp index b5eb8a3e64789dde7936bb6ef8851c45f1ff9d7c..b90182609112abf4a17807293a2994c26918e738 100644 --- a/fast-fourier-transforms/forward.lisp +++ b/fast-fourier-transforms/forward.lisp @@ -136,16 +136,19 @@ (export 'forward-fourier-transform) (defun forward-fourier-transform (vector &rest args - &key half-complex decimation-in-frequency (stride 1) &allow-other-keys) + &key half-complex decimation-in-frequency (stride 1) non-radix-2 &allow-other-keys) "Perform a forward fast Fourier transform on the given vector. If the length of the vector is not a power of 2, and the user has a suitable wavetable and/or workspace, these can be supplied as keyword arguments. If the (real) vector is in half-complex form, - then the key argument :half-complex should be non-NIL." + then the key argument :half-complex should be non-NIL. If the + length of the vector is a power of 2, use of a non-radix-2 transform + can be forced." (let ((pass-on-args (copy-list args))) (remf pass-on-args :half-complex) (remf pass-on-args :decimation-in-frequency) - (if (power-of-2-p (floor (size vector) stride)) + (remf pass-on-args :non-radix-2) + (if (and (not non-radix-2) (power-of-2-p (floor (size vector) stride))) (if half-complex (apply 'forward-fourier-transform-halfcomplex-radix2 vector pass-on-args) diff --git a/fast-fourier-transforms/inverse.lisp b/fast-fourier-transforms/inverse.lisp index 8ec19f30abf79dabc52ac482425b8a02e0a79e1a..77545e27e3449e921c1369d917eadeba57e28456 100644 --- a/fast-fourier-transforms/inverse.lisp +++ b/fast-fourier-transforms/inverse.lisp @@ -126,16 +126,18 @@ (export 'inverse-fourier-transform) (defun inverse-fourier-transform (vector &rest args - &key decimation-in-frequency (stride 1) &allow-other-keys) + &key decimation-in-frequency (stride 1) non-radix-2 &allow-other-keys) "Perform a inverse fast Fourier transform on the given vector. If the length of the vector is not a power of 2, and the user has a suitable wavetable and/or workspace, these can be supplied as keyword arguments. If the vector is real, it is assumed to be - in half-complex form." + in half-complex form. If the length of the vector is a power of 2, + use of a non-radix-2 transform can be forced." (let ((pass-on-args (copy-list args))) (remf pass-on-args :half-complex) (remf pass-on-args :decimation-in-frequency) - (if (power-of-2-p (floor (size vector) stride)) + (remf pass-on-args :non-radix-2) + (if (and (not non-radix-2) (power-of-2-p (floor (size vector) stride))) (if (subtypep (element-type vector) 'real) (apply 'inverse-fourier-transform-halfcomplex-radix2 vector pass-on-args) diff --git a/tests/fast-fourier-transform.lisp b/tests/fast-fourier-transform.lisp index 40209dacb0127737508134a4d0f9cb71f336ace7..dac5cda25708d500cfc06727e4576a1520da12f2 100644 --- a/tests/fast-fourier-transform.lisp +++ b/tests/fast-fourier-transform.lisp @@ -91,15 +91,21 @@ (eval-when (:compile-toplevel :load-toplevel :execute) (defun fft-test-forms (size stride) `((fft-real-result-check - (test-fft-noise 'double-float ,size :stride ,stride) double-float ,stride) + (test-fft-noise 'double-float ,size :stride ,stride :non-radix-2 t) double-float ,stride) (fft-real-result-check - (test-fft-noise 'single-float ,size :stride ,stride) single-float ,stride) + (test-fft-noise 'single-float ,size :stride ,stride :non-radix-2 t) single-float ,stride) (fft-complex-result-check (test-fft-noise - '(complex double-float) ,size :stride ,stride) double-float ,stride) + '(complex double-float) ,size :stride ,stride :non-radix-2 t) double-float ,stride) (fft-complex-result-check (test-fft-noise - '(complex single-float) ,size :stride ,stride) single-float ,stride)))) + '(complex single-float) ,size :stride ,stride :non-radix-2 t) single-float ,stride) + (when (power-of-2-p (floor ,size ,stride)) + (fft-real-result-check + (test-fft-noise 'double-float ,size :stride ,stride :non-radix-2 nil) double-float ,stride) + (fft-complex-result-check + (test-fft-noise + '(complex double-float) ,size :stride ,stride :non-radix-2 nil) double-float ,stride))))) (defmacro all-fft-test-forms (size-max stride-max &optional additional-single-stride) `(lisp-unit:define-test fast-fourier-transform @@ -113,5 +119,5 @@ (fft-test-forms size 1))))) -(all-fft-test-forms 20 3) +;(all-fft-test-forms 20 3) ;;; Tests commented out because they come out not so good.