- Increase dynamic range of filter coefficients from 28bit to 30bit.

This cause dramatic effect in overall precision and conversion quality
  by pushing down most aliasing artifacts around -180 dB.

  Spectrogram analysis/comparison:

  	http://people.freebsd.org/~ariff/z_comparison/z_28vs30/

- Guard against possible 64bit overflow during accumulation process by
  slightly normalize and saturate sample and coefficient multiplication,
  possible during extreme 32bit downsampling (eg. 380KHz -> 8KHz) with
  custom preset that require more than ~7000 taps filter (which is
  overkill).

- Add knobs through FEEDER_RATE_PRESETS to set dynamic range of filter
  coefficients/accumulator and prefered polynomial interpolator:

  	COEFFICIENT_BIT:X
	(where 1 <= X <= 30, default: 30)

	ACCUMULATOR_BIT:X
	(where 32 <= X <=64, default: 58)

	INTERPOLATOR:I
	(where I = ZOH, LINEAR, QUADRATIC, HERMITE, BSPLINE,
 	           OPT32X, OPT16X, OPT8X, OPT4X, OPT2X)

Approved by:	re (kib)
This commit is contained in:
Ariff Abdullah 2009-07-05 18:15:06 +00:00
parent 7634012302
commit 96831ae59b
Notes: svn2git 2020-12-20 02:59:44 +00:00
svn path=/head/; revision=195378
2 changed files with 120 additions and 13 deletions

View File

@ -527,7 +527,7 @@ z_feed_linear_##SIGN##BIT##ENDIAN(struct z_info *info, uint8_t *dst) \
z &= Z_MASK; \
coeff = Z_COEFF_INTERPOLATE(z, z_coeff[c], z_dcoeff[c]); \
x = _PCM_READ_##SIGN##BIT##_##ENDIAN(p); \
v += (intpcm64_t)x * coeff; \
v += Z_NORM_##BIT((intpcm64_t)x * coeff); \
z += info->z_dy; \
p adv##= info->channels * PCM_##BIT##_BPS
@ -582,7 +582,7 @@ z_feed_sinc_##SIGN##BIT##ENDIAN(struct z_info *info, uint8_t *dst) \
if (info->z_scale != Z_ONE) \
v = Z_SCALE_##BIT(v, info->z_scale); \
else \
v >>= Z_COEFF_SHIFT; \
v >>= Z_COEFF_SHIFT - Z_GUARD_BIT_##BIT; \
Z_CLIP_CHECK(v, BIT); \
_PCM_WRITE_##SIGN##BIT##_##ENDIAN(dst, Z_CLAMP(v, BIT)); \
} while (ch != 0); \
@ -610,18 +610,18 @@ z_feed_sinc_polyphase_##SIGN##BIT##ENDIAN(struct z_info *info, uint8_t *dst) \
((info->z_alpha * info->z_size) << 1); \
for (i = info->z_size; i != 0; i--) { \
x = _PCM_READ_##SIGN##BIT##_##ENDIAN(p); \
v += (intpcm64_t)x * *z_pcoeff; \
v += Z_NORM_##BIT((intpcm64_t)x * *z_pcoeff); \
z_pcoeff++; \
p += info->channels * PCM_##BIT##_BPS; \
x = _PCM_READ_##SIGN##BIT##_##ENDIAN(p); \
v += (intpcm64_t)x * *z_pcoeff; \
v += Z_NORM_##BIT((intpcm64_t)x * *z_pcoeff); \
z_pcoeff++; \
p += info->channels * PCM_##BIT##_BPS; \
} \
if (info->z_scale != Z_ONE) \
v = Z_SCALE_##BIT(v, info->z_scale); \
else \
v >>= Z_COEFF_SHIFT; \
v >>= Z_COEFF_SHIFT - Z_GUARD_BIT_##BIT; \
Z_CLIP_CHECK(v, BIT); \
_PCM_WRITE_##SIGN##BIT##_##ENDIAN(dst, Z_CLAMP(v, BIT)); \
} while (ch != 0); \
@ -882,7 +882,7 @@ z_coeff_interpolate(int32_t z, int32_t *z_coeff)
zoo2 = z_coeff[2] - z_coeff[-1];
zoo3 = z_coeff[3] - z_coeff[-2];
zoc0 = (((0x1ac2260dLL * zoe1)) >> 30) +
zoc0 = ((0x1ac2260dLL * zoe1) >> 30) +
((0x0526cdcaLL * zoe2) >> 30) + ((0x00170c29LL * zoe3) >> 30);
zoc1 = ((0x14f8a49aLL * zoo1) >> 30) +
((0x0d6d1109LL * zoo2) >> 30) + ((0x008cd4dcLL * zoo3) >> 30);
@ -911,7 +911,7 @@ z_coeff_interpolate(int32_t z, int32_t *z_coeff)
zoo2 = z_coeff[2] - z_coeff[-1];
zoo3 = z_coeff[3] - z_coeff[-2];
zoc0 = (((0x1ac2260dLL * zoe1)) >> 30) +
zoc0 = ((0x1ac2260dLL * zoe1) >> 30) +
((0x0526cdcaLL * zoe2) >> 30) + ((0x00170c29LL * zoe3) >> 30);
zoc1 = ((0x14f8a49aLL * zoo1) >> 30) +
((0x0d6d1109LL * zoo2) >> 30) + ((0x008cd4dcLL * zoo3) >> 30);
@ -940,7 +940,7 @@ z_coeff_interpolate(int32_t z, int32_t *z_coeff)
zoo2 = z_coeff[2] - z_coeff[-1];
zoo3 = z_coeff[3] - z_coeff[-2];
zoc0 = (((0x1aa9b47dLL * zoe1)) >> 30) +
zoc0 = ((0x1aa9b47dLL * zoe1) >> 30) +
((0x053d9944LL * zoe2) >> 30) + ((0x0018b23fLL * zoe3) >> 30);
zoc1 = ((0x14a104d1LL * zoo1) >> 30) +
((0x0d7d2504LL * zoo2) >> 30) + ((0x0094b599LL * zoo3) >> 30);
@ -969,7 +969,7 @@ z_coeff_interpolate(int32_t z, int32_t *z_coeff)
zoo2 = z_coeff[2] - z_coeff[-1];
zoo3 = z_coeff[3] - z_coeff[-2];
zoc0 = (((0x1a8eda43LL * zoe1)) >> 30) +
zoc0 = ((0x1a8eda43LL * zoe1) >> 30) +
((0x0556ee38LL * zoe2) >> 30) + ((0x001a3784LL * zoe3) >> 30);
zoc1 = ((0x143d863eLL * zoo1) >> 30) +
((0x0d910e36LL * zoo2) >> 30) + ((0x009ca889LL * zoo3) >> 30);
@ -998,7 +998,7 @@ z_coeff_interpolate(int32_t z, int32_t *z_coeff)
zoo2 = z_coeff[2] - z_coeff[-1];
zoo3 = z_coeff[3] - z_coeff[-2];
zoc0 = (((0x19edb6fdLL * zoe1)) >> 30) +
zoc0 = ((0x19edb6fdLL * zoe1) >> 30) +
((0x05ebd062LL * zoe2) >> 30) + ((0x00267881LL * zoe3) >> 30);
zoc1 = ((0x1223af76LL * zoo1) >> 30) +
((0x0de3dd6bLL * zoo2) >> 30) + ((0x00d683cdLL * zoo3) >> 30);

