media/libopus/celt/pitch.c

Thu, 15 Jan 2015 15:59:08 +0100

author
Michael Schloh von Bennewitz <michael@schloh.com>
date
Thu, 15 Jan 2015 15:59:08 +0100
branch
TOR_BUG_9701
changeset 10
ac0c01689b40
permissions
-rw-r--r--

Implement a real Private Browsing Mode condition by changing the API/ABI;
This solves Tor bug #9701, complying with disk avoidance documented in
https://www.torproject.org/projects/torbrowser/design/#disk-avoidance.

michael@0 1 /* Copyright (c) 2007-2008 CSIRO
michael@0 2 Copyright (c) 2007-2009 Xiph.Org Foundation
michael@0 3 Written by Jean-Marc Valin */
michael@0 4 /**
michael@0 5 @file pitch.c
michael@0 6 @brief Pitch analysis
michael@0 7 */
michael@0 8
michael@0 9 /*
michael@0 10 Redistribution and use in source and binary forms, with or without
michael@0 11 modification, are permitted provided that the following conditions
michael@0 12 are met:
michael@0 13
michael@0 14 - Redistributions of source code must retain the above copyright
michael@0 15 notice, this list of conditions and the following disclaimer.
michael@0 16
michael@0 17 - Redistributions in binary form must reproduce the above copyright
michael@0 18 notice, this list of conditions and the following disclaimer in the
michael@0 19 documentation and/or other materials provided with the distribution.
michael@0 20
michael@0 21 THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
michael@0 22 ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
michael@0 23 LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
michael@0 24 A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER
michael@0 25 OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
michael@0 26 EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
michael@0 27 PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
michael@0 28 PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
michael@0 29 LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
michael@0 30 NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
michael@0 31 SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
michael@0 32 */
michael@0 33
michael@0 34 #ifdef HAVE_CONFIG_H
michael@0 35 #include "config.h"
michael@0 36 #endif
michael@0 37
michael@0 38 #include "pitch.h"
michael@0 39 #include "os_support.h"
michael@0 40 #include "modes.h"
michael@0 41 #include "stack_alloc.h"
michael@0 42 #include "mathops.h"
michael@0 43 #include "celt_lpc.h"
michael@0 44
michael@0 45 static void find_best_pitch(opus_val32 *xcorr, opus_val16 *y, int len,
michael@0 46 int max_pitch, int *best_pitch
michael@0 47 #ifdef FIXED_POINT
michael@0 48 , int yshift, opus_val32 maxcorr
michael@0 49 #endif
michael@0 50 )
michael@0 51 {
michael@0 52 int i, j;
michael@0 53 opus_val32 Syy=1;
michael@0 54 opus_val16 best_num[2];
michael@0 55 opus_val32 best_den[2];
michael@0 56 #ifdef FIXED_POINT
michael@0 57 int xshift;
michael@0 58
michael@0 59 xshift = celt_ilog2(maxcorr)-14;
michael@0 60 #endif
michael@0 61
michael@0 62 best_num[0] = -1;
michael@0 63 best_num[1] = -1;
michael@0 64 best_den[0] = 0;
michael@0 65 best_den[1] = 0;
michael@0 66 best_pitch[0] = 0;
michael@0 67 best_pitch[1] = 1;
michael@0 68 for (j=0;j<len;j++)
michael@0 69 Syy = ADD32(Syy, SHR32(MULT16_16(y[j],y[j]), yshift));
michael@0 70 for (i=0;i<max_pitch;i++)
michael@0 71 {
michael@0 72 if (xcorr[i]>0)
michael@0 73 {
michael@0 74 opus_val16 num;
michael@0 75 opus_val32 xcorr16;
michael@0 76 xcorr16 = EXTRACT16(VSHR32(xcorr[i], xshift));
michael@0 77 #ifndef FIXED_POINT
michael@0 78 /* Considering the range of xcorr16, this should avoid both underflows
michael@0 79 and overflows (inf) when squaring xcorr16 */
michael@0 80 xcorr16 *= 1e-12f;
michael@0 81 #endif
michael@0 82 num = MULT16_16_Q15(xcorr16,xcorr16);
michael@0 83 if (MULT16_32_Q15(num,best_den[1]) > MULT16_32_Q15(best_num[1],Syy))
michael@0 84 {
michael@0 85 if (MULT16_32_Q15(num,best_den[0]) > MULT16_32_Q15(best_num[0],Syy))
michael@0 86 {
michael@0 87 best_num[1] = best_num[0];
michael@0 88 best_den[1] = best_den[0];
michael@0 89 best_pitch[1] = best_pitch[0];
michael@0 90 best_num[0] = num;
michael@0 91 best_den[0] = Syy;
michael@0 92 best_pitch[0] = i;
michael@0 93 } else {
michael@0 94 best_num[1] = num;
michael@0 95 best_den[1] = Syy;
michael@0 96 best_pitch[1] = i;
michael@0 97 }
michael@0 98 }
michael@0 99 }
michael@0 100 Syy += SHR32(MULT16_16(y[i+len],y[i+len]),yshift) - SHR32(MULT16_16(y[i],y[i]),yshift);
michael@0 101 Syy = MAX32(1, Syy);
michael@0 102 }
michael@0 103 }
michael@0 104
michael@0 105 static void celt_fir5(const opus_val16 *x,
michael@0 106 const opus_val16 *num,
michael@0 107 opus_val16 *y,
michael@0 108 int N,
michael@0 109 opus_val16 *mem)
michael@0 110 {
michael@0 111 int i;
michael@0 112 opus_val16 num0, num1, num2, num3, num4;
michael@0 113 opus_val32 mem0, mem1, mem2, mem3, mem4;
michael@0 114 num0=num[0];
michael@0 115 num1=num[1];
michael@0 116 num2=num[2];
michael@0 117 num3=num[3];
michael@0 118 num4=num[4];
michael@0 119 mem0=mem[0];
michael@0 120 mem1=mem[1];
michael@0 121 mem2=mem[2];
michael@0 122 mem3=mem[3];
michael@0 123 mem4=mem[4];
michael@0 124 for (i=0;i<N;i++)
michael@0 125 {
michael@0 126 opus_val32 sum = SHL32(EXTEND32(x[i]), SIG_SHIFT);
michael@0 127 sum = MAC16_16(sum,num0,mem0);
michael@0 128 sum = MAC16_16(sum,num1,mem1);
michael@0 129 sum = MAC16_16(sum,num2,mem2);
michael@0 130 sum = MAC16_16(sum,num3,mem3);
michael@0 131 sum = MAC16_16(sum,num4,mem4);
michael@0 132 mem4 = mem3;
michael@0 133 mem3 = mem2;
michael@0 134 mem2 = mem1;
michael@0 135 mem1 = mem0;
michael@0 136 mem0 = x[i];
michael@0 137 y[i] = ROUND16(sum, SIG_SHIFT);
michael@0 138 }
michael@0 139 mem[0]=mem0;
michael@0 140 mem[1]=mem1;
michael@0 141 mem[2]=mem2;
michael@0 142 mem[3]=mem3;
michael@0 143 mem[4]=mem4;
michael@0 144 }
michael@0 145
michael@0 146
michael@0 147 void pitch_downsample(celt_sig * OPUS_RESTRICT x[], opus_val16 * OPUS_RESTRICT x_lp,
michael@0 148 int len, int C, int arch)
michael@0 149 {
michael@0 150 int i;
michael@0 151 opus_val32 ac[5];
michael@0 152 opus_val16 tmp=Q15ONE;
michael@0 153 opus_val16 lpc[4], mem[5]={0,0,0,0,0};
michael@0 154 opus_val16 lpc2[5];
michael@0 155 opus_val16 c1 = QCONST16(.8f,15);
michael@0 156 #ifdef FIXED_POINT
michael@0 157 int shift;
michael@0 158 opus_val32 maxabs = celt_maxabs32(x[0], len);
michael@0 159 if (C==2)
michael@0 160 {
michael@0 161 opus_val32 maxabs_1 = celt_maxabs32(x[1], len);
michael@0 162 maxabs = MAX32(maxabs, maxabs_1);
michael@0 163 }
michael@0 164 if (maxabs<1)
michael@0 165 maxabs=1;
michael@0 166 shift = celt_ilog2(maxabs)-10;
michael@0 167 if (shift<0)
michael@0 168 shift=0;
michael@0 169 if (C==2)
michael@0 170 shift++;
michael@0 171 #endif
michael@0 172 for (i=1;i<len>>1;i++)
michael@0 173 x_lp[i] = SHR32(HALF32(HALF32(x[0][(2*i-1)]+x[0][(2*i+1)])+x[0][2*i]), shift);
michael@0 174 x_lp[0] = SHR32(HALF32(HALF32(x[0][1])+x[0][0]), shift);
michael@0 175 if (C==2)
michael@0 176 {
michael@0 177 for (i=1;i<len>>1;i++)
michael@0 178 x_lp[i] += SHR32(HALF32(HALF32(x[1][(2*i-1)]+x[1][(2*i+1)])+x[1][2*i]), shift);
michael@0 179 x_lp[0] += SHR32(HALF32(HALF32(x[1][1])+x[1][0]), shift);
michael@0 180 }
michael@0 181
michael@0 182 _celt_autocorr(x_lp, ac, NULL, 0,
michael@0 183 4, len>>1, arch);
michael@0 184
michael@0 185 /* Noise floor -40 dB */
michael@0 186 #ifdef FIXED_POINT
michael@0 187 ac[0] += SHR32(ac[0],13);
michael@0 188 #else
michael@0 189 ac[0] *= 1.0001f;
michael@0 190 #endif
michael@0 191 /* Lag windowing */
michael@0 192 for (i=1;i<=4;i++)
michael@0 193 {
michael@0 194 /*ac[i] *= exp(-.5*(2*M_PI*.002*i)*(2*M_PI*.002*i));*/
michael@0 195 #ifdef FIXED_POINT
michael@0 196 ac[i] -= MULT16_32_Q15(2*i*i, ac[i]);
michael@0 197 #else
michael@0 198 ac[i] -= ac[i]*(.008f*i)*(.008f*i);
michael@0 199 #endif
michael@0 200 }
michael@0 201
michael@0 202 _celt_lpc(lpc, ac, 4);
michael@0 203 for (i=0;i<4;i++)
michael@0 204 {
michael@0 205 tmp = MULT16_16_Q15(QCONST16(.9f,15), tmp);
michael@0 206 lpc[i] = MULT16_16_Q15(lpc[i], tmp);
michael@0 207 }
michael@0 208 /* Add a zero */
michael@0 209 lpc2[0] = lpc[0] + QCONST16(.8f,SIG_SHIFT);
michael@0 210 lpc2[1] = lpc[1] + MULT16_16_Q15(c1,lpc[0]);
michael@0 211 lpc2[2] = lpc[2] + MULT16_16_Q15(c1,lpc[1]);
michael@0 212 lpc2[3] = lpc[3] + MULT16_16_Q15(c1,lpc[2]);
michael@0 213 lpc2[4] = MULT16_16_Q15(c1,lpc[3]);
michael@0 214 celt_fir5(x_lp, lpc2, x_lp, len>>1, mem);
michael@0 215 }
michael@0 216
michael@0 217 #if 0 /* This is a simple version of the pitch correlation that should work
michael@0 218 well on DSPs like Blackfin and TI C5x/C6x */
michael@0 219
michael@0 220 #ifdef FIXED_POINT
michael@0 221 opus_val32
michael@0 222 #else
michael@0 223 void
michael@0 224 #endif
michael@0 225 celt_pitch_xcorr(opus_val16 *x, opus_val16 *y, opus_val32 *xcorr, int len, int max_pitch)
michael@0 226 {
michael@0 227 int i, j;
michael@0 228 #ifdef FIXED_POINT
michael@0 229 opus_val32 maxcorr=1;
michael@0 230 #endif
michael@0 231 for (i=0;i<max_pitch;i++)
michael@0 232 {
michael@0 233 opus_val32 sum = 0;
michael@0 234 for (j=0;j<len;j++)
michael@0 235 sum = MAC16_16(sum, x[j],y[i+j]);
michael@0 236 xcorr[i] = sum;
michael@0 237 #ifdef FIXED_POINT
michael@0 238 maxcorr = MAX32(maxcorr, sum);
michael@0 239 #endif
michael@0 240 }
michael@0 241 #ifdef FIXED_POINT
michael@0 242 return maxcorr;
michael@0 243 #endif
michael@0 244 }
michael@0 245
michael@0 246 #else /* Unrolled version of the pitch correlation -- runs faster on x86 and ARM */
michael@0 247
michael@0 248 #ifdef FIXED_POINT
michael@0 249 opus_val32
michael@0 250 #else
michael@0 251 void
michael@0 252 #endif
michael@0 253 celt_pitch_xcorr_c(const opus_val16 *_x, const opus_val16 *_y, opus_val32 *xcorr, int len, int max_pitch)
michael@0 254 {
michael@0 255 int i,j;
michael@0 256 /*The EDSP version requires that max_pitch is at least 1, and that _x is
michael@0 257 32-bit aligned.
michael@0 258 Since it's hard to put asserts in assembly, put them here.*/
michael@0 259 celt_assert(max_pitch>0);
michael@0 260 celt_assert((((unsigned char *)_x-(unsigned char *)NULL)&3)==0);
michael@0 261 #ifdef FIXED_POINT
michael@0 262 opus_val32 maxcorr=1;
michael@0 263 #endif
michael@0 264 for (i=0;i<max_pitch-3;i+=4)
michael@0 265 {
michael@0 266 opus_val32 sum[4]={0,0,0,0};
michael@0 267 xcorr_kernel(_x, _y+i, sum, len);
michael@0 268 xcorr[i]=sum[0];
michael@0 269 xcorr[i+1]=sum[1];
michael@0 270 xcorr[i+2]=sum[2];
michael@0 271 xcorr[i+3]=sum[3];
michael@0 272 #ifdef FIXED_POINT
michael@0 273 sum[0] = MAX32(sum[0], sum[1]);
michael@0 274 sum[2] = MAX32(sum[2], sum[3]);
michael@0 275 sum[0] = MAX32(sum[0], sum[2]);
michael@0 276 maxcorr = MAX32(maxcorr, sum[0]);
michael@0 277 #endif
michael@0 278 }
michael@0 279 /* In case max_pitch isn't a multiple of 4, do non-unrolled version. */
michael@0 280 for (;i<max_pitch;i++)
michael@0 281 {
michael@0 282 opus_val32 sum = 0;
michael@0 283 for (j=0;j<len;j++)
michael@0 284 sum = MAC16_16(sum, _x[j],_y[i+j]);
michael@0 285 xcorr[i] = sum;
michael@0 286 #ifdef FIXED_POINT
michael@0 287 maxcorr = MAX32(maxcorr, sum);
michael@0 288 #endif
michael@0 289 }
michael@0 290 #ifdef FIXED_POINT
michael@0 291 return maxcorr;
michael@0 292 #endif
michael@0 293 }
michael@0 294
michael@0 295 #endif
michael@0 296 void pitch_search(const opus_val16 * OPUS_RESTRICT x_lp, opus_val16 * OPUS_RESTRICT y,
michael@0 297 int len, int max_pitch, int *pitch, int arch)
michael@0 298 {
michael@0 299 int i, j;
michael@0 300 int lag;
michael@0 301 int best_pitch[2]={0,0};
michael@0 302 VARDECL(opus_val16, x_lp4);
michael@0 303 VARDECL(opus_val16, y_lp4);
michael@0 304 VARDECL(opus_val32, xcorr);
michael@0 305 #ifdef FIXED_POINT
michael@0 306 opus_val32 maxcorr;
michael@0 307 opus_val32 xmax, ymax;
michael@0 308 int shift=0;
michael@0 309 #endif
michael@0 310 int offset;
michael@0 311
michael@0 312 SAVE_STACK;
michael@0 313
michael@0 314 celt_assert(len>0);
michael@0 315 celt_assert(max_pitch>0);
michael@0 316 lag = len+max_pitch;
michael@0 317
michael@0 318 ALLOC(x_lp4, len>>2, opus_val16);
michael@0 319 ALLOC(y_lp4, lag>>2, opus_val16);
michael@0 320 ALLOC(xcorr, max_pitch>>1, opus_val32);
michael@0 321
michael@0 322 /* Downsample by 2 again */
michael@0 323 for (j=0;j<len>>2;j++)
michael@0 324 x_lp4[j] = x_lp[2*j];
michael@0 325 for (j=0;j<lag>>2;j++)
michael@0 326 y_lp4[j] = y[2*j];
michael@0 327
michael@0 328 #ifdef FIXED_POINT
michael@0 329 xmax = celt_maxabs16(x_lp4, len>>2);
michael@0 330 ymax = celt_maxabs16(y_lp4, lag>>2);
michael@0 331 shift = celt_ilog2(MAX32(1, MAX32(xmax, ymax)))-11;
michael@0 332 if (shift>0)
michael@0 333 {
michael@0 334 for (j=0;j<len>>2;j++)
michael@0 335 x_lp4[j] = SHR16(x_lp4[j], shift);
michael@0 336 for (j=0;j<lag>>2;j++)
michael@0 337 y_lp4[j] = SHR16(y_lp4[j], shift);
michael@0 338 /* Use double the shift for a MAC */
michael@0 339 shift *= 2;
michael@0 340 } else {
michael@0 341 shift = 0;
michael@0 342 }
michael@0 343 #endif
michael@0 344
michael@0 345 /* Coarse search with 4x decimation */
michael@0 346
michael@0 347 #ifdef FIXED_POINT
michael@0 348 maxcorr =
michael@0 349 #endif
michael@0 350 celt_pitch_xcorr(x_lp4, y_lp4, xcorr, len>>2, max_pitch>>2, arch);
michael@0 351
michael@0 352 find_best_pitch(xcorr, y_lp4, len>>2, max_pitch>>2, best_pitch
michael@0 353 #ifdef FIXED_POINT
michael@0 354 , 0, maxcorr
michael@0 355 #endif
michael@0 356 );
michael@0 357
michael@0 358 /* Finer search with 2x decimation */
michael@0 359 #ifdef FIXED_POINT
michael@0 360 maxcorr=1;
michael@0 361 #endif
michael@0 362 for (i=0;i<max_pitch>>1;i++)
michael@0 363 {
michael@0 364 opus_val32 sum=0;
michael@0 365 xcorr[i] = 0;
michael@0 366 if (abs(i-2*best_pitch[0])>2 && abs(i-2*best_pitch[1])>2)
michael@0 367 continue;
michael@0 368 for (j=0;j<len>>1;j++)
michael@0 369 sum += SHR32(MULT16_16(x_lp[j],y[i+j]), shift);
michael@0 370 xcorr[i] = MAX32(-1, sum);
michael@0 371 #ifdef FIXED_POINT
michael@0 372 maxcorr = MAX32(maxcorr, sum);
michael@0 373 #endif
michael@0 374 }
michael@0 375 find_best_pitch(xcorr, y, len>>1, max_pitch>>1, best_pitch
michael@0 376 #ifdef FIXED_POINT
michael@0 377 , shift+1, maxcorr
michael@0 378 #endif
michael@0 379 );
michael@0 380
michael@0 381 /* Refine by pseudo-interpolation */
michael@0 382 if (best_pitch[0]>0 && best_pitch[0]<(max_pitch>>1)-1)
michael@0 383 {
michael@0 384 opus_val32 a, b, c;
michael@0 385 a = xcorr[best_pitch[0]-1];
michael@0 386 b = xcorr[best_pitch[0]];
michael@0 387 c = xcorr[best_pitch[0]+1];
michael@0 388 if ((c-a) > MULT16_32_Q15(QCONST16(.7f,15),b-a))
michael@0 389 offset = 1;
michael@0 390 else if ((a-c) > MULT16_32_Q15(QCONST16(.7f,15),b-c))
michael@0 391 offset = -1;
michael@0 392 else
michael@0 393 offset = 0;
michael@0 394 } else {
michael@0 395 offset = 0;
michael@0 396 }
michael@0 397 *pitch = 2*best_pitch[0]-offset;
michael@0 398
michael@0 399 RESTORE_STACK;
michael@0 400 }
michael@0 401
michael@0 402 static const int second_check[16] = {0, 0, 3, 2, 3, 2, 5, 2, 3, 2, 3, 2, 5, 2, 3, 2};
michael@0 403 opus_val16 remove_doubling(opus_val16 *x, int maxperiod, int minperiod,
michael@0 404 int N, int *T0_, int prev_period, opus_val16 prev_gain)
michael@0 405 {
michael@0 406 int k, i, T, T0;
michael@0 407 opus_val16 g, g0;
michael@0 408 opus_val16 pg;
michael@0 409 opus_val32 xy,xx,yy,xy2;
michael@0 410 opus_val32 xcorr[3];
michael@0 411 opus_val32 best_xy, best_yy;
michael@0 412 int offset;
michael@0 413 int minperiod0;
michael@0 414 VARDECL(opus_val32, yy_lookup);
michael@0 415 SAVE_STACK;
michael@0 416
michael@0 417 minperiod0 = minperiod;
michael@0 418 maxperiod /= 2;
michael@0 419 minperiod /= 2;
michael@0 420 *T0_ /= 2;
michael@0 421 prev_period /= 2;
michael@0 422 N /= 2;
michael@0 423 x += maxperiod;
michael@0 424 if (*T0_>=maxperiod)
michael@0 425 *T0_=maxperiod-1;
michael@0 426
michael@0 427 T = T0 = *T0_;
michael@0 428 ALLOC(yy_lookup, maxperiod+1, opus_val32);
michael@0 429 dual_inner_prod(x, x, x-T0, N, &xx, &xy);
michael@0 430 yy_lookup[0] = xx;
michael@0 431 yy=xx;
michael@0 432 for (i=1;i<=maxperiod;i++)
michael@0 433 {
michael@0 434 yy = yy+MULT16_16(x[-i],x[-i])-MULT16_16(x[N-i],x[N-i]);
michael@0 435 yy_lookup[i] = MAX32(0, yy);
michael@0 436 }
michael@0 437 yy = yy_lookup[T0];
michael@0 438 best_xy = xy;
michael@0 439 best_yy = yy;
michael@0 440 #ifdef FIXED_POINT
michael@0 441 {
michael@0 442 opus_val32 x2y2;
michael@0 443 int sh, t;
michael@0 444 x2y2 = 1+HALF32(MULT32_32_Q31(xx,yy));
michael@0 445 sh = celt_ilog2(x2y2)>>1;
michael@0 446 t = VSHR32(x2y2, 2*(sh-7));
michael@0 447 g = g0 = VSHR32(MULT16_32_Q15(celt_rsqrt_norm(t), xy),sh+1);
michael@0 448 }
michael@0 449 #else
michael@0 450 g = g0 = xy/celt_sqrt(1+xx*yy);
michael@0 451 #endif
michael@0 452 /* Look for any pitch at T/k */
michael@0 453 for (k=2;k<=15;k++)
michael@0 454 {
michael@0 455 int T1, T1b;
michael@0 456 opus_val16 g1;
michael@0 457 opus_val16 cont=0;
michael@0 458 opus_val16 thresh;
michael@0 459 T1 = (2*T0+k)/(2*k);
michael@0 460 if (T1 < minperiod)
michael@0 461 break;
michael@0 462 /* Look for another strong correlation at T1b */
michael@0 463 if (k==2)
michael@0 464 {
michael@0 465 if (T1+T0>maxperiod)
michael@0 466 T1b = T0;
michael@0 467 else
michael@0 468 T1b = T0+T1;
michael@0 469 } else
michael@0 470 {
michael@0 471 T1b = (2*second_check[k]*T0+k)/(2*k);
michael@0 472 }
michael@0 473 dual_inner_prod(x, &x[-T1], &x[-T1b], N, &xy, &xy2);
michael@0 474 xy += xy2;
michael@0 475 yy = yy_lookup[T1] + yy_lookup[T1b];
michael@0 476 #ifdef FIXED_POINT
michael@0 477 {
michael@0 478 opus_val32 x2y2;
michael@0 479 int sh, t;
michael@0 480 x2y2 = 1+MULT32_32_Q31(xx,yy);
michael@0 481 sh = celt_ilog2(x2y2)>>1;
michael@0 482 t = VSHR32(x2y2, 2*(sh-7));
michael@0 483 g1 = VSHR32(MULT16_32_Q15(celt_rsqrt_norm(t), xy),sh+1);
michael@0 484 }
michael@0 485 #else
michael@0 486 g1 = xy/celt_sqrt(1+2.f*xx*1.f*yy);
michael@0 487 #endif
michael@0 488 if (abs(T1-prev_period)<=1)
michael@0 489 cont = prev_gain;
michael@0 490 else if (abs(T1-prev_period)<=2 && 5*k*k < T0)
michael@0 491 cont = HALF32(prev_gain);
michael@0 492 else
michael@0 493 cont = 0;
michael@0 494 thresh = MAX16(QCONST16(.3f,15), MULT16_16_Q15(QCONST16(.7f,15),g0)-cont);
michael@0 495 /* Bias against very high pitch (very short period) to avoid false-positives
michael@0 496 due to short-term correlation */
michael@0 497 if (T1<3*minperiod)
michael@0 498 thresh = MAX16(QCONST16(.4f,15), MULT16_16_Q15(QCONST16(.85f,15),g0)-cont);
michael@0 499 else if (T1<2*minperiod)
michael@0 500 thresh = MAX16(QCONST16(.5f,15), MULT16_16_Q15(QCONST16(.9f,15),g0)-cont);
michael@0 501 if (g1 > thresh)
michael@0 502 {
michael@0 503 best_xy = xy;
michael@0 504 best_yy = yy;
michael@0 505 T = T1;
michael@0 506 g = g1;
michael@0 507 }
michael@0 508 }
michael@0 509 best_xy = MAX32(0, best_xy);
michael@0 510 if (best_yy <= best_xy)
michael@0 511 pg = Q15ONE;
michael@0 512 else
michael@0 513 pg = SHR32(frac_div32(best_xy,best_yy+1),16);
michael@0 514
michael@0 515 for (k=0;k<3;k++)
michael@0 516 {
michael@0 517 int T1 = T+k-1;
michael@0 518 xy = 0;
michael@0 519 for (i=0;i<N;i++)
michael@0 520 xy = MAC16_16(xy, x[i], x[i-T1]);
michael@0 521 xcorr[k] = xy;
michael@0 522 }
michael@0 523 if ((xcorr[2]-xcorr[0]) > MULT16_32_Q15(QCONST16(.7f,15),xcorr[1]-xcorr[0]))
michael@0 524 offset = 1;
michael@0 525 else if ((xcorr[0]-xcorr[2]) > MULT16_32_Q15(QCONST16(.7f,15),xcorr[1]-xcorr[2]))
michael@0 526 offset = -1;
michael@0 527 else
michael@0 528 offset = 0;
michael@0 529 if (pg > g)
michael@0 530 pg = g;
michael@0 531 *T0_ = 2*T+offset;
michael@0 532
michael@0 533 if (*T0_<minperiod0)
michael@0 534 *T0_=minperiod0;
michael@0 535 RESTORE_STACK;
michael@0 536 return pg;
michael@0 537 }

mercurial