From 19ae9fc9446067cb119b66a77c13a90e396f5ef1 Mon Sep 17 00:00:00 2001 From: Jean-Marc Valin Date: Thu, 13 Mar 2008 11:18:15 +1100 Subject: [PATCH] fixed-point: simplifying the arithmetic in alg_quant() --- libcelt/fixed_generic.h | 2 ++ libcelt/mathops.h | 6 ++---- libcelt/vq.c | 20 ++++++++++---------- 3 files changed, 14 insertions(+), 14 deletions(-) diff --git a/libcelt/fixed_generic.h b/libcelt/fixed_generic.h index f0490a12..44b1f795 100644 --- a/libcelt/fixed_generic.h +++ b/libcelt/fixed_generic.h @@ -41,6 +41,8 @@ #define MULT32_32_Q31(a,b) ADD32(ADD32(SHL(MULT16_16(SHR((a),16),SHR((b),16)),1), SHR(MULT16_16SU(SHR((a),16),((b)&0x0000ffff)),15)), SHR(MULT16_16SU(SHR((b),16),((a)&0x0000ffff)),15)) +#define MULT32_32_Q32(a,b) ADD32(ADD32(MULT16_16(SHR((a),16),SHR((b),16)), SHR(MULT16_16SU(SHR((a),16),((b)&0x0000ffff)),16)), SHR(MULT16_16SU(SHR((b),16),((a)&0x0000ffff)),16)) + #define QCONST16(x,bits) ((celt_word16_t)(.5+(x)*(((celt_word32_t)1)<<(bits)))) #define QCONST32(x,bits) ((celt_word32_t)(.5+(x)*(((celt_word32_t)1)<<(bits)))) diff --git a/libcelt/mathops.h b/libcelt/mathops.h index 94cd2ba0..726d7188 100644 --- a/libcelt/mathops.h +++ b/libcelt/mathops.h @@ -164,7 +164,7 @@ static inline celt_word32_t celt_exp2(celt_word16_t x) return VSHR32(EXTEND32(frac), -integer-2); } -static inline celt_word32_t celt_rcp(celt_word16_t x) +static inline celt_word32_t celt_rcp(celt_word32_t x) { int i, neg=0; celt_word16_t n, frac; @@ -174,15 +174,13 @@ static inline celt_word32_t celt_rcp(celt_word16_t x) neg = 1; x = NEG16(x); } - if (x==0) - return 0; i = celt_ilog2(x); n = VSHR32(x,i-16)-SHL32(EXTEND32(3),15); frac = ADD16(C[0], MULT16_16_Q15(n, ADD16(C[1], MULT16_16_Q15(n, ADD16(C[2], MULT16_16_Q15(n, ADD16(C[3], MULT16_16_Q15(n, (C[4]))))))))); if (neg) frac = -frac; - return SHL32(EXTEND32(frac),16-i); + return VSHR32(EXTEND32(frac),i-16); } #endif /* FIXED_POINT */ diff --git a/libcelt/vq.c b/libcelt/vq.c index f265a886..7c89d02d 100644 --- a/libcelt/vq.c +++ b/libcelt/vq.c @@ -218,24 +218,24 @@ void alg_quant(celt_norm_t *X, celt_mask_t *W, int N, int K, const celt_norm_t * if (iy[m][j]*sign < 0) continue; - spj = MULT16_16_P14(s, P[j]); - aspj = MULT16_16_P15(alpha, spj); + spj = MULT16_16_Q14(s, P[j]); + aspj = MULT16_16_Q15(alpha, spj); /* Updating the sums of the new pulse(s) */ - Rxy = xy[m] + MULT16_16(s,X[j]) - MULT16_16(MULT16_16_P15(alpha,spj),Rxp); + Rxy = xy[m] + MULT16_16(s,X[j]) - MULT16_16(MULT16_16_Q15(alpha,spj),Rxp); Ryy = yy[m] + 2*MULT16_16(s,y[m][j]) + MULT16_16(s,s) +MULT16_16(aspj,MULT16_16_Q14(aspj,Rpp)) - 2*MULT16_32_Q14(aspj,yp[m]) - 2*MULT16_16(s,MULT16_16_Q14(aspj,P[j])); Ryp = yp[m] + MULT16_16(spj, SUB16(QCONST16(1.f,14),MULT16_16_Q15(alpha,Rpp))); /* Compute the gain such that ||p + g*y|| = 1 */ - g = MULT32_32_Q31( - SHL32(celt_sqrt(MULT16_16(ROUND(Ryp,14),ROUND(Ryp,14)) + Ryy - - MULT16_16(ROUND(Ryy,14),Rpp)) - - ROUND(Ryp,14), 14), - celt_rcp(ROUND(Ryy,14))); + g = MULT16_32_Q15( + celt_sqrt(MULT16_16(ROUND(Ryp,14),ROUND(Ryp,14)) + Ryy - + MULT16_16(ROUND(Ryy,14),Rpp)) + - ROUND(Ryp,14), + celt_rcp(SHR32(Ryy,12))); /* Knowing that gain, what's the error: (x-g*y)^2 (result is negated and we discard x^2 because it's constant) */ /*score = 2.f*g*Rxy - 1.f*g*g*Ryy*NORM_SCALING_1;*/ - score = 2*MULT16_32_Q14(ROUND(Rxy,14),g) - - MULT16_32_Q14(EXTRACT16(MULT16_32_Q14(ROUND(Ryy,14),g)),g); + score = 2*MULT16_32_Q14(ROUND(Rxy,14),g) + - MULT16_32_Q14(EXTRACT16(MULT16_32_Q14(ROUND(Ryy,14),g)),g); if (score>nbest[Lupdate-1]->score) {