resample.c 20.3 KB
Newer Older
Michael Niedermayer's avatar
Michael Niedermayer committed
1 2
/*
 * audio resampling
3
 * Copyright (c) 2004-2012 Michael Niedermayer <michaelni@gmx.at>
4
 * bessel function: Copyright (c) 2006 Xiaogang Zhang
Michael Niedermayer's avatar
Michael Niedermayer committed
5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28
 *
 * This file is part of FFmpeg.
 *
 * FFmpeg is free software; you can redistribute it and/or
 * modify it under the terms of the GNU Lesser General Public
 * License as published by the Free Software Foundation; either
 * version 2.1 of the License, or (at your option) any later version.
 *
 * FFmpeg is distributed in the hope that it will be useful,
 * but WITHOUT ANY WARRANTY; without even the implied warranty of
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
 * Lesser General Public License for more details.
 *
 * You should have received a copy of the GNU Lesser General Public
 * License along with FFmpeg; if not, write to the Free Software
 * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
 */

/**
 * @file
 * audio resampling
 * @author Michael Niedermayer <michaelni@gmx.at>
 */

29
#include "libavutil/avassert.h"
30
#include "resample.h"
Michael Niedermayer's avatar
Michael Niedermayer committed
31

32 33 34 35 36 37 38 39 40 41
static inline double eval_poly(const double *coeff, int size, double x) {
    double sum = coeff[size-1];
    int i;
    for (i = size-2; i >= 0; --i) {
        sum *= x;
        sum += coeff[i];
    }
    return sum;
}

Michael Niedermayer's avatar
Michael Niedermayer committed
42 43
/**
 * 0th order modified bessel function of the first kind.
44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69
 * Algorithm taken from the Boost project, source:
 * https://searchcode.com/codesearch/view/14918379/
 * Use, modification and distribution are subject to the
 * Boost Software License, Version 1.0 (see notice below).
 * Boost Software License - Version 1.0 - August 17th, 2003
Permission is hereby granted, free of charge, to any person or organization
obtaining a copy of the software and accompanying documentation covered by
this license (the "Software") to use, reproduce, display, distribute,
execute, and transmit the Software, and to prepare derivative works of the
Software, and to permit third-parties to whom the Software is furnished to
do so, all subject to the following:

The copyright notices in the Software and this entire statement, including
the above license grant, this restriction and the following disclaimer,
must be included in all copies of the Software, in whole or in part, and
all derivative works of the Software, unless such copies or derivative
works are solely in the form of machine-executable object code generated by
a source language processor.

THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
FITNESS FOR A PARTICULAR PURPOSE, TITLE AND NON-INFRINGEMENT. IN NO EVENT
SHALL THE COPYRIGHT HOLDERS OR ANYONE DISTRIBUTING THE SOFTWARE BE LIABLE
FOR ANY DAMAGES OR OTHER LIABILITY, WHETHER IN CONTRACT, TORT OR OTHERWISE,
ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
DEALINGS IN THE SOFTWARE.
Michael Niedermayer's avatar
Michael Niedermayer committed
70 71
 */

