/*
 * SDSPIDCT  IDCT S-function.
 * DSP Blockset S-Function to perform an IDCT complex data.
 *
 *  Copyright 1995-2000 The MathWorks, Inc.
 *  $Revision  $  $Date: 2000/05/05 20:16:01 $
 */

#define S_FUNCTION_NAME sdspidct

#include <math.h>
#include <string.h>

#include "dsp_sim.h"

static void mdlInitializeSizes(SimStruct *S)
{
  ssSetNumInputs(        S, DYNAMICALLY_SIZED);
  ssSetNumOutputs(       S, DYNAMICALLY_SIZED);
  ssSetDirectFeedThrough(S, 1);
  ssSetNumSampleTimes(   S, 1);
  ssSetNumRWork(         S, DYNAMICALLY_SIZED);
  ssSetOptions(          S, SS_OPTION_EXCEPTION_FREE_CODE |
	                    	SS_OPTION_USING_ssGetUPtrs);  
}


static void mdlInitializeSampleTimes(SimStruct *S)
{

  ssSetSampleTimeEvent(S, 0, INHERITED_SAMPLE_TIME);
  ssSetOffsetTimeEvent(S, 0, 0.0);
}


static void mdlInitializeConditions(real_T *x0, SimStruct *S)
{
 
    int_T numinputs = ssGetNumInputs(S);   /* real inputs */

    #ifdef MATLAB_MEX_FILE
      /*
       * Check to make sure input size is a power of 2:
       */

    if (numinputs != -1) {
   	int_T dummy;
	if (frexp((real_T)numinputs, &dummy) != 0.5) {
            THROW_ERROR(S,"Width of input to IDCT must be a power of 2");
	}
	if (numinputs == 1) {
    	    mexWarnMsgTxt("Computing the IDCT of a scalar input.\n");
 	}
    }
    #endif /* MATLAB_MEX_FILE */

  /*  Compute complex weights  */
  {
    const int_T  N	  = numinputs/2; 		 /* number of complex inputs */
    const real_T piN2 = 2.0*atan(1.0)/ N;  	 /* pi / (2*N) */
    const real_T fac  = sqrt(2.0*N);
    real_T      *ww   = ssGetRWork(S);
    int_T        i;
    
    for (i = 0; i < N; i++) {
      ww[i]     =  cos(i*piN2) * fac;
      ww[i + N] =  sin(i*piN2) * fac;
    }
    
    ww[0] *= sqrt(2.0);

  }
}


static void mdlOutputs(real_T *y, const real_T *x, const real_T *u, 
                       SimStruct *S, int_T tid)
{
  if (ssIsMajorTimeStep(S)) {
    const int_T numinputs = ssGetNumInputs(S);  /* real + imag inputs */
    const int_T	N         = numinputs/2;        /* number of complex input */
    real_T     *yr        = y;
    real_T     *yi        = y+numinputs;
    real_T     *wr        = ssGetRWork(S);
    real_T     *wi        = wr+N;
    UPtrsType   uptr      = ssGetUPtrs(S);
	UPtrsType   ur        = uptr;
	UPtrsType   ui        = ur+N;
	int_T       i;
    /*
     * NOTE: We leave a "cushion" after the real and imaginary parts:
     *       [ [Real (N/2)]  [Extra (N/2)] [Imag (N/2)] [Extra (N/2)] ]
     */

    /* Multiply by weights */

    for (i = 0; i< N; i++) {
		*yr++ = **(ur+i) * *(wr+i) - **(ui+i) * *(wi+i);
		*yi++ = **(ur+i) * *(wi+i) + **(ui+i) * *(wr+i);
    }

	*yr++=0;
	*yi++=0;
	
	/* Equivalent MATLAB code: 
	 * yy(n+2:2*n,:) = -j*W(2:n,:).*flipud(bb(2:n,:));
	 */
	ur = uptr + N;			/* Pointing to the last real input */
	ui = uptr + numinputs;  /* Pointing to the last imag input */

	for(i = 1; i < N; i++) {
		*yr++ = *(wr+i) * **(ui-i) + *(wi+i) * **(ur-i);
		*yi++ = *(wi+i) * **(ui-i) - *(wr+i) * **(ur-i); 
	}
  }
}


static void mdlUpdate(real_T *x, const real_T *u, SimStruct *S, int_T tid)
{
}


static void mdlDerivatives(real_T *dx, const real_T *x, const real_T *u, 
                           SimStruct *S, int_T tid)
{
}


static void mdlTerminate(SimStruct *S)
{
}

#ifdef MATLAB_MEX_FILE

#define MDL_GET_INPUT_PORT_WIDTH
/* Function: mdlGetInputPortWidth =============================================
 *
 */
static int_T mdlGetInputPortWidth(SimStruct *S, int_T outputPortWidth)
{
  if (outputPortWidth % 2 != 0) {
    ssSetErrorStatus(S,"Invalid output port width");
    return(1);
  }
  return(outputPortWidth/2);
}

#define MDL_GET_OUTPUT_PORT_WIDTH
/* Function: mdlGetOutputPortWidth ============================================
 *
 */
static int_T mdlGetOutputPortWidth(SimStruct *S, int_T inputPortWidth)
{
    return(inputPortWidth*2);
}

#endif

/* Function: mdlSetWorkWidths =================================================
 * Abstract:
 *
 */
#ifdef MATLAB_MEX_FILE
#define MDL_SET_WORK_WIDTHS
static void mdlSetWorkWidths(SimStruct *S)
{
  ssSetNumRWork(S, (int_T)(ssGetNumInputs(S)));
}
#endif


#include "dsp_trailer.c"
 
/* [EOF]: sdspidct.c */