viterbi.cpp

#include <stdio.h>
typedef struct{
	char model_name;
	int state_num = 3;
	int observ_num = 8;
	double initial[3];			//initial prob.
	double transition[3][3];	//transition prob.
	double observation[8][3];	//observation prob.
} HMM;

void viterbi(HMM*, int, int*, int, double*);

int main(){

	const int model_num = 3;
	HMM hmm[model_num];
	
	/*  part for load model
	    part for load model
	    part for load model   */
	const int seqlen= 500;
	int seq[seqlen];
	
	/*	load observation sequence
		load observation sequence
		load observation sequence   */
	double p[model_num]={0};
	
	viterbi(hmm, seqlen, seq, model_num, p);
	
	double maxp={0};
	int max_model;
	for(int i=0;i<model_num;i++)
	{
		if(p[i]>maxp)
		{
			maxp=p[i];
			max_model= i;
		}
	}

	printf("seq is generated by model %s", hmm[max_model].model_name);	
	return 0;
}	
void viterbi(HMM* hmm, int seqlen, int* seq, int model_num, double* p) {
	

	for(int i=0;i<model_num;i++)
	{
		double viterbi[hmm[model_num].state_num][seqlen];
		for(int j=0;j<hmm[model_num].state_num;j++)
		{
			viterbi[j][0]=hmm[i].initial[j]*hmm[i].observation[seq[0]][j];
		}
		for(int k=1;k<seqlen;k++)
		{
			for(int l=0;l<hmm[model_num].state_num;l++)
			{
				double max=0;
				for(int m=0;m<hmm[model_num].state_num;m++)
				{
					if(viterbi[m][k-1]*hmm[i].transition[m][l]>max)
					{
						max=viterbi[m][k-1]*hmm[i].transition[m][l];
					}
				}
				viterbi[l][k]=max*hmm[i].observation[seq[k]][l];
			}
		}
		for(int a=0;a<hmm[model_num].state_num;a++)
		{
			if(viterbi[a][seqlen-1]>p[i])
			{
				p[i]=viterbi[a][seqlen-1];
			}
		}
	}
	//return p;
	
}