2
|
1
|
|
2 /******************************************************************
|
|
3
|
|
4 iLBC Speech Coder ANSI-C Source Code
|
|
5
|
|
6 lsf.c
|
|
7
|
|
8 Copyright (C) The Internet Society (2004).
|
|
9 All Rights Reserved.
|
|
10
|
|
11 ******************************************************************/
|
|
12
|
|
13 #include <string.h>
|
|
14
|
|
15
|
|
16
|
|
17
|
|
18
|
|
19 #include <math.h>
|
|
20
|
|
21 #include "iLBC_define.h"
|
|
22
|
|
23 /*----------------------------------------------------------------*
|
|
24 * conversion from lpc coefficients to lsf coefficients
|
|
25 *---------------------------------------------------------------*/
|
|
26
|
|
27 void a2lsf(float *freq, /* (o) lsf coefficients */
|
|
28 float *a /* (i) lpc coefficients */
|
|
29 )
|
|
30 {
|
|
31 float steps[LSF_NUMBER_OF_STEPS] =
|
|
32 { (float) 0.00635, (float) 0.003175, (float) 0.0015875,
|
|
33 (float) 0.00079375
|
|
34 };
|
|
35 float step;
|
|
36 int step_idx;
|
|
37 int lsp_index;
|
|
38 float p[LPC_HALFORDER];
|
|
39 float q[LPC_HALFORDER];
|
|
40 float p_pre[LPC_HALFORDER];
|
|
41 float q_pre[LPC_HALFORDER];
|
|
42 float old_p, old_q, *old;
|
|
43 float *pq_coef;
|
|
44 float omega, old_omega;
|
|
45 int i;
|
|
46 float hlp, hlp1, hlp2, hlp3, hlp4, hlp5;
|
|
47
|
|
48 for (i = 0; i < LPC_HALFORDER; i++) {
|
|
49 p[i] = (float) -1.0 * (a[i + 1] + a[LPC_FILTERORDER - i]);
|
|
50 q[i] = a[LPC_FILTERORDER - i] - a[i + 1];
|
|
51 }
|
|
52
|
|
53 p_pre[0] = (float) -1.0 - p[0];
|
|
54 p_pre[1] = -p_pre[0] - p[1];
|
|
55 p_pre[2] = -p_pre[1] - p[2];
|
|
56 p_pre[3] = -p_pre[2] - p[3];
|
|
57 p_pre[4] = -p_pre[3] - p[4];
|
|
58 p_pre[4] = p_pre[4] / 2;
|
|
59
|
|
60 q_pre[0] = (float) 1.0 - q[0];
|
|
61 q_pre[1] = q_pre[0] - q[1];
|
|
62 q_pre[2] = q_pre[1] - q[2];
|
|
63 q_pre[3] = q_pre[2] - q[3];
|
|
64 q_pre[4] = q_pre[3] - q[4];
|
|
65 q_pre[4] = q_pre[4] / 2;
|
|
66
|
|
67 omega = 0.0;
|
|
68
|
|
69
|
|
70
|
|
71
|
|
72
|
|
73 old_omega = 0.0;
|
|
74
|
|
75 old_p = FLOAT_MAX;
|
|
76 old_q = FLOAT_MAX;
|
|
77
|
|
78 /* Here we loop through lsp_index to find all the
|
|
79 LPC_FILTERORDER roots for omega. */
|
|
80
|
|
81 for (lsp_index = 0; lsp_index < LPC_FILTERORDER; lsp_index++) {
|
|
82
|
|
83 /* Depending on lsp_index being even or odd, we
|
|
84 alternatively solve the roots for the two LSP equations. */
|
|
85
|
|
86
|
|
87 if ((lsp_index & 0x1) == 0) {
|
|
88 pq_coef = p_pre;
|
|
89 old = &old_p;
|
|
90 } else {
|
|
91 pq_coef = q_pre;
|
|
92 old = &old_q;
|
|
93 }
|
|
94
|
|
95 /* Start with low resolution grid */
|
|
96
|
|
97 for (step_idx = 0, step = steps[step_idx];
|
|
98 step_idx < LSF_NUMBER_OF_STEPS;) {
|
|
99
|
|
100 /* cos(10piw) + pq(0)cos(8piw) + pq(1)cos(6piw) +
|
|
101 pq(2)cos(4piw) + pq(3)cod(2piw) + pq(4) */
|
|
102
|
|
103 hlp = (float) cos(omega * TWO_PI);
|
|
104 hlp1 = (float) 2.0 *hlp + pq_coef[0];
|
|
105 hlp2 = (float) 2.0 *hlp * hlp1 - (float) 1.0 + pq_coef[1];
|
|
106 hlp3 = (float) 2.0 *hlp * hlp2 - hlp1 + pq_coef[2];
|
|
107 hlp4 = (float) 2.0 *hlp * hlp3 - hlp2 + pq_coef[3];
|
|
108 hlp5 = hlp * hlp4 - hlp3 + pq_coef[4];
|
|
109
|
|
110
|
|
111 if (((hlp5 * (*old)) <= 0.0) || (omega >= 0.5)) {
|
|
112
|
|
113 if (step_idx == (LSF_NUMBER_OF_STEPS - 1)) {
|
|
114
|
|
115 if (fabs(hlp5) >= fabs(*old)) {
|
|
116 freq[lsp_index] = omega - step;
|
|
117 } else {
|
|
118 freq[lsp_index] = omega;
|
|
119 }
|
|
120
|
|
121
|
|
122
|
|
123
|
|
124
|
|
125
|
|
126
|
|
127 if ((*old) >= 0.0) {
|
|
128 *old = (float) -1.0 * FLOAT_MAX;
|
|
129 } else {
|
|
130 *old = FLOAT_MAX;
|
|
131 }
|
|
132
|
|
133 omega = old_omega;
|
|
134 step_idx = 0;
|
|
135
|
|
136 step_idx = LSF_NUMBER_OF_STEPS;
|
|
137 } else {
|
|
138
|
|
139 if (step_idx == 0) {
|
|
140 old_omega = omega;
|
|
141 }
|
|
142
|
|
143 step_idx++;
|
|
144 omega -= steps[step_idx];
|
|
145
|
|
146 /* Go back one grid step */
|
|
147
|
|
148 step = steps[step_idx];
|
|
149 }
|
|
150 } else {
|
|
151
|
|
152 /* increment omega until they are of different sign,
|
|
153 and we know there is at least one root between omega
|
|
154 and old_omega */
|
|
155 *old = hlp5;
|
|
156 omega += step;
|
|
157 }
|
|
158 }
|
|
159 }
|
|
160
|
|
161 for (i = 0; i < LPC_FILTERORDER; i++) {
|
|
162 freq[i] = freq[i] * TWO_PI;
|
|
163 }
|
|
164 }
|
|
165
|
|
166 /*----------------------------------------------------------------*
|
|
167 * conversion from lsf coefficients to lpc coefficients
|
|
168 *---------------------------------------------------------------*/
|
|
169
|
|
170 void lsf2a(float *a_coef, /* (o) lpc coefficients */
|
|
171 float *freq /* (i) lsf coefficients */
|
|
172 )
|
|
173 {
|
|
174 int i, j;
|
|
175 float hlp;
|
|
176 float p[LPC_HALFORDER], q[LPC_HALFORDER];
|
|
177 float a[LPC_HALFORDER + 1], a1[LPC_HALFORDER], a2[LPC_HALFORDER];
|
|
178 float b[LPC_HALFORDER + 1], b1[LPC_HALFORDER], b2[LPC_HALFORDER];
|
|
179
|
|
180 for (i = 0; i < LPC_FILTERORDER; i++) {
|
|
181 freq[i] = freq[i] * PI2;
|
|
182 }
|
|
183
|
|
184 /* Check input for ill-conditioned cases. This part is not
|
|
185 found in the TIA standard. It involves the following 2 IF
|
|
186 blocks. If "freq" is judged ill-conditioned, then we first
|
|
187 modify freq[0] and freq[LPC_HALFORDER-1] (normally
|
|
188 LPC_HALFORDER = 10 for LPC applications), then we adjust
|
|
189 the other "freq" values slightly */
|
|
190
|
|
191
|
|
192 if ((freq[0] <= 0.0) || (freq[LPC_FILTERORDER - 1] >= 0.5)) {
|
|
193
|
|
194
|
|
195 if (freq[0] <= 0.0) {
|
|
196 freq[0] = (float) 0.022;
|
|
197 }
|
|
198
|
|
199
|
|
200 if (freq[LPC_FILTERORDER - 1] >= 0.5) {
|
|
201 freq[LPC_FILTERORDER - 1] = (float) 0.499;
|
|
202 }
|
|
203
|
|
204 hlp = (freq[LPC_FILTERORDER - 1] - freq[0]) /
|
|
205 (float) (LPC_FILTERORDER - 1);
|
|
206
|
|
207 for (i = 1; i < LPC_FILTERORDER; i++) {
|
|
208 freq[i] = freq[i - 1] + hlp;
|
|
209 }
|
|
210 }
|
|
211
|
|
212 memset(a1, 0, LPC_HALFORDER * sizeof(float));
|
|
213 memset(a2, 0, LPC_HALFORDER * sizeof(float));
|
|
214 memset(b1, 0, LPC_HALFORDER * sizeof(float));
|
|
215 memset(b2, 0, LPC_HALFORDER * sizeof(float));
|
|
216 memset(a, 0, (LPC_HALFORDER + 1) * sizeof(float));
|
|
217 memset(b, 0, (LPC_HALFORDER + 1) * sizeof(float));
|
|
218
|
|
219
|
|
220
|
|
221
|
|
222
|
|
223
|
|
224 /* p[i] and q[i] compute cos(2*pi*omega_{2j}) and
|
|
225 cos(2*pi*omega_{2j-1} in eqs. 4.2.2.2-1 and 4.2.2.2-2.
|
|
226 Note that for this code p[i] specifies the coefficients
|
|
227 used in .Q_A(z) while q[i] specifies the coefficients used
|
|
228 in .P_A(z) */
|
|
229
|
|
230 for (i = 0; i < LPC_HALFORDER; i++) {
|
|
231 p[i] = (float) cos(TWO_PI * freq[2 * i]);
|
|
232 q[i] = (float) cos(TWO_PI * freq[2 * i + 1]);
|
|
233 }
|
|
234
|
|
235 a[0] = 0.25;
|
|
236 b[0] = 0.25;
|
|
237
|
|
238 for (i = 0; i < LPC_HALFORDER; i++) {
|
|
239 a[i + 1] = a[i] - 2 * p[i] * a1[i] + a2[i];
|
|
240 b[i + 1] = b[i] - 2 * q[i] * b1[i] + b2[i];
|
|
241 a2[i] = a1[i];
|
|
242 a1[i] = a[i];
|
|
243 b2[i] = b1[i];
|
|
244 b1[i] = b[i];
|
|
245 }
|
|
246
|
|
247 for (j = 0; j < LPC_FILTERORDER; j++) {
|
|
248
|
|
249 if (j == 0) {
|
|
250 a[0] = 0.25;
|
|
251 b[0] = -0.25;
|
|
252 } else {
|
|
253 a[0] = b[0] = 0.0;
|
|
254 }
|
|
255
|
|
256 for (i = 0; i < LPC_HALFORDER; i++) {
|
|
257 a[i + 1] = a[i] - 2 * p[i] * a1[i] + a2[i];
|
|
258 b[i + 1] = b[i] - 2 * q[i] * b1[i] + b2[i];
|
|
259 a2[i] = a1[i];
|
|
260 a1[i] = a[i];
|
|
261 b2[i] = b1[i];
|
|
262 b1[i] = b[i];
|
|
263 }
|
|
264
|
|
265 a_coef[j + 1] = 2 * (a[LPC_HALFORDER] + b[LPC_HALFORDER]);
|
|
266 }
|
|
267
|
|
268 a_coef[0] = 1.0;
|
|
269 }
|