View Raw SPL
/*****************************************************************************
*                                                                            *
*   STEPZ.SPL    Copyright (C) 2015 DSP Development Corporation              *
*                               All Rights Reserved                          *
*                                                                            *
*   Author:      Randy Race                                                  *
*                                                                            *
*   Synopsis:    Step response of a Z transform                              *
*                                                                            *
*   Revisions:   17 Jun 2015  RRR  Creation from impz.spl                    *
*                                                                            *
*****************************************************************************/


#if @HELP_STEPZ

    STEPZ

    Purpose: Calculates the step response of a Z transform.

    Syntax:  STEPZ(b, a, N, Fs)
             (h, t) = STEPZ(b, a, N, Fs)

                 b - A series, the numerator coefficients in ascending
                     powers of z^(-1).

                 a - A series, the denominator coefficients in ascending
                     powers of z^(-1).

                 N - Optional. An integer specifying the length of the
                     output series. If not specifed or empty, the length
                     is determined by allowing the impulse response to decay
                     to -120dB or to display 5 periods if the output is not
                     significantly damped.

                Fs - Optional. A real, the sample rate of the output.
                     Defaults to rate(b) if rate(b) == rate(a), else
                     defaults to 1.0.

    Alternative Syntax:

             STEPZ(c, N, Fs)
             (h, t) = STEPZ(c, N, Fs)

                 c - A series, the system coefficients in cascaded
                     bi-quad form. If c contains 2 columns, c is
                     assumed to be in direct form, where the first
                     column is b and the second column is a.

                 N - Optional. An integer specifying the length of the
                     output series. If not specifed or empty, the length
                     is determined by allowing the impulse respone to decay
                     to -120dB or to display 5 periods if the output is not
                     significantly damped.

                Fs - Optional. A real, the sample rate of the output.
                     Defaults to rate(c).


    Returns: A series, the step response of the system.

             (h, t) = STEPZ(b, a) returns the step response and time
             values as two separate series.

    Example:
                             z           1
             Given H(z) = ------- = -----------
                          z - 0.5   1 - 0.5z^-1

                                            1                1
             The step response is H(z) * -------- = --------------------
                                         1 - z^-1   1 - 1.5z^-1 + 0.5z-2

                      2            -1
             =    --------  +  -----------
                  1 - z^-1     1 - 0.5z^-1


             or s[n] = 2 * (1)^n - (0.5)^n = 2 - 0.5^n

             n = 0..4;

             W1: stepz({1}, {1, -0.5}, 5)
             W2: impz({1}, {1, -1.5, 0.5}, 5)
             W3: 2 - 0.5^n

             W1 == W2 == W3 == {1, 1.5, 1.75, 1.875}

             W1 computes the step response directly from the system
             coefficients. W2 computes the step response from the transformed
             system coefficients. W3 computes the analytical step response.
             All three series are identical.

    Example:
                                 1 - 2z^-1
             Given H(z) = ----------------------
                          1 - 0.7z^-1 + 0.1z^- 2


             (r, p, k) = residuez({1, -2}, {1, -0.7, 0.1})

             r == {-5, 6}
             p == {.5, .2}
             k == {}

             Representing the partial fraction expansion:

                         -5               6
             H(z) = ------------- + -------------
                    (1 - 0.5z^-1)   (1 - 0.2z^-1)

             or h[n] = -5 * (0.5^n) + 6 * (0.2^n)

             n = ones(21, 1);
             W1: -5 * (0.5^n) + 6 * (0.2^n);stem
             W2: stepz({1, -2}, {1, -0.7, 0.1}, 21)
             W3: w1 - w2

             W1 contains 21 samples of the analytical step response
             as determined by the partial fraction expansion method.

             W2 displays 21 samples of the calculated response and
             W3 confirms that the difference is negligible.

    Remarks:
             The input series represent the terms of the rational polynomial
             H(z) = b(z)/a(z) where:

             M = length(b) and N = length(a):

                     b(z)    b[1] + b[2] z^(-1) + ... + b[M] z^(-M+1)
             H(z) = ------ = ----------------------------------------
                     a(z)    a[1] + a[2] z^(-1) + ... + a[N] z^(-N+1)


             If a[1] does not equal 1, the numerator and denominator terms
             are normalized by dividing each coefficient by a[1].

             If there are no output arguments, the result is displayed in
             the current window as a stem plot.

             For STEPZ(c, N, Fs), if the input c is a single column, the
             coefficients are assumed to be in cascaded bi-quad form.
             This is the output format of IIR filters designed by
             DADiSP/Filters.  If c contains 2 columns, the coefficients
             are assumed to be in direct form, where the first column
             is b and the second column is a.

    See Also:
             Filteq
             Impz
             Residuez
             Zfreq

#endif



/* Z transform step response */
stepz(b, a, N, Fs)
{
        local h, t;

        (h, t) = impz(b, a, N, Fs);

        h = cumsum(h);

        if (outargc > 1)
        {
                /* return step response and time values */
                return(h, xvals(h));
        }
        else if (outargc < 1)
        {
                /* default to stem plot */
                setplotstyle(h, 10);
        }
        
        return(h);
}