HMM forward algorithm - c source code

This is HMM forward algorithm source code.
I made this code today during about half the day.

I wouldn't introduce HMM and forward algorithm on this page. Because there are so many good material in the web.

The source code need Txt file included A, B, pi and observation values.
After load the txt file, the code return likelihood value.

Refer tp the real excute screen.

txt file include...
number of state, number of observation, number of observation sequence
state matrix by 1 line
observation matrix by 1 line
observation sequence

This picture describes the result of forward algorithm.

C source code

#include < stdio.h>


double HMM_forwardAlgorithm(double* &pi, double** &A, double ** &B, int* &obS, int sN, int oN, int obs_N);

void main()
{
 ///////////////////////////////////////////////
 //Input HMM ramda
 int stateN, observeN, ob_seq_N;
 double *pi;
 double **A;
 double **B;
 int *obS;
 ////////////////////////////////////////////////

 FILE * fp;
 fp = fopen("HMM_Ramda2.txt", "r");
 fscanf(fp,"%d %d %d",&stateN, &observeN, &ob_seq_N);
 /////////////////////////////////////////////////
 //alloc buffer
 pi = new double[stateN];
 A = new double*[stateN];
 B = new double*[stateN];
 obS = new int[ob_seq_N];

 for(int i=0; i< stateN; ++i)
 {
  A[i] = new double[stateN];
  B[i] = new double[observeN];
 }
 /////////////////////////////////////////////////
 //read Data 
 //pi
 for(int i=0; i< stateN; ++i)
  fscanf(fp, "%lf", &(pi[i]) );
 //A
 for(int i=0; i< stateN; ++i)
  for(int j=0; j< stateN; ++j)
   fscanf(fp, "%lf", &(A[i][j]) );
 //B
 for(int i=0; i< stateN; ++i)
  for(int j=0; j< observeN; ++j)
   fscanf(fp, "%lf", &(B[i][j]) );

 //observe sequence
 for(int i=0; i< ob_seq_N; ++i)
  fscanf(fp, "%d", &(obS[i]));

 ///////////////////////////////////////////////////
 HMM_forwardAlgorithm(pi, A, B, obS, stateN, observeN, ob_seq_N);
 

 ////////////////////////////////////////////////
 fclose(fp);

 for(int i=0; i< stateN; ++i)
 {
  delete[] A[i];
  delete[] B[i];
 }
 delete[] A;
 delete[] B;
 delete[] pi;
 delete[] obS; 

}



double HMM_forwardAlgorithm(double* &pi, double** &A, double ** &B, int* &obS, int sN, int oN, int obs_N)
{
 int i,j,t;
 //make forward matrix
 double** fMtx = new double *[sN];
 for(i=0; i< sN; ++i)
  fMtx[i] = new double[obs_N];

 //first step
 for(i=0; i< sN; ++i)
 {
  fMtx[i][0] = pi[i] * B[i][ obS[0]-1 ]; 
 }

 //routine
 double sum;
 for(t=1; t< obs_N; ++t)
 {  
  for(j=0; j< sN; ++j)
  {   
   sum=0;
   for(i=0; i< sN; ++i)
   {
    sum += (fMtx[i][t-1]*A[i][j]);    
   }
   fMtx[j][t] = sum * B[j][ obS[t]-1];
  }
 }


 //final
 sum=0; 
 for(j=0; j< sN; ++j)
 {
  sum += (fMtx[j][obs_N-1]);
 } 

 ////////////////////////////////////////////////////////
 //report
 printf("state N=%d / Observe N=%d\n\n", sN, oN);

 printf("Initial Matrix\n");
 for(int i=0; i< sN; ++i)
  printf("%lf ", pi[i] );
 printf("\n\n");

 printf("A matrix\n");
 for(int i=0; i< sN; ++i)
 {
  for(int j=0; j< sN; ++j)
  {
   printf("%lf ", A[i][j]);
  }
  printf("\n");
 }
 printf("\n");
 printf("B matrix\n");
 for(int i=0; i< sN; ++i)
 {
  for(int j=0; j< oN; ++j)
  {
   printf("%lf ", B[i][j]);
  }
  printf("\n");
 }
 printf("\n");
 printf("Observation sequence\n");
 for(int i=0; i< obs_N; ++i)
  printf("%d ", obS[i]);

 printf("\n\n");


 printf("** Trellis matrix **\n");
 for(j=0; j< sN; ++j)
 {
  for(t=0; t< obs_N; ++t)
  {
   printf("%lf  ", fMtx[j][t]);

   if(t!=obs_N-1)
    printf("->  ");
  }
  printf("\n");
 }


 printf("\nLikelihood = %lf \n", sum);

 for(i=0; i< sN; ++i)
  delete[] fMtx[i];
 delete[] fMtx;

 return sum;
}