ref: dce8fe1d141dd3912143feecdaebcd87374ef8e3
parent: 002ad34d50dbcbf7f410921d0a3f7aa14da4a02b
author: Rob Sykes <robs@users.sourceforge.net>
date: Tue Apr 17 15:41:10 EDT 2012
rate speed-ups: down-sample by 3/4, down-sample by > 4
--- a/src/effects_i_dsp.c
+++ b/src/effects_i_dsp.c
@@ -271,9 +271,10 @@
return h;
}
-int lsx_lpf_num_taps(double att, double tr_bw, int k)
+void lsx_kaiser_params(double att, double tr_bw, double * beta, int * num_taps)
{ /* TODO this could be cleaner, esp. for k != 0 */int n;
+ *beta = *beta < 0? lsx_kaiser_beta(att) : *beta;
if (att <= 80)
n = .25 / M_PI * (att - 7.95) / (2.285 * tr_bw) + .5;
else {@@ -280,33 +281,30 @@
double n160 = (.0425* att - 1.4) / tr_bw; /* Half order for att = 160 */
n = n160 * (16.556 / (att - 39.6) + .8625) + .5; /* For att [80,160) */
}
- return k? 2 * n * k - 1 : 2 * (n + (n & 1)) + 1; /* =1 %4 (0 phase 1/2 band) */
+ *num_taps = *num_taps? *num_taps : 2 * n;
}
double * lsx_design_lpf(
double Fp, /* End of pass-band; ~= 0.01dB point */
- double Fc, /* Start of stop-band */
+ double Fs, /* Start of stop-band */
double Fn, /* Nyquist freq; e.g. 0.5, 1, PI */
sox_bool allow_aliasing,
double att, /* Stop-band attenuation in dB */
- int * num_taps, /* (Single phase.) 0: value will be estimated */
- int k, /* Number of phases; 0 for single-phase */
- double beta)
+ int * num_taps, /* 0: value will be estimated */
+ int k, /* >0: number of phases; <0: num_taps ≡ 1 (mod -k) */
+ double beta) /* <0: value will be estimated */
{double tr_bw;
-
+ int n = *num_taps, phases = max(k, 1), modulo = max(-k, 1);
if (allow_aliasing)
- Fc += (Fc - Fp) * LSX_TO_3dB;
- Fp /= Fn, Fc /= Fn; /* Normalise to Fn = 1 */
- tr_bw = LSX_TO_6dB * (Fc-Fp); /* Transition band-width: 6dB to stop points */
-
- if (beta < 0)
- beta = lsx_kaiser_beta(att);
- if (!*num_taps)
- *num_taps = lsx_lpf_num_taps(att, tr_bw, k);
- k = max(k, 1);
- lsx_debug("design_lpf Fp=%g tr_bw=%g Fc=%g", Fp, tr_bw, Fc);- return lsx_make_lpf(*num_taps, (Fc - tr_bw) / k, beta, (double)k, sox_true);
+ Fs += (Fs - Fp) * LSX_TO_3dB;
+ Fp /= Fn, Fs /= Fn; /* Normalise to Fn = 1 */
+ tr_bw = LSX_TO_6dB * (Fs - Fp); /* Transition band-width: 6dB to stop points */
+ tr_bw /= phases, Fs /= phases;
+ lsx_kaiser_params(att, tr_bw, &beta, num_taps);
+ if (!n)
+ *num_taps = phases > 1? *num_taps / phases * phases + phases - 1 : (*num_taps + modulo - 2) / modulo * modulo + 1;
+ return lsx_make_lpf(*num_taps, Fs - tr_bw, beta, (double)phases, sox_true);
}
static double safe_log(double x)
--- a/src/rate.c
+++ b/src/rate.c
@@ -160,7 +160,7 @@
lsx_safe_rdft(f->dft_length, 1, output);
}
output[0] *= f->coefs[0];
- if (p->step.parts.integer) {+ if (p->step.parts.integer > 0) {output[1] *= f->coefs[1];
for (i = 2; i < f->dft_length; i += 2) {tmp = output[i];
@@ -177,25 +177,26 @@
else fifo_trim_by(output_fifo, overlap);
}
else { /* F-domain */- for (i = 2; i < (f->dft_length >> 1); i += 2) {+ int m = -p->step.parts.integer;
+ for (i = 2; i < (f->dft_length >> m); i += 2) {tmp = output[i];
output[i ] = f->coefs[i ] * tmp - f->coefs[i+1] * output[i+1];
output[i+1] = f->coefs[i+1] * tmp + f->coefs[i ] * output[i+1];
}
output[1] = f->coefs[i] * output[i] - f->coefs[i+1] * output[i+1];
- lsx_safe_rdft(f->dft_length >> 1, -1, output);
- fifo_trim_by(output_fifo, (f->dft_length + overlap) >> 1);
+ lsx_safe_rdft(f->dft_length >> m, -1, output);
+ fifo_trim_by(output_fifo, (((1 << m) - 1) * f->dft_length + overlap) >> m);
}
}
}
-static void setup_dft_stage(rate_shared_t * shared, int which, stage_t * stage, int L, int M)
+static void setup_dft_stage(rate_shared_t * shared, int which, stage_t * stage, int L, int M, sox_bool allow_aliasing)
{stage->fn = dft_stage_fn;
stage->preload = shared->dft_filter[which].post_peak / L;
stage->remL = shared->dft_filter[which].post_peak % L;
stage->L = L;
- stage->step.parts.integer = M;
+ stage->step.parts.integer = abs(3-M) == 1 && !allow_aliasing? -M/2 : M;
stage->which = which;
}
@@ -217,20 +218,12 @@
f->post_peak = num_taps / 2;
}
else {- double * h2 = lsx_design_lpf(Fp, Fc, Fn, allow_aliasing, att, &num_taps, 0, -1.);
+ int k = 4 << (phase == 50 && multiplier == 4 && Fn == 4);
+ double * h2 = lsx_design_lpf(Fp, Fc, Fn, allow_aliasing, att, &num_taps, -k, -1.);
if (phase != 50)
lsx_fir_to_phase(&h2, &num_taps, &f->post_peak, phase);
- else {- if (Fn == 4 && ((num_taps - 1) & 4)) { /* preserve phase */- double * h3 = calloc(num_taps + 4, sizeof(*h3));
- memcpy(h3 + 2, h2, num_taps * sizeof(*h3));
- free(h2);
- h2 = h3;
- num_taps += 4;
- }
- f->post_peak = num_taps / 2;
- }
+ else f->post_peak = num_taps / 2;
dft_length = lsx_set_dft_length(num_taps);
f->coefs = calloc(dft_length, sizeof(*f->coefs));
@@ -313,7 +306,7 @@
if (fracL > 2) {int L = fracL, M = realM;
for (i = level + 1; i && !(L & 1); L >>= 1, --i);
- if (L < 3 && (M <<= i) < 7) {+ if (((M <<= i) < 7 && L < 3) || M == 4) {preL = L, preM = M, realM = fracL = 1, level = 0, upsample = sox_true;
break;
}
@@ -395,7 +388,7 @@
if (preL * preM != 1) {init_dft_filter(shared, 0, 0, 0, bw, 1., (double)max(preL, preM), att, preL, phase, allow_aliasing);
- setup_dft_stage(shared, 0, &pre_stage, preL, preM == 2 && !allow_aliasing? 0 : preM);
+ setup_dft_stage(shared, 0, &pre_stage, preL, preM, allow_aliasing);
--p->input_stage_num;
}
else if (level && have_frac_stage && (1 - pass) / (1 - bw) > 2)
@@ -405,7 +398,7 @@
init_dft_filter(shared, 1, 0, 0,
bw * (upsample? factor * postL / postM : 1),
1., (double)(upsample? postL : postM), att, postL, phase, allow_aliasing);
- setup_dft_stage(shared, 1, &post_stage, postL, postM == 2 && !allow_aliasing? 0 : postM);
+ setup_dft_stage(shared, 1, &post_stage, postL, postM, allow_aliasing);
++p->output_stage_num;
}
}
@@ -413,12 +406,12 @@
stage_t * s = &p->stages[i];
if (i >= 0 && i < level - have_frac_stage) {s->fn = half_sample_25;
- s->pre_post = 2 * (array_length(half_fir_coefs_25) - 1);
+ s->pre_post = 4 * array_length(half_fir_coefs_25);
s->preload = s->pre = s->pre_post >> 1;
}
else if (level && i == level - 1) {if (shared->dft_filter[0].num_taps)
- setup_dft_stage(shared, 0, s, 1, (int)allow_aliasing << 1);
+ setup_dft_stage(shared, 0, s, 1, 2, allow_aliasing);
else *s = post_stage;
}
fifo_create(&s->fifo, (int)sizeof(sample_t));
--- a/src/rate_filters.h
+++ b/src/rate_filters.h
@@ -18,14 +18,10 @@
/* Generated by m4 */
static const sample_t half_fir_coefs_25[] = {- 4.9866643051942178e-001, 3.1333582318860204e-001, 1.2567743716165585e-003,
- -9.2035726038137103e-002, -1.0507348255277846e-003, 4.2764945027796687e-002,
- 7.7661461450703555e-004, -2.0673365323361139e-002, -5.0429677622613805e-004,
- 9.4223774565849357e-003, 2.8491539998284476e-004, -3.8562347294894628e-003,
- -1.3803431143314762e-004, 1.3634218103234187e-003, 5.6110366313398705e-005,
- -3.9872042837864422e-004, -1.8501044952475473e-005, 9.0580351350892191e-005,
- 4.6764104835321042e-006, -1.4284332593063177e-005, -8.1340436298087893e-007,
- 1.1833367010222812e-006, 7.3979325233687461e-008,
+ 3.1358440327836512e-01, -9.2701477245364594e-02, 4.3647483867630447e-02,
+ -2.1545228788689186e-02, 1.0119340890565588e-02, -4.3181204279612133e-03,
+ 1.6176661095102525e-03, -5.1348399782997947e-04, 1.3185858795078468e-04,
+ -2.5493512192147390e-05, 3.2461554199264636e-06, -1.9450196215470593e-07,
};
static const sample_t half_fir_coefs_low[] = {4.2759802493108773e-001, 3.0939308096100709e-001, 6.9285325719540158e-002,
@@ -45,10 +41,27 @@
-6.0893901283459912e-006, 1.4040363940567877e-005, 4.9834316581482789e-006,
};
#define FUNCTION half_sample_25
-#define CONVOLVE _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _
+#define CONVOLVE _ _ _ _ _ _ _ _ _ _ _ _
#define COEFS half_fir_coefs_25
-assert_static(!((array_length(COEFS)- 1) & 1), HALF_FIR_LENGTH_25 );
-#include "rate_half_fir.h"
+#define _ sum += (input[-(2*j +1)] + input[(2*j +1)]) * COEFS[j], ++j;
+static void FUNCTION(stage_t * p, fifo_t * output_fifo)
+{+ sample_t const * input = stage_read_p(p);
+ int i, num_out = (stage_occupancy(p) + 1) / 2;
+ sample_t * output = fifo_reserve(output_fifo, num_out);
+
+ for (i = 0; i < num_out; ++i, input += 2) {+ int j = 0;
+ sample_t sum = input[0] * .5;
+ CONVOLVE
+ output[i] = sum;
+ }
+ fifo_read(&p->fifo, 2 * num_out, NULL);
+}
+#undef _
+#undef COEFS
+#undef CONVOLVE
+#undef FUNCTION
#define FUNCTION half_sample_low
#define CONVOLVE _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _
#define COEFS half_fir_coefs_low
--- a/src/sinc.c
+++ b/src/sinc.c
@@ -80,14 +80,15 @@
static double * lpf(double Fn, double Fc, double tbw, int * num_taps, double att, double * beta, sox_bool round)
{+ int n = *num_taps;
if ((Fc /= Fn) <= 0 || Fc >= 1) {*num_taps = 0;
return NULL;
}
att = att? att : 120;
- *beta = *beta < 0? lsx_kaiser_beta(att) : *beta;
- if (!*num_taps) {- int n = lsx_lpf_num_taps(att, (tbw? tbw / Fn : .05) * .5, 0);
+ lsx_kaiser_params(att, (tbw? tbw / Fn : .05) * .5, beta, num_taps);
+ if (!n) {+ n = *num_taps;
*num_taps = range_limit(n, 11, 32767);
if (round)
*num_taps = 1 + 2 * (int)((int)((*num_taps / 2) * Fc + .5) / Fc + .5);
--- a/src/sox_i.h
+++ b/src/sox_i.h
@@ -101,16 +101,16 @@
double lsx_kaiser_beta(double att);
void lsx_apply_kaiser(double h[], const int num_points, double beta);
double * lsx_make_lpf(int num_taps, double Fc, double beta, double scale, sox_bool dc_norm);
-int lsx_lpf_num_taps(double att, double tr_bw, int k);
+void lsx_kaiser_params(double att, double tr_bw, double * beta, int * num_taps);
double * lsx_design_lpf(
double Fp, /* End of pass-band; ~= 0.01dB point */
- double Fc, /* Start of stop-band */
+ double Fs, /* Start of stop-band */
double Fn, /* Nyquist freq; e.g. 0.5, 1, PI */
sox_bool allow_aliasing,
double att, /* Stop-band attenuation in dB */
- int * num_taps, /* (Single phase.) 0: value will be estimated */
- int k, /* Number of phases; 0 for single-phase */
- double beta);
+ int * num_taps, /* 0: value will be estimated */
+ int k, /* >0: number of phases; <0: num_taps ≡ 1 (mod -k) */
+ double beta); /* <0: value will be estimated */
void lsx_fir_to_phase(double * * h, int * len,
int * post_len, double phase0);
#define LSX_TO_6dB .5869