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 */ 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 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 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 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 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 */ 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 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
