diff --git a/celt/bands.c b/celt/bands.c index aea37837..13ed576a 100644 --- a/celt/bands.c +++ b/celt/bands.c @@ -106,18 +106,35 @@ oac_int16 oaci_bitexact_cos(oac_int16 x) { return 1 + x2; } -int oaci_bitexact_log2tan(int isin, int icos) { - int lc; - int ls; - lc = EC_ILOG(icos); - ls = EC_ILOG(isin); - icos <<= 15 - lc; - isin <<= 15 - ls; - return (ls - lc)*(1<<11) - + FRAC_MUL16(isin, FRAC_MUL16(isin, -2597) + 7932) - - FRAC_MUL16(icos, FRAC_MUL16(icos, -2597) + 7932); +/* This is a log2() approximation designed to be bit-exact on any platform. Bit exactness + with this approximation is important because it has an impact on the bit allocation */ +static oac_int16 oaci_bitexact_log2(oac_int32 x) { + int i; + oac_int16 n, frac; + /* -0.41509302963303146, 0.9609890551383969, -0.31836011537636605, + 0.15530808010959576, -0.08556153059057618 */ + static const oac_int16 C[5] = {-6801 + (1<<(14 - 11)), 15746, -5217, 2545, -1401}; + if (x == 0) + return -32767; + i = EC_ILOG(x)-1; + n = (x << 15 >> i) - 32768 - 16384; + frac = (C[0] + FRAC_MUL16(n, (C[1] + FRAC_MUL16(n, (C[2] + FRAC_MUL16(n, (C[3] + FRAC_MUL16(n, C[4])))))))); + return ((i - 12) << 11) + (frac >> (14 - 11)); +} + +/* This is a log2(tan(x)) approximation designed to be bit-exact on any platform. Bit exactness + with this approximation is important because it has an impact on the bit allocation */ +int oaci_bitexact_delta(int itheta) { + int sign = 1; + if (itheta > 8192) { + itheta = 16384 - itheta; + sign = -1; + } + /* log2(tan(pi/4*x)) ~= log2(x) -.343 - .061*x + .404*x.^2. */ + return sign*(oaci_bitexact_log2(itheta) - 703 - FRAC_MUL16(itheta, (500 - FRAC_MUL16(itheta, 13238)))); } + #ifdef FIXED_POINT /* Compute the amplitude (sqrt energy) in each of the bands */ void oaci_compute_band_energies(const CELTMode *m, const celt_sig *X, celt_ener *bandE, int end, int C, int LM, int arch) { @@ -624,7 +641,7 @@ void oaci_haar1(celt_norm *X, int N0, int stride) { } } -static int oaci_compute_qn(int N, int b, int offset, int pulse_cap, int stereo) { +static int oaci_compute_qn(int N, oac_int32 b, int offset, int pulse_cap, int stereo) { static const oac_int16 exp2_table8[8] = {16384, 17866, 19483, 21247, 23170, 25267, 27554, 30048}; int qn, qb; @@ -669,22 +686,19 @@ struct band_ctx { struct split_ctx { int inv; - int imid; - int iside; - int delta; + oac_int32 delta; int itheta; int itheta_q30; }; static void oaci_compute_theta(struct band_ctx *ctx, struct split_ctx *sctx, - celt_norm *X, celt_norm *Y, int N, int *b, int B, int B0, + celt_norm *X, celt_norm *Y, int N, oac_int32 *b, int B, int B0, int LM, int stereo, int *fill) { int qn; int itheta = 0; int itheta_q30 = 0; int delta; - int imid, iside; int qalloc; int pulse_cap; int offset; @@ -728,9 +742,7 @@ static void oaci_compute_theta(struct band_ctx *ctx, struct split_ctx *sctx, to inject noise on one side. If so, make sure the energy of that side is zero. */ int unquantized = oaci_celt_udiv((oac_int32)itheta*16384, qn); - imid = oaci_bitexact_cos((oac_int16)unquantized); - iside = oaci_bitexact_cos((oac_int16)(16384 - unquantized)); - delta = FRAC_MUL16((N - 1)<<7, oaci_bitexact_log2tan(iside, imid)); + delta = FRAC_MUL16((N - 1)<<7, oaci_bitexact_delta(unquantized)); if (delta > *b) itheta = qn; else if (delta < -*b) @@ -843,26 +855,18 @@ static void oaci_compute_theta(struct band_ctx *ctx, struct split_ctx *sctx, *b -= qalloc; if (itheta == 0) { - imid = 32767; - iside = 0; *fill &= (1<inv = inv; - sctx->imid = imid; - sctx->iside = iside; sctx->delta = delta; sctx->itheta = itheta; sctx->itheta_q30 = itheta_q30; @@ -909,7 +913,6 @@ static unsigned oaci_quant_partition(struct band_ctx *ctx, celt_norm *X, const unsigned char *cache; int q; int curr_bits; - int imid = 0, iside = 0; int B0 = B; oac_val32 mid = 0, side = 0; unsigned cm = 0; @@ -944,12 +947,8 @@ static unsigned oaci_quant_partition(struct band_ctx *ctx, celt_norm *X, B = (B + 1)>>1; oaci_compute_theta(ctx, &sctx, X, Y, N, &b, B, B0, LM, 0, &fill); - imid = sctx.imid; - iside = sctx.iside; delta = sctx.delta; itheta = sctx.itheta; - (void)imid; - (void)iside; #ifdef FIXED_POINT mid = oaci_celt_cos_norm32(sctx.itheta_q30); side = oaci_celt_cos_norm32((1<<30) - sctx.itheta_q30); @@ -1124,7 +1123,7 @@ unsigned oaci_cubic_quant_partition(struct band_ctx *ctx, celt_norm *X, int N, i /* This function is responsible for encoding and decoding a band for the mono case. */ static unsigned oaci_quant_band(struct band_ctx *ctx, celt_norm *X, - int N, int b, int B, celt_norm *lowband, + int N, oac_int32 b, int B, celt_norm *lowband, int LM, celt_norm *lowband_out, oac_val32 gain, celt_norm *lowband_scratch, int fill) { int N0 = N; @@ -1245,10 +1244,9 @@ static unsigned oaci_quant_band(struct band_ctx *ctx, celt_norm *X, /* This function is responsible for encoding and decoding a band for the stereo case. */ static unsigned oaci_quant_band_stereo(struct band_ctx *ctx, celt_norm *X, celt_norm *Y, - int N, int b, int B, celt_norm *lowband, + int N, oac_int32 b, int B, celt_norm *lowband, int LM, celt_norm *lowband_out, celt_norm *lowband_scratch, int fill) { - int imid = 0, iside = 0; int inv = 0; oac_val32 mid = 0, side = 0; unsigned cm = 0; @@ -1277,12 +1275,8 @@ static unsigned oaci_quant_band_stereo(struct band_ctx *ctx, celt_norm *X, celt_ } oaci_compute_theta(ctx, &sctx, X, Y, N, &b, B, B, LM, 1, &fill); inv = sctx.inv; - imid = sctx.imid; - iside = sctx.iside; delta = sctx.delta; itheta = sctx.itheta; - (void)imid; - (void)iside; #ifdef FIXED_POINT mid = oaci_celt_cos_norm32(sctx.itheta_q30); side = oaci_celt_cos_norm32((1<<30) - sctx.itheta_q30); @@ -1482,7 +1476,7 @@ void oaci_quant_all_bands(int encode, const CELTMode *m, int start, int end, ctx.avoid_split_noise = B > 1; for (i = start; i < end; i++) { oac_int32 tell; - int b; + oac_int32 b; int N; oac_int32 curr_balance; int effective_lowband = -1; @@ -1511,7 +1505,7 @@ void oaci_quant_all_bands(int encode, const CELTMode *m, int start, int end, ctx.total_bits = total_bits; if (i <= codedBands - 1) { curr_balance = oaci_celt_sudiv(balance, IMIN(3, codedBands - i)); - b = IMAX(0, IMIN(16383, IMIN(remaining_bits, pulses[i] + curr_balance))); + b = IMAX(0, IMIN(32767, IMIN(remaining_bits, pulses[i] + curr_balance))); } else { b = 0; } diff --git a/celt/bands.h b/celt/bands.h index 8d9cb90a..829c83df 100644 --- a/celt/bands.h +++ b/celt/bands.h @@ -69,8 +69,7 @@ #include "entdec.h" #include "rate.h" -oac_int16 oaci_bitexact_cos(oac_int16 x); -int oaci_bitexact_log2tan(int isin, int icos); +int oaci_bitexact_delta(int itheta); /** Compute the amplitude (sqrt energy) in each of the bands * @param m Mode data diff --git a/celt/tests/test_unit_mathops.c b/celt/tests/test_unit_mathops.c index e69b09cf..2d170e50 100644 --- a/celt/tests/test_unit_mathops.c +++ b/celt/tests/test_unit_mathops.c @@ -119,48 +119,23 @@ void testsqrt(void) { } } -void testbitexactcos(void) { - int i; - oac_int32 min_d, max_d, last, chk; - chk = max_d = 0; - last = min_d = 32767; - for (i = 64; i <= 16320; i++) { - oac_int32 d; - oac_int32 q = oaci_bitexact_cos(i); - chk ^= q*i; - d = last - q; - if (d > max_d) max_d = d; - if (d < min_d) min_d = d; - last = q; - } - if ((chk != 89408644) || (max_d != 5) || (min_d != 0) || (oaci_bitexact_cos(64) != 32767) - || (oaci_bitexact_cos(16320) != 200) || (oaci_bitexact_cos(8192) != 23171)) { - fprintf (stderr, "oaci_bitexact_cos failed\n"); - ret = 1; - } +#ifndef DISABLE_FLOAT_API +static int oaci_delta_ground_truth(int itheta) { + return 2048*log2(tan(itheta/8192.*M_PI/4)); } +#endif void testbitexactlog2tan(void) { int i, fail; - oac_int32 min_d, max_d, last, chk; - fail = chk = max_d = 0; - last = min_d = 15059; - for (i = 64; i < 8193; i++) { - oac_int32 d; - oac_int32 mid = oaci_bitexact_cos(i); - oac_int32 side = oaci_bitexact_cos(16384 - i); - oac_int32 q = oaci_bitexact_log2tan(mid, side); - chk ^= q*i; - d = last - q; - if (q != -1*oaci_bitexact_log2tan(side, mid)) - fail = 1; - if (d > max_d) max_d = d; - if (d < min_d) min_d = d; - last = q; - } - if ((chk != 15821257) || (max_d != 61) || (min_d != -2) || fail - || (oaci_bitexact_log2tan(32767, 200) != 15059) || (oaci_bitexact_log2tan(30274, 12540) != 2611) - || (oaci_bitexact_log2tan(23171, 23171) != 0)) { + fail = 0; + for (i = 1; i <= 8192; i++) { + oac_int32 q = oaci_bitexact_delta(i); +#ifndef DISABLE_FLOAT_API + if (abs(q - oaci_delta_ground_truth(i)) > 16) fail=1; +#endif + if (q != -1*oaci_bitexact_delta(16384-i)) fail = 1; + } + if (fail) { fprintf (stderr, "oaci_bitexact_log2tan failed\n"); ret = 1; } @@ -722,7 +697,6 @@ int main(void) { int i; int use_ref_impl[2] = { 0, 1 }; - testbitexactcos(); testbitexactlog2tan(); testdiv(); testsqrt();