LCOV - code coverage report
Current view: top level - media/libspeex_resampler/src - resample.c (source / functions) Hit Total Coverage
Test: output.info Lines: 0 538 0.0 %
Date: 2017-07-14 16:53:18 Functions: 0 39 0.0 %
Legend: Lines: hit not hit

          Line data    Source code
       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://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           0 : static void *speex_alloc (int size) {return calloc(size,1);}
      69           0 : static void *speex_realloc (void *ptr, int size) {return realloc(ptr, size);}
      70           0 : static void speex_free (void *ptr) {free(ptr);}
      71             : #include "speex_resampler.h"
      72             : #include "arch.h"
      73             : #else /* OUTSIDE_SPEEX */
      74             : 
      75             : #include "speex/speex_resampler.h"
      76             : #include "arch.h"
      77             : #include "os_support.h"
      78             : #endif /* OUTSIDE_SPEEX */
      79             : 
      80             : #include "stack_alloc.h"
      81             : #include <math.h>
      82             : #include <limits.h>
      83             : 
      84             : #ifndef M_PI
      85             : #define M_PI 3.14159265358979323846
      86             : #endif
      87             : 
      88             : #define IMAX(a,b) ((a) > (b) ? (a) : (b))
      89             : #define IMIN(a,b) ((a) < (b) ? (a) : (b))
      90             : 
      91             : #ifndef NULL
      92             : #define NULL 0
      93             : #endif
      94             : 
      95             : #ifndef UINT32_MAX
      96             : #define UINT32_MAX 4294967296U
      97             : #endif
      98             : 
      99             : #include "simd_detect.h"
     100             : 
     101             : /* Numer of elements to allocate on the stack */
     102             : #ifdef VAR_ARRAYS
     103             : #define FIXED_STACK_ALLOC 8192
     104             : #else
     105             : #define FIXED_STACK_ALLOC 1024
     106             : #endif
     107             : 
     108             : typedef int (*resampler_basic_func)(SpeexResamplerState *, spx_uint32_t , const spx_word16_t *, spx_uint32_t *, spx_word16_t *, spx_uint32_t *);
     109             : 
     110             : struct SpeexResamplerState_ {
     111             :    spx_uint32_t in_rate;
     112             :    spx_uint32_t out_rate;
     113             :    spx_uint32_t num_rate;
     114             :    spx_uint32_t den_rate;
     115             : 
     116             :    int    quality;
     117             :    spx_uint32_t nb_channels;
     118             :    spx_uint32_t filt_len;
     119             :    spx_uint32_t mem_alloc_size;
     120             :    spx_uint32_t buffer_size;
     121             :    int          int_advance;
     122             :    int          frac_advance;
     123             :    float  cutoff;
     124             :    spx_uint32_t oversample;
     125             :    int          initialised;
     126             :    int          started;
     127             : 
     128             :    /* These are per-channel */
     129             :    spx_int32_t  *last_sample;
     130             :    spx_uint32_t *samp_frac_num;
     131             :    spx_uint32_t *magic_samples;
     132             : 
     133             :    spx_word16_t *mem;
     134             :    spx_word16_t *sinc_table;
     135             :    spx_uint32_t sinc_table_length;
     136             :    resampler_basic_func resampler_ptr;
     137             : 
     138             :    int    in_stride;
     139             :    int    out_stride;
     140             : } ;
     141             : 
     142             : static const double kaiser12_table[68] = {
     143             :    0.99859849, 1.00000000, 0.99859849, 0.99440475, 0.98745105, 0.97779076,
     144             :    0.96549770, 0.95066529, 0.93340547, 0.91384741, 0.89213598, 0.86843014,
     145             :    0.84290116, 0.81573067, 0.78710866, 0.75723148, 0.72629970, 0.69451601,
     146             :    0.66208321, 0.62920216, 0.59606986, 0.56287762, 0.52980938, 0.49704014,
     147             :    0.46473455, 0.43304576, 0.40211431, 0.37206735, 0.34301800, 0.31506490,
     148             :    0.28829195, 0.26276832, 0.23854851, 0.21567274, 0.19416736, 0.17404546,
     149             :    0.15530766, 0.13794294, 0.12192957, 0.10723616, 0.09382272, 0.08164178,
     150             :    0.07063950, 0.06075685, 0.05193064, 0.04409466, 0.03718069, 0.03111947,
     151             :    0.02584161, 0.02127838, 0.01736250, 0.01402878, 0.01121463, 0.00886058,
     152             :    0.00691064, 0.00531256, 0.00401805, 0.00298291, 0.00216702, 0.00153438,
     153             :    0.00105297, 0.00069463, 0.00043489, 0.00025272, 0.00013031, 0.0000527734,
     154             :    0.00001000, 0.00000000};
     155             : /*
     156             : static const double kaiser12_table[36] = {
     157             :    0.99440475, 1.00000000, 0.99440475, 0.97779076, 0.95066529, 0.91384741,
     158             :    0.86843014, 0.81573067, 0.75723148, 0.69451601, 0.62920216, 0.56287762,
     159             :    0.49704014, 0.43304576, 0.37206735, 0.31506490, 0.26276832, 0.21567274,
     160             :    0.17404546, 0.13794294, 0.10723616, 0.08164178, 0.06075685, 0.04409466,
     161             :    0.03111947, 0.02127838, 0.01402878, 0.00886058, 0.00531256, 0.00298291,
     162             :    0.00153438, 0.00069463, 0.00025272, 0.0000527734, 0.00000500, 0.00000000};
     163             : */
     164             : static const double kaiser10_table[36] = {
     165             :    0.99537781, 1.00000000, 0.99537781, 0.98162644, 0.95908712, 0.92831446,
     166             :    0.89005583, 0.84522401, 0.79486424, 0.74011713, 0.68217934, 0.62226347,
     167             :    0.56155915, 0.50119680, 0.44221549, 0.38553619, 0.33194107, 0.28205962,
     168             :    0.23636152, 0.19515633, 0.15859932, 0.12670280, 0.09935205, 0.07632451,
     169             :    0.05731132, 0.04193980, 0.02979584, 0.02044510, 0.01345224, 0.00839739,
     170             :    0.00488951, 0.00257636, 0.00115101, 0.00035515, 0.00000000, 0.00000000};
     171             : 
     172             : static const double kaiser8_table[36] = {
     173             :    0.99635258, 1.00000000, 0.99635258, 0.98548012, 0.96759014, 0.94302200,
     174             :    0.91223751, 0.87580811, 0.83439927, 0.78875245, 0.73966538, 0.68797126,
     175             :    0.63451750, 0.58014482, 0.52566725, 0.47185369, 0.41941150, 0.36897272,
     176             :    0.32108304, 0.27619388, 0.23465776, 0.19672670, 0.16255380, 0.13219758,
     177             :    0.10562887, 0.08273982, 0.06335451, 0.04724088, 0.03412321, 0.02369490,
     178             :    0.01563093, 0.00959968, 0.00527363, 0.00233883, 0.00050000, 0.00000000};
     179             : 
     180             : static const double kaiser6_table[36] = {
     181             :    0.99733006, 1.00000000, 0.99733006, 0.98935595, 0.97618418, 0.95799003,
     182             :    0.93501423, 0.90755855, 0.87598009, 0.84068475, 0.80211977, 0.76076565,
     183             :    0.71712752, 0.67172623, 0.62508937, 0.57774224, 0.53019925, 0.48295561,
     184             :    0.43647969, 0.39120616, 0.34752997, 0.30580127, 0.26632152, 0.22934058,
     185             :    0.19505503, 0.16360756, 0.13508755, 0.10953262, 0.08693120, 0.06722600,
     186             :    0.05031820, 0.03607231, 0.02432151, 0.01487334, 0.00752000, 0.00000000};
     187             : 
     188             : struct FuncDef {
     189             :    const double *table;
     190             :    int oversample;
     191             : };
     192             : 
     193             : static const struct FuncDef _KAISER12 = {kaiser12_table, 64};
     194             : #define KAISER12 (&_KAISER12)
     195             : /*static struct FuncDef _KAISER12 = {kaiser12_table, 32};
     196             : #define KAISER12 (&_KAISER12)*/
     197             : static const struct FuncDef _KAISER10 = {kaiser10_table, 32};
     198             : #define KAISER10 (&_KAISER10)
     199             : static const struct FuncDef _KAISER8 = {kaiser8_table, 32};
     200             : #define KAISER8 (&_KAISER8)
     201             : static const struct FuncDef _KAISER6 = {kaiser6_table, 32};
     202             : #define KAISER6 (&_KAISER6)
     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           0 : static double compute_func(float x, const struct FuncDef *func)
     237             : {
     238             :    float y, frac;
     239             :    double interp[4];
     240             :    int ind;
     241           0 :    y = x*func->oversample;
     242           0 :    ind = (int)floor(y);
     243           0 :    frac = (y-ind);
     244             :    /* CSE with handle the repeated powers */
     245           0 :    interp[3] =  -0.1666666667*frac + 0.1666666667*(frac*frac*frac);
     246           0 :    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           0 :    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           0 :    interp[1] = 1.f-interp[3]-interp[2]-interp[0];
     251             : 
     252             :    /*sum = frac*accum[1] + (1-frac)*accum[2];*/
     253           0 :    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           0 : static spx_word16_t sinc(float cutoff, float x, int N, const struct FuncDef *window_func)
     285             : {
     286             :    /*fprintf (stderr, "%f ", x);*/
     287           0 :    float xx = x * cutoff;
     288           0 :    if (fabs(x)<1e-6)
     289           0 :       return cutoff;
     290           0 :    else if (fabs(x) > .5*N)
     291           0 :       return 0;
     292             :    /*FIXME: Can it really be any slower than this? */
     293           0 :    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           0 : 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           0 :    interp[0] =  -0.16667f*frac + 0.16667f*frac*frac*frac;
     319           0 :    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           0 :    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           0 :    interp[2] = 1.-interp[0]-interp[1]-interp[3];
     324           0 : }
     325             : #endif
     326             : 
     327           0 : 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           0 :    const int N = st->filt_len;
     330           0 :    int out_sample = 0;
     331           0 :    int last_sample = st->last_sample[channel_index];
     332           0 :    spx_uint32_t samp_frac_num = st->samp_frac_num[channel_index];
     333           0 :    const spx_word16_t *sinc_table = st->sinc_table;
     334           0 :    const int out_stride = st->out_stride;
     335           0 :    const int int_advance = st->int_advance;
     336           0 :    const int frac_advance = st->frac_advance;
     337           0 :    const spx_uint32_t den_rate = st->den_rate;
     338             :    spx_word32_t sum;
     339             : 
     340           0 :    while (!(last_sample >= (spx_int32_t)*in_len || out_sample >= (spx_int32_t)*out_len))
     341             :    {
     342           0 :       const spx_word16_t *sinct = & sinc_table[samp_frac_num*N];
     343           0 :       const spx_word16_t *iptr = & in[last_sample];
     344             : 
     345             : #ifdef OVERRIDE_INNER_PRODUCT_SINGLE
     346           0 :       if (!moz_speex_have_single_simd()) {
     347             : #endif
     348             :       int j;
     349           0 :       sum = 0;
     350           0 :       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           0 :       sum = SATURATE32PSHR(sum, 15, 32767);
     365             : #ifdef OVERRIDE_INNER_PRODUCT_SINGLE
     366             :       } else {
     367           0 :       sum = inner_product_single(sinct, iptr, N);
     368             :       }
     369             : #endif
     370             : 
     371           0 :       out[out_stride * out_sample++] = sum;
     372           0 :       last_sample += int_advance;
     373           0 :       samp_frac_num += frac_advance;
     374           0 :       if (samp_frac_num >= den_rate)
     375             :       {
     376           0 :          samp_frac_num -= den_rate;
     377           0 :          last_sample++;
     378             :       }
     379             :    }
     380             : 
     381           0 :    st->last_sample[channel_index] = last_sample;
     382           0 :    st->samp_frac_num[channel_index] = samp_frac_num;
     383           0 :    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           0 : 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           0 :    const int N = st->filt_len;
     392           0 :    int out_sample = 0;
     393           0 :    int last_sample = st->last_sample[channel_index];
     394           0 :    spx_uint32_t samp_frac_num = st->samp_frac_num[channel_index];
     395           0 :    const spx_word16_t *sinc_table = st->sinc_table;
     396           0 :    const int out_stride = st->out_stride;
     397           0 :    const int int_advance = st->int_advance;
     398           0 :    const int frac_advance = st->frac_advance;
     399           0 :    const spx_uint32_t den_rate = st->den_rate;
     400             :    double sum;
     401             : 
     402           0 :    while (!(last_sample >= (spx_int32_t)*in_len || out_sample >= (spx_int32_t)*out_len))
     403             :    {
     404           0 :       const spx_word16_t *sinct = & sinc_table[samp_frac_num*N];
     405           0 :       const spx_word16_t *iptr = & in[last_sample];
     406             : 
     407             : #ifdef OVERRIDE_INNER_PRODUCT_DOUBLE
     408           0 :       if(moz_speex_have_double_simd()) {
     409             : #endif
     410             :       int j;
     411           0 :       double accum[4] = {0,0,0,0};
     412             : 
     413           0 :       for(j=0;j<N;j+=4) {
     414           0 :         accum[0] += sinct[j]*iptr[j];
     415           0 :         accum[1] += sinct[j+1]*iptr[j+1];
     416           0 :         accum[2] += sinct[j+2]*iptr[j+2];
     417           0 :         accum[3] += sinct[j+3]*iptr[j+3];
     418             :       }
     419           0 :       sum = accum[0] + accum[1] + accum[2] + accum[3];
     420             : #ifdef OVERRIDE_INNER_PRODUCT_DOUBLE
     421             :       } else {
     422           0 :       sum = inner_product_double(sinct, iptr, N);
     423             :       }
     424             : #endif
     425             : 
     426           0 :       out[out_stride * out_sample++] = PSHR32(sum, 15);
     427           0 :       last_sample += int_advance;
     428           0 :       samp_frac_num += frac_advance;
     429           0 :       if (samp_frac_num >= den_rate)
     430             :       {
     431           0 :          samp_frac_num -= den_rate;
     432           0 :          last_sample++;
     433             :       }
     434             :    }
     435             : 
     436           0 :    st->last_sample[channel_index] = last_sample;
     437           0 :    st->samp_frac_num[channel_index] = samp_frac_num;
     438           0 :    return out_sample;
     439             : }
     440             : #endif
     441             : 
     442           0 : 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           0 :    const int N = st->filt_len;
     445           0 :    int out_sample = 0;
     446           0 :    int last_sample = st->last_sample[channel_index];
     447           0 :    spx_uint32_t samp_frac_num = st->samp_frac_num[channel_index];
     448           0 :    const int out_stride = st->out_stride;
     449           0 :    const int int_advance = st->int_advance;
     450           0 :    const int frac_advance = st->frac_advance;
     451           0 :    const spx_uint32_t den_rate = st->den_rate;
     452             :    spx_word32_t sum;
     453             : 
     454           0 :    while (!(last_sample >= (spx_int32_t)*in_len || out_sample >= (spx_int32_t)*out_len))
     455             :    {
     456           0 :       const spx_word16_t *iptr = & in[last_sample];
     457             : 
     458           0 :       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           0 :       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           0 :       if (!moz_speex_have_single_simd()) {
     469             : #endif
     470             :       int j;
     471           0 :       spx_word32_t accum[4] = {0,0,0,0};
     472             : 
     473           0 :       for(j=0;j<N;j++) {
     474           0 :         const spx_word16_t curr_in=iptr[j];
     475           0 :         accum[0] += MULT16_16(curr_in,st->sinc_table[4+(j+1)*st->oversample-offset-2]);
     476           0 :         accum[1] += MULT16_16(curr_in,st->sinc_table[4+(j+1)*st->oversample-offset-1]);
     477           0 :         accum[2] += MULT16_16(curr_in,st->sinc_table[4+(j+1)*st->oversample-offset]);
     478           0 :         accum[3] += MULT16_16(curr_in,st->sinc_table[4+(j+1)*st->oversample-offset+1]);
     479             :       }
     480             : 
     481           0 :       cubic_coef(frac, interp);
     482           0 :       sum = MULT16_32_Q15(interp[0],SHR32(accum[0], 1)) + MULT16_32_Q15(interp[1],SHR32(accum[1], 1)) + MULT16_32_Q15(interp[2],SHR32(accum[2], 1)) + MULT16_32_Q15(interp[3],SHR32(accum[3], 1));
     483           0 :       sum = SATURATE32PSHR(sum, 15, 32767);
     484             : #ifdef OVERRIDE_INTERPOLATE_PRODUCT_SINGLE
     485             :       } else {
     486           0 :       cubic_coef(frac, interp);
     487           0 :       sum = interpolate_product_single(iptr, st->sinc_table + st->oversample + 4 - offset - 2, N, st->oversample, interp);
     488             :       }
     489             : #endif
     490             : 
     491           0 :       out[out_stride * out_sample++] = sum;
     492           0 :       last_sample += int_advance;
     493           0 :       samp_frac_num += frac_advance;
     494           0 :       if (samp_frac_num >= den_rate)
     495             :       {
     496           0 :          samp_frac_num -= den_rate;
     497           0 :          last_sample++;
     498             :       }
     499             :    }
     500             : 
     501           0 :    st->last_sample[channel_index] = last_sample;
     502           0 :    st->samp_frac_num[channel_index] = samp_frac_num;
     503           0 :    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           0 : 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           0 :    const int N = st->filt_len;
     512           0 :    int out_sample = 0;
     513           0 :    int last_sample = st->last_sample[channel_index];
     514           0 :    spx_uint32_t samp_frac_num = st->samp_frac_num[channel_index];
     515           0 :    const int out_stride = st->out_stride;
     516           0 :    const int int_advance = st->int_advance;
     517           0 :    const int frac_advance = st->frac_advance;
     518           0 :    const spx_uint32_t den_rate = st->den_rate;
     519             :    spx_word32_t sum;
     520             : 
     521           0 :    while (!(last_sample >= (spx_int32_t)*in_len || out_sample >= (spx_int32_t)*out_len))
     522             :    {
     523           0 :       const spx_word16_t *iptr = & in[last_sample];
     524             : 
     525           0 :       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           0 :       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           0 :       if (!moz_speex_have_double_simd()) {
     536             : #endif
     537             :       int j;
     538           0 :       double accum[4] = {0,0,0,0};
     539             : 
     540           0 :       for(j=0;j<N;j++) {
     541           0 :         const double curr_in=iptr[j];
     542           0 :         accum[0] += MULT16_16(curr_in,st->sinc_table[4+(j+1)*st->oversample-offset-2]);
     543           0 :         accum[1] += MULT16_16(curr_in,st->sinc_table[4+(j+1)*st->oversample-offset-1]);
     544           0 :         accum[2] += MULT16_16(curr_in,st->sinc_table[4+(j+1)*st->oversample-offset]);
     545           0 :         accum[3] += MULT16_16(curr_in,st->sinc_table[4+(j+1)*st->oversample-offset+1]);
     546             :       }
     547             : 
     548           0 :       cubic_coef(frac, interp);
     549           0 :       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           0 :       cubic_coef(frac, interp);
     553           0 :       sum = interpolate_product_double(iptr, st->sinc_table + st->oversample + 4 - offset - 2, N, st->oversample, interp);
     554             :       }
     555             : #endif
     556             : 
     557           0 :       out[out_stride * out_sample++] = PSHR32(sum,15);
     558           0 :       last_sample += int_advance;
     559           0 :       samp_frac_num += frac_advance;
     560           0 :       if (samp_frac_num >= den_rate)
     561             :       {
     562           0 :          samp_frac_num -= den_rate;
     563           0 :          last_sample++;
     564             :       }
     565             :    }
     566             : 
     567           0 :    st->last_sample[channel_index] = last_sample;
     568           0 :    st->samp_frac_num[channel_index] = samp_frac_num;
     569           0 :    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           0 : 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           0 :    int out_sample = 0;
     580           0 :    int last_sample = st->last_sample[channel_index];
     581           0 :    spx_uint32_t samp_frac_num = st->samp_frac_num[channel_index];
     582           0 :    const int out_stride = st->out_stride;
     583           0 :    const int int_advance = st->int_advance;
     584           0 :    const int frac_advance = st->frac_advance;
     585           0 :    const spx_uint32_t den_rate = st->den_rate;
     586             : 
     587           0 :    while (!(last_sample >= (spx_int32_t)*in_len || out_sample >= (spx_int32_t)*out_len))
     588             :    {
     589           0 :       out[out_stride * out_sample++] = 0;
     590           0 :       last_sample += int_advance;
     591           0 :       samp_frac_num += frac_advance;
     592           0 :       if (samp_frac_num >= den_rate)
     593             :       {
     594           0 :          samp_frac_num -= den_rate;
     595           0 :          last_sample++;
     596             :       }
     597             :    }
     598             : 
     599           0 :    st->last_sample[channel_index] = last_sample;
     600           0 :    st->samp_frac_num[channel_index] = samp_frac_num;
     601           0 :    return out_sample;
     602             : }
     603             : 
     604           0 : static int _muldiv(spx_uint32_t *result, spx_uint32_t value, spx_uint32_t mul, spx_uint32_t div)
     605             : {
     606             :    speex_assert(result);
     607           0 :    spx_uint32_t major = value / div;
     608           0 :    spx_uint32_t remainder = value % div;
     609             :    /* TODO: Could use 64 bits operation to check for overflow. But only guaranteed in C99+ */
     610           0 :    if (remainder > UINT32_MAX / mul || major > UINT32_MAX / mul
     611           0 :        || major * mul > UINT32_MAX - remainder * mul / div)
     612           0 :       return RESAMPLER_ERR_OVERFLOW;
     613           0 :    *result = remainder * mul / div + major * mul;
     614           0 :    return RESAMPLER_ERR_SUCCESS;
     615             : }
     616             : 
     617           0 : static int update_filter(SpeexResamplerState *st)
     618             : {
     619           0 :    spx_uint32_t old_length = st->filt_len;
     620           0 :    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           0 :    st->int_advance = st->num_rate/st->den_rate;
     626           0 :    st->frac_advance = st->num_rate%st->den_rate;
     627           0 :    st->oversample = quality_map[st->quality].oversample;
     628           0 :    st->filt_len = quality_map[st->quality].base_length;
     629             : 
     630           0 :    if (st->num_rate > st->den_rate)
     631             :    {
     632             :       /* down-sampling */
     633           0 :       st->cutoff = quality_map[st->quality].downsample_bandwidth * st->den_rate / st->num_rate;
     634           0 :       if (_muldiv(&st->filt_len,st->filt_len,st->num_rate,st->den_rate) != RESAMPLER_ERR_SUCCESS)
     635           0 :          goto fail;
     636             :       /* Round up to make sure we have a multiple of 8 for SSE */
     637           0 :       st->filt_len = ((st->filt_len-1)&(~0x7))+8;
     638           0 :       if (2*st->den_rate < st->num_rate)
     639           0 :          st->oversample >>= 1;
     640           0 :       if (4*st->den_rate < st->num_rate)
     641           0 :          st->oversample >>= 1;
     642           0 :       if (8*st->den_rate < st->num_rate)
     643           0 :          st->oversample >>= 1;
     644           0 :       if (16*st->den_rate < st->num_rate)
     645           0 :          st->oversample >>= 1;
     646           0 :       if (st->oversample < 1)
     647           0 :          st->oversample = 1;
     648             :    } else {
     649             :       /* up-sampling */
     650           0 :       st->cutoff = quality_map[st->quality].upsample_bandwidth;
     651             :    }
     652             : 
     653           0 :    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           0 :       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           0 :                 && INT_MAX/sizeof(spx_word16_t)/st->den_rate >= st->filt_len;
     663           0 :    if (use_direct)
     664             :    {
     665           0 :       min_sinc_table_length = st->filt_len*st->den_rate;
     666             :    } else {
     667           0 :       if ((INT_MAX/sizeof(spx_word16_t)-8)/st->oversample < st->filt_len)
     668           0 :          goto fail;
     669             : 
     670           0 :       min_sinc_table_length = st->filt_len*st->oversample+8;
     671             :    }
     672           0 :    if (st->sinc_table_length < min_sinc_table_length)
     673             :    {
     674           0 :       spx_word16_t *sinc_table = (spx_word16_t *)speex_realloc(st->sinc_table,min_sinc_table_length*sizeof(spx_word16_t));
     675           0 :       if (!sinc_table)
     676           0 :          goto fail;
     677             : 
     678           0 :       st->sinc_table = sinc_table;
     679           0 :       st->sinc_table_length = min_sinc_table_length;
     680             :    }
     681           0 :    if (use_direct)
     682             :    {
     683             :       spx_uint32_t i;
     684           0 :       for (i=0;i<st->den_rate;i++)
     685             :       {
     686             :          spx_int32_t j;
     687           0 :          for (j=0;j<st->filt_len;j++)
     688             :          {
     689           0 :             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           0 :       if (st->quality>8)
     696           0 :          st->resampler_ptr = resampler_basic_direct_double;
     697             :       else
     698           0 :          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           0 :       for (i=-4;i<(spx_int32_t)(st->oversample*st->filt_len+4);i++)
     704           0 :          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           0 :       if (st->quality>8)
     709           0 :          st->resampler_ptr = resampler_basic_interpolate_double;
     710             :       else
     711           0 :          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           0 :    min_alloc_size = st->filt_len-1 + st->buffer_size;
     723           0 :    if (min_alloc_size > st->mem_alloc_size)
     724             :    {
     725             :       spx_word16_t *mem;
     726           0 :       if (INT_MAX/sizeof(spx_word16_t)/st->nb_channels < min_alloc_size)
     727           0 :           goto fail;
     728           0 :       else if (!(mem = (spx_word16_t*)speex_realloc(st->mem, st->nb_channels*min_alloc_size * sizeof(*mem))))
     729           0 :           goto fail;
     730             : 
     731           0 :       st->mem = mem;
     732           0 :       st->mem_alloc_size = min_alloc_size;
     733             :    }
     734           0 :    if (!st->started)
     735             :    {
     736             :       spx_uint32_t i;
     737           0 :       for (i=0;i<st->nb_channels*st->mem_alloc_size;i++)
     738           0 :          st->mem[i] = 0;
     739             :       /*speex_warning("reinit filter");*/
     740           0 :    } 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           0 :       for (i=st->nb_channels;i--;)
     746             :       {
     747             :          spx_uint32_t j;
     748           0 :          spx_uint32_t olen = old_length;
     749             :          /*if (st->magic_samples[i])*/
     750             :          {
     751             :             /* Try and remove the magic samples as if nothing had happened */
     752             : 
     753             :             /* FIXME: This is wrong but for now we need it to avoid going over the array bounds */
     754           0 :             olen = old_length + 2*st->magic_samples[i];
     755           0 :             for (j=old_length-1+st->magic_samples[i];j--;)
     756           0 :                st->mem[i*st->mem_alloc_size+j+st->magic_samples[i]] = st->mem[i*old_alloc_size+j];
     757           0 :             for (j=0;j<st->magic_samples[i];j++)
     758           0 :                st->mem[i*st->mem_alloc_size+j] = 0;
     759           0 :             st->magic_samples[i] = 0;
     760             :          }
     761           0 :          if (st->filt_len > olen)
     762             :          {
     763             :             /* If the new filter length is still bigger than the "augmented" length */
     764             :             /* Copy data going backward */
     765           0 :             for (j=0;j<olen-1;j++)
     766           0 :                st->mem[i*st->mem_alloc_size+(st->filt_len-2-j)] = st->mem[i*st->mem_alloc_size+(olen-2-j)];
     767             :             /* Then put zeros for lack of anything better */
     768           0 :             for (;j<st->filt_len-1;j++)
     769           0 :                st->mem[i*st->mem_alloc_size+(st->filt_len-2-j)] = 0;
     770             :             /* Adjust last_sample */
     771           0 :             st->last_sample[i] += (st->filt_len - olen)/2;
     772             :          } else {
     773             :             /* Put back some of the magic! */
     774           0 :             st->magic_samples[i] = (olen - st->filt_len)/2;
     775           0 :             for (j=0;j<st->filt_len-1+st->magic_samples[i];j++)
     776           0 :                st->mem[i*st->mem_alloc_size+j] = st->mem[i*st->mem_alloc_size+j+st->magic_samples[i]];
     777             :          }
     778             :       }
     779           0 :    } else if (st->filt_len < old_length)
     780             :    {
     781             :       spx_uint32_t i;
     782             :       /* Reduce filter length, this a bit tricky. We need to store some of the memory as "magic"
     783             :          samples so they can be used directly as input the next time(s) */
     784           0 :       for (i=0;i<st->nb_channels;i++)
     785             :       {
     786             :          spx_uint32_t j;
     787           0 :          spx_uint32_t old_magic = st->magic_samples[i];
     788           0 :          st->magic_samples[i] = (old_length - st->filt_len)/2;
     789             :          /* We must copy some of the memory that's no longer used */
     790             :          /* Copy data going backward */
     791           0 :          for (j=0;j<st->filt_len-1+st->magic_samples[i]+old_magic;j++)
     792           0 :             st->mem[i*st->mem_alloc_size+j] = st->mem[i*st->mem_alloc_size+j+st->magic_samples[i]];
     793           0 :          st->magic_samples[i] += old_magic;
     794             :       }
     795             :    }
     796           0 :    return RESAMPLER_ERR_SUCCESS;
     797             : 
     798             : fail:
     799           0 :    st->resampler_ptr = resampler_basic_zero;
     800             :    /* st->mem may still contain consumed input samples for the filter.
     801             :       Restore filt_len so that filt_len - 1 still points to the position after
     802             :       the last of these samples. */
     803           0 :    st->filt_len = old_length;
     804           0 :    return RESAMPLER_ERR_ALLOC_FAILED;
     805             : }
     806             : 
     807           0 : EXPORT SpeexResamplerState *speex_resampler_init(spx_uint32_t nb_channels, spx_uint32_t in_rate, spx_uint32_t out_rate, int quality, int *err)
     808             : {
     809           0 :    return speex_resampler_init_frac(nb_channels, in_rate, out_rate, in_rate, out_rate, quality, err);
     810             : }
     811             : 
     812           0 : 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)
     813             : {
     814             :    spx_uint32_t i;
     815             :    SpeexResamplerState *st;
     816             :    int filter_err;
     817             : 
     818           0 :    if (nb_channels == 0 || ratio_num == 0 || ratio_den == 0 || quality > 10 || quality < 0)
     819             :    {
     820           0 :       if (err)
     821           0 :          *err = RESAMPLER_ERR_INVALID_ARG;
     822           0 :       return NULL;
     823             :    }
     824           0 :    st = (SpeexResamplerState *)speex_alloc(sizeof(SpeexResamplerState));
     825           0 :    if (!st)
     826             :    {
     827           0 :       if (err)
     828           0 :          *err = RESAMPLER_ERR_ALLOC_FAILED;
     829           0 :       return NULL;
     830             :    }
     831           0 :    st->initialised = 0;
     832           0 :    st->started = 0;
     833           0 :    st->in_rate = 0;
     834           0 :    st->out_rate = 0;
     835           0 :    st->num_rate = 0;
     836           0 :    st->den_rate = 0;
     837           0 :    st->quality = -1;
     838           0 :    st->sinc_table_length = 0;
     839           0 :    st->mem_alloc_size = 0;
     840           0 :    st->filt_len = 0;
     841           0 :    st->mem = 0;
     842           0 :    st->resampler_ptr = 0;
     843             : 
     844           0 :    st->cutoff = 1.f;
     845           0 :    st->nb_channels = nb_channels;
     846           0 :    st->in_stride = 1;
     847           0 :    st->out_stride = 1;
     848             : 
     849           0 :    st->buffer_size = 160;
     850             : 
     851             :    /* Per channel data */
     852           0 :    if (!(st->last_sample = (spx_int32_t*)speex_alloc(nb_channels*sizeof(spx_int32_t))))
     853           0 :       goto fail;
     854           0 :    if (!(st->magic_samples = (spx_uint32_t*)speex_alloc(nb_channels*sizeof(spx_uint32_t))))
     855           0 :       goto fail;
     856           0 :    if (!(st->samp_frac_num = (spx_uint32_t*)speex_alloc(nb_channels*sizeof(spx_uint32_t))))
     857           0 :       goto fail;
     858             : 
     859           0 :    speex_resampler_set_quality(st, quality);
     860           0 :    speex_resampler_set_rate_frac(st, ratio_num, ratio_den, in_rate, out_rate);
     861             : 
     862           0 :    filter_err = update_filter(st);
     863           0 :    if (filter_err == RESAMPLER_ERR_SUCCESS)
     864             :    {
     865           0 :       st->initialised = 1;
     866             :    } else {
     867           0 :       speex_resampler_destroy(st);
     868           0 :       st = NULL;
     869             :    }
     870           0 :    if (err)
     871           0 :       *err = filter_err;
     872             : 
     873           0 :    return st;
     874             : 
     875             : fail:
     876           0 :    if (err)
     877           0 :       *err = RESAMPLER_ERR_ALLOC_FAILED;
     878           0 :    speex_resampler_destroy(st);
     879           0 :    return NULL;
     880             : }
     881             : 
     882           0 : EXPORT void speex_resampler_destroy(SpeexResamplerState *st)
     883             : {
     884           0 :    speex_free(st->mem);
     885           0 :    speex_free(st->sinc_table);
     886           0 :    speex_free(st->last_sample);
     887           0 :    speex_free(st->magic_samples);
     888           0 :    speex_free(st->samp_frac_num);
     889           0 :    speex_free(st);
     890           0 : }
     891             : 
     892           0 : 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)
     893             : {
     894           0 :    int j=0;
     895           0 :    const int N = st->filt_len;
     896           0 :    int out_sample = 0;
     897           0 :    spx_word16_t *mem = st->mem + channel_index * st->mem_alloc_size;
     898             :    spx_uint32_t ilen;
     899             : 
     900           0 :    st->started = 1;
     901             : 
     902             :    /* Call the right resampler through the function ptr */
     903           0 :    out_sample = st->resampler_ptr(st, channel_index, mem, in_len, out, out_len);
     904             : 
     905           0 :    if (st->last_sample[channel_index] < (spx_int32_t)*in_len)
     906           0 :       *in_len = st->last_sample[channel_index];
     907           0 :    *out_len = out_sample;
     908           0 :    st->last_sample[channel_index] -= *in_len;
     909             : 
     910           0 :    ilen = *in_len;
     911             : 
     912           0 :    for(j=0;j<N-1;++j)
     913           0 :      mem[j] = mem[j+ilen];
     914             : 
     915           0 :    return RESAMPLER_ERR_SUCCESS;
     916             : }
     917             : 
     918           0 : static int speex_resampler_magic(SpeexResamplerState *st, spx_uint32_t channel_index, spx_word16_t **out, spx_uint32_t out_len) {
     919           0 :    spx_uint32_t tmp_in_len = st->magic_samples[channel_index];
     920           0 :    spx_word16_t *mem = st->mem + channel_index * st->mem_alloc_size;
     921           0 :    const int N = st->filt_len;
     922             : 
     923           0 :    speex_resampler_process_native(st, channel_index, &tmp_in_len, *out, &out_len);
     924             : 
     925           0 :    st->magic_samples[channel_index] -= tmp_in_len;
     926             : 
     927             :    /* If we couldn't process all "magic" input samples, save the rest for next time */
     928           0 :    if (st->magic_samples[channel_index])
     929             :    {
     930             :       spx_uint32_t i;
     931           0 :       for (i=0;i<st->magic_samples[channel_index];i++)
     932           0 :          mem[N-1+i]=mem[N-1+i+tmp_in_len];
     933             :    }
     934           0 :    *out += out_len*st->out_stride;
     935           0 :    return out_len;
     936             : }
     937             : 
     938             : #ifdef FIXED_POINT
     939             : 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)
     940             : #else
     941           0 : 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)
     942             : #endif
     943             : {
     944             :    int j;
     945           0 :    spx_uint32_t ilen = *in_len;
     946           0 :    spx_uint32_t olen = *out_len;
     947           0 :    spx_word16_t *x = st->mem + channel_index * st->mem_alloc_size;
     948           0 :    const int filt_offs = st->filt_len - 1;
     949           0 :    const spx_uint32_t xlen = st->mem_alloc_size - filt_offs;
     950           0 :    const int istride = st->in_stride;
     951             : 
     952           0 :    if (st->magic_samples[channel_index])
     953           0 :       olen -= speex_resampler_magic(st, channel_index, &out, olen);
     954           0 :    if (! st->magic_samples[channel_index]) {
     955           0 :       while (ilen && olen) {
     956           0 :         spx_uint32_t ichunk = (ilen > xlen) ? xlen : ilen;
     957           0 :         spx_uint32_t ochunk = olen;
     958             : 
     959           0 :         if (in) {
     960           0 :            for(j=0;j<ichunk;++j)
     961           0 :               x[j+filt_offs]=in[j*istride];
     962             :         } else {
     963           0 :           for(j=0;j<ichunk;++j)
     964           0 :             x[j+filt_offs]=0;
     965             :         }
     966           0 :         speex_resampler_process_native(st, channel_index, &ichunk, out, &ochunk);
     967           0 :         ilen -= ichunk;
     968           0 :         olen -= ochunk;
     969           0 :         out += ochunk * st->out_stride;
     970           0 :         if (in)
     971           0 :            in += ichunk * istride;
     972             :       }
     973             :    }
     974           0 :    *in_len -= ilen;
     975           0 :    *out_len -= olen;
     976           0 :    return st->resampler_ptr == resampler_basic_zero ? RESAMPLER_ERR_ALLOC_FAILED : RESAMPLER_ERR_SUCCESS;
     977             : }
     978             : 
     979             : #ifdef FIXED_POINT
     980             : 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)
     981             : #else
     982           0 : 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)
     983             : #endif
     984             : {
     985             :    int j;
     986           0 :    const int istride_save = st->in_stride;
     987           0 :    const int ostride_save = st->out_stride;
     988           0 :    spx_uint32_t ilen = *in_len;
     989           0 :    spx_uint32_t olen = *out_len;
     990           0 :    spx_word16_t *x = st->mem + channel_index * st->mem_alloc_size;
     991           0 :    const spx_uint32_t xlen = st->mem_alloc_size - (st->filt_len - 1);
     992             : #ifdef VAR_ARRAYS
     993             :    const unsigned int ylen = (olen < FIXED_STACK_ALLOC) ? olen : FIXED_STACK_ALLOC;
     994             :    VARDECL(spx_word16_t *ystack);
     995             :    ALLOC(ystack, ylen, spx_word16_t);
     996             : #else
     997           0 :    const unsigned int ylen = FIXED_STACK_ALLOC;
     998             :    spx_word16_t ystack[FIXED_STACK_ALLOC];
     999             : #endif
    1000             : 
    1001           0 :    st->out_stride = 1;
    1002             : 
    1003           0 :    while (ilen && olen) {
    1004           0 :      spx_word16_t *y = ystack;
    1005           0 :      spx_uint32_t ichunk = (ilen > xlen) ? xlen : ilen;
    1006           0 :      spx_uint32_t ochunk = (olen > ylen) ? ylen : olen;
    1007           0 :      spx_uint32_t omagic = 0;
    1008             : 
    1009           0 :      if (st->magic_samples[channel_index]) {
    1010           0 :        omagic = speex_resampler_magic(st, channel_index, &y, ochunk);
    1011           0 :        ochunk -= omagic;
    1012           0 :        olen -= omagic;
    1013             :      }
    1014           0 :      if (! st->magic_samples[channel_index]) {
    1015           0 :        if (in) {
    1016           0 :          for(j=0;j<ichunk;++j)
    1017             : #ifdef FIXED_POINT
    1018             :            x[j+st->filt_len-1]=WORD2INT(in[j*istride_save]);
    1019             : #else
    1020           0 :            x[j+st->filt_len-1]=in[j*istride_save];
    1021             : #endif
    1022             :        } else {
    1023           0 :          for(j=0;j<ichunk;++j)
    1024           0 :            x[j+st->filt_len-1]=0;
    1025             :        }
    1026             : 
    1027           0 :        speex_resampler_process_native(st, channel_index, &ichunk, y, &ochunk);
    1028             :      } else {
    1029           0 :        ichunk = 0;
    1030           0 :        ochunk = 0;
    1031             :      }
    1032             : 
    1033           0 :      for (j=0;j<ochunk+omagic;++j)
    1034             : #ifdef FIXED_POINT
    1035             :         out[j*ostride_save] = ystack[j];
    1036             : #else
    1037           0 :         out[j*ostride_save] = WORD2INT(ystack[j]);
    1038             : #endif
    1039             : 
    1040           0 :      ilen -= ichunk;
    1041           0 :      olen -= ochunk;
    1042           0 :      out += (ochunk+omagic) * ostride_save;
    1043           0 :      if (in)
    1044           0 :        in += ichunk * istride_save;
    1045             :    }
    1046           0 :    st->out_stride = ostride_save;
    1047           0 :    *in_len -= ilen;
    1048           0 :    *out_len -= olen;
    1049             : 
    1050           0 :    return st->resampler_ptr == resampler_basic_zero ? RESAMPLER_ERR_ALLOC_FAILED : RESAMPLER_ERR_SUCCESS;
    1051             : }
    1052             : 
    1053           0 : EXPORT int speex_resampler_process_interleaved_float(SpeexResamplerState *st, const float *in, spx_uint32_t *in_len, float *out, spx_uint32_t *out_len)
    1054             : {
    1055             :    spx_uint32_t i;
    1056             :    int istride_save, ostride_save;
    1057           0 :    spx_uint32_t bak_out_len = *out_len;
    1058           0 :    spx_uint32_t bak_in_len = *in_len;
    1059           0 :    istride_save = st->in_stride;
    1060           0 :    ostride_save = st->out_stride;
    1061           0 :    st->in_stride = st->out_stride = st->nb_channels;
    1062           0 :    for (i=0;i<st->nb_channels;i++)
    1063             :    {
    1064           0 :       *out_len = bak_out_len;
    1065           0 :       *in_len = bak_in_len;
    1066           0 :       if (in != NULL)
    1067           0 :          speex_resampler_process_float(st, i, in+i, in_len, out+i, out_len);
    1068             :       else
    1069           0 :          speex_resampler_process_float(st, i, NULL, in_len, out+i, out_len);
    1070             :    }
    1071           0 :    st->in_stride = istride_save;
    1072           0 :    st->out_stride = ostride_save;
    1073           0 :    return st->resampler_ptr == resampler_basic_zero ? RESAMPLER_ERR_ALLOC_FAILED : RESAMPLER_ERR_SUCCESS;
    1074             : }
    1075             : 
    1076           0 : 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)
    1077             : {
    1078             :    spx_uint32_t i;
    1079             :    int istride_save, ostride_save;
    1080           0 :    spx_uint32_t bak_out_len = *out_len;
    1081           0 :    spx_uint32_t bak_in_len = *in_len;
    1082           0 :    istride_save = st->in_stride;
    1083           0 :    ostride_save = st->out_stride;
    1084           0 :    st->in_stride = st->out_stride = st->nb_channels;
    1085           0 :    for (i=0;i<st->nb_channels;i++)
    1086             :    {
    1087           0 :       *out_len = bak_out_len;
    1088           0 :       *in_len = bak_in_len;
    1089           0 :       if (in != NULL)
    1090           0 :          speex_resampler_process_int(st, i, in+i, in_len, out+i, out_len);
    1091             :       else
    1092           0 :          speex_resampler_process_int(st, i, NULL, in_len, out+i, out_len);
    1093             :    }
    1094           0 :    st->in_stride = istride_save;
    1095           0 :    st->out_stride = ostride_save;
    1096           0 :    return st->resampler_ptr == resampler_basic_zero ? RESAMPLER_ERR_ALLOC_FAILED : RESAMPLER_ERR_SUCCESS;
    1097             : }
    1098             : 
    1099           0 : EXPORT int speex_resampler_set_rate(SpeexResamplerState *st, spx_uint32_t in_rate, spx_uint32_t out_rate)
    1100             : {
    1101           0 :    return speex_resampler_set_rate_frac(st, in_rate, out_rate, in_rate, out_rate);
    1102             : }
    1103             : 
    1104           0 : EXPORT void speex_resampler_get_rate(SpeexResamplerState *st, spx_uint32_t *in_rate, spx_uint32_t *out_rate)
    1105             : {
    1106           0 :    *in_rate = st->in_rate;
    1107           0 :    *out_rate = st->out_rate;
    1108           0 : }
    1109             : 
    1110           0 : static inline spx_uint32_t _gcd(spx_uint32_t a, spx_uint32_t b)
    1111             : {
    1112           0 :    while (b != 0)
    1113             :    {
    1114           0 :       spx_uint32_t temp = a;
    1115             : 
    1116           0 :       a = b;
    1117           0 :       b = temp % b;
    1118             :    }
    1119           0 :    return a;
    1120             : }
    1121             : 
    1122           0 : 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)
    1123             : {
    1124             :    spx_uint32_t fact;
    1125             :    spx_uint32_t old_den;
    1126             :    spx_uint32_t i;
    1127             : 
    1128           0 :    if (ratio_num == 0 || ratio_den == 0)
    1129           0 :       return RESAMPLER_ERR_INVALID_ARG;
    1130             : 
    1131           0 :    if (st->in_rate == in_rate && st->out_rate == out_rate && st->num_rate == ratio_num && st->den_rate == ratio_den)
    1132           0 :       return RESAMPLER_ERR_SUCCESS;
    1133             : 
    1134           0 :    old_den = st->den_rate;
    1135           0 :    st->in_rate = in_rate;
    1136           0 :    st->out_rate = out_rate;
    1137           0 :    st->num_rate = ratio_num;
    1138           0 :    st->den_rate = ratio_den;
    1139             : 
    1140           0 :    fact = _gcd (st->num_rate, st->den_rate);
    1141             : 
    1142           0 :    st->num_rate /= fact;
    1143           0 :    st->den_rate /= fact;
    1144             : 
    1145           0 :    if (old_den > 0)
    1146             :    {
    1147           0 :       for (i=0;i<st->nb_channels;i++)
    1148             :       {
    1149           0 :          if (_muldiv(&st->samp_frac_num[i],st->samp_frac_num[i],st->den_rate,old_den) != RESAMPLER_ERR_SUCCESS) {
    1150           0 :            st->samp_frac_num[i] = st->den_rate-1;
    1151             :          }
    1152             :          /* Safety net */
    1153           0 :          if (st->samp_frac_num[i] >= st->den_rate)
    1154           0 :             st->samp_frac_num[i] = st->den_rate-1;
    1155             :       }
    1156             :    }
    1157             : 
    1158           0 :    if (st->initialised)
    1159           0 :       return update_filter(st);
    1160           0 :    return RESAMPLER_ERR_SUCCESS;
    1161             : }
    1162             : 
    1163           0 : EXPORT void speex_resampler_get_ratio(SpeexResamplerState *st, spx_uint32_t *ratio_num, spx_uint32_t *ratio_den)
    1164             : {
    1165           0 :    *ratio_num = st->num_rate;
    1166           0 :    *ratio_den = st->den_rate;
    1167           0 : }
    1168             : 
    1169           0 : EXPORT int speex_resampler_set_quality(SpeexResamplerState *st, int quality)
    1170             : {
    1171           0 :    if (quality > 10 || quality < 0)
    1172           0 :       return RESAMPLER_ERR_INVALID_ARG;
    1173           0 :    if (st->quality == quality)
    1174           0 :       return RESAMPLER_ERR_SUCCESS;
    1175           0 :    st->quality = quality;
    1176           0 :    if (st->initialised)
    1177           0 :       return update_filter(st);
    1178           0 :    return RESAMPLER_ERR_SUCCESS;
    1179             : }
    1180             : 
    1181           0 : EXPORT void speex_resampler_get_quality(SpeexResamplerState *st, int *quality)
    1182             : {
    1183           0 :    *quality = st->quality;
    1184           0 : }
    1185             : 
    1186           0 : EXPORT void speex_resampler_set_input_stride(SpeexResamplerState *st, spx_uint32_t stride)
    1187             : {
    1188           0 :    st->in_stride = stride;
    1189           0 : }
    1190             : 
    1191           0 : EXPORT void speex_resampler_get_input_stride(SpeexResamplerState *st, spx_uint32_t *stride)
    1192             : {
    1193           0 :    *stride = st->in_stride;
    1194           0 : }
    1195             : 
    1196           0 : EXPORT void speex_resampler_set_output_stride(SpeexResamplerState *st, spx_uint32_t stride)
    1197             : {
    1198           0 :    st->out_stride = stride;
    1199           0 : }
    1200             : 
    1201           0 : EXPORT void speex_resampler_get_output_stride(SpeexResamplerState *st, spx_uint32_t *stride)
    1202             : {
    1203           0 :    *stride = st->out_stride;
    1204           0 : }
    1205             : 
    1206           0 : EXPORT int speex_resampler_get_input_latency(SpeexResamplerState *st)
    1207             : {
    1208           0 :   return st->filt_len / 2;
    1209             : }
    1210             : 
    1211           0 : EXPORT int speex_resampler_get_output_latency(SpeexResamplerState *st)
    1212             : {
    1213           0 :   return ((st->filt_len / 2) * st->den_rate + (st->num_rate >> 1)) / st->num_rate;
    1214             : }
    1215             : 
    1216           0 : EXPORT int speex_resampler_skip_zeros(SpeexResamplerState *st)
    1217             : {
    1218             :    spx_uint32_t i;
    1219           0 :    for (i=0;i<st->nb_channels;i++)
    1220           0 :       st->last_sample[i] = st->filt_len/2;
    1221           0 :    return RESAMPLER_ERR_SUCCESS;
    1222             : }
    1223             : 
    1224           0 : EXPORT int speex_resampler_set_skip_frac_num(SpeexResamplerState *st, spx_uint32_t skip_frac_num)
    1225             : {
    1226             :    spx_uint32_t i;
    1227           0 :    spx_uint32_t last_sample = skip_frac_num / st->den_rate;
    1228           0 :    spx_uint32_t samp_frac_num = skip_frac_num % st->den_rate;
    1229           0 :    for (i=0;i<st->nb_channels;i++) {
    1230           0 :       st->last_sample[i] = last_sample;
    1231           0 :       st->samp_frac_num[i] = samp_frac_num;
    1232             :    }
    1233           0 :    return RESAMPLER_ERR_SUCCESS;
    1234             : }
    1235             : 
    1236           0 : EXPORT int speex_resampler_reset_mem(SpeexResamplerState *st)
    1237             : {
    1238             :    spx_uint32_t i;
    1239           0 :    for (i=0;i<st->nb_channels;i++)
    1240             :    {
    1241           0 :       st->last_sample[i] = 0;
    1242           0 :       st->magic_samples[i] = 0;
    1243           0 :       st->samp_frac_num[i] = 0;
    1244             :    }
    1245           0 :    for (i=0;i<st->nb_channels*(st->filt_len-1);i++)
    1246           0 :       st->mem[i] = 0;
    1247           0 :    return RESAMPLER_ERR_SUCCESS;
    1248             : }
    1249             : 
    1250           0 : EXPORT const char *speex_resampler_strerror(int err)
    1251             : {
    1252           0 :    switch (err)
    1253             :    {
    1254             :       case RESAMPLER_ERR_SUCCESS:
    1255           0 :          return "Success.";
    1256             :       case RESAMPLER_ERR_ALLOC_FAILED:
    1257           0 :          return "Memory allocation failed.";
    1258             :       case RESAMPLER_ERR_BAD_STATE:
    1259           0 :          return "Bad resampler state.";
    1260             :       case RESAMPLER_ERR_INVALID_ARG:
    1261           0 :          return "Invalid argument.";
    1262             :       case RESAMPLER_ERR_PTR_OVERLAP:
    1263           0 :          return "Input and output buffers overlap.";
    1264             :       default:
    1265           0 :          return "Unknown error. Bad error code or strange version mismatch.";
    1266             :    }
    1267             : }

Generated by: LCOV version 1.13