72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132
static double bessel(double x) {
// Modified Bessel function of the first kind of order zero
// minimax rational approximations on intervals, see
// Blair and Edwards, Chalk River Report AECL-4928, 1974
    static const double p1[] = {
        -2.2335582639474375249e+15,
        -5.5050369673018427753e+14,
        -3.2940087627407749166e+13,
        -8.4925101247114157499e+11,
        -1.1912746104985237192e+10,
        -1.0313066708737980747e+08,
        -5.9545626019847898221e+05,
        -2.4125195876041896775e+03,
        -7.0935347449210549190e+00,
        -1.5453977791786851041e-02,
        -2.5172644670688975051e-05,
        -3.0517226450451067446e-08,
        -2.6843448573468483278e-11,
        -1.5982226675653184646e-14,
        -5.2487866627945699800e-18,
    };
    static const double q1[] = {
        -2.2335582639474375245e+15,
         7.8858692566751002988e+12,
        -1.2207067397808979846e+10,
         1.0377081058062166144e+07,
        -4.8527560179962773045e+03,
         1.0L,
    };
    static const double p2[] = {
        -2.2210262233306573296e-04,
         1.3067392038106924055e-02,
        -4.4700805721174453923e-01,
         5.5674518371240761397e+00,
        -2.3517945679239481621e+01,
         3.1611322818701131207e+01,
        -9.6090021968656180000e+00,
    };
    static const double q2[] = {
        -5.5194330231005480228e-04,
         3.2547697594819615062e-02,
        -1.1151759188741312645e+00,
         1.3982595353892851542e+01,
        -6.0228002066743340583e+01,
         8.5539563258012929600e+01,
        -3.1446690275135491500e+01,
        1.0L,
    };
    double y, r, factor;
    if (x == 0)
        return 1.0;
    x = fabs(x);
    if (x <= 15) {
        y = x * x;
        return eval_poly(p1, FF_ARRAY_ELEMS(p1), y) / eval_poly(q1, FF_ARRAY_ELEMS(q1), y);
    }
    else {
        y = 1 / x - 1.0 / 15;
        r = eval_poly(p2, FF_ARRAY_ELEMS(p2), y) / eval_poly(q2, FF_ARRAY_ELEMS(q2), y);
        factor = exp(x) / sqrt(x);
        return factor * r;
Michael Niedermayer's avatar
Michael Niedermayer committed
133 134 135 136 137 138 139
    }
}

/**
 * builds a polyphase filterbank.
 * @param factor resampling factor
 * @param scale wanted sum of coefficients for each filter
140 141
 * @param filter_type  filter type
 * @param kaiser_beta  kaiser window beta
Michael Niedermayer's avatar
Michael Niedermayer committed
142 143
 * @return 0 on success, negative on error
 */
