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