1 /*
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.
5 */
6
7 #include <stdio.h>
8 #include <assert.h>
9
10 #include "gsm610_priv.h"
11
12 /* 4.2.13 .. 4.2.17 RPE ENCODING SECTION
13 */
14
15 /* 4.2.13 */
16
Weighting_filter(register int16_t * e,int16_t * x)17 static void Weighting_filter (
18 register int16_t * e, /* signal [-5..0.39.44] IN */
19 int16_t * x /* signal [0..39] OUT */
20 )
21 /*
22 * The coefficients of the weighting filter are stored in a table
23 * (see table 4.4). The following scaling is used:
24 *
25 * H[0..10] = integer(real_H [0..10] * 8192) ;
26 */
27 {
28 /* int16_t wt [50] ; */
29
30 register int32_t L_result ;
31 register int k /* , i */ ;
32
33 /* Initialization of a temporary working array wt[0...49]
34 */
35
36 /* for (k = 0 ; k <= 4 ; k++) wt[k] = 0 ;
37 * for (k = 5 ; k <= 44 ; k++) wt[k] = *e++;
38 * for (k = 45 ; k <= 49 ; k++) wt[k] = 0 ;
39 *
40 * (e[-5..-1] and e[40..44] are allocated by the caller,
41 * are initially zero and are not written anywhere.)
42 */
43 e -= 5 ;
44
45 /* Compute the signal x[0..39]
46 */
47 for (k = 0 ; k <= 39 ; k++)
48 { L_result = 8192 >> 1 ;
49
50 /* for (i = 0 ; i <= 10 ; i++) {
51 * L_temp = GSM_L_MULT(wt[k+i], gsm_H[i]) ;
52 * L_result = GSM_L_ADD(L_result, L_temp) ;
53 * }
54 */
55
56 #undef STEP
57 #define STEP(i, H) (e [k + i] * (int32_t) H)
58
59 /* Every one of these multiplications is done twice --
60 * but I don't see an elegant way to optimize this.
61 * Do you?
62 */
63
64 #ifdef STUPID_COMPILER
65 L_result += STEP (0, -134) ;
66 L_result += STEP (1, -374) ;
67 /* + STEP (2, 0) */
68 L_result += STEP (3, 2054) ;
69 L_result += STEP (4, 5741) ;
70 L_result += STEP (5, 8192) ;
71 L_result += STEP (6, 5741) ;
72 L_result += STEP (7, 2054) ;
73 /* + STEP (8, 0) */
74 L_result += STEP (9, -374) ;
75 L_result += STEP (10, -134) ;
76 #else
77 L_result += STEP (0, -134)
78 + STEP (1, -374)
79 /* + STEP (2, 0) */
80 + STEP (3, 2054)
81 + STEP (4, 5741)
82 + STEP (5, 8192)
83 + STEP (6, 5741)
84 + STEP (7, 2054)
85 /* + STEP (8, 0) */
86 + STEP (9, -374)
87 + STEP (10, -134) ;
88 #endif
89
90 /* L_result = GSM_L_ADD(L_result, L_result) ; (* scaling(x2) *)
91 * L_result = GSM_L_ADD(L_result, L_result) ; (* scaling(x4) *)
92 *
93 * x[k] = SASR(L_result, 16) ;
94 */
95
96 /* 2 adds vs. >>16 => 14, minus one shift to compensate for
97 * those we lost when replacing L_MULT by '*'.
98 */
99
100 L_result = SASR_L (L_result, 13) ;
101 x [k] = (L_result < MIN_WORD ? MIN_WORD
102 : (L_result > MAX_WORD ? MAX_WORD : L_result)) ;
103 }
104 }
105
106 /* 4.2.14 */
107
RPE_grid_selection(int16_t * x,int16_t * xM,int16_t * Mc_out)108 static void RPE_grid_selection (
109 int16_t * x, /* [0..39] IN */
110 int16_t * xM, /* [0..12] OUT */
111 int16_t * Mc_out /* OUT */
112 )
113 /*
114 * The signal x[0..39] is used to select the RPE grid which is
115 * represented by Mc.
116 */
117 {
118 register int i ;
119 register int32_t L_result, L_temp ;
120 int32_t EM ; /* xxx should be L_EM? */
121 int16_t Mc ;
122
123 int32_t L_common_0_3 ;
124
125 EM = 0 ;
126 Mc = 0 ;
127
128 /* for (m = 0 ; m <= 3 ; m++) {
129 * L_result = 0 ;
130 *
131 *
132 * for (i = 0 ; i <= 12 ; i++) {
133 *
134 * temp1 = SASR_W (x[m + 3*i], 2) ;
135 *
136 * assert (temp1 != MIN_WORD) ;
137 *
138 * L_temp = GSM_L_MULT(temp1, temp1) ;
139 * L_result = GSM_L_ADD(L_temp, L_result) ;
140 * }
141 *
142 * if (L_result > EM) {
143 * Mc = m ;
144 * EM = L_result ;
145 * }
146 * }
147 */
148
149 #undef STEP
150 #define STEP(m, i) L_temp = SASR_W (x [m + 3 * i], 2) ; \
151 L_result += L_temp * L_temp ;
152
153 /* common part of 0 and 3 */
154
155 L_result = 0 ;
156 STEP (0, 1) ; STEP (0, 2) ; STEP (0, 3) ; STEP (0, 4) ;
157 STEP (0, 5) ; STEP (0, 6) ; STEP (0, 7) ; STEP (0, 8) ;
158 STEP (0, 9) ; STEP (0, 10) ; STEP (0, 11) ; STEP (0, 12) ;
159 L_common_0_3 = L_result ;
160
161 /* i = 0 */
162
163 STEP (0, 0) ;
164 L_result <<= 1 ; /* implicit in L_MULT */
165 EM = L_result ;
166
167 /* i = 1 */
168
169 L_result = 0 ;
170 STEP (1, 0) ;
171 STEP (1, 1) ; STEP (1, 2) ; STEP (1, 3) ; STEP (1, 4) ;
172 STEP (1, 5) ; STEP (1, 6) ; STEP (1, 7) ; STEP (1, 8) ;
173 STEP (1, 9) ; STEP (1, 10) ; STEP (1, 11) ; STEP (1, 12) ;
174 L_result <<= 1 ;
175 if (L_result > EM)
176 { Mc = 1 ;
177 EM = L_result ;
178 }
179
180 /* i = 2 */
181
182 L_result = 0 ;
183 STEP (2, 0) ;
184 STEP (2, 1) ; STEP (2, 2) ; STEP (2, 3) ; STEP (2, 4) ;
185 STEP (2, 5) ; STEP (2, 6) ; STEP (2, 7) ; STEP (2, 8) ;
186 STEP (2, 9) ; STEP (2, 10) ; STEP (2, 11) ; STEP (2, 12) ;
187 L_result <<= 1 ;
188 if (L_result > EM)
189 { Mc = 2 ;
190 EM = L_result ;
191 }
192
193 /* i = 3 */
194
195 L_result = L_common_0_3 ;
196 STEP (3, 12) ;
197 L_result <<= 1 ;
198 if (L_result > EM)
199 { Mc = 3 ;
200 EM = L_result ;
201 }
202
203 /* Down-sampling by a factor 3 to get the selected xM [0..12]
204 * RPE sequence.
205 */
206 for (i = 0 ; i <= 12 ; i ++) xM [i] = x [Mc + 3 * i] ;
207 *Mc_out = Mc ;
208 }
209
210 /* 4.12.15 */
211
APCM_quantization_xmaxc_to_exp_mant(int16_t xmaxc,int16_t * expon_out,int16_t * mant_out)212 static void APCM_quantization_xmaxc_to_exp_mant (
213 int16_t xmaxc, /* IN */
214 int16_t * expon_out, /* OUT */
215 int16_t * mant_out) /* OUT */
216 {
217 int16_t expon, mant ;
218
219 /* Compute expononent and mantissa of the decoded version of xmaxc
220 */
221
222 expon = 0 ;
223 if (xmaxc > 15) expon = SASR_W (xmaxc, 3) - 1 ;
224 mant = xmaxc - (expon << 3) ;
225
226 if (mant == 0)
227 { expon = -4 ;
228 mant = 7 ;
229 }
230 else
231 { while (mant <= 7)
232 { mant = mant << 1 | 1 ;
233 expon-- ;
234 }
235 mant -= 8 ;
236 }
237
238 assert (expon >= -4 && expon <= 6) ;
239 assert (mant >= 0 && mant <= 7) ;
240
241 *expon_out = expon ;
242 *mant_out = mant ;
243 }
244
APCM_quantization(int16_t * xM,int16_t * xMc,int16_t * mant_out,int16_t * expon_out,int16_t * xmaxc_out)245 static void APCM_quantization (
246 int16_t * xM, /* [0..12] IN */
247 int16_t * xMc, /* [0..12] OUT */
248 int16_t * mant_out, /* OUT */
249 int16_t * expon_out, /* OUT */
250 int16_t * xmaxc_out /* OUT */
251 )
252 {
253 int i, itest ;
254
255 int16_t xmax, xmaxc, temp, temp1, temp2 ;
256 int16_t expon, mant ;
257
258
259 /* Find the maximum absolute value xmax of xM [0..12].
260 */
261
262 xmax = 0 ;
263 for (i = 0 ; i <= 12 ; i++)
264 { temp = xM [i] ;
265 temp = GSM_ABS (temp) ;
266 if (temp > xmax) xmax = temp ;
267 }
268
269 /* Qantizing and coding of xmax to get xmaxc.
270 */
271
272 expon = 0 ;
273 temp = SASR_W (xmax, 9) ;
274 itest = 0 ;
275
276 for (i = 0 ; i <= 5 ; i++)
277 { itest |= (temp <= 0) ;
278 temp = SASR_W (temp, 1) ;
279
280 assert (expon <= 5) ;
281 if (itest == 0) expon++ ; /* expon = add (expon, 1) */
282 }
283
284 assert (expon <= 6 && expon >= 0) ;
285 temp = expon + 5 ;
286
287 assert (temp <= 11 && temp >= 0) ;
288 xmaxc = gsm_add (SASR_W (xmax, temp), (int16_t) (expon << 3)) ;
289
290 /* Quantizing and coding of the xM [0..12] RPE sequence
291 * to get the xMc [0..12]
292 */
293
294 APCM_quantization_xmaxc_to_exp_mant (xmaxc, &expon, &mant) ;
295
296 /* This computation uses the fact that the decoded version of xmaxc
297 * can be calculated by using the expononent and the mantissa part of
298 * xmaxc (logarithmic table).
299 * So, this method avoids any division and uses only a scaling
300 * of the RPE samples by a function of the expononent. A direct
301 * multiplication by the inverse of the mantissa (NRFAC[0..7]
302 * found in table 4.5) gives the 3 bit coded version xMc [0..12]
303 * of the RPE samples.
304 */
305
306
307 /* Direct computation of xMc [0..12] using table 4.5
308 */
309
310 assert (expon <= 4096 && expon >= -4096) ;
311 assert (mant >= 0 && mant <= 7) ;
312
313 temp1 = 6 - expon ; /* normalization by the expononent */
314 temp2 = gsm_NRFAC [mant] ; /* inverse mantissa */
315
316 for (i = 0 ; i <= 12 ; i++)
317 { assert (temp1 >= 0 && temp1 < 16) ;
318
319 temp = arith_shift_left (xM [i], temp1) ;
320 temp = GSM_MULT (temp, temp2) ;
321 temp = SASR_W (temp, 12) ;
322 xMc [i] = temp + 4 ; /* see note below */
323 }
324
325 /* NOTE: This equation is used to make all the xMc [i] positive.
326 */
327
328 *mant_out = mant ;
329 *expon_out = expon ;
330 *xmaxc_out = xmaxc ;
331 }
332
333 /* 4.2.16 */
334
APCM_inverse_quantization(register int16_t * xMc,int16_t mant,int16_t expon,register int16_t * xMp)335 static void APCM_inverse_quantization (
336 register int16_t * xMc, /* [0..12] IN */
337 int16_t mant,
338 int16_t expon,
339 register int16_t * xMp) /* [0..12] OUT */
340 /*
341 * This part is for decoding the RPE sequence of coded xMc [0..12]
342 * samples to obtain the xMp[0..12] array. Table 4.6 is used to get
343 * the mantissa of xmaxc (FAC[0..7]).
344 */
345 {
346 int i ;
347 int16_t temp, temp1, temp2, temp3 ;
348
349 assert (mant >= 0 && mant <= 7) ;
350
351 temp1 = gsm_FAC [mant] ; /* see 4.2-15 for mant */
352 temp2 = gsm_sub (6, expon) ; /* see 4.2-15 for exp */
353 temp3 = gsm_asl (1, gsm_sub (temp2, 1)) ;
354
355 for (i = 13 ; i-- ;)
356 { assert (*xMc <= 7 && *xMc >= 0) ; /* 3 bit unsigned */
357
358 /* temp = gsm_sub (*xMc++ << 1, 7) ; */
359 temp = (*xMc++ << 1) - 7 ; /* restore sign */
360 assert (temp <= 7 && temp >= -7) ; /* 4 bit signed */
361
362 temp = arith_shift_left (temp, 12) ; /* 16 bit signed */
363 temp = GSM_MULT_R (temp1, temp) ;
364 temp = GSM_ADD (temp, temp3) ;
365 *xMp++ = gsm_asr (temp, temp2) ;
366 }
367 }
368
369 /* 4.2.17 */
370
RPE_grid_positioning(int16_t Mc,register int16_t * xMp,register int16_t * ep)371 static void RPE_grid_positioning (
372 int16_t Mc, /* grid position IN */
373 register int16_t * xMp, /* [0..12] IN */
374 register int16_t * ep /* [0..39] OUT */
375 )
376 /*
377 * This procedure computes the reconstructed long term residual signal
378 * ep[0..39] for the LTP analysis filter. The inputs are the Mc
379 * which is the grid position selection and the xMp[0..12] decoded
380 * RPE samples which are upsampled by a factor of 3 by inserting zero
381 * values.
382 */
383 {
384 int i = 13 ;
385
386 assert (0 <= Mc && Mc <= 3) ;
387
388 switch (Mc)
389 { case 3: *ep++ = 0 ;
390 /* Falls through. */
391 case 2: do
392 { *ep++ = 0 ;
393 /* Falls through. */
394 case 1: *ep++ = 0 ;
395 /* Falls through. */
396 case 0: *ep++ = *xMp++ ;
397 } while (--i) ;
398 }
399 while (++Mc < 4) *ep++ = 0 ;
400 }
401
402 /* 4.2.18 */
403
404 /* This procedure adds the reconstructed long term residual signal
405 * ep[0..39] to the estimated signal dpp[0..39] from the long term
406 * analysis filter to compute the reconstructed short term residual
407 * signal dp[-40..-1] ; also the reconstructed short term residual
408 * array dp[-120..-41] is updated.
409 */
410
411 #if 0 /* Has been inlined in code.c */
412 void Gsm_Update_of_reconstructed_short_time_residual_signal (
413 int16_t * dpp, /* [0...39] IN */
414 int16_t * ep, /* [0...39] IN */
415 int16_t * dp) /* [-120...-1] IN/OUT */
416 {
417 int k ;
418
419 for (k = 0 ; k <= 79 ; k++)
420 dp [-120 + k] = dp [-80 + k] ;
421
422 for (k = 0 ; k <= 39 ; k++)
423 dp [-40 + k] = gsm_add (ep [k], dpp [k]) ;
424 }
425 #endif /* Has been inlined in code.c */
426
Gsm_RPE_Encoding(int16_t * e,int16_t * xmaxc,int16_t * Mc,int16_t * xMc)427 void Gsm_RPE_Encoding (
428 int16_t * e, /* -5..-1][0..39][40..44 IN/OUT */
429 int16_t * xmaxc, /* OUT */
430 int16_t * Mc, /* OUT */
431 int16_t * xMc) /* [0..12] OUT */
432 {
433 int16_t x [40] ;
434 int16_t xM [13], xMp [13] ;
435 int16_t mant, expon ;
436
437 Weighting_filter (e, x) ;
438 RPE_grid_selection (x, xM, Mc) ;
439
440 APCM_quantization (xM, xMc, &mant, &expon, xmaxc) ;
441 APCM_inverse_quantization (xMc, mant, expon, xMp) ;
442
443 RPE_grid_positioning (*Mc, xMp, e) ;
444
445 }
446
Gsm_RPE_Decoding(int16_t xmaxcr,int16_t Mcr,int16_t * xMcr,int16_t * erp)447 void Gsm_RPE_Decoding (
448 int16_t xmaxcr,
449 int16_t Mcr,
450 int16_t * xMcr, /* [0..12], 3 bits IN */
451 int16_t * erp /* [0..39] OUT */
452 )
453 {
454 int16_t expon, mant ;
455 int16_t xMp [13] ;
456
457 APCM_quantization_xmaxc_to_exp_mant (xmaxcr, &expon, &mant) ;
458 APCM_inverse_quantization (xMcr, mant, expon, xMp) ;
459 RPE_grid_positioning (Mcr, xMp, erp) ;
460 }
461