Mercurial > hg > peckfft

comparison _peck_fft_guts.h @ 6:fee54f1878f7

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
kill FIXED_POINT stuff, simplify
author Peter Meerwald <p.meerwald@bct-electronic.com>
date Fri, 16 Sep 2011 15:18:28 +0200
parents 2d6c49fcafcb
children 655dc5c14169
comparison
equal deleted inserted replaced
5:c7237a7544eb 6:fee54f1878f7
31 int factors[2*MAXFACTORS]; 31 int factors[2*MAXFACTORS];
32 peck_fft_cpx twiddles[1]; 32 peck_fft_cpx twiddles[1];
33 }; 33 };
34 34
35 /* 35 /*
36 Explanation of macros dealing with complex math: 36 * Explanation of macros dealing with complex math:
37 * C_MUL(m,a,b) : m = a*b
38 * C_SUB( res, a,b) : res = a - b
39 * C_SUBFROM( res , a) : res -= a
40 * C_ADDTO( res , a) : res += a
41 */
37 42
38 C_MUL(m,a,b) : m = a*b 43 #define S_MUL(a, b) ((a) * (b))
39 C_FIXDIV( c , div ) : if a fixed point impl., c /= div. noop otherwise 44 #define C_MUL(m, a, b) \
40 C_SUB( res, a,b) : res = a - b 45 do { \
41 C_SUBFROM( res , a) : res -= a 46 (m).r = (a).r*(b).r - (a).i*(b).i; \
42 C_ADDTO( res , a) : res += a 47 (m).i = (a).r*(b).i + (a).i*(b).r; \
43 * */ 48 } while(0)
44 #ifdef FIXED_POINT 49 #define C_MULBYSCALAR(c, s ) \
45 #if (FIXED_POINT==32) 50 do { \
46 # define FRACBITS 31 51 (c).r *= (s); \
47 # define SAMPPROD int64_t 52 (c).i *= (s); \
48 #define SAMP_MAX 2147483647 53 } while(0)
49 #else 54 #define C_ADD(res, a, b) \
50 # define FRACBITS 15 55 do { \
51 # define SAMPPROD int32_t 56 (res).r = (a).r + (b).r; \
52 #define SAMP_MAX 32767 57 (res).i = (a).i + (b).i; \
53 #endif 58 } while(0)
59 #define C_SUB(res, a, b) \
60 do { \
61 (res).r = (a).r - (b).r; \
62 (res).i = (a).i - (b).i; \
63 } while(0)
64 #define C_ADDTO(res, a) \
65 do { \
66 (res).r += (a).r; \
67 (res).i += (a).i; \
68 } while(0)
69 #define C_SUBFROM(res, a) \
70 do { \
71 (res).r -= (a).r; \
72 (res).i -= (a).i; \
73 } while(0)
54 74
55 #define SAMP_MIN -SAMP_MAX 75 #if USE_SIMD == SIMD_SSE2
56
57 #if defined(CHECK_OVERFLOW)
58 # define CHECK_OVERFLOW_OP(a,op,b) \
59 if ( (SAMPPROD)(a) op (SAMPPROD)(b) > SAMP_MAX || (SAMPPROD)(a) op (SAMPPROD)(b) < SAMP_MIN ) { \
60 fprintf(stderr,"WARNING:overflow @ " __FILE__ "(%d): (%d " #op" %d) = %ld\n",__LINE__,(a),(b),(SAMPPROD)(a) op (SAMPPROD)(b) ); }
61 #endif
62
63
64 # define smul(a,b) ( (SAMPPROD)(a)*(b) )
65 # define sround( x ) (peck_fft_scalar)( ( (x) + (1<<(FRACBITS-1)) ) >> FRACBITS )
66
67 # define S_MUL(a,b) sround( smul(a,b) )
68
69 # define C_MUL(m,a,b) \
70 do{ (m).r = sround( smul((a).r,(b).r) - smul((a).i,(b).i) ); \
71 (m).i = sround( smul((a).r,(b).i) + smul((a).i,(b).r) ); }while(0)
72
73 # define DIVSCALAR(x,k) \
74 (x) = sround( smul( x, SAMP_MAX/k ) )
75
76 # define C_FIXDIV(c,div) \
77 do { DIVSCALAR( (c).r , div); \
78 DIVSCALAR( (c).i , div); }while (0)
79
80 # define C_MULBYSCALAR( c, s ) \
81 do{ (c).r = sround( smul( (c).r , s ) ) ;\
82 (c).i = sround( smul( (c).i , s ) ) ; }while(0)
83
84 #else /* not FIXED_POINT*/
85
86 # define S_MUL(a,b) ( (a)*(b) )
87 #define C_MUL(m,a,b) \
88 do{ (m).r = (a).r*(b).r - (a).i*(b).i;\
89 (m).i = (a).r*(b).i + (a).i*(b).r; }while(0)
90 # define C_FIXDIV(c,div) /* NOOP */
91 # define C_MULBYSCALAR( c, s ) \
92 do{ (c).r *= (s);\
93 (c).i *= (s); }while(0)
94 #endif
95
96 #ifndef CHECK_OVERFLOW_OP
97 # define CHECK_OVERFLOW_OP(a,op,b) /* noop */
98 #endif
99
100 #define C_ADD( res, a,b)\
101 do { \
102 CHECK_OVERFLOW_OP((a).r,+,(b).r)\
103 CHECK_OVERFLOW_OP((a).i,+,(b).i)\
104 (res).r=(a).r+(b).r; (res).i=(a).i+(b).i; \
105 }while(0)
106 #define C_SUB( res, a,b)\
107 do { \
108 CHECK_OVERFLOW_OP((a).r,-,(b).r)\
109 CHECK_OVERFLOW_OP((a).i,-,(b).i)\
110 (res).r=(a).r-(b).r; (res).i=(a).i-(b).i; \
111 }while(0)
112 #define C_ADDTO( res , a)\
113 do { \
114 CHECK_OVERFLOW_OP((res).r,+,(a).r)\
115 CHECK_OVERFLOW_OP((res).i,+,(a).i)\
116 (res).r += (a).r; (res).i += (a).i;\
117 }while(0)
118
119 #define C_SUBFROM( res , a)\
120 do {\
121 CHECK_OVERFLOW_OP((res).r,-,(a).r)\
122 CHECK_OVERFLOW_OP((res).i,-,(a).i)\
123 (res).r -= (a).r; (res).i -= (a).i; \
124 }while(0)
125
126
127 #ifdef FIXED_POINT
128 #define PECK_FFT_COS(phase) floorf(0.5f+SAMP_MAX * cosf(phase))
129 #define PECK_FFT_SIN(phase) floorf(0.5f+SAMP_MAX * sinf(phase))
130 #define HALF_OF(x) ((x)>>1)
131 #elif USE_SIMD == SIMD_SSE2
132 #define PECK_FFT_COS(phase) _mm_set1_ps(cosf(phase)) 76 #define PECK_FFT_COS(phase) _mm_set1_ps(cosf(phase))
133 #define PECK_FFT_SIN(phase) _mm_set1_ps(sinf(phase)) 77 #define PECK_FFT_SIN(phase) _mm_set1_ps(sinf(phase))
134 #define HALF_OF(x) ((x)*_mm_set1_ps(0.5f)) 78 #define HALF_OF(x) ((x)*_mm_set1_ps(0.5f))
135 #elif USE_SIMD == SIMD_NEON4 79 #elif USE_SIMD == SIMD_NEON4
136 #define PECK_FFT_COS(phase) vdupq_n_f32(cosf(phase)) 80 #define PECK_FFT_COS(phase) vdupq_n_f32(cosf(phase))
145 #define PECK_FFT_SIN(phase) (peck_fft_scalar) sinf(phase) 89 #define PECK_FFT_SIN(phase) (peck_fft_scalar) sinf(phase)
146 #define HALF_OF(x) ((x)*0.5f) 90 #define HALF_OF(x) ((x)*0.5f)
147 #endif 91 #endif
148 92
149 #define kf_cexp(x,phase) \ 93 #define kf_cexp(x,phase) \
150 do{ \ 94 do { \
151 (x)->r = PECK_FFT_COS(phase);\ 95 (x)->r = PECK_FFT_COS(phase); \
152 (x)->i = PECK_FFT_SIN(phase);\ 96 (x)->i = PECK_FFT_SIN(phase); \
153 }while(0) 97 } while(0)
154
155
156 /* a debugging function */
157 #define pcpx(c)\
158 fprintf(stderr,"%g + %gi\n",(double)((c)->r),(double)((c)->i) )
159
160 98
161 #ifdef PECK_FFT_USE_ALLOCA 99 #ifdef PECK_FFT_USE_ALLOCA
162 // define this to allow use of alloca instead of malloc for temporary buffers 100 // define this to allow use of alloca instead of malloc for temporary buffers
163 // Temporary buffers are used in two case: 101 // Temporary buffers are used in two case:
164 // 1. FFT sizes that have "bad" factors. i.e. not 2,3 and 5 102 // 1. FFT sizes that have "bad" factors. i.e. not 2,3 and 5

Repositories maintained by Peter Meerwald, pmeerw@pmeerw.net.