144
static int build_filter(ResampleContext *c, void *filter, double factor, int tap_count, int alloc, int phase_count, int scale,
145
                        int filter_type, double kaiser_beta){
Michael Niedermayer's avatar
Michael Niedermayer committed
146
    int ph, i;
147
    double x, y, w, t, s;
148
    double *tab = av_malloc_array(tap_count+1,  sizeof(*tab));
149
    double *sin_lut = av_malloc_array(phase_count / 2 + 1, sizeof(*sin_lut));
Michael Niedermayer's avatar
Michael Niedermayer committed
150 151
    const int center= (tap_count-1)/2;

152 153
    if (!tab || !sin_lut)
        goto fail;
Michael Niedermayer's avatar
Michael Niedermayer committed
154 155 156 157 158

    /* if upsampling, only need to interpolate, no filter */
    if (factor > 1.0)
        factor = 1.0;

159
    av_assert0(phase_count == 1 || phase_count % 2 == 0);
160 161 162 163 164

    if (factor == 1.0) {
        for (ph = 0; ph <= phase_count / 2; ph++)
            sin_lut[ph] = sin(M_PI * ph / phase_count);
    }
165
    for(ph = 0; ph <= phase_count / 2; ph++) {
Michael Niedermayer's avatar
Michael Niedermayer committed
166
        double norm = 0;
167
        s = sin_lut[ph];
168
        for(i=0;i<=tap_count;i++) {
Michael Niedermayer's avatar
Michael Niedermayer committed
169 170
            x = M_PI * ((double)(i - center) - (double)ph / phase_count) * factor;
            if (x == 0) y = 1.0;
171 172 173 174
            else if (factor == 1.0)
                y = s / x;
            else
                y = sin(x) / x;
175 176
            switch(filter_type){
            case SWR_FILTER_TYPE_CUBIC:{
Michael Niedermayer's avatar
Michael Niedermayer committed
177 178 179 180 181
                const float d= -0.5; //first order derivative = -0.5
                x = fabs(((double)(i - center) - (double)ph / phase_count) * factor);
                if(x<1.0) y= 1 - 3*x*x + 2*x*x*x + d*(            -x*x + x*x*x);
                else      y=                       d*(-4 + 8*x - 5*x*x + x*x*x);
                break;}
182
            case SWR_FILTER_TYPE_BLACKMAN_NUTTALL:
183 184
                w = 2.0*x / (factor*tap_count);
                t = -cos(w);
185
                y *= 0.3635819 - 0.4891775 * t + 0.1365995 * (2*t*t-1) - 0.0106411 * (4*t*t*t - 3*t);
Michael Niedermayer's avatar
Michael Niedermayer committed
186
                break;
187
            case SWR_FILTER_TYPE_KAISER:
Michael Niedermayer's avatar
Michael Niedermayer committed
188
                w = 2.0*x / (factor*tap_count*M_PI);
189
                y *= bessel(kaiser_beta*sqrt(FFMAX(1-w*w, 0)));
Michael Niedermayer's avatar
Michael Niedermayer committed
190
                break;
191 192
            default:
                av_assert0(0);
Michael Niedermayer's avatar
Michael Niedermayer committed
193 194 195
            }

            tab[i] = y;
196
            s = -s;
197 198
            if (i < tap_count)
                norm += y;
Michael Niedermayer's avatar
Michael Niedermayer committed
199 200 201
        }

        /* normalize so that an uniform color remains the same */
202
        switch(c->format){
203
        case AV_SAMPLE_FMT_S16P:
204
            for(i=0;i<tap_count;i++)
205
                ((int16_t*)filter)[ph * alloc + i] = av_clip(lrintf(tab[i] * scale / norm), INT16_MIN, INT16_MAX);
206 207 208 209 210 211 212 213 214
            if (tap_count % 2 == 0) {
                for (i = 0; i < tap_count; i++)
                    ((int16_t*)filter)[(phase_count-ph) * alloc + tap_count-1-i] = ((int16_t*)filter)[ph * alloc + i];
            }
            else {
                for (i = 1; i <= tap_count; i++)
                    ((int16_t*)filter)[(phase_count-ph) * alloc + tap_count-i] =
                        av_clip(lrintf(tab[i] * scale / (norm - tab[0] + tab[tap_count])), INT16_MIN, INT16_MAX);
            }
215
            break;
216
        case AV_SAMPLE_FMT_S32P:
217
            for(i=0;i<tap_count;i++)
218
                ((int32_t*)filter)[ph * alloc + i] = av_clipl_int32(llrint(tab[i] * scale / norm));
219 220 221 222 223 224 225 226 227
            if (tap_count % 2 == 0) {
                for (i = 0; i < tap_count; i++)
                    ((int32_t*)filter)[(phase_count-ph) * alloc + tap_count-1-i] = ((int32_t*)filter)[ph * alloc + i];
            }
            else {
                for (i = 1; i <= tap_count; i++)
                    ((int32_t*)filter)[(phase_count-ph) * alloc + tap_count-i] =
                        av_clipl_int32(llrint(tab[i] * scale / (norm - tab[0] + tab[tap_count])));
            }
228
            break;
229
        case AV_SAMPLE_FMT_FLTP:
230
            for(i=0;i<tap_count;i++)
231
                ((float*)filter)[ph * alloc + i] = tab[i] * scale / norm;
232 233 234 235 236 237 238 239
            if (tap_count % 2 == 0) {
                for (i = 0; i < tap_count; i++)
                    ((float*)filter)[(phase_count-ph) * alloc + tap_count-1-i] = ((float*)filter)[ph * alloc + i];
            }
            else {
                for (i = 1; i <= tap_count; i++)
                    ((float*)filter)[(phase_count-ph) * alloc + tap_count-i] = tab[i] * scale / (norm - tab[0] + tab[tap_count]);
            }
240
            break;
241
        case AV_SAMPLE_FMT_DBLP:
242
            for(i=0;i<tap_count;i++)
243
                ((double*)filter)[ph * alloc + i] = tab[i] * scale / norm;
244 245 246 247 248 249 250 251
            if (tap_count % 2 == 0) {
                for (i = 0; i < tap_count; i++)
                    ((double*)filter)[(phase_count-ph) * alloc + tap_count-1-i] = ((double*)filter)[ph * alloc + i];
            }
            else {
                for (i = 1; i <= tap_count; i++)
                    ((double*)filter)[(phase_count-ph) * alloc + tap_count-i] = tab[i] * scale / (norm - tab[0] + tab[tap_count]);
            }
252
            break;
Michael Niedermayer's avatar
Michael Niedermayer committed
253 254 255 256 257 258 259 260 261 262 263 264 265 266 267 268 269 270 271 272 273 274 275 276 277 278 279 280 281 282 283 284 285 286 287 288 289 290 291
        }
    }
#if 0
    {
#define LEN 1024
        int j,k;
        double sine[LEN + tap_count];
        double filtered[LEN];
        double maxff=-2, minff=2, maxsf=-2, minsf=2;
        for(i=0; i<LEN; i++){
            double ss=0, sf=0, ff=0;
            for(j=0; j<LEN+tap_count; j++)
                sine[j]= cos(i*j*M_PI/LEN);
            for(j=0; j<LEN; j++){
                double sum=0;
                ph=0;
                for(k=0; k<tap_count; k++)
                    sum += filter[ph * tap_count + k] * sine[k+j];
                filtered[j]= sum / (1<<FILTER_SHIFT);
                ss+= sine[j + center] * sine[j + center];
                ff+= filtered[j] * filtered[j];
                sf+= sine[j + center] * filtered[j];
            }
            ss= sqrt(2*ss/LEN);
            ff= sqrt(2*ff/LEN);
            sf= 2*sf/LEN;
            maxff= FFMAX(maxff, ff);
            minff= FFMIN(minff, ff);
            maxsf= FFMAX(maxsf, sf);
            minsf= FFMIN(minsf, sf);
            if(i%11==0){
                av_log(NULL, AV_LOG_ERROR, "i:%4d ss:%f ff:%13.6e-%13.6e sf:%13.6e-%13.6e\n", i, ss, maxff, minff, maxsf, minsf);
                minff=minsf= 2;
                maxff=maxsf= -2;
            }
        }
    }
#endif

292
fail:
Michael Niedermayer's avatar
Michael Niedermayer committed
293
    av_free(tab);
294
    av_free(sin_lut);
Michael Niedermayer's avatar
Michael Niedermayer committed
295 296 297
    return 0;
}

