Skip to content

GitLab

Commit b9d6f6cd authored by Ivica Bukvic's avatar Ivica Bukvic
Browse files

improved fftw3 support

parent 5302598f
Hide whitespace changes
Inline Side-by-side
pd/src/configure.in
View file @ b9d6f6cd
......@@ -47,8 +47,10 @@ AC_ARG_ENABLE(setuid, [ --enable-setuid install as setuid (linux)],
setuid=$enableval)
AC_ARG_ENABLE(fftw, [ --enable-fftw use FFTW package],
fftw=$enableval)
AC_ARG_ENABLE(fat, [ --disable-fat build fat binary on Mac OS X],
AC_ARG_ENABLE(fat, [ --disable-fat build fat binary on Mac OS X],
fat=$enableval, fat="yes")
AC_ARG_ENABLE(tk, [ --disable-tk build without tcl/tk-GUI],
tk=$enableval)
dnl Checks for programs.
AC_PROG_CC
......@@ -96,8 +98,8 @@ AC_CHECK_LIB(pthread, pthread_create,PDLIB="$PDLIB -lpthread",
dnl Check for fftw package
if test x$fftw = "xyes";
then
AC_CHECK_LIB(fftw, fftw_one,PDLIB="$PDLIB -lfftw",
echo "fftw package not found - using built-in FFT"; fftw=no)
AC_CHECK_LIB(fftw3f, fftwf_version,PDLIB="$PDLIB -lfftw3f",
echo "fftw3 package not found - using built-in FFT"; fftw=no)
fi
dnl look for tcl 8.x... do I really have to go through all this!?
......@@ -176,13 +178,10 @@ else
GUISRC=
fi
# leave the executable as 'pd' on Mac OS X and Windows, but call it
# 'pdl2ork' on GNU/Linux
PDEXEC=pd
PDEXEC=pd-l2ork
if test `uname -s` = Linux;
if test `uname -s` = Linux -o `uname -s` = "GNU/kFreeBSD" -o `uname -s` = "GNU";
then
PDEXEC=pd-l2ork
dnl Ckecking for ALSA
echo .................... alsa= $alsa
dnl This should be fixed so Pd can use ALSA shared libraries where appropriate.
......@@ -206,7 +205,16 @@ dnl This should be fixed so Pd can use ALSA shared libraries where appropriate.
CPPFLAGS="-DDL_OPEN -DPA_USE_OSS -DUNIX -DUNISTD\
-DUSEAPI_OSS \
-fno-strict-aliasing"
SYSSRC="s_midi_oss.c s_audio_oss.c"
dnl No OSS on hurd
if test `uname -s` = "GNU";
then
SYSSRC="s_midi_dummy.c"
else
SYSSRC="s_midi_oss.c s_audio_oss.c"
CPPFLAGS=$CPPFLAGS" -DPA_USE_OSS -DUSEAPI_OSS"
fi
if test x$alsa = "xyes";
then
SYSSRC=$SYSSRC" s_audio_alsa.c s_audio_alsamm.c s_midi_alsa.c"
......@@ -223,7 +231,6 @@ dnl This should be fixed so Pd can use ALSA shared libraries where appropriate.
-I../portmidi/pm_common \
-I../portmidi/pm_linux"
SYSSRC="s_audio_pa.c \
s_audio_pablio.c \
s_audio_paring.c \
../portaudio/src/common/pa_allocation.c \
../portaudio/src/common/pa_converters.c \
......@@ -281,7 +288,7 @@ then
-framework AudioUnit -framework AudioToolbox \
-framework Carbon -framework CoreMIDI"
EXT=pd_darwin
CPPFLAGS="-DDL_OPEN -DMACOSX -DUNISTD -I/usr/X11R6/include \
CPPFLAGS="-DHAVE_LIBDL -DMACOSX -DHAVE_UNISTD_H -I/usr/X11R6/include \
-I../portaudio/include -I../portaudio/src/common \
-I../portaudio/src/os/mac_osx/ \
-I../portmidi/pm_common -I../portmidi/pm_mac \
......@@ -303,7 +310,6 @@ then
EXTERNTARGET=d_ppc
fi
SYSSRC="s_midi_pm.c s_audio_pa.c \
s_audio_pablio.c \
s_audio_paring.c \
../portaudio/src/common/pa_allocation.c \
../portaudio/src/common/pa_converters.c \
......@@ -323,8 +329,6 @@ then
../portaudio/src/hostapi/coreaudio/pa_mac_core_utilities.c \
../portmidi/pm_mac/pmmac.c \
../portmidi/pm_mac/pmmacosxcm.c \
../portmidi/pm_mac/finddefault.c \
../portmidi/pm_mac/readbinaryplist.c \
../portmidi/pm_common/pmutil.c \
../portmidi/pm_common/portmidi.c \
../portmidi/porttime/ptmacosx_cf.c "
......@@ -357,7 +361,7 @@ then
OSNUMBER=2
if test x$jack = "xyes";
then
LDFLAGS=$LDFLAGS" -weak_framework Jackmp"
LDFLAGS=$LDFLAGS" -weak_framework Jackmp"
fi
if test x$jack = "xrun";
then
......@@ -377,7 +381,7 @@ then
-mwindows -mms-bitfields "$MORECFLAGS
PDLIB=$PDLIB" -lwsock32 -lwinmm -lole32 -lstdc++"
SYSSRC="s_audio_pa.c s_audio_pablio.c s_audio_paring.c \
SYSSRC="s_audio_pa.c s_audio_paring.c \
s_audio_mmio.c s_midi_mmio.c \
../portaudio/src/common/pa_allocation.c \
../portaudio/src/common/pa_converters.c \
......@@ -418,7 +422,7 @@ fi
if test x$fftw = "xyes";
then
SYSSRC=$SYSSRC" d_fft_fftw.c d_fftroutine.c"
LDFLAGS=$LDFLAGS" -lfftw"
LDFLAGS=$LDFLAGS" -lfftw3"
else
SYSSRC=$SYSSRC" d_fft_mayer.c d_fftroutine.c"
fi
......
pd/src/d_fft_fftw.c
View file @ b9d6f6cd
......@@ -4,31 +4,32 @@
/* --------- Pd interface to FFTW library; imitate Mayer API ---------- */
#include "m_pd.h"
#ifdef MSW
#include <malloc.h>
#else
#include <alloca.h>
#endif
#error oops -- I'm talking to the old fftw. Apparently it's still changing.
/* changes and additions for FFTW3 by Thomas Grill */
#include <fftw.h>
#include "m_pd.h"
#include <fftw3.h>
int ilog2(int n);
#define MINFFT 5
#define MINFFT 0
#define MAXFFT 30
/* from the FFTW website:
fftw_complex in[N], out[N];
#include <fftw3.h>
...
{
fftw_complex *in, *out;
fftw_plan p;
...
p = fftw_create_plan(N, FFTW_FORWARD, FFTW_ESTIMATE);
in = (fftw_complex*) fftw_malloc(sizeof(fftw_complex) * N);
out = (fftw_complex*) fftw_malloc(sizeof(fftw_complex) * N);
p = fftw_plan_dft_1d(N, in, out, FFTW_FORWARD, FFTW_ESTIMATE);
...
fftw_one(p, in, out);
fftw_execute(p);
...
fftw_destroy_plan(p);
fftw_free(in); fftw_free(out);
}
FFTW_FORWARD or FFTW_BACKWARD, and indicates the direction of the transform you
are interested in. Alternatively, you can use the sign of the exponent in the
......@@ -37,57 +38,126 @@ respectively. The flags argument is either FFTW_MEASURE
*/
static fftw_plan fftw_pvec[2 * (MAXFFT+1 - MINFFT)];
/* complex stuff */
typedef struct {
fftwf_plan plan;
fftwf_complex *in,*out;
} cfftw_info;
static cfftw_info cfftw_fwd[MAXFFT+1 - MINFFT],cfftw_bwd[MAXFFT+1 - MINFFT];
static cfftw_info *cfftw_getplan(int n,int fwd)
{
cfftw_info *info;
int logn = ilog2(n);
if (logn < MINFFT || logn > MAXFFT)
return (0);
info = (fwd?cfftw_fwd:cfftw_bwd)+(logn-MINFFT);
if (!info->plan)
{
info->in = (fftwf_complex*) fftwf_malloc(sizeof(fftwf_complex) * n);
info->out = (fftwf_complex*) fftwf_malloc(sizeof(fftwf_complex) * n);
info->plan = fftwf_plan_dft_1d(n, info->in, info->out, fwd?FFTW_FORWARD:FFTW_BACKWARD, FFTW_MEASURE);
}
return info;
}
static fftw_plan fftw_getplan(int n, int dir)
/* real stuff */
typedef struct {
fftwf_plan plan;
float *in,*out;
} rfftw_info;
static rfftw_info rfftw_fwd[MAXFFT+1 - MINFFT],rfftw_bwd[MAXFFT+1 - MINFFT];
static rfftw_info *rfftw_getplan(int n,int fwd)
{
logn = ilog2(n);
rfftw_info *info;
int logn = ilog2(n);
if (logn < MINFFT || logn > MAXFFT)
return (0);
int indx = 2*(logn-MINFFT) + inverse);
if (!fftw_pvec[indx]
fftw_pvec[indx] = fftw_create_plan(N, dir, FFTW_MEASURE);
return (fftw_pvec[indx]);
info = (fwd?rfftw_fwd:rfftw_bwd)+(logn-MINFFT);
if (!info->plan)
{
info->in = (float*) fftwf_malloc(sizeof(float) * n);
info->out = (float*) fftwf_malloc(sizeof(float) * n);
info->plan = fftwf_plan_r2r_1d(n, info->in, info->out, fwd?FFTW_R2HC:FFTW_HC2R, FFTW_MEASURE);
}
return info;
}
EXTERN void mayer_fht(float *fz, int n)
{
post("FHT: not yet implemented");
}
static void mayer_dofft(int n, float *fz1, float *fz2, int dir)
static void mayer_do_cfft(int n, float *fz1, float *fz2, int fwd)
{
float *inbuf, *outbuf, *fp1, *fp2, *fp3;
int i;
fftw_plan p = fftw_getplan(n, dir);
inbuf = alloca(n * (4 * sizeof(float)));
outbuf = inbuf + 2*n;
float *fz;
cfftw_info *p = cfftw_getplan(n, fwd);
if (!p)
return;
for (i = 0, fp1 = fz1, fp2 = fz2, fp3 = inbuf; i < n; i++)
fp3[0] = *fp1++, fp3[1] = *fp2++, fp3 += 2;
fftw_one(p, inbuf, outbuf);
for (i = 0, fp1 = fz1, fp2 = fz2, fp3 = outbuf; i < n; i++)
*fp1++ = fp3[0], *fp2++ = fp3[1], fp3 += 2;
for (i = 0, fz = (float *)p->in; i < n; i++)
fz[i*2] = fz1[i], fz[i*2+1] = fz2[i];
fftwf_execute(p->plan);
for (i = 0, fz = (float *)p->out; i < n; i++)
fz1[i] = fz[i*2], fz2[i] = fz[i*2+1];
}
EXTERN void mayer_fft(int n, float *fz1, float *fz2)
{
mayer_dofft(n, fz1, fz2, FFTW_FORWARD);
mayer_do_cfft(n, fz1, fz2, 1);
}
EXTERN void mayer_ifft(int n, float *fz1, float *fz2)
{
mayer_dofft(n, fz1, fz2, FFTW_BACKWARD);
mayer_do_cfft(n, fz1, fz2, 0);
}
/*
in the following the sign flips are done to
be compatible with the mayer_fft implementation,
but it's probably the mayer_fft that should be corrected...
*/
EXTERN void mayer_realfft(int n, float *fz)
{
post("rfft: not yet implemented");
int i;
rfftw_info *p = rfftw_getplan(n, 1);
if (!p)
return;
for (i = 0; i < n; i++)
p->in[i] = fz[i];
fftwf_execute(p->plan);
for (i = 0; i < n/2+1; i++)
fz[i] = p->out[i];
for (; i < n; i++)
fz[i] = -p->out[i];
}
EXTERN void mayer_realifft(int n, float *fz)
{
post("rifft: not yet implemented");
int i;
rfftw_info *p = rfftw_getplan(n, 0);
if (!p)
return;
for (i = 0; i < n/2+1; i++)
p->in[i] = fz[i];
for (; i < n; i++)
p->in[i] = -fz[i];
fftwf_execute(p->plan);
for (i = 0; i < n; i++)
fz[i] = p->out[i];
}
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