www.pudn.com > EMGMM.rar > guassian.cpp

//#include "stdafx.h"  
#include  
#include   
#include  
#include  
#include   
#include  
#include "gaussian.h" 
//#include "KMEANS.h" 
using namespace std; 
//#define n (pow(2,31)-1) 
 
#define PI 3.1415926535897532 
//using namespace std; 
 
void Gaussian::Initialize(char *fname) 
//int pattern,int dimension,double **sample,int mix_num) 
/* Error handler */ 
{ 
	m_kmeans.LoadPatterns(fname); 
	m_kmeans.InitClusters(); 
	m_kmeans.RunKMeans(); 
 
	NumPatterns=m_kmeans.NumPatterns; 
	SizeVector=m_kmeans.SizeVector; 
	NumMixture=m_kmeans.NumClusters; 
 
	samplelist=matrix(NumPatterns,SizeVector); 
	mean = vect(NumMixture*SizeVector); 
	devia= vect(NumMixture*SizeVector); 
	a= vect(NumMixture); 
 
	clu_idx=(int *) malloc (NumPatterns*sizeof(int)); 
	num_samples_clu=(int *) malloc (NumMixture*sizeof(int)); 
 
	int i,j; 
	for (i=0;iNumMembers ; i++) 
		{ 
			Uni_mean[j] += data[Clus->Member[i]*SizeVector+j]; 
		} 
		Uni_mean[j] /= (double)Clus->NumMembers; 
	} 
	 
	/* Center the column vectors. */ 
	 
	for (i = 0; i < SizeVector; i++) 
    { 
		Uni_devia[i]=0; 
		for (j = 0; j < Clus->NumMembers; j++) 
        { 
			Uni_devia[i]+=(data[Clus->Member[j]*SizeVector+i]-Uni_mean[i])*(data[Clus->Member[j]*SizeVector+i]-Uni_mean[i]); 
        } 
		Uni_devia[i]/=(double)Clus->NumMembers; 
    } 
	 
	return; 
	 
} 
 
void Gaussian::erhand(char err_msg[]) 
/* Error handler */ 
{ 
    fprintf(stderr,"Run-time error:\n"); 
    fprintf(stderr,"%s\n", err_msg); 
    fprintf(stderr,"Exiting to system.\n"); 
    exit(1); 
} 
 
double ** Gaussian::matrix(int n,int m) 
/* Allocate a double matrix with range [1..n][1..m]. */ 
{ 
    int i; 
    double **mat; 
	 
    /* Allocate pointers to rows. */ 
    mat = (double **) malloc((unsigned) (n)*sizeof(double*)); 
    if (!mat) erhand("Allocation failure 1 in matrix()."); 
	 
	for (i = 0; i < n; i++) 
	{ 
        mat[i] = (double *) malloc((unsigned) (m)*sizeof(double)); 
        if (!mat[i]) erhand("Allocation failure 2 in matrix()."); 
	} 
	 
	/* Return pointer to array of pointers to rows. */ 
	return mat; 
	 
} 
 
/**  Allocation of vector storage  ***********************************/ 
 
double * Gaussian::vect(int n) 
/* Allocates a double vector with range [1..n]. */ 
{ 
	 
    double *v; 
	 
    v = (double *) malloc ((unsigned) n*sizeof(double)); 
    if (!v) erhand("Allocation failure in vect()."); 
    return v; 
	 
} 
 
double * vecth(int n) 
/* Allocates a double vector with range [1..n]. */ 
{ 
	 
    double *v; 
	 
    v = (double *) malloc ((unsigned) n*sizeof(double)); 
    return v; 
	 
} 
 
void Gaussian::GaussEM(double **samples, int num_samples, int vec_size, \ 
		int num_cluster, int *cluster_index,int* num_samples_in_clusters,\ 
		double *means,double* delta,double* weight,int maxIter) 
{ 
	int j,k,x;//j:cluster;   x:sample;     k: dimension 
	//E step 
	double** p_x_j = new double*[num_samples];//P(sample_x|Cluster_j) 
	for(x=0;x0.0); 
			} 
			for(x=0;x=0.0); 
			} 
		} 
		 
		for(x=0;x=0.0); 
			} 
		} 
		 
		//M step 
		for(j=0;j0.0); 
			p_j[j] = sum_pjx/num_samples;//各类的概率 
			//update means vector 
			for(k=0;k=0.0&&means[j*vec_size+k]<100.0); 
			} 
			weight[j]=p_j[j]; 
 
			//update delta vector 
			for(k=0;k 0.0) 
            {  
                l_prob = ComputeUniModeGauss(vect, mu + m*vect_size,  
					inv_var + m * vect_size, 
					log_var_val[m],  
					vect_size);  
				 
                prob = prob + weight[m]*exp((double)l_prob); 
            } 
        }  
        prob = log(prob);  
    }                         
    return (double)prob;    
} 
 
void Gaussian::Classifier(double *vector) 
{ 
	int i,j; 
	double* inv_var; 
	double* log_var_val; 
	double temp; 
	inv_var=vect(NumMixture*SizeVector); 
	log_var_val=vect(NumMixture); 
	for (i=0;i