1 /*
2 ** Copyright 2003-2010, VisualOn, Inc.
3 **
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
7 **
8 ** http://www.apache.org/licenses/LICENSE-2.0
9 **
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 express or implied.
13 ** See the License for the specific language governing permissions and
14 ** limitations under the License.
15 */
16
17 /***********************************************************************
18 * File: isp_az.c *
19 * *
20 * Description:Compute the LPC coefficients from isp (order=M) *
21 * *
22 ************************************************************************/
23
24 #include "typedef.h"
25 #include "basic_op.h"
26 #include "oper_32b.h"
27 #include "cnst.h"
28
29 #define NC (M/2)
30 #define NC16k (M16k/2)
31
32 /* local function */
33
34 static void Get_isp_pol(Word16 * isp, Word32 * f, Word16 n);
35 static void Get_isp_pol_16kHz(Word16 * isp, Word32 * f, Word16 n);
36
Isp_Az(Word16 isp[],Word16 a[],Word16 m,Word16 adaptive_scaling)37 void Isp_Az(
38 Word16 isp[], /* (i) Q15 : Immittance spectral pairs */
39 Word16 a[], /* (o) Q12 : predictor coefficients (order = M) */
40 Word16 m,
41 Word16 adaptive_scaling /* (i) 0 : adaptive scaling disabled */
42 /* 1 : adaptive scaling enabled */
43 )
44 {
45 Word32 i, j;
46 Word16 hi, lo;
47 Word32 f1[NC16k + 1], f2[NC16k];
48 Word16 nc;
49 Word32 t0;
50 Word16 q, q_sug;
51 Word32 tmax;
52
53 nc = (m >> 1);
54 if(nc > 8)
55 {
56 Get_isp_pol_16kHz(&isp[0], f1, nc);
57 for (i = 0; i <= nc; i++)
58 {
59 f1[i] = f1[i] << 2;
60 }
61 } else
62 Get_isp_pol(&isp[0], f1, nc);
63
64 if (nc > 8)
65 {
66 Get_isp_pol_16kHz(&isp[1], f2, (nc - 1));
67 for (i = 0; i <= nc - 1; i++)
68 {
69 f2[i] = f2[i] << 2;
70 }
71 } else
72 Get_isp_pol(&isp[1], f2, (nc - 1));
73
74 /*-----------------------------------------------------*
75 * Multiply F2(z) by (1 - z^-2) *
76 *-----------------------------------------------------*/
77
78 for (i = (nc - 1); i > 1; i--)
79 {
80 f2[i] = vo_L_sub(f2[i], f2[i - 2]); /* f2[i] -= f2[i-2]; */
81 }
82
83 /*----------------------------------------------------------*
84 * Scale F1(z) by (1+isp[m-1]) and F2(z) by (1-isp[m-1]) *
85 *----------------------------------------------------------*/
86
87 for (i = 0; i < nc; i++)
88 {
89 /* f1[i] *= (1.0 + isp[M-1]); */
90
91 hi = f1[i] >> 16;
92 lo = (f1[i] & 0xffff)>>1;
93
94 t0 = Mpy_32_16(hi, lo, isp[m - 1]);
95 f1[i] = vo_L_add(f1[i], t0);
96
97 /* f2[i] *= (1.0 - isp[M-1]); */
98
99 hi = f2[i] >> 16;
100 lo = (f2[i] & 0xffff)>>1;
101 t0 = Mpy_32_16(hi, lo, isp[m - 1]);
102 f2[i] = vo_L_sub(f2[i], t0);
103 }
104
105 /*-----------------------------------------------------*
106 * A(z) = (F1(z)+F2(z))/2 *
107 * F1(z) is symmetric and F2(z) is antisymmetric *
108 *-----------------------------------------------------*/
109
110 /* a[0] = 1.0; */
111 a[0] = 4096;
112 tmax = 1;
113 for (i = 1, j = m - 1; i < nc; i++, j--)
114 {
115 /* a[i] = 0.5*(f1[i] + f2[i]); */
116
117 t0 = vo_L_add(f1[i], f2[i]); /* f1[i] + f2[i] */
118 tmax |= L_abs(t0);
119 a[i] = (Word16)(vo_L_shr_r(t0, 12)); /* from Q23 to Q12 and * 0.5 */
120
121 /* a[j] = 0.5*(f1[i] - f2[i]); */
122
123 t0 = vo_L_sub(f1[i], f2[i]); /* f1[i] - f2[i] */
124 tmax |= L_abs(t0);
125 a[j] = (Word16)(vo_L_shr_r(t0, 12)); /* from Q23 to Q12 and * 0.5 */
126 }
127
128 /* rescale data if overflow has occured and reprocess the loop */
129 if(adaptive_scaling == 1)
130 q = 4 - norm_l(tmax); /* adaptive scaling enabled */
131 else
132 q = 0; /* adaptive scaling disabled */
133
134 if (q > 0)
135 {
136 q_sug = (12 + q);
137 for (i = 1, j = m - 1; i < nc; i++, j--)
138 {
139 /* a[i] = 0.5*(f1[i] + f2[i]); */
140 t0 = vo_L_add(f1[i], f2[i]); /* f1[i] + f2[i] */
141 a[i] = (Word16)(vo_L_shr_r(t0, q_sug)); /* from Q23 to Q12 and * 0.5 */
142
143 /* a[j] = 0.5*(f1[i] - f2[i]); */
144 t0 = vo_L_sub(f1[i], f2[i]); /* f1[i] - f2[i] */
145 a[j] = (Word16)(vo_L_shr_r(t0, q_sug)); /* from Q23 to Q12 and * 0.5 */
146 }
147 a[0] = shr(a[0], q);
148 }
149 else
150 {
151 q_sug = 12;
152 q = 0;
153 }
154 /* a[NC] = 0.5*f1[NC]*(1.0 + isp[M-1]); */
155 hi = f1[nc] >> 16;
156 lo = (f1[nc] & 0xffff)>>1;
157 t0 = Mpy_32_16(hi, lo, isp[m - 1]);
158 t0 = vo_L_add(f1[nc], t0);
159 a[nc] = (Word16)(L_shr_r(t0, q_sug)); /* from Q23 to Q12 and * 0.5 */
160 /* a[m] = isp[m-1]; */
161
162 a[m] = vo_shr_r(isp[m - 1], (3 + q)); /* from Q15 to Q12 */
163 return;
164 }
165
166 /*-----------------------------------------------------------*
167 * procedure Get_isp_pol: *
168 * ~~~~~~~~~~~ *
169 * Find the polynomial F1(z) or F2(z) from the ISPs. *
170 * This is performed by expanding the product polynomials: *
171 * *
172 * F1(z) = product ( 1 - 2 isp_i z^-1 + z^-2 ) *
173 * i=0,2,4,6,8 *
174 * F2(z) = product ( 1 - 2 isp_i z^-1 + z^-2 ) *
175 * i=1,3,5,7 *
176 * *
177 * where isp_i are the ISPs in the cosine domain. *
178 *-----------------------------------------------------------*
179 * *
180 * Parameters: *
181 * isp[] : isp vector (cosine domaine) in Q15 *
182 * f[] : the coefficients of F1 or F2 in Q23 *
183 * n : == NC for F1(z); == NC-1 for F2(z) *
184 *-----------------------------------------------------------*/
185
Get_isp_pol(Word16 * isp,Word32 * f,Word16 n)186 static void Get_isp_pol(Word16 * isp, Word32 * f, Word16 n)
187 {
188 Word16 hi, lo;
189 Word32 i, j, t0;
190 /* All computation in Q23 */
191
192 f[0] = vo_L_mult(4096, 1024); /* f[0] = 1.0; in Q23 */
193 f[1] = vo_L_mult(isp[0], -256); /* f[1] = -2.0*isp[0] in Q23 */
194
195 f += 2; /* Advance f pointer */
196 isp += 2; /* Advance isp pointer */
197 for (i = 2; i <= n; i++)
198 {
199 *f = f[-2];
200 for (j = 1; j < i; j++, f--)
201 {
202 hi = f[-1]>>16;
203 lo = (f[-1] & 0xffff)>>1;
204
205 t0 = Mpy_32_16(hi, lo, *isp); /* t0 = f[-1] * isp */
206 t0 = t0 << 1;
207 *f = vo_L_sub(*f, t0); /* *f -= t0 */
208 *f = vo_L_add(*f, f[-2]); /* *f += f[-2] */
209 }
210 *f -= (*isp << 9); /* *f -= isp<<8 */
211 f += i; /* Advance f pointer */
212 isp += 2; /* Advance isp pointer */
213 }
214 return;
215 }
216
Get_isp_pol_16kHz(Word16 * isp,Word32 * f,Word16 n)217 static void Get_isp_pol_16kHz(Word16 * isp, Word32 * f, Word16 n)
218 {
219 Word16 hi, lo;
220 Word32 i, j, t0;
221
222 /* All computation in Q23 */
223 f[0] = L_mult(4096, 256); /* f[0] = 1.0; in Q23 */
224 f[1] = L_mult(isp[0], -64); /* f[1] = -2.0*isp[0] in Q23 */
225
226 f += 2; /* Advance f pointer */
227 isp += 2; /* Advance isp pointer */
228
229 for (i = 2; i <= n; i++)
230 {
231 *f = f[-2];
232 for (j = 1; j < i; j++, f--)
233 {
234 VO_L_Extract(f[-1], &hi, &lo);
235 t0 = Mpy_32_16(hi, lo, *isp); /* t0 = f[-1] * isp */
236 t0 = L_shl2(t0, 1);
237 *f = L_sub(*f, t0); /* *f -= t0 */
238 *f = L_add(*f, f[-2]); /* *f += f[-2] */
239 }
240 *f = L_msu(*f, *isp, 64); /* *f -= isp<<8 */
241 f += i; /* Advance f pointer */
242 isp += 2; /* Advance isp pointer */
243 }
244 return;
245 }
246
247
248