diff options
-rw-r--r-- | include/ffts.h | 5 | ||||
-rw-r--r-- | src/Makefile.am | 2 | ||||
-rw-r--r-- | src/Makefile.in | 11 | ||||
-rw-r--r-- | src/ffts.h | 3 | ||||
-rw-r--r-- | src/ffts_real.c | 117 | ||||
-rw-r--r-- | src/ffts_real.h | 51 |
6 files changed, 181 insertions, 8 deletions
diff --git a/include/ffts.h b/include/ffts.h index f2194dc..1888982 100644 --- a/include/ffts.h +++ b/include/ffts.h @@ -50,7 +50,10 @@ ffts_plan_t *ffts_init_1d(size_t N, int sign); ffts_plan_t *ffts_init_2d(size_t N1, size_t N2, int sign); ffts_plan_t *ffts_init_nd(int rank, size_t *Ns, int sign); -void ffts_execute(ffts_plan_t * , const void * , const void * ); +ffts_plan_t *ffts_init_1d_real(size_t N, int sign); + + +void ffts_execute(ffts_plan_t * , const void * , void * ); void ffts_free(ffts_plan_t *); #ifdef __cplusplus diff --git a/src/Makefile.am b/src/Makefile.am index 80c4dcb..ce73303 100644 --- a/src/Makefile.am +++ b/src/Makefile.am @@ -2,7 +2,7 @@ lib_LTLIBRARIES = libffts.la -libffts_la_SOURCES = ffts.c ffts_nd.c patterns.c codegen.c +libffts_la_SOURCES = ffts.c ffts_nd.c ffts_real.c patterns.c codegen.c libffts_includedir=$(includedir)/ffts libffts_include_HEADERS = ../include/ffts.h diff --git a/src/Makefile.in b/src/Makefile.in index ee42f18..0c3789b 100644 --- a/src/Makefile.in +++ b/src/Makefile.in @@ -95,12 +95,12 @@ am__installdirs = "$(DESTDIR)$(libdir)" \ "$(DESTDIR)$(libffts_includedir)" LTLIBRARIES = $(lib_LTLIBRARIES) libffts_la_LIBADD = -am__libffts_la_SOURCES_DIST = ffts.c ffts_nd.c patterns.c codegen.c \ - neon.s sse.s +am__libffts_la_SOURCES_DIST = ffts.c ffts_nd.c ffts_real.c patterns.c \ + codegen.c neon.s sse.s @HAVE_NEON_TRUE@am__objects_1 = neon.lo @HAVE_NEON_FALSE@@HAVE_SSE_TRUE@am__objects_2 = sse.lo -am_libffts_la_OBJECTS = ffts.lo ffts_nd.lo patterns.lo codegen.lo \ - $(am__objects_1) $(am__objects_2) +am_libffts_la_OBJECTS = ffts.lo ffts_nd.lo ffts_real.lo patterns.lo \ + codegen.lo $(am__objects_1) $(am__objects_2) libffts_la_OBJECTS = $(am_libffts_la_OBJECTS) DEFAULT_INCLUDES = -I.@am__isrc@ -I$(top_builddir) depcomp = $(SHELL) $(top_srcdir)/depcomp @@ -251,7 +251,7 @@ top_build_prefix = @top_build_prefix@ top_builddir = @top_builddir@ top_srcdir = @top_srcdir@ lib_LTLIBRARIES = libffts.la -libffts_la_SOURCES = ffts.c ffts_nd.c patterns.c codegen.c \ +libffts_la_SOURCES = ffts.c ffts_nd.c ffts_real.c patterns.c codegen.c \ $(am__append_1) $(am__append_2) libffts_includedir = $(includedir)/ffts libffts_include_HEADERS = ../include/ffts.h @@ -335,6 +335,7 @@ distclean-compile: @AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/codegen.Plo@am__quote@ @AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/ffts.Plo@am__quote@ @AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/ffts_nd.Plo@am__quote@ +@AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/ffts_real.Plo@am__quote@ @AMDEP_TRUE@@am__include@ @am__quote@./$(DEPDIR)/patterns.Plo@am__quote@ .c.o: @@ -83,7 +83,8 @@ struct _ffts_plan_t { void *transpose_buf; void (*destroy)(struct _ffts_plan_t *); - + + float *A, *B; }; typedef struct _ffts_plan_t ffts_plan_t; diff --git a/src/ffts_real.c b/src/ffts_real.c new file mode 100644 index 0000000..5f3bb02 --- /dev/null +++ b/src/ffts_real.c @@ -0,0 +1,117 @@ +/* + + This file is part of FFTS -- The Fastest Fourier Transform in the South + + Copyright (c) 2012, Anthony M. Blake <amb@anthonix.com> + Copyright (c) 2012, The University of Waikato + + All rights reserved. + + Redistribution and use in source and binary forms, with or without + modification, are permitted provided that the following conditions are met: + * Redistributions of source code must retain the above copyright + notice, this list of conditions and the following disclaimer. + * Redistributions in binary form must reproduce the above copyright + notice, this list of conditions and the following disclaimer in the + documentation and/or other materials provided with the distribution. + * Neither the name of the organization nor the + names of its contributors may be used to endorse or promote products + derived from this software without specific prior written permission. + + THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND + ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED + WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE + DISCLAIMED. IN NO EVENT SHALL ANTHONY M. BLAKE BE LIABLE FOR ANY + DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES + (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; + LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND + ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT + (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS + SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. + +*/ + +#include "ffts_real.h" + +void ffts_free_1d_real(ffts_plan_t *p) { + ffts_free(p->plans[0]); + free(p->A); + free(p->B); + free(p->plans); + free(p->buf); + free(p); +} + +void ffts_execute_1d_real(ffts_plan_t *p, const void *vin, void *vout) { + float *out = (float *)vout; + float *buf = (float *)p->buf; + float *A = p->A; + float *B = p->B; + + p->plans[0]->transform(p->plans[0], vin, buf); + + size_t N = p->N; + buf[N] = buf[0]; + buf[N+1] = buf[1]; + + size_t i; + for(i=0;i<N/2;i++) { + out[2*i] = buf[2*i]*A[2*i] - buf[2*i+1]*A[2*i+1] + buf[N-2*i]*B[2*i] + buf[N-2*i+1]*B[2*i+1]; + out[2*i+1] = buf[2*i+1]*A[2*i] + buf[2*i]*A[2*i+1] + buf[N-2*i]*B[2*i+1] - buf[N-2*i+1]*B[2*i]; + out[2*N - 2*i] = out[2*i]; + out[2*N - 2*i + 1] = -out[2*i+1]; + } + + out[N] = buf[0] - buf[1]; + out[N+1] = 0.0f; + +} + +void ffts_execute_1d_real_inv(ffts_plan_t *p, const void *vin, void *vout) { + float *out = (float *)vout; + float *in = (float *)vin; + float *buf = (float *)p->buf; + float *A = p->A; + float *B = p->B; + size_t N = p->N; + + size_t i; + for(i=0;i<N/2;i++) { + buf[2*i] = in[2*i]*A[2*i] + in[2*i+1]*A[2*i+1] + in[N-2*i]*B[2*i] - in[N-2*i+1]*B[2*i+1]; + buf[2*i+1] = in[2*i+1]*A[2*i] - in[2*i]*A[2*i+1] - in[N-2*i]*B[2*i+1] - in[N-2*i+1]*B[2*i]; + } + + p->plans[0]->transform(p->plans[0], buf, out); + +} + +ffts_plan_t *ffts_init_1d_real(size_t N, int sign) { + ffts_plan_t *p = malloc(sizeof(ffts_plan_t)); + + if(sign < 0) p->transform = &ffts_execute_1d_real; + else p->transform = &ffts_execute_1d_real_inv; + + p->destroy = &ffts_free_1d_real; + p->N = N; + p->rank = 1; + p->plans = malloc(sizeof(ffts_plan_t **) * 1); + + p->plans[0] = ffts_init_1d(N/2, sign); + + p->buf = valloc(sizeof(float) * 2 * (N/2) + 1); + + p->A = valloc(sizeof(float) * N); + p->B = valloc(sizeof(float) * N); + + int i; + for (i = 0; i < N/2; i++) { + p->A[2 * i] = 0.5 * (1.0 - sin (2 * PI / (double) N * (double) i)); + p->A[2 * i + 1] = 0.5 * (-1.0 * cos (2 * PI / (double) N * (double) i)); + p->B[2 * i] = 0.5 * (1.0 + sin (2 * PI / (double) N * (double) i)); + p->B[2 * i + 1] = 0.5 * (1.0 * cos (2 * PI / (double) N * (double) i)); + } + + return p; +} + + diff --git a/src/ffts_real.h b/src/ffts_real.h new file mode 100644 index 0000000..06e15bb --- /dev/null +++ b/src/ffts_real.h @@ -0,0 +1,51 @@ +/* + + This file is part of FFTS -- The Fastest Fourier Transform in the South + + Copyright (c) 2012, Anthony M. Blake <amb@anthonix.com> + Copyright (c) 2012, The University of Waikato + + All rights reserved. + + Redistribution and use in source and binary forms, with or without + modification, are permitted provided that the following conditions are met: + * Redistributions of source code must retain the above copyright + notice, this list of conditions and the following disclaimer. + * Redistributions in binary form must reproduce the above copyright + notice, this list of conditions and the following disclaimer in the + documentation and/or other materials provided with the distribution. + * Neither the name of the organization nor the + names of its contributors may be used to endorse or promote products + derived from this software without specific prior written permission. + + THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" AND + ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED + WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE + DISCLAIMED. IN NO EVENT SHALL ANTHONY M. BLAKE BE LIABLE FOR ANY + DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES + (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; + LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND + ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT + (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS + SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. + +*/ + +#ifndef __FFTS_REAL_H__ +#define __FFTS_REAL_H__ + +#include <stdint.h> +#include <stddef.h> +#include <stdio.h> + +#include "ffts.h" + +#ifdef __ARM_NEON__ + #include <arm_neon.h> +#else + #include <xmmintrin.h> +#endif + + +#endif + |