2
* Copyright 1992 by Jutta Degener and Carsten Bormann, Technische
3
* Universitaet Berlin. See the accompanying file "COPYRIGHT" for
4
* details. THERE IS ABSOLUTELY NO WARRANTY FOR THIS SOFTWARE.
7
/* $Header: /tmp_amd/presto/export/kbs/jutta/src/gsm/RCS/lpc.c,v 1.5 1994/12/30 23:14:54 jutta Exp $ */
21
* 4.2.4 .. 4.2.7 LPC ANALYSIS SECTION
27
static void Autocorrelation P2((s, L_ACF),
28
word * s, /* [0..159] IN/OUT */
29
longword * L_ACF) /* [0..8] OUT */
31
* The goal is to compute the array L_ACF[k]. The signal s[i] must
32
* be scaled in order to avoid an overflow situation.
37
word temp, smax, scalauto;
43
/* Dynamic scaling of the array s[0..159]
46
/* Search for the maximum.
49
for (k = 0; k <= 159; k++) {
50
temp = GSM_ABS( s[k] );
51
if (temp > smax) smax = temp;
54
/* Computation of the scaling factor.
56
if (smax == 0) scalauto = 0;
59
scalauto = 4 - gsm_norm( (longword)smax << 16 );/* sub(4,..) */
62
/* Scaling of the array s[0...159]
69
case n: for (k = 0; k <= 159; k++) \
70
float_s[k] = (float) \
71
(s[k] = GSM_MULT_R(s[k], 16384 >> (n-1)));\
75
case n: for (k = 0; k <= 159; k++) \
76
s[k] = GSM_MULT_R( s[k], 16384 >> (n-1) );\
78
# endif /* USE_FLOAT_MUL */
89
else for (k = 0; k <= 159; k++) float_s[k] = (float) s[k];
92
/* Compute the L_ACF[..].
96
register float * sp = float_s;
97
register float sl = *sp;
99
# define STEP(k) L_ACF[k] += (longword)(sl * sp[ -(k) ]);
104
# define STEP(k) L_ACF[k] += ((longword)sl * sp[ -(k) ]);
107
# define NEXTI sl = *++sp
110
for (k = 9; k--; L_ACF[k] = 0) ;
116
STEP(0); STEP(1); STEP(2);
118
STEP(0); STEP(1); STEP(2); STEP(3);
120
STEP(0); STEP(1); STEP(2); STEP(3); STEP(4);
122
STEP(0); STEP(1); STEP(2); STEP(3); STEP(4); STEP(5);
124
STEP(0); STEP(1); STEP(2); STEP(3); STEP(4); STEP(5); STEP(6);
126
STEP(0); STEP(1); STEP(2); STEP(3); STEP(4); STEP(5); STEP(6); STEP(7);
128
for (i = 8; i <= 159; i++) {
133
STEP(1); STEP(2); STEP(3); STEP(4);
134
STEP(5); STEP(6); STEP(7); STEP(8);
137
for (k = 9; k--; L_ACF[k] <<= 1) ;
140
/* Rescaling of the array s[0..159]
143
assert(scalauto <= 4);
144
for (k = 160; k--; *s++ <<= scalauto) ;
148
#if defined(USE_FLOAT_MUL) && defined(FAST)
150
static void Fast_Autocorrelation P2((s, L_ACF),
151
word * s, /* [0..159] IN/OUT */
152
longword * L_ACF) /* [0..8] OUT */
159
register float *sf = s_f;
161
for (i = 0; i < 160; ++i) sf[i] = s[i];
162
for (k = 0; k <= 8; k++) {
163
register float L_temp2 = 0;
164
register float *sfl = sf - k;
165
for (i = k; i < 160; ++i) L_temp2 += sf[i] * sfl[i];
166
f_L_ACF[k] = L_temp2;
168
scale = MAX_LONGWORD / f_L_ACF[0];
170
for (k = 0; k <= 8; k++) {
171
L_ACF[k] = f_L_ACF[k] * scale;
174
#endif /* defined (USE_FLOAT_MUL) && defined (FAST) */
178
static void Reflection_coefficients P2( (L_ACF, r),
179
longword * L_ACF, /* 0...8 IN */
180
register word * r /* 0...7 OUT */
183
register int i, m, n;
185
register longword ltmp;
186
word ACF[9]; /* 0..8 */
187
word P[ 9]; /* 0..8 */
188
word K[ 9]; /* 2..8 */
190
/* Schur recursion with 16 bits arithmetic.
194
for (i = 8; i--; *r++ = 0) ;
198
assert( L_ACF[0] != 0 );
199
temp = gsm_norm( L_ACF[0] );
201
assert(temp >= 0 && temp < 32);
204
for (i = 0; i <= 8; i++) ACF[i] = SASR( L_ACF[i] << temp, 16 );
206
/* Initialize array P[..] and K[..] for the recursion.
209
for (i = 1; i <= 7; i++) K[ i ] = ACF[ i ];
210
for (i = 0; i <= 8; i++) P[ i ] = ACF[ i ];
212
/* Compute reflection coefficients
214
for (n = 1; n <= 8; n++, r++) {
217
temp = GSM_ABS(temp);
219
for (i = n; i <= 8; i++) *r++ = 0;
223
*r = gsm_div( temp, P[0] );
226
if (P[1] > 0) *r = -*r; /* r[n] = sub(0, r[n]) */
227
assert (*r != MIN_WORD);
232
temp = GSM_MULT_R( P[1], *r );
233
P[0] = GSM_ADD( P[0], temp );
235
for (m = 1; m <= 8 - n; m++) {
236
temp = GSM_MULT_R( K[ m ], *r );
237
P[m] = GSM_ADD( P[ m+1 ], temp );
239
temp = GSM_MULT_R( P[ m+1 ], *r );
240
K[m] = GSM_ADD( K[ m ], temp );
247
static void Transformation_to_Log_Area_Ratios P1((r),
248
register word * r /* 0..7 IN/OUT */
251
* The following scaling for r[..] and LAR[..] has been used:
253
* r[..] = integer( real_r[..]*32768. ); -1 <= real_r < 1.
254
* LAR[..] = integer( real_LAR[..] * 16384 );
255
* with -1.625 <= real_LAR <= 1.625
262
/* Computation of the LAR[0..7] from the r[0..7]
264
for (i = 1; i <= 8; i++, r++) {
267
temp = GSM_ABS(temp);
272
} else if (temp < 31130) {
273
assert( temp >= 11059 );
276
assert( temp >= 26112 );
281
*r = *r < 0 ? -temp : temp;
282
assert( *r != MIN_WORD );
288
static void Quantization_and_coding P1((LAR),
289
register word * LAR /* [0..7] IN/OUT */
296
/* This procedure needs four tables; the following equations
297
* give the optimum scaling for the constants:
299
* A[0..7] = integer( real_A[0..7] * 1024 )
300
* B[0..7] = integer( real_B[0..7] * 512 )
301
* MAC[0..7] = maximum of the LARc[0..7]
302
* MIC[0..7] = minimum of the LARc[0..7]
306
# define STEP( A, B, MAC, MIC ) \
307
temp = GSM_MULT( A, *LAR ); \
308
temp = GSM_ADD( temp, B ); \
309
temp = GSM_ADD( temp, 256 ); \
310
temp = SASR( temp, 9 ); \
311
*LAR = temp>MAC ? MAC - MIC : (temp<MIC ? 0 : temp - MIC); \
314
STEP( 20480, 0, 31, -32 );
315
STEP( 20480, 0, 31, -32 );
316
STEP( 20480, 2048, 15, -16 );
317
STEP( 20480, -2560, 15, -16 );
319
STEP( 13964, 94, 7, -8 );
320
STEP( 15360, -1792, 7, -8 );
321
STEP( 8534, -341, 3, -4 );
322
STEP( 9036, -1144, 3, -4 );
327
void Gsm_LPC_Analysis P3((S, s,LARc),
329
word * s, /* 0..159 signals IN/OUT */
330
word * LARc) /* 0..7 LARc's OUT */
334
#if defined(USE_FLOAT_MUL) && defined(FAST)
335
if (S->fast) Fast_Autocorrelation (s, L_ACF );
338
Autocorrelation (s, L_ACF );
339
Reflection_coefficients (L_ACF, LARc );
340
Transformation_to_Log_Area_Ratios (LARc);
341
Quantization_and_coding (LARc);