298
static ResampleContext *resample_init(ResampleContext *c, int out_rate, int in_rate, int filter_size, int phase_shift, int linear,
299
                                    double cutoff0, enum AVSampleFormat format, enum SwrFilterType filter_type, double kaiser_beta,
300 301
                                    double precision, int cheby)
{
302
    double cutoff = cutoff0? cutoff0 : 0.97;
Michael Niedermayer's avatar
Michael Niedermayer committed
303 304 305
    double factor= FFMIN(out_rate * cutoff / in_rate, 1.0);
    int phase_count= 1<<phase_shift;

306
    if (!c || c->phase_shift != phase_shift || c->linear!=linear || c->factor != factor
307 308
           || c->filter_length != FFMAX((int)ceil(filter_size/factor), 1) || c->format != format
           || c->filter_type != filter_type || c->kaiser_beta != kaiser_beta) {
309
        c = av_mallocz(sizeof(*c));
310 311 312
        if (!c)
            return NULL;

313 314
        c->format= format;

315 316
        c->felem_size= av_get_bytes_per_sample(c->format);

317
        switch(c->format){
318
        case AV_SAMPLE_FMT_S16P:
319 320
            c->filter_shift = 15;
            break;
321
        case AV_SAMPLE_FMT_S32P:
322 323
            c->filter_shift = 30;
            break;
324 325
        case AV_SAMPLE_FMT_FLTP:
        case AV_SAMPLE_FMT_DBLP:
326 327
            c->filter_shift = 0;
            break;
328 329
        default:
            av_log(NULL, AV_LOG_ERROR, "Unsupported sample format\n");
330
            av_assert0(0);
331 332
        }

333 334 335 336 337
        if (filter_size/factor > INT32_MAX/256) {
            av_log(NULL, AV_LOG_ERROR, "Filter length too large\n");
            goto error;
        }

338 339 340 341 342
        c->phase_shift   = phase_shift;
        c->phase_mask    = phase_count - 1;
        c->linear        = linear;
        c->factor        = factor;
        c->filter_length = FFMAX((int)ceil(filter_size/factor), 1);
343
        c->filter_alloc  = FFALIGN(c->filter_length, 8);
344
        c->filter_bank   = av_calloc(c->filter_alloc, (phase_count+1)*c->felem_size);
345 346
        c->filter_type   = filter_type;
        c->kaiser_beta   = kaiser_beta;
347 348
        if (!c->filter_bank)
            goto error;
349
        if (build_filter(c, (void*)c->filter_bank, factor, c->filter_length, c->filter_alloc, phase_count, 1<<c->filter_shift, filter_type, kaiser_beta))
350
            goto error;
351 352
        memcpy(c->filter_bank + (c->filter_alloc*phase_count+1)*c->felem_size, c->filter_bank, (c->filter_alloc-1)*c->felem_size);
        memcpy(c->filter_bank + (c->filter_alloc*phase_count  )*c->felem_size, c->filter_bank + (c->filter_alloc - 1)*c->felem_size, c->felem_size);
Michael Niedermayer's avatar
Michael Niedermayer committed
353 354 355
    }

    c->compensation_distance= 0;
356 357
    if(!av_reduce(&c->src_incr, &c->dst_incr, out_rate, in_rate * (int64_t)phase_count, INT32_MAX/2))
        goto error;
358 359 360
    c->ideal_dst_incr = c->dst_incr;
    c->dst_incr_div   = c->dst_incr / c->src_incr;
    c->dst_incr_mod   = c->dst_incr % c->src_incr;
361

Michael Niedermayer's avatar
Michael Niedermayer committed
362 363 364
    c->index= -phase_count*((c->filter_length-1)/2);
    c->frac= 0;

365
    swri_resample_dsp_init(c);
366

Michael Niedermayer's avatar
Michael Niedermayer committed
367 368
    return c;
error:
369
    av_freep(&c->filter_bank);
Michael Niedermayer's avatar
Michael Niedermayer committed
370 371 372 373
    av_free(c);
    return NULL;
}

