/* libFLAC - Free Lossless Audio Codec library * Copyright (C) 2000-2009  Josh Coalson * Copyright (C) 2011-2016  Xiph.Org Foundation * * Redistribution and use in source and binary forms, with or without * modification, are permitted provided that the following conditions * are met: * * - Redistributions of source code must retain the above copyright * notice, this list of conditions and the following disclaimer. * * - Redistributions in binary form must reproduce the above copyright * notice, this list of conditions and the following disclaimer in the * documentation and/or other materials provided with the distribution. * * - Neither the name of the Xiph.org Foundation nor the names of its * contributors may be used to endorse or promote products derived from * this software without specific prior written permission. * * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS * ``AS IS'' AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT * LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR * A PARTICULAR PURPOSE ARE DISCLAIMED.  IN NO EVENT SHALL THE FOUNDATION OR * CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, * EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, * PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR * PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF * LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING * NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS * SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. */ #include <math.h>#include "FLAC/format.h"#include "share/compat.h"#include "private/bitmath.h"#include "private/lpc.h"#include "private/macros.h" /* OPT: #undef'ing this may improve the speed on some architectures */#define FLAC__LPC_UNROLLED_FILTER_LOOPS #if defined(_MSC_VER) && (_MSC_VER < 1800)#include <float.h>static inline long int lround(double x) {	return (long)(x + _copysign(0.5, x));}#endif void FLAC__lpc_window_data(const FLAC__int32 in[], const FLAC__real window[], FLAC__real out[], uint32_t data_len){	uint32_t i;	for(i = 0; i < data_len; i++)		out[i] = in[i] * window[i];} void FLAC__lpc_compute_autocorrelation(const FLAC__real data[], uint32_t data_len, uint32_t lag, FLAC__real autoc[]){	/*	 * this version tends to run faster because of better data locality	 * ('data_len' is usually much larger than 'lag')	 */	FLAC__real d;	uint32_t sample, coeff;	const uint32_t limit = data_len - lag; 	for(coeff = 0; coeff < lag; coeff++)		autoc[coeff] = 0.0;	for(sample = 0; sample <= limit; sample++) {		d = data[sample];		for(coeff = 0; coeff < lag; coeff++)			autoc[coeff] += d * data[sample+coeff];	}	for(; sample < data_len; sample++) {		d = data[sample];		for(coeff = 0; coeff < data_len - sample; coeff++)			autoc[coeff] += d * data[sample+coeff];	}} void FLAC__lpc_compute_lp_coefficients(const FLAC__real autoc[], uint32_t *max_order, FLAC__real lp_coeff[][FLAC__MAX_LPC_ORDER], double error[]){	uint32_t i, j;	double r, err, lpc[FLAC__MAX_LPC_ORDER]; 	err = autoc[0]; 	for(i = 0; i < *max_order; i++) {		/* Sum up this iteration's reflection coefficient. */		r = -autoc[i+1];		for(j = 0; j < i; j++)			r -= lpc[j] * autoc[i-j];		r /= err; 		/* Update LPC coefficients and total error. */		lpc[i]=r;		for(j = 0; j < (i>>1); j++) {			double tmp = lpc[j];			lpc[j] += r * lpc[i-1-j];			lpc[i-1-j] += r * tmp;		}		if(i & 1)			lpc[j] += lpc[j] * r; 		err *= (1.0 - r * r); 		/* save this order */		for(j = 0; j <= i; j++)			lp_coeff[i][j] = (FLAC__real)(-lpc[j]); /* negate FIR filter coeff to get predictor coeff */		error[i] = err; 		/* see SF bug https://sourceforge.net/p/flac/bugs/234/ */		if(err == 0.0) {			*max_order = i+1;			return;		}	}} int FLAC__lpc_quantize_coefficients(const FLAC__real lp_coeff[], uint32_t order, uint32_t precision, FLAC__int32 qlp_coeff[], int *shift){	uint32_t i;	double cmax;	FLAC__int32 qmax, qmin; 	/* drop one bit for the sign; from here on out we consider only |lp_coeff[i]| */	precision--;	qmax = 1 << precision;	qmin = -qmax;	qmax--; 	/* calc cmax = max( |lp_coeff[i]| ) */	cmax = 0.0;	for(i = 0; i < order; i++) {		const double d = fabs(lp_coeff[i]);		if(d > cmax)			cmax = d;	} 	if(cmax <= 0.0) {		/* => coefficients are all 0, which means our constant-detect didn't work */		return 2;	}	else {		const int max_shiftlimit = (1 << (FLAC__SUBFRAME_LPC_QLP_SHIFT_LEN-1)) - 1;		const int min_shiftlimit = -max_shiftlimit - 1;		int log2cmax; 		(void)frexp(cmax, &log2cmax);		log2cmax--;		*shift = (int)precision - log2cmax - 1; 		if(*shift > max_shiftlimit)			*shift = max_shiftlimit;		else if(*shift < min_shiftlimit)			return 1;	} 	if(*shift >= 0) {		double error = 0.0;		FLAC__int32 q;		for(i = 0; i < order; i++) {			error += lp_coeff[i] * (1 << *shift);			q = lround(error); 			if(q > qmax)				q = qmax;			else if(q < qmin)				q = qmin;			error -= q;			qlp_coeff[i] = q;		}	}	/* negative shift is very rare but due to design flaw, negative shift is	 * not allowed in the decoder, so it must be handled specially by scaling	 * down coeffs	 */	else {		const int nshift = -(*shift);		double error = 0.0;		FLAC__int32 q; 		for(i = 0; i < order; i++) {			error += lp_coeff[i] / (1 << nshift);			q = lround(error);			if(q > qmax)				q = qmax;			else if(q < qmin)				q = qmin;			error -= q;			qlp_coeff[i] = q;		}		*shift = 0;	} 	return 0;} #if defined(_MSC_VER)// silence MSVC warnings about __restrict modifier#pragma warning ( disable : 4028 )#endif 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)#if !defined(FLAC__LPC_UNROLLED_FILTER_LOOPS){	FLAC__int64 sumo;	uint32_t i, j;	FLAC__int32 sum;	const FLAC__int32 *history; 	FLAC__ASSERT(order > 0); 	for(i = 0; i < data_len; i++) {		sumo = 0;		sum = 0;		history = data;		for(j = 0; j < order; j++) {			sum += qlp_coeff[j] * (*(--history));			sumo += (FLAC__int64)qlp_coeff[j] * (FLAC__int64)(*history);			if(sumo > 2147483647ll || sumo < -2147483648ll)				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);		}		*(residual++) = *(data++) - (sum >> lp_quantization);	}}#else /* fully unrolled version for normal use */{	int i;	FLAC__int32 sum; 	/*	 * We do unique versions up to 12th order since that's the subset limit.	 * Also they are roughly ordered to match frequency of occurrence to	 * minimize branching.	 */	if(order <= 12) {		if(order > 8) {			if(order > 10) {				if(order == 12) {					for(i = 0; i < (int)data_len; i++) {						sum = 0;						sum += qlp_coeff[11] * data[i-12];						sum += qlp_coeff[10] * data[i-11];						sum += qlp_coeff[9] * data[i-10];						sum += qlp_coeff[8] * data[i-9];						sum += qlp_coeff[7] * data[i-8];						sum += qlp_coeff[6] * data[i-7];						sum += qlp_coeff[5] * data[i-6];						sum += qlp_coeff[4] * data[i-5];						sum += qlp_coeff[3] * data[i-4];						sum += qlp_coeff[2] * data[i-3];						sum += qlp_coeff[1] * data[i-2];						sum += qlp_coeff[0] * data[i-1];						residual[i] = data[i] - (sum >> lp_quantization);					}				}				else { /* order == 11 */					for(i = 0; i < (int)data_len; i++) {						sum = 0;						sum += qlp_coeff[10] * data[i-11];						sum += qlp_coeff[9] * data[i-10];						sum += qlp_coeff[8] * data[i-9];						sum += qlp_coeff[7] * data[i-8];						sum += qlp_coeff[6] * data[i-7];						sum += qlp_coeff[5] * data[i-6];						sum += qlp_coeff[4] * data[i-5];						sum += qlp_coeff[3] * data[i-4];						sum += qlp_coeff[2] * data[i-3];						sum += qlp_coeff[1] * data[i-2];						sum += qlp_coeff[0] * data[i-1];						residual[i] = data[i] - (sum >> lp_quantization);					}				}			}			else {				if(order == 10) {					for(i = 0; i < (int)data_len; i++) {						sum = 0;						sum += qlp_coeff[9] * data[i-10];						sum += qlp_coeff[8] * data[i-9];						sum += qlp_coeff[7] * data[i-8];						sum += qlp_coeff[6] * data[i-7];						sum += qlp_coeff[5] * data[i-6];						sum += qlp_coeff[4] * data[i-5];						sum += qlp_coeff[3] * data[i-4];						sum += qlp_coeff[2] * data[i-3];						sum += qlp_coeff[1] * data[i-2];						sum += qlp_coeff[0] * data[i-1];						residual[i] = data[i] - (sum >> lp_quantization);					}				}				else { /* order == 9 */					for(i = 0; i < (int)data_len; i++) {						sum = 0;						sum += qlp_coeff[8] * data[i-9];						sum += qlp_coeff[7] * data[i-8];						sum += qlp_coeff[6] * data[i-7];						sum += qlp_coeff[5] * data[i-6];						sum += qlp_coeff[4] * data[i-5];						sum += qlp_coeff[3] * data[i-4];						sum += qlp_coeff[2] * data[i-3];						sum += qlp_coeff[1] * data[i-2];						sum += qlp_coeff[0] * data[i-1];						residual[i] = data[i] - (sum >> lp_quantization);					}				}			}		}		else if(order > 4) {			if(order > 6) {				if(order == 8) {					for(i = 0; i < (int)data_len; i++) {						sum = 0;						sum += qlp_coeff[7] * data[i-8];						sum += qlp_coeff[6] * data[i-7];						sum += qlp_coeff[5] * data[i-6];						sum += qlp_coeff[4] * data[i-5];						sum += qlp_coeff[3] * data[i-4];						sum += qlp_coeff[2] * data[i-3];						sum += qlp_coeff[1] * data[i-2];						sum += qlp_coeff[0] * data[i-1];						residual[i] = data[i] - (sum >> lp_quantization);					}				}				else { /* order == 7 */					for(i = 0; i < (int)data_len; i++) {						sum = 0;						sum += qlp_coeff[6] * data[i-7];						sum += qlp_coeff[5] * data[i-6];						sum += qlp_coeff[4] * data[i-5];						sum += qlp_coeff[3] * data[i-4];						sum += qlp_coeff[2] * data[i-3];						sum += qlp_coeff[1] * data[i-2];						sum += qlp_coeff[0] * data[i-1];						residual[i] = data[i] - (sum >> lp_quantization);					}				}			}			else {				if(order == 6) {					for(i = 0; i < (int)data_len; i++) {						sum = 0;						sum += qlp_coeff[5] * data[i-6];						sum += qlp_coeff[4] * data[i-5];						sum += qlp_coeff[3] * data[i-4];						sum += qlp_coeff[2] * data[i-3];						sum += qlp_coeff[1] * data[i-2];						sum += qlp_coeff[0] * data[i-1];						residual[i] = data[i] - (sum >> lp_quantization);					}				}				else { /* order == 5 */					for(i = 0; i < (int)data_len; i++) {						sum = 0;						sum += qlp_coeff[4] * data[i-5];						sum += qlp_coeff[3] * data[i-4];						sum += qlp_coeff[2] * data[i-3];						sum += qlp_coeff[1] * data[i-2];						sum += qlp_coeff[0] * data[i-1];						residual[i] = data[i] - (sum >> lp_quantization);					}				}			}		}		else {			if(order > 2) {				if(order == 4) {					for(i = 0; i < (int)data_len; i++) {						sum = 0;						sum += qlp_coeff[3] * data[i-4];						sum += qlp_coeff[2] * data[i-3];						sum += qlp_coeff[1] * data[i-2];						sum += qlp_coeff[0] * data[i-1];						residual[i] = data[i] - (sum >> lp_quantization);					}				}				else { /* order == 3 */					for(i = 0; i < (int)data_len; i++) {						sum = 0;						sum += qlp_coeff[2] * data[i-3];						sum += qlp_coeff[1] * data[i-2];						sum += qlp_coeff[0] * data[i-1];						residual[i] = data[i] - (sum >> lp_quantization);					}				}			}			else {				if(order == 2) {					for(i = 0; i < (int)data_len; i++) {						sum = 0;						sum += qlp_coeff[1] * data[i-2];						sum += qlp_coeff[0] * data[i-1];						residual[i] = data[i] - (sum >> lp_quantization);					}				}				else { /* order == 1 */					for(i = 0; i < (int)data_len; i++)						residual[i] = data[i] - ((qlp_coeff[0] * data[i-1]) >> lp_quantization);				}			}		}	}	else { /* order > 12 */		for(i = 0; i < (int)data_len; i++) {			sum = 0;			switch(order) {				default:				case 32: sum += qlp_coeff[31] * data[i-32]; /* Falls through. */				case 31: sum += qlp_coeff[30] * data[i-31]; /* Falls through. */				case 30: sum += qlp_coeff[29] * data[i-30]; /* Falls through. */				case 29: sum += qlp_coeff[28] * data[i-29]; /* Falls through. */				case 28: sum += qlp_coeff[27] * data[i-28]; /* Falls through. */				case 27: sum += qlp_coeff[26] * data[i-27]; /* Falls through. */				case 26: sum += qlp_coeff[25] * data[i-26]; /* Falls through. */				case 25: sum += qlp_coeff[24] * data[i-25]; /* Falls through. */				case 24: sum += qlp_coeff[23] * data[i-24]; /* Falls through. */				case 23: sum += qlp_coeff[22] * data[i-23]; /* Falls through. */				case 22: sum += qlp_coeff[21] * data[i-22]; /* Falls through. */				case 21: sum += qlp_coeff[20] * data[i-21]; /* Falls through. */				case 20: sum += qlp_coeff[19] * data[i-20]; /* Falls through. */				case 19: sum += qlp_coeff[18] * data[i-19]; /* Falls through. */				case 18: sum += qlp_coeff[17] * data[i-18]; /* Falls through. */				case 17: sum += qlp_coeff[16] * data[i-17]; /* Falls through. */				case 16: sum += qlp_coeff[15] * data[i-16]; /* Falls through. */				case 15: sum += qlp_coeff[14] * data[i-15]; /* Falls through. */				case 14: sum += qlp_coeff[13] * data[i-14]; /* Falls through. */				case 13: sum += qlp_coeff[12] * data[i-13];				         sum += qlp_coeff[11] * data[i-12];				         sum += qlp_coeff[10] * data[i-11];				         sum += qlp_coeff[ 9] * data[i-10];				         sum += qlp_coeff[ 8] * data[i- 9];				         sum += qlp_coeff[ 7] * data[i- 8];				         sum += qlp_coeff[ 6] * data[i- 7];				         sum += qlp_coeff[ 5] * data[i- 6];				         sum += qlp_coeff[ 4] * data[i- 5];				         sum += qlp_coeff[ 3] * data[i- 4];				         sum += qlp_coeff[ 2] * data[i- 3];				         sum += qlp_coeff[ 1] * data[i- 2];				         sum += qlp_coeff[ 0] * data[i- 1];			}			residual[i] = data[i] - (sum >> lp_quantization);		}	}}#endif 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)#if !defined(FLAC__LPC_UNROLLED_FILTER_LOOPS){	uint32_t i, j;	FLAC__int64 sum;	const FLAC__int32 *history; 	FLAC__ASSERT(order > 0); 	for(i = 0; i < data_len; i++) {		sum = 0;		history = data;		for(j = 0; j < order; j++)			sum += (FLAC__int64)qlp_coeff[j] * (FLAC__int64)(*(--history));		if(FLAC__bitmath_silog2(sum >> lp_quantization) > 32) {			fprintf(stderr,"FLAC__lpc_compute_residual_from_qlp_coefficients_wide: OVERFLOW, i=%u, sum=%" PRId64 "\n", i, (sum >> lp_quantization));			break;		}		if(FLAC__bitmath_silog2((FLAC__int64)(*data) - (sum >> lp_quantization)) > 32) {			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)));			break;		}		*(residual++) = *(data++) - (FLAC__int32)(sum >> lp_quantization);	}}#else /* fully unrolled version for normal use */{	int i;	FLAC__int64 sum; 	/*	 * We do unique versions up to 12th order since that's the subset limit.	 * Also they are roughly ordered to match frequency of occurrence to	 * minimize branching.	 */	if(order <= 12) {		if(order > 8) {			if(order > 10) {				if(order == 12) {					for(i = 0; i < (int)data_len; i++) {						sum = 0;						sum += qlp_coeff[11] * (FLAC__int64)data[i-12];						sum += qlp_coeff[10] * (FLAC__int64)data[i-11];						sum += qlp_coeff[9] * (FLAC__int64)data[i-10];						sum += qlp_coeff[8] * (FLAC__int64)data[i-9];						sum += qlp_coeff[7] * (FLAC__int64)data[i-8];						sum += qlp_coeff[6] * (FLAC__int64)data[i-7];						sum += qlp_coeff[5] * (FLAC__int64)data[i-6];						sum += qlp_coeff[4] * (FLAC__int64)data[i-5];						sum += qlp_coeff[3] * (FLAC__int64)data[i-4];						sum += qlp_coeff[2] * (FLAC__int64)data[i-3];						sum += qlp_coeff[1] * (FLAC__int64)data[i-2];						sum += qlp_coeff[0] * (FLAC__int64)data[i-1];

Lines 1–500 of 1233. Load more