annotate lls.c @ 395:da2c1190902f libavutil

Fix typo. Patch by Jeremy Mordkoff: jlm zazzletech com
author benoit
date Fri, 19 Oct 2007 15:39:02 +0000
parents f13e5473611e
children f9a4c04ebb0e
Ignore whitespace changes - Everywhere: Within whitespace: At end of lines:
rev   line source
76
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
1 /*
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
2 * linear least squares model
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
3 *
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
4 * Copyright (c) 2006 Michael Niedermayer <michaelni@gmx.at>
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
5 *
116
d76a36742464 Change license headers to say 'FFmpeg' instead of 'this program/this library'
diego
parents: 90
diff changeset
6 * This file is part of FFmpeg.
d76a36742464 Change license headers to say 'FFmpeg' instead of 'this program/this library'
diego
parents: 90
diff changeset
7 *
d76a36742464 Change license headers to say 'FFmpeg' instead of 'this program/this library'
diego
parents: 90
diff changeset
8 * FFmpeg is free software; you can redistribute it and/or
76
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
9 * modify it under the terms of the GNU Lesser General Public
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
10 * License as published by the Free Software Foundation; either
116
d76a36742464 Change license headers to say 'FFmpeg' instead of 'this program/this library'
diego
parents: 90
diff changeset
11 * version 2.1 of the License, or (at your option) any later version.
76
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
12 *
116
d76a36742464 Change license headers to say 'FFmpeg' instead of 'this program/this library'
diego
parents: 90
diff changeset
13 * FFmpeg is distributed in the hope that it will be useful,
76
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
14 * but WITHOUT ANY WARRANTY; without even the implied warranty of
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
15 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
16 * Lesser General Public License for more details.
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
17 *
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
18 * You should have received a copy of the GNU Lesser General Public
116
d76a36742464 Change license headers to say 'FFmpeg' instead of 'this program/this library'
diego
parents: 90
diff changeset
19 * License along with FFmpeg; if not, write to the Free Software
358
f13e5473611e license header consistency cosmetics
diego
parents: 116
diff changeset
20 * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
76
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
21 */
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
22
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
23 /**
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
24 * @file lls.c
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
25 * linear least squares model
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
26 */
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
27
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
28 #include <math.h>
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
29 #include <string.h>
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
30
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
31 #include "lls.h"
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
32
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
33 #ifdef TEST
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
34 #define av_log(a,b,...) printf(__VA_ARGS__)
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
35 #endif
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
36
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
37 void av_init_lls(LLSModel *m, int indep_count){
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
38 memset(m, 0, sizeof(LLSModel));
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
39
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
40 m->indep_count= indep_count;
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
41 }
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
42
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
43 void av_update_lls(LLSModel *m, double *var, double decay){
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
44 int i,j;
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
45
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
46 for(i=0; i<=m->indep_count; i++){
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
47 for(j=i; j<=m->indep_count; j++){
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
48 m->covariance[i][j] *= decay;
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
49 m->covariance[i][j] += var[i]*var[j];
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
50 }
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
51 }
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
52 }
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
53
79
adbb5540fa47 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
michael
parents: 78
diff changeset
54 void av_solve_lls(LLSModel *m, double threshold, int min_order){
76
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
55 int i,j,k;
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
56 double (*factor)[MAX_VARS+1]= &m->covariance[1][0];
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
57 double (*covar )[MAX_VARS+1]= &m->covariance[1][1];
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
58 double *covar_y = m->covariance[0];
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
59 int count= m->indep_count;
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
60
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
61 for(i=0; i<count; i++){
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
62 for(j=i; j<count; j++){
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
63 double sum= covar[i][j];
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
64
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
65 for(k=i-1; k>=0; k--)
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
66 sum -= factor[i][k]*factor[j][k];
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
67
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
68 if(i==j){
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
69 if(sum < threshold)
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
70 sum= 1.0;
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
71 factor[i][i]= sqrt(sum);
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
72 }else
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
73 factor[j][i]= sum / factor[i][i];
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
74 }
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
75 }
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
76 for(i=0; i<count; i++){
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
77 double sum= covar_y[i+1];
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
78 for(k=i-1; k>=0; k--)
79
adbb5540fa47 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
michael
parents: 78
diff changeset
79 sum -= factor[i][k]*m->coeff[0][k];
adbb5540fa47 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
michael
parents: 78
diff changeset
80 m->coeff[0][i]= sum / factor[i][i];
76
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
81 }
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
82
79
adbb5540fa47 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
michael
parents: 78
diff changeset
83 for(j=count-1; j>=min_order; j--){
adbb5540fa47 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
michael
parents: 78
diff changeset
84 for(i=j; i>=0; i--){
adbb5540fa47 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
michael
parents: 78
diff changeset
85 double sum= m->coeff[0][i];
adbb5540fa47 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
michael
parents: 78
diff changeset
86 for(k=i+1; k<=j; k++)
adbb5540fa47 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
michael
parents: 78
diff changeset
87 sum -= factor[k][i]*m->coeff[j][k];
adbb5540fa47 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
michael
parents: 78
diff changeset
88 m->coeff[j][i]= sum / factor[i][i];
adbb5540fa47 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
michael
parents: 78
diff changeset
89 }
adbb5540fa47 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
michael
parents: 78
diff changeset
90
adbb5540fa47 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
michael
parents: 78
diff changeset
91 m->variance[j]= covar_y[0];
adbb5540fa47 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
michael
parents: 78
diff changeset
92 for(i=0; i<=j; i++){
adbb5540fa47 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
michael
parents: 78
diff changeset
93 double sum= m->coeff[j][i]*covar[i][i] - 2*covar_y[i+1];
adbb5540fa47 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
michael
parents: 78
diff changeset
94 for(k=0; k<i; k++)
adbb5540fa47 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
michael
parents: 78
diff changeset
95 sum += 2*m->coeff[j][k]*covar[k][i];
adbb5540fa47 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
michael
parents: 78
diff changeset
96 m->variance[j] += m->coeff[j][i]*sum;
adbb5540fa47 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
michael
parents: 78
diff changeset
97 }
76
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
98 }
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
99 }
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
100
79
adbb5540fa47 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
michael
parents: 78
diff changeset
101 double av_evaluate_lls(LLSModel *m, double *param, int order){
76
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
102 int i;
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
103 double out= 0;
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
104
79
adbb5540fa47 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
michael
parents: 78
diff changeset
105 for(i=0; i<=order; i++)
adbb5540fa47 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
michael
parents: 78
diff changeset
106 out+= param[i]*m->coeff[order][i];
76
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
107
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
108 return out;
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
109 }
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
110
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
111 #ifdef TEST
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
112
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
113 #include <stdlib.h>
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
114 #include <stdio.h>
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
115
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
116 int main(){
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
117 LLSModel m;
79
adbb5540fa47 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
michael
parents: 78
diff changeset
118 int i, order;
76
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
119
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
120 av_init_lls(&m, 3);
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
121
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
122 for(i=0; i<100; i++){
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
123 double var[4];
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
124 double eval, variance;
79
adbb5540fa47 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
michael
parents: 78
diff changeset
125 #if 0
76
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
126 var[1] = rand() / (double)RAND_MAX;
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
127 var[2] = rand() / (double)RAND_MAX;
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
128 var[3] = rand() / (double)RAND_MAX;
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
129
79
adbb5540fa47 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
michael
parents: 78
diff changeset
130 var[2]= var[1] + var[3]/2;
76
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
131
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
132 var[0] = var[1] + var[2] + var[3] + var[1]*var[2]/100;
79
adbb5540fa47 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
michael
parents: 78
diff changeset
133 #else
adbb5540fa47 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
michael
parents: 78
diff changeset
134 var[0] = (rand() / (double)RAND_MAX - 0.5)*2;
adbb5540fa47 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
michael
parents: 78
diff changeset
135 var[1] = var[0] + rand() / (double)RAND_MAX - 0.5;
adbb5540fa47 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
michael
parents: 78
diff changeset
136 var[2] = var[1] + rand() / (double)RAND_MAX - 0.5;
adbb5540fa47 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
michael
parents: 78
diff changeset
137 var[3] = var[2] + rand() / (double)RAND_MAX - 0.5;
adbb5540fa47 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
michael
parents: 78
diff changeset
138 #endif
76
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
139 av_update_lls(&m, var, 0.99);
79
adbb5540fa47 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
michael
parents: 78
diff changeset
140 av_solve_lls(&m, 0.001, 0);
adbb5540fa47 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
michael
parents: 78
diff changeset
141 for(order=0; order<3; order++){
adbb5540fa47 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
michael
parents: 78
diff changeset
142 eval= av_evaluate_lls(&m, var+1, order);
adbb5540fa47 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
michael
parents: 78
diff changeset
143 av_log(NULL, AV_LOG_DEBUG, "real:%f order:%d pred:%f var:%f coeffs:%f %f %f\n",
adbb5540fa47 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
michael
parents: 78
diff changeset
144 var[0], order, eval, sqrt(m.variance[order] / (i+1)),
adbb5540fa47 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
michael
parents: 78
diff changeset
145 m.coeff[order][0], m.coeff[order][1], m.coeff[order][2]);
adbb5540fa47 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
michael
parents: 78
diff changeset
146 }
76
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
147 }
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
148 return 0;
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
149 }
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
150
8c75234388b5 linear least squares solver using cholesky factorization
michael
parents:
diff changeset
151 #endif