/*****************************************************************************
*                                                                            *
*   POLYSTAB.SPL  Copyright (C) 2009 DSP Development Corporation             *
*                               All Rights Reserved                          *
*                                                                            *
*   Author:       Randy Race                                                 *
*                                                                            *
*   Synopsis:     Reflects poles outside the unit circle to the inside       *
*                                                                            *
*   Revisions:    23 Jul 2009  RRR  Creation                                 *
*                                                                            *
*****************************************************************************/


#if @HELP_POLYSTAB

    POLYSTAB

    Purpose: Reflects polynomial roots outside the unit circle to the inside.

    Syntax:  POLYSTAB(a)

                a - A series, the polynomial coefficients in descending order.

    Returns: A series, the coefficients of the reflected polynomial.

    Example:
             W1: {1, -2, 5}
             W2: polystab(W1)
             W3: zplane({1}, W1)
             W4: zplane({1}, W2)

             W1 represents the coefficients of the polynomial
             
               1 - 2*z^-1 + 5*z^-2

             The roots of the polynomial are 1 + 2i and 1 - 2i, outside the
             unit circle. If these coefficients represent the denominator
             polynomial of a digital system, the system is unstable.

             W2 == {1, -.4, .2} the coefficients of the reflected polynomial.
             The roots of this polynomial are 0.2 + 4i, 0.2 - 4i, inside
             the unit circle. If the coefficients represent the denominator
             polynomial of a digital system, the system is stable.

             W3 and W4 graphically show the pole locations with respect
             to the unit circle.

    Remarks:
             POLYSTAB replaces the roots of a polynomial that lie outside
             the unit circle with roots that are inside the unit circle.

             Typically, the input series represents the coefficients of
             the following digital system:

               a[1] + a[2]*z^-1 + a[3]*z^-2 + ... + a[n+1]*z^-n

    See Also:
             poly
             roots
#endif



/* reflect polynomial roots outside of the unit circle to the inside */
ITERATE polystab(a)
{
        local r, idx, ind, b;

        if (isempty(a))     return(a);
        if (length(a) == 1) return(a);

        r   = roots(a);
        idx = find(mag(r) > 1);

        r[idx] = 1 / conj(r[idx]);

        ind = find(a != 0);
        b   = a[ind[1]] * poly(r);

        /* return real if input real */
        if (isreal(a))
        {
                b = real(b);
        }
        
        return(b);
}