From 408ec4e2a6c3fb40e14ac4f0fb2fb9e40ff3e6a3 Mon Sep 17 00:00:00 2001 From: Michael Niedermayer Date: Sat, 15 Jul 2006 23:43:38 +0000 Subject: calculate all coefficients for several orders during cholesky factorization, the resulting coefficients are not strictly optimal though as there is a small difference in the autocorrelation matrixes which is ignored for the smaller orders Originally committed as revision 5758 to svn://svn.ffmpeg.org/ffmpeg/trunk --- libavutil/lls.c | 64 ++++++++++++++++++++++++++++++++------------------------- 1 file changed, 36 insertions(+), 28 deletions(-) (limited to 'libavutil/lls.c') diff --git a/libavutil/lls.c b/libavutil/lls.c index 0556d8c80f..50a5003763 100644 --- a/libavutil/lls.c +++ b/libavutil/lls.c @@ -49,12 +49,11 @@ void av_update_lls(LLSModel *m, double *var, double decay){ } } -double av_solve_lls(LLSModel *m, double threshold){ +void av_solve_lls(LLSModel *m, double threshold, int min_order){ int i,j,k; double (*factor)[MAX_VARS+1]= &m->covariance[1][0]; double (*covar )[MAX_VARS+1]= &m->covariance[1][1]; double *covar_y = m->covariance[0]; - double variance; int count= m->indep_count; for(i=0; i=0; k--) - sum -= factor[i][k]*m->coeff[k]; - m->coeff[i]= sum / factor[i][i]; + sum -= factor[i][k]*m->coeff[0][k]; + m->coeff[0][i]= sum / factor[i][i]; } - for(i=count-1; i>=0; i--){ - double sum= m->coeff[i]; - for(k=i+1; kcoeff[k]; - m->coeff[i]= sum / factor[i][i]; - } + for(j=count-1; j>=min_order; j--){ + for(i=j; i>=0; i--){ + double sum= m->coeff[0][i]; + for(k=i+1; k<=j; k++) + sum -= factor[k][i]*m->coeff[j][k]; + m->coeff[j][i]= sum / factor[i][i]; + } - variance= covar_y[0]; - for(i=0; icoeff[i]*covar[i][i] - 2*covar_y[i+1]; - for(j=0; jcoeff[j]*covar[j][i]; - variance += m->coeff[i]*sum; + m->variance[j]= covar_y[0]; + for(i=0; i<=j; i++){ + double sum= m->coeff[j][i]*covar[i][i] - 2*covar_y[i+1]; + for(k=0; kcoeff[j][k]*covar[k][i]; + m->variance[j] += m->coeff[j][i]*sum; + } } - return variance; } -double av_evaluate_lls(LLSModel *m, double *param){ +double av_evaluate_lls(LLSModel *m, double *param, int order){ int i; double out= 0; - for(i=0; iindep_count; i++) - out+= param[i]*m->coeff[i]; + for(i=0; i<=order; i++) + out+= param[i]*m->coeff[order][i]; return out; } @@ -113,27 +113,35 @@ double av_evaluate_lls(LLSModel *m, double *param){ int main(){ LLSModel m; - int i; + int i, order; av_init_lls(&m, 3); for(i=0; i<100; i++){ double var[4]; double eval, variance; +#if 0 var[1] = rand() / (double)RAND_MAX; var[2] = rand() / (double)RAND_MAX; var[3] = rand() / (double)RAND_MAX; - var[2]= var[1] + var[3]; + var[2]= var[1] + var[3]/2; var[0] = var[1] + var[2] + var[3] + var[1]*var[2]/100; - - eval= av_evaluate_lls(&m, var+1); +#else + var[0] = (rand() / (double)RAND_MAX - 0.5)*2; + var[1] = var[0] + rand() / (double)RAND_MAX - 0.5; + var[2] = var[1] + rand() / (double)RAND_MAX - 0.5; + var[3] = var[2] + rand() / (double)RAND_MAX - 0.5; +#endif av_update_lls(&m, var, 0.99); - variance= av_solve_lls(&m, 0.001); - av_log(NULL, AV_LOG_DEBUG, "real:%f pred:%f var:%f coeffs:%f %f %f\n", - var[0], eval, sqrt(variance / (i+1)), - m.coeff[0], m.coeff[1], m.coeff[2]); + av_solve_lls(&m, 0.001, 0); + for(order=0; order<3; order++){ + eval= av_evaluate_lls(&m, var+1, order); + av_log(NULL, AV_LOG_DEBUG, "real:%f order:%d pred:%f var:%f coeffs:%f %f %f\n", + var[0], order, eval, sqrt(m.variance[order] / (i+1)), + m.coeff[order][0], m.coeff[order][1], m.coeff[order][2]); + } } return 0; } -- cgit v1.2.3