/*----------------------------MegaWave2 module----------------------------*/
/* mwcommand
name = {sample_gm2};
version = {"1.0"};
author = {"Lionel Moisan"};
function = {"Draw samples from a Gaussian Mixture"};
usage = {
 'n':[n=100]->n    "number of samples (default: 100)",
 's':[s=1.]->s     "standart deviation (default 1.)",
 c->c              "centroids (Flist, arbitrary dimension)",
 out<-sample_gm2   "output Flist"
        };
*/

#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include <time.h>
#include "mw.h"

double loi_normale()
{
  double a,b;

  a = drand48();
  b = drand48();
  return (sqrt(-2.*log(a))*cos(2.*M_PI*b));
}

Flist sample_gm2(c,n,s)
     Flist c;
     int *n;
     float *s;
{
  Flist out;
  int i,j,k;

  /*** Initialize random seed  ***/
  srand48( (long int) time (NULL) );

  out = mw_change_flist(NULL,*n,*n,c->dim);
  for (k=0;k<*n;k++) {
    i = (int)(drand48()*(double)c->size);
    if (i>=c->size) i=c->size-1;
    for (j=0;j<c->dim;j++) {
      out->values[j+c->dim*k]
	= c->values[j+c->dim*i] + *s*(float)loi_normale();
    }
  }
    
  return(out);
}