77 static const double p1[] = {
78 -2.2335582639474375249e+15,
79 -5.5050369673018427753e+14,
80 -3.2940087627407749166e+13,
81 -8.4925101247114157499e+11,
82 -1.1912746104985237192e+10,
83 -1.0313066708737980747e+08,
84 -5.9545626019847898221e+05,
85 -2.4125195876041896775e+03,
86 -7.0935347449210549190e+00,
87 -1.5453977791786851041e-02,
88 -2.5172644670688975051e-05,
89 -3.0517226450451067446e-08,
90 -2.6843448573468483278e-11,
91 -1.5982226675653184646e-14,
92 -5.2487866627945699800e-18,
94 static const double q1[] = {
95 -2.2335582639474375245e+15,
96 7.8858692566751002988e+12,
97 -1.2207067397808979846e+10,
98 1.0377081058062166144e+07,
99 -4.8527560179962773045e+03,
102 static const double p2[] = {
103 -2.2210262233306573296e-04,
104 1.3067392038106924055e-02,
105 -4.4700805721174453923e-01,
106 5.5674518371240761397e+00,
107 -2.3517945679239481621e+01,
108 3.1611322818701131207e+01,
109 -9.6090021968656180000e+00,
111 static const double q2[] = {
112 -5.5194330231005480228e-04,
113 3.2547697594819615062e-02,
114 -1.1151759188741312645e+00,
115 1.3982595353892851542e+01,
116 -6.0228002066743340583e+01,
117 8.5539563258012929600e+01,
118 -3.1446690275135491500e+01,
130 y = 1 / x - 1.0 / 15;
148 int ph_nb = phase_count % 2 ? phase_count : phase_count / 2 + 1;
149 double x, y,
w, t,
s;
152 const int center= (tap_count-1)/2;
156 if (!
tab || !sin_lut)
159 av_assert0(tap_count == 1 || tap_count % 2 == 0);
166 for (ph = 0; ph < ph_nb; ph++)
167 sin_lut[ph] = sin(
M_PI * ph / phase_count) * (center & 1 ? 1 : -1);
169 for(ph = 0; ph < ph_nb; ph++) {
171 for(
i=0;
i<tap_count;
i++) {
181 x =
fabs(((
double)(
i - center) - (
double)ph / phase_count) *
factor);
182 if(x<1.0) y= 1 - 3*x*x + 2*x*x*x +
d*( -x*x + x*x*x);
183 else y=
d*(-4 + 8*x - 5*x*x + x*x*x);
186 w = 2.0*x / (
factor*tap_count);
188 y *= 0.3635819 - 0.4891775 * t + 0.1365995 * (2*t*t-1) - 0.0106411 * (4*t*t*t - 3*t);
207 for(
i=0;
i<tap_count;
i++)
209 if (phase_count % 2)
break;
210 for (
i = 0;
i < tap_count;
i++)
211 ((int16_t*)
filter)[(phase_count-ph) * alloc + tap_count-1-
i] = ((int16_t*)
filter)[ph * alloc +
i];
214 for(
i=0;
i<tap_count;
i++)
216 if (phase_count % 2)
break;
217 for (
i = 0;
i < tap_count;
i++)
221 for(
i=0;
i<tap_count;
i++)
223 if (phase_count % 2)
break;
224 for (
i = 0;
i < tap_count;
i++)
225 ((
float*)
filter)[(phase_count-ph) * alloc + tap_count-1-
i] = ((
float*)
filter)[ph * alloc +
i];
228 for(
i=0;
i<tap_count;
i++)
230 if (phase_count % 2)
break;
231 for (
i = 0;
i < tap_count;
i++)
232 ((
double*)
filter)[(phase_count-ph) * alloc + tap_count-1-
i] = ((
double*)
filter)[ph * alloc +
i];
240 double sine[
LEN + tap_count];
241 double filtered[
LEN];
242 double maxff=-2, minff=2, maxsf=-2, minsf=2;
244 double ss=0, sf=0, ff=0;
245 for(j=0; j<
LEN+tap_count; j++)
247 for(j=0; j<
LEN; j++){
250 for(k=0; k<tap_count; k++)
251 sum +=
filter[ph * tap_count + k] * sine[k+j];
252 filtered[j]= sum / (1<<FILTER_SHIFT);
253 ss+= sine[j + center] * sine[j + center];
254 ff+= filtered[j] * filtered[j];
255 sf+= sine[j + center] * filtered[j];
260 maxff=
FFMAX(maxff, ff);
261 minff=
FFMIN(minff, ff);
262 maxsf=
FFMAX(maxsf, sf);
263 minsf=
FFMIN(minsf, sf);
290 double precision,
int cheby,
int exact_rational)
292 double cutoff = cutoff0? cutoff0 : 0.97;
293 double factor=
FFMIN(out_rate * cutoff / in_rate, 1.0);
294 int phase_count= 1<<phase_shift;
295 int phase_count_compensation = phase_count;
298 if (filter_length > 1)
299 filter_length =
FFALIGN(filter_length, 2);
301 if (exact_rational) {
302 int phase_count_exact, phase_count_exact_den;
304 av_reduce(&phase_count_exact, &phase_count_exact_den, out_rate, in_rate, INT_MAX);
305 if (phase_count_exact <= phase_count) {
306 phase_count_compensation = phase_count_exact * (phase_count / phase_count_exact);
307 phase_count = phase_count_exact;
311 if (!
c ||
c->phase_count != phase_count ||
c->linear!=
linear ||
c->factor !=
factor
312 ||
c->filter_length != filter_length ||
c->format !=
format
313 ||
c->filter_type != filter_type ||
c->kaiser_beta !=
kaiser_beta) {
325 c->filter_shift = 15;
328 c->filter_shift = 30;
339 if (filter_size/
factor > INT32_MAX/256) {
344 c->phase_count = phase_count;
347 c->filter_length = filter_length;
348 c->filter_alloc =
FFALIGN(
c->filter_length, 8);
349 c->filter_bank =
av_calloc(
c->filter_alloc, (phase_count+1)*
c->felem_size);
350 c->filter_type = filter_type;
352 c->phase_count_compensation = phase_count_compensation;
357 memcpy(
c->filter_bank + (
c->filter_alloc*phase_count+1)*
c->felem_size,
c->filter_bank, (
c->filter_alloc-1)*
c->felem_size);
358 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);
361 c->compensation_distance= 0;
362 if(!
av_reduce(&
c->src_incr, &
c->dst_incr, out_rate, in_rate * (int64_t)phase_count, INT32_MAX/2))
364 while (
c->dst_incr < (1<<20) &&
c->src_incr < (1<<20)) {
368 c->ideal_dst_incr =
c->dst_incr;
369 c->dst_incr_div =
c->dst_incr /
c->src_incr;
370 c->dst_incr_mod =
c->dst_incr %
c->src_incr;
372 c->index= -phase_count*((
c->filter_length-1)/2);
386 uint8_t *new_filter_bank;
387 int new_src_incr, new_dst_incr;
388 int phase_count =
c->phase_count_compensation;
391 if (phase_count ==
c->phase_count)
396 new_filter_bank =
av_calloc(
c->filter_alloc, (phase_count + 1) *
c->felem_size);
397 if (!new_filter_bank)
401 phase_count, 1 <<
c->filter_shift,
c->filter_type,
c->kaiser_beta);
406 memcpy(new_filter_bank + (
c->filter_alloc*phase_count+1)*
c->felem_size, new_filter_bank, (
c->filter_alloc-1)*
c->felem_size);
407 memcpy(new_filter_bank + (
c->filter_alloc*phase_count )*
c->felem_size, new_filter_bank + (
c->filter_alloc - 1)*
c->felem_size,
c->felem_size);
409 if (!
av_reduce(&new_src_incr, &new_dst_incr,
c->src_incr,
410 c->dst_incr * (int64_t)(phase_count/
c->phase_count), INT32_MAX/2))
416 c->src_incr = new_src_incr;
417 c->dst_incr = new_dst_incr;
418 while (
c->dst_incr < (1<<20) &&
c->src_incr < (1<<20)) {
422 c->ideal_dst_incr =
c->dst_incr;
423 c->dst_incr_div =
c->dst_incr /
c->src_incr;
424 c->dst_incr_mod =
c->dst_incr %
c->src_incr;
425 c->index *= phase_count /
c->phase_count;
426 c->phase_count = phase_count;
428 c->filter_bank = new_filter_bank;
435 if (compensation_distance && sample_delta) {
441 c->compensation_distance= compensation_distance;
442 if (compensation_distance)
443 c->dst_incr =
c->ideal_dst_incr -
c->ideal_dst_incr * (int64_t)sample_delta / compensation_distance;
445 c->dst_incr =
c->ideal_dst_incr;
447 c->dst_incr_div =
c->dst_incr /
c->src_incr;
448 c->dst_incr_mod =
c->dst_incr %
c->src_incr;
455 int64_t max_src_size = (INT64_MAX/2 /
c->phase_count) /
c->src_incr;
457 if (
c->compensation_distance)
458 dst_size =
FFMIN(dst_size,
c->compensation_distance);
459 src_size =
FFMIN(src_size, max_src_size);
463 if (
c->filter_length == 1 &&
c->phase_count == 1) {
464 int64_t index2= (1LL<<32)*
c->frac/
c->src_incr + (1LL<<32)*
c->index;
465 int64_t incr= (1LL<<32) *
c->dst_incr /
c->src_incr;
466 int new_size = (src_size * (int64_t)
c->src_incr -
c->frac +
c->dst_incr - 1) /
c->dst_incr;
468 dst_size =
FFMAX(
FFMIN(dst_size, new_size), 0);
471 c->dsp.resample_one(dst->
ch[
i],
src->ch[
i], dst_size, index2, incr);
473 c->index += dst_size *
c->dst_incr_div;
474 c->index += (
c->frac + dst_size * (int64_t)
c->dst_incr_mod) /
c->src_incr;
476 *consumed =
c->index;
477 c->frac = (
c->frac + dst_size * (int64_t)
c->dst_incr_mod) %
c->src_incr;
483 int64_t end_index = (1LL + src_size -
c->filter_length) *
c->phase_count;
484 int64_t delta_frac = (end_index -
c->index) *
c->src_incr -
c->frac;
485 int delta_n = (delta_frac +
c->dst_incr - 1) /
c->dst_incr;
487 const void *
src,
int n,
int update_ctx);
489 dst_size =
FFMAX(
FFMIN(dst_size, delta_n), 0);
493 resample_func = (
c->linear && (
c->frac ||
c->dst_incr_mod)) ?
494 c->dsp.resample_linear :
c->dsp.resample_common;
495 for (
i = 0;
i < dst->ch_count;
i++)
496 *consumed = resample_func(
c, dst->ch[
i],
src->ch[
i], dst_size,
i+1 == dst->ch_count);
500 if (
c->compensation_distance) {
501 c->compensation_distance -= dst_size;
502 if (!
c->compensation_distance) {
503 c->dst_incr =
c->ideal_dst_incr;
504 c->dst_incr_div =
c->dst_incr /
c->src_incr;
505 c->dst_incr_mod =
c->dst_incr %
c->src_incr;
514 int64_t num =
s->in_buffer_count - (
c->filter_length-1)/2;
515 num *=
c->phase_count;
519 return av_rescale(num,
base,
s->in_sample_rate*(int64_t)
c->src_incr *
c->phase_count);
527 int64_t num =
s->in_buffer_count + 2LL + in_samples;
528 num *=
c->phase_count;
532 if (
c->compensation_distance) {
536 num =
FFMAX(num, (num *
c->ideal_dst_incr - 1) /
c->dst_incr + 1);
545 int reflection = (
FFMIN(
s->in_buffer_count,
c->filter_length) + 1) / 2;
550 for(
i=0;
i<
a->ch_count;
i++){
551 for(j=0; j<reflection; j++){
552 memcpy(
a->ch[
i] + (
s->in_buffer_index+
s->in_buffer_count+j )*
a->bps,
553 a->ch[
i] + (
s->in_buffer_index+
s->in_buffer_count-j-1)*
a->bps,
a->bps);
556 s->in_buffer_count += reflection;
562 int in_count,
int *out_idx,
int *out_sz)
564 int n, ch, num =
FFMIN(in_count + *out_sz,
c->filter_length + 1), res;
573 for (n = *out_sz; n < num; n++) {
574 for (ch = 0; ch <
src->ch_count; ch++) {
575 memcpy(dst->
ch[ch] + ((
c->filter_length + n) *
c->felem_size),
576 src->ch[ch] + ((n - *out_sz) *
c->felem_size),
c->felem_size);
583 *out_idx =
c->filter_length;
588 for (n = 1; n <=
c->filter_length; n++) {
589 for (ch = 0; ch <
src->ch_count; ch++) {
590 memcpy(dst->
ch[ch] + ((
c->filter_length - n) *
c->felem_size),
591 dst->
ch[ch] + ((
c->filter_length + n) *
c->felem_size),
597 *out_idx =
c->filter_length;
598 while (
c->index < 0) {
600 c->index +=
c->phase_count;
602 *out_sz =
FFMAX(*out_sz +
c->filter_length,
603 1 +
c->filter_length * 2) - *out_idx;
605 return FFMAX(res, 0);