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