## ffmpeg / libavutil / lls.c @ d71ad089

History | View | Annotate | Download (3.74 KB)

1 | 82ab5ad7 | Michael Niedermayer | ```
/*
``` |
---|---|---|---|

2 | ```
* linear least squares model
``` |
||

3 | ```
*
``` |
||

4 | ```
* Copyright (c) 2006 Michael Niedermayer <michaelni@gmx.at>
``` |
||

5 | ```
*
``` |
||

6 | b78e7197 | Diego Biurrun | ```
* This file is part of FFmpeg.
``` |

7 | ```
*
``` |
||

8 | ```
* FFmpeg is free software; you can redistribute it and/or
``` |
||

9 | 82ab5ad7 | Michael Niedermayer | ```
* modify it under the terms of the GNU Lesser General Public
``` |

10 | ```
* License as published by the Free Software Foundation; either
``` |
||

11 | b78e7197 | Diego Biurrun | ```
* version 2.1 of the License, or (at your option) any later version.
``` |

12 | 82ab5ad7 | Michael Niedermayer | ```
*
``` |

13 | b78e7197 | Diego Biurrun | ```
* FFmpeg is distributed in the hope that it will be useful,
``` |

14 | 82ab5ad7 | Michael Niedermayer | ```
* but WITHOUT ANY WARRANTY; without even the implied warranty of
``` |

15 | ```
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
``` |
||

16 | ```
* Lesser General Public License for more details.
``` |
||

17 | ```
*
``` |
||

18 | ```
* You should have received a copy of the GNU Lesser General Public
``` |
||

19 | b78e7197 | Diego Biurrun | ```
* License along with FFmpeg; if not, write to the Free Software
``` |

20 | e5a389a1 | Diego Biurrun | ```
* Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
``` |

21 | 82ab5ad7 | Michael Niedermayer | ```
*/
``` |

22 | |||

23 | ```
/**
``` |
||

24 | bad5537e | Diego Biurrun | ```
* @file libavutil/lls.c
``` |

25 | 82ab5ad7 | Michael Niedermayer | ```
* linear least squares model
``` |

26 | ```
*/
``` |
||

27 | |||

28 | #include <math.h> |
||

29 | #include <string.h> |
||

30 | |||

31 | #include "lls.h" |
||

32 | |||

33 | void av_init_lls(LLSModel *m, int indep_count){ |
||

34 | memset(m, 0, sizeof(LLSModel)); |
||

35 | |||

36 | m->indep_count= indep_count; |
||

37 | } |
||

38 | |||

39 | void av_update_lls(LLSModel *m, double *var, double decay){ |
||

40 | ```
int i,j;
``` |
||

41 | |||

42 | for(i=0; i<=m->indep_count; i++){ |
||

43 | ```
for(j=i; j<=m->indep_count; j++){
``` |
||

44 | m->covariance[i][j] *= decay; |
||

45 | m->covariance[i][j] += var[i]*var[j]; |
||

46 | } |
||

47 | } |
||

48 | } |
||

49 | |||

50 | 408ec4e2 | Michael Niedermayer | void av_solve_lls(LLSModel *m, double threshold, int min_order){ |

51 | 82ab5ad7 | Michael Niedermayer | ```
int i,j,k;
``` |

52 | b8fe8ab0 | Michael Niedermayer | double (*factor)[MAX_VARS+1]= (void*)&m->covariance[1][0]; |

53 | double (*covar )[MAX_VARS+1]= (void*)&m->covariance[1][1]; |
||

54 | 82ab5ad7 | Michael Niedermayer | double *covar_y = m->covariance[0]; |

55 | ```
int count= m->indep_count;
``` |
||

56 | |||

57 | for(i=0; i<count; i++){ |
||

58 | ```
for(j=i; j<count; j++){
``` |
||

59 | ```
double sum= covar[i][j];
``` |
||

60 | |||

61 | for(k=i-1; k>=0; k--) |
||

62 | sum -= factor[i][k]*factor[j][k]; |
||

63 | |||

64 | ```
if(i==j){
``` |
||

65 | ```
if(sum < threshold)
``` |
||

66 | sum= 1.0; |
||

67 | factor[i][i]= sqrt(sum); |
||

68 | ```
}else
``` |
||

69 | factor[j][i]= sum / factor[i][i]; |
||

70 | } |
||

71 | } |
||

