/* Change this to reflect the apropriate header file */
#include <math.h> 
#include "c:\matlab6p1\extern\include\mex.h"
// #include "c:\matlabr12\extern\include\mex.h"
// #include "c:\matlabr11\extern\include\mex.h"

/*
* garchcore.c - 
* This is a helper function and is part of the UCSD_GARCH toolbox 
* You can compile it and should work on any platform.
* 
* Author: Kevin Sheppard
* kksheppard@ucsd.edu
* Revision: 2    Date: 12/31/2001
*/
void makeh(double *data, double *omega, double *alpha, double *beta, int m, int T, int k, int dimP, int dimQ, double *covEst, double *Ht)
{
	int i, j, ii, jj, index;
	for (i=0; i<k; i++)
	{
		index=i*T;
		for (j=0; j<m; j++) 
		{
		*(Ht+j+index)=*(covEst+i);
		}
	}
	
	for (i=m; i<T; i++)  // loop over all T
	{
		for (j=0;j<k;j++) // loop over all elements of H
		{
			*(Ht+i+(j*T))=*(omega+j);
		}

		index=0;
		for (ii=0;ii<dimP;ii++) // loop over all dims of P
		{
			for (j=0;j<k;j++) // loop over all elements of H
			{
				for (jj=0;jj<k;jj++)
				{
					*(Ht+i+j*T)+=*(data+(i-1-ii+jj*T))* *(alpha+index);
					index++;
				}
			}
		}
		index=0;
		for (ii=0;ii<dimP;ii++) // loop over all dims of P
		{
			for (j=0;j<k;j++) // loop over all elements of H
			{
				for (jj=0;jj<k;jj++)
				{
					*(Ht+i+j*T)+=*(Ht+(i-1-ii+jj*T)) * *(beta+index);
					index++;
				}
			}
		}	
	}
}





/* The gateway routine */
void mexFunction( int nlhs, mxArray *plhs[],
				 int nrhs, const mxArray *prhs[])
{
	double *data, *omega, *alpha, *beta, *Ht, *covEst;
	int m,k,t,dimP,dimQ;
	int     mrows,ncols;
	
	/*  Check for proper number of arguments. */
	if(nrhs!=10) 
		mexErrMsgTxt("Ten inputs required.");
	if(nlhs!=1) 
		mexErrMsgTxt("One output required.");
	
	/*  Get the scalar inputs */
	m = (int)mxGetScalar(prhs[4]);
	k = (int)mxGetScalar(prhs[5]);
	t = (int)mxGetScalar(prhs[6]);
	dimP = (int)mxGetScalar(prhs[7]);
	dimQ = (int)mxGetScalar(prhs[8]);
    
	/*  Create a pointer to the input matrices . */
	data           = mxGetPr(prhs[0]);
	omega     = mxGetPr(prhs[1]);
	alpha     = mxGetPr(prhs[2]);
	beta     = mxGetPr(prhs[3]);
	covEst         = mxGetPr(prhs[9]);
	
	/*  Get the dimensions of the matrix input ht to make an output matrix. */
	mrows = mxGetM(prhs[0]);
	ncols = mxGetN(prhs[0]);
	
	/*  Set the output pointer to the output matrix. */
	plhs[0] = mxCreateDoubleMatrix(mrows,ncols, mxREAL);
	
	/*  Create a C pointer to a copy of the output matrix. */
	Ht = mxGetPr(plhs[0]);
	
	/*  Call the C subroutine. */
	makeh(data, omega, alpha, beta, m, mrows, k, dimP, dimQ, covEst, Ht);
	
}