git.sesse.net Git - ffmpeg/blob - libavutil/tests/pca.c

   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
  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 }