Line data Source code
1 : /*
2 : * Copyright (c) 2011 The WebRTC project authors. All Rights Reserved.
3 : *
4 : * Use of this source code is governed by a BSD-style license
5 : * that can be found in the LICENSE file in the root of the source
6 : * tree. An additional intellectual property rights grant can be found
7 : * in the file PATENTS. All contributing project authors may
8 : * be found in the AUTHORS file in the root of the source tree.
9 : */
10 :
11 :
12 : /*
13 : * This file contains the function WebRtcSpl_ComplexFFT().
14 : * The description header can be found in signal_processing_library.h
15 : *
16 : */
17 :
18 : #include "webrtc/common_audio/signal_processing/complex_fft_tables.h"
19 : #include "webrtc/common_audio/signal_processing/include/signal_processing_library.h"
20 :
21 : #define CFFTSFT 14
22 : #define CFFTRND 1
23 : #define CFFTRND2 16384
24 :
25 : #define CIFFTSFT 14
26 : #define CIFFTRND 1
27 :
28 :
29 0 : int WebRtcSpl_ComplexFFT(int16_t frfi[], int stages, int mode)
30 : {
31 : int i, j, l, k, istep, n, m;
32 : int16_t wr, wi;
33 : int32_t tr32, ti32, qr32, qi32;
34 :
35 : /* The 1024-value is a constant given from the size of kSinTable1024[],
36 : * and should not be changed depending on the input parameter 'stages'
37 : */
38 0 : n = 1 << stages;
39 0 : if (n > 1024)
40 0 : return -1;
41 :
42 0 : l = 1;
43 0 : k = 10 - 1; /* Constant for given kSinTable1024[]. Do not change
44 : depending on the input parameter 'stages' */
45 :
46 0 : if (mode == 0)
47 : {
48 : // mode==0: Low-complexity and Low-accuracy mode
49 0 : while (l < n)
50 : {
51 0 : istep = l << 1;
52 :
53 0 : for (m = 0; m < l; ++m)
54 : {
55 0 : j = m << k;
56 :
57 : /* The 256-value is a constant given as 1/4 of the size of
58 : * kSinTable1024[], and should not be changed depending on the input
59 : * parameter 'stages'. It will result in 0 <= j < N_SINE_WAVE/2
60 : */
61 0 : wr = kSinTable1024[j + 256];
62 0 : wi = -kSinTable1024[j];
63 :
64 0 : for (i = m; i < n; i += istep)
65 : {
66 0 : j = i + l;
67 :
68 0 : tr32 = (wr * frfi[2 * j] - wi * frfi[2 * j + 1]) >> 15;
69 :
70 0 : ti32 = (wr * frfi[2 * j + 1] + wi * frfi[2 * j]) >> 15;
71 :
72 0 : qr32 = (int32_t)frfi[2 * i];
73 0 : qi32 = (int32_t)frfi[2 * i + 1];
74 0 : frfi[2 * j] = (int16_t)((qr32 - tr32) >> 1);
75 0 : frfi[2 * j + 1] = (int16_t)((qi32 - ti32) >> 1);
76 0 : frfi[2 * i] = (int16_t)((qr32 + tr32) >> 1);
77 0 : frfi[2 * i + 1] = (int16_t)((qi32 + ti32) >> 1);
78 : }
79 : }
80 :
81 0 : --k;
82 0 : l = istep;
83 :
84 : }
85 :
86 : } else
87 : {
88 : // mode==1: High-complexity and High-accuracy mode
89 0 : while (l < n)
90 : {
91 0 : istep = l << 1;
92 :
93 0 : for (m = 0; m < l; ++m)
94 : {
95 0 : j = m << k;
96 :
97 : /* The 256-value is a constant given as 1/4 of the size of
98 : * kSinTable1024[], and should not be changed depending on the input
99 : * parameter 'stages'. It will result in 0 <= j < N_SINE_WAVE/2
100 : */
101 0 : wr = kSinTable1024[j + 256];
102 0 : wi = -kSinTable1024[j];
103 :
104 : #ifdef WEBRTC_ARCH_ARM_V7
105 : int32_t wri = 0;
106 : __asm __volatile("pkhbt %0, %1, %2, lsl #16" : "=r"(wri) :
107 : "r"((int32_t)wr), "r"((int32_t)wi));
108 : #endif
109 :
110 0 : for (i = m; i < n; i += istep)
111 : {
112 0 : j = i + l;
113 :
114 : #ifdef WEBRTC_ARCH_ARM_V7
115 : register int32_t frfi_r;
116 : __asm __volatile(
117 : "pkhbt %[frfi_r], %[frfi_even], %[frfi_odd],"
118 : " lsl #16\n\t"
119 : "smlsd %[tr32], %[wri], %[frfi_r], %[cfftrnd]\n\t"
120 : "smladx %[ti32], %[wri], %[frfi_r], %[cfftrnd]\n\t"
121 : :[frfi_r]"=&r"(frfi_r),
122 : [tr32]"=&r"(tr32),
123 : [ti32]"=r"(ti32)
124 : :[frfi_even]"r"((int32_t)frfi[2*j]),
125 : [frfi_odd]"r"((int32_t)frfi[2*j +1]),
126 : [wri]"r"(wri),
127 : [cfftrnd]"r"(CFFTRND));
128 : #else
129 0 : tr32 = wr * frfi[2 * j] - wi * frfi[2 * j + 1] + CFFTRND;
130 :
131 0 : ti32 = wr * frfi[2 * j + 1] + wi * frfi[2 * j] + CFFTRND;
132 : #endif
133 :
134 0 : tr32 >>= 15 - CFFTSFT;
135 0 : ti32 >>= 15 - CFFTSFT;
136 :
137 0 : qr32 = ((int32_t)frfi[2 * i]) << CFFTSFT;
138 0 : qi32 = ((int32_t)frfi[2 * i + 1]) << CFFTSFT;
139 :
140 0 : frfi[2 * j] = (int16_t)(
141 0 : (qr32 - tr32 + CFFTRND2) >> (1 + CFFTSFT));
142 0 : frfi[2 * j + 1] = (int16_t)(
143 0 : (qi32 - ti32 + CFFTRND2) >> (1 + CFFTSFT));
144 0 : frfi[2 * i] = (int16_t)(
145 0 : (qr32 + tr32 + CFFTRND2) >> (1 + CFFTSFT));
146 0 : frfi[2 * i + 1] = (int16_t)(
147 0 : (qi32 + ti32 + CFFTRND2) >> (1 + CFFTSFT));
148 : }
149 : }
150 :
151 0 : --k;
152 0 : l = istep;
153 : }
154 : }
155 0 : return 0;
156 : }
157 :
158 0 : int WebRtcSpl_ComplexIFFT(int16_t frfi[], int stages, int mode)
159 : {
160 : size_t i, j, l, istep, n, m;
161 : int k, scale, shift;
162 : int16_t wr, wi;
163 : int32_t tr32, ti32, qr32, qi32;
164 : int32_t tmp32, round2;
165 :
166 : /* The 1024-value is a constant given from the size of kSinTable1024[],
167 : * and should not be changed depending on the input parameter 'stages'
168 : */
169 0 : n = 1 << stages;
170 0 : if (n > 1024)
171 0 : return -1;
172 :
173 0 : scale = 0;
174 :
175 0 : l = 1;
176 0 : k = 10 - 1; /* Constant for given kSinTable1024[]. Do not change
177 : depending on the input parameter 'stages' */
178 :
179 0 : while (l < n)
180 : {
181 : // variable scaling, depending upon data
182 0 : shift = 0;
183 0 : round2 = 8192;
184 :
185 0 : tmp32 = WebRtcSpl_MaxAbsValueW16(frfi, 2 * n);
186 0 : if (tmp32 > 13573)
187 : {
188 0 : shift++;
189 0 : scale++;
190 0 : round2 <<= 1;
191 : }
192 0 : if (tmp32 > 27146)
193 : {
194 0 : shift++;
195 0 : scale++;
196 0 : round2 <<= 1;
197 : }
198 :
199 0 : istep = l << 1;
200 :
201 0 : if (mode == 0)
202 : {
203 : // mode==0: Low-complexity and Low-accuracy mode
204 0 : for (m = 0; m < l; ++m)
205 : {
206 0 : j = m << k;
207 :
208 : /* The 256-value is a constant given as 1/4 of the size of
209 : * kSinTable1024[], and should not be changed depending on the input
210 : * parameter 'stages'. It will result in 0 <= j < N_SINE_WAVE/2
211 : */
212 0 : wr = kSinTable1024[j + 256];
213 0 : wi = kSinTable1024[j];
214 :
215 0 : for (i = m; i < n; i += istep)
216 : {
217 0 : j = i + l;
218 :
219 0 : tr32 = (wr * frfi[2 * j] - wi * frfi[2 * j + 1]) >> 15;
220 :
221 0 : ti32 = (wr * frfi[2 * j + 1] + wi * frfi[2 * j]) >> 15;
222 :
223 0 : qr32 = (int32_t)frfi[2 * i];
224 0 : qi32 = (int32_t)frfi[2 * i + 1];
225 0 : frfi[2 * j] = (int16_t)((qr32 - tr32) >> shift);
226 0 : frfi[2 * j + 1] = (int16_t)((qi32 - ti32) >> shift);
227 0 : frfi[2 * i] = (int16_t)((qr32 + tr32) >> shift);
228 0 : frfi[2 * i + 1] = (int16_t)((qi32 + ti32) >> shift);
229 : }
230 : }
231 : } else
232 : {
233 : // mode==1: High-complexity and High-accuracy mode
234 :
235 0 : for (m = 0; m < l; ++m)
236 : {
237 0 : j = m << k;
238 :
239 : /* The 256-value is a constant given as 1/4 of the size of
240 : * kSinTable1024[], and should not be changed depending on the input
241 : * parameter 'stages'. It will result in 0 <= j < N_SINE_WAVE/2
242 : */
243 0 : wr = kSinTable1024[j + 256];
244 0 : wi = kSinTable1024[j];
245 :
246 : #ifdef WEBRTC_ARCH_ARM_V7
247 : int32_t wri = 0;
248 : __asm __volatile("pkhbt %0, %1, %2, lsl #16" : "=r"(wri) :
249 : "r"((int32_t)wr), "r"((int32_t)wi));
250 : #endif
251 :
252 0 : for (i = m; i < n; i += istep)
253 : {
254 0 : j = i + l;
255 :
256 : #ifdef WEBRTC_ARCH_ARM_V7
257 : register int32_t frfi_r;
258 : __asm __volatile(
259 : "pkhbt %[frfi_r], %[frfi_even], %[frfi_odd], lsl #16\n\t"
260 : "smlsd %[tr32], %[wri], %[frfi_r], %[cifftrnd]\n\t"
261 : "smladx %[ti32], %[wri], %[frfi_r], %[cifftrnd]\n\t"
262 : :[frfi_r]"=&r"(frfi_r),
263 : [tr32]"=&r"(tr32),
264 : [ti32]"=r"(ti32)
265 : :[frfi_even]"r"((int32_t)frfi[2*j]),
266 : [frfi_odd]"r"((int32_t)frfi[2*j +1]),
267 : [wri]"r"(wri),
268 : [cifftrnd]"r"(CIFFTRND)
269 : );
270 : #else
271 :
272 0 : tr32 = wr * frfi[2 * j] - wi * frfi[2 * j + 1] + CIFFTRND;
273 :
274 0 : ti32 = wr * frfi[2 * j + 1] + wi * frfi[2 * j] + CIFFTRND;
275 : #endif
276 0 : tr32 >>= 15 - CIFFTSFT;
277 0 : ti32 >>= 15 - CIFFTSFT;
278 :
279 0 : qr32 = ((int32_t)frfi[2 * i]) << CIFFTSFT;
280 0 : qi32 = ((int32_t)frfi[2 * i + 1]) << CIFFTSFT;
281 :
282 0 : frfi[2 * j] = (int16_t)(
283 0 : (qr32 - tr32 + round2) >> (shift + CIFFTSFT));
284 0 : frfi[2 * j + 1] = (int16_t)(
285 0 : (qi32 - ti32 + round2) >> (shift + CIFFTSFT));
286 0 : frfi[2 * i] = (int16_t)(
287 0 : (qr32 + tr32 + round2) >> (shift + CIFFTSFT));
288 0 : frfi[2 * i + 1] = (int16_t)(
289 0 : (qi32 + ti32 + round2) >> (shift + CIFFTSFT));
290 : }
291 : }
292 :
293 : }
294 0 : --k;
295 0 : l = istep;
296 : }
297 0 : return scale;
298 : }
|