1 /* ------------------------------------------------------------------
2 * Copyright (C) 1998-2009 PacketVideo
4 * Licensed under the Apache License, Version 2.0 (the "License");
5 * you may not use this file except in compliance with the License.
6 * You may obtain a copy of the License at
8 * http://www.apache.org/licenses/LICENSE-2.0
10 * Unless required by applicable law or agreed to in writing, software
11 * distributed under the License is distributed on an "AS IS" BASIS,
12 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either
14 * See the License for the specific language governing permissions
15 * and limitations under the License.
16 * -------------------------------------------------------------------
18 /****************************************************************************************
19 Portions of this file are derived from the following 3GPP standard:
22 ANSI-C code for the Adaptive Multi-Rate (AMR) speech codec
23 Available from http://www.3gpp.org
25 (C) 2004, 3GPP Organizational Partners (ARIB, ATIS, CCSA, ETSI, TTA, TTC)
26 Permission to distribute, modify and use this file under the standard license
27 terms listed above has been obtained from the copyright holder.
28 ****************************************************************************************/
36 ------------------------------------------------------------------------------
39 These modules compute the LSPs from the LP coefficients.
41 ------------------------------------------------------------------------------
45 /*----------------------------------------------------------------------------
47 ----------------------------------------------------------------------------*/
52 /*----------------------------------------------------------------------------
54 ; Define module specific macros here
55 ----------------------------------------------------------------------------*/
58 /*----------------------------------------------------------------------------
60 ; Include all pre-processor statements here. Include conditional
61 ; compile variables also.
62 ----------------------------------------------------------------------------*/
63 #define NC M/2 /* M = LPC order, NC = M/2 */
65 /*----------------------------------------------------------------------------
66 ; LOCAL FUNCTION DEFINITIONS
67 ; Function Prototype declaration
68 ----------------------------------------------------------------------------*/
71 /*----------------------------------------------------------------------------
72 ; LOCAL VARIABLE DEFINITIONS
73 ; Variable declaration - defined here and used outside this module
74 ----------------------------------------------------------------------------*/
77 ------------------------------------------------------------------------------
79 ------------------------------------------------------------------------------
80 INPUT AND OUTPUT DEFINITIONS
83 x = input value (Word16)
84 f = polynomial (Word16)
85 n = polynomial order (Word16)
87 pOverflow = pointer to overflow (Flag)
90 pOverflow -> 1 if the operations in the function resulted in saturation.
93 cheb = Chebyshev polynomial for the input value x.(Word16)
95 Global Variables Used:
98 Local Variables Needed:
101 ------------------------------------------------------------------------------
104 This module evaluates the Chebyshev polynomial series.
105 - The polynomial order is n = m/2 = 5
106 - The polynomial F(z) (F1(z) or F2(z)) is given by
107 F(w) = 2 exp(-j5w) C(x)
109 C(x) = T_n(x) + f(1)T_n-1(x) + ... +f(n-1)T_1(x) + f(n)/2
110 and T_m(x) = cos(mw) is the mth order Chebyshev
111 polynomial ( x=cos(w) )
112 - C(x) for the input x is returned.
114 ------------------------------------------------------------------------------
119 ------------------------------------------------------------------------------
122 az_lsp.c, UMTS GSM AMR speech codec, R99 - Version 3.2.0, March 2, 2001
124 ------------------------------------------------------------------------------
127 static Word16 Chebps (Word16 x,
132 Word16 b0_h, b0_l, b1_h, b1_l, b2_h, b2_l;
135 // The reference ETSI code uses a global flag for Overflow. However, in the
136 // actual implementation a pointer to Overflow flag is passed in as a
137 // parameter to the function. This pointer is passed into all the basic math
140 b2_h = 256; // b2 = 1.0
143 t0 = L_mult (x, 512); // 2*x
144 t0 = L_mac (t0, f[1], 8192); // + f[1]
145 L_Extract (t0, &b1_h, &b1_l); // b1 = 2*x + f[1]
147 for (i = 2; i < n; i++)
149 t0 = Mpy_32_16 (b1_h, b1_l, x); // t0 = 2.0*x*b1
151 t0 = L_mac (t0, b2_h, (Word16) 0x8000); // t0 = 2.0*x*b1 - b2
152 t0 = L_msu (t0, b2_l, 1);
153 t0 = L_mac (t0, f[i], 8192); // t0 = 2.0*x*b1 - b2 + f[i]
155 L_Extract (t0, &b0_h, &b0_l); // b0 = 2.0*x*b1 - b2 + f[i]
157 b2_l = b1_l; // b2 = b1;
159 b1_l = b0_l; // b1 = b0;
163 t0 = Mpy_32_16 (b1_h, b1_l, x); // t0 = x*b1;
164 t0 = L_mac (t0, b2_h, (Word16) 0x8000); // t0 = x*b1 - b2
165 t0 = L_msu (t0, b2_l, 1);
166 t0 = L_mac (t0, f[i], 4096); // t0 = x*b1 - b2 + f[i]/2
170 cheb = extract_h (t0);
175 ------------------------------------------------------------------------------
177 [State any special notes, constraints or cautions for users of this function]
179 ------------------------------------------------------------------------------
182 static Word16 Chebps(Word16 x,
183 Word16 f[], /* (n) */
195 OSCL_UNUSED_ARG(pOverflow);
199 L_temp = 0x01000000L;
201 t0 = ((Word32) x << 10) + ((Word32) * (p_f++) << 14);
203 /* b1 = t0 = 2*x + f[1] */
205 b1_h = (Word16)(t0 >> 16);
206 b1_l = (Word16)((t0 >> 1) - (b1_h << 15));
209 for (i = 2; i < n; i++)
212 t0 = ((Word32) b1_h * x);
213 t0 += ((Word32) b1_l * x) >> 15;
216 /* t0 = 2.0*x*b1 - b2 */
219 /* t0 = 2.0*x*b1 - b2 + f[i] */
220 t0 += (Word32) * (p_f++) << 14;
222 L_temp = ((Word32) b1_h << 16) + ((Word32) b1_l << 1);
224 /* b0 = 2.0*x*b1 - b2 + f[i]*/
225 b1_h = (Word16)(t0 >> 16);
226 b1_l = (Word16)((t0 >> 1) - (b1_h << 15));
231 t0 = ((Word32) b1_h * x);
232 t0 += ((Word32) b1_l * x) >> 15;
239 /* t0 = x*b1 - b2 + f[i]/2 */
240 t0 += (Word32) * (p_f) << 13;
243 if ((UWord32)(t0 + 33554432) < 67108863)
245 cheb = (Word16)(t0 >> 10);
249 if (t0 > (Word32) 0x01ffffffL)
265 ------------------------------------------------------------------------------
266 FUNCTION NAME: Az_lsp
267 ------------------------------------------------------------------------------
268 INPUT AND OUTPUT DEFINITIONS FOR Az_lsp
271 a = predictor coefficients (Word16)
272 lsp = line spectral pairs (Word16)
273 old_lsp = old line spectral pairs (Word16)
275 pOverflow = pointer to overflow (Flag)
278 pOverflow -> 1 if the operations in the function resulted in saturation.
283 Global Variables Used:
286 Local Variables Needed:
289 ------------------------------------------------------------------------------
292 This function computes the LSPs from the LP coefficients.
294 The sum and difference filters are computed and divided by 1+z^{-1} and
295 1-z^{-1}, respectively.
297 f1[i] = a[i] + a[11-i] - f1[i-1] ; i=1,...,5
298 f2[i] = a[i] - a[11-i] + f2[i-1] ; i=1,...,5
300 The roots of F1(z) and F2(z) are found using Chebyshev polynomial evaluation.
301 The polynomials are evaluated at 60 points regularly spaced in the
302 frequency domain. The sign change interval is subdivided 4 times to better
303 track the root. The LSPs are found in the cosine domain [1,-1].
305 If less than 10 roots are found, the LSPs from the past frame are used.
307 ------------------------------------------------------------------------------
312 ------------------------------------------------------------------------------
315 az_lsp.c, UMTS GSM AMR speech codec, R99 - Version 3.2.0, March 2, 2001
317 ------------------------------------------------------------------------------
321 Word16 a[], // (i) : predictor coefficients (MP1)
322 Word16 lsp[], // (o) : line spectral pairs (M)
323 Word16 old_lsp[] // (i) : old lsp[] (in case not found 10 roots) (M)
327 Word16 xlow, ylow, xhigh, yhigh, xmid, ymid, xint;
328 Word16 x, y, sign, exp;
330 Word16 f1[M / 2 + 1], f2[M / 2 + 1];
333 *-------------------------------------------------------------*
334 * find the sum and diff. pol. F1(z) and F2(z) *
335 * F1(z) <--- F1(z)/(1+z**-1) & F2(z) <--- F2(z)/(1-z**-1) *
340 * for (i = 0; i< NC; i++) *
342 * f1[i+1] = a[i+1] + a[M-i] - f1[i] ; *
343 * f2[i+1] = a[i+1] - a[M-i] + f2[i] ; *
345 *-------------------------------------------------------------*
347 f1[0] = 1024; // f1[0] = 1.0
348 f2[0] = 1024; // f2[0] = 1.0
350 // The reference ETSI code uses a global flag for Overflow. However, in the
351 // actual implementation a pointer to Overflow flag is passed in as a
352 // parameter to the function. This pointer is passed into all the basic math
355 for (i = 0; i < NC; i++)
357 t0 = L_mult (a[i + 1], 8192); // x = (a[i+1] + a[M-i]) >> 2
358 t0 = L_mac (t0, a[M - i], 8192);
360 // f1[i+1] = a[i+1] + a[M-i] - f1[i]
361 f1[i + 1] = sub (x, f1[i]);
363 t0 = L_mult (a[i + 1], 8192); // x = (a[i+1] - a[M-i]) >> 2
364 t0 = L_msu (t0, a[M - i], 8192);
366 // f2[i+1] = a[i+1] - a[M-i] + f2[i]
367 f2[i + 1] = add (x, f2[i]);
370 *-------------------------------------------------------------*
371 * find the LSPs using the Chebychev pol. evaluation *
372 *-------------------------------------------------------------*
374 nf = 0; // number of found frequencies
375 ip = 0; // indicator for f1 or f2
380 ylow = Chebps (xlow, coef, NC);
383 // while ( (nf < M) && (j < grid_points) )
384 while ((sub (nf, M) < 0) && (sub (j, grid_points) < 0))
390 ylow = Chebps (xlow, coef, NC);
392 if (L_mult (ylow, yhigh) <= (Word32) 0L)
395 // divide 4 times the interval
397 for (i = 0; i < 4; i++)
399 // xmid = (xlow + xhigh)/2
400 xmid = add (shr (xlow, 1), shr (xhigh, 1));
401 ymid = Chebps (xmid, coef, NC);
403 if (L_mult (ylow, ymid) <= (Word32) 0L)
415 *-------------------------------------------------------------*
416 * Linear interpolation *
417 * xint = xlow - ylow*(xhigh-xlow)/(yhigh-ylow); *
418 *-------------------------------------------------------------*
420 x = sub (xhigh, xlow);
421 y = sub (yhigh, ylow);
433 y = div_s ((Word16) 16383, y);
435 t0 = L_shr (t0, sub (20, exp));
436 y = extract_l (t0); // y= (xhigh-xlow)/(yhigh-ylow)
441 t0 = L_mult (ylow, y);
443 xint = sub (xlow, extract_l (t0)); // xint = xlow - ylow*y
460 ylow = Chebps (xlow, coef, NC);
465 // Check if M roots found
469 for (i = 0; i < M; i++)
478 ------------------------------------------------------------------------------
480 [State any special notes, constraints or cautions for users of this function]
482 ------------------------------------------------------------------------------
486 Word16 a[], /* (i) : predictor coefficients (MP1) */
487 Word16 lsp[], /* (o) : line spectral pairs (M) */
488 Word16 old_lsp[], /* (i) : old lsp[] (in case not found 10 roots) (M) */
489 Flag *pOverflow /* (i/o): overflow flag */
515 /*-------------------------------------------------------------*
516 * find the sum and diff. pol. F1(z) and F2(z) *
517 * F1(z) <--- F1(z)/(1+z**-1) & F2(z) <--- F2(z)/(1-z**-1) *
522 * for (i = 0; i< NC; i++) *
524 * f1[i+1] = a[i+1] + a[M-i] - f1[i] ; *
525 * f2[i+1] = a[i+1] - a[M-i] + f2[i] ; *
527 *-------------------------------------------------------------*/
529 *p_f1 = 1024; /* f1[0] = 1.0 */
530 *p_f2 = 1024; /* f2[0] = 1.0 */
532 for (i = 0; i < NC; i++)
534 L_temp1 = (Word32) * (a + i + 1);
535 L_temp2 = (Word32) * (a + M - i);
536 /* x = (a[i+1] + a[M-i]) >> 2 */
537 x = (Word16)((L_temp1 + L_temp2) >> 2);
538 /* y = (a[i+1] - a[M-i]) >> 2 */
539 y = (Word16)((L_temp1 - L_temp2) >> 2);
540 /* f1[i+1] = a[i+1] + a[M-i] - f1[i] */
543 /* f2[i+1] = a[i+1] - a[M-i] + f2[i] */
548 /*-------------------------------------------------------------*
549 * find the LSPs using the Chebychev pol. evaluation *
550 *-------------------------------------------------------------*/
552 nf = 0; /* number of found frequencies */
553 ip = 0; /* indicator for f1 or f2 */
558 ylow = Chebps(xlow, coef, NC, pOverflow);
562 while ((nf < M) && (j < grid_points))
568 ylow = Chebps(xlow, coef, NC, pOverflow);
570 if (((Word32)ylow*yhigh) <= 0)
572 /* divide 4 times the interval */
573 for (i = 4; i != 0; i--)
575 /* xmid = (xlow + xhigh)/2 */
580 ymid = Chebps(xmid, coef, NC, pOverflow);
582 if (((Word32)ylow*ymid) <= 0)
594 /*-------------------------------------------------------------*
595 * Linear interpolation *
596 * xint = xlow - ylow*(xhigh-xlow)/(yhigh-ylow); *
597 *-------------------------------------------------------------*/
612 y = div_s((Word16) 16383, y);
614 y = ((Word32)x * y) >> (19 - exp);
621 /* xint = xlow - ylow*y */
622 xint = xlow - (((Word32) ylow * y) >> 10);
640 ylow = Chebps(xlow, coef, NC, pOverflow);
645 /* Check if M roots found */
649 for (i = NC; i != 0 ; i--)
658 Word16 Chebps_Wrapper(Word16 x,
659 Word16 f[], /* (n) */
663 return Chebps(x, f, n, pOverflow);