diff options
Diffstat (limited to 'gr-vocoder/lib/gsm/rpe.c')
-rw-r--r-- | gr-vocoder/lib/gsm/rpe.c | 488 |
1 files changed, 488 insertions, 0 deletions
diff --git a/gr-vocoder/lib/gsm/rpe.c b/gr-vocoder/lib/gsm/rpe.c new file mode 100644 index 000000000..8a6b81fae --- /dev/null +++ b/gr-vocoder/lib/gsm/rpe.c @@ -0,0 +1,488 @@ +/* + * Copyright 1992 by Jutta Degener and Carsten Bormann, Technische + * Universitaet Berlin. See the accompanying file "COPYRIGHT" for + * details. THERE IS ABSOLUTELY NO WARRANTY FOR THIS SOFTWARE. + */ + +/* $Header$ */ + +#include <stdio.h> +#include <assert.h> + +#include "private.h" + +#include "gsm.h" +#include "proto.h" + +/* 4.2.13 .. 4.2.17 RPE ENCODING SECTION + */ + +/* 4.2.13 */ + +static void Weighting_filter P2((e, x), + register word * e, /* signal [-5..0.39.44] IN */ + word * x /* signal [0..39] OUT */ +) +/* + * The coefficients of the weighting filter are stored in a table + * (see table 4.4). The following scaling is used: + * + * H[0..10] = integer( real_H[ 0..10] * 8192 ); + */ +{ + /* word wt[ 50 ]; */ + + register longword L_result; + register int k /* , i */ ; + + /* Initialization of a temporary working array wt[0...49] + */ + + /* for (k = 0; k <= 4; k++) wt[k] = 0; + * for (k = 5; k <= 44; k++) wt[k] = *e++; + * for (k = 45; k <= 49; k++) wt[k] = 0; + * + * (e[-5..-1] and e[40..44] are allocated by the caller, + * are initially zero and are not written anywhere.) + */ + e -= 5; + + /* Compute the signal x[0..39] + */ + for (k = 0; k <= 39; k++) { + + L_result = 8192 >> 1; + + /* for (i = 0; i <= 10; i++) { + * L_temp = GSM_L_MULT( wt[k+i], gsm_H[i] ); + * L_result = GSM_L_ADD( L_result, L_temp ); + * } + */ + +#undef STEP +#define STEP( i, H ) (e[ k + i ] * (longword)H) + + /* Every one of these multiplications is done twice -- + * but I don't see an elegant way to optimize this. + * Do you? + */ + +#ifdef STUPID_COMPILER + L_result += STEP( 0, -134 ) ; + L_result += STEP( 1, -374 ) ; + /* + STEP( 2, 0 ) */ + L_result += STEP( 3, 2054 ) ; + L_result += STEP( 4, 5741 ) ; + L_result += STEP( 5, 8192 ) ; + L_result += STEP( 6, 5741 ) ; + L_result += STEP( 7, 2054 ) ; + /* + STEP( 8, 0 ) */ + L_result += STEP( 9, -374 ) ; + L_result += STEP( 10, -134 ) ; +#else + L_result += + STEP( 0, -134 ) + + STEP( 1, -374 ) + /* + STEP( 2, 0 ) */ + + STEP( 3, 2054 ) + + STEP( 4, 5741 ) + + STEP( 5, 8192 ) + + STEP( 6, 5741 ) + + STEP( 7, 2054 ) + /* + STEP( 8, 0 ) */ + + STEP( 9, -374 ) + + STEP(10, -134 ) + ; +#endif + + /* L_result = GSM_L_ADD( L_result, L_result ); (* scaling(x2) *) + * L_result = GSM_L_ADD( L_result, L_result ); (* scaling(x4) *) + * + * x[k] = SASR( L_result, 16 ); + */ + + /* 2 adds vs. >>16 => 14, minus one shift to compensate for + * those we lost when replacing L_MULT by '*'. + */ + + L_result = SASR( L_result, 13 ); + x[k] = ( L_result < MIN_WORD ? MIN_WORD + : (L_result > MAX_WORD ? MAX_WORD : L_result )); + } +} + +/* 4.2.14 */ + +static void RPE_grid_selection P3((x,xM,Mc_out), + word * x, /* [0..39] IN */ + word * xM, /* [0..12] OUT */ + word * Mc_out /* OUT */ +) +/* + * The signal x[0..39] is used to select the RPE grid which is + * represented by Mc. + */ +{ + /* register word temp1; */ + register int /* m, */ i; + register longword L_result, L_temp; + longword EM; /* xxx should be L_EM? */ + word Mc; + + longword L_common_0_3; + + EM = 0; + Mc = 0; + + /* for (m = 0; m <= 3; m++) { + * L_result = 0; + * + * + * for (i = 0; i <= 12; i++) { + * + * temp1 = SASR( x[m + 3*i], 2 ); + * + * assert(temp1 != MIN_WORD); + * + * L_temp = GSM_L_MULT( temp1, temp1 ); + * L_result = GSM_L_ADD( L_temp, L_result ); + * } + * + * if (L_result > EM) { + * Mc = m; + * EM = L_result; + * } + * } + */ + +#undef STEP +#define STEP( m, i ) L_temp = SASR( x[m + 3 * i], 2 ); \ + L_result += L_temp * L_temp; + + /* common part of 0 and 3 */ + + L_result = 0; + STEP( 0, 1 ); STEP( 0, 2 ); STEP( 0, 3 ); STEP( 0, 4 ); + STEP( 0, 5 ); STEP( 0, 6 ); STEP( 0, 7 ); STEP( 0, 8 ); + STEP( 0, 9 ); STEP( 0, 10); STEP( 0, 11); STEP( 0, 12); + L_common_0_3 = L_result; + + /* i = 0 */ + + STEP( 0, 0 ); + L_result <<= 1; /* implicit in L_MULT */ + EM = L_result; + + /* i = 1 */ + + L_result = 0; + STEP( 1, 0 ); + STEP( 1, 1 ); STEP( 1, 2 ); STEP( 1, 3 ); STEP( 1, 4 ); + STEP( 1, 5 ); STEP( 1, 6 ); STEP( 1, 7 ); STEP( 1, 8 ); + STEP( 1, 9 ); STEP( 1, 10); STEP( 1, 11); STEP( 1, 12); + L_result <<= 1; + if (L_result > EM) { + Mc = 1; + EM = L_result; + } + + /* i = 2 */ + + L_result = 0; + STEP( 2, 0 ); + STEP( 2, 1 ); STEP( 2, 2 ); STEP( 2, 3 ); STEP( 2, 4 ); + STEP( 2, 5 ); STEP( 2, 6 ); STEP( 2, 7 ); STEP( 2, 8 ); + STEP( 2, 9 ); STEP( 2, 10); STEP( 2, 11); STEP( 2, 12); + L_result <<= 1; + if (L_result > EM) { + Mc = 2; + EM = L_result; + } + + /* i = 3 */ + + L_result = L_common_0_3; + STEP( 3, 12 ); + L_result <<= 1; + if (L_result > EM) { + Mc = 3; + EM = L_result; + } + + /**/ + + /* Down-sampling by a factor 3 to get the selected xM[0..12] + * RPE sequence. + */ + for (i = 0; i <= 12; i ++) xM[i] = x[Mc + 3*i]; + *Mc_out = Mc; +} + +/* 4.12.15 */ + +static void APCM_quantization_xmaxc_to_exp_mant P3((xmaxc,exp_out,mant_out), + word xmaxc, /* IN */ + word * exp_out, /* OUT */ + word * mant_out ) /* OUT */ +{ + word exp, mant; + + /* Compute exponent and mantissa of the decoded version of xmaxc + */ + + exp = 0; + if (xmaxc > 15) exp = SASR(xmaxc, 3) - 1; + mant = xmaxc - (exp << 3); + + if (mant == 0) { + exp = -4; + mant = 7; + } + else { + while (mant <= 7) { + mant = mant << 1 | 1; + exp--; + } + mant -= 8; + } + + assert( exp >= -4 && exp <= 6 ); + assert( mant >= 0 && mant <= 7 ); + + *exp_out = exp; + *mant_out = mant; +} + +static void APCM_quantization P5((xM,xMc,mant_out,exp_out,xmaxc_out), + word * xM, /* [0..12] IN */ + + word * xMc, /* [0..12] OUT */ + word * mant_out, /* OUT */ + word * exp_out, /* OUT */ + word * xmaxc_out /* OUT */ +) +{ + int i, itest; + + word xmax, xmaxc, temp, temp1, temp2; + word exp, mant; + + + /* Find the maximum absolute value xmax of xM[0..12]. + */ + + xmax = 0; + for (i = 0; i <= 12; i++) { + temp = xM[i]; + temp = GSM_ABS(temp); + if (temp > xmax) xmax = temp; + } + + /* Qantizing and coding of xmax to get xmaxc. + */ + + exp = 0; + temp = SASR( xmax, 9 ); + itest = 0; + + for (i = 0; i <= 5; i++) { + + itest |= (temp <= 0); + temp = SASR( temp, 1 ); + + assert(exp <= 5); + if (itest == 0) exp++; /* exp = add (exp, 1) */ + } + + assert(exp <= 6 && exp >= 0); + temp = exp + 5; + + assert(temp <= 11 && temp >= 0); + xmaxc = gsm_add( SASR(xmax, temp), exp << 3 ); + + /* Quantizing and coding of the xM[0..12] RPE sequence + * to get the xMc[0..12] + */ + + APCM_quantization_xmaxc_to_exp_mant( xmaxc, &exp, &mant ); + + /* This computation uses the fact that the decoded version of xmaxc + * can be calculated by using the exponent and the mantissa part of + * xmaxc (logarithmic table). + * So, this method avoids any division and uses only a scaling + * of the RPE samples by a function of the exponent. A direct + * multiplication by the inverse of the mantissa (NRFAC[0..7] + * found in table 4.5) gives the 3 bit coded version xMc[0..12] + * of the RPE samples. + */ + + + /* Direct computation of xMc[0..12] using table 4.5 + */ + + assert( exp <= 4096 && exp >= -4096); + assert( mant >= 0 && mant <= 7 ); + + temp1 = 6 - exp; /* normalization by the exponent */ + temp2 = gsm_NRFAC[ mant ]; /* inverse mantissa */ + + for (i = 0; i <= 12; i++) { + + assert(temp1 >= 0 && temp1 < 16); + + temp = xM[i] << temp1; + temp = GSM_MULT( temp, temp2 ); + temp = SASR(temp, 12); + xMc[i] = temp + 4; /* see note below */ + } + + /* NOTE: This equation is used to make all the xMc[i] positive. + */ + + *mant_out = mant; + *exp_out = exp; + *xmaxc_out = xmaxc; +} + +/* 4.2.16 */ + +static void APCM_inverse_quantization P4((xMc,mant,exp,xMp), + register word * xMc, /* [0..12] IN */ + word mant, + word exp, + register word * xMp) /* [0..12] OUT */ +/* + * This part is for decoding the RPE sequence of coded xMc[0..12] + * samples to obtain the xMp[0..12] array. Table 4.6 is used to get + * the mantissa of xmaxc (FAC[0..7]). + */ +{ + int i; + word temp, temp1, temp2, temp3; + longword ltmp; + + assert( mant >= 0 && mant <= 7 ); + + temp1 = gsm_FAC[ mant ]; /* see 4.2-15 for mant */ + temp2 = gsm_sub( 6, exp ); /* see 4.2-15 for exp */ + temp3 = gsm_asl( 1, gsm_sub( temp2, 1 )); + + for (i = 13; i--;) { + + assert( *xMc <= 7 && *xMc >= 0 ); /* 3 bit unsigned */ + + /* temp = gsm_sub( *xMc++ << 1, 7 ); */ + temp = (*xMc++ << 1) - 7; /* restore sign */ + assert( temp <= 7 && temp >= -7 ); /* 4 bit signed */ + + temp <<= 12; /* 16 bit signed */ + temp = GSM_MULT_R( temp1, temp ); + temp = GSM_ADD( temp, temp3 ); + *xMp++ = gsm_asr( temp, temp2 ); + } +} + +/* 4.2.17 */ + +static void RPE_grid_positioning P3((Mc,xMp,ep), + word Mc, /* grid position IN */ + register word * xMp, /* [0..12] IN */ + register word * ep /* [0..39] OUT */ +) +/* + * This procedure computes the reconstructed long term residual signal + * ep[0..39] for the LTP analysis filter. The inputs are the Mc + * which is the grid position selection and the xMp[0..12] decoded + * RPE samples which are upsampled by a factor of 3 by inserting zero + * values. + */ +{ + int i = 13; + + assert(0 <= Mc && Mc <= 3); + + switch (Mc) { + case 3: *ep++ = 0; + case 2: do { + *ep++ = 0; + case 1: *ep++ = 0; + case 0: *ep++ = *xMp++; + } while (--i); + } + while (++Mc < 4) *ep++ = 0; + + /* + + int i, k; + for (k = 0; k <= 39; k++) ep[k] = 0; + for (i = 0; i <= 12; i++) { + ep[ Mc + (3*i) ] = xMp[i]; + } + */ +} + +/* 4.2.18 */ + +/* This procedure adds the reconstructed long term residual signal + * ep[0..39] to the estimated signal dpp[0..39] from the long term + * analysis filter to compute the reconstructed short term residual + * signal dp[-40..-1]; also the reconstructed short term residual + * array dp[-120..-41] is updated. + */ + +#if 0 /* Has been inlined in code.c */ +void Gsm_Update_of_reconstructed_short_time_residual_signal P3((dpp, ep, dp), + word * dpp, /* [0...39] IN */ + word * ep, /* [0...39] IN */ + word * dp) /* [-120...-1] IN/OUT */ +{ + int k; + + for (k = 0; k <= 79; k++) + dp[ -120 + k ] = dp[ -80 + k ]; + + for (k = 0; k <= 39; k++) + dp[ -40 + k ] = gsm_add( ep[k], dpp[k] ); +} +#endif /* Has been inlined in code.c */ + +void Gsm_RPE_Encoding P5((S,e,xmaxc,Mc,xMc), + + struct gsm_state * S, + + word * e, /* -5..-1][0..39][40..44 IN/OUT */ + word * xmaxc, /* OUT */ + word * Mc, /* OUT */ + word * xMc) /* [0..12] OUT */ +{ + word x[40]; + word xM[13], xMp[13]; + word mant, exp; + + Weighting_filter(e, x); + RPE_grid_selection(x, xM, Mc); + + APCM_quantization( xM, xMc, &mant, &exp, xmaxc); + APCM_inverse_quantization( xMc, mant, exp, xMp); + + RPE_grid_positioning( *Mc, xMp, e ); + +} + +void Gsm_RPE_Decoding P5((S, xmaxcr, Mcr, xMcr, erp), + struct gsm_state * S, + + word xmaxcr, + word Mcr, + word * xMcr, /* [0..12], 3 bits IN */ + word * erp /* [0..39] OUT */ +) +{ + word exp, mant; + word xMp[ 13 ]; + + APCM_quantization_xmaxc_to_exp_mant( xmaxcr, &exp, &mant ); + APCM_inverse_quantization( xMcr, mant, exp, xMp ); + RPE_grid_positioning( Mcr, xMp, erp ); + +} |