ffmpeg / libavcodec / lpc.c @ 81fc8a63
History  View  Annotate  Download (5.75 KB)
1 
/**


2 
* LPC utility code

3 
* Copyright (c) 2006 Justin Ruggles <jruggle@earthlink.net>

4 
*

5 
* This file is part of FFmpeg.

6 
*

7 
* FFmpeg is free software; you can redistribute it and/or

8 
* modify it under the terms of the GNU Lesser General Public

9 
* License as published by the Free Software Foundation; either

10 
* version 2.1 of the License, or (at your option) any later version.

11 
*

12 
* FFmpeg is distributed in the hope that it will be useful,

13 
* but WITHOUT ANY WARRANTY; without even the implied warranty of

14 
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU

15 
* Lesser General Public License for more details.

16 
*

17 
* You should have received a copy of the GNU Lesser General Public

18 
* License along with FFmpeg; if not, write to the Free Software

19 
* Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 021101301 USA

20 
*/

21  
22 
#include "libavutil/lls.h" 
23 
#include "dsputil.h" 
24 
#include "lpc.h" 
25  
26  
27 
/**

28 
* LevinsonDurbin recursion.

29 
* Produces LPC coefficients from autocorrelation data.

30 
*/

31 
static void compute_lpc_coefs(const double *autoc, int max_order, 
32 
double lpc[][MAX_LPC_ORDER], double *ref) 
33 
{ 
34 
int i, j, i2;

35 
double r, err, tmp;

36 
double lpc_tmp[MAX_LPC_ORDER];

37  
38 
for(i=0; i<max_order; i++) lpc_tmp[i] = 0; 
39 
err = autoc[0];

40  
41 
for(i=0; i<max_order; i++) { 
42 
r = autoc[i+1];

43 
for(j=0; j<i; j++) { 
44 
r = lpc_tmp[j] * autoc[ij]; 
45 
} 
46 
r /= err; 
47 
ref[i] = fabs(r); 
48  
49 
err *= 1.0  (r * r); 
50  
51 
i2 = (i >> 1);

52 
lpc_tmp[i] = r; 
53 
for(j=0; j<i2; j++) { 
54 
tmp = lpc_tmp[j]; 
55 
lpc_tmp[j] += r * lpc_tmp[i1j];

56 
lpc_tmp[i1j] += r * tmp;

57 
} 
58 
if(i & 1) { 
59 
lpc_tmp[j] += lpc_tmp[j] * r; 
60 
} 
61  
62 
for(j=0; j<=i; j++) { 
63 
lpc[i][j] = lpc_tmp[j]; 
64 
} 
65 
} 
66 
} 
67  
68 
/**

69 
* Quantize LPC coefficients

70 
*/

71 
static void quantize_lpc_coefs(double *lpc_in, int order, int precision, 
72 
int32_t *lpc_out, int *shift, int max_shift, int zero_shift) 
73 
{ 
74 
int i;

75 
double cmax, error;

76 
int32_t qmax; 
77 
int sh;

78  
79 
/* define maximum levels */

80 
qmax = (1 << (precision  1))  1; 
81  
82 
/* find maximum coefficient value */

83 
cmax = 0.0; 
84 
for(i=0; i<order; i++) { 
85 
cmax= FFMAX(cmax, fabs(lpc_in[i])); 
86 
} 
87  
88 
/* if maximum value quantizes to zero, return all zeros */

89 
if(cmax * (1 << max_shift) < 1.0) { 
90 
*shift = zero_shift; 
91 
memset(lpc_out, 0, sizeof(int32_t) * order); 
92 
return;

93 
} 
94  
95 
/* calculate level shift which scales max coeff to available bits */

96 
sh = max_shift; 
97 
while((cmax * (1 << sh) > qmax) && (sh > 0)) { 
98 
sh; 
99 
} 
100  
101 
/* since negative shift values are unsupported in decoder, scale down

102 
coefficients instead */

103 
if(sh == 0 && cmax > qmax) { 
104 
double scale = ((double)qmax) / cmax; 
105 
for(i=0; i<order; i++) { 
106 
lpc_in[i] *= scale; 
107 
} 
108 
} 
109  
110 
/* output quantized coefficients and level shift */

111 
error=0;

112 
for(i=0; i<order; i++) { 
113 
error += lpc_in[i] * (1 << sh);

114 
lpc_out[i] = av_clip(lrintf(error), qmax, qmax); 
115 
error = lpc_out[i]; 
116 
} 
117 
*shift = sh; 
118 
} 
119  
120 
static int estimate_best_order(double *ref, int min_order, int max_order) 
121 
{ 
122 
int i, est;

123  
124 
est = min_order; 
125 
for(i=max_order1; i>=min_order1; i) { 
126 
if(ref[i] > 0.10) { 
127 
est = i+1;

128 
break;

129 
} 
130 
} 
131 
return est;

132 
} 
133  
134 
/**

135 
* Calculate LPC coefficients for multiple orders

136 
*/

