1 /*
2 * principal component analysis (PCA)
3 * Copyright (c) 2004 Michael Niedermayer <michaelni@gmx.at>
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 02110-1301 USA
20 */
21
22 #include "libavutil/pca.c"
23 #include "libavutil/lfg.h"
24
25 #undef printf
26 #include <stdio.h>
27 #include <stdlib.h>
28
main(void)29 int main(void){
30 PCA *pca;
31 int i, j, k;
32 #define LEN 8
33 double eigenvector[LEN*LEN];
34 double eigenvalue[LEN];
35 AVLFG prng;
36
37 av_lfg_init(&prng, 1);
38
39 pca= ff_pca_init(LEN);
40
41 for(i=0; i<9000000; i++){
42 double v[2*LEN+100];
43 double sum=0;
44 int pos = av_lfg_get(&prng) % LEN;
45 int v2 = av_lfg_get(&prng) % 101 - 50;
46 v[0] = av_lfg_get(&prng) % 101 - 50;
47 for(j=1; j<8; j++){
48 if(j<=pos) v[j]= v[0];
49 else v[j]= v2;
50 sum += v[j];
51 }
52 /* for(j=0; j<LEN; j++){
53 v[j] -= v[pos];
54 }*/
55 // sum += av_lfg_get(&prng) % 10;
56 /* for(j=0; j<LEN; j++){
57 v[j] -= sum/LEN;
58 }*/
59 // lbt1(v+100,v+100,LEN);
60 ff_pca_add(pca, v);
61 }
62
63
64 ff_pca(pca, eigenvector, eigenvalue);
65 for(i=0; i<LEN; i++){
66 pca->count= 1;
67 pca->mean[i]= 0;
68
69 // (0.5^|x|)^2 = 0.5^2|x| = 0.25^|x|
70
71
72 // pca.covariance[i + i*LEN]= pow(0.5, fabs
73 for(j=i; j<LEN; j++){
74 printf("%f ", pca->covariance[i + j*LEN]);
75 }
76 printf("\n");
77 }
78
79 for(i=0; i<LEN; i++){
80 double v[LEN];
81 double error=0;
82 memset(v, 0, sizeof(v));
83 for(j=0; j<LEN; j++){
84 for(k=0; k<LEN; k++){
85 v[j] += pca->covariance[FFMIN(k,j) + FFMAX(k,j)*LEN] * eigenvector[i + k*LEN];
86 }
87 v[j] /= eigenvalue[i];
88 error += fabs(v[j] - eigenvector[i + j*LEN]);
89 }
90 printf("%f ", error);
91 }
92 printf("\n");
93
94 for(i=0; i<LEN; i++){
95 for(j=0; j<LEN; j++){
96 printf("%9.6f ", eigenvector[i + j*LEN]);
97 }
98 printf(" %9.1f %f\n", eigenvalue[i], eigenvalue[i]/eigenvalue[0]);
99 }
100
101 return 0;
102 }
103