moveint.c

数字通信第四版原书的例程

/*
 * MOVEINT   A SIMULINK moving integration
 *
 *           Syntax:  [sys, x0] = sbuffer(t, x, u, flag, In_length, Time_window, Sample_time, Method)
 *
 *      In-length: the vector length of the input vector.
 *      Time_window: int_t^(t+T) f(x) dx, the scale T, window size of the integration.
 *      Sample_time: sampling time of the integration. This is a discrete-time block.
 *      Method: integration method. Could be choose from: forward, backward, first_order.
 *
 * Wes Wang Dec. 8, 1994
 * Copyright (c) 1994-96 by The MathWorks, Inc.
 * All Rights Reserved
 * $Revision: 1.1 $  $Date: 1996/04/01 19:03:36 $
 */

#define S_FUNCTION_NAME moveint

#ifdef MATLAB_MEX_FILE
#include <stdio.h>    /* needed for declaration of sprintf */
#include "mex.h"      /* needed for declaration of mexErrMsgTxt */
#endif

/*
 * need to include simstruc.h for the definition of the SimStruct and
 * its associated macro definitions.
 */

#include "simstruc.h"

/*
 * Defines for easy access of the input parameters
 */

#define NUM_ARGS   4
#define IN_LENGTH      ssGetArg(S,0)
#define TIME_WINDOW        ssGetArg(S,1)
#define SAMPLE_TIME        ssGetArg(S,2)
#define METHOD     ssGetArg(S,3)

/*
 * mdlInitializeSizes - called to initialize the sizes array stored in
 *                      the SimStruct.  The sizes array defines the
 *                      characteristics (number of inputs, outputs,
 *                      states, etc.) of the S-Function.
 */

static void mdlInitializeSizes(S)
    SimStruct *S;
{
    /*
     * Set-up size information.
     */ 
    
    if (ssGetNumArgs(S) == NUM_ARGS) {
      int in_length, num_points;
    char method[2];
/*    char err_msg[256];
 */
    
    if ((mxGetN(IN_LENGTH) != 1) || (mxGetM(IN_LENGTH) != 1)) {
#ifdef MATLAB_MEX_FILE
           mexErrMsgTxt("The inout vector length must be a scalar");
#endif
     }

       if ((mxGetN(TIME_WINDOW) != 1) || (mxGetM(TIME_WINDOW) != 1)) {
#ifdef MATLAB_MEX_FILE
       mexErrMsgTxt("The integration time window must be a scalar");
#endif
 }

       if ((mxGetN(SAMPLE_TIME)*mxGetM(SAMPLE_TIME) != 1) &&
       (mxGetN(SAMPLE_TIME)*mxGetM(SAMPLE_TIME) != 2)) {
#ifdef MATLAB_MEX_FILE
         mexErrMsgTxt("The sample time must be a scalar or a vector of length 2");
#endif
     }

       if (mxGetM(METHOD) != 1) {
#ifdef MATLAB_MEX_FILE
            mexErrMsgTxt("Integration method must be a string");
#endif
  }
       
        in_length  = mxGetPr(IN_LENGTH)[0];
     num_points = (int)(mxGetPr(TIME_WINDOW)[0] / mxGetPr(SAMPLE_TIME)[0] + 0.49) + 1;
       if (num_points <= 2) {
#ifdef MATLAB_MEX_FILE
        mexErrMsgTxt("The integration time window must be greater than sample time.");
#endif
        } else if (in_length < 1) {
#ifdef MATLAB_MEX_FILE
           mexErrMsgTxt("The input vector length must be a positive number.");
#endif
   }
       method[0] = (char)mxGetPr(METHOD)[0];
   method[1] = (char)mxGetPr(METHOD)[1];   
        
   ssSetNumContStates(    S, 0);
   ssSetNumDiscStates(    S, 0);
   ssSetNumInputs(        S, in_length);
   ssSetNumOutputs(       S, in_length);
   if ((method[0] == 'f') && (method[1] == 'o'))
       ssSetDirectFeedThrough(S, 0);
   else
       ssSetDirectFeedThrough(S, 1);
   ssSetNumInputArgs(     S, NUM_ARGS);
   ssSetNumSampleTimes(   S, 1);
   ssSetNumRWork(         S, in_length*(num_points + 1) + 1);
      /* storage: t, points, last_integration_value */
   ssSetNumIWork(         S, 2);
   /* Total_points=num_points, current points */
   ssSetNumPWork(         S, 0);
/*    sprintf(err_msg, "Initialization num_points = %i, in_length = %i\n",num_points, in_length);
    mexPrintf(err_msg);
 */
    } else {
#ifdef MATLAB_MEX_FILE
      char err_msg[256];
      sprintf(err_msg, "Wrong number of input arguments passed to S-function MEX-file.\n"
             "%d input arguments were passed in when expecting %d input arguments.\n", ssGetNumArgs(S) + 4, NUM_ARGS + 4);
   mexErrMsgTxt(err_msg);
#endif
    }
}

/*
 * mdlInitializeSampleTimes - initializes the array of sample times stored in
 *                            the SimStruct associated with this S-Function.
 */

static void mdlInitializeSampleTimes(S)
    SimStruct *S;
{
    double sampleTime, offsetTime;

    /*
     * Note, blocks that are continuous in nature should have a single
     * sample time of 0.0.
     */

    sampleTime = mxGetPr(SAMPLE_TIME)[0];
    if ((mxGetN(SAMPLE_TIME) * mxGetM(SAMPLE_TIME)) == 2)
    offsetTime = mxGetPr(SAMPLE_TIME)[1];
    else
   offsetTime = 0.;
    
    ssSetSampleTimeEvent(S, 0, sampleTime);
    ssSetOffsetTimeEvent(S, 0, offsetTime);
}

