1 /*
2  *  Copyright (c) 2011 The WebRTC project authors. All Rights Reserved.
3  *
4  *  Use of this source code is governed by a BSD-style license
5  *  that can be found in the LICENSE file in the root of the source
6  *  tree. An additional intellectual property rights grant can be found
7  *  in the file PATENTS.  All contributing project authors may
8  *  be found in the AUTHORS file in the root of the source tree.
9  */
10 
11 /******************************************************************
12 
13  iLBC Speech Coder ANSI-C Source Code
14 
15  WebRtcIlbcfix_Chebyshev.c
16 
17 ******************************************************************/
18 
19 #include "modules/audio_coding/codecs/ilbc/chebyshev.h"
20 
21 #include "modules/audio_coding/codecs/ilbc/constants.h"
22 #include "modules/audio_coding/codecs/ilbc/defines.h"
23 
24 /*------------------------------------------------------------------*
25  *  Calculate the Chevyshev polynomial series
26  *  F(w) = 2*exp(-j5w)*C(x)
27  *   C(x) = (T_0(x) + f(1)T_1(x) + ... + f(4)T_1(x) + f(5)/2)
28  *   T_i(x) is the i:th order Chebyshev polynomial
29  *------------------------------------------------------------------*/
30 
WebRtcIlbcfix_Chebyshev(int16_t x,int16_t * f)31 int16_t WebRtcIlbcfix_Chebyshev(
32     /* (o) Result of C(x) */
33     int16_t x,  /* (i) Value to the Chevyshev polynomial */
34     int16_t *f  /* (i) The coefficients in the polynomial */
35                                       ) {
36   int16_t b1_high, b1_low; /* Use the high, low format to increase the accuracy */
37   int32_t b2;
38   int32_t tmp1W32;
39   int32_t tmp2W32;
40   int i;
41 
42   b2 = (int32_t)0x1000000; /* b2 = 1.0 (Q23) */
43   /* Calculate b1 = 2*x + f[1] */
44   tmp1W32 = (x << 10) + (f[1] << 14);
45 
46   for (i = 2; i < 5; i++) {
47     tmp2W32 = tmp1W32;
48 
49     /* Split b1 (in tmp1W32) into a high and low part */
50     b1_high = (int16_t)(tmp1W32 >> 16);
51     b1_low = (int16_t)((tmp1W32 - ((int32_t)b1_high << 16)) >> 1);
52 
53     /* Calculate 2*x*b1-b2+f[i] */
54     tmp1W32 = ((b1_high * x + ((b1_low * x) >> 15)) << 2) - b2 + (f[i] << 14);
55 
56     /* Update b2 for next round */
57     b2 = tmp2W32;
58   }
59 
60   /* Split b1 (in tmp1W32) into a high and low part */
61   b1_high = (int16_t)(tmp1W32 >> 16);
62   b1_low = (int16_t)((tmp1W32 - ((int32_t)b1_high << 16)) >> 1);
63 
64   /* tmp1W32 = x*b1 - b2 + f[i]/2 */
65   tmp1W32 = ((b1_high * x) << 1) + (((b1_low * x) >> 15) << 1) -
66       b2 + (f[i] << 13);
67 
68   /* Handle overflows and set to maximum or minimum int16_t instead */
69   if (tmp1W32>((int32_t)33553408)) {
70     return(WEBRTC_SPL_WORD16_MAX);
71   } else if (tmp1W32<((int32_t)-33554432)) {
72     return(WEBRTC_SPL_WORD16_MIN);
73   } else {
74     return (int16_t)(tmp1W32 >> 10);
75   }
76 }
77