137 
int ff_lpc_calc_coefs(DSPContext *s,

138 
const int32_t *samples, int blocksize, int min_order, 
139 
int max_order, int precision, 
140 
int32_t coefs[][MAX_LPC_ORDER], int *shift, int use_lpc, 
141 
int omethod, int max_shift, int zero_shift) 
142 
{ 
143 
double autoc[MAX_LPC_ORDER+1]; 
144 
double ref[MAX_LPC_ORDER];

145 
double lpc[MAX_LPC_ORDER][MAX_LPC_ORDER];

146 
int i, j, pass;

147 
int opt_order;

148  
149 
assert(max_order >= MIN_LPC_ORDER && max_order <= MAX_LPC_ORDER); 
150  
151 
if(use_lpc == 1){ 
152 
s>flac_compute_autocorr(samples, blocksize, max_order, autoc); 
153  
154 
compute_lpc_coefs(autoc, max_order, lpc, ref); 
155 
}else{

156 
LLSModel m[2];

157 
double var[MAX_LPC_ORDER+1], weight; 
158  
159 
for(pass=0; pass<use_lpc1; pass++){ 
160 
av_init_lls(&m[pass&1], max_order);

161  
162 
weight=0;

163 
for(i=max_order; i<blocksize; i++){

164 
for(j=0; j<=max_order; j++) 
165 
var[j]= samples[ij]; 
166  
167 
if(pass){

168 
double eval, inv, rinv;

169 
eval= av_evaluate_lls(&m[(pass1)&1], var+1, max_order1); 
170 
eval= (512>>pass) + fabs(eval  var[0]); 
171 
inv = 1/eval;

172 
rinv = sqrt(inv); 
173 
for(j=0; j<=max_order; j++) 
174 
var[j] *= rinv; 
175 
weight += inv; 
176 
}else

177 
weight++; 
178  
179 
av_update_lls(&m[pass&1], var, 1.0); 
180 
} 
181 
av_solve_lls(&m[pass&1], 0.001, 0); 
182 
} 
183  
184 
for(i=0; i<max_order; i++){ 
185 
for(j=0; j<max_order; j++) 
186 
lpc[i][j]= m[(pass1)&1].coeff[i][j]; 
187 
ref[i]= sqrt(m[(pass1)&1].variance[i] / weight) * (blocksize  max_order) / 4000; 
188 
} 
189 
for(i=max_order1; i>0; i) 
190 
ref[i] = ref[i1]  ref[i];

191 
} 
192 
opt_order = max_order; 
193  
194 
if(omethod == ORDER_METHOD_EST) {

195 
opt_order = estimate_best_order(ref, min_order, max_order); 
196 
i = opt_order1;

197 
quantize_lpc_coefs(lpc[i], i+1, precision, coefs[i], &shift[i], max_shift, zero_shift);

198 
} else {

199 
for(i=min_order1; i<max_order; i++) { 
200 
quantize_lpc_coefs(lpc[i], i+1, precision, coefs[i], &shift[i], max_shift, zero_shift);

201 
} 
202 
} 
203  
204 
return opt_order;

205 
} 