1 /*
2 * Copyright (c) 2012 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 * lpc_masking_model.c
13 *
14 * LPC analysis and filtering functions
15 *
16 */
17
18 #include "lpc_masking_model.h"
19
20 #include <limits.h> /* For LLONG_MAX and LLONG_MIN. */
21 #include "codec.h"
22 #include "entropy_coding.h"
23 #include "settings.h"
24
25 /* The conversion is implemented by the step-down algorithm */
WebRtcSpl_AToK_JSK(int16_t * a16,int16_t useOrder,int16_t * k16)26 void WebRtcSpl_AToK_JSK(
27 int16_t *a16, /* Q11 */
28 int16_t useOrder,
29 int16_t *k16 /* Q15 */
30 )
31 {
32 int m, k;
33 int32_t tmp32[MAX_AR_MODEL_ORDER];
34 int32_t tmp32b;
35 int32_t tmp_inv_denum32;
36 int16_t tmp_inv_denum16;
37
38 k16[useOrder-1]= WEBRTC_SPL_LSHIFT_W16(a16[useOrder], 4); //Q11<<4 => Q15
39
40 for (m=useOrder-1; m>0; m--) {
41 tmp_inv_denum32 = ((int32_t) 1073741823) - WEBRTC_SPL_MUL_16_16(k16[m], k16[m]); // (1 - k^2) in Q30
42 tmp_inv_denum16 = (int16_t) WEBRTC_SPL_RSHIFT_W32(tmp_inv_denum32, 15); // (1 - k^2) in Q15
43
44 for (k=1; k<=m; k++) {
45 tmp32b = WEBRTC_SPL_LSHIFT_W32((int32_t)a16[k], 16) -
46 WEBRTC_SPL_LSHIFT_W32(WEBRTC_SPL_MUL_16_16(k16[m], a16[m-k+1]), 1);
47
48 tmp32[k] = WebRtcSpl_DivW32W16(tmp32b, tmp_inv_denum16); //Q27/Q15 = Q12
49 }
50
51 for (k=1; k<m; k++) {
52 a16[k] = (int16_t) WEBRTC_SPL_RSHIFT_W32(tmp32[k], 1); //Q12>>1 => Q11
53 }
54
55 tmp32[m] = WEBRTC_SPL_SAT(4092, tmp32[m], -4092);
56 k16[m-1] = (int16_t) WEBRTC_SPL_LSHIFT_W32(tmp32[m], 3); //Q12<<3 => Q15
57 }
58
59 return;
60 }
61
62
63
64
65
WebRtcSpl_LevinsonW32_JSK(int32_t * R,int16_t * A,int16_t * K,int16_t order)66 int16_t WebRtcSpl_LevinsonW32_JSK(
67 int32_t *R, /* (i) Autocorrelation of length >= order+1 */
68 int16_t *A, /* (o) A[0..order] LPC coefficients (Q11) */
69 int16_t *K, /* (o) K[0...order-1] Reflection coefficients (Q15) */
70 int16_t order /* (i) filter order */
71 ) {
72 int16_t i, j;
73 int16_t R_hi[LEVINSON_MAX_ORDER+1], R_low[LEVINSON_MAX_ORDER+1];
74 /* Aurocorr coefficients in high precision */
75 int16_t A_hi[LEVINSON_MAX_ORDER+1], A_low[LEVINSON_MAX_ORDER+1];
76 /* LPC coefficients in high precicion */
77 int16_t A_upd_hi[LEVINSON_MAX_ORDER+1], A_upd_low[LEVINSON_MAX_ORDER+1];
78 /* LPC coefficients for next iteration */
79 int16_t K_hi, K_low; /* reflection coefficient in high precision */
80 int16_t Alpha_hi, Alpha_low, Alpha_exp; /* Prediction gain Alpha in high precision
81 and with scale factor */
82 int16_t tmp_hi, tmp_low;
83 int32_t temp1W32, temp2W32, temp3W32;
84 int16_t norm;
85
86 /* Normalize the autocorrelation R[0]...R[order+1] */
87
88 norm = WebRtcSpl_NormW32(R[0]);
89
90 for (i=order;i>=0;i--) {
91 temp1W32 = WEBRTC_SPL_LSHIFT_W32(R[i], norm);
92 /* Put R in hi and low format */
93 R_hi[i] = (int16_t) WEBRTC_SPL_RSHIFT_W32(temp1W32, 16);
94 R_low[i] = (int16_t)WEBRTC_SPL_RSHIFT_W32((temp1W32 - WEBRTC_SPL_LSHIFT_W32((int32_t)R_hi[i], 16)), 1);
95 }
96
97 /* K = A[1] = -R[1] / R[0] */
98
99 temp2W32 = WEBRTC_SPL_LSHIFT_W32((int32_t)R_hi[1],16) +
100 WEBRTC_SPL_LSHIFT_W32((int32_t)R_low[1],1); /* R[1] in Q31 */
101 temp3W32 = WEBRTC_SPL_ABS_W32(temp2W32); /* abs R[1] */
102 temp1W32 = WebRtcSpl_DivW32HiLow(temp3W32, R_hi[0], R_low[0]); /* abs(R[1])/R[0] in Q31 */
103 /* Put back the sign on R[1] */
104 if (temp2W32 > 0) {
105 temp1W32 = -temp1W32;
106 }
107
108 /* Put K in hi and low format */
109 K_hi = (int16_t) WEBRTC_SPL_RSHIFT_W32(temp1W32, 16);
110 K_low = (int16_t)WEBRTC_SPL_RSHIFT_W32((temp1W32 - WEBRTC_SPL_LSHIFT_W32((int32_t)K_hi, 16)), 1);
111
112 /* Store first reflection coefficient */
113 K[0] = K_hi;
114
115 temp1W32 = WEBRTC_SPL_RSHIFT_W32(temp1W32, 4); /* A[1] in Q27 */
116
117 /* Put A[1] in hi and low format */
118 A_hi[1] = (int16_t) WEBRTC_SPL_RSHIFT_W32(temp1W32, 16);
119 A_low[1] = (int16_t)WEBRTC_SPL_RSHIFT_W32((temp1W32 - WEBRTC_SPL_LSHIFT_W32((int32_t)A_hi[1], 16)), 1);
120
121 /* Alpha = R[0] * (1-K^2) */
122
123 temp1W32 = WEBRTC_SPL_LSHIFT_W32((WEBRTC_SPL_RSHIFT_W32(WEBRTC_SPL_MUL_16_16(K_hi, K_low), 14) +
124 WEBRTC_SPL_MUL_16_16(K_hi, K_hi)), 1); /* temp1W32 = k^2 in Q31 */
125
126 temp1W32 = WEBRTC_SPL_ABS_W32(temp1W32); /* Guard against <0 */
127 temp1W32 = (int32_t)0x7fffffffL - temp1W32; /* temp1W32 = (1 - K[0]*K[0]) in Q31 */
128
129 /* Store temp1W32 = 1 - K[0]*K[0] on hi and low format */
130 tmp_hi = (int16_t) WEBRTC_SPL_RSHIFT_W32(temp1W32, 16);
131 tmp_low = (int16_t)WEBRTC_SPL_RSHIFT_W32((temp1W32 - WEBRTC_SPL_LSHIFT_W32((int32_t)tmp_hi, 16)), 1);
132
133 /* Calculate Alpha in Q31 */
134 temp1W32 = WEBRTC_SPL_LSHIFT_W32((WEBRTC_SPL_MUL_16_16(R_hi[0], tmp_hi) +
135 WEBRTC_SPL_RSHIFT_W32(WEBRTC_SPL_MUL_16_16(R_hi[0], tmp_low), 15) +
136 WEBRTC_SPL_RSHIFT_W32(WEBRTC_SPL_MUL_16_16(R_low[0], tmp_hi), 15) ), 1);
137
138 /* Normalize Alpha and put it in hi and low format */
139
140 Alpha_exp = WebRtcSpl_NormW32(temp1W32);
141 temp1W32 = WEBRTC_SPL_LSHIFT_W32(temp1W32, Alpha_exp);
142 Alpha_hi = (int16_t) WEBRTC_SPL_RSHIFT_W32(temp1W32, 16);
143 Alpha_low = (int16_t)WEBRTC_SPL_RSHIFT_W32((temp1W32 - WEBRTC_SPL_LSHIFT_W32((int32_t)Alpha_hi, 16)), 1);
144
145 /* Perform the iterative calculations in the
146 Levinson Durbin algorithm */
147
148 for (i=2; i<=order; i++)
149 {
150
151 /* ----
152 \
153 temp1W32 = R[i] + > R[j]*A[i-j]
154 /
155 ----
156 j=1..i-1
157 */
158
159 temp1W32 = 0;
160
161 for(j=1; j<i; j++) {
162 /* temp1W32 is in Q31 */
163 temp1W32 += (WEBRTC_SPL_LSHIFT_W32(WEBRTC_SPL_MUL_16_16(R_hi[j], A_hi[i-j]), 1) +
164 WEBRTC_SPL_LSHIFT_W32(( WEBRTC_SPL_RSHIFT_W32(WEBRTC_SPL_MUL_16_16(R_hi[j], A_low[i-j]), 15) +
165 WEBRTC_SPL_RSHIFT_W32(WEBRTC_SPL_MUL_16_16(R_low[j], A_hi[i-j]), 15) ), 1));
166 }
167
168 temp1W32 = WEBRTC_SPL_LSHIFT_W32(temp1W32, 4);
169 temp1W32 += (WEBRTC_SPL_LSHIFT_W32((int32_t)R_hi[i], 16) +
170 WEBRTC_SPL_LSHIFT_W32((int32_t)R_low[i], 1));
171
172 /* K = -temp1W32 / Alpha */
173 temp2W32 = WEBRTC_SPL_ABS_W32(temp1W32); /* abs(temp1W32) */
174 temp3W32 = WebRtcSpl_DivW32HiLow(temp2W32, Alpha_hi, Alpha_low); /* abs(temp1W32)/Alpha */
175
176 /* Put the sign of temp1W32 back again */
177 if (temp1W32 > 0) {
178 temp3W32 = -temp3W32;
179 }
180
181 /* Use the Alpha shifts from earlier to denormalize */
182 norm = WebRtcSpl_NormW32(temp3W32);
183 if ((Alpha_exp <= norm)||(temp3W32==0)) {
184 temp3W32 = WEBRTC_SPL_LSHIFT_W32(temp3W32, Alpha_exp);
185 } else {
186 if (temp3W32 > 0)
187 {
188 temp3W32 = (int32_t)0x7fffffffL;
189 } else
190 {
191 temp3W32 = (int32_t)0x80000000L;
192 }
193 }
194
195 /* Put K on hi and low format */
196 K_hi = (int16_t) WEBRTC_SPL_RSHIFT_W32(temp3W32, 16);
197 K_low = (int16_t)WEBRTC_SPL_RSHIFT_W32((temp3W32 - WEBRTC_SPL_LSHIFT_W32((int32_t)K_hi, 16)), 1);
198
199 /* Store Reflection coefficient in Q15 */
200 K[i-1] = K_hi;
201
202 /* Test for unstable filter. If unstable return 0 and let the
203 user decide what to do in that case
204 */
205
206 if ((int32_t)WEBRTC_SPL_ABS_W16(K_hi) > (int32_t)32740) {
207 return(-i); /* Unstable filter */
208 }
209
210 /*
211 Compute updated LPC coefficient: Anew[i]
212 Anew[j]= A[j] + K*A[i-j] for j=1..i-1
213 Anew[i]= K
214 */
215
216 for(j=1; j<i; j++)
217 {
218 temp1W32 = WEBRTC_SPL_LSHIFT_W32((int32_t)A_hi[j],16) +
219 WEBRTC_SPL_LSHIFT_W32((int32_t)A_low[j],1); /* temp1W32 = A[j] in Q27 */
220
221 temp1W32 += WEBRTC_SPL_LSHIFT_W32(( WEBRTC_SPL_MUL_16_16(K_hi, A_hi[i-j]) +
222 WEBRTC_SPL_RSHIFT_W32(WEBRTC_SPL_MUL_16_16(K_hi, A_low[i-j]), 15) +
223 WEBRTC_SPL_RSHIFT_W32(WEBRTC_SPL_MUL_16_16(K_low, A_hi[i-j]), 15) ), 1); /* temp1W32 += K*A[i-j] in Q27 */
224
225 /* Put Anew in hi and low format */
226 A_upd_hi[j] = (int16_t) WEBRTC_SPL_RSHIFT_W32(temp1W32, 16);
227 A_upd_low[j] = (int16_t)WEBRTC_SPL_RSHIFT_W32((temp1W32 - WEBRTC_SPL_LSHIFT_W32((int32_t)A_upd_hi[j], 16)), 1);
228 }
229
230 temp3W32 = WEBRTC_SPL_RSHIFT_W32(temp3W32, 4); /* temp3W32 = K in Q27 (Convert from Q31 to Q27) */
231
232 /* Store Anew in hi and low format */
233 A_upd_hi[i] = (int16_t) WEBRTC_SPL_RSHIFT_W32(temp3W32, 16);
234 A_upd_low[i] = (int16_t)WEBRTC_SPL_RSHIFT_W32((temp3W32 - WEBRTC_SPL_LSHIFT_W32((int32_t)A_upd_hi[i], 16)), 1);
235
236 /* Alpha = Alpha * (1-K^2) */
237
238 temp1W32 = WEBRTC_SPL_LSHIFT_W32((WEBRTC_SPL_RSHIFT_W32(WEBRTC_SPL_MUL_16_16(K_hi, K_low), 14) +
239 WEBRTC_SPL_MUL_16_16(K_hi, K_hi)), 1); /* K*K in Q31 */
240
241 temp1W32 = WEBRTC_SPL_ABS_W32(temp1W32); /* Guard against <0 */
242 temp1W32 = (int32_t)0x7fffffffL - temp1W32; /* 1 - K*K in Q31 */
243
244 /* Convert 1- K^2 in hi and low format */
245 tmp_hi = (int16_t) WEBRTC_SPL_RSHIFT_W32(temp1W32, 16);
246 tmp_low = (int16_t)WEBRTC_SPL_RSHIFT_W32((temp1W32 - WEBRTC_SPL_LSHIFT_W32((int32_t)tmp_hi, 16)), 1);
247
248 /* Calculate Alpha = Alpha * (1-K^2) in Q31 */
249 temp1W32 = WEBRTC_SPL_LSHIFT_W32(( WEBRTC_SPL_MUL_16_16(Alpha_hi, tmp_hi) +
250 WEBRTC_SPL_RSHIFT_W32(WEBRTC_SPL_MUL_16_16(Alpha_hi, tmp_low), 15) +
251 WEBRTC_SPL_RSHIFT_W32(WEBRTC_SPL_MUL_16_16(Alpha_low, tmp_hi), 15)), 1);
252
253 /* Normalize Alpha and store it on hi and low format */
254
255 norm = WebRtcSpl_NormW32(temp1W32);
256 temp1W32 = WEBRTC_SPL_LSHIFT_W32(temp1W32, norm);
257
258 Alpha_hi = (int16_t) WEBRTC_SPL_RSHIFT_W32(temp1W32, 16);
259 Alpha_low = (int16_t)WEBRTC_SPL_RSHIFT_W32((temp1W32 - WEBRTC_SPL_LSHIFT_W32((int32_t)Alpha_hi, 16)), 1);
260
261 /* Update the total nomalization of Alpha */
262 Alpha_exp = Alpha_exp + norm;
263
264 /* Update A[] */
265
266 for(j=1; j<=i; j++)
267 {
268 A_hi[j] =A_upd_hi[j];
269 A_low[j] =A_upd_low[j];
270 }
271 }
272
273 /*
274 Set A[0] to 1.0 and store the A[i] i=1...order in Q12
275 (Convert from Q27 and use rounding)
276 */
277
278 A[0] = 2048;
279
280 for(i=1; i<=order; i++) {
281 /* temp1W32 in Q27 */
282 temp1W32 = WEBRTC_SPL_LSHIFT_W32((int32_t)A_hi[i], 16) +
283 WEBRTC_SPL_LSHIFT_W32((int32_t)A_low[i], 1);
284 /* Round and store upper word */
285 A[i] = (int16_t)WEBRTC_SPL_RSHIFT_W32(temp1W32+(int32_t)32768, 16);
286 }
287 return(1); /* Stable filters */
288 }
289
290
291
292
293
294 /* window */
295 /* Matlab generation of floating point code:
296 * t = (1:256)/257; r = 1-(1-t).^.45; w = sin(r*pi).^3; w = w/sum(w); plot((1:256)/8, w); grid;
297 * for k=1:16, fprintf(1, '%.8f, ', w(k*16 + (-15:0))); fprintf(1, '\n'); end
298 * All values are multiplyed with 2^21 in fixed point code.
299 */
300 static const int16_t kWindowAutocorr[WINLEN] = {
301 0, 0, 0, 0, 0, 1, 1, 2, 2, 3, 5, 6,
302 8, 10, 12, 14, 17, 20, 24, 28, 33, 38, 43, 49,
303 56, 63, 71, 79, 88, 98, 108, 119, 131, 143, 157, 171,
304 186, 202, 219, 237, 256, 275, 296, 318, 341, 365, 390, 416,
305 444, 472, 502, 533, 566, 600, 635, 671, 709, 748, 789, 831,
306 875, 920, 967, 1015, 1065, 1116, 1170, 1224, 1281, 1339, 1399, 1461,
307 1525, 1590, 1657, 1726, 1797, 1870, 1945, 2021, 2100, 2181, 2263, 2348,
308 2434, 2523, 2614, 2706, 2801, 2898, 2997, 3099, 3202, 3307, 3415, 3525,
309 3637, 3751, 3867, 3986, 4106, 4229, 4354, 4481, 4611, 4742, 4876, 5012,
310 5150, 5291, 5433, 5578, 5725, 5874, 6025, 6178, 6333, 6490, 6650, 6811,
311 6974, 7140, 7307, 7476, 7647, 7820, 7995, 8171, 8349, 8529, 8711, 8894,
312 9079, 9265, 9453, 9642, 9833, 10024, 10217, 10412, 10607, 10803, 11000, 11199,
313 11398, 11597, 11797, 11998, 12200, 12401, 12603, 12805, 13008, 13210, 13412, 13614,
314 13815, 14016, 14216, 14416, 14615, 14813, 15009, 15205, 15399, 15591, 15782, 15971,
315 16157, 16342, 16524, 16704, 16881, 17056, 17227, 17395, 17559, 17720, 17877, 18030,
316 18179, 18323, 18462, 18597, 18727, 18851, 18970, 19082, 19189, 19290, 19384, 19471,
317 19551, 19623, 19689, 19746, 19795, 19835, 19867, 19890, 19904, 19908, 19902, 19886,
318 19860, 19823, 19775, 19715, 19644, 19561, 19465, 19357, 19237, 19102, 18955, 18793,
319 18618, 18428, 18223, 18004, 17769, 17518, 17252, 16970, 16672, 16357, 16025, 15677,
320 15311, 14929, 14529, 14111, 13677, 13225, 12755, 12268, 11764, 11243, 10706, 10152,
321 9583, 8998, 8399, 7787, 7162, 6527, 5883, 5231, 4576, 3919, 3265, 2620,
322 1990, 1386, 825, 333
323 };
324
325
326 /* By using a hearing threshold level in dB of -28 dB (higher value gives more noise),
327 the H_T_H (in float) can be calculated as:
328 H_T_H = pow(10.0, 0.05 * (-28.0)) = 0.039810717055350
329 In Q19, H_T_H becomes round(0.039810717055350*2^19) ~= 20872, i.e.
330 H_T_H = 20872/524288.0, and H_T_HQ19 = 20872;
331 */
332
333
334 /* The bandwidth expansion vectors are created from:
335 kPolyVecLo=[0.900000,0.810000,0.729000,0.656100,0.590490,0.531441,0.478297,0.430467,0.387420,0.348678,0.313811,0.282430];
336 kPolyVecHi=[0.800000,0.640000,0.512000,0.409600,0.327680,0.262144];
337 round(kPolyVecLo*32768)
338 round(kPolyVecHi*32768)
339 */
340 static const int16_t kPolyVecLo[12] = {
341 29491, 26542, 23888, 21499, 19349, 17414, 15673, 14106, 12695, 11425, 10283, 9255
342 };
343 static const int16_t kPolyVecHi[6] = {
344 26214, 20972, 16777, 13422, 10737, 8590
345 };
346
log2_Q8_LPC(uint32_t x)347 static __inline int32_t log2_Q8_LPC( uint32_t x ) {
348
349 int32_t zeros, lg2;
350 int16_t frac;
351
352 zeros=WebRtcSpl_NormU32(x);
353 frac=(int16_t)WEBRTC_SPL_RSHIFT_W32(((uint32_t)WEBRTC_SPL_LSHIFT_W32(x, zeros)&0x7FFFFFFF), 23);
354
355 /* log2(x) */
356
357 lg2= (WEBRTC_SPL_LSHIFT_W16((31-zeros), 8)+frac);
358 return lg2;
359
360 }
361
362 static const int16_t kMulPitchGain = -25; /* 200/256 in Q5 */
363 static const int16_t kChngFactor = 3523; /* log10(2)*10/4*0.4/1.4=log10(2)/1.4= 0.2150 in Q14 */
364 static const int16_t kExp2 = 11819; /* 1/log(2) */
365 const int kShiftLowerBand = 11; /* Shift value for lower band in Q domain. */
366 const int kShiftHigherBand = 12; /* Shift value for higher band in Q domain. */
367
WebRtcIsacfix_GetVars(const int16_t * input,const int16_t * pitchGains_Q12,uint32_t * oldEnergy,int16_t * varscale)368 void WebRtcIsacfix_GetVars(const int16_t *input, const int16_t *pitchGains_Q12,
369 uint32_t *oldEnergy, int16_t *varscale)
370 {
371 int k;
372 uint32_t nrgQ[4];
373 int16_t nrgQlog[4];
374 int16_t tmp16, chng1, chng2, chng3, chng4, tmp, chngQ, oldNrgQlog, pgQ, pg3;
375 int32_t expPg32;
376 int16_t expPg, divVal;
377 int16_t tmp16_1, tmp16_2;
378
379 /* Calculate energies of first and second frame halfs */
380 nrgQ[0]=0;
381 for (k = QLOOKAHEAD/2; k < (FRAMESAMPLES/4 + QLOOKAHEAD) / 2; k++) {
382 nrgQ[0] +=WEBRTC_SPL_MUL_16_16(input[k],input[k]);
383 }
384 nrgQ[1]=0;
385 for ( ; k < (FRAMESAMPLES/2 + QLOOKAHEAD) / 2; k++) {
386 nrgQ[1] +=WEBRTC_SPL_MUL_16_16(input[k],input[k]);
387 }
388 nrgQ[2]=0;
389 for ( ; k < (WEBRTC_SPL_MUL_16_16(FRAMESAMPLES, 3)/4 + QLOOKAHEAD) / 2; k++) {
390 nrgQ[2] +=WEBRTC_SPL_MUL_16_16(input[k],input[k]);
391 }
392 nrgQ[3]=0;
393 for ( ; k < (FRAMESAMPLES + QLOOKAHEAD) / 2; k++) {
394 nrgQ[3] +=WEBRTC_SPL_MUL_16_16(input[k],input[k]);
395 }
396
397 for ( k=0; k<4; k++) {
398 nrgQlog[k] = (int16_t)log2_Q8_LPC(nrgQ[k]); /* log2(nrgQ) */
399 }
400 oldNrgQlog = (int16_t)log2_Q8_LPC(*oldEnergy);
401
402 /* Calculate average level change */
403 chng1 = WEBRTC_SPL_ABS_W16(nrgQlog[3]-nrgQlog[2]);
404 chng2 = WEBRTC_SPL_ABS_W16(nrgQlog[2]-nrgQlog[1]);
405 chng3 = WEBRTC_SPL_ABS_W16(nrgQlog[1]-nrgQlog[0]);
406 chng4 = WEBRTC_SPL_ABS_W16(nrgQlog[0]-oldNrgQlog);
407 tmp = chng1+chng2+chng3+chng4;
408 chngQ = (int16_t)WEBRTC_SPL_MUL_16_16_RSFT(tmp, kChngFactor, 10); /* Q12 */
409 chngQ += 2926; /* + 1.0/1.4 in Q12 */
410
411 /* Find average pitch gain */
412 pgQ = 0;
413 for (k=0; k<4; k++)
414 {
415 pgQ += pitchGains_Q12[k];
416 }
417
418 pg3 = (int16_t)WEBRTC_SPL_MUL_16_16_RSFT(pgQ, pgQ,11); /* pgQ in Q(12+2)=Q14. Q14*Q14>>11 => Q17 */
419 pg3 = (int16_t)WEBRTC_SPL_MUL_16_16_RSFT(pgQ, pg3,13); /* Q17*Q14>>13 =>Q18 */
420 pg3 = (int16_t)WEBRTC_SPL_MUL_16_16_RSFT(pg3, kMulPitchGain ,5); /* Q10 kMulPitchGain = -25 = -200 in Q-3. */
421
422 tmp16=(int16_t)WEBRTC_SPL_MUL_16_16_RSFT_WITH_ROUND(kExp2,pg3,13);/* Q13*Q10>>13 => Q10*/
423 if (tmp16<0) {
424 tmp16_2 = (0x0400 | (tmp16 & 0x03FF));
425 tmp16_1 = (WEBRTC_SPL_RSHIFT_W16((uint16_t)(tmp16 ^ 0xFFFF), 10)-3); /* Gives result in Q14 */
426 if (tmp16_1<0)
427 expPg=(int16_t) -WEBRTC_SPL_LSHIFT_W16(tmp16_2, -tmp16_1);
428 else
429 expPg=(int16_t) -WEBRTC_SPL_RSHIFT_W16(tmp16_2, tmp16_1);
430 } else
431 expPg = (int16_t) -16384; /* 1 in Q14, since 2^0=1 */
432
433 expPg32 = (int32_t)WEBRTC_SPL_LSHIFT_W16((int32_t)expPg, 8); /* Q22 */
434 divVal = WebRtcSpl_DivW32W16ResW16(expPg32, chngQ); /* Q22/Q12=Q10 */
435
436 tmp16=(int16_t)WEBRTC_SPL_MUL_16_16_RSFT_WITH_ROUND(kExp2,divVal,13);/* Q13*Q10>>13 => Q10*/
437 if (tmp16<0) {
438 tmp16_2 = (0x0400 | (tmp16 & 0x03FF));
439 tmp16_1 = (WEBRTC_SPL_RSHIFT_W16((uint16_t)(tmp16 ^ 0xFFFF), 10)-3); /* Gives result in Q14 */
440 if (tmp16_1<0)
441 expPg=(int16_t) WEBRTC_SPL_LSHIFT_W16(tmp16_2, -tmp16_1);
442 else
443 expPg=(int16_t) WEBRTC_SPL_RSHIFT_W16(tmp16_2, tmp16_1);
444 } else
445 expPg = (int16_t) 16384; /* 1 in Q14, since 2^0=1 */
446
447 *varscale = expPg-1;
448 *oldEnergy = nrgQ[3];
449 }
450
451
452
exp2_Q10_T(int16_t x)453 static __inline int16_t exp2_Q10_T(int16_t x) { // Both in and out in Q10
454
455 int16_t tmp16_1, tmp16_2;
456
457 tmp16_2=(int16_t)(0x0400|(x&0x03FF));
458 tmp16_1=-(int16_t)WEBRTC_SPL_RSHIFT_W16(x,10);
459 if(tmp16_1>0)
460 return (int16_t) WEBRTC_SPL_RSHIFT_W16(tmp16_2, tmp16_1);
461 else
462 return (int16_t) WEBRTC_SPL_LSHIFT_W16(tmp16_2, -tmp16_1);
463
464 }
465
466
467 // Declare function pointers.
468 AutocorrFix WebRtcIsacfix_AutocorrFix;
469 CalculateResidualEnergy WebRtcIsacfix_CalculateResidualEnergy;
470
471 /* This routine calculates the residual energy for LPC.
472 * Formula as shown in comments inside.
473 */
WebRtcIsacfix_CalculateResidualEnergyC(int lpc_order,int32_t q_val_corr,int q_val_polynomial,int16_t * a_polynomial,int32_t * corr_coeffs,int * q_val_residual_energy)474 int32_t WebRtcIsacfix_CalculateResidualEnergyC(int lpc_order,
475 int32_t q_val_corr,
476 int q_val_polynomial,
477 int16_t* a_polynomial,
478 int32_t* corr_coeffs,
479 int* q_val_residual_energy) {
480 int i = 0, j = 0;
481 int shift_internal = 0, shift_norm = 0;
482 int32_t tmp32 = 0, word32_high = 0, word32_low = 0, residual_energy = 0;
483 int64_t sum64 = 0, sum64_tmp = 0;
484
485 for (i = 0; i <= lpc_order; i++) {
486 for (j = i; j <= lpc_order; j++) {
487 /* For the case of i == 0: residual_energy +=
488 * a_polynomial[j] * corr_coeffs[i] * a_polynomial[j - i];
489 * For the case of i != 0: residual_energy +=
490 * a_polynomial[j] * corr_coeffs[i] * a_polynomial[j - i] * 2;
491 */
492
493 tmp32 = WEBRTC_SPL_MUL_16_16(a_polynomial[j], a_polynomial[j - i]);
494 /* tmp32 in Q(q_val_polynomial * 2). */
495 if (i != 0) {
496 tmp32 <<= 1;
497 }
498 sum64_tmp = (int64_t)tmp32 * (int64_t)corr_coeffs[i];
499 sum64_tmp >>= shift_internal;
500
501 /* Test overflow and sum the result. */
502 if(((sum64_tmp > 0 && sum64 > 0) && (LLONG_MAX - sum64 < sum64_tmp)) ||
503 ((sum64_tmp < 0 && sum64 < 0) && (LLONG_MIN - sum64 > sum64_tmp))) {
504 /* Shift right for overflow. */
505 shift_internal += 1;
506 sum64 >>= 1;
507 sum64 += sum64_tmp >> 1;
508 } else {
509 sum64 += sum64_tmp;
510 }
511 }
512 }
513
514 word32_high = (int32_t)(sum64 >> 32);
515 word32_low = (int32_t)sum64;
516
517 // Calculate the value of shifting (shift_norm) for the 64-bit sum.
518 if(word32_high != 0) {
519 shift_norm = 32 - WebRtcSpl_NormW32(word32_high);
520 residual_energy = (int32_t)(sum64 >> shift_norm);
521 } else {
522 if((word32_low & 0x80000000) != 0) {
523 shift_norm = 1;
524 residual_energy = (uint32_t)word32_low >> 1;
525 } else {
526 shift_norm = WebRtcSpl_NormW32(word32_low);
527 residual_energy = word32_low << shift_norm;
528 shift_norm = -shift_norm;
529 }
530 }
531
532 /* Q(q_val_polynomial * 2) * Q(q_val_corr) >> shift_internal >> shift_norm
533 * = Q(q_val_corr - shift_internal - shift_norm + q_val_polynomial * 2)
534 */
535 *q_val_residual_energy = q_val_corr - shift_internal - shift_norm
536 + q_val_polynomial * 2;
537
538 return residual_energy;
539 }
540
WebRtcIsacfix_GetLpcCoef(int16_t * inLoQ0,int16_t * inHiQ0,MaskFiltstr_enc * maskdata,int16_t snrQ10,const int16_t * pitchGains_Q12,int32_t * gain_lo_hiQ17,int16_t * lo_coeffQ15,int16_t * hi_coeffQ15)541 void WebRtcIsacfix_GetLpcCoef(int16_t *inLoQ0,
542 int16_t *inHiQ0,
543 MaskFiltstr_enc *maskdata,
544 int16_t snrQ10,
545 const int16_t *pitchGains_Q12,
546 int32_t *gain_lo_hiQ17,
547 int16_t *lo_coeffQ15,
548 int16_t *hi_coeffQ15)
549 {
550 int k, n, ii;
551 int pos1, pos2;
552 int sh_lo, sh_hi, sh, ssh, shMem;
553 int16_t varscaleQ14;
554
555 int16_t tmpQQlo, tmpQQhi;
556 int32_t tmp32;
557 int16_t tmp16,tmp16b;
558
559 int16_t polyHI[ORDERHI+1];
560 int16_t rcQ15_lo[ORDERLO], rcQ15_hi[ORDERHI];
561
562
563 int16_t DataLoQ6[WINLEN], DataHiQ6[WINLEN];
564 int32_t corrloQQ[ORDERLO+2];
565 int32_t corrhiQQ[ORDERHI+1];
566 int32_t corrlo2QQ[ORDERLO+1];
567 int16_t scale;
568 int16_t QdomLO, QdomHI, newQdomHI, newQdomLO;
569
570 int32_t res_nrgQQ;
571 int32_t sqrt_nrg;
572
573 /* less-noise-at-low-frequencies factor */
574 int16_t aaQ14;
575
576 /* Multiplication with 1/sqrt(12) ~= 0.28901734104046 can be done by convertion to
577 Q15, i.e. round(0.28901734104046*32768) = 9471, and use 9471/32768.0 ~= 0.289032
578 */
579 int16_t snrq;
580 int shft;
581
582 int16_t tmp16a;
583 int32_t tmp32a, tmp32b, tmp32c;
584
585 int16_t a_LOQ11[ORDERLO+1];
586 int16_t k_vecloQ15[ORDERLO];
587 int16_t a_HIQ12[ORDERHI+1];
588 int16_t k_vechiQ15[ORDERHI];
589
590 int16_t stab;
591
592 snrq=snrQ10;
593
594 /* SNR= C * 2 ^ (D * snrq) ; C=0.289, D=0.05*log2(10)=0.166 (~=172 in Q10)*/
595 tmp16 = (int16_t) WEBRTC_SPL_MUL_16_16_RSFT(snrq, 172, 10); // Q10
596 tmp16b = exp2_Q10_T(tmp16); // Q10
597 snrq = (int16_t) WEBRTC_SPL_MUL_16_16_RSFT(tmp16b, 285, 10); // Q10
598
599 /* change quallevel depending on pitch gains and level fluctuations */
600 WebRtcIsacfix_GetVars(inLoQ0, pitchGains_Q12, &(maskdata->OldEnergy), &varscaleQ14);
601
602 /* less-noise-at-low-frequencies factor */
603 /* Calculation of 0.35 * (0.5 + 0.5 * varscale) in fixpoint:
604 With 0.35 in Q16 (0.35 ~= 22938/65536.0 = 0.3500061) and varscaleQ14 in Q14,
605 we get Q16*Q14>>16 = Q14
606 */
607 aaQ14 = (int16_t) WEBRTC_SPL_RSHIFT_W32(
608 (WEBRTC_SPL_MUL_16_16(22938, (8192 + WEBRTC_SPL_RSHIFT_W32(varscaleQ14, 1)))
609 + ((int32_t)32768)), 16);
610
611 /* Calculate tmp = (1.0 + aa*aa); in Q12 */
612 tmp16 = (int16_t) WEBRTC_SPL_MUL_16_16_RSFT(aaQ14, aaQ14, 15); //Q14*Q14>>15 = Q13
613 tmpQQlo = 4096 + WEBRTC_SPL_RSHIFT_W16(tmp16, 1); // Q12 + Q13>>1 = Q12
614
615 /* Calculate tmp = (1.0+aa) * (1.0+aa); */
616 tmp16 = 8192 + WEBRTC_SPL_RSHIFT_W16(aaQ14, 1); // 1+a in Q13
617 tmpQQhi = (int16_t) WEBRTC_SPL_MUL_16_16_RSFT(tmp16, tmp16, 14); //Q13*Q13>>14 = Q12
618
619 /* replace data in buffer by new look-ahead data */
620 for (pos1 = 0; pos1 < QLOOKAHEAD; pos1++) {
621 maskdata->DataBufferLoQ0[pos1 + WINLEN - QLOOKAHEAD] = inLoQ0[pos1];
622 }
623
624 for (k = 0; k < SUBFRAMES; k++) {
625
626 /* Update input buffer and multiply signal with window */
627 for (pos1 = 0; pos1 < WINLEN - UPDATE/2; pos1++) {
628 maskdata->DataBufferLoQ0[pos1] = maskdata->DataBufferLoQ0[pos1 + UPDATE/2];
629 maskdata->DataBufferHiQ0[pos1] = maskdata->DataBufferHiQ0[pos1 + UPDATE/2];
630 DataLoQ6[pos1] = (int16_t) WEBRTC_SPL_MUL_16_16_RSFT(
631 maskdata->DataBufferLoQ0[pos1], kWindowAutocorr[pos1], 15); // Q0*Q21>>15 = Q6
632 DataHiQ6[pos1] = (int16_t) WEBRTC_SPL_MUL_16_16_RSFT(
633 maskdata->DataBufferHiQ0[pos1], kWindowAutocorr[pos1], 15); // Q0*Q21>>15 = Q6
634 }
635 pos2 = (int16_t)(WEBRTC_SPL_MUL_16_16(k, UPDATE)/2);
636 for (n = 0; n < UPDATE/2; n++, pos1++) {
637 maskdata->DataBufferLoQ0[pos1] = inLoQ0[QLOOKAHEAD + pos2];
638 maskdata->DataBufferHiQ0[pos1] = inHiQ0[pos2++];
639 DataLoQ6[pos1] = (int16_t) WEBRTC_SPL_MUL_16_16_RSFT(
640 maskdata->DataBufferLoQ0[pos1], kWindowAutocorr[pos1], 15); // Q0*Q21>>15 = Q6
641 DataHiQ6[pos1] = (int16_t) WEBRTC_SPL_MUL_16_16_RSFT(
642 maskdata->DataBufferHiQ0[pos1], kWindowAutocorr[pos1], 15); // Q0*Q21>>15 = Q6
643 }
644
645 /* Get correlation coefficients */
646 /* The highest absolute value measured inside DataLo in the test set
647 For DataHi, corresponding value was 160.
648
649 This means that it should be possible to represent the input values
650 to WebRtcSpl_AutoCorrelation() as Q6 values (since 307*2^6 =
651 19648). Of course, Q0 will also work, but due to the low energy in
652 DataLo and DataHi, the outputted autocorrelation will be more accurate
653 and mimic the floating point code better, by being in an high as possible
654 Q-domain.
655 */
656
657 WebRtcIsacfix_AutocorrFix(corrloQQ,DataLoQ6,WINLEN, ORDERLO+1, &scale);
658 QdomLO = 12-scale; // QdomLO is the Q-domain of corrloQQ
659 sh_lo = WebRtcSpl_NormW32(corrloQQ[0]);
660 QdomLO += sh_lo;
661 for (ii=0; ii<ORDERLO+2; ii++) {
662 corrloQQ[ii] = WEBRTC_SPL_LSHIFT_W32(corrloQQ[ii], sh_lo);
663 }
664 /* It is investigated whether it was possible to use 16 bits for the
665 32-bit vector corrloQQ, but it didn't work. */
666
667 WebRtcIsacfix_AutocorrFix(corrhiQQ,DataHiQ6,WINLEN, ORDERHI, &scale);
668
669 QdomHI = 12-scale; // QdomHI is the Q-domain of corrhiQQ
670 sh_hi = WebRtcSpl_NormW32(corrhiQQ[0]);
671 QdomHI += sh_hi;
672 for (ii=0; ii<ORDERHI+1; ii++) {
673 corrhiQQ[ii] = WEBRTC_SPL_LSHIFT_W32(corrhiQQ[ii], sh_hi);
674 }
675
676 /* less noise for lower frequencies, by filtering/scaling autocorrelation sequences */
677
678 /* Calculate corrlo2[0] = tmpQQlo * corrlo[0] - 2.0*tmpQQlo * corrlo[1];*/
679 corrlo2QQ[0] = WEBRTC_SPL_RSHIFT_W32(WEBRTC_SPL_MUL_16_32_RSFT16(tmpQQlo, corrloQQ[0]), 1)- // Q(12+QdomLO-16)>>1 = Q(QdomLO-5)
680 WEBRTC_SPL_RSHIFT_W32(WEBRTC_SPL_MUL_16_32_RSFT16(aaQ14, corrloQQ[1]), 2); // 2*Q(14+QdomLO-16)>>3 = Q(QdomLO-2)>>2 = Q(QdomLO-5)
681
682 /* Calculate corrlo2[n] = tmpQQlo * corrlo[n] - tmpQQlo * (corrlo[n-1] + corrlo[n+1]);*/
683 for (n = 1; n <= ORDERLO; n++) {
684
685 tmp32 = WEBRTC_SPL_RSHIFT_W32(corrloQQ[n-1], 1) + WEBRTC_SPL_RSHIFT_W32(corrloQQ[n+1], 1); // Q(QdomLO-1)
686 corrlo2QQ[n] = WEBRTC_SPL_RSHIFT_W32(WEBRTC_SPL_MUL_16_32_RSFT16(tmpQQlo, corrloQQ[n]), 1)- // Q(12+QdomLO-16)>>1 = Q(QdomLO-5)
687 WEBRTC_SPL_RSHIFT_W32(WEBRTC_SPL_MUL_16_32_RSFT16(aaQ14, tmp32), 2); // Q(14+QdomLO-1-16)>>2 = Q(QdomLO-3)>>2 = Q(QdomLO-5)
688
689 }
690 QdomLO -= 5;
691
692 /* Calculate corrhi[n] = tmpQQhi * corrhi[n]; */
693 for (n = 0; n <= ORDERHI; n++) {
694 corrhiQQ[n] = WEBRTC_SPL_MUL_16_32_RSFT16(tmpQQhi, corrhiQQ[n]); // Q(12+QdomHI-16) = Q(QdomHI-4)
695 }
696 QdomHI -= 4;
697
698 /* add white noise floor */
699 /* corrlo2QQ is in Q(QdomLO) and corrhiQQ is in Q(QdomHI) */
700 /* Calculate corrlo2[0] += 9.5367431640625e-7; and
701 corrhi[0] += 9.5367431640625e-7, where the constant is 1/2^20 */
702
703 tmp32 = WEBRTC_SPL_SHIFT_W32((int32_t) 1, QdomLO-20);
704 corrlo2QQ[0] += tmp32;
705 tmp32 = WEBRTC_SPL_SHIFT_W32((int32_t) 1, QdomHI-20);
706 corrhiQQ[0] += tmp32;
707
708 /* corrlo2QQ is in Q(QdomLO) and corrhiQQ is in Q(QdomHI) before the following
709 code segment, where we want to make sure we get a 1-bit margin */
710 for (n = 0; n <= ORDERLO; n++) {
711 corrlo2QQ[n] = WEBRTC_SPL_RSHIFT_W32(corrlo2QQ[n], 1); // Make sure we have a 1-bit margin
712 }
713 QdomLO -= 1; // Now, corrlo2QQ is in Q(QdomLO), with a 1-bit margin
714
715 for (n = 0; n <= ORDERHI; n++) {
716 corrhiQQ[n] = WEBRTC_SPL_RSHIFT_W32(corrhiQQ[n], 1); // Make sure we have a 1-bit margin
717 }
718 QdomHI -= 1; // Now, corrhiQQ is in Q(QdomHI), with a 1-bit margin
719
720
721 newQdomLO = QdomLO;
722
723 for (n = 0; n <= ORDERLO; n++) {
724 int32_t tmp, tmpB, tmpCorr;
725 int16_t alpha=328; //0.01 in Q15
726 int16_t beta=324; //(1-0.01)*0.01=0.0099 in Q15
727 int16_t gamma=32440; //(1-0.01)=0.99 in Q15
728
729 if (maskdata->CorrBufLoQQ[n] != 0) {
730 shMem=WebRtcSpl_NormW32(maskdata->CorrBufLoQQ[n]);
731 sh = QdomLO - maskdata->CorrBufLoQdom[n];
732 if (sh<=shMem) {
733 tmp = WEBRTC_SPL_SHIFT_W32(maskdata->CorrBufLoQQ[n], sh); // Get CorrBufLoQQ to same domain as corrlo2
734 tmp = WEBRTC_SPL_MUL_16_32_RSFT15(alpha, tmp);
735 } else if ((sh-shMem)<7){
736 tmp = WEBRTC_SPL_SHIFT_W32(maskdata->CorrBufLoQQ[n], shMem); // Shift up CorrBufLoQQ as much as possible
737 tmp = WEBRTC_SPL_MUL_16_32_RSFT15(WEBRTC_SPL_LSHIFT_W16(alpha, (sh-shMem)), tmp); // Shift alpha the number of times required to get tmp in QdomLO
738 } else {
739 tmp = WEBRTC_SPL_SHIFT_W32(maskdata->CorrBufLoQQ[n], shMem); // Shift up CorrBufHiQQ as much as possible
740 tmp = WEBRTC_SPL_MUL_16_32_RSFT15(WEBRTC_SPL_LSHIFT_W16(alpha, 6), tmp); // Shift alpha as much as possible without overflow the number of times required to get tmp in QdomHI
741 tmpCorr = WEBRTC_SPL_RSHIFT_W32(corrloQQ[n], sh-shMem-6);
742 tmp = tmp + tmpCorr;
743 maskdata->CorrBufLoQQ[n] = tmp;
744 newQdomLO = QdomLO-(sh-shMem-6);
745 maskdata->CorrBufLoQdom[n] = newQdomLO;
746 }
747 } else
748 tmp = 0;
749
750 tmp = tmp + corrlo2QQ[n];
751
752 maskdata->CorrBufLoQQ[n] = tmp;
753 maskdata->CorrBufLoQdom[n] = QdomLO;
754
755 tmp=WEBRTC_SPL_MUL_16_32_RSFT15(beta, tmp);
756 tmpB=WEBRTC_SPL_MUL_16_32_RSFT15(gamma, corrlo2QQ[n]);
757 corrlo2QQ[n] = tmp + tmpB;
758 }
759 if( newQdomLO!=QdomLO) {
760 for (n = 0; n <= ORDERLO; n++) {
761 if (maskdata->CorrBufLoQdom[n] != newQdomLO)
762 corrloQQ[n] = WEBRTC_SPL_RSHIFT_W32(corrloQQ[n], maskdata->CorrBufLoQdom[n]-newQdomLO);
763 }
764 QdomLO = newQdomLO;
765 }
766
767
768 newQdomHI = QdomHI;
769
770 for (n = 0; n <= ORDERHI; n++) {
771 int32_t tmp, tmpB, tmpCorr;
772 int16_t alpha=328; //0.01 in Q15
773 int16_t beta=324; //(1-0.01)*0.01=0.0099 in Q15
774 int16_t gamma=32440; //(1-0.01)=0.99 in Q1
775 if (maskdata->CorrBufHiQQ[n] != 0) {
776 shMem=WebRtcSpl_NormW32(maskdata->CorrBufHiQQ[n]);
777 sh = QdomHI - maskdata->CorrBufHiQdom[n];
778 if (sh<=shMem) {
779 tmp = WEBRTC_SPL_SHIFT_W32(maskdata->CorrBufHiQQ[n], sh); // Get CorrBufHiQQ to same domain as corrhi
780 tmp = WEBRTC_SPL_MUL_16_32_RSFT15(alpha, tmp);
781 tmpCorr = corrhiQQ[n];
782 tmp = tmp + tmpCorr;
783 maskdata->CorrBufHiQQ[n] = tmp;
784 maskdata->CorrBufHiQdom[n] = QdomHI;
785 } else if ((sh-shMem)<7) {
786 tmp = WEBRTC_SPL_SHIFT_W32(maskdata->CorrBufHiQQ[n], shMem); // Shift up CorrBufHiQQ as much as possible
787 tmp = WEBRTC_SPL_MUL_16_32_RSFT15(WEBRTC_SPL_LSHIFT_W16(alpha, (sh-shMem)), tmp); // Shift alpha the number of times required to get tmp in QdomHI
788 tmpCorr = corrhiQQ[n];
789 tmp = tmp + tmpCorr;
790 maskdata->CorrBufHiQQ[n] = tmp;
791 maskdata->CorrBufHiQdom[n] = QdomHI;
792 } else {
793 tmp = WEBRTC_SPL_SHIFT_W32(maskdata->CorrBufHiQQ[n], shMem); // Shift up CorrBufHiQQ as much as possible
794 tmp = WEBRTC_SPL_MUL_16_32_RSFT15(WEBRTC_SPL_LSHIFT_W16(alpha, 6), tmp); // Shift alpha as much as possible without overflow the number of times required to get tmp in QdomHI
795 tmpCorr = WEBRTC_SPL_RSHIFT_W32(corrhiQQ[n], sh-shMem-6);
796 tmp = tmp + tmpCorr;
797 maskdata->CorrBufHiQQ[n] = tmp;
798 newQdomHI = QdomHI-(sh-shMem-6);
799 maskdata->CorrBufHiQdom[n] = newQdomHI;
800 }
801 } else {
802 tmp = corrhiQQ[n];
803 tmpCorr = tmp;
804 maskdata->CorrBufHiQQ[n] = tmp;
805 maskdata->CorrBufHiQdom[n] = QdomHI;
806 }
807
808 tmp=WEBRTC_SPL_MUL_16_32_RSFT15(beta, tmp);
809 tmpB=WEBRTC_SPL_MUL_16_32_RSFT15(gamma, tmpCorr);
810 corrhiQQ[n] = tmp + tmpB;
811 }
812
813 if( newQdomHI!=QdomHI) {
814 for (n = 0; n <= ORDERHI; n++) {
815 if (maskdata->CorrBufHiQdom[n] != newQdomHI)
816 corrhiQQ[n] = WEBRTC_SPL_RSHIFT_W32(corrhiQQ[n], maskdata->CorrBufHiQdom[n]-newQdomHI);
817 }
818 QdomHI = newQdomHI;
819 }
820
821 stab=WebRtcSpl_LevinsonW32_JSK(corrlo2QQ, a_LOQ11, k_vecloQ15, ORDERLO);
822
823 if (stab<0) { // If unstable use lower order
824 a_LOQ11[0]=2048;
825 for (n = 1; n <= ORDERLO; n++) {
826 a_LOQ11[n]=0;
827 }
828
829 stab=WebRtcSpl_LevinsonW32_JSK(corrlo2QQ, a_LOQ11, k_vecloQ15, 8);
830 }
831
832
833 WebRtcSpl_LevinsonDurbin(corrhiQQ, a_HIQ12, k_vechiQ15, ORDERHI);
834
835 /* bandwidth expansion */
836 for (n = 1; n <= ORDERLO; n++) {
837 a_LOQ11[n] = (int16_t) ((WEBRTC_SPL_MUL_16_16(
838 kPolyVecLo[n-1], a_LOQ11[n]) + ((int32_t) (1 << 14))) >> 15);
839 }
840
841
842 polyHI[0] = a_HIQ12[0];
843 for (n = 1; n <= ORDERHI; n++) {
844 a_HIQ12[n] = (int16_t) ((WEBRTC_SPL_MUL_16_16(
845 kPolyVecHi[n-1], a_HIQ12[n]) + ((int32_t) (1 << 14))) >> 15);
846 polyHI[n] = a_HIQ12[n];
847 }
848
849 /* Normalize the corrlo2 vector */
850 sh = WebRtcSpl_NormW32(corrlo2QQ[0]);
851 for (n = 0; n <= ORDERLO; n++) {
852 corrlo2QQ[n] = WEBRTC_SPL_LSHIFT_W32(corrlo2QQ[n], sh);
853 }
854 QdomLO += sh; /* Now, corrlo2QQ is still in Q(QdomLO) */
855
856
857 /* residual energy */
858
859 sh_lo = 31;
860 res_nrgQQ = WebRtcIsacfix_CalculateResidualEnergy(ORDERLO, QdomLO,
861 kShiftLowerBand, a_LOQ11, corrlo2QQ, &sh_lo);
862
863 /* Convert to reflection coefficients */
864 WebRtcSpl_AToK_JSK(a_LOQ11, ORDERLO, rcQ15_lo);
865
866 if (sh_lo & 0x0001) {
867 res_nrgQQ=WEBRTC_SPL_RSHIFT_W32(res_nrgQQ, 1);
868 sh_lo-=1;
869 }
870
871
872 if( res_nrgQQ > 0 )
873 {
874 sqrt_nrg=WebRtcSpl_Sqrt(res_nrgQQ);
875
876 /* add hearing threshold and compute the gain */
877 /* lo_coeff = varscale * S_N_R / (sqrt_nrg + varscale * H_T_H); */
878
879
880 //tmp32a=WEBRTC_SPL_MUL_16_16_RSFT(varscaleQ14, H_T_HQ19, 17); // Q14
881 tmp32a=WEBRTC_SPL_RSHIFT_W32((int32_t) varscaleQ14,1); // H_T_HQ19=65536 (16-17=-1) ssh= WEBRTC_SPL_RSHIFT_W16(sh_lo, 1); // sqrt_nrg is in Qssh
882 ssh= WEBRTC_SPL_RSHIFT_W16(sh_lo, 1); // sqrt_nrg is in Qssh
883 sh = ssh - 14;
884 tmp32b = WEBRTC_SPL_SHIFT_W32(tmp32a, sh); // Q14->Qssh
885 tmp32c = sqrt_nrg + tmp32b; // Qssh (denominator)
886 tmp32a = WEBRTC_SPL_MUL_16_16_RSFT(varscaleQ14, snrq, 0); //Q24 (numerator)
887
888 sh = WebRtcSpl_NormW32(tmp32c);
889 shft = 16 - sh;
890 tmp16a = (int16_t) WEBRTC_SPL_SHIFT_W32(tmp32c, -shft); // Q(ssh-shft) (denominator)
891
892 tmp32b = WebRtcSpl_DivW32W16(tmp32a, tmp16a); // Q(24-ssh+shft)
893 sh = ssh-shft-7;
894 *gain_lo_hiQ17 = WEBRTC_SPL_SHIFT_W32(tmp32b, sh); // Gains in Q17
895 }
896 else
897 {
898 *gain_lo_hiQ17 = 100; //(int32_t)WEBRTC_SPL_LSHIFT_W32( (int32_t)1, 17); // Gains in Q17
899 }
900 gain_lo_hiQ17++;
901
902 /* copy coefficients to output array */
903 for (n = 0; n < ORDERLO; n++) {
904 *lo_coeffQ15 = (int16_t) (rcQ15_lo[n]);
905 lo_coeffQ15++;
906 }
907 /* residual energy */
908 sh_hi = 31;
909 res_nrgQQ = WebRtcIsacfix_CalculateResidualEnergy(ORDERHI, QdomHI,
910 kShiftHigherBand, a_HIQ12, corrhiQQ, &sh_hi);
911
912 /* Convert to reflection coefficients */
913 WebRtcSpl_LpcToReflCoef(polyHI, ORDERHI, rcQ15_hi);
914
915 if (sh_hi & 0x0001) {
916 res_nrgQQ=WEBRTC_SPL_RSHIFT_W32(res_nrgQQ, 1);
917 sh_hi-=1;
918 }
919
920
921 if( res_nrgQQ > 0 )
922 {
923 sqrt_nrg=WebRtcSpl_Sqrt(res_nrgQQ);
924
925
926 /* add hearing threshold and compute the gain */
927 /* hi_coeff = varscale * S_N_R / (sqrt_nrg + varscale * H_T_H); */
928
929 //tmp32a=WEBRTC_SPL_MUL_16_16_RSFT(varscaleQ14, H_T_HQ19, 17); // Q14
930 tmp32a=WEBRTC_SPL_RSHIFT_W32((int32_t) varscaleQ14,1); // H_T_HQ19=65536 (16-17=-1)
931
932 ssh= WEBRTC_SPL_RSHIFT_W32(sh_hi, 1); // sqrt_nrg is in Qssh
933 sh = ssh - 14;
934 tmp32b = WEBRTC_SPL_SHIFT_W32(tmp32a, sh); // Q14->Qssh
935 tmp32c = sqrt_nrg + tmp32b; // Qssh (denominator)
936 tmp32a = WEBRTC_SPL_MUL_16_16_RSFT(varscaleQ14, snrq, 0); //Q24 (numerator)
937
938 sh = WebRtcSpl_NormW32(tmp32c);
939 shft = 16 - sh;
940 tmp16a = (int16_t) WEBRTC_SPL_SHIFT_W32(tmp32c, -shft); // Q(ssh-shft) (denominator)
941
942 tmp32b = WebRtcSpl_DivW32W16(tmp32a, tmp16a); // Q(24-ssh+shft)
943 sh = ssh-shft-7;
944 *gain_lo_hiQ17 = WEBRTC_SPL_SHIFT_W32(tmp32b, sh); // Gains in Q17
945 }
946 else
947 {
948 *gain_lo_hiQ17 = 100; //(int32_t)WEBRTC_SPL_LSHIFT_W32( (int32_t)1, 17); // Gains in Q17
949 }
950 gain_lo_hiQ17++;
951
952
953 /* copy coefficients to output array */
954 for (n = 0; n < ORDERHI; n++) {
955 *hi_coeffQ15 = rcQ15_hi[n];
956 hi_coeffQ15++;
957 }
958 }
959 }
960