374
static void resample_free(ResampleContext **c){
Michael Niedermayer's avatar
Michael Niedermayer committed
375 376 377 378 379 380
    if(!*c)
        return;
    av_freep(&(*c)->filter_bank);
    av_freep(c);
}

381
static int set_compensation(ResampleContext *c, int sample_delta, int compensation_distance){
Michael Niedermayer's avatar
Michael Niedermayer committed
382
    c->compensation_distance= compensation_distance;
383 384 385 386
    if (compensation_distance)
        c->dst_incr = c->ideal_dst_incr - c->ideal_dst_incr * (int64_t)sample_delta / compensation_distance;
    else
        c->dst_incr = c->ideal_dst_incr;
387 388 389 390

    c->dst_incr_div   = c->dst_incr / c->src_incr;
    c->dst_incr_mod   = c->dst_incr % c->src_incr;

391
    return 0;
Michael Niedermayer's avatar
Michael Niedermayer committed
392 393
}

394 395 396 397 398 399 400 401 402 403 404 405
static int swri_resample(ResampleContext *c,
                         uint8_t *dst, const uint8_t *src, int *consumed,
                         int src_size, int dst_size, int update_ctx)
{
    if (c->filter_length == 1 && c->phase_shift == 0) {
        int index= c->index;
        int frac= c->frac;
        int64_t index2= (1LL<<32)*c->frac/c->src_incr + (1LL<<32)*index;
        int64_t incr= (1LL<<32) * c->dst_incr / c->src_incr;
        int new_size = (src_size * (int64_t)c->src_incr - frac + c->dst_incr - 1) / c->dst_incr;

        dst_size= FFMIN(dst_size, new_size);
406
        c->dsp.resample_one(dst, src, dst_size, index2, incr);
407

408 409
        index += dst_size * c->dst_incr_div;
        index += (frac + dst_size * (int64_t)c->dst_incr_mod) / c->src_incr;
410 411 412
        av_assert2(index >= 0);
        *consumed= index;
        if (update_ctx) {
413
            c->frac   = (frac + dst_size * (int64_t)c->dst_incr_mod) % c->src_incr;
414 415 416 417 418 419 420 421
            c->index = 0;
        }
    } else {
        int64_t end_index = (1LL + src_size - c->filter_length) << c->phase_shift;
        int64_t delta_frac = (end_index - c->index) * c->src_incr - c->frac;
        int delta_n = (delta_frac + c->dst_incr - 1) / c->dst_incr;

        dst_size = FFMIN(dst_size, delta_n);
422
        if (dst_size > 0) {
423
            *consumed = c->dsp.resample(c, dst, src, dst_size, update_ctx);
424
        } else {
425
            *consumed = 0;
426 427
        }
    }
428

429 430
    return dst_size;
}
Michael Niedermayer's avatar
Michael Niedermayer committed
431