View File

@ -386,6 +386,27 @@ function filter_parse(s, a, i, attn, alen)
return (-1);
}
if (alen > 0 && a[1] == "COEFFICIENT_BIT") {
if (alen != 2)
return (-1);
init_coeff_bit(floor(a[2]));
return (-1);
}
if (alen > 0 && a[1] == "ACCUMULATOR_BIT") {
if (alen != 2)
return (-1);
init_accum_bit(floor(a[2]));
return (-1);
}
if (alen > 0 && a[1] == "INTERPOLATOR") {
if (alen != 2)
return (-1);
init_coeff_interpolator(toupper(a[2]));
return (-1);
}
if (alen == 1 || alen == 2) {
if (a[1] == "NYQUIST_HOVER") {
i = 1.0 * a[2];
@ -448,7 +469,12 @@ function filter_parse(s, a, i, attn, alen)
function genscale(bit, s1, s2, scale)
{
s1 = Z_COEFF_SHIFT - (32 - bit);
if ((bit + Z_COEFF_SHIFT) > Z_ACCUMULATOR_BIT)
s1 = Z_COEFF_SHIFT - \
(32 - (Z_ACCUMULATOR_BIT - Z_COEFF_SHIFT));
else
s1 = Z_COEFF_SHIFT - (32 - bit);
s2 = Z_SHIFT + (32 - bit);
if (s1 == 0)
@ -527,6 +553,62 @@ function init_drift(drift, xdrift)
Z_MASK = Z_ONE - 1;
}
function init_coeff_bit(cbit, xcbit)
{
xcbit = floor(cbit);
if (Z_COEFF_SHIFT != 0) {
if (xcbit != Z_COEFF_SHIFT)
printf("#error Z_COEFF_SHIFT reinitialize!\n");
return;
}
#
# Initialize dynamic range of coefficients.
#
if (xcbit < 1)
xcbit = 1;
else if (xcbit > 30)
xcbit = 30;
Z_COEFF_SHIFT = xcbit;
Z_COEFF_ONE = shl(1, Z_COEFF_SHIFT);
}
function init_accum_bit(accbit, xaccbit)
{
xaccbit = floor(accbit);
if (Z_ACCUMULATOR_BIT != 0) {
if (xaccbit != Z_ACCUMULATOR_BIT)
printf("#error Z_ACCUMULATOR_BIT reinitialize!\n");
return;
}
#
# Initialize dynamic range of accumulator.
#
if (xaccbit > 64)
xaccbit = 64;
else if (xaccbit < 32)
xaccbit = 32;
Z_ACCUMULATOR_BIT = xaccbit;
}
function init_coeff_interpolator(interp)
{
#
# Validate interpolator type.
#
if (interp == "ZOH" || interp == "LINEAR" || \
interp == "QUADRATIC" || interp == "HERMITE" || \
interp == "BSPLINE" || interp == "OPT32X" || \
interp == "OPT16X" || interp == "OPT8X" || \
interp == "OPT4X" || interp == "OPT2X")
Z_COEFF_INTERPOLATOR = interp;
}
BEGIN {
I0_EPSILON = 1e-21;
M_PI = atan2(0.0, -1.0);
@ -536,11 +618,17 @@ BEGIN {
Z_COEFF_OFFSET = 5;
Z_ACCUMULATOR_BIT_DEFAULT = 58;
Z_ACCUMULATOR_BIT = 0;
Z_FULL_SHIFT = 30;
Z_FULL_ONE = shl(1, Z_FULL_SHIFT);
Z_COEFF_SHIFT = 28;
Z_COEFF_ONE = shl(1, Z_COEFF_SHIFT);
Z_COEFF_SHIFT_DEFAULT = 30;
Z_COEFF_SHIFT = 0;
Z_COEFF_ONE = 0;
Z_COEFF_INTERPOLATOR = 0;
Z_INTERP_COEFF_SHIFT = 24;
Z_INTERP_COEFF_ONE = shl(1, Z_INTERP_COEFF_SHIFT);
@ -620,6 +708,10 @@ BEGIN {
rolloff = Popts["rolloff"];
if (Z_DRIFT_SHIFT == -1)
init_drift(Z_DRIFT_SHIFT_DEFAULT);
if (Z_COEFF_SHIFT == 0)
init_coeff_bit(Z_COEFF_SHIFT_DEFAULT);
if (Z_ACCUMULATOR_BIT == 0)
init_accum_bit(Z_ACCUMULATOR_BIT_DEFAULT);
ztab[imp["quality"] - 2] = \
mkfilter(imp, nmult, rolloff, beta, Z_DRIFT_ONE);
imp["quality"]++;
@ -751,6 +843,18 @@ BEGIN {
genscale(24);
genscale(32);
printf("\n");
printf("#define Z_ACCUMULATOR_BIT\t%d\n\n", Z_ACCUMULATOR_BIT)
for (i = 8; i <= 32; i += 8) {
gbit = ((i + Z_COEFF_SHIFT) > Z_ACCUMULATOR_BIT) ? \
(i - (Z_ACCUMULATOR_BIT - Z_COEFF_SHIFT)) : 0;
printf("#define Z_GUARD_BIT_%d\t\t%d\n", i, gbit);
if (gbit == 0)
printf("#define Z_NORM_%d(v)\t\t(v)\n\n", i);
else
printf("#define Z_NORM_%d(v)\t\t" \
"((v) >> Z_GUARD_BIT_%d)\n\n", i, i);
}
printf("\n");
printf("#define Z_LINEAR_FULL_ONE\t0x%08xU\n", Z_LINEAR_FULL_ONE);
printf("#define Z_LINEAR_SHIFT\t\t%d\n", Z_LINEAR_SHIFT);
printf("#define Z_LINEAR_UNSHIFT\t%d\n", Z_LINEAR_UNSHIFT);
@ -775,6 +879,9 @@ BEGIN {
printf("\n");
printf("#define Z_QUALITY_MIN\t\t0\n");
printf("#define Z_QUALITY_MAX\t\t%d\n", length(ztab) + 1);
if (Z_COEFF_INTERPOLATOR != 0)
printf("\n#define Z_COEFF_INTERP_%s\t1\n", \
Z_COEFF_INTERPOLATOR);
printf("\n/*\n * smallest: %.32f\n * largest: %.32f\n *\n", \
smallest, largest);
printf(" * z_unshift=%d, z_interp_shift=%d\n *\n", \