File indexing completed on 2024-05-19 04:54:31
0001 /* 0002 SPDX-FileCopyrightText: 2003-2010 Mark Borgerding <Mark@Borgerding.net> 0003 SPDX-License-Identifier: BSD-3-Clause 0004 */ 0005 0006 #pragma once 0007 0008 /* kiss_fft.h 0009 defines kiss_fft_scalar as either short or a float type 0010 and defines 0011 typedef struct { kiss_fft_scalar r; kiss_fft_scalar i; }kiss_fft_cpx; */ 0012 #include "kiss_fft.h" 0013 #include <limits.h> 0014 0015 #define MAXFACTORS 32 0016 /* e.g. an fft of length 128 has 4 factors 0017 as far as kissfft is concerned 0018 4*4*4*2 0019 */ 0020 0021 struct kiss_fft_state 0022 { 0023 int nfft; 0024 int inverse; 0025 int factors[2 * MAXFACTORS]; 0026 kiss_fft_cpx twiddles[1]; 0027 }; 0028 0029 /* 0030 Explanation of macros dealing with complex math: 0031 0032 C_MUL(m,a,b) : m = a*b 0033 C_FIXDIV( c , div ) : if a fixed point impl., c /= div. noop otherwise 0034 C_SUB( res, a,b) : res = a - b 0035 C_SUBFROM( res , a) : res -= a 0036 C_ADDTO( res , a) : res += a 0037 * */ 0038 #ifdef FIXED_POINT 0039 #if (FIXED_POINT == 32) 0040 #define FRACBITS 31 0041 #define SAMPPROD int64_t 0042 #define SAMP_MAX 2147483647 0043 #else 0044 #define FRACBITS 15 0045 #define SAMPPROD int32_t 0046 #define SAMP_MAX 32767 0047 #endif 0048 0049 #define SAMP_MIN -SAMP_MAX 0050 0051 #if defined(CHECK_OVERFLOW) 0052 #define CHECK_OVERFLOW_OP(a, op, b) \ 0053 if ((SAMPPROD)(a)op(SAMPPROD)(b) > SAMP_MAX || (SAMPPROD)(a)op(SAMPPROD)(b) < SAMP_MIN) { \ 0054 fprintf(stderr, "WARNING:overflow @ " __FILE__ "(%d): (%d " #op " %d) = %ld\n", __LINE__, (a), (b), (SAMPPROD)(a)op(SAMPPROD)(b)); \ 0055 } 0056 #endif 0057 0058 #define smul(a, b) ((SAMPPROD)(a) * (b)) 0059 #define sround(x) (kiss_fft_scalar)(((x) + (1 << (FRACBITS - 1))) >> FRACBITS) 0060 0061 #define S_MUL(a, b) sround(smul(a, b)) 0062 0063 #define C_MUL(m, a, b) \ 0064 do { \ 0065 (m).r = sround(smul((a).r, (b).r) - smul((a).i, (b).i)); \ 0066 (m).i = sround(smul((a).r, (b).i) + smul((a).i, (b).r)); \ 0067 } while (0) 0068 0069 #define DIVSCALAR(x, k) (x) = sround(smul(x, SAMP_MAX / k)) 0070 0071 #define C_FIXDIV(c, div) \ 0072 do { \ 0073 DIVSCALAR((c).r, div); \ 0074 DIVSCALAR((c).i, div); \ 0075 } while (0) 0076 0077 #define C_MULBYSCALAR(c, s) \ 0078 do { \ 0079 (c).r = sround(smul((c).r, s)); \ 0080 (c).i = sround(smul((c).i, s)); \ 0081 } while (0) 0082 0083 #else /* not FIXED_POINT*/ 0084 0085 #define S_MUL(a, b) ((a) * (b)) 0086 #define C_MUL(m, a, b) \ 0087 do { \ 0088 (m).r = (a).r * (b).r - (a).i * (b).i; \ 0089 (m).i = (a).r * (b).i + (a).i * (b).r; \ 0090 } while (0) 0091 #define C_FIXDIV(c, div) /* NOOP */ 0092 #define C_MULBYSCALAR(c, s) \ 0093 do { \ 0094 (c).r *= (s); \ 0095 (c).i *= (s); \ 0096 } while (0) 0097 #endif 0098 0099 #ifndef CHECK_OVERFLOW_OP 0100 #define CHECK_OVERFLOW_OP(a, op, b) /* noop */ 0101 #endif 0102 0103 #define C_ADD(res, a, b) \ 0104 do { \ 0105 CHECK_OVERFLOW_OP((a).r, +, (b).r) \ 0106 CHECK_OVERFLOW_OP((a).i, +, (b).i) \ 0107 (res).r = (a).r + (b).r; \ 0108 (res).i = (a).i + (b).i; \ 0109 } while (0) 0110 #define C_SUB(res, a, b) \ 0111 do { \ 0112 CHECK_OVERFLOW_OP((a).r, -, (b).r) \ 0113 CHECK_OVERFLOW_OP((a).i, -, (b).i) \ 0114 (res).r = (a).r - (b).r; \ 0115 (res).i = (a).i - (b).i; \ 0116 } while (0) 0117 #define C_ADDTO(res, a) \ 0118 do { \ 0119 CHECK_OVERFLOW_OP((res).r, +, (a).r) \ 0120 CHECK_OVERFLOW_OP((res).i, +, (a).i) \ 0121 (res).r += (a).r; \ 0122 (res).i += (a).i; \ 0123 } while (0) 0124 0125 #define C_SUBFROM(res, a) \ 0126 do { \ 0127 CHECK_OVERFLOW_OP((res).r, -, (a).r) \ 0128 CHECK_OVERFLOW_OP((res).i, -, (a).i) \ 0129 (res).r -= (a).r; \ 0130 (res).i -= (a).i; \ 0131 } while (0) 0132 0133 #ifdef FIXED_POINT 0134 #define KISS_FFT_COS(phase) floor(.5 + SAMP_MAX * cos(phase)) 0135 #define KISS_FFT_SIN(phase) floor(.5 + SAMP_MAX * sin(phase)) 0136 #define HALF_OF(x) ((x) >> 1) 0137 #elif defined(USE_SIMD) 0138 #define KISS_FFT_COS(phase) _mm_set1_ps(cos(phase)) 0139 #define KISS_FFT_SIN(phase) _mm_set1_ps(sin(phase)) 0140 #define HALF_OF(x) ((x)*_mm_set1_ps(.5)) 0141 #else 0142 #define KISS_FFT_COS(phase) (kiss_fft_scalar) cos(phase) 0143 #define KISS_FFT_SIN(phase) (kiss_fft_scalar) sin(phase) 0144 #define HALF_OF(x) ((x)*.5) 0145 #endif 0146 0147 #define kf_cexp(x, phase) \ 0148 do { \ 0149 (x)->r = KISS_FFT_COS(phase); \ 0150 (x)->i = KISS_FFT_SIN(phase); \ 0151 } while (0) 0152 0153 /* a debugging function */ 0154 #define pcpx(c) fprintf(stderr, "%g + %gi\n", (double)((c)->r), (double)((c)->i)) 0155 0156 #ifdef KISS_FFT_USE_ALLOCA 0157 // define this to allow use of alloca instead of malloc for temporary buffers 0158 // Temporary buffers are used in two case: 0159 // 1. FFT sizes that have "bad" factors. i.e. not 2,3 and 5 0160 // 2. "in-place" FFTs. Notice the quotes, since kissfft does not really do an in-place transform. 0161 #include <alloca.h> 0162 #define KISS_FFT_TMP_ALLOC(nbytes) alloca(nbytes) 0163 #define KISS_FFT_TMP_FREE(ptr) 0164 #else 0165 #define KISS_FFT_TMP_ALLOC(nbytes) KISS_FFT_MALLOC(nbytes) 0166 #define KISS_FFT_TMP_FREE(ptr) KISS_FFT_FREE(ptr) 0167 #endif