432
static int multiple_resample(ResampleContext *c, AudioData *dst, int dst_size, AudioData *src, int src_size, int *consumed){
Michael Niedermayer's avatar
Michael Niedermayer committed
433
    int i, ret= -1;
434
    int av_unused mm_flags = av_get_cpu_flags();
435 436
    int need_emms = c->format == AV_SAMPLE_FMT_S16P && ARCH_X86_32 &&
                    (mm_flags & (AV_CPU_FLAG_MMX2 | AV_CPU_FLAG_SSE2)) == AV_CPU_FLAG_MMX2;
437
    int64_t max_src_size = (INT64_MAX >> (c->phase_shift+1)) / c->src_incr;
Michael Niedermayer's avatar
Michael Niedermayer committed
438

439 440
    if (c->compensation_distance)
        dst_size = FFMIN(dst_size, c->compensation_distance);
441
    src_size = FFMIN(src_size, max_src_size);
442

Michael Niedermayer's avatar
Michael Niedermayer committed
443
    for(i=0; i<dst->ch_count; i++){
444 445
        ret= swri_resample(c, dst->ch[i], src->ch[i],
                           consumed, src_size, dst_size, i+1==dst->ch_count);
Michael Niedermayer's avatar
Michael Niedermayer committed
446
    }
447 448
    if(need_emms)
        emms_c();
449 450 451

    if (c->compensation_distance) {
        c->compensation_distance -= ret;
452 453 454 455 456
        if (!c->compensation_distance) {
            c->dst_incr     = c->ideal_dst_incr;
            c->dst_incr_div = c->dst_incr / c->src_incr;
            c->dst_incr_mod = c->dst_incr % c->src_incr;
        }
457 458
    }

Michael Niedermayer's avatar
Michael Niedermayer committed
459 460
    return ret;
}
461

462
static int64_t get_delay(struct SwrContext *s, int64_t base){
463
    ResampleContext *c = s->resample;
464
    int64_t num = s->in_buffer_count - (c->filter_length-1)/2;
465
    num *= 1 << c->phase_shift;
466 467 468 469
    num -= c->index;
    num *= c->src_incr;
    num -= c->frac;
    return av_rescale(num, base, s->in_sample_rate*(int64_t)c->src_incr << c->phase_shift);
470
}
471

