1/* libFLAC - Free Lossless Audio Codec library 2 * Copyright (C) 2000,2001,2002,2003,2004,2005,2006,2007 Josh Coalson 3 * 4 * Redistribution and use in source and binary forms, with or without 5 * modification, are permitted provided that the following conditions 6 * are met: 7 * 8 * - Redistributions of source code must retain the above copyright 9 * notice, this list of conditions and the following disclaimer. 10 * 11 * - Redistributions in binary form must reproduce the above copyright 12 * notice, this list of conditions and the following disclaimer in the 13 * documentation and/or other materials provided with the distribution. 14 * 15 * - Neither the name of the Xiph.org Foundation nor the names of its 16 * contributors may be used to endorse or promote products derived from 17 * this software without specific prior written permission. 18 * 19 * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS 20 * ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT 21 * LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR 22 * A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE FOUNDATION OR 23 * CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, 24 * EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, 25 * PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR 26 * PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF 27 * LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING 28 * NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS 29 * SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. 30 */ 31 32#if HAVE_CONFIG_H 33# include <config.h> 34#endif 35 36#include <math.h> 37#include "FLAC/assert.h" 38#include "FLAC/format.h" 39#include "private/bitmath.h" 40#include "private/lpc.h" 41#if defined DEBUG || defined FLAC__OVERFLOW_DETECT || defined FLAC__OVERFLOW_DETECT_VERBOSE 42#include <stdio.h> 43#endif 44 45#ifndef FLAC__INTEGER_ONLY_LIBRARY 46 47#ifndef M_LN2 48/* math.h in VC++ doesn't seem to have this (how Microsoft is that?) */ 49#define M_LN2 0.69314718055994530942 50#endif 51 52/* OPT: #undef'ing this may improve the speed on some architectures */ 53#define FLAC__LPC_UNROLLED_FILTER_LOOPS 54 55 56void FLAC__lpc_window_data(const FLAC__int32 in[], const FLAC__real window[], FLAC__real out[], unsigned data_len) 57{ 58 unsigned i; 59 for(i = 0; i < data_len; i++) 60 out[i] = in[i] * window[i]; 61} 62 63void FLAC__lpc_compute_autocorrelation(const FLAC__real data[], unsigned data_len, unsigned lag, FLAC__real autoc[]) 64{ 65 /* a readable, but slower, version */ 66#if 0 67 FLAC__real d; 68 unsigned i; 69 70 FLAC__ASSERT(lag > 0); 71 FLAC__ASSERT(lag <= data_len); 72 73 /* 74 * Technically we should subtract the mean first like so: 75 * for(i = 0; i < data_len; i++) 76 * data[i] -= mean; 77 * but it appears not to make enough of a difference to matter, and 78 * most signals are already closely centered around zero 79 */ 80 while(lag--) { 81 for(i = lag, d = 0.0; i < data_len; i++) 82 d += data[i] * data[i - lag]; 83 autoc[lag] = d; 84 } 85#endif 86 87 /* 88 * this version tends to run faster because of better data locality 89 * ('data_len' is usually much larger than 'lag') 90 */ 91 FLAC__real d; 92 unsigned sample, coeff; 93 const unsigned limit = data_len - lag; 94 95 FLAC__ASSERT(lag > 0); 96 FLAC__ASSERT(lag <= data_len); 97 98 for(coeff = 0; coeff < lag; coeff++) 99 autoc[coeff] = 0.0; 100 for(sample = 0; sample <= limit; sample++) { 101 d = data[sample]; 102 for(coeff = 0; coeff < lag; coeff++) 103 autoc[coeff] += d * data[sample+coeff]; 104 } 105 for(; sample < data_len; sample++) { 106 d = data[sample]; 107 for(coeff = 0; coeff < data_len - sample; coeff++) 108 autoc[coeff] += d * data[sample+coeff]; 109 } 110} 111 112void FLAC__lpc_compute_lp_coefficients(const FLAC__real autoc[], unsigned *max_order, FLAC__real lp_coeff[][FLAC__MAX_LPC_ORDER], FLAC__double error[]) 113{ 114 unsigned i, j; 115 FLAC__double r, err, ref[FLAC__MAX_LPC_ORDER], lpc[FLAC__MAX_LPC_ORDER]; 116 117 FLAC__ASSERT(0 != max_order); 118 FLAC__ASSERT(0 < *max_order); 119 FLAC__ASSERT(*max_order <= FLAC__MAX_LPC_ORDER); 120 FLAC__ASSERT(autoc[0] != 0.0); 121 122 err = autoc[0]; 123 124 for(i = 0; i < *max_order; i++) { 125 /* Sum up this iteration's reflection coefficient. */ 126 r = -autoc[i+1]; 127 for(j = 0; j < i; j++) 128 r -= lpc[j] * autoc[i-j]; 129 ref[i] = (r/=err); 130 131 /* Update LPC coefficients and total error. */ 132 lpc[i]=r; 133 for(j = 0; j < (i>>1); j++) { 134 FLAC__double tmp = lpc[j]; 135 lpc[j] += r * lpc[i-1-j]; 136 lpc[i-1-j] += r * tmp; 137 } 138 if(i & 1) 139 lpc[j] += lpc[j] * r; 140 141 err *= (1.0 - r * r); 142 143 /* save this order */ 144 for(j = 0; j <= i; j++) 145 lp_coeff[i][j] = (FLAC__real)(-lpc[j]); /* negate FIR filter coeff to get predictor coeff */ 146 error[i] = err; 147 148 /* see SF bug #1601812 http://sourceforge.net/tracker/index.php?func=detail&aid=1601812&group_id=13478&atid=113478 */ 149 if(err == 0.0) { 150 *max_order = i+1; 151 return; 152 } 153 } 154} 155 156int FLAC__lpc_quantize_coefficients(const FLAC__real lp_coeff[], unsigned order, unsigned precision, FLAC__int32 qlp_coeff[], int *shift) 157{ 158 unsigned i; 159 FLAC__double cmax; 160 FLAC__int32 qmax, qmin; 161 162 FLAC__ASSERT(precision > 0); 163 FLAC__ASSERT(precision >= FLAC__MIN_QLP_COEFF_PRECISION); 164 165 /* drop one bit for the sign; from here on out we consider only |lp_coeff[i]| */ 166 precision--; 167 qmax = 1 << precision; 168 qmin = -qmax; 169 qmax--; 170 171 /* calc cmax = max( |lp_coeff[i]| ) */ 172 cmax = 0.0; 173 for(i = 0; i < order; i++) { 174 const FLAC__double d = fabs(lp_coeff[i]); 175 if(d > cmax) 176 cmax = d; 177 } 178 179 if(cmax <= 0.0) { 180 /* => coefficients are all 0, which means our constant-detect didn't work */ 181 return 2; 182 } 183 else { 184 const int max_shiftlimit = (1 << (FLAC__SUBFRAME_LPC_QLP_SHIFT_LEN-1)) - 1; 185 const int min_shiftlimit = -max_shiftlimit - 1; 186 int log2cmax; 187 188 (void)frexp(cmax, &log2cmax); 189 log2cmax--; 190 *shift = (int)precision - log2cmax - 1; 191 192 if(*shift > max_shiftlimit) 193 *shift = max_shiftlimit; 194 else if(*shift < min_shiftlimit) 195 return 1; 196 } 197 198 if(*shift >= 0) { 199 FLAC__double error = 0.0; 200 FLAC__int32 q; 201 for(i = 0; i < order; i++) { 202 error += lp_coeff[i] * (1 << *shift); 203#if 1 /* unfortunately lround() is C99 */ 204 if(error >= 0.0) 205 q = (FLAC__int32)(error + 0.5); 206 else 207 q = (FLAC__int32)(error - 0.5); 208#else 209 q = lround(error); 210#endif 211#ifdef FLAC__OVERFLOW_DETECT 212 if(q > qmax+1) /* we expect q==qmax+1 occasionally due to rounding */ 213 fprintf(stderr,"FLAC__lpc_quantize_coefficients: quantizer overflow: q>qmax %d>%d shift=%d cmax=%f precision=%u lpc[%u]=%f\n",q,qmax,*shift,cmax,precision+1,i,lp_coeff[i]); 214 else if(q < qmin) 215 fprintf(stderr,"FLAC__lpc_quantize_coefficients: quantizer overflow: q<qmin %d<%d shift=%d cmax=%f precision=%u lpc[%u]=%f\n",q,qmin,*shift,cmax,precision+1,i,lp_coeff[i]); 216#endif 217 if(q > qmax) 218 q = qmax; 219 else if(q < qmin) 220 q = qmin; 221 error -= q; 222 qlp_coeff[i] = q; 223 } 224 } 225 /* negative shift is very rare but due to design flaw, negative shift is 226 * a NOP in the decoder, so it must be handled specially by scaling down 227 * coeffs 228 */ 229 else { 230 const int nshift = -(*shift); 231 FLAC__double error = 0.0; 232 FLAC__int32 q; 233#ifdef DEBUG 234 fprintf(stderr,"FLAC__lpc_quantize_coefficients: negative shift=%d order=%u cmax=%f\n", *shift, order, cmax); 235#endif 236 for(i = 0; i < order; i++) { 237 error += lp_coeff[i] / (1 << nshift); 238#if 1 /* unfortunately lround() is C99 */ 239 if(error >= 0.0) 240 q = (FLAC__int32)(error + 0.5); 241 else 242 q = (FLAC__int32)(error - 0.5); 243#else 244 q = lround(error); 245#endif 246#ifdef FLAC__OVERFLOW_DETECT 247 if(q > qmax+1) /* we expect q==qmax+1 occasionally due to rounding */ 248 fprintf(stderr,"FLAC__lpc_quantize_coefficients: quantizer overflow: q>qmax %d>%d shift=%d cmax=%f precision=%u lpc[%u]=%f\n",q,qmax,*shift,cmax,precision+1,i,lp_coeff[i]); 249 else if(q < qmin) 250 fprintf(stderr,"FLAC__lpc_quantize_coefficients: quantizer overflow: q<qmin %d<%d shift=%d cmax=%f precision=%u lpc[%u]=%f\n",q,qmin,*shift,cmax,precision+1,i,lp_coeff[i]); 251#endif 252 if(q > qmax) 253 q = qmax; 254 else if(q < qmin) 255 q = qmin; 256 error -= q; 257 qlp_coeff[i] = q; 258 } 259 *shift = 0; 260 } 261 262 return 0; 263} 264 265void FLAC__lpc_compute_residual_from_qlp_coefficients(const FLAC__int32 *data, unsigned data_len, const FLAC__int32 qlp_coeff[], unsigned order, int lp_quantization, FLAC__int32 residual[]) 266#if defined(FLAC__OVERFLOW_DETECT) || !defined(FLAC__LPC_UNROLLED_FILTER_LOOPS) 267{ 268 FLAC__int64 sumo; 269 unsigned i, j; 270 FLAC__int32 sum; 271 const FLAC__int32 *history; 272 273#ifdef FLAC__OVERFLOW_DETECT_VERBOSE 274 fprintf(stderr,"FLAC__lpc_compute_residual_from_qlp_coefficients: data_len=%d, order=%u, lpq=%d",data_len,order,lp_quantization); 275 for(i=0;i<order;i++) 276 fprintf(stderr,", q[%u]=%d",i,qlp_coeff[i]); 277 fprintf(stderr,"\n"); 278#endif 279 FLAC__ASSERT(order > 0); 280 281 for(i = 0; i < data_len; i++) { 282 sumo = 0; 283 sum = 0; 284 history = data; 285 for(j = 0; j < order; j++) { 286 sum += qlp_coeff[j] * (*(--history)); 287 sumo += (FLAC__int64)qlp_coeff[j] * (FLAC__int64)(*history); 288#if defined _MSC_VER 289 if(sumo > 2147483647I64 || sumo < -2147483648I64) 290 fprintf(stderr,"FLAC__lpc_compute_residual_from_qlp_coefficients: OVERFLOW, i=%u, j=%u, c=%d, d=%d, sumo=%I64d\n",i,j,qlp_coeff[j],*history,sumo); 291#else 292 if(sumo > 2147483647ll || sumo < -2147483648ll) 293 fprintf(stderr,"FLAC__lpc_compute_residual_from_qlp_coefficients: OVERFLOW, i=%u, j=%u, c=%d, d=%d, sumo=%lld\n",i,j,qlp_coeff[j],*history,(long long)sumo); 294#endif 295 } 296 *(residual++) = *(data++) - (sum >> lp_quantization); 297 } 298 299 /* Here's a slower but clearer version: 300 for(i = 0; i < data_len; i++) { 301 sum = 0; 302 for(j = 0; j < order; j++) 303 sum += qlp_coeff[j] * data[i-j-1]; 304 residual[i] = data[i] - (sum >> lp_quantization); 305 } 306 */ 307} 308#else /* fully unrolled version for normal use */ 309{ 310 int i; 311 FLAC__int32 sum; 312 313 FLAC__ASSERT(order > 0); 314 FLAC__ASSERT(order <= 32); 315 316 /* 317 * We do unique versions up to 12th order since that's the subset limit. 318 * Also they are roughly ordered to match frequency of occurrence to 319 * minimize branching. 320 */ 321 if(order <= 12) { 322 if(order > 8) { 323 if(order > 10) { 324 if(order == 12) { 325 for(i = 0; i < (int)data_len; i++) { 326 sum = 0; 327 sum += qlp_coeff[11] * data[i-12]; 328 sum += qlp_coeff[10] * data[i-11]; 329 sum += qlp_coeff[9] * data[i-10]; 330 sum += qlp_coeff[8] * data[i-9]; 331 sum += qlp_coeff[7] * data[i-8]; 332 sum += qlp_coeff[6] * data[i-7]; 333 sum += qlp_coeff[5] * data[i-6]; 334 sum += qlp_coeff[4] * data[i-5]; 335 sum += qlp_coeff[3] * data[i-4]; 336 sum += qlp_coeff[2] * data[i-3]; 337 sum += qlp_coeff[1] * data[i-2]; 338 sum += qlp_coeff[0] * data[i-1]; 339 residual[i] = data[i] - (sum >> lp_quantization); 340 } 341 } 342 else { /* order == 11 */ 343 for(i = 0; i < (int)data_len; i++) { 344 sum = 0; 345 sum += qlp_coeff[10] * data[i-11]; 346 sum += qlp_coeff[9] * data[i-10]; 347 sum += qlp_coeff[8] * data[i-9]; 348 sum += qlp_coeff[7] * data[i-8]; 349 sum += qlp_coeff[6] * data[i-7]; 350 sum += qlp_coeff[5] * data[i-6]; 351 sum += qlp_coeff[4] * data[i-5]; 352 sum += qlp_coeff[3] * data[i-4]; 353 sum += qlp_coeff[2] * data[i-3]; 354 sum += qlp_coeff[1] * data[i-2]; 355 sum += qlp_coeff[0] * data[i-1]; 356 residual[i] = data[i] - (sum >> lp_quantization); 357 } 358 } 359 } 360 else { 361 if(order == 10) { 362 for(i = 0; i < (int)data_len; i++) { 363 sum = 0; 364 sum += qlp_coeff[9] * data[i-10]; 365 sum += qlp_coeff[8] * data[i-9]; 366 sum += qlp_coeff[7] * data[i-8]; 367 sum += qlp_coeff[6] * data[i-7]; 368 sum += qlp_coeff[5] * data[i-6]; 369 sum += qlp_coeff[4] * data[i-5]; 370 sum += qlp_coeff[3] * data[i-4]; 371 sum += qlp_coeff[2] * data[i-3]; 372 sum += qlp_coeff[1] * data[i-2]; 373 sum += qlp_coeff[0] * data[i-1]; 374 residual[i] = data[i] - (sum >> lp_quantization); 375 } 376 } 377 else { /* order == 9 */ 378 for(i = 0; i < (int)data_len; i++) { 379 sum = 0; 380 sum += qlp_coeff[8] * data[i-9]; 381 sum += qlp_coeff[7] * data[i-8]; 382 sum += qlp_coeff[6] * data[i-7]; 383 sum += qlp_coeff[5] * data[i-6]; 384 sum += qlp_coeff[4] * data[i-5]; 385 sum += qlp_coeff[3] * data[i-4]; 386 sum += qlp_coeff[2] * data[i-3]; 387 sum += qlp_coeff[1] * data[i-2]; 388 sum += qlp_coeff[0] * data[i-1]; 389 residual[i] = data[i] - (sum >> lp_quantization); 390 } 391 } 392 } 393 } 394 else if(order > 4) { 395 if(order > 6) { 396 if(order == 8) { 397 for(i = 0; i < (int)data_len; i++) { 398 sum = 0; 399 sum += qlp_coeff[7] * data[i-8]; 400 sum += qlp_coeff[6] * data[i-7]; 401 sum += qlp_coeff[5] * data[i-6]; 402 sum += qlp_coeff[4] * data[i-5]; 403 sum += qlp_coeff[3] * data[i-4]; 404 sum += qlp_coeff[2] * data[i-3]; 405 sum += qlp_coeff[1] * data[i-2]; 406 sum += qlp_coeff[0] * data[i-1]; 407 residual[i] = data[i] - (sum >> lp_quantization); 408 } 409 } 410 else { /* order == 7 */ 411 for(i = 0; i < (int)data_len; i++) { 412 sum = 0; 413 sum += qlp_coeff[6] * data[i-7]; 414 sum += qlp_coeff[5] * data[i-6]; 415 sum += qlp_coeff[4] * data[i-5]; 416 sum += qlp_coeff[3] * data[i-4]; 417 sum += qlp_coeff[2] * data[i-3]; 418 sum += qlp_coeff[1] * data[i-2]; 419 sum += qlp_coeff[0] * data[i-1]; 420 residual[i] = data[i] - (sum >> lp_quantization); 421 } 422 } 423 } 424 else { 425 if(order == 6) { 426 for(i = 0; i < (int)data_len; i++) { 427 sum = 0; 428 sum += qlp_coeff[5] * data[i-6]; 429 sum += qlp_coeff[4] * data[i-5]; 430 sum += qlp_coeff[3] * data[i-4]; 431 sum += qlp_coeff[2] * data[i-3]; 432 sum += qlp_coeff[1] * data[i-2]; 433 sum += qlp_coeff[0] * data[i-1]; 434 residual[i] = data[i] - (sum >> lp_quantization); 435 } 436 } 437 else { /* order == 5 */ 438 for(i = 0; i < (int)data_len; i++) { 439 sum = 0; 440 sum += qlp_coeff[4] * data[i-5]; 441 sum += qlp_coeff[3] * data[i-4]; 442 sum += qlp_coeff[2] * data[i-3]; 443 sum += qlp_coeff[1] * data[i-2]; 444 sum += qlp_coeff[0] * data[i-1]; 445 residual[i] = data[i] - (sum >> lp_quantization); 446 } 447 } 448 } 449 } 450 else { 451 if(order > 2) { 452 if(order == 4) { 453 for(i = 0; i < (int)data_len; i++) { 454 sum = 0; 455 sum += qlp_coeff[3] * data[i-4]; 456 sum += qlp_coeff[2] * data[i-3]; 457 sum += qlp_coeff[1] * data[i-2]; 458 sum += qlp_coeff[0] * data[i-1]; 459 residual[i] = data[i] - (sum >> lp_quantization); 460 } 461 } 462 else { /* order == 3 */ 463 for(i = 0; i < (int)data_len; i++) { 464 sum = 0; 465 sum += qlp_coeff[2] * data[i-3]; 466 sum += qlp_coeff[1] * data[i-2]; 467 sum += qlp_coeff[0] * data[i-1]; 468 residual[i] = data[i] - (sum >> lp_quantization); 469 } 470 } 471 } 472 else { 473 if(order == 2) { 474 for(i = 0; i < (int)data_len; i++) { 475 sum = 0; 476 sum += qlp_coeff[1] * data[i-2]; 477 sum += qlp_coeff[0] * data[i-1]; 478 residual[i] = data[i] - (sum >> lp_quantization); 479 } 480 } 481 else { /* order == 1 */ 482 for(i = 0; i < (int)data_len; i++) 483 residual[i] = data[i] - ((qlp_coeff[0] * data[i-1]) >> lp_quantization); 484 } 485 } 486 } 487 } 488 else { /* order > 12 */ 489 for(i = 0; i < (int)data_len; i++) { 490 sum = 0; 491 switch(order) { 492 case 32: sum += qlp_coeff[31] * data[i-32]; 493 case 31: sum += qlp_coeff[30] * data[i-31]; 494 case 30: sum += qlp_coeff[29] * data[i-30]; 495 case 29: sum += qlp_coeff[28] * data[i-29]; 496 case 28: sum += qlp_coeff[27] * data[i-28]; 497 case 27: sum += qlp_coeff[26] * data[i-27]; 498 case 26: sum += qlp_coeff[25] * data[i-26]; 499 case 25: sum += qlp_coeff[24] * data[i-25]; 500 case 24: sum += qlp_coeff[23] * data[i-24]; 501 case 23: sum += qlp_coeff[22] * data[i-23]; 502 case 22: sum += qlp_coeff[21] * data[i-22]; 503 case 21: sum += qlp_coeff[20] * data[i-21]; 504 case 20: sum += qlp_coeff[19] * data[i-20]; 505 case 19: sum += qlp_coeff[18] * data[i-19]; 506 case 18: sum += qlp_coeff[17] * data[i-18]; 507 case 17: sum += qlp_coeff[16] * data[i-17]; 508 case 16: sum += qlp_coeff[15] * data[i-16]; 509 case 15: sum += qlp_coeff[14] * data[i-15]; 510 case 14: sum += qlp_coeff[13] * data[i-14]; 511 case 13: sum += qlp_coeff[12] * data[i-13]; 512 sum += qlp_coeff[11] * data[i-12]; 513 sum += qlp_coeff[10] * data[i-11]; 514 sum += qlp_coeff[ 9] * data[i-10]; 515 sum += qlp_coeff[ 8] * data[i- 9]; 516 sum += qlp_coeff[ 7] * data[i- 8]; 517 sum += qlp_coeff[ 6] * data[i- 7]; 518 sum += qlp_coeff[ 5] * data[i- 6]; 519 sum += qlp_coeff[ 4] * data[i- 5]; 520 sum += qlp_coeff[ 3] * data[i- 4]; 521 sum += qlp_coeff[ 2] * data[i- 3]; 522 sum += qlp_coeff[ 1] * data[i- 2]; 523 sum += qlp_coeff[ 0] * data[i- 1]; 524 } 525 residual[i] = data[i] - (sum >> lp_quantization); 526 } 527 } 528} 529#endif 530 531void FLAC__lpc_compute_residual_from_qlp_coefficients_wide(const FLAC__int32 *data, unsigned data_len, const FLAC__int32 qlp_coeff[], unsigned order, int lp_quantization, FLAC__int32 residual[]) 532#if defined(FLAC__OVERFLOW_DETECT) || !defined(FLAC__LPC_UNROLLED_FILTER_LOOPS) 533{ 534 unsigned i, j; 535 FLAC__int64 sum; 536 const FLAC__int32 *history; 537 538#ifdef FLAC__OVERFLOW_DETECT_VERBOSE 539 fprintf(stderr,"FLAC__lpc_compute_residual_from_qlp_coefficients_wide: data_len=%d, order=%u, lpq=%d",data_len,order,lp_quantization); 540 for(i=0;i<order;i++) 541 fprintf(stderr,", q[%u]=%d",i,qlp_coeff[i]); 542 fprintf(stderr,"\n"); 543#endif 544 FLAC__ASSERT(order > 0); 545 546 for(i = 0; i < data_len; i++) { 547 sum = 0; 548 history = data; 549 for(j = 0; j < order; j++) 550 sum += (FLAC__int64)qlp_coeff[j] * (FLAC__int64)(*(--history)); 551 if(FLAC__bitmath_silog2_wide(sum >> lp_quantization) > 32) { 552#if defined _MSC_VER 553 fprintf(stderr,"FLAC__lpc_compute_residual_from_qlp_coefficients_wide: OVERFLOW, i=%u, sum=%I64d\n", i, sum >> lp_quantization); 554#else 555 fprintf(stderr,"FLAC__lpc_compute_residual_from_qlp_coefficients_wide: OVERFLOW, i=%u, sum=%lld\n", i, (long long)(sum >> lp_quantization)); 556#endif 557 break; 558 } 559 if(FLAC__bitmath_silog2_wide((FLAC__int64)(*data) - (sum >> lp_quantization)) > 32) { 560#if defined _MSC_VER 561 fprintf(stderr,"FLAC__lpc_compute_residual_from_qlp_coefficients_wide: OVERFLOW, i=%u, data=%d, sum=%I64d, residual=%I64d\n", i, *data, sum >> lp_quantization, (FLAC__int64)(*data) - (sum >> lp_quantization)); 562#else 563 fprintf(stderr,"FLAC__lpc_compute_residual_from_qlp_coefficients_wide: OVERFLOW, i=%u, data=%d, sum=%lld, residual=%lld\n", i, *data, (long long)(sum >> lp_quantization), (long long)((FLAC__int64)(*data) - (sum >> lp_quantization))); 564#endif 565 break; 566 } 567 *(residual++) = *(data++) - (FLAC__int32)(sum >> lp_quantization); 568 } 569} 570#else /* fully unrolled version for normal use */ 571{ 572 int i; 573 FLAC__int64 sum; 574 575 FLAC__ASSERT(order > 0); 576 FLAC__ASSERT(order <= 32); 577 578 /* 579 * We do unique versions up to 12th order since that's the subset limit. 580 * Also they are roughly ordered to match frequency of occurrence to 581 * minimize branching. 582 */ 583 if(order <= 12) { 584 if(order > 8) { 585 if(order > 10) { 586 if(order == 12) { 587 for(i = 0; i < (int)data_len; i++) { 588 sum = 0; 589 sum += qlp_coeff[11] * (FLAC__int64)data[i-12]; 590 sum += qlp_coeff[10] * (FLAC__int64)data[i-11]; 591 sum += qlp_coeff[9] * (FLAC__int64)data[i-10]; 592 sum += qlp_coeff[8] * (FLAC__int64)data[i-9]; 593 sum += qlp_coeff[7] * (FLAC__int64)data[i-8]; 594 sum += qlp_coeff[6] * (FLAC__int64)data[i-7]; 595 sum += qlp_coeff[5] * (FLAC__int64)data[i-6]; 596 sum += qlp_coeff[4] * (FLAC__int64)data[i-5]; 597 sum += qlp_coeff[3] * (FLAC__int64)data[i-4]; 598 sum += qlp_coeff[2] * (FLAC__int64)data[i-3]; 599 sum += qlp_coeff[1] * (FLAC__int64)data[i-2]; 600 sum += qlp_coeff[0] * (FLAC__int64)data[i-1]; 601 residual[i] = data[i] - (FLAC__int32)(sum >> lp_quantization); 602 } 603 } 604 else { /* order == 11 */ 605 for(i = 0; i < (int)data_len; i++) { 606 sum = 0; 607 sum += qlp_coeff[10] * (FLAC__int64)data[i-11]; 608 sum += qlp_coeff[9] * (FLAC__int64)data[i-10]; 609 sum += qlp_coeff[8] * (FLAC__int64)data[i-9]; 610 sum += qlp_coeff[7] * (FLAC__int64)data[i-8]; 611 sum += qlp_coeff[6] * (FLAC__int64)data[i-7]; 612 sum += qlp_coeff[5] * (FLAC__int64)data[i-6]; 613 sum += qlp_coeff[4] * (FLAC__int64)data[i-5]; 614 sum += qlp_coeff[3] * (FLAC__int64)data[i-4]; 615 sum += qlp_coeff[2] * (FLAC__int64)data[i-3]; 616 sum += qlp_coeff[1] * (FLAC__int64)data[i-2]; 617 sum += qlp_coeff[0] * (FLAC__int64)data[i-1]; 618 residual[i] = data[i] - (FLAC__int32)(sum >> lp_quantization); 619 } 620 } 621 } 622 else { 623 if(order == 10) { 624 for(i = 0; i < (int)data_len; i++) { 625 sum = 0; 626 sum += qlp_coeff[9] * (FLAC__int64)data[i-10]; 627 sum += qlp_coeff[8] * (FLAC__int64)data[i-9]; 628 sum += qlp_coeff[7] * (FLAC__int64)data[i-8]; 629 sum += qlp_coeff[6] * (FLAC__int64)data[i-7]; 630 sum += qlp_coeff[5] * (FLAC__int64)data[i-6]; 631 sum += qlp_coeff[4] * (FLAC__int64)data[i-5]; 632 sum += qlp_coeff[3] * (FLAC__int64)data[i-4]; 633 sum += qlp_coeff[2] * (FLAC__int64)data[i-3]; 634 sum += qlp_coeff[1] * (FLAC__int64)data[i-2]; 635 sum += qlp_coeff[0] * (FLAC__int64)data[i-1]; 636 residual[i] = data[i] - (FLAC__int32)(sum >> lp_quantization); 637 } 638 } 639 else { /* order == 9 */ 640 for(i = 0; i < (int)data_len; i++) { 641 sum = 0; 642 sum += qlp_coeff[8] * (FLAC__int64)data[i-9]; 643 sum += qlp_coeff[7] * (FLAC__int64)data[i-8]; 644 sum += qlp_coeff[6] * (FLAC__int64)data[i-7]; 645 sum += qlp_coeff[5] * (FLAC__int64)data[i-6]; 646 sum += qlp_coeff[4] * (FLAC__int64)data[i-5]; 647 sum += qlp_coeff[3] * (FLAC__int64)data[i-4]; 648 sum += qlp_coeff[2] * (FLAC__int64)data[i-3]; 649 sum += qlp_coeff[1] * (FLAC__int64)data[i-2]; 650 sum += qlp_coeff[0] * (FLAC__int64)data[i-1]; 651 residual[i] = data[i] - (FLAC__int32)(sum >> lp_quantization); 652 } 653 } 654 } 655 } 656 else if(order > 4) { 657 if(order > 6) { 658 if(order == 8) { 659 for(i = 0; i < (int)data_len; i++) { 660 sum = 0; 661 sum += qlp_coeff[7] * (FLAC__int64)data[i-8]; 662 sum += qlp_coeff[6] * (FLAC__int64)data[i-7]; 663 sum += qlp_coeff[5] * (FLAC__int64)data[i-6]; 664 sum += qlp_coeff[4] * (FLAC__int64)data[i-5]; 665 sum += qlp_coeff[3] * (FLAC__int64)data[i-4]; 666 sum += qlp_coeff[2] * (FLAC__int64)data[i-3]; 667 sum += qlp_coeff[1] * (FLAC__int64)data[i-2]; 668 sum += qlp_coeff[0] * (FLAC__int64)data[i-1]; 669 residual[i] = data[i] - (FLAC__int32)(sum >> lp_quantization); 670 } 671 } 672 else { /* order == 7 */ 673 for(i = 0; i < (int)data_len; i++) { 674 sum = 0; 675 sum += qlp_coeff[6] * (FLAC__int64)data[i-7]; 676 sum += qlp_coeff[5] * (FLAC__int64)data[i-6]; 677 sum += qlp_coeff[4] * (FLAC__int64)data[i-5]; 678 sum += qlp_coeff[3] * (FLAC__int64)data[i-4]; 679 sum += qlp_coeff[2] * (FLAC__int64)data[i-3]; 680 sum += qlp_coeff[1] * (FLAC__int64)data[i-2]; 681 sum += qlp_coeff[0] * (FLAC__int64)data[i-1]; 682 residual[i] = data[i] - (FLAC__int32)(sum >> lp_quantization); 683 } 684 } 685 } 686 else { 687 if(order == 6) { 688 for(i = 0; i < (int)data_len; i++) { 689 sum = 0; 690 sum += qlp_coeff[5] * (FLAC__int64)data[i-6]; 691 sum += qlp_coeff[4] * (FLAC__int64)data[i-5]; 692 sum += qlp_coeff[3] * (FLAC__int64)data[i-4]; 693 sum += qlp_coeff[2] * (FLAC__int64)data[i-3]; 694 sum += qlp_coeff[1] * (FLAC__int64)data[i-2]; 695 sum += qlp_coeff[0] * (FLAC__int64)data[i-1]; 696 residual[i] = data[i] - (FLAC__int32)(sum >> lp_quantization); 697 } 698 } 699 else { /* order == 5 */ 700 for(i = 0; i < (int)data_len; i++) { 701 sum = 0; 702 sum += qlp_coeff[4] * (FLAC__int64)data[i-5]; 703 sum += qlp_coeff[3] * (FLAC__int64)data[i-4]; 704 sum += qlp_coeff[2] * (FLAC__int64)data[i-3]; 705 sum += qlp_coeff[1] * (FLAC__int64)data[i-2]; 706 sum += qlp_coeff[0] * (FLAC__int64)data[i-1]; 707 residual[i] = data[i] - (FLAC__int32)(sum >> lp_quantization); 708 } 709 } 710 } 711 } 712 else { 713 if(order > 2) { 714 if(order == 4) { 715 for(i = 0; i < (int)data_len; i++) { 716 sum = 0; 717 sum += qlp_coeff[3] * (FLAC__int64)data[i-4]; 718 sum += qlp_coeff[2] * (FLAC__int64)data[i-3]; 719 sum += qlp_coeff[1] * (FLAC__int64)data[i-2]; 720 sum += qlp_coeff[0] * (FLAC__int64)data[i-1]; 721 residual[i] = data[i] - (FLAC__int32)(sum >> lp_quantization); 722 } 723 } 724 else { /* order == 3 */ 725 for(i = 0; i < (int)data_len; i++) { 726 sum = 0; 727 sum += qlp_coeff[2] * (FLAC__int64)data[i-3]; 728 sum += qlp_coeff[1] * (FLAC__int64)data[i-2]; 729 sum += qlp_coeff[0] * (FLAC__int64)data[i-1]; 730 residual[i] = data[i] - (FLAC__int32)(sum >> lp_quantization); 731 } 732 } 733 } 734 else { 735 if(order == 2) { 736 for(i = 0; i < (int)data_len; i++) { 737 sum = 0; 738 sum += qlp_coeff[1] * (FLAC__int64)data[i-2]; 739 sum += qlp_coeff[0] * (FLAC__int64)data[i-1]; 740 residual[i] = data[i] - (FLAC__int32)(sum >> lp_quantization); 741 } 742 } 743 else { /* order == 1 */ 744 for(i = 0; i < (int)data_len; i++) 745 residual[i] = data[i] - (FLAC__int32)((qlp_coeff[0] * (FLAC__int64)data[i-1]) >> lp_quantization); 746 } 747 } 748 } 749 } 750 else { /* order > 12 */ 751 for(i = 0; i < (int)data_len; i++) { 752 sum = 0; 753 switch(order) { 754 case 32: sum += qlp_coeff[31] * (FLAC__int64)data[i-32]; 755 case 31: sum += qlp_coeff[30] * (FLAC__int64)data[i-31]; 756 case 30: sum += qlp_coeff[29] * (FLAC__int64)data[i-30]; 757 case 29: sum += qlp_coeff[28] * (FLAC__int64)data[i-29]; 758 case 28: sum += qlp_coeff[27] * (FLAC__int64)data[i-28]; 759 case 27: sum += qlp_coeff[26] * (FLAC__int64)data[i-27]; 760 case 26: sum += qlp_coeff[25] * (FLAC__int64)data[i-26]; 761 case 25: sum += qlp_coeff[24] * (FLAC__int64)data[i-25]; 762 case 24: sum += qlp_coeff[23] * (FLAC__int64)data[i-24]; 763 case 23: sum += qlp_coeff[22] * (FLAC__int64)data[i-23]; 764 case 22: sum += qlp_coeff[21] * (FLAC__int64)data[i-22]; 765 case 21: sum += qlp_coeff[20] * (FLAC__int64)data[i-21]; 766 case 20: sum += qlp_coeff[19] * (FLAC__int64)data[i-20]; 767 case 19: sum += qlp_coeff[18] * (FLAC__int64)data[i-19]; 768 case 18: sum += qlp_coeff[17] * (FLAC__int64)data[i-18]; 769 case 17: sum += qlp_coeff[16] * (FLAC__int64)data[i-17]; 770 case 16: sum += qlp_coeff[15] * (FLAC__int64)data[i-16]; 771 case 15: sum += qlp_coeff[14] * (FLAC__int64)data[i-15]; 772 case 14: sum += qlp_coeff[13] * (FLAC__int64)data[i-14]; 773 case 13: sum += qlp_coeff[12] * (FLAC__int64)data[i-13]; 774 sum += qlp_coeff[11] * (FLAC__int64)data[i-12]; 775 sum += qlp_coeff[10] * (FLAC__int64)data[i-11]; 776 sum += qlp_coeff[ 9] * (FLAC__int64)data[i-10]; 777 sum += qlp_coeff[ 8] * (FLAC__int64)data[i- 9]; 778 sum += qlp_coeff[ 7] * (FLAC__int64)data[i- 8]; 779 sum += qlp_coeff[ 6] * (FLAC__int64)data[i- 7]; 780 sum += qlp_coeff[ 5] * (FLAC__int64)data[i- 6]; 781 sum += qlp_coeff[ 4] * (FLAC__int64)data[i- 5]; 782 sum += qlp_coeff[ 3] * (FLAC__int64)data[i- 4]; 783 sum += qlp_coeff[ 2] * (FLAC__int64)data[i- 3]; 784 sum += qlp_coeff[ 1] * (FLAC__int64)data[i- 2]; 785 sum += qlp_coeff[ 0] * (FLAC__int64)data[i- 1]; 786 } 787 residual[i] = data[i] - (FLAC__int32)(sum >> lp_quantization); 788 } 789 } 790} 791#endif 792 793#endif /* !defined FLAC__INTEGER_ONLY_LIBRARY */ 794 795void FLAC__lpc_restore_signal(const FLAC__int32 residual[], unsigned data_len, const FLAC__int32 qlp_coeff[], unsigned order, int lp_quantization, FLAC__int32 data[]) 796#if defined(FLAC__OVERFLOW_DETECT) || !defined(FLAC__LPC_UNROLLED_FILTER_LOOPS) 797{ 798 FLAC__int64 sumo; 799 unsigned i, j; 800 FLAC__int32 sum; 801 const FLAC__int32 *r = residual, *history; 802 803#ifdef FLAC__OVERFLOW_DETECT_VERBOSE 804 fprintf(stderr,"FLAC__lpc_restore_signal: data_len=%d, order=%u, lpq=%d",data_len,order,lp_quantization); 805 for(i=0;i<order;i++) 806 fprintf(stderr,", q[%u]=%d",i,qlp_coeff[i]); 807 fprintf(stderr,"\n"); 808#endif 809 FLAC__ASSERT(order > 0); 810 811 for(i = 0; i < data_len; i++) { 812 sumo = 0; 813 sum = 0; 814 history = data; 815 for(j = 0; j < order; j++) { 816 sum += qlp_coeff[j] * (*(--history)); 817 sumo += (FLAC__int64)qlp_coeff[j] * (FLAC__int64)(*history); 818#if defined _MSC_VER 819 if(sumo > 2147483647I64 || sumo < -2147483648I64) 820 fprintf(stderr,"FLAC__lpc_restore_signal: OVERFLOW, i=%u, j=%u, c=%d, d=%d, sumo=%I64d\n",i,j,qlp_coeff[j],*history,sumo); 821#else 822 if(sumo > 2147483647ll || sumo < -2147483648ll) 823 fprintf(stderr,"FLAC__lpc_restore_signal: OVERFLOW, i=%u, j=%u, c=%d, d=%d, sumo=%lld\n",i,j,qlp_coeff[j],*history,(long long)sumo); 824#endif 825 } 826 *(data++) = *(r++) + (sum >> lp_quantization); 827 } 828 829 /* Here's a slower but clearer version: 830 for(i = 0; i < data_len; i++) { 831 sum = 0; 832 for(j = 0; j < order; j++) 833 sum += qlp_coeff[j] * data[i-j-1]; 834 data[i] = residual[i] + (sum >> lp_quantization); 835 } 836 */ 837} 838#else /* fully unrolled version for normal use */ 839{ 840 int i; 841 FLAC__int32 sum; 842 843 FLAC__ASSERT(order > 0); 844 FLAC__ASSERT(order <= 32); 845 846 /* 847 * We do unique versions up to 12th order since that's the subset limit. 848 * Also they are roughly ordered to match frequency of occurrence to 849 * minimize branching. 850 */ 851 if(order <= 12) { 852 if(order > 8) { 853 if(order > 10) { 854 if(order == 12) { 855 for(i = 0; i < (int)data_len; i++) { 856 sum = 0; 857 sum += qlp_coeff[11] * data[i-12]; 858 sum += qlp_coeff[10] * data[i-11]; 859 sum += qlp_coeff[9] * data[i-10]; 860 sum += qlp_coeff[8] * data[i-9]; 861 sum += qlp_coeff[7] * data[i-8]; 862 sum += qlp_coeff[6] * data[i-7]; 863 sum += qlp_coeff[5] * data[i-6]; 864 sum += qlp_coeff[4] * data[i-5]; 865 sum += qlp_coeff[3] * data[i-4]; 866 sum += qlp_coeff[2] * data[i-3]; 867 sum += qlp_coeff[1] * data[i-2]; 868 sum += qlp_coeff[0] * data[i-1]; 869 data[i] = residual[i] + (sum >> lp_quantization); 870 } 871 } 872 else { /* order == 11 */ 873 for(i = 0; i < (int)data_len; i++) { 874 sum = 0; 875 sum += qlp_coeff[10] * data[i-11]; 876 sum += qlp_coeff[9] * data[i-10]; 877 sum += qlp_coeff[8] * data[i-9]; 878 sum += qlp_coeff[7] * data[i-8]; 879 sum += qlp_coeff[6] * data[i-7]; 880 sum += qlp_coeff[5] * data[i-6]; 881 sum += qlp_coeff[4] * data[i-5]; 882 sum += qlp_coeff[3] * data[i-4]; 883 sum += qlp_coeff[2] * data[i-3]; 884 sum += qlp_coeff[1] * data[i-2]; 885 sum += qlp_coeff[0] * data[i-1]; 886 data[i] = residual[i] + (sum >> lp_quantization); 887 } 888 } 889 } 890 else { 891 if(order == 10) { 892 for(i = 0; i < (int)data_len; i++) { 893 sum = 0; 894 sum += qlp_coeff[9] * data[i-10]; 895 sum += qlp_coeff[8] * data[i-9]; 896 sum += qlp_coeff[7] * data[i-8]; 897 sum += qlp_coeff[6] * data[i-7]; 898 sum += qlp_coeff[5] * data[i-6]; 899 sum += qlp_coeff[4] * data[i-5]; 900 sum += qlp_coeff[3] * data[i-4]; 901 sum += qlp_coeff[2] * data[i-3]; 902 sum += qlp_coeff[1] * data[i-2]; 903 sum += qlp_coeff[0] * data[i-1]; 904 data[i] = residual[i] + (sum >> lp_quantization); 905 } 906 } 907 else { /* order == 9 */ 908 for(i = 0; i < (int)data_len; i++) { 909 sum = 0; 910 sum += qlp_coeff[8] * data[i-9]; 911 sum += qlp_coeff[7] * data[i-8]; 912 sum += qlp_coeff[6] * data[i-7]; 913 sum += qlp_coeff[5] * data[i-6]; 914 sum += qlp_coeff[4] * data[i-5]; 915 sum += qlp_coeff[3] * data[i-4]; 916 sum += qlp_coeff[2] * data[i-3]; 917 sum += qlp_coeff[1] * data[i-2]; 918 sum += qlp_coeff[0] * data[i-1]; 919 data[i] = residual[i] + (sum >> lp_quantization); 920 } 921 } 922 } 923 } 924 else if(order > 4) { 925 if(order > 6) { 926 if(order == 8) { 927 for(i = 0; i < (int)data_len; i++) { 928 sum = 0; 929 sum += qlp_coeff[7] * data[i-8]; 930 sum += qlp_coeff[6] * data[i-7]; 931 sum += qlp_coeff[5] * data[i-6]; 932 sum += qlp_coeff[4] * data[i-5]; 933 sum += qlp_coeff[3] * data[i-4]; 934 sum += qlp_coeff[2] * data[i-3]; 935 sum += qlp_coeff[1] * data[i-2]; 936 sum += qlp_coeff[0] * data[i-1]; 937 data[i] = residual[i] + (sum >> lp_quantization); 938 } 939 } 940 else { /* order == 7 */ 941 for(i = 0; i < (int)data_len; i++) { 942 sum = 0; 943 sum += qlp_coeff[6] * data[i-7]; 944 sum += qlp_coeff[5] * data[i-6]; 945 sum += qlp_coeff[4] * data[i-5]; 946 sum += qlp_coeff[3] * data[i-4]; 947 sum += qlp_coeff[2] * data[i-3]; 948 sum += qlp_coeff[1] * data[i-2]; 949 sum += qlp_coeff[0] * data[i-1]; 950 data[i] = residual[i] + (sum >> lp_quantization); 951 } 952 } 953 } 954 else { 955 if(order == 6) { 956 for(i = 0; i < (int)data_len; i++) { 957 sum = 0; 958 sum += qlp_coeff[5] * data[i-6]; 959 sum += qlp_coeff[4] * data[i-5]; 960 sum += qlp_coeff[3] * data[i-4]; 961 sum += qlp_coeff[2] * data[i-3]; 962 sum += qlp_coeff[1] * data[i-2]; 963 sum += qlp_coeff[0] * data[i-1]; 964 data[i] = residual[i] + (sum >> lp_quantization); 965 } 966 } 967 else { /* order == 5 */ 968 for(i = 0; i < (int)data_len; i++) { 969 sum = 0; 970 sum += qlp_coeff[4] * data[i-5]; 971 sum += qlp_coeff[3] * data[i-4]; 972 sum += qlp_coeff[2] * data[i-3]; 973 sum += qlp_coeff[1] * data[i-2]; 974 sum += qlp_coeff[0] * data[i-1]; 975 data[i] = residual[i] + (sum >> lp_quantization); 976 } 977 } 978 } 979 } 980 else { 981 if(order > 2) { 982 if(order == 4) { 983 for(i = 0; i < (int)data_len; i++) { 984 sum = 0; 985 sum += qlp_coeff[3] * data[i-4]; 986 sum += qlp_coeff[2] * data[i-3]; 987 sum += qlp_coeff[1] * data[i-2]; 988 sum += qlp_coeff[0] * data[i-1]; 989 data[i] = residual[i] + (sum >> lp_quantization); 990 } 991 } 992 else { /* order == 3 */ 993 for(i = 0; i < (int)data_len; i++) { 994 sum = 0; 995 sum += qlp_coeff[2] * data[i-3]; 996 sum += qlp_coeff[1] * data[i-2]; 997 sum += qlp_coeff[0] * data[i-1]; 998 data[i] = residual[i] + (sum >> lp_quantization); 999 } 1000 } 1001 } 1002 else { 1003 if(order == 2) { 1004 for(i = 0; i < (int)data_len; i++) { 1005 sum = 0; 1006 sum += qlp_coeff[1] * data[i-2]; 1007 sum += qlp_coeff[0] * data[i-1]; 1008 data[i] = residual[i] + (sum >> lp_quantization); 1009 } 1010 } 1011 else { /* order == 1 */ 1012 for(i = 0; i < (int)data_len; i++) 1013 data[i] = residual[i] + ((qlp_coeff[0] * data[i-1]) >> lp_quantization); 1014 } 1015 } 1016 } 1017 } 1018 else { /* order > 12 */ 1019 for(i = 0; i < (int)data_len; i++) { 1020 sum = 0; 1021 switch(order) { 1022 case 32: sum += qlp_coeff[31] * data[i-32]; 1023 case 31: sum += qlp_coeff[30] * data[i-31]; 1024 case 30: sum += qlp_coeff[29] * data[i-30]; 1025 case 29: sum += qlp_coeff[28] * data[i-29]; 1026 case 28: sum += qlp_coeff[27] * data[i-28]; 1027 case 27: sum += qlp_coeff[26] * data[i-27]; 1028 case 26: sum += qlp_coeff[25] * data[i-26]; 1029 case 25: sum += qlp_coeff[24] * data[i-25]; 1030 case 24: sum += qlp_coeff[23] * data[i-24]; 1031 case 23: sum += qlp_coeff[22] * data[i-23]; 1032 case 22: sum += qlp_coeff[21] * data[i-22]; 1033 case 21: sum += qlp_coeff[20] * data[i-21]; 1034 case 20: sum += qlp_coeff[19] * data[i-20]; 1035 case 19: sum += qlp_coeff[18] * data[i-19]; 1036 case 18: sum += qlp_coeff[17] * data[i-18]; 1037 case 17: sum += qlp_coeff[16] * data[i-17]; 1038 case 16: sum += qlp_coeff[15] * data[i-16]; 1039 case 15: sum += qlp_coeff[14] * data[i-15]; 1040 case 14: sum += qlp_coeff[13] * data[i-14]; 1041 case 13: sum += qlp_coeff[12] * data[i-13]; 1042 sum += qlp_coeff[11] * data[i-12]; 1043 sum += qlp_coeff[10] * data[i-11]; 1044 sum += qlp_coeff[ 9] * data[i-10]; 1045 sum += qlp_coeff[ 8] * data[i- 9]; 1046 sum += qlp_coeff[ 7] * data[i- 8]; 1047 sum += qlp_coeff[ 6] * data[i- 7]; 1048 sum += qlp_coeff[ 5] * data[i- 6]; 1049 sum += qlp_coeff[ 4] * data[i- 5]; 1050 sum += qlp_coeff[ 3] * data[i- 4]; 1051 sum += qlp_coeff[ 2] * data[i- 3]; 1052 sum += qlp_coeff[ 1] * data[i- 2]; 1053 sum += qlp_coeff[ 0] * data[i- 1]; 1054 } 1055 data[i] = residual[i] + (sum >> lp_quantization); 1056 } 1057 } 1058} 1059#endif 1060 1061void FLAC__lpc_restore_signal_wide(const FLAC__int32 residual[], unsigned data_len, const FLAC__int32 qlp_coeff[], unsigned order, int lp_quantization, FLAC__int32 data[]) 1062#if defined(FLAC__OVERFLOW_DETECT) || !defined(FLAC__LPC_UNROLLED_FILTER_LOOPS) 1063{ 1064 unsigned i, j; 1065 FLAC__int64 sum; 1066 const FLAC__int32 *r = residual, *history; 1067 1068#ifdef FLAC__OVERFLOW_DETECT_VERBOSE 1069 fprintf(stderr,"FLAC__lpc_restore_signal_wide: data_len=%d, order=%u, lpq=%d",data_len,order,lp_quantization); 1070 for(i=0;i<order;i++) 1071 fprintf(stderr,", q[%u]=%d",i,qlp_coeff[i]); 1072 fprintf(stderr,"\n"); 1073#endif 1074 FLAC__ASSERT(order > 0); 1075 1076 for(i = 0; i < data_len; i++) { 1077 sum = 0; 1078 history = data; 1079 for(j = 0; j < order; j++) 1080 sum += (FLAC__int64)qlp_coeff[j] * (FLAC__int64)(*(--history)); 1081 if(FLAC__bitmath_silog2_wide(sum >> lp_quantization) > 32) { 1082#ifdef _MSC_VER 1083 fprintf(stderr,"FLAC__lpc_restore_signal_wide: OVERFLOW, i=%u, sum=%I64d\n", i, sum >> lp_quantization); 1084#else 1085 fprintf(stderr,"FLAC__lpc_restore_signal_wide: OVERFLOW, i=%u, sum=%lld\n", i, (long long)(sum >> lp_quantization)); 1086#endif 1087 break; 1088 } 1089 if(FLAC__bitmath_silog2_wide((FLAC__int64)(*r) + (sum >> lp_quantization)) > 32) { 1090#ifdef _MSC_VER 1091 fprintf(stderr,"FLAC__lpc_restore_signal_wide: OVERFLOW, i=%u, residual=%d, sum=%I64d, data=%I64d\n", i, *r, sum >> lp_quantization, (FLAC__int64)(*r) + (sum >> lp_quantization)); 1092#else 1093 fprintf(stderr,"FLAC__lpc_restore_signal_wide: OVERFLOW, i=%u, residual=%d, sum=%lld, data=%lld\n", i, *r, (long long)(sum >> lp_quantization), (long long)((FLAC__int64)(*r) + (sum >> lp_quantization))); 1094#endif 1095 break; 1096 } 1097 *(data++) = *(r++) + (FLAC__int32)(sum >> lp_quantization); 1098 } 1099} 1100#else /* fully unrolled version for normal use */ 1101{ 1102 int i; 1103 FLAC__int64 sum; 1104 1105 FLAC__ASSERT(order > 0); 1106 FLAC__ASSERT(order <= 32); 1107 1108 /* 1109 * We do unique versions up to 12th order since that's the subset limit. 1110 * Also they are roughly ordered to match frequency of occurrence to 1111 * minimize branching. 1112 */ 1113 if(order <= 12) { 1114 if(order > 8) { 1115 if(order > 10) { 1116 if(order == 12) { 1117 for(i = 0; i < (int)data_len; i++) { 1118 sum = 0; 1119 sum += qlp_coeff[11] * (FLAC__int64)data[i-12]; 1120 sum += qlp_coeff[10] * (FLAC__int64)data[i-11]; 1121 sum += qlp_coeff[9] * (FLAC__int64)data[i-10]; 1122 sum += qlp_coeff[8] * (FLAC__int64)data[i-9]; 1123 sum += qlp_coeff[7] * (FLAC__int64)data[i-8]; 1124 sum += qlp_coeff[6] * (FLAC__int64)data[i-7]; 1125 sum += qlp_coeff[5] * (FLAC__int64)data[i-6]; 1126 sum += qlp_coeff[4] * (FLAC__int64)data[i-5]; 1127 sum += qlp_coeff[3] * (FLAC__int64)data[i-4]; 1128 sum += qlp_coeff[2] * (FLAC__int64)data[i-3]; 1129 sum += qlp_coeff[1] * (FLAC__int64)data[i-2]; 1130 sum += qlp_coeff[0] * (FLAC__int64)data[i-1]; 1131 data[i] = residual[i] + (FLAC__int32)(sum >> lp_quantization); 1132 } 1133 } 1134 else { /* order == 11 */ 1135 for(i = 0; i < (int)data_len; i++) { 1136 sum = 0; 1137 sum += qlp_coeff[10] * (FLAC__int64)data[i-11]; 1138 sum += qlp_coeff[9] * (FLAC__int64)data[i-10]; 1139 sum += qlp_coeff[8] * (FLAC__int64)data[i-9]; 1140 sum += qlp_coeff[7] * (FLAC__int64)data[i-8]; 1141 sum += qlp_coeff[6] * (FLAC__int64)data[i-7]; 1142 sum += qlp_coeff[5] * (FLAC__int64)data[i-6]; 1143 sum += qlp_coeff[4] * (FLAC__int64)data[i-5]; 1144 sum += qlp_coeff[3] * (FLAC__int64)data[i-4]; 1145 sum += qlp_coeff[2] * (FLAC__int64)data[i-3]; 1146 sum += qlp_coeff[1] * (FLAC__int64)data[i-2]; 1147 sum += qlp_coeff[0] * (FLAC__int64)data[i-1]; 1148 data[i] = residual[i] + (FLAC__int32)(sum >> lp_quantization); 1149 } 1150 } 1151 } 1152 else { 1153 if(order == 10) { 1154 for(i = 0; i < (int)data_len; i++) { 1155 sum = 0; 1156 sum += qlp_coeff[9] * (FLAC__int64)data[i-10]; 1157 sum += qlp_coeff[8] * (FLAC__int64)data[i-9]; 1158 sum += qlp_coeff[7] * (FLAC__int64)data[i-8]; 1159 sum += qlp_coeff[6] * (FLAC__int64)data[i-7]; 1160 sum += qlp_coeff[5] * (FLAC__int64)data[i-6]; 1161 sum += qlp_coeff[4] * (FLAC__int64)data[i-5]; 1162 sum += qlp_coeff[3] * (FLAC__int64)data[i-4]; 1163 sum += qlp_coeff[2] * (FLAC__int64)data[i-3]; 1164 sum += qlp_coeff[1] * (FLAC__int64)data[i-2]; 1165 sum += qlp_coeff[0] * (FLAC__int64)data[i-1]; 1166 data[i] = residual[i] + (FLAC__int32)(sum >> lp_quantization); 1167 } 1168 } 1169 else { /* order == 9 */ 1170 for(i = 0; i < (int)data_len; i++) { 1171 sum = 0; 1172 sum += qlp_coeff[8] * (FLAC__int64)data[i-9]; 1173 sum += qlp_coeff[7] * (FLAC__int64)data[i-8]; 1174 sum += qlp_coeff[6] * (FLAC__int64)data[i-7]; 1175 sum += qlp_coeff[5] * (FLAC__int64)data[i-6]; 1176 sum += qlp_coeff[4] * (FLAC__int64)data[i-5]; 1177 sum += qlp_coeff[3] * (FLAC__int64)data[i-4]; 1178 sum += qlp_coeff[2] * (FLAC__int64)data[i-3]; 1179 sum += qlp_coeff[1] * (FLAC__int64)data[i-2]; 1180 sum += qlp_coeff[0] * (FLAC__int64)data[i-1]; 1181 data[i] = residual[i] + (FLAC__int32)(sum >> lp_quantization); 1182 } 1183 } 1184 } 1185 } 1186 else if(order > 4) { 1187 if(order > 6) { 1188 if(order == 8) { 1189 for(i = 0; i < (int)data_len; i++) { 1190 sum = 0; 1191 sum += qlp_coeff[7] * (FLAC__int64)data[i-8]; 1192 sum += qlp_coeff[6] * (FLAC__int64)data[i-7]; 1193 sum += qlp_coeff[5] * (FLAC__int64)data[i-6]; 1194 sum += qlp_coeff[4] * (FLAC__int64)data[i-5]; 1195 sum += qlp_coeff[3] * (FLAC__int64)data[i-4]; 1196 sum += qlp_coeff[2] * (FLAC__int64)data[i-3]; 1197 sum += qlp_coeff[1] * (FLAC__int64)data[i-2]; 1198 sum += qlp_coeff[0] * (FLAC__int64)data[i-1]; 1199 data[i] = residual[i] + (FLAC__int32)(sum >> lp_quantization); 1200 } 1201 } 1202 else { /* order == 7 */ 1203 for(i = 0; i < (int)data_len; i++) { 1204 sum = 0; 1205 sum += qlp_coeff[6] * (FLAC__int64)data[i-7]; 1206 sum += qlp_coeff[5] * (FLAC__int64)data[i-6]; 1207 sum += qlp_coeff[4] * (FLAC__int64)data[i-5]; 1208 sum += qlp_coeff[3] * (FLAC__int64)data[i-4]; 1209 sum += qlp_coeff[2] * (FLAC__int64)data[i-3]; 1210 sum += qlp_coeff[1] * (FLAC__int64)data[i-2]; 1211 sum += qlp_coeff[0] * (FLAC__int64)data[i-1]; 1212 data[i] = residual[i] + (FLAC__int32)(sum >> lp_quantization); 1213 } 1214 } 1215 } 1216 else { 1217 if(order == 6) { 1218 for(i = 0; i < (int)data_len; i++) { 1219 sum = 0; 1220 sum += qlp_coeff[5] * (FLAC__int64)data[i-6]; 1221 sum += qlp_coeff[4] * (FLAC__int64)data[i-5]; 1222 sum += qlp_coeff[3] * (FLAC__int64)data[i-4]; 1223 sum += qlp_coeff[2] * (FLAC__int64)data[i-3]; 1224 sum += qlp_coeff[1] * (FLAC__int64)data[i-2]; 1225 sum += qlp_coeff[0] * (FLAC__int64)data[i-1]; 1226 data[i] = residual[i] + (FLAC__int32)(sum >> lp_quantization); 1227 } 1228 } 1229 else { /* order == 5 */ 1230 for(i = 0; i < (int)data_len; i++) { 1231 sum = 0; 1232 sum += qlp_coeff[4] * (FLAC__int64)data[i-5]; 1233 sum += qlp_coeff[3] * (FLAC__int64)data[i-4]; 1234 sum += qlp_coeff[2] * (FLAC__int64)data[i-3]; 1235 sum += qlp_coeff[1] * (FLAC__int64)data[i-2]; 1236 sum += qlp_coeff[0] * (FLAC__int64)data[i-1]; 1237 data[i] = residual[i] + (FLAC__int32)(sum >> lp_quantization); 1238 } 1239 } 1240 } 1241 } 1242 else { 1243 if(order > 2) { 1244 if(order == 4) { 1245 for(i = 0; i < (int)data_len; i++) { 1246 sum = 0; 1247 sum += qlp_coeff[3] * (FLAC__int64)data[i-4]; 1248 sum += qlp_coeff[2] * (FLAC__int64)data[i-3]; 1249 sum += qlp_coeff[1] * (FLAC__int64)data[i-2]; 1250 sum += qlp_coeff[0] * (FLAC__int64)data[i-1]; 1251 data[i] = residual[i] + (FLAC__int32)(sum >> lp_quantization); 1252 } 1253 } 1254 else { /* order == 3 */ 1255 for(i = 0; i < (int)data_len; i++) { 1256 sum = 0; 1257 sum += qlp_coeff[2] * (FLAC__int64)data[i-3]; 1258 sum += qlp_coeff[1] * (FLAC__int64)data[i-2]; 1259 sum += qlp_coeff[0] * (FLAC__int64)data[i-1]; 1260 data[i] = residual[i] + (FLAC__int32)(sum >> lp_quantization); 1261 } 1262 } 1263 } 1264 else { 1265 if(order == 2) { 1266 for(i = 0; i < (int)data_len; i++) { 1267 sum = 0; 1268 sum += qlp_coeff[1] * (FLAC__int64)data[i-2]; 1269 sum += qlp_coeff[0] * (FLAC__int64)data[i-1]; 1270 data[i] = residual[i] + (FLAC__int32)(sum >> lp_quantization); 1271 } 1272 } 1273 else { /* order == 1 */ 1274 for(i = 0; i < (int)data_len; i++) 1275 data[i] = residual[i] + (FLAC__int32)((qlp_coeff[0] * (FLAC__int64)data[i-1]) >> lp_quantization); 1276 } 1277 } 1278 } 1279 } 1280 else { /* order > 12 */ 1281 for(i = 0; i < (int)data_len; i++) { 1282 sum = 0; 1283 switch(order) { 1284 case 32: sum += qlp_coeff[31] * (FLAC__int64)data[i-32]; 1285 case 31: sum += qlp_coeff[30] * (FLAC__int64)data[i-31]; 1286 case 30: sum += qlp_coeff[29] * (FLAC__int64)data[i-30]; 1287 case 29: sum += qlp_coeff[28] * (FLAC__int64)data[i-29]; 1288 case 28: sum += qlp_coeff[27] * (FLAC__int64)data[i-28]; 1289 case 27: sum += qlp_coeff[26] * (FLAC__int64)data[i-27]; 1290 case 26: sum += qlp_coeff[25] * (FLAC__int64)data[i-26]; 1291 case 25: sum += qlp_coeff[24] * (FLAC__int64)data[i-25]; 1292 case 24: sum += qlp_coeff[23] * (FLAC__int64)data[i-24]; 1293 case 23: sum += qlp_coeff[22] * (FLAC__int64)data[i-23]; 1294 case 22: sum += qlp_coeff[21] * (FLAC__int64)data[i-22]; 1295 case 21: sum += qlp_coeff[20] * (FLAC__int64)data[i-21]; 1296 case 20: sum += qlp_coeff[19] * (FLAC__int64)data[i-20]; 1297 case 19: sum += qlp_coeff[18] * (FLAC__int64)data[i-19]; 1298 case 18: sum += qlp_coeff[17] * (FLAC__int64)data[i-18]; 1299 case 17: sum += qlp_coeff[16] * (FLAC__int64)data[i-17]; 1300 case 16: sum += qlp_coeff[15] * (FLAC__int64)data[i-16]; 1301 case 15: sum += qlp_coeff[14] * (FLAC__int64)data[i-15]; 1302 case 14: sum += qlp_coeff[13] * (FLAC__int64)data[i-14]; 1303 case 13: sum += qlp_coeff[12] * (FLAC__int64)data[i-13]; 1304 sum += qlp_coeff[11] * (FLAC__int64)data[i-12]; 1305 sum += qlp_coeff[10] * (FLAC__int64)data[i-11]; 1306 sum += qlp_coeff[ 9] * (FLAC__int64)data[i-10]; 1307 sum += qlp_coeff[ 8] * (FLAC__int64)data[i- 9]; 1308 sum += qlp_coeff[ 7] * (FLAC__int64)data[i- 8]; 1309 sum += qlp_coeff[ 6] * (FLAC__int64)data[i- 7]; 1310 sum += qlp_coeff[ 5] * (FLAC__int64)data[i- 6]; 1311 sum += qlp_coeff[ 4] * (FLAC__int64)data[i- 5]; 1312 sum += qlp_coeff[ 3] * (FLAC__int64)data[i- 4]; 1313 sum += qlp_coeff[ 2] * (FLAC__int64)data[i- 3]; 1314 sum += qlp_coeff[ 1] * (FLAC__int64)data[i- 2]; 1315 sum += qlp_coeff[ 0] * (FLAC__int64)data[i- 1]; 1316 } 1317 data[i] = residual[i] + (FLAC__int32)(sum >> lp_quantization); 1318 } 1319 } 1320} 1321#endif 1322 1323#ifndef FLAC__INTEGER_ONLY_LIBRARY 1324 1325FLAC__double FLAC__lpc_compute_expected_bits_per_residual_sample(FLAC__double lpc_error, unsigned total_samples) 1326{ 1327 FLAC__double error_scale; 1328 1329 FLAC__ASSERT(total_samples > 0); 1330 1331 error_scale = 0.5 * M_LN2 * M_LN2 / (FLAC__double)total_samples; 1332 1333 return FLAC__lpc_compute_expected_bits_per_residual_sample_with_error_scale(lpc_error, error_scale); 1334} 1335 1336FLAC__double FLAC__lpc_compute_expected_bits_per_residual_sample_with_error_scale(FLAC__double lpc_error, FLAC__double error_scale) 1337{ 1338 if(lpc_error > 0.0) { 1339 FLAC__double bps = (FLAC__double)0.5 * log(error_scale * lpc_error) / M_LN2; 1340 if(bps >= 0.0) 1341 return bps; 1342 else 1343 return 0.0; 1344 } 1345 else if(lpc_error < 0.0) { /* error should not be negative but can happen due to inadequate floating-point resolution */ 1346 return 1e32; 1347 } 1348 else { 1349 return 0.0; 1350 } 1351} 1352 1353unsigned FLAC__lpc_compute_best_order(const FLAC__double lpc_error[], unsigned max_order, unsigned total_samples, unsigned overhead_bits_per_order) 1354{ 1355 unsigned order, index, best_index; /* 'index' the index into lpc_error; index==order-1 since lpc_error[0] is for order==1, lpc_error[1] is for order==2, etc */ 1356 FLAC__double bits, best_bits, error_scale; 1357 1358 FLAC__ASSERT(max_order > 0); 1359 FLAC__ASSERT(total_samples > 0); 1360 1361 error_scale = 0.5 * M_LN2 * M_LN2 / (FLAC__double)total_samples; 1362 1363 best_index = 0; 1364 best_bits = (unsigned)(-1); 1365 1366 for(index = 0, order = 1; index < max_order; index++, order++) { 1367 bits = FLAC__lpc_compute_expected_bits_per_residual_sample_with_error_scale(lpc_error[index], error_scale) * (FLAC__double)(total_samples - order) + (FLAC__double)(order * overhead_bits_per_order); 1368 if(bits < best_bits) { 1369 best_index = index; 1370 best_bits = bits; 1371 } 1372 } 1373 1374 return best_index+1; /* +1 since index of lpc_error[] is order-1 */ 1375} 1376 1377#endif /* !defined FLAC__INTEGER_ONLY_LIBRARY */ 1378