1 /* Copyright (C) 2007-2008 Jean-Marc Valin
2 Copyright (C) 2008 Thorvald Natvig
3
4 File: resample.c
5 Arbitrary resampling code
6
7 Redistribution and use in source and binary forms, with or without
8 modification, are permitted provided that the following conditions are
9 met:
10
11 1. Redistributions of source code must retain the above copyright notice,
12 this list of conditions and the following disclaimer.
13
14 2. Redistributions in binary form must reproduce the above copyright
15 notice, this list of conditions and the following disclaimer in the
16 documentation and/or other materials provided with the distribution.
17
18 3. The name of the author may not be used to endorse or promote products
19 derived from this software without specific prior written permission.
20
21 THIS SOFTWARE IS PROVIDED BY THE AUTHOR ``AS IS'' AND ANY EXPRESS OR
22 IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES
23 OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
24 DISCLAIMED. IN NO EVENT SHALL THE AUTHOR BE LIABLE FOR ANY DIRECT,
25 INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES
26 (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR
27 SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION)
28 HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT,
29 STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN
30 ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
31 POSSIBILITY OF SUCH DAMAGE.
32 */
33
34 /*
35 The design goals of this code are:
36 - Very fast algorithm
37 - SIMD-friendly algorithm
38 - Low memory requirement
39 - Good *perceptual* quality (and not best SNR)
40
41 Warning: This resampler is relatively new. Although I think I got rid of
42 all the major bugs and I don't expect the API to change anymore, there
43 may be something I've missed. So use with caution.
44
45 This algorithm is based on this original resampling algorithm:
46 Smith, Julius O. Digital Audio Resampling Home Page
47 Center for Computer Research in Music and Acoustics (CCRMA),
48 Stanford University, 2007.
49 Web published at http://www-ccrma.stanford.edu/~jos/resample/.
50
51 There is one main difference, though. This resampler uses cubic
52 interpolation instead of linear interpolation in the above paper. This
53 makes the table much smaller and makes it possible to compute that table
54 on a per-stream basis. In turn, being able to tweak the table for each
55 stream makes it possible to both reduce complexity on simple ratios
56 (e.g. 2/3), and get rid of the rounding operations in the inner loop.
57 The latter both reduces CPU time and makes the algorithm more SIMD-friendly.
58 */
59
60 #include <config.h>
61
62 #define OUTSIDE_SPEEX 1
63 #define FLOATING_POINT 1
64
65 #ifdef OUTSIDE_SPEEX
66 #include <glib.h>
speex_alloc(size_t size)67 static void *speex_alloc (size_t size) {return g_malloc0(size);}
speex_realloc(void * ptr,size_t size)68 static void *speex_realloc (void *ptr, size_t size) {return g_realloc(ptr, size);}
speex_free(void * ptr)69 static void speex_free (void *ptr) {g_free(ptr);}
70 #include "speex_resampler.h"
71 #include "arch.h"
72 #else /* OUTSIDE_SPEEX */
73
74 #include "speex/speex_resampler.h"
75 #include "arch.h"
76 #include "os_support.h"
77 #endif /* OUTSIDE_SPEEX */
78
79 #include "stack_alloc.h"
80 #include <math.h>
81 #include <limits.h>
82
83 #include "ws_attributes.h"
84
85 #ifndef M_PI
86 #define M_PI 3.14159265358979323846
87 #endif
88
89 #define IMAX(a,b) ((a) > (b) ? (a) : (b))
90 #define IMIN(a,b) ((a) < (b) ? (a) : (b))
91
92 #ifndef NULL
93 #define NULL 0
94 #endif
95
96 #ifdef _USE_SSE
97 #include "resample_sse.h"
98 #endif
99
100 #ifdef _USE_NEON
101 #include "resample_neon.h"
102 #endif
103
104 /* Numer of elements to allocate on the stack */
105 #ifdef VAR_ARRAYS
106 #define FIXED_STACK_ALLOC 8192
107 #else
108 #define FIXED_STACK_ALLOC 1024
109 #endif
110
111 typedef int (*resampler_basic_func)(SpeexResamplerState *, spx_uint32_t , const spx_word16_t *, spx_uint32_t *, spx_word16_t *, spx_uint32_t *);
112
113 struct SpeexResamplerState_ {
114 spx_uint32_t in_rate;
115 spx_uint32_t out_rate;
116 spx_uint32_t num_rate;
117 spx_uint32_t den_rate;
118
119 int quality;
120 spx_uint32_t nb_channels;
121 spx_uint32_t filt_len;
122 spx_uint32_t mem_alloc_size;
123 spx_uint32_t buffer_size;
124 int int_advance;
125 int frac_advance;
126 float cutoff;
127 spx_uint32_t oversample;
128 int initialised;
129 int started;
130
131 /* These are per-channel */
132 spx_int32_t *last_sample;
133 spx_uint32_t *samp_frac_num;
134 spx_uint32_t *magic_samples;
135
136 spx_word16_t *mem;
137 spx_word16_t *sinc_table;
138 spx_uint32_t sinc_table_length;
139 resampler_basic_func resampler_ptr;
140
141 int in_stride;
142 int out_stride;
143 } ;
144
145 static const double kaiser12_table[68] = {
146 0.99859849, 1.00000000, 0.99859849, 0.99440475, 0.98745105, 0.97779076,
147 0.96549770, 0.95066529, 0.93340547, 0.91384741, 0.89213598, 0.86843014,
148 0.84290116, 0.81573067, 0.78710866, 0.75723148, 0.72629970, 0.69451601,
149 0.66208321, 0.62920216, 0.59606986, 0.56287762, 0.52980938, 0.49704014,
150 0.46473455, 0.43304576, 0.40211431, 0.37206735, 0.34301800, 0.31506490,
151 0.28829195, 0.26276832, 0.23854851, 0.21567274, 0.19416736, 0.17404546,
152 0.15530766, 0.13794294, 0.12192957, 0.10723616, 0.09382272, 0.08164178,
153 0.07063950, 0.06075685, 0.05193064, 0.04409466, 0.03718069, 0.03111947,
154 0.02584161, 0.02127838, 0.01736250, 0.01402878, 0.01121463, 0.00886058,
155 0.00691064, 0.00531256, 0.00401805, 0.00298291, 0.00216702, 0.00153438,
156 0.00105297, 0.00069463, 0.00043489, 0.00025272, 0.00013031, 0.0000527734,
157 0.00001000, 0.00000000};
158 /*
159 static const double kaiser12_table[36] = {
160 0.99440475, 1.00000000, 0.99440475, 0.97779076, 0.95066529, 0.91384741,
161 0.86843014, 0.81573067, 0.75723148, 0.69451601, 0.62920216, 0.56287762,
162 0.49704014, 0.43304576, 0.37206735, 0.31506490, 0.26276832, 0.21567274,
163 0.17404546, 0.13794294, 0.10723616, 0.08164178, 0.06075685, 0.04409466,
164 0.03111947, 0.02127838, 0.01402878, 0.00886058, 0.00531256, 0.00298291,
165 0.00153438, 0.00069463, 0.00025272, 0.0000527734, 0.00000500, 0.00000000};
166 */
167 static const double kaiser10_table[36] = {
168 0.99537781, 1.00000000, 0.99537781, 0.98162644, 0.95908712, 0.92831446,
169 0.89005583, 0.84522401, 0.79486424, 0.74011713, 0.68217934, 0.62226347,
170 0.56155915, 0.50119680, 0.44221549, 0.38553619, 0.33194107, 0.28205962,
171 0.23636152, 0.19515633, 0.15859932, 0.12670280, 0.09935205, 0.07632451,
172 0.05731132, 0.04193980, 0.02979584, 0.02044510, 0.01345224, 0.00839739,
173 0.00488951, 0.00257636, 0.00115101, 0.00035515, 0.00000000, 0.00000000};
174
175 static const double kaiser8_table[36] = {
176 0.99635258, 1.00000000, 0.99635258, 0.98548012, 0.96759014, 0.94302200,
177 0.91223751, 0.87580811, 0.83439927, 0.78875245, 0.73966538, 0.68797126,
178 0.63451750, 0.58014482, 0.52566725, 0.47185369, 0.41941150, 0.36897272,
179 0.32108304, 0.27619388, 0.23465776, 0.19672670, 0.16255380, 0.13219758,
180 0.10562887, 0.08273982, 0.06335451, 0.04724088, 0.03412321, 0.02369490,
181 0.01563093, 0.00959968, 0.00527363, 0.00233883, 0.00050000, 0.00000000};
182
183 static const double kaiser6_table[36] = {
184 0.99733006, 1.00000000, 0.99733006, 0.98935595, 0.97618418, 0.95799003,
185 0.93501423, 0.90755855, 0.87598009, 0.84068475, 0.80211977, 0.76076565,
186 0.71712752, 0.67172623, 0.62508937, 0.57774224, 0.53019925, 0.48295561,
187 0.43647969, 0.39120616, 0.34752997, 0.30580127, 0.26632152, 0.22934058,
188 0.19505503, 0.16360756, 0.13508755, 0.10953262, 0.08693120, 0.06722600,
189 0.05031820, 0.03607231, 0.02432151, 0.01487334, 0.00752000, 0.00000000};
190
191 struct FuncDef {
192 const double *table;
193 int oversample;
194 };
195
196 static const struct FuncDef _KAISER12 = {kaiser12_table, 64};
197 #define KAISER12 (&_KAISER12)
198 /*static struct FuncDef _KAISER12 = {kaiser12_table, 32};
199 #define KAISER12 (&_KAISER12)*/
200 static const struct FuncDef _KAISER10 = {kaiser10_table, 32};
201 #define KAISER10 (&_KAISER10)
202 static const struct FuncDef _KAISER8 = {kaiser8_table, 32};
203 #define KAISER8 (&_KAISER8)
204 static const struct FuncDef _KAISER6 = {kaiser6_table, 32};
205 #define KAISER6 (&_KAISER6)
206
207 struct QualityMapping {
208 int base_length;
209 int oversample;
210 float downsample_bandwidth;
211 float upsample_bandwidth;
212 const struct FuncDef *window_func;
213 };
214
215
216 /* This table maps conversion quality to internal parameters. There are two
217 reasons that explain why the up-sampling bandwidth is larger than the
218 down-sampling bandwidth:
219 1) When up-sampling, we can assume that the spectrum is already attenuated
220 close to the Nyquist rate (from an A/D or a previous resampling filter)
221 2) Any aliasing that occurs very close to the Nyquist rate will be masked
222 by the sinusoids/noise just below the Nyquist rate (guaranteed only for
223 up-sampling).
224 */
225 static const struct QualityMapping quality_map[11] = {
226 { 8, 4, 0.830f, 0.860f, KAISER6 }, /* Q0 */
227 { 16, 4, 0.850f, 0.880f, KAISER6 }, /* Q1 */
228 { 32, 4, 0.882f, 0.910f, KAISER6 }, /* Q2 */ /* 82.3% cutoff ( ~60 dB stop) 6 */
229 { 48, 8, 0.895f, 0.917f, KAISER8 }, /* Q3 */ /* 84.9% cutoff ( ~80 dB stop) 8 */
230 { 64, 8, 0.921f, 0.940f, KAISER8 }, /* Q4 */ /* 88.7% cutoff ( ~80 dB stop) 8 */
231 { 80, 16, 0.922f, 0.940f, KAISER10}, /* Q5 */ /* 89.1% cutoff (~100 dB stop) 10 */
232 { 96, 16, 0.940f, 0.945f, KAISER10}, /* Q6 */ /* 91.5% cutoff (~100 dB stop) 10 */
233 {128, 16, 0.950f, 0.950f, KAISER10}, /* Q7 */ /* 93.1% cutoff (~100 dB stop) 10 */
234 {160, 16, 0.960f, 0.960f, KAISER10}, /* Q8 */ /* 94.5% cutoff (~100 dB stop) 10 */
235 {192, 32, 0.968f, 0.968f, KAISER12}, /* Q9 */ /* 95.5% cutoff (~100 dB stop) 10 */
236 {256, 32, 0.975f, 0.975f, KAISER12}, /* Q10 */ /* 96.6% cutoff (~100 dB stop) 10 */
237 };
238 /*8,24,40,56,80,104,128,160,200,256,320*/
compute_func(float x,const struct FuncDef * func)239 static double compute_func(float x, const struct FuncDef *func)
240 {
241 float y, frac;
242 double interp[4];
243 int ind;
244 y = x*func->oversample;
245 ind = (int)floor(y);
246 frac = (y-ind);
247 /* CSE with handle the repeated powers */
248 interp[3] = -0.1666666667*frac + 0.1666666667*(frac*frac*frac);
249 interp[2] = frac + 0.5*(frac*frac) - 0.5*(frac*frac*frac);
250 /*interp[2] = 1.f - 0.5f*frac - frac*frac + 0.5f*frac*frac*frac;*/
251 interp[0] = -0.3333333333*frac + 0.5*(frac*frac) - 0.1666666667*(frac*frac*frac);
252 /* Just to make sure we don't have rounding problems */
253 interp[1] = 1.f-interp[3]-interp[2]-interp[0];
254
255 /*sum = frac*accum[1] + (1-frac)*accum[2];*/
256 return interp[0]*func->table[ind] + interp[1]*func->table[ind+1] + interp[2]*func->table[ind+2] + interp[3]*func->table[ind+3];
257 }
258
259 #if 0
260 #include <stdio.h>
261 int main(int argc, char **argv)
262 {
263 int i;
264 for (i=0;i<256;i++)
265 {
266 printf ("%f\n", compute_func(i/256., KAISER12));
267 }
268 return 0;
269 }
270 #endif
271
272 DIAG_OFF_CLANG(self-assign) /* SATURATE32PSHR */
273
274 #ifdef FIXED_POINT
275 /* The slow way of computing a sinc for the table. Should improve that some day */
sinc(float cutoff,float x,int N,const struct FuncDef * window_func)276 static spx_word16_t sinc(float cutoff, float x, int N, const struct FuncDef *window_func)
277 {
278 /*fprintf (stderr, "%f ", x);*/
279 float xx = x * cutoff;
280 if (fabs(x)<1e-6f)
281 return WORD2INT(32768.*cutoff);
282 else if (fabs(x) > .5f*N)
283 return 0;
284 /*FIXME: Can it really be any slower than this? */
285 return WORD2INT(32768.*cutoff*sin(M_PI*xx)/(M_PI*xx) * compute_func(fabs(2.*x/N), window_func));
286 }
287 #else
288 /* The slow way of computing a sinc for the table. Should improve that some day */
289 static spx_word16_t sinc(float cutoff, float x, int N, const struct FuncDef *window_func)
290 {
291 /*fprintf (stderr, "%f ", x);*/
292 float xx = x * cutoff;
293 if (fabs(x)<1e-6)
294 return cutoff;
295 else if (fabs(x) > .5*N)
296 return 0;
297 /*FIXME: Can it really be any slower than this? */
298 return (spx_word16_t)(cutoff*sin(M_PI*xx)/(M_PI*xx) * compute_func((float)fabs(2.*x/N), window_func));
299 }
300 #endif
301
302 #ifdef FIXED_POINT
cubic_coef(spx_word16_t x,spx_word16_t interp[4])303 static void cubic_coef(spx_word16_t x, spx_word16_t interp[4])
304 {
305 /* Compute interpolation coefficients. I'm not sure whether this corresponds to cubic interpolation
306 but I know it's MMSE-optimal on a sinc */
307 spx_word16_t x2, x3;
308 x2 = MULT16_16_P15(x, x);
309 x3 = MULT16_16_P15(x, x2);
310 interp[0] = PSHR32(MULT16_16(QCONST16(-0.16667f, 15),x) + MULT16_16(QCONST16(0.16667f, 15),x3),15);
311 interp[1] = EXTRACT16(EXTEND32(x) + SHR32(SUB32(EXTEND32(x2),EXTEND32(x3)),1));
312 interp[3] = PSHR32(MULT16_16(QCONST16(-0.33333f, 15),x) + MULT16_16(QCONST16(.5f,15),x2) - MULT16_16(QCONST16(0.16667f, 15),x3),15);
313 /* Just to make sure we don't have rounding problems */
314 interp[2] = Q15_ONE-interp[0]-interp[1]-interp[3];
315 if (interp[2]<32767)
316 interp[2]+=1;
317 }
318 #else
cubic_coef(spx_word16_t frac,spx_word16_t interp[4])319 static void cubic_coef(spx_word16_t frac, spx_word16_t interp[4])
320 {
321 /* Compute interpolation coefficients. I'm not sure whether this corresponds to cubic interpolation
322 but I know it's MMSE-optimal on a sinc */
323 interp[0] = -0.16667f*frac + 0.16667f*frac*frac*frac;
324 interp[1] = frac + 0.5f*frac*frac - 0.5f*frac*frac*frac;
325 /*interp[2] = 1.f - 0.5f*frac - frac*frac + 0.5f*frac*frac*frac;*/
326 interp[3] = -0.33333f*frac + 0.5f*frac*frac - 0.16667f*frac*frac*frac;
327 /* Just to make sure we don't have rounding problems */
328 interp[2] = (spx_word16_t)(1.-interp[0]-interp[1]-interp[3]);
329 }
330 #endif
331
resampler_basic_direct_single(SpeexResamplerState * st,spx_uint32_t channel_index,const spx_word16_t * in,spx_uint32_t * in_len,spx_word16_t * out,spx_uint32_t * out_len)332 static int resampler_basic_direct_single(SpeexResamplerState *st, spx_uint32_t channel_index, const spx_word16_t *in, spx_uint32_t *in_len, spx_word16_t *out, spx_uint32_t *out_len)
333 {
334 const int N = st->filt_len;
335 int out_sample = 0;
336 int last_sample = st->last_sample[channel_index];
337 spx_uint32_t samp_frac_num = st->samp_frac_num[channel_index];
338 const spx_word16_t *sinc_table = st->sinc_table;
339 const int out_stride = st->out_stride;
340 const int int_advance = st->int_advance;
341 const int frac_advance = st->frac_advance;
342 const spx_uint32_t den_rate = st->den_rate;
343 spx_word32_t sum;
344
345 while (!(last_sample >= (spx_int32_t)*in_len || out_sample >= (spx_int32_t)*out_len))
346 {
347 const spx_word16_t *sinct = & sinc_table[samp_frac_num*N];
348 const spx_word16_t *iptr = & in[last_sample];
349
350 #ifndef OVERRIDE_INNER_PRODUCT_SINGLE
351 int j;
352 sum = 0;
353 for(j=0;j<N;j++) sum += MULT16_16(sinct[j], iptr[j]);
354
355 /* This code is slower on most DSPs which have only 2 accumulators.
356 Plus this this forces truncation to 32 bits and you lose the HW guard bits.
357 I think we can trust the compiler and let it vectorize and/or unroll itself.
358 spx_word32_t accum[4] = {0,0,0,0};
359 for(j=0;j<N;j+=4) {
360 accum[0] += MULT16_16(sinct[j], iptr[j]);
361 accum[1] += MULT16_16(sinct[j+1], iptr[j+1]);
362 accum[2] += MULT16_16(sinct[j+2], iptr[j+2]);
363 accum[3] += MULT16_16(sinct[j+3], iptr[j+3]);
364 }
365 sum = accum[0] + accum[1] + accum[2] + accum[3];
366 */
367 sum = SATURATE32PSHR(sum, 15, 32767);
368 #else
369 sum = inner_product_single(sinct, iptr, N);
370 #endif
371
372 out[out_stride * out_sample++] = sum;
373 last_sample += int_advance;
374 samp_frac_num += frac_advance;
375 if (samp_frac_num >= den_rate)
376 {
377 samp_frac_num -= den_rate;
378 last_sample++;
379 }
380 }
381
382 st->last_sample[channel_index] = last_sample;
383 st->samp_frac_num[channel_index] = samp_frac_num;
384 return out_sample;
385 }
386
387 #ifdef FIXED_POINT
388 #else
389 /* This is the same as the previous function, except with a double-precision accumulator */
resampler_basic_direct_double(SpeexResamplerState * st,spx_uint32_t channel_index,const spx_word16_t * in,spx_uint32_t * in_len,spx_word16_t * out,spx_uint32_t * out_len)390 static int resampler_basic_direct_double(SpeexResamplerState *st, spx_uint32_t channel_index, const spx_word16_t *in, spx_uint32_t *in_len, spx_word16_t *out, spx_uint32_t *out_len)
391 {
392 const int N = st->filt_len;
393 int out_sample = 0;
394 int last_sample = st->last_sample[channel_index];
395 spx_uint32_t samp_frac_num = st->samp_frac_num[channel_index];
396 const spx_word16_t *sinc_table = st->sinc_table;
397 const int out_stride = st->out_stride;
398 const int int_advance = st->int_advance;
399 const int frac_advance = st->frac_advance;
400 const spx_uint32_t den_rate = st->den_rate;
401 double sum;
402
403 while (!(last_sample >= (spx_int32_t)*in_len || out_sample >= (spx_int32_t)*out_len))
404 {
405 const spx_word16_t *sinct = & sinc_table[samp_frac_num*N];
406 const spx_word16_t *iptr = & in[last_sample];
407
408 #ifndef OVERRIDE_INNER_PRODUCT_DOUBLE
409 int j;
410 double accum[4] = {0,0,0,0};
411
412 for(j=0;j<N;j+=4) {
413 accum[0] += sinct[j]*iptr[j];
414 accum[1] += sinct[j+1]*iptr[j+1];
415 accum[2] += sinct[j+2]*iptr[j+2];
416 accum[3] += sinct[j+3]*iptr[j+3];
417 }
418 sum = accum[0] + accum[1] + accum[2] + accum[3];
419 #else
420 sum = inner_product_double(sinct, iptr, N);
421 #endif
422
423 out[out_stride * out_sample++] = PSHR32(sum, 15);
424 last_sample += int_advance;
425 samp_frac_num += frac_advance;
426 if (samp_frac_num >= den_rate)
427 {
428 samp_frac_num -= den_rate;
429 last_sample++;
430 }
431 }
432
433 st->last_sample[channel_index] = last_sample;
434 st->samp_frac_num[channel_index] = samp_frac_num;
435 return out_sample;
436 }
437 #endif
438
resampler_basic_interpolate_single(SpeexResamplerState * st,spx_uint32_t channel_index,const spx_word16_t * in,spx_uint32_t * in_len,spx_word16_t * out,spx_uint32_t * out_len)439 static int resampler_basic_interpolate_single(SpeexResamplerState *st, spx_uint32_t channel_index, const spx_word16_t *in, spx_uint32_t *in_len, spx_word16_t *out, spx_uint32_t *out_len)
440 {
441 const int N = st->filt_len;
442 int out_sample = 0;
443 int last_sample = st->last_sample[channel_index];
444 spx_uint32_t samp_frac_num = st->samp_frac_num[channel_index];
445 const int out_stride = st->out_stride;
446 const int int_advance = st->int_advance;
447 const int frac_advance = st->frac_advance;
448 const spx_uint32_t den_rate = st->den_rate;
449 spx_word32_t sum;
450
451 while (!(last_sample >= (spx_int32_t)*in_len || out_sample >= (spx_int32_t)*out_len))
452 {
453 const spx_word16_t *iptr = & in[last_sample];
454
455 const int offset = samp_frac_num*st->oversample/st->den_rate;
456 #ifdef FIXED_POINT
457 const spx_word16_t frac = PDIV32(SHL32((samp_frac_num*st->oversample) % st->den_rate,15),st->den_rate);
458 #else
459 const spx_word16_t frac = ((float)((samp_frac_num*st->oversample) % st->den_rate))/st->den_rate;
460 #endif
461 spx_word16_t interp[4];
462
463
464 #ifndef OVERRIDE_INTERPOLATE_PRODUCT_SINGLE
465 int j;
466 spx_word32_t accum[4] = {0,0,0,0};
467
468 for(j=0;j<N;j++) {
469 const spx_word16_t curr_in=iptr[j];
470 accum[0] += MULT16_16(curr_in,st->sinc_table[4+(j+1)*st->oversample-offset-2]);
471 accum[1] += MULT16_16(curr_in,st->sinc_table[4+(j+1)*st->oversample-offset-1]);
472 accum[2] += MULT16_16(curr_in,st->sinc_table[4+(j+1)*st->oversample-offset]);
473 accum[3] += MULT16_16(curr_in,st->sinc_table[4+(j+1)*st->oversample-offset+1]);
474 }
475
476 cubic_coef(frac, interp);
477 sum = MULT16_32_Q15(interp[0],SHR32(accum[0], 1)) + MULT16_32_Q15(interp[1],SHR32(accum[1], 1)) + MULT16_32_Q15(interp[2],SHR32(accum[2], 1)) + MULT16_32_Q15(interp[3],SHR32(accum[3], 1));
478 sum = SATURATE32PSHR(sum, 15, 32767);
479 #else
480 cubic_coef(frac, interp);
481 sum = interpolate_product_single(iptr, st->sinc_table + st->oversample + 4 - offset - 2, N, st->oversample, interp);
482 #endif
483
484 out[out_stride * out_sample++] = sum;
485 last_sample += int_advance;
486 samp_frac_num += frac_advance;
487 if (samp_frac_num >= den_rate)
488 {
489 samp_frac_num -= den_rate;
490 last_sample++;
491 }
492 }
493
494 st->last_sample[channel_index] = last_sample;
495 st->samp_frac_num[channel_index] = samp_frac_num;
496 return out_sample;
497 }
498
499 #ifdef FIXED_POINT
500 #else
501 /* This is the same as the previous function, except with a double-precision accumulator */
resampler_basic_interpolate_double(SpeexResamplerState * st,spx_uint32_t channel_index,const spx_word16_t * in,spx_uint32_t * in_len,spx_word16_t * out,spx_uint32_t * out_len)502 static int resampler_basic_interpolate_double(SpeexResamplerState *st, spx_uint32_t channel_index, const spx_word16_t *in, spx_uint32_t *in_len, spx_word16_t *out, spx_uint32_t *out_len)
503 {
504 const int N = st->filt_len;
505 int out_sample = 0;
506 int last_sample = st->last_sample[channel_index];
507 spx_uint32_t samp_frac_num = st->samp_frac_num[channel_index];
508 const int out_stride = st->out_stride;
509 const int int_advance = st->int_advance;
510 const int frac_advance = st->frac_advance;
511 const spx_uint32_t den_rate = st->den_rate;
512 spx_word32_t sum;
513
514 while (!(last_sample >= (spx_int32_t)*in_len || out_sample >= (spx_int32_t)*out_len))
515 {
516 const spx_word16_t *iptr = & in[last_sample];
517
518 const int offset = samp_frac_num*st->oversample/st->den_rate;
519 #ifdef FIXED_POINT
520 const spx_word16_t frac = PDIV32(SHL32((samp_frac_num*st->oversample) % st->den_rate,15),st->den_rate);
521 #else
522 const spx_word16_t frac = ((float)((samp_frac_num*st->oversample) % st->den_rate))/st->den_rate;
523 #endif
524 spx_word16_t interp[4];
525
526
527 #ifndef OVERRIDE_INTERPOLATE_PRODUCT_DOUBLE
528 int j;
529 double accum[4] = {0,0,0,0};
530
531 for(j=0;j<N;j++) {
532 const double curr_in=iptr[j];
533 accum[0] += MULT16_16(curr_in,st->sinc_table[4+(j+1)*st->oversample-offset-2]);
534 accum[1] += MULT16_16(curr_in,st->sinc_table[4+(j+1)*st->oversample-offset-1]);
535 accum[2] += MULT16_16(curr_in,st->sinc_table[4+(j+1)*st->oversample-offset]);
536 accum[3] += MULT16_16(curr_in,st->sinc_table[4+(j+1)*st->oversample-offset+1]);
537 }
538
539 cubic_coef(frac, interp);
540 sum = MULT16_32_Q15(interp[0],accum[0]) + MULT16_32_Q15(interp[1],accum[1]) + MULT16_32_Q15(interp[2],accum[2]) + MULT16_32_Q15(interp[3],accum[3]);
541 #else
542 cubic_coef(frac, interp);
543 sum = interpolate_product_double(iptr, st->sinc_table + st->oversample + 4 - offset - 2, N, st->oversample, interp);
544 #endif
545
546 out[out_stride * out_sample++] = PSHR32(sum,15);
547 last_sample += int_advance;
548 samp_frac_num += frac_advance;
549 if (samp_frac_num >= den_rate)
550 {
551 samp_frac_num -= den_rate;
552 last_sample++;
553 }
554 }
555
556 st->last_sample[channel_index] = last_sample;
557 st->samp_frac_num[channel_index] = samp_frac_num;
558 return out_sample;
559 }
560 #endif
561
562 /* This resampler is used to produce zero output in situations where memory
563 for the filter could not be allocated. The expected numbers of input and
564 output samples are still processed so that callers failing to check error
565 codes are not surprised, possibly getting into infinite loops. */
resampler_basic_zero(SpeexResamplerState * st,spx_uint32_t channel_index,const spx_word16_t * in _U_,spx_uint32_t * in_len,spx_word16_t * out,spx_uint32_t * out_len)566 static int resampler_basic_zero(SpeexResamplerState *st, spx_uint32_t channel_index, const spx_word16_t *in _U_, spx_uint32_t *in_len, spx_word16_t *out, spx_uint32_t *out_len)
567 {
568 int out_sample = 0;
569 int last_sample = st->last_sample[channel_index];
570 spx_uint32_t samp_frac_num = st->samp_frac_num[channel_index];
571 const int out_stride = st->out_stride;
572 const int int_advance = st->int_advance;
573 const int frac_advance = st->frac_advance;
574 const spx_uint32_t den_rate = st->den_rate;
575
576 while (!(last_sample >= (spx_int32_t)*in_len || out_sample >= (spx_int32_t)*out_len))
577 {
578 out[out_stride * out_sample++] = 0;
579 last_sample += int_advance;
580 samp_frac_num += frac_advance;
581 if (samp_frac_num >= den_rate)
582 {
583 samp_frac_num -= den_rate;
584 last_sample++;
585 }
586 }
587
588 st->last_sample[channel_index] = last_sample;
589 st->samp_frac_num[channel_index] = samp_frac_num;
590 return out_sample;
591 }
592
update_filter(SpeexResamplerState * st)593 static int update_filter(SpeexResamplerState *st)
594 {
595 spx_uint32_t old_length = st->filt_len;
596 spx_uint32_t old_alloc_size = st->mem_alloc_size;
597 int use_direct;
598 spx_uint32_t min_sinc_table_length;
599 spx_uint32_t min_alloc_size;
600
601 st->int_advance = st->num_rate/st->den_rate;
602 st->frac_advance = st->num_rate%st->den_rate;
603 st->oversample = quality_map[st->quality].oversample;
604 st->filt_len = quality_map[st->quality].base_length;
605
606 if (st->num_rate > st->den_rate)
607 {
608 /* down-sampling */
609 st->cutoff = quality_map[st->quality].downsample_bandwidth * st->den_rate / st->num_rate;
610 /* FIXME: divide the numerator and denominator by a certain amount if they're too large */
611 st->filt_len = st->filt_len*st->num_rate / st->den_rate;
612 /* Round up to make sure we have a multiple of 8 for SSE */
613 st->filt_len = ((st->filt_len-1)&(~0x7))+8;
614 if (2*st->den_rate < st->num_rate)
615 st->oversample >>= 1;
616 if (4*st->den_rate < st->num_rate)
617 st->oversample >>= 1;
618 if (8*st->den_rate < st->num_rate)
619 st->oversample >>= 1;
620 if (16*st->den_rate < st->num_rate)
621 st->oversample >>= 1;
622 if (st->oversample < 1)
623 st->oversample = 1;
624 } else {
625 /* up-sampling */
626 st->cutoff = quality_map[st->quality].upsample_bandwidth;
627 }
628
629 /* Choose the resampling type that requires the least amount of memory */
630 #ifdef RESAMPLE_FULL_SINC_TABLE
631 use_direct = 1;
632 if (INT_MAX/sizeof(spx_word16_t)/st->den_rate < st->filt_len)
633 goto fail;
634 #else
635 use_direct = st->filt_len*st->den_rate <= st->filt_len*st->oversample+8
636 && INT_MAX/(spx_uint32_t)sizeof(spx_word16_t)/st->den_rate >= st->filt_len;
637 #endif
638 if (use_direct)
639 {
640 min_sinc_table_length = st->filt_len*st->den_rate;
641 } else {
642 if ((INT_MAX/(spx_uint32_t)sizeof(spx_word16_t)-8)/st->oversample < st->filt_len)
643 goto fail;
644
645 min_sinc_table_length = st->filt_len*st->oversample+8;
646 }
647 if (st->sinc_table_length < min_sinc_table_length)
648 {
649 spx_word16_t *sinc_table = (spx_word16_t *)speex_realloc(st->sinc_table,min_sinc_table_length*sizeof(spx_word16_t));
650 if (!sinc_table)
651 goto fail;
652
653 st->sinc_table = sinc_table;
654 st->sinc_table_length = min_sinc_table_length;
655 }
656 if (use_direct)
657 {
658 spx_uint32_t i;
659 for (i=0;i<st->den_rate;i++)
660 {
661 spx_uint32_t j;
662 for (j=0;j<st->filt_len;j++)
663 {
664 st->sinc_table[i*st->filt_len+j] = sinc(st->cutoff,((j-(spx_int32_t)st->filt_len/2+1)-((float)i)/st->den_rate), st->filt_len, quality_map[st->quality].window_func);
665 }
666 }
667 #ifdef FIXED_POINT
668 st->resampler_ptr = resampler_basic_direct_single;
669 #else
670 if (st->quality>8)
671 st->resampler_ptr = resampler_basic_direct_double;
672 else
673 st->resampler_ptr = resampler_basic_direct_single;
674 #endif
675 /*fprintf (stderr, "resampler uses direct sinc table and normalised cutoff %f\n", cutoff);*/
676 } else {
677 spx_int32_t i;
678 for (i=-4;i<(spx_int32_t)(st->oversample*st->filt_len+4);i++)
679 st->sinc_table[i+4] = sinc(st->cutoff,(i/(float)st->oversample - st->filt_len/2), st->filt_len, quality_map[st->quality].window_func);
680 #ifdef FIXED_POINT
681 st->resampler_ptr = resampler_basic_interpolate_single;
682 #else
683 if (st->quality>8)
684 st->resampler_ptr = resampler_basic_interpolate_double;
685 else
686 st->resampler_ptr = resampler_basic_interpolate_single;
687 #endif
688 /*fprintf (stderr, "resampler uses interpolated sinc table and normalised cutoff %f\n", cutoff);*/
689 }
690
691 /* Here's the place where we update the filter memory to take into account
692 the change in filter length. It's probably the messiest part of the code
693 due to handling of lots of corner cases. */
694
695 /* Adding buffer_size to filt_len won't overflow here because filt_len
696 could be multiplied by sizeof(spx_word16_t) above. */
697 min_alloc_size = st->filt_len-1 + st->buffer_size;
698 if (min_alloc_size > st->mem_alloc_size)
699 {
700 spx_word16_t *mem;
701 if (INT_MAX/(spx_uint32_t)sizeof(spx_word16_t)/st->nb_channels < min_alloc_size)
702 goto fail;
703 else if (!(mem = (spx_word16_t*)speex_realloc(st->mem, st->nb_channels*min_alloc_size * sizeof(*mem))))
704 goto fail;
705
706 st->mem = mem;
707 st->mem_alloc_size = min_alloc_size;
708 }
709 if (!st->started)
710 {
711 spx_uint32_t i;
712 for (i=0;i<st->nb_channels*st->mem_alloc_size;i++)
713 st->mem[i] = 0;
714 /*speex_warning("reinit filter");*/
715 } else if (st->filt_len > old_length)
716 {
717 spx_uint32_t i;
718 /* Increase the filter length */
719 /*speex_warning("increase filter size");*/
720 for (i=st->nb_channels;i--;)
721 {
722 spx_uint32_t j;
723 spx_uint32_t olen = old_length;
724 /*if (st->magic_samples[i])*/
725 {
726 /* Try and remove the magic samples as if nothing had happened */
727
728 /* FIXME: This is wrong but for now we need it to avoid going over the array bounds */
729 olen = old_length + 2*st->magic_samples[i];
730 for (j=old_length-1+st->magic_samples[i];j--;)
731 st->mem[i*st->mem_alloc_size+j+st->magic_samples[i]] = st->mem[i*old_alloc_size+j];
732 for (j=0;j<st->magic_samples[i];j++)
733 st->mem[i*st->mem_alloc_size+j] = 0;
734 st->magic_samples[i] = 0;
735 }
736 if (st->filt_len > olen)
737 {
738 /* If the new filter length is still bigger than the "augmented" length */
739 /* Copy data going backward */
740 for (j=0;j<olen-1;j++)
741 st->mem[i*st->mem_alloc_size+(st->filt_len-2-j)] = st->mem[i*st->mem_alloc_size+(olen-2-j)];
742 /* Then put zeros for lack of anything better */
743 for (;j<st->filt_len-1;j++)
744 st->mem[i*st->mem_alloc_size+(st->filt_len-2-j)] = 0;
745 /* Adjust last_sample */
746 st->last_sample[i] += (st->filt_len - olen)/2;
747 } else {
748 /* Put back some of the magic! */
749 st->magic_samples[i] = (olen - st->filt_len)/2;
750 for (j=0;j<st->filt_len-1+st->magic_samples[i];j++)
751 st->mem[i*st->mem_alloc_size+j] = st->mem[i*st->mem_alloc_size+j+st->magic_samples[i]];
752 }
753 }
754 } else if (st->filt_len < old_length)
755 {
756 spx_uint32_t i;
757 /* Reduce filter length, this a bit tricky. We need to store some of the memory as "magic"
758 samples so they can be used directly as input the next time(s) */
759 for (i=0;i<st->nb_channels;i++)
760 {
761 spx_uint32_t j;
762 spx_uint32_t old_magic = st->magic_samples[i];
763 st->magic_samples[i] = (old_length - st->filt_len)/2;
764 /* We must copy some of the memory that's no longer used */
765 /* Copy data going backward */
766 for (j=0;j<st->filt_len-1+st->magic_samples[i]+old_magic;j++)
767 st->mem[i*st->mem_alloc_size+j] = st->mem[i*st->mem_alloc_size+j+st->magic_samples[i]];
768 st->magic_samples[i] += old_magic;
769 }
770 }
771 return RESAMPLER_ERR_SUCCESS;
772
773 fail:
774 st->resampler_ptr = resampler_basic_zero;
775 /* st->mem may still contain consumed input samples for the filter.
776 Restore filt_len so that filt_len - 1 still points to the position after
777 the last of these samples. */
778 st->filt_len = old_length;
779 return RESAMPLER_ERR_ALLOC_FAILED;
780 }
781
speex_resampler_init(spx_uint32_t nb_channels,spx_uint32_t in_rate,spx_uint32_t out_rate,int quality,int * err)782 EXPORT SpeexResamplerState *speex_resampler_init(spx_uint32_t nb_channels, spx_uint32_t in_rate, spx_uint32_t out_rate, int quality, int *err)
783 {
784 return speex_resampler_init_frac(nb_channels, in_rate, out_rate, in_rate, out_rate, quality, err);
785 }
786
speex_resampler_init_frac(spx_uint32_t nb_channels,spx_uint32_t ratio_num,spx_uint32_t ratio_den,spx_uint32_t in_rate,spx_uint32_t out_rate,int quality,int * err)787 EXPORT SpeexResamplerState *speex_resampler_init_frac(spx_uint32_t nb_channels, spx_uint32_t ratio_num, spx_uint32_t ratio_den, spx_uint32_t in_rate, spx_uint32_t out_rate, int quality, int *err)
788 {
789 spx_uint32_t i;
790 SpeexResamplerState *st;
791 int filter_err;
792
793 if (nb_channels == 0 || ratio_num == 0 || ratio_den == 0 || quality > 10 || quality < 0)
794 {
795 if (err)
796 *err = RESAMPLER_ERR_INVALID_ARG;
797 return NULL;
798 }
799 st = (SpeexResamplerState *)speex_alloc(sizeof(SpeexResamplerState));
800 st->initialised = 0;
801 st->started = 0;
802 st->in_rate = 0;
803 st->out_rate = 0;
804 st->num_rate = 0;
805 st->den_rate = 0;
806 st->quality = -1;
807 st->sinc_table_length = 0;
808 st->mem_alloc_size = 0;
809 st->filt_len = 0;
810 st->mem = 0;
811 st->resampler_ptr = 0;
812
813 st->cutoff = 1.f;
814 st->nb_channels = nb_channels;
815 st->in_stride = 1;
816 st->out_stride = 1;
817
818 st->buffer_size = 160;
819
820 /* Per channel data */
821 st->last_sample = (spx_int32_t*)speex_alloc(nb_channels*sizeof(spx_int32_t));
822 st->magic_samples = (spx_uint32_t*)speex_alloc(nb_channels*sizeof(spx_uint32_t));
823 st->samp_frac_num = (spx_uint32_t*)speex_alloc(nb_channels*sizeof(spx_uint32_t));
824 for (i=0;i<nb_channels;i++)
825 {
826 st->last_sample[i] = 0;
827 st->magic_samples[i] = 0;
828 st->samp_frac_num[i] = 0;
829 }
830
831 speex_resampler_set_quality(st, quality);
832 speex_resampler_set_rate_frac(st, ratio_num, ratio_den, in_rate, out_rate);
833
834 filter_err = update_filter(st);
835 if (filter_err == RESAMPLER_ERR_SUCCESS)
836 {
837 st->initialised = 1;
838 } else {
839 speex_resampler_destroy(st);
840 st = NULL;
841 }
842 if (err)
843 *err = filter_err;
844
845 return st;
846 }
847
speex_resampler_destroy(SpeexResamplerState * st)848 EXPORT void speex_resampler_destroy(SpeexResamplerState *st)
849 {
850 speex_free(st->mem);
851 speex_free(st->sinc_table);
852 speex_free(st->last_sample);
853 speex_free(st->magic_samples);
854 speex_free(st->samp_frac_num);
855 speex_free(st);
856 }
857
speex_resampler_process_native(SpeexResamplerState * st,spx_uint32_t channel_index,spx_uint32_t * in_len,spx_word16_t * out,spx_uint32_t * out_len)858 static int speex_resampler_process_native(SpeexResamplerState *st, spx_uint32_t channel_index, spx_uint32_t *in_len, spx_word16_t *out, spx_uint32_t *out_len)
859 {
860 int j=0;
861 const int N = st->filt_len;
862 int out_sample = 0;
863 spx_word16_t *mem = st->mem + channel_index * st->mem_alloc_size;
864 spx_uint32_t ilen;
865
866 st->started = 1;
867
868 /* Call the right resampler through the function ptr */
869 out_sample = st->resampler_ptr(st, channel_index, mem, in_len, out, out_len);
870
871 if (st->last_sample[channel_index] < (spx_int32_t)*in_len)
872 *in_len = st->last_sample[channel_index];
873 *out_len = out_sample;
874 st->last_sample[channel_index] -= *in_len;
875
876 ilen = *in_len;
877
878 for(j=0;j<N-1;++j)
879 mem[j] = mem[j+ilen];
880
881 return RESAMPLER_ERR_SUCCESS;
882 }
883
speex_resampler_magic(SpeexResamplerState * st,spx_uint32_t channel_index,spx_word16_t ** out,spx_uint32_t out_len)884 static int speex_resampler_magic(SpeexResamplerState *st, spx_uint32_t channel_index, spx_word16_t **out, spx_uint32_t out_len) {
885 spx_uint32_t tmp_in_len = st->magic_samples[channel_index];
886 spx_word16_t *mem = st->mem + channel_index * st->mem_alloc_size;
887 const int N = st->filt_len;
888
889 speex_resampler_process_native(st, channel_index, &tmp_in_len, *out, &out_len);
890
891 st->magic_samples[channel_index] -= tmp_in_len;
892
893 /* If we couldn't process all "magic" input samples, save the rest for next time */
894 if (st->magic_samples[channel_index])
895 {
896 spx_uint32_t i;
897 for (i=0;i<st->magic_samples[channel_index];i++)
898 mem[N-1+i]=mem[N-1+i+tmp_in_len];
899 }
900 *out += out_len*st->out_stride;
901 return out_len;
902 }
903
904 #ifdef FIXED_POINT
speex_resampler_process_int(SpeexResamplerState * st,spx_uint32_t channel_index,const spx_int16_t * in,spx_uint32_t * in_len,spx_int16_t * out,spx_uint32_t * out_len)905 EXPORT int speex_resampler_process_int(SpeexResamplerState *st, spx_uint32_t channel_index, const spx_int16_t *in, spx_uint32_t *in_len, spx_int16_t *out, spx_uint32_t *out_len)
906 #else
907 EXPORT int speex_resampler_process_float(SpeexResamplerState *st, spx_uint32_t channel_index, const float *in, spx_uint32_t *in_len, float *out, spx_uint32_t *out_len)
908 #endif
909 {
910 spx_uint32_t j;
911 spx_uint32_t ilen = *in_len;
912 spx_uint32_t olen = *out_len;
913 spx_word16_t *x = st->mem + channel_index * st->mem_alloc_size;
914 const int filt_offs = st->filt_len - 1;
915 const spx_uint32_t xlen = st->mem_alloc_size - filt_offs;
916 const int istride = st->in_stride;
917
918 if (st->magic_samples[channel_index])
919 olen -= speex_resampler_magic(st, channel_index, &out, olen);
920 if (! st->magic_samples[channel_index]) {
921 while (ilen && olen) {
922 spx_uint32_t ichunk = (ilen > xlen) ? xlen : ilen;
923 spx_uint32_t ochunk = olen;
924
925 if (in) {
926 for(j=0;j<ichunk;++j)
927 x[j+filt_offs]=in[j*istride];
928 } else {
929 for(j=0;j<ichunk;++j)
930 x[j+filt_offs]=0;
931 }
932 speex_resampler_process_native(st, channel_index, &ichunk, out, &ochunk);
933 ilen -= ichunk;
934 olen -= ochunk;
935 out += ochunk * st->out_stride;
936 if (in)
937 in += ichunk * istride;
938 }
939 }
940 *in_len -= ilen;
941 *out_len -= olen;
942 return st->resampler_ptr == resampler_basic_zero ? RESAMPLER_ERR_ALLOC_FAILED : RESAMPLER_ERR_SUCCESS;
943 }
944
945 #ifdef FIXED_POINT
speex_resampler_process_float(SpeexResamplerState * st,spx_uint32_t channel_index,const float * in,spx_uint32_t * in_len,float * out,spx_uint32_t * out_len)946 EXPORT int speex_resampler_process_float(SpeexResamplerState *st, spx_uint32_t channel_index, const float *in, spx_uint32_t *in_len, float *out, spx_uint32_t *out_len)
947 #else
948 EXPORT int speex_resampler_process_int(SpeexResamplerState *st, spx_uint32_t channel_index, const spx_int16_t *in, spx_uint32_t *in_len, spx_int16_t *out, spx_uint32_t *out_len)
949 #endif
950 {
951 spx_uint32_t j;
952 const int istride_save = st->in_stride;
953 const int ostride_save = st->out_stride;
954 spx_uint32_t ilen = *in_len;
955 spx_uint32_t olen = *out_len;
956 spx_word16_t *x = st->mem + channel_index * st->mem_alloc_size;
957 const spx_uint32_t xlen = st->mem_alloc_size - (st->filt_len - 1);
958 #ifdef VAR_ARRAYS
959 const unsigned int ylen = (olen < FIXED_STACK_ALLOC) ? olen : FIXED_STACK_ALLOC;
960 VARDECL(spx_word16_t *ystack);
961 ALLOC(ystack, ylen, spx_word16_t);
962 #else
963 const unsigned int ylen = FIXED_STACK_ALLOC;
964 spx_word16_t ystack[FIXED_STACK_ALLOC];
965 #endif
966
967 st->out_stride = 1;
968
969 while (ilen && olen) {
970 spx_word16_t *y = ystack;
971 spx_uint32_t ichunk = (ilen > xlen) ? xlen : ilen;
972 spx_uint32_t ochunk = (olen > ylen) ? ylen : olen;
973 spx_uint32_t omagic = 0;
974
975 if (st->magic_samples[channel_index]) {
976 omagic = speex_resampler_magic(st, channel_index, &y, ochunk);
977 ochunk -= omagic;
978 olen -= omagic;
979 }
980 if (! st->magic_samples[channel_index]) {
981 if (in) {
982 for(j=0;j<ichunk;++j)
983 #ifdef FIXED_POINT
984 x[j+st->filt_len-1]=WORD2INT(in[j*istride_save]);
985 #else
986 x[j+st->filt_len-1]=in[j*istride_save];
987 #endif
988 } else {
989 for(j=0;j<ichunk;++j)
990 x[j+st->filt_len-1]=0;
991 }
992
993 speex_resampler_process_native(st, channel_index, &ichunk, y, &ochunk);
994 } else {
995 ichunk = 0;
996 ochunk = 0;
997 }
998
999 for (j=0;j<ochunk+omagic;++j)
1000 #ifdef FIXED_POINT
1001 out[j*ostride_save] = ystack[j];
1002 #else
1003 out[j*ostride_save] = WORD2INT(ystack[j]);
1004 #endif
1005
1006 ilen -= ichunk;
1007 olen -= ochunk;
1008 out += (ochunk+omagic) * ostride_save;
1009 if (in)
1010 in += ichunk * istride_save;
1011 }
1012 st->out_stride = ostride_save;
1013 *in_len -= ilen;
1014 *out_len -= olen;
1015
1016 return st->resampler_ptr == resampler_basic_zero ? RESAMPLER_ERR_ALLOC_FAILED : RESAMPLER_ERR_SUCCESS;
1017 }
1018
speex_resampler_process_interleaved_float(SpeexResamplerState * st,const float * in,spx_uint32_t * in_len,float * out,spx_uint32_t * out_len)1019 EXPORT int speex_resampler_process_interleaved_float(SpeexResamplerState *st, const float *in, spx_uint32_t *in_len, float *out, spx_uint32_t *out_len)
1020 {
1021 spx_uint32_t i;
1022 int istride_save, ostride_save;
1023 spx_uint32_t bak_out_len = *out_len;
1024 spx_uint32_t bak_in_len = *in_len;
1025 istride_save = st->in_stride;
1026 ostride_save = st->out_stride;
1027 st->in_stride = st->out_stride = st->nb_channels;
1028 for (i=0;i<st->nb_channels;i++)
1029 {
1030 *out_len = bak_out_len;
1031 *in_len = bak_in_len;
1032 if (in != NULL)
1033 speex_resampler_process_float(st, i, in+i, in_len, out+i, out_len);
1034 else
1035 speex_resampler_process_float(st, i, NULL, in_len, out+i, out_len);
1036 }
1037 st->in_stride = istride_save;
1038 st->out_stride = ostride_save;
1039 return st->resampler_ptr == resampler_basic_zero ? RESAMPLER_ERR_ALLOC_FAILED : RESAMPLER_ERR_SUCCESS;
1040 }
1041
speex_resampler_process_interleaved_int(SpeexResamplerState * st,const spx_int16_t * in,spx_uint32_t * in_len,spx_int16_t * out,spx_uint32_t * out_len)1042 EXPORT int speex_resampler_process_interleaved_int(SpeexResamplerState *st, const spx_int16_t *in, spx_uint32_t *in_len, spx_int16_t *out, spx_uint32_t *out_len)
1043 {
1044 spx_uint32_t i;
1045 int istride_save, ostride_save;
1046 spx_uint32_t bak_out_len = *out_len;
1047 spx_uint32_t bak_in_len = *in_len;
1048 istride_save = st->in_stride;
1049 ostride_save = st->out_stride;
1050 st->in_stride = st->out_stride = st->nb_channels;
1051 for (i=0;i<st->nb_channels;i++)
1052 {
1053 *out_len = bak_out_len;
1054 *in_len = bak_in_len;
1055 if (in != NULL)
1056 speex_resampler_process_int(st, i, in+i, in_len, out+i, out_len);
1057 else
1058 speex_resampler_process_int(st, i, NULL, in_len, out+i, out_len);
1059 }
1060 st->in_stride = istride_save;
1061 st->out_stride = ostride_save;
1062 return st->resampler_ptr == resampler_basic_zero ? RESAMPLER_ERR_ALLOC_FAILED : RESAMPLER_ERR_SUCCESS;
1063 }
1064
speex_resampler_set_rate(SpeexResamplerState * st,spx_uint32_t in_rate,spx_uint32_t out_rate)1065 EXPORT int speex_resampler_set_rate(SpeexResamplerState *st, spx_uint32_t in_rate, spx_uint32_t out_rate)
1066 {
1067 return speex_resampler_set_rate_frac(st, in_rate, out_rate, in_rate, out_rate);
1068 }
1069
speex_resampler_get_rate(SpeexResamplerState * st,spx_uint32_t * in_rate,spx_uint32_t * out_rate)1070 EXPORT void speex_resampler_get_rate(SpeexResamplerState *st, spx_uint32_t *in_rate, spx_uint32_t *out_rate)
1071 {
1072 *in_rate = st->in_rate;
1073 *out_rate = st->out_rate;
1074 }
1075
speex_resampler_set_rate_frac(SpeexResamplerState * st,spx_uint32_t ratio_num,spx_uint32_t ratio_den,spx_uint32_t in_rate,spx_uint32_t out_rate)1076 EXPORT int speex_resampler_set_rate_frac(SpeexResamplerState *st, spx_uint32_t ratio_num, spx_uint32_t ratio_den, spx_uint32_t in_rate, spx_uint32_t out_rate)
1077 {
1078 spx_uint32_t fact;
1079 spx_uint32_t old_den;
1080 spx_uint32_t i;
1081
1082 if (ratio_num == 0 || ratio_den == 0)
1083 return RESAMPLER_ERR_INVALID_ARG;
1084
1085 if (st->in_rate == in_rate && st->out_rate == out_rate && st->num_rate == ratio_num && st->den_rate == ratio_den)
1086 return RESAMPLER_ERR_SUCCESS;
1087
1088 old_den = st->den_rate;
1089 st->in_rate = in_rate;
1090 st->out_rate = out_rate;
1091 st->num_rate = ratio_num;
1092 st->den_rate = ratio_den;
1093 /* FIXME: This is terribly inefficient, but who cares (at least for now)? */
1094 for (fact=2;fact<=IMIN(st->num_rate, st->den_rate);fact++)
1095 {
1096 while ((st->num_rate % fact == 0) && (st->den_rate % fact == 0))
1097 {
1098 st->num_rate /= fact;
1099 st->den_rate /= fact;
1100 }
1101 }
1102
1103 if (old_den > 0)
1104 {
1105 for (i=0;i<st->nb_channels;i++)
1106 {
1107 st->samp_frac_num[i]=st->samp_frac_num[i]*st->den_rate/old_den;
1108 /* Safety net */
1109 if (st->samp_frac_num[i] >= st->den_rate)
1110 st->samp_frac_num[i] = st->den_rate-1;
1111 }
1112 }
1113
1114 if (st->initialised)
1115 return update_filter(st);
1116 return RESAMPLER_ERR_SUCCESS;
1117 }
1118
speex_resampler_get_ratio(SpeexResamplerState * st,spx_uint32_t * ratio_num,spx_uint32_t * ratio_den)1119 EXPORT void speex_resampler_get_ratio(SpeexResamplerState *st, spx_uint32_t *ratio_num, spx_uint32_t *ratio_den)
1120 {
1121 *ratio_num = st->num_rate;
1122 *ratio_den = st->den_rate;
1123 }
1124
speex_resampler_set_quality(SpeexResamplerState * st,int quality)1125 EXPORT int speex_resampler_set_quality(SpeexResamplerState *st, int quality)
1126 {
1127 if (quality > 10 || quality < 0)
1128 return RESAMPLER_ERR_INVALID_ARG;
1129 if (st->quality == quality)
1130 return RESAMPLER_ERR_SUCCESS;
1131 st->quality = quality;
1132 if (st->initialised)
1133 return update_filter(st);
1134 return RESAMPLER_ERR_SUCCESS;
1135 }
1136
speex_resampler_get_quality(SpeexResamplerState * st,int * quality)1137 EXPORT void speex_resampler_get_quality(SpeexResamplerState *st, int *quality)
1138 {
1139 *quality = st->quality;
1140 }
1141
speex_resampler_set_input_stride(SpeexResamplerState * st,spx_uint32_t stride)1142 EXPORT void speex_resampler_set_input_stride(SpeexResamplerState *st, spx_uint32_t stride)
1143 {
1144 st->in_stride = stride;
1145 }
1146
speex_resampler_get_input_stride(SpeexResamplerState * st,spx_uint32_t * stride)1147 EXPORT void speex_resampler_get_input_stride(SpeexResamplerState *st, spx_uint32_t *stride)
1148 {
1149 *stride = st->in_stride;
1150 }
1151
speex_resampler_set_output_stride(SpeexResamplerState * st,spx_uint32_t stride)1152 EXPORT void speex_resampler_set_output_stride(SpeexResamplerState *st, spx_uint32_t stride)
1153 {
1154 st->out_stride = stride;
1155 }
1156
speex_resampler_get_output_stride(SpeexResamplerState * st,spx_uint32_t * stride)1157 EXPORT void speex_resampler_get_output_stride(SpeexResamplerState *st, spx_uint32_t *stride)
1158 {
1159 *stride = st->out_stride;
1160 }
1161
speex_resampler_get_input_latency(SpeexResamplerState * st)1162 EXPORT int speex_resampler_get_input_latency(SpeexResamplerState *st)
1163 {
1164 return st->filt_len / 2;
1165 }
1166
speex_resampler_get_output_latency(SpeexResamplerState * st)1167 EXPORT int speex_resampler_get_output_latency(SpeexResamplerState *st)
1168 {
1169 return ((st->filt_len / 2) * st->den_rate + (st->num_rate >> 1)) / st->num_rate;
1170 }
1171
speex_resampler_skip_zeros(SpeexResamplerState * st)1172 EXPORT int speex_resampler_skip_zeros(SpeexResamplerState *st)
1173 {
1174 spx_uint32_t i;
1175 for (i=0;i<st->nb_channels;i++)
1176 st->last_sample[i] = st->filt_len/2;
1177 return RESAMPLER_ERR_SUCCESS;
1178 }
1179
speex_resampler_reset_mem(SpeexResamplerState * st)1180 EXPORT int speex_resampler_reset_mem(SpeexResamplerState *st)
1181 {
1182 spx_uint32_t i;
1183 for (i=0;i<st->nb_channels;i++)
1184 {
1185 st->last_sample[i] = 0;
1186 st->magic_samples[i] = 0;
1187 st->samp_frac_num[i] = 0;
1188 }
1189 for (i=0;i<st->nb_channels*(st->filt_len-1);i++)
1190 st->mem[i] = 0;
1191 return RESAMPLER_ERR_SUCCESS;
1192 }
1193
speex_resampler_strerror(int err)1194 EXPORT const char *speex_resampler_strerror(int err)
1195 {
1196 switch (err)
1197 {
1198 case RESAMPLER_ERR_SUCCESS:
1199 return "Success.";
1200 case RESAMPLER_ERR_ALLOC_FAILED:
1201 return "Memory allocation failed.";
1202 case RESAMPLER_ERR_BAD_STATE:
1203 return "Bad resampler state.";
1204 case RESAMPLER_ERR_INVALID_ARG:
1205 return "Invalid argument.";
1206 case RESAMPLER_ERR_PTR_OVERLAP:
1207 return "Input and output buffers overlap.";
1208 default:
1209 return "Unknown error. Bad error code or strange version mismatch.";
1210 }
1211 }
1212