Josh Coalson | 26560dd | 2001-02-08 00:38:41 +0000 | [diff] [blame] | 1 | /* libFLAC - Free Lossless Audio Codec library |
Josh Coalson | 305ae2e | 2002-01-26 17:36:39 +0000 | [diff] [blame] | 2 | * Copyright (C) 2000,2001,2002 Josh Coalson |
Josh Coalson | bb7f6b9 | 2000-12-10 04:09:52 +0000 | [diff] [blame] | 3 | * |
| 4 | * This library is free software; you can redistribute it and/or |
| 5 | * modify it under the terms of the GNU Library General Public |
| 6 | * License as published by the Free Software Foundation; either |
| 7 | * version 2 of the License, or (at your option) any later version. |
| 8 | * |
| 9 | * This library is distributed in the hope that it will be useful, |
| 10 | * but WITHOUT ANY WARRANTY; without even the implied warranty of |
| 11 | * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU |
| 12 | * Library General Public License for more details. |
| 13 | * |
| 14 | * You should have received a copy of the GNU Library General Public |
| 15 | * License along with this library; if not, write to the |
| 16 | * Free Software Foundation, Inc., 59 Temple Place - Suite 330, |
| 17 | * Boston, MA 02111-1307, USA. |
| 18 | */ |
| 19 | |
Josh Coalson | bb7f6b9 | 2000-12-10 04:09:52 +0000 | [diff] [blame] | 20 | #include <math.h> |
Josh Coalson | 1b68982 | 2001-05-31 20:11:02 +0000 | [diff] [blame] | 21 | #include "FLAC/assert.h" |
Josh Coalson | bb7f6b9 | 2000-12-10 04:09:52 +0000 | [diff] [blame] | 22 | #include "FLAC/format.h" |
| 23 | #include "private/lpc.h" |
Josh Coalson | 03ed88b | 2002-05-17 06:19:28 +0000 | [diff] [blame^] | 24 | #if defined DEBUG || defined FLAC__OVERFLOW_DETECT || defined FLAC__OVERFLOW_DETECT_VERBOSE |
| 25 | #include <stdio.h> |
| 26 | #endif |
Josh Coalson | bb7f6b9 | 2000-12-10 04:09:52 +0000 | [diff] [blame] | 27 | |
| 28 | #ifndef M_LN2 |
| 29 | /* math.h in VC++ doesn't seem to have this (how Microsoft is that?) */ |
| 30 | #define M_LN2 0.69314718055994530942 |
| 31 | #endif |
| 32 | |
Josh Coalson | 77e3f31 | 2001-06-23 03:03:24 +0000 | [diff] [blame] | 33 | void FLAC__lpc_compute_autocorrelation(const FLAC__real data[], unsigned data_len, unsigned lag, FLAC__real autoc[]) |
Josh Coalson | bb7f6b9 | 2000-12-10 04:09:52 +0000 | [diff] [blame] | 34 | { |
Josh Coalson | c0785cd | 2001-05-10 19:29:41 +0000 | [diff] [blame] | 35 | /* a readable, but slower, version */ |
| 36 | #if 0 |
Josh Coalson | 77e3f31 | 2001-06-23 03:03:24 +0000 | [diff] [blame] | 37 | FLAC__real d; |
Josh Coalson | bb7f6b9 | 2000-12-10 04:09:52 +0000 | [diff] [blame] | 38 | unsigned i; |
| 39 | |
Josh Coalson | 1b68982 | 2001-05-31 20:11:02 +0000 | [diff] [blame] | 40 | FLAC__ASSERT(lag > 0); |
| 41 | FLAC__ASSERT(lag <= data_len); |
Josh Coalson | bb7f6b9 | 2000-12-10 04:09:52 +0000 | [diff] [blame] | 42 | |
| 43 | while(lag--) { |
| 44 | for(i = lag, d = 0.0; i < data_len; i++) |
| 45 | d += data[i] * data[i - lag]; |
| 46 | autoc[lag] = d; |
| 47 | } |
Josh Coalson | c0785cd | 2001-05-10 19:29:41 +0000 | [diff] [blame] | 48 | #endif |
| 49 | |
| 50 | /* |
| 51 | * this version tends to run faster because of better data locality |
| 52 | * ('data_len' is usually much larger than 'lag') |
| 53 | */ |
Josh Coalson | 77e3f31 | 2001-06-23 03:03:24 +0000 | [diff] [blame] | 54 | FLAC__real d; |
Josh Coalson | c0785cd | 2001-05-10 19:29:41 +0000 | [diff] [blame] | 55 | unsigned sample, coeff; |
| 56 | const unsigned limit = data_len - lag; |
| 57 | |
Josh Coalson | 1b68982 | 2001-05-31 20:11:02 +0000 | [diff] [blame] | 58 | FLAC__ASSERT(lag > 0); |
| 59 | FLAC__ASSERT(lag <= data_len); |
Josh Coalson | c0785cd | 2001-05-10 19:29:41 +0000 | [diff] [blame] | 60 | |
| 61 | for(coeff = 0; coeff < lag; coeff++) |
| 62 | autoc[coeff] = 0.0; |
Josh Coalson | 376807d | 2001-05-16 19:23:35 +0000 | [diff] [blame] | 63 | for(sample = 0; sample <= limit; sample++) { |
Josh Coalson | c0785cd | 2001-05-10 19:29:41 +0000 | [diff] [blame] | 64 | d = data[sample]; |
| 65 | for(coeff = 0; coeff < lag; coeff++) |
| 66 | autoc[coeff] += d * data[sample+coeff]; |
| 67 | } |
Josh Coalson | 376807d | 2001-05-16 19:23:35 +0000 | [diff] [blame] | 68 | for(; sample < data_len; sample++) { |
Josh Coalson | c0785cd | 2001-05-10 19:29:41 +0000 | [diff] [blame] | 69 | d = data[sample]; |
| 70 | for(coeff = 0; coeff < data_len - sample; coeff++) |
| 71 | autoc[coeff] += d * data[sample+coeff]; |
| 72 | } |
Josh Coalson | bb7f6b9 | 2000-12-10 04:09:52 +0000 | [diff] [blame] | 73 | } |
| 74 | |
Josh Coalson | 77e3f31 | 2001-06-23 03:03:24 +0000 | [diff] [blame] | 75 | void FLAC__lpc_compute_lp_coefficients(const FLAC__real autoc[], unsigned max_order, FLAC__real lp_coeff[][FLAC__MAX_LPC_ORDER], FLAC__real error[]) |
Josh Coalson | bb7f6b9 | 2000-12-10 04:09:52 +0000 | [diff] [blame] | 76 | { |
| 77 | unsigned i, j; |
Josh Coalson | b35bebd | 2001-07-03 04:37:18 +0000 | [diff] [blame] | 78 | double r, err, ref[FLAC__MAX_LPC_ORDER], lpc[FLAC__MAX_LPC_ORDER]; |
Josh Coalson | bb7f6b9 | 2000-12-10 04:09:52 +0000 | [diff] [blame] | 79 | |
Josh Coalson | 1b68982 | 2001-05-31 20:11:02 +0000 | [diff] [blame] | 80 | FLAC__ASSERT(0 < max_order); |
| 81 | FLAC__ASSERT(max_order <= FLAC__MAX_LPC_ORDER); |
| 82 | FLAC__ASSERT(autoc[0] != 0.0); |
Josh Coalson | bb7f6b9 | 2000-12-10 04:09:52 +0000 | [diff] [blame] | 83 | |
| 84 | err = autoc[0]; |
| 85 | |
| 86 | for(i = 0; i < max_order; i++) { |
| 87 | /* Sum up this iteration's reflection coefficient. */ |
Josh Coalson | a1b53c4 | 2001-05-24 19:27:08 +0000 | [diff] [blame] | 88 | r = -autoc[i+1]; |
Josh Coalson | bb7f6b9 | 2000-12-10 04:09:52 +0000 | [diff] [blame] | 89 | for(j = 0; j < i; j++) |
| 90 | r -= lpc[j] * autoc[i-j]; |
| 91 | ref[i] = (r/=err); |
| 92 | |
| 93 | /* Update LPC coefficients and total error. */ |
| 94 | lpc[i]=r; |
| 95 | for(j = 0; j < (i>>1); j++) { |
Josh Coalson | b35bebd | 2001-07-03 04:37:18 +0000 | [diff] [blame] | 96 | double tmp = lpc[j]; |
Josh Coalson | bb7f6b9 | 2000-12-10 04:09:52 +0000 | [diff] [blame] | 97 | lpc[j] += r * lpc[i-1-j]; |
| 98 | lpc[i-1-j] += r * tmp; |
| 99 | } |
| 100 | if(i & 1) |
| 101 | lpc[j] += lpc[j] * r; |
| 102 | |
| 103 | err *= (1.0 - r * r); |
| 104 | |
| 105 | /* save this order */ |
| 106 | for(j = 0; j <= i; j++) |
Josh Coalson | b35bebd | 2001-07-03 04:37:18 +0000 | [diff] [blame] | 107 | lp_coeff[i][j] = (FLAC__real)(-lpc[j]); /* negate FIR filter coeff to get predictor coeff */ |
| 108 | error[i] = (FLAC__real)err; |
Josh Coalson | bb7f6b9 | 2000-12-10 04:09:52 +0000 | [diff] [blame] | 109 | } |
| 110 | } |
| 111 | |
Josh Coalson | 77e3f31 | 2001-06-23 03:03:24 +0000 | [diff] [blame] | 112 | int FLAC__lpc_quantize_coefficients(const FLAC__real lp_coeff[], unsigned order, unsigned precision, unsigned bits_per_sample, FLAC__int32 qlp_coeff[], int *shift) |
Josh Coalson | bb7f6b9 | 2000-12-10 04:09:52 +0000 | [diff] [blame] | 113 | { |
| 114 | unsigned i; |
Josh Coalson | b35bebd | 2001-07-03 04:37:18 +0000 | [diff] [blame] | 115 | double d, cmax = -1e32; |
Josh Coalson | 5975e8a | 2001-07-03 04:10:21 +0000 | [diff] [blame] | 116 | FLAC__int32 qmax, qmin; |
| 117 | const int max_shiftlimit = (1 << (FLAC__SUBFRAME_LPC_QLP_SHIFT_LEN-1)) - 1; |
| 118 | const int min_shiftlimit = -max_shiftlimit - 1; |
Josh Coalson | c625f90 | 2001-02-28 23:56:03 +0000 | [diff] [blame] | 119 | |
Josh Coalson | 1b68982 | 2001-05-31 20:11:02 +0000 | [diff] [blame] | 120 | FLAC__ASSERT(bits_per_sample > 0); |
Josh Coalson | 77e3f31 | 2001-06-23 03:03:24 +0000 | [diff] [blame] | 121 | FLAC__ASSERT(bits_per_sample <= sizeof(FLAC__int32)*8); |
Josh Coalson | 1b68982 | 2001-05-31 20:11:02 +0000 | [diff] [blame] | 122 | FLAC__ASSERT(precision > 0); |
| 123 | FLAC__ASSERT(precision >= FLAC__MIN_QLP_COEFF_PRECISION); |
Josh Coalson | 77e3f31 | 2001-06-23 03:03:24 +0000 | [diff] [blame] | 124 | FLAC__ASSERT(precision + bits_per_sample < sizeof(FLAC__int32)*8); |
Josh Coalson | c625f90 | 2001-02-28 23:56:03 +0000 | [diff] [blame] | 125 | #ifdef NDEBUG |
| 126 | (void)bits_per_sample; /* silence compiler warning about unused parameter */ |
| 127 | #endif |
| 128 | |
| 129 | /* drop one bit for the sign; from here on out we consider only |lp_coeff[i]| */ |
| 130 | precision--; |
Josh Coalson | 5975e8a | 2001-07-03 04:10:21 +0000 | [diff] [blame] | 131 | qmax = 1 << precision; |
| 132 | qmin = -qmax; |
| 133 | qmax--; |
Josh Coalson | c625f90 | 2001-02-28 23:56:03 +0000 | [diff] [blame] | 134 | |
| 135 | for(i = 0; i < order; i++) { |
| 136 | if(lp_coeff[i] == 0.0) |
| 137 | continue; |
Josh Coalson | f52360a | 2001-08-13 23:10:06 +0000 | [diff] [blame] | 138 | d = fabs(lp_coeff[i]); |
Josh Coalson | c625f90 | 2001-02-28 23:56:03 +0000 | [diff] [blame] | 139 | if(d > cmax) |
| 140 | cmax = d; |
Josh Coalson | c625f90 | 2001-02-28 23:56:03 +0000 | [diff] [blame] | 141 | } |
Josh Coalson | 5975e8a | 2001-07-03 04:10:21 +0000 | [diff] [blame] | 142 | redo_it: |
Josh Coalson | d3ed49f | 2001-07-18 23:47:19 +0000 | [diff] [blame] | 143 | if(cmax <= 0.0) { |
Josh Coalson | 456db82 | 2001-03-01 19:14:05 +0000 | [diff] [blame] | 144 | /* => coefficients are all 0, which means our constant-detect didn't work */ |
Josh Coalson | c625f90 | 2001-02-28 23:56:03 +0000 | [diff] [blame] | 145 | return 2; |
| 146 | } |
| 147 | else { |
Josh Coalson | 8c6f90f | 2001-07-19 17:07:13 +0000 | [diff] [blame] | 148 | int log2cmax; |
Josh Coalson | c625f90 | 2001-02-28 23:56:03 +0000 | [diff] [blame] | 149 | |
Josh Coalson | 8c6f90f | 2001-07-19 17:07:13 +0000 | [diff] [blame] | 150 | (void)frexp(cmax, &log2cmax); |
| 151 | log2cmax--; |
| 152 | *shift = (int)precision - log2cmax - 1; |
Josh Coalson | c625f90 | 2001-02-28 23:56:03 +0000 | [diff] [blame] | 153 | |
| 154 | if(*shift < min_shiftlimit || *shift > max_shiftlimit) { |
Josh Coalson | c625f90 | 2001-02-28 23:56:03 +0000 | [diff] [blame] | 155 | return 1; |
| 156 | } |
| 157 | } |
| 158 | |
Josh Coalson | 4e6b3ac | 2001-07-06 00:37:57 +0000 | [diff] [blame] | 159 | if(*shift >= 0) { |
| 160 | for(i = 0; i < order; i++) { |
| 161 | qlp_coeff[i] = (FLAC__int32)floor((double)lp_coeff[i] * (double)(1 << *shift)); |
Josh Coalson | 5975e8a | 2001-07-03 04:10:21 +0000 | [diff] [blame] | 162 | |
Josh Coalson | 8c6f90f | 2001-07-19 17:07:13 +0000 | [diff] [blame] | 163 | /* double-check the result */ |
Josh Coalson | 4e6b3ac | 2001-07-06 00:37:57 +0000 | [diff] [blame] | 164 | if(qlp_coeff[i] > qmax || qlp_coeff[i] < qmin) { |
Josh Coalson | 5975e8a | 2001-07-03 04:10:21 +0000 | [diff] [blame] | 165 | #ifdef FLAC__OVERFLOW_DETECT |
Josh Coalson | 4e6b3ac | 2001-07-06 00:37:57 +0000 | [diff] [blame] | 166 | fprintf(stderr,"FLAC__lpc_quantize_coefficients: compensating for overflow, qlp_coeff[%u]=%d, lp_coeff[%u]=%f, cmax=%f, precision=%u, shift=%d, q=%f, f(q)=%f\n", i, qlp_coeff[i], i, lp_coeff[i], cmax, precision, *shift, (double)lp_coeff[i] * (double)(1 << *shift), floor((double)lp_coeff[i] * (double)(1 << *shift))); |
Josh Coalson | 5975e8a | 2001-07-03 04:10:21 +0000 | [diff] [blame] | 167 | #endif |
Josh Coalson | 4e6b3ac | 2001-07-06 00:37:57 +0000 | [diff] [blame] | 168 | cmax *= 2.0; |
| 169 | goto redo_it; |
Josh Coalson | 5975e8a | 2001-07-03 04:10:21 +0000 | [diff] [blame] | 170 | } |
| 171 | } |
Josh Coalson | bb7f6b9 | 2000-12-10 04:09:52 +0000 | [diff] [blame] | 172 | } |
Josh Coalson | 4e6b3ac | 2001-07-06 00:37:57 +0000 | [diff] [blame] | 173 | else { /* (*shift < 0) */ |
| 174 | const int nshift = -(*shift); |
Josh Coalson | 3dbbf94 | 2001-07-12 21:27:40 +0000 | [diff] [blame] | 175 | #ifdef DEBUG |
Josh Coalson | 4e6b3ac | 2001-07-06 00:37:57 +0000 | [diff] [blame] | 176 | fprintf(stderr,"FLAC__lpc_quantize_coefficients: negative shift = %d\n", *shift); |
Josh Coalson | 3dbbf94 | 2001-07-12 21:27:40 +0000 | [diff] [blame] | 177 | #endif |
Josh Coalson | 4e6b3ac | 2001-07-06 00:37:57 +0000 | [diff] [blame] | 178 | for(i = 0; i < order; i++) { |
| 179 | qlp_coeff[i] = (FLAC__int32)floor((double)lp_coeff[i] / (double)(1 << nshift)); |
| 180 | |
Josh Coalson | 8c6f90f | 2001-07-19 17:07:13 +0000 | [diff] [blame] | 181 | /* double-check the result */ |
Josh Coalson | 4e6b3ac | 2001-07-06 00:37:57 +0000 | [diff] [blame] | 182 | if(qlp_coeff[i] > qmax || qlp_coeff[i] < qmin) { |
| 183 | #ifdef FLAC__OVERFLOW_DETECT |
| 184 | fprintf(stderr,"FLAC__lpc_quantize_coefficients: compensating for overflow, qlp_coeff[%u]=%d, lp_coeff[%u]=%f, cmax=%f, precision=%u, shift=%d, q=%f, f(q)=%f\n", i, qlp_coeff[i], i, lp_coeff[i], cmax, precision, *shift, (double)lp_coeff[i] / (double)(1 << nshift), floor((double)lp_coeff[i] / (double)(1 << nshift))); |
| 185 | #endif |
| 186 | cmax *= 2.0; |
| 187 | goto redo_it; |
| 188 | } |
| 189 | } |
| 190 | } |
| 191 | |
Josh Coalson | bb7f6b9 | 2000-12-10 04:09:52 +0000 | [diff] [blame] | 192 | return 0; |
| 193 | } |
| 194 | |
Josh Coalson | 77e3f31 | 2001-06-23 03:03:24 +0000 | [diff] [blame] | 195 | void 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[]) |
Josh Coalson | bb7f6b9 | 2000-12-10 04:09:52 +0000 | [diff] [blame] | 196 | { |
Josh Coalson | bb6712e | 2001-04-24 22:54:07 +0000 | [diff] [blame] | 197 | #ifdef FLAC__OVERFLOW_DETECT |
Josh Coalson | 77e3f31 | 2001-06-23 03:03:24 +0000 | [diff] [blame] | 198 | FLAC__int64 sumo; |
Josh Coalson | bb7f6b9 | 2000-12-10 04:09:52 +0000 | [diff] [blame] | 199 | #endif |
| 200 | unsigned i, j; |
Josh Coalson | 77e3f31 | 2001-06-23 03:03:24 +0000 | [diff] [blame] | 201 | FLAC__int32 sum; |
| 202 | const FLAC__int32 *history; |
Josh Coalson | bb7f6b9 | 2000-12-10 04:09:52 +0000 | [diff] [blame] | 203 | |
Josh Coalson | bb6712e | 2001-04-24 22:54:07 +0000 | [diff] [blame] | 204 | #ifdef FLAC__OVERFLOW_DETECT_VERBOSE |
Josh Coalson | bb7f6b9 | 2000-12-10 04:09:52 +0000 | [diff] [blame] | 205 | fprintf(stderr,"FLAC__lpc_compute_residual_from_qlp_coefficients: data_len=%d, order=%u, lpq=%d",data_len,order,lp_quantization); |
| 206 | for(i=0;i<order;i++) |
| 207 | fprintf(stderr,", q[%u]=%d",i,qlp_coeff[i]); |
| 208 | fprintf(stderr,"\n"); |
| 209 | #endif |
Josh Coalson | 1b68982 | 2001-05-31 20:11:02 +0000 | [diff] [blame] | 210 | FLAC__ASSERT(order > 0); |
Josh Coalson | bb7f6b9 | 2000-12-10 04:09:52 +0000 | [diff] [blame] | 211 | |
| 212 | for(i = 0; i < data_len; i++) { |
Josh Coalson | bb6712e | 2001-04-24 22:54:07 +0000 | [diff] [blame] | 213 | #ifdef FLAC__OVERFLOW_DETECT |
Josh Coalson | bb7f6b9 | 2000-12-10 04:09:52 +0000 | [diff] [blame] | 214 | sumo = 0; |
| 215 | #endif |
| 216 | sum = 0; |
| 217 | history = data; |
| 218 | for(j = 0; j < order; j++) { |
| 219 | sum += qlp_coeff[j] * (*(--history)); |
Josh Coalson | bb6712e | 2001-04-24 22:54:07 +0000 | [diff] [blame] | 220 | #ifdef FLAC__OVERFLOW_DETECT |
Josh Coalson | 77e3f31 | 2001-06-23 03:03:24 +0000 | [diff] [blame] | 221 | sumo += (FLAC__int64)qlp_coeff[j] * (FLAC__int64)(*history); |
Josh Coalson | 40333b1 | 2001-11-13 21:37:04 +0000 | [diff] [blame] | 222 | #if defined _MSC_VER || defined __MINGW32__ /* don't know how to do 64-bit literals in VC++ */ |
Josh Coalson | d37acf4 | 2001-07-09 18:22:46 +0000 | [diff] [blame] | 223 | if(sumo < 0) sumo = -sumo; |
| 224 | if(sumo > 2147483647) |
| 225 | #else |
| 226 | if(sumo > 2147483647ll || sumo < -2147483648ll) |
| 227 | #endif |
| 228 | { |
Josh Coalson | b35bebd | 2001-07-03 04:37:18 +0000 | [diff] [blame] | 229 | 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,sumo); |
Josh Coalson | 9f77a19 | 2001-02-08 00:26:45 +0000 | [diff] [blame] | 230 | } |
Josh Coalson | bb7f6b9 | 2000-12-10 04:09:52 +0000 | [diff] [blame] | 231 | #endif |
| 232 | } |
| 233 | *(residual++) = *(data++) - (sum >> lp_quantization); |
| 234 | } |
| 235 | |
Josh Coalson | 9f77a19 | 2001-02-08 00:26:45 +0000 | [diff] [blame] | 236 | /* Here's a slower but clearer version: |
Josh Coalson | bb7f6b9 | 2000-12-10 04:09:52 +0000 | [diff] [blame] | 237 | for(i = 0; i < data_len; i++) { |
| 238 | sum = 0; |
Josh Coalson | bb7f6b9 | 2000-12-10 04:09:52 +0000 | [diff] [blame] | 239 | for(j = 0; j < order; j++) |
Josh Coalson | 9f77a19 | 2001-02-08 00:26:45 +0000 | [diff] [blame] | 240 | sum += qlp_coeff[j] * data[i-j-1]; |
Josh Coalson | bb7f6b9 | 2000-12-10 04:09:52 +0000 | [diff] [blame] | 241 | residual[i] = data[i] - (sum >> lp_quantization); |
| 242 | } |
| 243 | */ |
| 244 | } |
| 245 | |
Josh Coalson | 77e3f31 | 2001-06-23 03:03:24 +0000 | [diff] [blame] | 246 | void FLAC__lpc_restore_signal(const FLAC__int32 residual[], unsigned data_len, const FLAC__int32 qlp_coeff[], unsigned order, int lp_quantization, FLAC__int32 data[]) |
Josh Coalson | bb7f6b9 | 2000-12-10 04:09:52 +0000 | [diff] [blame] | 247 | { |
Josh Coalson | bb6712e | 2001-04-24 22:54:07 +0000 | [diff] [blame] | 248 | #ifdef FLAC__OVERFLOW_DETECT |
Josh Coalson | 77e3f31 | 2001-06-23 03:03:24 +0000 | [diff] [blame] | 249 | FLAC__int64 sumo; |
Josh Coalson | bb7f6b9 | 2000-12-10 04:09:52 +0000 | [diff] [blame] | 250 | #endif |
| 251 | unsigned i, j; |
Josh Coalson | 77e3f31 | 2001-06-23 03:03:24 +0000 | [diff] [blame] | 252 | FLAC__int32 sum; |
| 253 | const FLAC__int32 *history; |
Josh Coalson | bb7f6b9 | 2000-12-10 04:09:52 +0000 | [diff] [blame] | 254 | |
Josh Coalson | bb6712e | 2001-04-24 22:54:07 +0000 | [diff] [blame] | 255 | #ifdef FLAC__OVERFLOW_DETECT_VERBOSE |
Josh Coalson | bb7f6b9 | 2000-12-10 04:09:52 +0000 | [diff] [blame] | 256 | fprintf(stderr,"FLAC__lpc_restore_signal: data_len=%d, order=%u, lpq=%d",data_len,order,lp_quantization); |
| 257 | for(i=0;i<order;i++) |
| 258 | fprintf(stderr,", q[%u]=%d",i,qlp_coeff[i]); |
| 259 | fprintf(stderr,"\n"); |
| 260 | #endif |
Josh Coalson | 1b68982 | 2001-05-31 20:11:02 +0000 | [diff] [blame] | 261 | FLAC__ASSERT(order > 0); |
Josh Coalson | bb7f6b9 | 2000-12-10 04:09:52 +0000 | [diff] [blame] | 262 | |
| 263 | for(i = 0; i < data_len; i++) { |
Josh Coalson | bb6712e | 2001-04-24 22:54:07 +0000 | [diff] [blame] | 264 | #ifdef FLAC__OVERFLOW_DETECT |
Josh Coalson | bb7f6b9 | 2000-12-10 04:09:52 +0000 | [diff] [blame] | 265 | sumo = 0; |
| 266 | #endif |
| 267 | sum = 0; |
Josh Coalson | 9f77a19 | 2001-02-08 00:26:45 +0000 | [diff] [blame] | 268 | history = data; |
Josh Coalson | bb7f6b9 | 2000-12-10 04:09:52 +0000 | [diff] [blame] | 269 | for(j = 0; j < order; j++) { |
| 270 | sum += qlp_coeff[j] * (*(--history)); |
Josh Coalson | bb6712e | 2001-04-24 22:54:07 +0000 | [diff] [blame] | 271 | #ifdef FLAC__OVERFLOW_DETECT |
Josh Coalson | 77e3f31 | 2001-06-23 03:03:24 +0000 | [diff] [blame] | 272 | sumo += (FLAC__int64)qlp_coeff[j] * (FLAC__int64)(*history); |
Josh Coalson | 40333b1 | 2001-11-13 21:37:04 +0000 | [diff] [blame] | 273 | #if defined _MSC_VER || defined __MINGW32__ /* don't know how to do 64-bit literals in VC++ */ |
Josh Coalson | d37acf4 | 2001-07-09 18:22:46 +0000 | [diff] [blame] | 274 | if(sumo < 0) sumo = -sumo; |
| 275 | if(sumo > 2147483647) |
| 276 | #else |
| 277 | if(sumo > 2147483647ll || sumo < -2147483648ll) |
| 278 | #endif |
| 279 | { |
Josh Coalson | b35bebd | 2001-07-03 04:37:18 +0000 | [diff] [blame] | 280 | fprintf(stderr,"FLAC__lpc_restore_signal: OVERFLOW, i=%u, j=%u, c=%d, d=%d, sumo=%lld\n",i,j,qlp_coeff[j],*history,sumo); |
Josh Coalson | 9f77a19 | 2001-02-08 00:26:45 +0000 | [diff] [blame] | 281 | } |
Josh Coalson | bb7f6b9 | 2000-12-10 04:09:52 +0000 | [diff] [blame] | 282 | #endif |
| 283 | } |
Josh Coalson | 9f77a19 | 2001-02-08 00:26:45 +0000 | [diff] [blame] | 284 | *(data++) = *(residual++) + (sum >> lp_quantization); |
| 285 | } |
| 286 | |
| 287 | /* Here's a slower but clearer version: |
| 288 | for(i = 0; i < data_len; i++) { |
| 289 | sum = 0; |
| 290 | for(j = 0; j < order; j++) |
| 291 | sum += qlp_coeff[j] * data[i-j-1]; |
Josh Coalson | bb7f6b9 | 2000-12-10 04:09:52 +0000 | [diff] [blame] | 292 | data[i] = residual[i] + (sum >> lp_quantization); |
| 293 | } |
Josh Coalson | 9f77a19 | 2001-02-08 00:26:45 +0000 | [diff] [blame] | 294 | */ |
Josh Coalson | bb7f6b9 | 2000-12-10 04:09:52 +0000 | [diff] [blame] | 295 | } |
| 296 | |
Josh Coalson | 77e3f31 | 2001-06-23 03:03:24 +0000 | [diff] [blame] | 297 | FLAC__real FLAC__lpc_compute_expected_bits_per_residual_sample(FLAC__real lpc_error, unsigned total_samples) |
Josh Coalson | bb7f6b9 | 2000-12-10 04:09:52 +0000 | [diff] [blame] | 298 | { |
Josh Coalson | b35bebd | 2001-07-03 04:37:18 +0000 | [diff] [blame] | 299 | double error_scale; |
Josh Coalson | bb7f6b9 | 2000-12-10 04:09:52 +0000 | [diff] [blame] | 300 | |
Josh Coalson | 1b68982 | 2001-05-31 20:11:02 +0000 | [diff] [blame] | 301 | FLAC__ASSERT(total_samples > 0); |
Josh Coalson | bb7f6b9 | 2000-12-10 04:09:52 +0000 | [diff] [blame] | 302 | |
Josh Coalson | 77e3f31 | 2001-06-23 03:03:24 +0000 | [diff] [blame] | 303 | error_scale = 0.5 * M_LN2 * M_LN2 / (FLAC__real)total_samples; |
Josh Coalson | bb7f6b9 | 2000-12-10 04:09:52 +0000 | [diff] [blame] | 304 | |
Josh Coalson | 9f77a19 | 2001-02-08 00:26:45 +0000 | [diff] [blame] | 305 | if(lpc_error > 0.0) { |
Josh Coalson | b35bebd | 2001-07-03 04:37:18 +0000 | [diff] [blame] | 306 | FLAC__real bps = (FLAC__real)((double)0.5 * log(error_scale * lpc_error) / M_LN2); |
Josh Coalson | 9f77a19 | 2001-02-08 00:26:45 +0000 | [diff] [blame] | 307 | if(bps >= 0.0) |
| 308 | return bps; |
| 309 | else |
| 310 | return 0.0; |
| 311 | } |
Josh Coalson | a1b53c4 | 2001-05-24 19:27:08 +0000 | [diff] [blame] | 312 | else if(lpc_error < 0.0) { /* error should not be negative but can happen due to inadequate float resolution */ |
Josh Coalson | b35bebd | 2001-07-03 04:37:18 +0000 | [diff] [blame] | 313 | return (FLAC__real)1e32; |
Josh Coalson | a1b53c4 | 2001-05-24 19:27:08 +0000 | [diff] [blame] | 314 | } |
| 315 | else { |
| 316 | return 0.0; |
| 317 | } |
| 318 | } |
| 319 | |
Josh Coalson | b35bebd | 2001-07-03 04:37:18 +0000 | [diff] [blame] | 320 | FLAC__real FLAC__lpc_compute_expected_bits_per_residual_sample_with_error_scale(FLAC__real lpc_error, double error_scale) |
Josh Coalson | a1b53c4 | 2001-05-24 19:27:08 +0000 | [diff] [blame] | 321 | { |
| 322 | if(lpc_error > 0.0) { |
Josh Coalson | b35bebd | 2001-07-03 04:37:18 +0000 | [diff] [blame] | 323 | FLAC__real bps = (FLAC__real)((double)0.5 * log(error_scale * lpc_error) / M_LN2); |
Josh Coalson | a1b53c4 | 2001-05-24 19:27:08 +0000 | [diff] [blame] | 324 | if(bps >= 0.0) |
| 325 | return bps; |
| 326 | else |
| 327 | return 0.0; |
| 328 | } |
| 329 | else if(lpc_error < 0.0) { /* error should not be negative but can happen due to inadequate float resolution */ |
Josh Coalson | b35bebd | 2001-07-03 04:37:18 +0000 | [diff] [blame] | 330 | return (FLAC__real)1e32; |
Josh Coalson | a1b53c4 | 2001-05-24 19:27:08 +0000 | [diff] [blame] | 331 | } |
Josh Coalson | 9f77a19 | 2001-02-08 00:26:45 +0000 | [diff] [blame] | 332 | else { |
Josh Coalson | bb7f6b9 | 2000-12-10 04:09:52 +0000 | [diff] [blame] | 333 | return 0.0; |
Josh Coalson | 9f77a19 | 2001-02-08 00:26:45 +0000 | [diff] [blame] | 334 | } |
Josh Coalson | bb7f6b9 | 2000-12-10 04:09:52 +0000 | [diff] [blame] | 335 | } |
| 336 | |
Josh Coalson | 77e3f31 | 2001-06-23 03:03:24 +0000 | [diff] [blame] | 337 | unsigned FLAC__lpc_compute_best_order(const FLAC__real lpc_error[], unsigned max_order, unsigned total_samples, unsigned bits_per_signal_sample) |
Josh Coalson | bb7f6b9 | 2000-12-10 04:09:52 +0000 | [diff] [blame] | 338 | { |
| 339 | unsigned order, best_order; |
Josh Coalson | b35bebd | 2001-07-03 04:37:18 +0000 | [diff] [blame] | 340 | FLAC__real best_bits, tmp_bits; |
| 341 | double error_scale; |
Josh Coalson | bb7f6b9 | 2000-12-10 04:09:52 +0000 | [diff] [blame] | 342 | |
Josh Coalson | 1b68982 | 2001-05-31 20:11:02 +0000 | [diff] [blame] | 343 | FLAC__ASSERT(max_order > 0); |
| 344 | FLAC__ASSERT(total_samples > 0); |
Josh Coalson | a1b53c4 | 2001-05-24 19:27:08 +0000 | [diff] [blame] | 345 | |
Josh Coalson | 77e3f31 | 2001-06-23 03:03:24 +0000 | [diff] [blame] | 346 | error_scale = 0.5 * M_LN2 * M_LN2 / (FLAC__real)total_samples; |
Josh Coalson | bb7f6b9 | 2000-12-10 04:09:52 +0000 | [diff] [blame] | 347 | |
| 348 | best_order = 0; |
Josh Coalson | 77e3f31 | 2001-06-23 03:03:24 +0000 | [diff] [blame] | 349 | best_bits = FLAC__lpc_compute_expected_bits_per_residual_sample_with_error_scale(lpc_error[0], error_scale) * (FLAC__real)total_samples; |
Josh Coalson | bb7f6b9 | 2000-12-10 04:09:52 +0000 | [diff] [blame] | 350 | |
| 351 | for(order = 1; order < max_order; order++) { |
Josh Coalson | 77e3f31 | 2001-06-23 03:03:24 +0000 | [diff] [blame] | 352 | tmp_bits = FLAC__lpc_compute_expected_bits_per_residual_sample_with_error_scale(lpc_error[order], error_scale) * (FLAC__real)(total_samples - order) + (FLAC__real)(order * bits_per_signal_sample); |
Josh Coalson | bb7f6b9 | 2000-12-10 04:09:52 +0000 | [diff] [blame] | 353 | if(tmp_bits < best_bits) { |
| 354 | best_order = order; |
| 355 | best_bits = tmp_bits; |
| 356 | } |
| 357 | } |
| 358 | |
| 359 | return best_order+1; /* +1 since index of lpc_error[] is order-1 */ |
| 360 | } |