/*
 * mdlInitializeConditions - initializes the states for the S-Function
 */

static void mdlInitializeConditions(x0, S)
    double *x0;
    SimStruct *S;
{
    int in_length = mxGetPr(IN_LENGTH)[0];
    int num_points = (int)(mxGetPr(TIME_WINDOW)[0] / mxGetPr(SAMPLE_TIME)[0] + 0.49) + 1;
    
    double *last_time  = ssGetRWork(S);
    double *storage_point = ssGetRWork(S) + 1;
    double *integ_value   = ssGetRWork(S) + in_length * num_points + 1;
    
    int *cur_point       = ssGetIWork(S);
    int i;
    
    /* 
     * Initialize all storage and integration values
     */

    for (i = 0; i < in_length * (num_points + 1); i++)
      *storage_point++ = 0.;
    *last_time = -1;
    /*
     * IWork
     */
    
    *cur_point++       = num_points;
    *cur_point         = 0;
}

/*
 * mdlOutputs - computes the outputs of the S-Function
 */

static void mdlOutputs(y, x, u, S, tid)
    double *y, *x, *u;
    SimStruct *S;
    int tid;
{
    int *num_points     = ssGetIWork(S);
    int *cur_point      = ssGetIWork(S)+1;
    int in_length       = mxGetPr(IN_LENGTH)[0];

    double *last_time = ssGetRWork(S);
    double *storage_point = ssGetRWork(S) + 1;
    double *integ_value   = ssGetRWork(S) + (in_length * (*num_points) + 1);
    
    double sample_time = mxGetPr(SAMPLE_TIME)[0];

    double current_time;
    int i;
    char method[2];

    current_time  = ssGetT(S);
     method[0] = (char)mxGetPr(METHOD)[0];
   method[1] = (char)mxGetPr(METHOD)[1];   

        if (*num_points > 0) {
    if ((current_time - *last_time) > sample_time/1000000.) {
/*        char err_msg[256];        */
        int increament, point_before, point_after;

        point_before = (*num_points + *cur_point - 1) % (*num_points);
        point_after  = (*cur_point + 1) % (*num_points);
/*    sprintf(err_msg, "Output: num_points = %i, in_length = %i\n",*num_points, in_length);
    mexPrintf(err_msg);
 */
        for (i=0; i<in_length; i++) {
            increament = (*num_points) * i;
/*    sprintf(err_msg, "cur_point = %i, point_before = %i, point_after %i, increament %i\n",*cur_point, point_before, point_after, increament);
    mexPrintf(err_msg);
 */
         if ((method[0] == 'f') && (method[1] == 'o')) {
/*    sprintf(err_msg, "integ_value = %f, storage_before = %f, current_point = %f \n", integ_value[i], storage_point[point_before+increament], storage_point[*cur_point+increament]);
    mexPrintf(err_msg);
 */
              integ_value[i] = integ_value[i]
                     + (storage_point[point_before + increament]
                        - storage_point[*cur_point + increament]) * sample_time;
             } else if ((method[0] == 'b') && (method[1] == 'a')) {
              integ_value[i] = integ_value[i]
                    + (u[i]
                        - storage_point[point_after + increament]) * sample_time;
         } else {
                    integ_value[i] = integ_value[i]
                    + (u[i]
                        + storage_point[point_before + increament]
                        - storage_point[(*cur_point) + increament]                    
                        - storage_point[point_after + increament]) * sample_time / 2;
             }
            storage_point[(*cur_point) + increament] = u[i];
        }
        *cur_point = point_after;
        *last_time = current_time;
    } else {
/*        char err_msg[256];        */
        for (i=0; i<in_length; i++) {
           int increament;
            increament = (*num_points) * i;
/*    sprintf(err_msg, "\n cur_point = %i, increament %i\n",*cur_point, increament);
    mexPrintf(err_msg);
 */
             if ((method[0] == 'f') && (method[1] == 'o')) {
         } else if ((method[0] == 'b') && (method[1] == 'a')) {
              integ_value[i] = integ_value[i]
                    + (u[i]
                        - storage_point[*cur_point + increament]) * sample_time;
          } else {
                    integ_value[i] = integ_value[i]
                    + (u[i]
                        - storage_point[(*cur_point) + increament]) * sample_time / 2;
            }
            storage_point[(*cur_point) + increament] = u[i];
        }
    }
     }
    for (i = 0; i<in_length; i++)
        y[i] = integ_value[i];
}

/*
 * mdlUpdate - computes the discrete states of the S-Function
 */

static void mdlUpdate(x, u, S, tid)
    double *x, *u;
    SimStruct *S;
    int tid;
{
}

/*
 * mdlDerivatives - computes the derivatives of the S-Function
 */

static void mdlDerivatives(dx, x, u, S, tid)
    double *dx, *x, *u;
    SimStruct *S;
    int tid;
{
}

/*
 * mdlTerminate - called at termination of model execution.
 */

static void mdlTerminate(S)
    SimStruct *S;
{
}

#ifdef   MATLAB_MEX_FILE    /* Is this file being compiled as a MEX-file? */
#include "simulink.c"      /* MEX-File interface mechanism */
#else
#include "cg_sfun.h"       /* Code generation registration function */
#endif