--- /dev/null
+all: sizes
+
+sizes: sizes.c
+ gcc -I/opt/local/include sizes.c -L/opt/local/lib -lfftw3 -lm -o sizes
+
+clean:
+ rm -f *.o *~ sizes
+
--- /dev/null
+#include <stdio.h>
+#include <fftw3.h>
+
+main () {
+ printf("Size of fftw_complex: %d\n", sizeof (fftw_complex));
+ printf("Size of fftw_plan: %d\n", sizeof (fftw_plan));
+}
--- /dev/null
+;;;; -*- Mode: Lisp; Syntax: ANSI-Common-Lisp; Base: 10 -*-
+;;;; *************************************************************************
+;;;; FILE IDENTIFICATION
+;;;;
+;;;; Name: fftw3-tests.asd
+;;;; Purpose: ASDF system definitionf for fftw3 testing package
+;;;; Author: Kevin M. Rosenberg
+;;;; Date Started: Jan 2009
+;;;;
+;;;; $Id$
+;;;; *************************************************************************
+
+(defpackage #:fftw3-tests-system
+ (:use #:asdf #:cl))
+(in-package #:fftw3-tests-system)
+
+(defsystem cl-fftw3-tests
+ :depends-on (rt cl-fftw3)
+ :components
+ ((:file "tests")))
+
+(defmethod perform ((o test-op) (c (eql (find-system 'cl-fftw3-tests))))
+ (or (funcall (intern (symbol-name '#:do-tests)
+ (find-package '#:regression-test)))
+ (error "test-op failed")))
+
--- /dev/null
+;;;; -*- Mode: Lisp; Syntax: ANSI-Common-Lisp; Base: 10 -*-
+;;;; *************************************************************************
+;;;; FILE IDENTIFICATION
+;;;;
+;;;; Name: cl-fftw3.asd
+;;;; Purpose: ASDF system definition for FFTW3 package
+;;;; Programmer: Kevin M. Rosenberg
+;;;; Date Started: Jan 2009
+;;;;
+;;;; $Id$
+;;;;
+;;;; This file, part of FFTW3, is Copyright (c) 2009 by Kevin M. Rosenberg
+;;;;
+;;;; FFTW3 users are granted the rights to distribute and use this software
+;;;; as governed by the terms of the Lisp Lesser GNU Public License
+;;;; (http://opensource.franz.com/preamble.html), also known as the LLGPL.
+;;;; *************************************************************************
+
+(in-package #:cl-user)
+(defpackage #:fftw3-system (:use #:asdf #:cl))
+(in-package #:fftw3-system)
+
+(defsystem cl-fftw3
+ :name "fftw3"
+ :author "Kevin M. Rosenberg <kevin@rosenberg.net>"
+ :maintainer "Kevin M. Rosenberg <kmr@debian.org>"
+ :licence "LLGPL"
+ :depends-on (:kmrcl :cffi)
+ :components
+ ((:file "package")
+ (:file "specials" :depends-on ("package"))
+ (:file "common" :depends-on ("specials"))
+ (:file "ffi" :depends-on ("common"))
+ (:file "transform" :depends-on ("ffi"))
+ (:file "wisdom" :depends-on ("ffi"))))
+
+(defmethod perform ((o test-op) (c (eql (find-system 'cl-fftw3))))
+ (operate 'load-op 'cl-fftw3-tests)
+ (operate 'test-op 'cl-fftw3-tests :force t))
--- /dev/null
+;;;; -*- Mode: Lisp; Syntax: ANSI-Common-Lisp; Base: 10 -*-
+;;;; *************************************************************************
+;;;; FILE IDENTIFICATION
+;;;;
+;;;; Name: common.lisp
+;;;; Purpose: Common functions for FFTW3 package
+;;;; Programmer: Kevin M. Rosenberg
+;;;; Date Started: Jan 2009
+;;;;
+;;;; $Id$
+;;;;
+;;;; This file, part of FFTW3, is Copyright (c) 2009 by Kevin M. Rosenberg
+;;;;
+;;;; FFTW3 users are granted the rights to distribute and use this software
+;;;; as governed by the terms of the Lisp Lesser GNU Public License
+;;;; (http://opensource.franz.com/preamble.html), also known as the LLGPL.
+;;;; *************************************************************************
+
+(in-package #:fftw3)
+
+(defclass fftw-multi ()
+ ((plan :initarg :plan :initform nil :accessor plan)))
+
+(defclass fftw-r2r-1d-multi (fftw-multi)
+ ((in-count :initarg :in-count :initform 0 :type fixnum :accessor in-count)
+ (out-c :initarg :out-c :initform nil :accessor out-c)
+ (in-cf :initarg :in-cf :initform nil :accessor in-cf)
+ (out-cf :initarg :out-cf :initform nil :accessor out-cf)))
+
+(defclass fftw-r2c-1d-multi (fftw-multi)
+ ((in-count :initarg :in-count :initform 0 :type fixnum :accessor in-count)
+ (out-count :initarg :out-count :initform 0 :type fixnum :accessor out-count)
+ (complex-output :initarg :complex-output :initform t :type boolean :accessor complex-output)
+ (in-cf :initarg :in-cf :initform nil :accessor in-cf)
+ (out-cf :initarg :out-cf :initform nil :accessor out-cf)
+ (out :initarg :out :initform nil :accessor out)))
+
+(defclass fftw-1d-multi (fftw-multi)
+ ((in-count :initarg :in-count :initform 0 :type fixnum :accessor in-count)
+ (out-count :initarg :out-count :initform 0 :type fixnum :accessor out-count)
+ (in-cf :initarg :in-cf :initform nil :accessor in-cf)
+ (out-cf :initarg :out-cf :initform nil :accessor out-cf)
+ (out :initarg :out :initform nil :accessor out)))
+
+(defclass fftw-c2r-1d-multi (fftw-multi)
+ ((plan-even :initarg :plan-even :initform nil :accessor plan-even)
+ (plan-odd :initarg :plan-odd :initform nil :accessor plan-odd)
+ (in-count :initarg :in-count :initform 0 :type fixnum :accessor in-count)
+ (out-count-odd :initarg :out-count-odd :initform 0 :type fixnum :accessor out-count-odd)
+ (out-count-even :initarg :out-count-even :initform 0 :type fixnum :accessor out-count-even)
+ (in-cf :initarg :in-cf :initform nil :accessor in-cf)
+ (out-cf :initarg :out-cf :initform nil :accessor out-cf)
+ (out-odd :initarg :out-odd :initform nil :accessor out-odd)
+ (out-even :initarg :out-even :initform nil :accessor out-even)))
+
+(defun normalize-input-range (in start count)
+ "Returns the start and count of an input range with range checking."
+ (declare #.*standard-optimize-settings*)
+ (unless (vectorp in)
+ (error "Input needs to be a vector."))
+ (when (or (not (integerp start)) (minusp start))
+ (error "Start must be an integer, got ~A." start))
+ (let ((len (length in)))
+ (declare (fixnum len))
+ (when (>= start len)
+ (error "Start (~D) greater than last element (~D)." start (1- len)))
+
+ (cond
+ ((null count)
+ (setq count (- len start)))
+ ((or (not (integerp count)) (minusp count))
+ (error "Count must be positive integer or nil, got ~A." count))
+ ((> (+ start count) len)
+ (error "Requesting element (~D) past end of vector (~D)." (+ start count) len)))
+ (values start count)))
+
+(defun complex-is-even-p (v)
+ "Returns T if first and last members have zero imaginary component."
+ ;; Do not use *standard-optimize-settings-here* -- they fail on lispworks 5.1.2
+ (declare (optimize (speed 3) (space 0)))
+ (and (zerop (imagpart (aref v 0)))
+ (zerop (imagpart (aref v (1- (length v)))))))
+
+(defun complex-to-mag-phase (cmplx &optional (elements-complex t))
+ (declare #.*standard-optimize-settings*)
+ (let* ((len (length cmplx))
+ (out-len (if elements-complex len (/ len 2)))
+ (mag (make-array out-len :element-type 'double-float))
+ (phase (make-array out-len :element-type 'double-float)))
+ (declare (fixnum len out-len))
+ (cond
+ (elements-complex
+ (dotimes (i out-len)
+ (declare (fixnum i))
+ (let* ((e (aref cmplx i)))
+ (setf (aref mag i) (abs e))
+ (setf (aref phase i) (phase e)))))
+ (t
+ (dotimes (i out-len)
+ (declare (fixnum i))
+ (let* ((base (+ i i))
+ (r (aref cmplx base))
+ (i (aref cmplx (1+ base))))
+ (setf (aref mag i) (sqrt (+ (* r r) (* i i))))
+ (setf (aref phase i) (atan i r))))))
+ (values mag phase)))
+
+(defun hc-to-mag-phase (hc)
+ "Turns half-complex into magniture and phase vectors."
+ (declare #.*standard-optimize-settings*)
+ (let* ((n (length hc))
+ (out-n (1+ (floor n 2)))
+ (mag (make-array out-n :element-type 'double-float))
+ (phase (make-array out-n :element-type 'double-float)))
+ (declare (fixnum n out-n))
+ (dotimes (i out-n)
+ (cond
+ ((or (zerop i)
+ (and (evenp n) (eql i (1- out-n))))
+ (setf (aref mag i) (aref hc (/ i 2)))
+ (setf (aref phase i) 0d0))
+ (t
+ (let ((re (aref hc i))
+ (im (aref hc (- n i))))
+ (setf (aref mag i) (sqrt (+ (* re re) (* im im))))
+ (setf (aref phase i) (atan im re))))))
+ (values mag phase)))
+
--- /dev/null
+;;;; -*- Mode: Lisp; Syntax: ANSI-Common-Lisp; Base: 10 -*-
+;;;; *************************************************************************
+;;;; FILE IDENTIFICATION
+;;;;
+;;;; Name: ffi.lisp
+;;;; Purpose: CFFI interface for FFTW3 package
+;;;; Programmer: Kevin M. Rosenberg
+;;;; Date Started: Mar 2009
+;;;;
+;;;; $Id$
+;;;;
+;;;; This file, part of FFTW3, is Copyright (c) 2009 by Kevin M. Rosenberg
+;;;;
+;;;; FFTW3 users are granted the rights to distribute and use this software
+;;;; as governed by the terms of the Lisp Lesser GNU Public License
+;;;; (http://opensource.franz.com/preamble.html), also known as the LLGPL.
+;;;; *************************************************************************
+
+(in-package #:fftw3)
+
+(cffi:define-foreign-library fftw3
+ ((or :darwin :macosx) (:or #P"/opt/local/lib/libfftw3.dylib"
+ #P"/usr/lib/libfftw3.dylib"))
+ (:linux (:or #P"/usr/lib/libfftw3.so"
+ #P"/usr/local/lib/libfftw3.so"))
+ (t (:default "libfftw3")))
+
+(cffi:use-foreign-library fftw3)
+
+(cffi:defcstruct fftw-complex-struct
+ "Complex number structure."
+ (re :double)
+ (im :double))
+
+(cffi:defctype fftw-plan :pointer)
+
+(declaim (inline fftw-plan-dft-1d))
+(cffi:defcfun ("fftw_plan_dft_1d" fftw-plan-dft-1d) fftw-plan
+ (n :int)
+ (in (:pointer fftw-complex-struct))
+ (out (:pointer fftw-complex-struct))
+ (sign :int)
+ (flags :uint))
+
+(declaim (inline fftw-plan-dft-2d))
+(cffi:defcfun ("fftw_plan_dft_2d" fftw-plan-dft-2d) fftw-plan
+ (n0 :int)
+ (n1 :int)
+ (in (:pointer fftw-complex-struct))
+ (out (:pointer fftw-complex-struct))
+ (sign :int)
+ (flags :uint))
+
+(declaim (inline fftw-plan-dft-3d))
+(cffi:defcfun ("fftw_plan_dft_3d" fftw-plan-dft-3d) fftw-plan
+ (n0 :int)
+ (n1 :int)
+ (n2 :int)
+ (in (:pointer fftw-complex-struct))
+ (out (:pointer fftw-complex-struct))
+ (sign :int)
+ (flags :uint))
+
+(declaim (inline fftw-plan-dft))
+(cffi:defcfun ("fftw_plan_dft" fftw-plan-dft) fftw-plan
+ (rank :int)
+ (n (:pointer :int))
+ (in (:pointer fftw-complex-struct))
+ (out (:pointer fftw-complex-struct))
+ (sign :int)
+ (flags :uint))
+
+(declaim (inline fftw-plan-r2r-1d))
+(cffi:defcfun ("fftw_plan_r2r_1d" fftw-plan-r2r-1d) fftw-plan
+ (n :int)
+ (in :pointer)
+ (out :pointer)
+ (kind :int)
+ (flags :uint))
+
+(declaim (inline fftw-plan-r2r-2d))
+(cffi:defcfun ("fftw_plan_r2r_2d" fftw-plan-r2r-2d) fftw-plan
+ (n0 :int)
+ (n1 :int)
+ (in :pointer)
+ (out :pointer)
+ (kind :int)
+ (flags :uint))
+
+(declaim (inline fftw-plan-r2r-3d))
+(cffi:defcfun ("fftw_plan_r2r_3d" fftw-plan-r2r-3d) fftw-plan
+ (n0 :int)
+ (n1 :int)
+ (n2 :int)
+ (in :pointer)
+ (out :pointer)
+ (kind :int)
+ (flags :uint))
+
+(declaim (inline fftw-plan-dft-r2c-1d))
+(cffi:defcfun ("fftw_plan_dft_r2c_1d" fftw-plan-dft-r2c-1d) fftw-plan
+ (n :int)
+ (in :pointer)
+ (out :pointer)
+ (flags :uint))
+
+(declaim (inline fftw-plan-dft-r2c-2d))
+(cffi:defcfun ("fftw_plan_dft_r2c_2d" fftw-plan-dft-r2c-2d) fftw-plan
+ (n0 :int)
+ (n1 :int)
+ (in :pointer)
+ (out :pointer)
+ (flags :uint))
+
+(declaim (inline fftw-plan-dft-r2c-3d))
+(cffi:defcfun ("fftw_plan_dft_r2c_3d" fftw-plan-dft-r2c-3d) fftw-plan
+ (n0 :int)
+ (n1 :int)
+ (n2 :int)
+ (in :pointer)
+ (out :pointer)
+ (flags :uint))
+
+(declaim (inline fftw-plan-dft-r2c))
+(cffi:defcfun ("fftw_plan_dft_r2c" fftw-plan-dft-r2c) fftw-plan
+ (rank :int)
+ (n (:pointer :int))
+ (in :pointer)
+ (out :pointer)
+ (flags :uint))
+
+(declaim (inline fftw-plan-dft-c2r-1d))
+(cffi:defcfun ("fftw_plan_dft_c2r_1d" fftw-plan-dft-c2r-1d) fftw-plan
+ (n :int)
+ (in :pointer)
+ (out :pointer)
+ (flags :uint))
+
+(declaim (inline fftw-plan-dft-c2r-2d))
+(cffi:defcfun ("fftw_plan_dft_c2r_2d" fftw-plan-dft-c2r-2d) fftw-plan
+ (n0 :int)
+ (n1 :int)
+ (in :pointer)
+ (out :pointer)
+ (flags :uint))
+
+(declaim (inline fftw-plan-dft-c2r-3d))
+(cffi:defcfun ("fftw_plan_dft_c2r_3d" fftw-plan-dft-c2r-3d) fftw-plan
+ (n0 :int)
+ (n1 :int)
+ (n2 :int)
+ (in :pointer)
+ (out :pointer)
+ (flags :uint))
+
+(declaim (inline fftw-plan-dft-c2r))
+(cffi:defcfun ("fftw_plan_dft_c2r" fftw-plan-dft-c2r) fftw-plan
+ (rank :int)
+ (n (:pointer :int))
+ (in :pointer)
+ (out :pointer)
+ (flags :uint))
+
+(declaim (inline fftw-plan-many-dft))
+(cffi:defcfun ("fftw_plan_many_dft" fftw-plan-many-dft) fftw-plan
+ (rank :int)
+ (n (:pointer :int))
+ (howmany :int)
+ (in (:pointer fftw-complex-struct))
+ (inembed (:pointer :int))
+ (istride :int)
+ (idist :int)
+ (out (:pointer fftw-complex-struct))
+ (onembed (:pointer :int))
+ (ostride :int)
+ (odist :int)
+ (sign :int)
+ (flags :uint))
+
+(declaim (inline fftw-plan-many-dft_r2c))
+(cffi:defcfun ("fftw_plan_many_dft_r2c" fftw-plan-many-dft-r2c) fftw-plan
+ (rank :int)
+ (n (:pointer :int))
+ (howmany :int)
+ (in (:pointer :double))
+ (inembed (:pointer :int))
+ (istride :int)
+ (idist :int)
+ (out (:pointer fftw-complex-struct))
+ (onembed (:pointer :int))
+ (ostride :int)
+ (odist :int)
+ (flags :uint))
+
+(declaim (inline fftw-plan-many-dft_c2r))
+(cffi:defcfun ("fftw_plan_many_dft_c2r" fftw-plan-many-dft-c2r) fftw-plan
+ (rank :int)
+ (n (:pointer :int))
+ (howmany :int)
+ (in (:pointer fftw-complex-struct))
+ (inembed (:pointer :int))
+ (istride :int)
+ (idist :int)
+ (out (:pointer :double))
+ (onembed (:pointer :int))
+ (ostride :int)
+ (odist :int)
+ (flags :uint))
+
+(declaim (inline fftw-plan-many-dft_r2r))
+(cffi:defcfun ("fftw_plan_many_dft_r2r" fftw-plan-many-dft-r2r) fftw-plan
+ (rank :int)
+ (n (:pointer :int))
+ (howmany :int)
+ (in (:pointer :double))
+ (inembed (:pointer :int))
+ (istride :int)
+ (idist :int)
+ (out (:pointer :double))
+ (onembed (:pointer :int))
+ (ostride :int)
+ (odist :int)
+ (kind (:pointer :int))
+ (flags :uint))
+
+
+(declaim (inline fftw-malloc))
+(cffi:defcfun ("fftw_malloc" fftw-malloc) (:pointer fftw-complex-struct)
+ (n :int))
+
+(declaim (inline fftw-execute))
+(cffi:defcfun ("fftw_execute" fftw-execute) :void
+ (p fftw-plan))
+
+(declaim (inline fftw-destroy-plan))
+(cffi:defcfun ("fftw_destroy_plan" fftw-destroy-plan) :void
+ (p fftw-plan))
+
+(declaim (inline fftw-free))
+(cffi:defcfun ("fftw_free" fftw-free) :void
+ (p :pointer))
+
+
+;;; Wisdom functions
+
+(cffi:defcfun ("fftw_import_wisdom_from_string" fftw-import-wisdom-from-string) :int
+ (input-string :string))
+
+(cffi:defcfun ("fftw_export_wisdom_to_string" fftw-export-wisdom-to-string) :string+ptr
+ )
+
+(cffi:defcfun free :void
+ (ptr :pointer))
--- /dev/null
+;;;; -*- Mode: Lisp; Syntax: ANSI-Common-Lisp; Base: 10 -*-
+;;;; *************************************************************************
+;;;; FILE IDENTIFICATION
+;;;;
+;;;; Name: package.lisp
+;;;; Purpose: Package definition for fftw3 package
+;;;; Programmer: Kevin M. Rosenberg
+;;;; Date Started: Jan 2009
+;;;;
+;;;; $Id$
+;;;;
+;;;; This file, part of FFTW3, is Copyright (c) 2009 by Kevin M. Rosenberg
+;;;;
+;;;; FFTW3 users are granted the rights to distribute and use this software
+;;;; as governed by the terms of the Lisp Lesser GNU Public License
+;;;; (http://opensource.franz.com/preamble.html), also known as the LLGPL.
+;;;; *************************************************************************
+
+(in-package #:cl-user)
+
+(defpackage #:cl-fftw3
+ (:nicknames #:fftw3)
+ (:use #:common-lisp)
+ (:export
+ ;; transform.lisp
+ #:fftw-r2hc-1d
+ #:make-fftw-r2r-1d-multi
+ #:fftw-r2r-1d-multi
+ #:hc-to-mag-phase
+ #:fftw-hc2r-1d
+ #:fftw-r2c-1d
+ #:make-fftw-r2c-1d-multi
+ #:make-fftw-c2r-1d-multi
+ #:destroy-fftw-multi
+ #:fftw-r2c-1d-multi
+ #:fftw-c2r-1d-multi
+ #:complex-to-mag-phase
+ #:fftw-c2r-1d
+ #:fftw-c-1d
+
+ ;; wisdom.lisp
+ #:import-user-wisdom
+ #:export-user-wisdom
+
+ ;; specials
+ #:+fftw-forward+
+ #:+fftw-backward+
+ #:+fftw-hc2r+
+ #:+fftw-r2hc+
+
+ #:+fftw-measure+
+ #:+fftw-destroy-input+
+ #:+fftw-unaligned+
+ #:+fftw-conservative-memory+
+ #:+fftw-exhaustive+
+ #:+fftw-preserve-input+
+ #:+fftw-patiento+
+ #:+fftw-estimate+
+ ))
--- /dev/null
+;;;; -*- Mode: Lisp; Syntax: ANSI-Common-Lisp; Base: 10 -*-
+;;;; *************************************************************************
+;;;; FILE IDENTIFICATION
+;;;;
+;;;; Name: specials.lisp
+;;;; Purpose: Special data declarations for FFTW3 package
+;;;; Programmer: Kevin M. Rosenberg
+;;;; Date Started: March 2009
+;;;;
+;;;; $Id$
+;;;;
+;;;; This file, part of FFTW3, is Copyright (c) 2009 by Kevin M. Rosenberg
+;;;;
+;;;; FFTW3 users are granted the rights to distribute and use this software
+;;;; as governed by the terms of the Lisp Lesser GNU Public License
+;;;; (http://opensource.franz.com/preamble.html), also known as the LLGPL.
+;;;; *************************************************************************
+
+(in-package #:fftw3)
+
+(defconstant +fftw-forward+ -1)
+(defconstant +fftw-backward+ +1)
+
+(defconstant +fftw-r2hc+ 0)
+(defconstant +fftw-hc2r+ 1)
+
+(defconstant +fftw-measure+ 0)
+(defconstant +fftw-destroy-input+ (ash 1 0))
+(defconstant +fftw-unaligned+ (ash 1 1))
+(defconstant +fftw-conservative-memory+ (ash 1 2))
+(defconstant +fftw-exhaustive+ (ash 1 3))
+(defconstant +fftw-preserve-input+ (ash 1 4))
+(defconstant +fftw-patiento+ (ash 1 5))
+(defconstant +fftw-estimate+ (ash 1 6))
+
+(defparameter *user-wisdom-file*
+ (concatenate 'string (namestring (user-homedir-pathname)) ".fftw3-wisdom"))
+
+(defvar *standard-optimize-settings*
+ '(optimize
+ speed
+ (safety 0)
+ (space 0)
+ (debug 1)
+ (compilation-speed 0)
+ #+:lispworks (hcl:fixnum-safety 0))
+ "The standard optimize settings used by most declaration expressions.")
--- /dev/null
+(defpackage #:cl-fftw3-tests
+ (:nicknames #:fftw3-tests)
+ (:use #:cl #:rtest #:cl-fftw3))
+
+(in-package #:cl-fftw3-tests)
+
+(defun double-float-close-p (d1 d2 &key (epsilon 1d-7))
+ (when (and (zerop d1) (zerop d2))
+ (return-from double-float-close-p t))
+
+ (let* ((delta (abs (- d1 d2)))
+ (eps (/ delta (max (abs d1) (abs d2)))))
+ (when (or (and (zerop (min (abs d1) (abs d2))) (< (max (abs d1) (abs d2)) epsilon))
+ (< eps epsilon))
+ t)))
+
+(defun vectors-close-p (v1 v2)
+ (unless (eql (length v1) (length v2))
+ (return-from vectors-close-p nil))
+
+ (dotimes (i (length v1))
+ (unless (and (double-float-close-p (realpart (aref v1 i)) (realpart (aref v2 i)))
+ (double-float-close-p (imagpart (aref v1 i)) (imagpart (aref v2 i))))
+ (return-from vectors-close-p nil)))
+ t)
+
+(rem-all-tests)
+
+(defparameter *a7* #(0 1 2 3 2 1 0))
+(defparameter *a8* #(0 1 2 3 3 2 1 0))
+(defparameter *ac7* #(#C(0 0) #C(1 1) #C(2 2) #C(3 3) #C(2 2) #C(1 1) #C(0 0)))
+(defparameter *ac8* #(#C(0 0) #C(1 1) #C(2 2) #C(3 3) #C(3 3) #C(2 2) #C(1 1) #C(0 0)))
+(defparameter *cc7* #(#C(0 1) #C(1 1) #C(2 2) #C(3 3) #C(2 2) #C(1 1) #C(0 0)))
+(defparameter *cc8* #(#C(0 1) #C(1 1) #C(2 2) #C(3 3) #C(3 3) #C(2 2) #C(1 1) #C(0 0)))
+
+(defparameter *b7* #(3 2 1 0 1 2 3))
+(defparameter *b8* #(3 2 1 0 0 1 2 3))
+(defparameter *bc7* #(#C(3 3) #C(2 2) #C(1 1) #C(0 0) #C(1 1) #C(2 2) #C(3 3)))
+(defparameter *bc8* #(#C(3 3) #C(2 2) #C(1 1) #C(0 0) #C(0 0) #C(1 1) #C(2 2) #C(3 3)))
+
+(defparameter *a7hc* (fftw-r2hc-1d *a7*))
+(defparameter *a8hc* (fftw-r2hc-1d *a8*))
+(defparameter *a7c* (fftw-r2c-1d *a7*))
+(defparameter *a8c* (fftw-r2c-1d *a8*))
+(defparameter *ac7c* (fftw-c-1d *a7*))
+(defparameter *ac8c* (fftw-c-1d *a8*))
+
+(defparameter *a7c-octave* #(#C(1.28571428571428580945d0 0d0)
+ #C(-0.64984533421747225912d0 -0.31294901910963018876d0)
+ #C(0.02743163880429938875d0 0.03439818705768063478d0)
+ #C(-0.02044344744397007946d0 -0.08956859554733599682d0)))
+
+(defparameter *a8c-octave* #(#C(1.50000000000000000000d0 0.00000000000000000000d0)
+ #C(-0.72855339059327373086d0 -0.30177669529663686543d0)
+ #C(0.00000000000000000000d0 0.00000000000000000000d0)
+ #C(-0.02144660940672621363d0 -0.05177669529663689318d0)
+ #C(0.00000000000000000000d0 0.00000000000000000000d0)))
+
+(defparameter *ac7c-octave* #(#C(1.28571428571428580945d0 1.28571428571428580945d0)
+ #C(-0.33689631510784201485d0 -0.96279435332710250339d0)
+ #C(-0.00696654825338124516d0 0.06182982586198002700d0)
+ #C(0.06912514810336592430d0 -0.11001204299130606934d0)
+ #C(-0.11001204299130606934d0 0.06912514810336592430d0)
+ #C(0.06182982586198002700d0 -0.00696654825338124516d0)
+ #C(-0.96279435332710250339d0 -0.33689631510784201485d0)))
+
+(defparameter *cc7c-octave* #(#C(1.28571428571428580945d0 1.42857142857142860315d0)
+ #C(-0.33689631510784201485d0 -0.81993721046995970969d0)
+ #C(-0.00696654825338124516d0 0.20468696871912286928d0)
+ #C(0.06912514810336592430d0 0.03284509986583677293d0)
+ #C(-0.11001204299130606934d0 0.21198229096050877351d0)
+ #C(0.06182982586198002700d0 0.13589059460376159971d0)
+ #C(-0.96279435332710250339d0 -0.19403917225069916563d0)))
+
+(defparameter *ac8c-octave* #(#C(1.50000000000000000000d0 1.50000000000000000000d0)
+ #C(-0.42677669529663686543d0 -1.03033008588991070731d0)
+ #C(0.00000000000000000000d0 0.00000000000000000000d0)
+ #C(0.03033008588991070731d0 -0.07322330470336310682d0)
+ #C(0.00000000000000000000d0 0.00000000000000000000d0)
+ #C(-0.07322330470336310682d0 0.03033008588991070731d0)
+ #C(0.00000000000000000000d0 0.00000000000000000000d0)
+ #C(-1.03033008588991070731d0 -0.42677669529663686543d0)))
+
+(defparameter *cc8c-octave* #(#C(1.50000000000000000000d0 1.62500000000000000000d0)
+ #C(-0.42677669529663686543d0 -0.90533008588991070731)
+ #C(0.00000000000000000000d0 0.12500000000000000000d0)
+ #C(0.03033008588991070731d0 0.05177669529663689318d0)
+ #C(0.00000000000000000000d0 0.12500000000000000000d0)
+ #C(-0.07322330470336310682d0 0.15533008588991070731d0)
+ #C(0.00000000000000000000d0 0.12500000000000000000d0)
+ #C(-1.03033008588991070731d0 -0.30177669529663686543d0)))
+
+(deftest :len.1 (length *a7*) 7)
+(deftest :len.2 (length *a8*) 8)
+(deftest :len.3 (length *a7hc*) 7)
+(deftest :len.4 (length *a8hc*) 8)
+(deftest :len.5 (length *a7c*) 4)
+(deftest :len.6 (length *a8c*) 5)
+(deftest :len.7 (length *b7*) 7)
+(deftest :len.8 (length *b8*) 8)
+
+(deftest :eql.1 (vectors-close-p *a7* (fftw-hc2r-1d *a7hc*)) t)
+(deftest :eql.2 (vectors-close-p *a8* (fftw-hc2r-1d *a8hc*)) t)
+(deftest :eql.3 (vectors-close-p *a7* (fftw-c2r-1d *a7c*)) t)
+(deftest :eql.4 (vectors-close-p *a8* (fftw-c2r-1d *a8c*)) t)
+(deftest :eql.5 (vectors-close-p *a7* (fftw-c-1d *ac7c* :normalize nil :direction +fftw-backward+)) t)
+(deftest :eql.6 (vectors-close-p *a8* (fftw-c-1d *ac8c* :normalize nil :direction +fftw-backward+)) t)
+
+(deftest :eqlc.1 (vectors-close-p
+ *ac7*
+ (fftw-c-1d (fftw-c-1d *ac7* :normalize t :direction +fftw-forward+)
+ :normalize nil :direction +fftw-backward+))
+ t)
+
+(deftest :multi.1
+ (let* ((multi-r2hc (make-fftw-r2r-1d-multi (length *a7*) +fftw-r2hc+))
+ (multi-hc2r (make-fftw-r2r-1d-multi (length *a7*) +fftw-hc2r+))
+ (a7hc (copy-seq (fftw-r2r-1d-multi multi-r2hc *a7*)))
+ (b7hc (copy-seq (fftw-r2r-1d-multi multi-r2hc *b7*)))
+ (a7r (copy-seq (fftw-r2r-1d-multi multi-hc2r a7hc :normalize nil)))
+ (b7r (copy-seq (fftw-r2r-1d-multi multi-hc2r b7hc :normalize nil)))
+ (eq (and (vectors-close-p *a7* a7r) (vectors-close-p *b7* b7r) t)))
+ (destroy-fftw-multi multi-r2hc)
+ (destroy-fftw-multi multi-hc2r)
+ eq)
+ t)
+
+(deftest :multi.2
+ (let* ((multi-r2hc (make-fftw-r2r-1d-multi (length *a8*) +fftw-r2hc+))
+ (multi-hc2r (make-fftw-r2r-1d-multi (length *a8*) +fftw-hc2r+))
+ (a8hc (copy-seq (fftw-r2r-1d-multi multi-r2hc *a8*)))
+ (b8hc (copy-seq (fftw-r2r-1d-multi multi-r2hc *b8*)))
+ (a8r (copy-seq (fftw-r2r-1d-multi multi-hc2r a8hc :normalize nil)))
+ (b8r (copy-seq (fftw-r2r-1d-multi multi-hc2r b8hc :normalize nil)))
+ (eq (and (vectors-close-p *a8* a8r) (vectors-close-p *b8* b8r) t)))
+ (destroy-fftw-multi multi-r2hc)
+ (destroy-fftw-multi multi-hc2r)
+ eq)
+ t)
+
+(deftest :multi.3
+ (let* ((multi-r2c (make-fftw-r2c-1d-multi (length *a7*)))
+ (a7c (copy-seq (fftw-r2c-1d-multi multi-r2c *a7*)))
+ (b7c (copy-seq (fftw-r2c-1d-multi multi-r2c *b7*)))
+ (a7 (fftw-c2r-1d a7c))
+ (b7 (fftw-c2r-1d b7c))
+ (eq (and (vectors-close-p *a7* a7)
+ (vectors-close-p *b7* b7))))
+ (destroy-fftw-multi multi-r2c)
+ eq)
+ t)
+
+(deftest :multi.4
+ (let* ((multi-r2c (make-fftw-r2c-1d-multi (length *a8*)))
+ (a8c (copy-seq (fftw-r2c-1d-multi multi-r2c *a8*)))
+ (b8c (copy-seq (fftw-r2c-1d-multi multi-r2c *b8*)))
+ (a8 (fftw-c2r-1d a8c))
+ (b8 (fftw-c2r-1d b8c))
+ (eq (and (vectors-close-p *a8* a8)
+ (vectors-close-p *b8* b8))))
+ (destroy-fftw-multi multi-r2c)
+ eq)
+ t)
+
+(deftest :multi.5
+ (let* ((multi-r2c (make-fftw-r2c-1d-multi (length *a7*)))
+ (a7c (copy-seq (fftw-r2c-1d-multi multi-r2c *a7*)))
+ (b7c (copy-seq (fftw-r2c-1d-multi multi-r2c *b7*)))
+ (multi-c2r (make-fftw-c2r-1d-multi (length a7c)))
+ (a7r (copy-seq (fftw-c2r-1d-multi multi-c2r a7c)))
+ (b7r (copy-seq (fftw-c2r-1d-multi multi-c2r b7c)))
+ (eq (and (vectors-close-p *a7* a7r) (vectors-close-p *b7* b7r) t)))
+ (destroy-fftw-multi multi-r2c)
+ (destroy-fftw-multi multi-c2r)
+ eq)
+ t)
+
+(deftest :multi.6
+ (let* ((multi-r2c (make-fftw-r2c-1d-multi (length *a8*)))
+ (a8c (copy-seq (fftw-r2c-1d-multi multi-r2c *a8*)))
+ (b8c (copy-seq (fftw-r2c-1d-multi multi-r2c *b8*)))
+ (multi-c2r (make-fftw-c2r-1d-multi (length a8c)))
+ (a8r (copy-seq (fftw-c2r-1d-multi multi-c2r a8c)))
+ (b8r (copy-seq (fftw-c2r-1d-multi multi-c2r b8c)))
+ (eq (and (vectors-close-p *a8* a8r) (vectors-close-p *b8* b8r) t)))
+ (destroy-fftw-multi multi-r2c)
+ (destroy-fftw-multi multi-c2r)
+ eq)
+ t)
+
+(deftest :c2c.1
+ (let* ((ac7c (fftw-c-1d *ac7*))
+ (bc7c (fftw-c-1d *bc7*))
+ (cc7c (fftw-c-1d *cc7*))
+ (ac7i (fftw-c-1d ac7c :direction +fftw-backward+ :normalize nil))
+ (bc7i (fftw-c-1d bc7c :direction +fftw-backward+ :normalize nil))
+ (cc7i (fftw-c-1d cc7c :direction +fftw-backward+ :normalize nil))
+ (eq (and (vectors-close-p *ac7* ac7i) (vectors-close-p *bc7* bc7i)
+ (vectors-close-p *cc7* cc7i))))
+ eq)
+ t)
+
+(deftest :c2c.2
+ (let* ((ac8c (fftw-c-1d *ac8*))
+ (bc8c (fftw-c-1d *bc8*))
+ (cc8c (fftw-c-1d *cc8*))
+ (ac8i (fftw-c-1d ac8c :direction +fftw-backward+ :normalize nil))
+ (bc8i (fftw-c-1d bc8c :direction +fftw-backward+ :normalize nil))
+ (cc8i (fftw-c-1d cc8c :direction +fftw-backward+ :normalize nil))
+ (eq (and (vectors-close-p *ac8* ac8i) (vectors-close-p *bc8* bc8i)
+ (vectors-close-p *cc8* cc8i))))
+ eq)
+ t)
+
+(deftest :r2c-octave.1
+ (vectors-close-p *a7c* *a7c-octave*) t)
+
+(deftest :r2c-octave.2
+ (vectors-close-p *a8c* *a8c-octave*) t)
+
+(deftest :c2c-octave.1
+ (vectors-close-p (fftw-c-1d *ac7*) *ac7c-octave*) t)
+
+(deftest :c2c-octave.2
+ (vectors-close-p (fftw-c-1d *cc7*) *cc7c-octave*) t)
+
+(deftest :c2c-octave.3
+ (vectors-close-p (fftw-c-1d *ac8*) *ac8c-octave*) t)
+
+(deftest :c2c-octave.4
+ (vectors-close-p (fftw-c-1d *cc8*) *cc8c-octave*) t)
--- /dev/null
+;;;; -*- Mode: Lisp; Syntax: ANSI-Common-Lisp; Base: 10 -*-
+;;;; *************************************************************************
+;;;; FILE IDENTIFICATION
+;;;;
+;;;; Name: transform.lisp
+;;;; Purpose: Transformation functions
+;;;; Programmer: Kevin M. Rosenberg
+;;;; Date Started: March 2009
+;;;;
+;;;; $Id$
+;;;;
+;;;; This file, part of FFTW3, is Copyright (c) 2009 by Kevin M. Rosenberg
+;;;;
+;;;; FFTW3 users are granted the rights to distribute and use this software
+;;;; as governed by the terms of the Lisp Lesser GNU Public License
+;;;; (http://opensource.franz.com/preamble.html), also known as the LLGPL.
+;;;; *************************************************************************
+
+(in-package #:fftw3)
+
+(defconstant +size-double+ 8)
+
+(defun fftw-r2hc-1d (in &key (start 0) (count nil) (normalize t))
+ (declare #.*standard-optimize-settings*)
+ (multiple-value-bind (start count) (normalize-input-range in start count)
+ (declare (fixnum start count))
+
+ (let ((in-pos start)
+ (out-c (make-array count :element-type 'double-float))
+ (in-cf (fftw-malloc (* count +size-double+)))
+ (out-cf (fftw-malloc (* count +size-double+))))
+ (declare (fixnum in-pos))
+ (dotimes (i count)
+ (declare (fixnum i))
+ (setf (cffi:mem-aref in-cf :double i) (coerce (aref in in-pos) 'double-float))
+ (incf in-pos))
+
+ (let ((p (fftw-plan-r2r-1d count in-cf out-cf +fftw-r2hc+
+ (logior +fftw-estimate+ +fftw-destroy-input+))))
+ (fftw-execute p)
+ (fftw-destroy-plan p))
+
+ (if normalize
+ (let ((factor (coerce (/ 1 count) 'double-float)))
+ (declare (double-float factor))
+ (dotimes (i count)
+ (declare (fixnum i))
+ (setf (aref out-c i) (* factor (cffi:mem-aref out-cf :double i)))))
+ (dotimes (i count)
+ (declare (fixnum i))
+ (setf (aref out-c i) (cffi:mem-aref out-cf :double i))))
+
+ (fftw-free in-cf)
+ (fftw-free out-cf)
+
+ out-c)))
+
+(defun make-fftw-r2r-1d-multi (count type)
+ (declare #.*standard-optimize-settings*)
+ (declare (fixnum count type))
+ (let* ((out-c (make-array count :element-type 'double-float))
+ (in-cf (fftw-malloc (* count +size-double+)))
+ (out-cf (fftw-malloc (* count +size-double+)))
+ (plan (fftw-plan-r2r-1d count in-cf out-cf type
+ (logior +fftw-estimate+ +fftw-destroy-input+))))
+ (make-instance 'fftw-r2r-1d-multi
+ :in-count count :out-c out-c :plan plan
+ :in-cf in-cf :out-cf out-cf)))
+
+(defun fftw-r2r-1d-multi (multi in &key (start 0) (count nil) (normalize t))
+ (declare #.*standard-optimize-settings*)
+ (multiple-value-bind (start count) (normalize-input-range in start count)
+ (declare (fixnum start count))
+ (unless (equal (in-count multi) count)
+ (error "Different plan and vector lengths."))
+
+ (let ((in-cf (in-cf multi))
+ (out-cf (out-cf multi))
+ (out-c (out-c multi))
+ (in-pos start))
+ (declare (fixnum in-pos))
+
+ (dotimes (i count)
+ (declare (fixnum i))
+ (setf (cffi:mem-aref in-cf :double i) (coerce (aref in in-pos) 'double-float))
+ (incf in-pos))
+
+ (fftw-execute (plan multi))
+ (if normalize
+ (let ((factor (coerce (/ 1 count) 'double-float)))
+ (declare (double-float factor))
+ (dotimes (i count)
+ (declare (fixnum i))
+ (setf (aref out-c i) (* factor (cffi:mem-aref out-cf :double i)))))
+ (dotimes (i count)
+ (declare (fixnum i))
+ (setf (aref out-c i) (cffi:mem-aref out-cf :double i))))
+
+ out-c)))
+
+
+(defun fftw-hc2r-1d (in)
+ (declare #.*standard-optimize-settings*)
+ (let* ((count (length in))
+ (out-c (make-array count :element-type 'double-float))
+ (in-cf (fftw-malloc (* count +size-double+)))
+ (out-cf (fftw-malloc (* count +size-double+))))
+ (declare (fixnum count))
+
+ (dotimes (i count)
+ (declare (fixnum i))
+ (setf (cffi:mem-aref in-cf :double i) (aref in i)))
+
+ (let ((p (fftw-plan-r2r-1d count in-cf out-cf +fftw-hc2r+
+ (logior +fftw-estimate+ +fftw-destroy-input+))))
+ (fftw-execute p)
+ (fftw-destroy-plan p))
+
+ (dotimes (i count)
+ (declare (fixnum i))
+ (setf (aref out-c i) (cffi:mem-aref out-cf :double i)))
+ (fftw-free in-cf)
+ (fftw-free out-cf)
+ out-c))
+
+(defun fftw-r2c-1d (in &key (start 0) (count nil) (normalize t))
+ (declare #.*standard-optimize-settings*)
+ (multiple-value-bind (start count) (normalize-input-range in start count)
+ (declare (fixnum start count))
+
+ (let* ((out-n (1+ (floor count 2)))
+ (in-cf (fftw-malloc (* count +size-double+)))
+ (out-cf (fftw-malloc (* 2 out-n +size-double+)))
+ (in-pos start))
+ (declare (fixnum in-pos out-n))
+
+ (dotimes (i count)
+ (declare (fixnum i))
+ (setf (cffi:mem-aref in-cf :double i) (coerce (aref in in-pos) 'double-float))
+ (incf in-pos))
+
+ (let ((p (fftw-plan-dft-r2c-1d count in-cf out-cf
+ (logior +fftw-estimate+ +fftw-destroy-input+))))
+ (fftw-execute p)
+ (fftw-destroy-plan p))
+
+ (let ((out (make-array out-n :element-type '(or rational complex)))
+ (opos 0))
+ (declare (fixnum opos))
+ (if normalize
+ (let ((factor (coerce (/ 1 count) 'double-float)))
+ (declare (double-float factor))
+ (dotimes (i out-n)
+ (declare (fixnum i))
+ (let ((impos (1+ opos)))
+ (declare (fixnum impos))
+ (setf (aref out i) (complex (* (cffi:mem-aref out-cf :double opos) factor)
+ (* (cffi:mem-aref out-cf :double impos) factor)))
+ (incf opos 2))))
+ (dotimes (i out-n)
+ (declare (fixnum i))
+ (let ((impos (1+ opos)))
+ (declare (fixnum impos))
+ (setf (aref out i) (complex (cffi:mem-aref out-cf :double opos)
+ (cffi:mem-aref out-cf :double impos)))
+ (incf opos 2))))
+ (fftw-free in-cf)
+ (fftw-free out-cf)
+ out))))
+
+(defun make-fftw-r2c-1d-multi (count &key (complex-output t))
+ (declare #.*standard-optimize-settings*)
+ (declare (fixnum count))
+ (let* ((out-n (1+ (floor count 2)))
+ (in-cf (fftw-malloc (* count +size-double+)))
+ (out-cf (fftw-malloc (* 2 out-n +size-double+)))
+ (plan (fftw-plan-dft-r2c-1d count in-cf out-cf
+ (logior +fftw-estimate+ +fftw-destroy-input+))))
+ (declare (fixnum out-n))
+
+ (unless complex-output
+ (incf out-n out-n))
+
+ (make-instance 'fftw-r2c-1d-multi
+ :in-count count :out-count out-n :plan plan
+ :in-cf in-cf
+ :out-cf out-cf
+ :complex-output complex-output
+ :out
+ (if complex-output
+ (make-array out-n :element-type '(or double-float complex) :initial-element 0d0)
+ (make-array out-n :element-type 'double-float)))))
+
+(defgeneric destroy-fftw-multi (multi))
+(defmethod destroy-fftw-multi ((multi fftw-multi))
+ (declare #.*standard-optimize-settings*)
+ (when (plan multi)
+ (fftw-destroy-plan (plan multi))
+ (setf (plan multi) nil))
+ (dolist (slot '(in-cf out-cf))
+ (when (slot-value multi slot)
+ (fftw-free (slot-value multi slot))
+ (setf (slot-value multi slot) nil)))
+ (dolist (slot '(plan-even plan-odd))
+ (when (and (slot-exists-p multi slot)
+ (slot-value multi slot))
+ (fftw-destroy-plan (slot-value multi slot))
+ (setf (slot-value multi slot) nil))))
+
+
+(defun fftw-r2c-1d-multi (multi in &key (start 0) (normalize t))
+ (declare #.*standard-optimize-settings*)
+ (multiple-value-bind (start count) (normalize-input-range in start (in-count multi))
+ (declare (fixnum start count))
+
+ (unless (equal (length in) (in-count multi))
+ (error "In count of multi plan doesn't equal length of in vector."))
+
+ (let ((in-cf (in-cf multi))
+ (out-cf (out-cf multi))
+ (out (out multi))
+ (in-pos start)
+ (opos 0)
+ (in-count (in-count multi))
+ (out-count (out-count multi)))
+ (declare (fixnum in-pos opos in-count out-count))
+
+ (dotimes (i in-count)
+ (declare (fixnum i))
+ (setf (cffi:mem-aref in-cf :double i) (coerce (aref in in-pos) 'double-float))
+ (incf in-pos))
+
+ (fftw-execute (plan multi))
+
+ (setf opos 0)
+ (cond
+ ((complex-output multi)
+ (if normalize
+ (let ((factor (coerce (/ 1 count) 'double-float)))
+ (declare (double-float factor))
+ (dotimes (i out-count)
+ (declare (fixnum i))
+ (let ((impos (1+ opos)))
+ (declare (fixnum impos))
+ (setf (aref out i) (complex (* (cffi:mem-aref out-cf :double opos) factor)
+ (* (cffi:mem-aref out-cf :double impos) factor)))
+ (incf opos 2))))
+ (dotimes (i out-count)
+ (declare (fixnum i))
+ (let ((impos (1+ opos)))
+ (declare (fixnum impos))
+ (setf (aref out i) (complex (cffi:mem-aref out-cf :double opos)
+ (cffi:mem-aref out-cf :double impos)))
+ (incf opos 2)))))
+ (t
+ (when normalize
+ (let ((o2 (* out-count 2))
+ (factor (coerce (/ 1 count) 'double-float)))
+ (declare (fixnum o2)
+ (double-float factor))
+ (dotimes (i o2)
+ (declare (fixnum i))
+ (setf (aref out i) (* factor (aref out i))))))))
+ out)))
+
+(defun fftw-c2r-1d (in)
+ (declare #.*standard-optimize-settings*)
+ (let* ((count (length in))
+ (out-n (if (complex-is-even-p in)
+ (* 2 (1- count))
+ (1+ (* 2 (1- count)))))
+ (out-c (make-array out-n :element-type 'double-float))
+ (in-cf (fftw-malloc (* 2 count +size-double+)))
+ (out-cf (fftw-malloc (* 2 count +size-double+)))
+ (pos 0))
+ (declare (fixnum count out-n pos))
+
+ (setq pos 0)
+ (dotimes (i count)
+ (declare (fixnum i))
+ (let ((c (aref in i)))
+ (declare (complex c))
+ (setf (cffi:mem-aref in-cf :double pos) (coerce (realpart c) 'double-float))
+ (incf pos)
+ (setf (cffi:mem-aref in-cf :double pos) (coerce (imagpart c) 'double-float))
+ (incf pos)))
+
+ (let ((p (fftw-plan-dft-c2r-1d out-n in-cf out-cf
+ (logior +fftw-estimate+ +fftw-destroy-input+))))
+ (fftw-execute p)
+ (fftw-destroy-plan p))
+ (dotimes (i out-n)
+ (declare (fixnum i))
+ (setf (aref out-c i) (cffi:mem-aref out-cf :double i)))
+
+ (fftw-free in-cf)
+ (fftw-free out-cf)
+ out-c))
+
+(defun make-fftw-c2r-1d-multi (count)
+ (declare #.*standard-optimize-settings*)
+ (declare (fixnum count))
+ (let* ((out-count-even (* 2 (1- count)))
+ (out-count-odd (1+ (* 2 (1- count))))
+ (in-cf (fftw-malloc (* 2 count +size-double+)))
+ (out-cf (fftw-malloc (* 2 (max out-count-odd out-count-even) +size-double+))) ;; use larger size
+ (plan-odd (fftw-plan-dft-c2r-1d out-count-odd in-cf out-cf
+ (logior +fftw-estimate+ +fftw-destroy-input+)))
+ (plan-even (fftw-plan-dft-c2r-1d out-count-even in-cf out-cf +fftw-measure+)))
+ (declare (fixnum out-count-even out-count-odd))
+
+ (make-instance 'fftw-c2r-1d-multi
+ :in-count count
+ :plan-even plan-even
+ :plan-odd plan-odd
+ :in-cf in-cf
+ :out-cf out-cf
+ :out-count-even out-count-even
+ :out-count-odd out-count-odd
+ :out-even (make-array out-count-even :element-type 'double-float)
+ :out-odd (make-array out-count-odd :element-type 'double-float))))
+
+(defun fftw-c2r-1d-multi (multi in)
+ (declare #.*standard-optimize-settings*)
+ (let* ((is-even (complex-is-even-p in))
+ (out (if is-even (out-even multi) (out-odd multi)))
+ (out-count (if is-even (out-count-even multi) (out-count-odd multi)))
+ (plan (if is-even (plan-even multi) (plan-odd multi)))
+ (count (in-count multi))
+ (in-cf (in-cf multi))
+ (out-cf (out-cf multi))
+ (pos 0))
+ (declare (fixnum count out-count pos))
+
+ (setq pos 0)
+ (dotimes (i count)
+ (declare (fixnum i))
+ (let ((c (aref in i)))
+ (setf (cffi:mem-aref in-cf :double pos) (coerce (realpart c) 'double-float))
+ (incf pos)
+ (setf (cffi:mem-aref in-cf :double pos) (coerce (imagpart c) 'double-float))
+ (incf pos)))
+
+ (fftw-execute plan)
+
+ (dotimes (i out-count)
+ (declare (fixnum i))
+ (setf (aref out i) (cffi:mem-aref out-cf :double i)))
+
+ out))
+
+(defun fftw-c-1d (in &key (direction +fftw-forward+)
+ (start 0) (count nil) (normalize t))
+ (declare #.*standard-optimize-settings*)
+ (multiple-value-bind (start count) (normalize-input-range in start count)
+ (declare (fixnum start count))
+
+ (let* ((out-n count)
+ (out (make-array out-n :element-type '(or rational complex)))
+ (in-cf (fftw-malloc (* 2 count +size-double+)))
+ (out-cf (fftw-malloc (* 2 out-n +size-double+)))
+ (opos 0)
+ (ipos 0))
+ (declare (fixnum out-n ipos opos))
+
+ (setf opos 0)
+ (setf ipos start)
+ (dotimes (i count)
+ (declare (fixnum i))
+ (let ((c (aref in ipos)))
+ (incf ipos)
+ (setf (cffi:mem-aref in-cf :double opos) (coerce (realpart c) 'double-float))
+ (incf opos)
+ (setf (cffi:mem-aref in-cf :double opos) (coerce (imagpart c) 'double-float))
+ (incf opos)))
+
+ (let ((p (fftw-plan-dft-1d count in-cf out-cf direction
+ (logior +fftw-estimate+ +fftw-destroy-input+))))
+ (fftw-execute p)
+ (fftw-destroy-plan p))
+
+ (setf opos 0)
+ (if normalize
+ (let ((factor (coerce (/ 1 count) 'double-float)))
+ (declare (double-float factor))
+ (dotimes (i out-n)
+ (declare (fixnum i))
+ (let ((impos (1+ opos)))
+ (declare (fixnum impos))
+ (setf (aref out i) (complex (* factor (cffi:mem-aref out-cf :double opos))
+ (* factor (cffi:mem-aref out-cf :double impos)))))
+ (incf opos 2)))
+ (dotimes (i out-n)
+ (declare (fixnum i))
+ (let ((impos (1+ opos)))
+ (declare (fixnum impos))
+ (setf (aref out i) (complex (cffi:mem-aref out-cf :double opos)
+ (cffi:mem-aref out-cf :double impos))))
+ (incf opos 2)))
+
+ (fftw-free in-cf)
+ (fftw-free out-cf)
+ out)))
--- /dev/null
+;;;; -*- Mode: Lisp; Syntax: ANSI-Common-Lisp; Base: 10 -*-
+;;;; *************************************************************************
+;;;; FILE IDENTIFICATION
+;;;;
+;;;; Name: wisdowm.lisp
+;;;; Purpose: Functions for handling FFTW wisdom
+;;;; Programmer: Kevin M. Rosenberg
+;;;; Date Started: March 2009
+;;;;
+;;;; $Id$
+;;;;
+;;;; This file, part of FFTW3, is Copyright (c) 2009 by Kevin M. Rosenberg
+;;;;
+;;;; FFTW3 users are granted the rights to distribute and use this software
+;;;; as governed by the terms of the Lisp Lesser GNU Public License
+;;;; (http://opensource.franz.com/preamble.html), also known as the LLGPL.
+;;;; *************************************************************************
+
+(in-package #:fftw3)
+
+(defun import-user-wisdom ()
+ (let ((file-contents (kmrcl:read-file-to-string *user-wisdom-file*)))
+ (when (and (stringp file-contents) (plusp (length file-contents)))
+ (fftw-import-wisdom-from-string file-contents))))
+
+(defun export-user-wisdom ()
+ (let ((str+ptr (fftw-export-wisdom-to-string)))
+ (when (probe-file *user-wisdom-file*)
+ (delete-file *user-wisdom-file*))
+ (with-open-file (out *user-wisdom-file* :direction :output :if-exists :overwrite
+ :if-does-not-exist :create)
+ (format out "~A" (first str+ptr)))
+ (free (second str+ptr))))
+
+
projects / cl-fftw3.git / commitdiff