ExampleRandom2.cpp

/* 
This file is part of the BIAS library (Basic ImageAlgorithmS).

Copyright (C) 2003-2009    (see file CONTACT for details)
  Multimediale Systeme der Informationsverarbeitung
  Institut fuer Informatik
  Christian-Albrechts-Universitaet Kiel


BIAS is free software; you can redistribute it and/or modify
it under the terms of the GNU Lesser General Public License as published by
the Free Software Foundation; either version 2.1 of the License, or
(at your option) any later version.

BIAS is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
GNU Lesser General Public License for more details.

You should have received a copy of the GNU Lesser General Public License
along with BIAS; if not, write to the Free Software
Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
*/

/**
   @example ExampleRandom2.cpp
   @relates Random
   @brief Example2 for random number generator
   @ingroup g_examples
   @author MIP 
*/
#include <cmath>
#include <cstring>
#include <Base/Math/Random.hh>

using namespace BIAS;
using namespace std;

#define COUNT 1000

int main()
{
  Random ran1, ran2;
  double num1, num2;
  double mean1 = 10.0, mean2 = -30.0;
  double sigma1 = 0.5, sigma2 = 10.0;
  double angle = 0.16;
  double dist[COUNT][2];
  double covm[2][2];
  double mean[2], var[2];

  // zero mean and var
  for (int i=0; i<2; i++){
    mean[i] = var[i] = 0.0;
    for (int j=0; j<2; j++)
      covm[i][j]=0.0;
  }

  // calculate a "rotated" 2dimensional distribution for PCA
  for (register unsigned int i=0; i<COUNT; i++){
    num1=ran1.GetNormalDistributed(mean1, sigma1);
    num2=ran2.GetNormalDistributed(mean2, sigma2);
    dist[i][0]=num1 * cos(angle) - num2 * sin(angle);
    dist[i][1]=num1 * sin(angle) + num2 * cos(angle);
    mean[0]+=dist[i][0];
    mean[1]+=dist[i][1];
    cout << dist[i][0] << "\t" << dist[i][1] << endl;
  }
  
  mean[0]/=COUNT;
  mean[1]/=COUNT;
  cerr << "mean:" << endl;
  cerr << mean[0] << "\t" << mean[1] << endl;

  // calculate covariance matrix
  for (register unsigned int i=0; i<COUNT; i++){
    covm[0][0]+=(dist[i][0]-mean[0])*(dist[i][0]-mean[0]);
    covm[1][1]+=(dist[i][1]-mean[1])*(dist[i][1]-mean[1]);
    covm[0][1]+=(dist[i][1]-mean[1])*(dist[i][0]-mean[0]);
  }
  covm[0][0]/=COUNT-1;
  covm[1][1]/=COUNT-1;
  covm[0][1]/=COUNT-1;
  covm[1][0]=covm[0][1];

  // output of the covariance matrix
  cerr << "covariance matrix:" << endl;
  for (int i=0; i<2; i++){
    for (int j=0; j<2; j++){
      cerr << covm[i][j] << "\t";
    }
    cerr << endl;
  }

  return 0;
}