Mercurial > hg > audiostuff

comparison intercom/ilbc/lsf.c @ 2:13be24d74cd2

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
import intercom-0.4.1
author Peter Meerwald <pmeerw@cosy.sbg.ac.at>
date Fri, 25 Jun 2010 09:57:52 +0200
parents
children
comparison
equal deleted inserted replaced
1:9cadc470e3da 2:13be24d74cd2
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 }

Repositories maintained by Peter Meerwald, pmeerw@pmeerw.net.