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