00001
00022 #include "libavutil/lls.h"
00023 #include "dsputil.h"
00024
00025 #define LPC_USE_DOUBLE
00026 #include "lpc.h"
00027
00028
00032 static void quantize_lpc_coefs(double *lpc_in, int order, int precision,
00033 int32_t *lpc_out, int *shift, int max_shift, int zero_shift)
00034 {
00035 int i;
00036 double cmax, error;
00037 int32_t qmax;
00038 int sh;
00039
00040
00041 qmax = (1 << (precision - 1)) - 1;
00042
00043
00044 cmax = 0.0;
00045 for(i=0; i<order; i++) {
00046 cmax= FFMAX(cmax, fabs(lpc_in[i]));
00047 }
00048
00049
00050 if(cmax * (1 << max_shift) < 1.0) {
00051 *shift = zero_shift;
00052 memset(lpc_out, 0, sizeof(int32_t) * order);
00053 return;
00054 }
00055
00056
00057 sh = max_shift;
00058 while((cmax * (1 << sh) > qmax) && (sh > 0)) {
00059 sh--;
00060 }
00061
00062
00063
00064 if(sh == 0 && cmax > qmax) {
00065 double scale = ((double)qmax) / cmax;
00066 for(i=0; i<order; i++) {
00067 lpc_in[i] *= scale;
00068 }
00069 }
00070
00071
00072 error=0;
00073 for(i=0; i<order; i++) {
00074 error -= lpc_in[i] * (1 << sh);
00075 lpc_out[i] = av_clip(lrintf(error), -qmax, qmax);
00076 error -= lpc_out[i];
00077 }
00078 *shift = sh;
00079 }
00080
00081 static int estimate_best_order(double *ref, int min_order, int max_order)
00082 {
00083 int i, est;
00084
00085 est = min_order;
00086 for(i=max_order-1; i>=min_order-1; i--) {
00087 if(ref[i] > 0.10) {
00088 est = i+1;
00089 break;
00090 }
00091 }
00092 return est;
00093 }
00094
00103 int ff_lpc_calc_coefs(DSPContext *s,
00104 const int32_t *samples, int blocksize, int min_order,
00105 int max_order, int precision,
00106 int32_t coefs[][MAX_LPC_ORDER], int *shift, int use_lpc,
00107 int omethod, int max_shift, int zero_shift)
00108 {
00109 double autoc[MAX_LPC_ORDER+1];
00110 double ref[MAX_LPC_ORDER];
00111 double lpc[MAX_LPC_ORDER][MAX_LPC_ORDER];
00112 int i, j, pass;
00113 int opt_order;
00114
00115 assert(max_order >= MIN_LPC_ORDER && max_order <= MAX_LPC_ORDER && use_lpc > 0);
00116
00117 if(use_lpc == 1){
00118 s->flac_compute_autocorr(samples, blocksize, max_order, autoc);
00119
00120 compute_lpc_coefs(autoc, max_order, &lpc[0][0], MAX_LPC_ORDER, 0, 1);
00121
00122 for(i=0; i<max_order; i++)
00123 ref[i] = fabs(lpc[i][i]);
00124 }else{
00125 LLSModel m[2];
00126 double var[MAX_LPC_ORDER+1], av_uninit(weight);
00127
00128 for(pass=0; pass<use_lpc-1; pass++){
00129 av_init_lls(&m[pass&1], max_order);
00130
00131 weight=0;
00132 for(i=max_order; i<blocksize; i++){
00133 for(j=0; j<=max_order; j++)
00134 var[j]= samples[i-j];
00135
00136 if(pass){
00137 double eval, inv, rinv;
00138 eval= av_evaluate_lls(&m[(pass-1)&1], var+1, max_order-1);
00139 eval= (512>>pass) + fabs(eval - var[0]);
00140 inv = 1/eval;
00141 rinv = sqrt(inv);
00142 for(j=0; j<=max_order; j++)
00143 var[j] *= rinv;
00144 weight += inv;
00145 }else
00146 weight++;
00147
00148 av_update_lls(&m[pass&1], var, 1.0);
00149 }
00150 av_solve_lls(&m[pass&1], 0.001, 0);
00151 }
00152
00153 for(i=0; i<max_order; i++){
00154 for(j=0; j<max_order; j++)
00155 lpc[i][j]=-m[(pass-1)&1].coeff[i][j];
00156 ref[i]= sqrt(m[(pass-1)&1].variance[i] / weight) * (blocksize - max_order) / 4000;
00157 }
00158 for(i=max_order-1; i>0; i--)
00159 ref[i] = ref[i-1] - ref[i];
00160 }
00161 opt_order = max_order;
00162
00163 if(omethod == ORDER_METHOD_EST) {
00164 opt_order = estimate_best_order(ref, min_order, max_order);
00165 i = opt_order-1;
00166 quantize_lpc_coefs(lpc[i], i+1, precision, coefs[i], &shift[i], max_shift, zero_shift);
00167 } else {
00168 for(i=min_order-1; i<max_order; i++) {
00169 quantize_lpc_coefs(lpc[i], i+1, precision, coefs[i], &shift[i], max_shift, zero_shift);
00170 }
00171 }
00172
00173 return opt_order;
00174 }