472 473
static int64_t get_out_samples(struct SwrContext *s, int in_samples) {
    ResampleContext *c = s->resample;
474
    // The + 2 are added to allow implementations to be slightly inaccurate, they should not be needed currently.
475
    // They also make it easier to proof that changes and optimizations do not
476
    // break the upper bound.
477 478 479 480 481 482 483 484 485 486 487 488 489 490
    int64_t num = s->in_buffer_count + 2LL + in_samples;
    num *= 1 << c->phase_shift;
    num -= c->index;
    num = av_rescale_rnd(num, s->out_sample_rate, ((int64_t)s->in_sample_rate) << c->phase_shift, AV_ROUND_UP) + 2;

    if (c->compensation_distance) {
        if (num > INT_MAX)
            return AVERROR(EINVAL);

        num = FFMAX(num, (num * c->ideal_dst_incr - 1) / c->dst_incr + 1);
    }
    return num;
}

491 492 493 494 495 496 497 498 499 500 501 502 503 504 505 506
static int resample_flush(struct SwrContext *s) {
    AudioData *a= &s->in_buffer;
    int i, j, ret;
    if((ret = swri_realloc_audio(a, s->in_buffer_index + 2*s->in_buffer_count)) < 0)
        return ret;
    av_assert0(a->planar);
    for(i=0; i<a->ch_count; i++){
        for(j=0; j<s->in_buffer_count; j++){
            memcpy(a->ch[i] + (s->in_buffer_index+s->in_buffer_count+j  )*a->bps,
                a->ch[i] + (s->in_buffer_index+s->in_buffer_count-j-1)*a->bps, a->bps);
        }
    }
    s->in_buffer_count += (s->in_buffer_count+1)/2;
    return 0;
}

507 508 509 510 511 512 513 514 515 516 517 518 519 520 521 522 523 524 525 526 527 528 529 530 531 532 533 534 535 536 537 538 539 540 541 542 543 544
// in fact the whole handle multiple ridiculously small buffers might need more thinking...
static int invert_initial_buffer(ResampleContext *c, AudioData *dst, const AudioData *src,
                                 int in_count, int *out_idx, int *out_sz)
{
    int n, ch, num = FFMIN(in_count + *out_sz, c->filter_length + 1), res;

    if (c->index >= 0)
        return 0;

    if ((res = swri_realloc_audio(dst, c->filter_length * 2 + 1)) < 0)
        return res;

    // copy
    for (n = *out_sz; n < num; n++) {
        for (ch = 0; ch < src->ch_count; ch++) {
            memcpy(dst->ch[ch] + ((c->filter_length + n) * c->felem_size),
                   src->ch[ch] + ((n - *out_sz) * c->felem_size), c->felem_size);
        }
    }

    // if not enough data is in, return and wait for more
    if (num < c->filter_length + 1) {
        *out_sz = num;
        *out_idx = c->filter_length;
        return INT_MAX;
    }

    // else invert
    for (n = 1; n <= c->filter_length; n++) {
        for (ch = 0; ch < src->ch_count; ch++) {
            memcpy(dst->ch[ch] + ((c->filter_length - n) * c->felem_size),
                   dst->ch[ch] + ((c->filter_length + n) * c->felem_size),
                   c->felem_size);
        }
    }

    res = num - *out_sz;
    *out_idx = c->filter_length + (c->index >> c->phase_shift);
545 546
    *out_sz = FFMAX(*out_sz + c->filter_length,
                    1 + c->filter_length * 2) - *out_idx;
547 548
    c->index &= c->phase_mask;

549
    return FFMAX(res, 0);
550 551
}

552 553 554 555
struct Resampler const swri_resampler={
  resample_init,
  resample_free,
  multiple_resample,
556
  resample_flush,
557 558
  set_compensation,
  get_delay,
559
  invert_initial_buffer,
560
  get_out_samples,
561
};