72 | for(i=0; i<count; i++){ |
||

73 | double sum= covar_y[i+1]; |
||

74 | for(k=i-1; k>=0; k--) |
||

75 | 408ec4e2 | Michael Niedermayer | ```
sum -= factor[i][k]*m->coeff[0][k];
``` |

76 | ```
m->coeff[0][i]= sum / factor[i][i];
``` |
||

77 | 82ab5ad7 | Michael Niedermayer | } |

78 | |||

79 | 408ec4e2 | Michael Niedermayer | for(j=count-1; j>=min_order; j--){ |

80 | for(i=j; i>=0; i--){ |
||

81 | double sum= m->coeff[0][i]; |
||

82 | for(k=i+1; k<=j; k++) |
||

83 | sum -= factor[k][i]*m->coeff[j][k]; |
||

84 | m->coeff[j][i]= sum / factor[i][i]; |
||

85 | } |
||

86 | 82ab5ad7 | Michael Niedermayer | |

87 | 408ec4e2 | Michael Niedermayer | ```
m->variance[j]= covar_y[0];
``` |

88 | for(i=0; i<=j; i++){ |
||

89 | double sum= m->coeff[j][i]*covar[i][i] - 2*covar_y[i+1]; |
||

90 | for(k=0; k<i; k++) |
||

91 | ```
sum += 2*m->coeff[j][k]*covar[k][i];
``` |
||

92 | m->variance[j] += m->coeff[j][i]*sum; |
||

93 | } |
||

94 | 82ab5ad7 | Michael Niedermayer | } |

95 | } |
||

96 | |||

97 | 408ec4e2 | Michael Niedermayer | double av_evaluate_lls(LLSModel *m, double *param, int order){ |

98 | 82ab5ad7 | Michael Niedermayer | ```
int i;
``` |

99 | double out= 0; |
||

100 | |||

101 | 408ec4e2 | Michael Niedermayer | for(i=0; i<=order; i++) |

102 | out+= param[i]*m->coeff[order][i]; |
||

103 | 82ab5ad7 | Michael Niedermayer | |

104 | ```
return out;
``` |
||

105 | } |
||

106 | |||

107 | ```
#ifdef TEST
``` |
||

108 | |||

109 | #include <stdlib.h> |
||

110 | #include <stdio.h> |
||

111 | |||

112 | f8a80fd6 | Diego Biurrun | int main(void){ |

113 | 82ab5ad7 | Michael Niedermayer | LLSModel m; |

114 | 408ec4e2 | Michael Niedermayer | ```
int i, order;
``` |

115 | 82ab5ad7 | Michael Niedermayer | |

116 | ```
av_init_lls(&m, 3);
``` |
||

117 | |||

118 | for(i=0; i<100; i++){ |
||

119 | double var[4]; |
||

120 | 38c162c1 | Diego Biurrun | ```
double eval;
``` |

121 | 408ec4e2 | Michael Niedermayer | var[0] = (rand() / (double)RAND_MAX - 0.5)*2; |

122 | var[1] = var[0] + rand() / (double)RAND_MAX - 0.5; |
||

123 | var[2] = var[1] + rand() / (double)RAND_MAX - 0.5; |
||

124 | var[3] = var[2] + rand() / (double)RAND_MAX - 0.5; |
||

125 | 82ab5ad7 | Michael Niedermayer | av_update_lls(&m, var, 0.99); |

126 | 408ec4e2 | Michael Niedermayer | av_solve_lls(&m, 0.001, 0); |

127 | for(order=0; order<3; order++){ |
||

128 | ```
eval= av_evaluate_lls(&m, var+1, order);
``` |
||

129 | 578f90a8 | Diego Biurrun | ```
printf("real:%9f order:%d pred:%9f var:%f coeffs:%f %9f %9f\n",
``` |

130 | 408ec4e2 | Michael Niedermayer | var[0], order, eval, sqrt(m.variance[order] / (i+1)), |

131 | m.coeff[order][0], m.coeff[order][1], m.coeff[order][2]); |
||

132 | } |
||

133 | 82ab5ad7 | Michael Niedermayer | } |

134 | return 0; |
||

135 | } |
||

136 | |||

137 | `#endif` |