1 /* 2 * Copyright (c) 2003-2010, Mark Borgerding. All rights reserved. 3 * This file is part of KISS FFT - https://github.com/mborgerding/kissfft 4 * 5 * SPDX-License-Identifier: BSD-3-Clause 6 * See COPYING file for more information. 7 */ 8 9 /* kiss_fft_f32.h 10 defines kiss_fft_f32_scalar as either short or a float type 11 and defines 12 typedef struct { kiss_fft_f32_scalar r; kiss_fft_f32_scalar i; }kiss_fft_f32_cpx; */ 13 #include "kiss_fft_f32.h" 14 #include <limits.h> 15 16 /* The 2*sizeof(size_t) alignment here is borrowed from 17 * GNU libc, so it should be good most everywhere. 18 * It is more conservative than is needed on some 64-bit 19 * platforms, but ia64 does require a 16-byte alignment. 20 * The SIMD extensions for x86 and ppc32 would want a 21 * larger alignment than this, but we don't need to 22 * do better than malloc. 23 * 24 * Borrowed from GLib's gobject/gtype.c 25 */ 26 #define STRUCT_ALIGNMENT (2 * sizeof (size_t)) 27 #define ALIGN_STRUCT(offset) \ 28 ((offset + (STRUCT_ALIGNMENT - 1)) & -STRUCT_ALIGNMENT) 29 30 #define MAXFACTORS 32 31 /* e.g. an fft of length 128 has 4 factors 32 as far as kissfft is concerned 33 4*4*4*2 34 */ 35 36 struct kiss_fft_f32_state{ 37 int nfft; 38 int inverse; 39 int factors[2*MAXFACTORS]; 40 kiss_fft_f32_cpx twiddles[1]; 41 }; 42 43 /* 44 Explanation of macros dealing with complex math: 45 46 C_MUL(m,a,b) : m = a*b 47 C_FIXDIV( c , div ) : if a fixed point impl., c /= div. noop otherwise 48 C_SUB( res, a,b) : res = a - b 49 C_SUBFROM( res , a) : res -= a 50 C_ADDTO( res , a) : res += a 51 * */ 52 #ifdef FIXED_POINT 53 #include <stdint.h> 54 #if (FIXED_POINT==32) 55 # define FRACBITS 31 56 # define SAMPPROD int64_t 57 #define SAMP_MAX INT32_MAX 58 #define SAMP_MIN INT32_MIN 59 #else 60 # define FRACBITS 15 61 # define SAMPPROD int32_t 62 #define SAMP_MAX INT16_MAX 63 #define SAMP_MIN INT16_MIN 64 #endif 65 66 #if defined(CHECK_OVERFLOW) 67 # define CHECK_OVERFLOW_OP(a,op,b) \ 68 if ( (SAMPPROD)(a) op (SAMPPROD)(b) > SAMP_MAX || (SAMPPROD)(a) op (SAMPPROD)(b) < SAMP_MIN ) { \ 69 g_critical("overflow @ " __FILE__ "(%d): (%d " #op" %d) = %ld",__LINE__,(a),(b),(SAMPPROD)(a) op (SAMPPROD)(b) ); } 70 #endif 71 72 73 # define smul(a,b) ( (SAMPPROD)(a)*(b) ) 74 # define sround( x ) (kiss_fft_f32_scalar)( ( (x) + (1<<(FRACBITS-1)) ) >> FRACBITS ) 75 76 # define S_MUL(a,b) sround( smul(a,b) ) 77 78 # define C_MUL(m,a,b) \ 79 do{ (m).r = sround( smul((a).r,(b).r) - smul((a).i,(b).i) ); \ 80 (m).i = sround( smul((a).r,(b).i) + smul((a).i,(b).r) ); }while(0) 81 82 # define DIVSCALAR(x,k) \ 83 (x) = sround( smul( x, SAMP_MAX/k ) ) 84 85 # define C_FIXDIV(c,div) \ 86 do { DIVSCALAR( (c).r , div); \ 87 DIVSCALAR( (c).i , div); }while (0) 88 89 # define C_MULBYSCALAR( c, s ) \ 90 do{ (c).r = sround( smul( (c).r , s ) ) ;\ 91 (c).i = sround( smul( (c).i , s ) ) ; }while(0) 92 93 #else /* not FIXED_POINT*/ 94 95 # define S_MUL(a,b) ( (a)*(b) ) 96 #define C_MUL(m,a,b) \ 97 do{ (m).r = (a).r*(b).r - (a).i*(b).i;\ 98 (m).i = (a).r*(b).i + (a).i*(b).r; }while(0) 99 # define C_FIXDIV(c,div) /* NOOP */ 100 # define C_MULBYSCALAR( c, s ) \ 101 do{ (c).r *= (s);\ 102 (c).i *= (s); }while(0) 103 #endif 104 105 #ifndef CHECK_OVERFLOW_OP 106 # define CHECK_OVERFLOW_OP(a,op,b) /* noop */ 107 #endif 108 109 #define C_ADD( res, a,b)\ 110 do { \ 111 CHECK_OVERFLOW_OP((a).r,+,(b).r)\ 112 CHECK_OVERFLOW_OP((a).i,+,(b).i)\ 113 (res).r=(a).r+(b).r; (res).i=(a).i+(b).i; \ 114 }while(0) 115 #define C_SUB( res, a,b)\ 116 do { \ 117 CHECK_OVERFLOW_OP((a).r,-,(b).r)\ 118 CHECK_OVERFLOW_OP((a).i,-,(b).i)\ 119 (res).r=(a).r-(b).r; (res).i=(a).i-(b).i; \ 120 }while(0) 121 #define C_ADDTO( res , a)\ 122 do { \ 123 CHECK_OVERFLOW_OP((res).r,+,(a).r)\ 124 CHECK_OVERFLOW_OP((res).i,+,(a).i)\ 125 (res).r += (a).r; (res).i += (a).i;\ 126 }while(0) 127 128 #define C_SUBFROM( res , a)\ 129 do {\ 130 CHECK_OVERFLOW_OP((res).r,-,(a).r)\ 131 CHECK_OVERFLOW_OP((res).i,-,(a).i)\ 132 (res).r -= (a).r; (res).i -= (a).i; \ 133 }while(0) 134 135 136 #ifdef FIXED_POINT 137 # define KISS_FFT_F32_COS(phase) floor(.5+SAMP_MAX * cos (phase)) 138 # define KISS_FFT_F32_SIN(phase) floor(.5+SAMP_MAX * sin (phase)) 139 # define HALF_OF(x) ((x)>>1) 140 #elif defined(USE_SIMD) 141 # define KISS_FFT_F32_COS(phase) _mm_set1_ps( cos(phase) ) 142 # define KISS_FFT_F32_SIN(phase) _mm_set1_ps( sin(phase) ) 143 # define HALF_OF(x) ((x)*_mm_set1_ps(.5)) 144 #else 145 # define KISS_FFT_F32_COS(phase) (kiss_fft_f32_scalar) cos(phase) 146 # define KISS_FFT_F32_SIN(phase) (kiss_fft_f32_scalar) sin(phase) 147 # define HALF_OF(x) ((x)*.5) 148 #endif 149 150 #define kf_cexp(x,phase) \ 151 do{ \ 152 (x)->r = KISS_FFT_F32_COS(phase);\ 153 (x)->i = KISS_FFT_F32_SIN(phase);\ 154 }while(0) 155 156 157 /* a debugging function */ 158 #define pcpx(c)\ 159 fprintf(stderr,"%g + %gi\n",(double)((c)->r),(double)((c)->i) ) 160 161 162 #ifdef KISS_FFT_F32_USE_ALLOCA 163 // define this to allow use of alloca instead of malloc for temporary buffers 164 // Temporary buffers are used in two case: 165 // 1. FFT sizes that have "bad" factors. i.e. not 2,3 and 5 166 // 2. "in-place" FFTs. Notice the quotes, since kissfft does not really do an in-place transform. 167 #include <alloca.h> 168 #define KISS_FFT_F32_TMP_ALLOC(nbytes) alloca(nbytes) 169 #define KISS_FFT_F32_TMP_FREE(ptr) 170 #else 171 #define KISS_FFT_F32_TMP_ALLOC(nbytes) KISS_FFT_F32_MALLOC(nbytes) 172 #define KISS_FFT_F32_TMP_FREE(ptr) KISS_FFT_F32_FREE(ptr) 173 #endif 174