tor-browser

The Tor Browser
git clone https://git.dasho.dev/tor-browser.git
Log | Files | Refs | README | LICENSE

resample.c (45045B)


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