View Raw SPL
/*****************************************************************************
*                                                                            *
*   ZP2SOS.SPL   Copyright (C) 2015 DSP Development Corporation              *
*                               All Rights Reserved                          *
*                                                                            *
*   Author:      Randy Race                                                  *
*                                                                            *
*   Synopsis:    Zeros, poles and gain to SOS form                           *
*                                                                            *
*   Revisions:   15 Jun 2015  RRR  Creation                                  *
*                                                                            *
*****************************************************************************/

#if @HELP_ZP2SOS

    ZP2SOS

    Purpose: Converts zeros, poles and gain to secord order section form.

    Syntax:  ZP2SOS(z, p, k, "flag")

             (sos, g) = ZP2SOS(z, p, k, "flag")

                   z - A series, the zeros.

                   p - A series, the poles.

                   k - Optional. A scalar, the gain. Defaults to 1.0.

              "flag" - Optional. A string, the stage ordering flag.

                          "up" : order stages such that the poles of each
                                 stage are closer to the unit circle than 
                                 the previous stage (default)

                        "down" : order stages such that the poles of each
                                 stage are closer to zero than the previous
                                 stage

    Returns: An Nx6 array, the system coefficients in second order section
             form.

             (sos, g) = ZP2SOS(z, p, k) returns the SOS coefficients and
             gain in two separate variables.

    Example:
             z = {0.0, 2.0};
             p = {0.5, 0.2};
             k = 1.0;

             sos = zp2sos(z, p, k);

             sos == {{1, -2, 0, 1, -0.7, 0.1}}

             The SOS coefficients represent the system:

                                 1 - 2z^(-1)
             H(z) = (1) * --------------------------
                          1 - 0.7 z^(-1) + 0.1 z^(-2)

    Remarks:

             ZP2SOS converts the zeros, poles and gain of a discrete
             system to second order section form where the coefficients
             represent the following Z transform:

                               -1       -2                -1       -2
                    b10 + b11 z  + b12 z       b20 + b21 z  + b22 z
         H(z) = G * ________________________ *  ________________________ ...
                               -1       -2                -1       -2
                      1 + a11 z  + a12 z         1 + a21 z  + a22 z


             The SOS form is an Nx6 array:

             {{b10, b11, b12, 1, a11, a12},
              {b20, b21, b22, 1, a21, a22}, 

                            ...

              {bN0, bN1, bN2, 1, aN1, aN2}}

             Each row of the SOS coefficients represents a second order stage.

             By default, the stages are ordered such that the poles of
             each stage are closer to the unit circle than the previous
             stage. The zeros of each stage are chosen to be closest to
             the poles of the same stage.

             If FLAG == "down", the stages are reversed.

             See SOS2ZP to convert SOS coefficients to zeros, poles and
             gain terms.

    See Also:
             Cas2sos
             Cas2zp
             Cascade
             Sos2zp
             Tf2sos
             Zp2cas
             Zp2sos
             Zp2tf
#endif


/* zeros, poles and gain to seond order section form */
zp2sos(z, p, k, flag)
{
        local cas, sos, g;

        /* convert to cascade form */
        cas = zp2cas(z, p, k, flag);

        /* convert to SOS form */
        (sos, g) = cas2sos(cas);

        if (outargc > 1)
        {
                return(sos, g);
        }
        else
        {
                /* combine gain to first numerator stage */
                sos[1, 1..3] *= g;

                return(sos);
        }
}