/*****************************************************************************
*                                                                            *
*   FFTSHIFT.SPL Copyright (C) 1997-2005 DSP Development Corporation         *
*                           All Rights Reserved                              *
*                                                                            *
*   Author:      Randy Race                                                  *
*                                                                            *
*   Synopsis:    Flips FFT so zero frequency is the middle point             *
*                                                                            *
*   Revisions:   18 Aug 1997  RRR  Creation                                  *
*                21 Jan 1998  RRR  Odd length fixup for rows                 *
*                 8 Jul 1999  RRR  don't reverse halves                      *
*                19 May 2000  RRR  optimize buffering for arrays             *
*                 9 Feb 2005  RRR  buffering now internally optimized        *
*                                                                            *
*****************************************************************************/


#if @HELP_FFTSHIFT

    FFTSHIFT

    Purpose: Shifts a 1D or 2D FFT so the 0 frequency is the midpoint

    Syntax:  FFTSHIFT(ser, dim)

               ser - A series or array.

               dim - Optional. An integer, the shift dimension.

                      -1: Shift rows and columns for arrays (default).
                       0: Do not shift quadrants for arrays (same as 1).
                       1: Shift rows only.
                       2: Shift columns only.

    Returns: A series or array.

    Example:
             W1: {1, 2, 3, 2, 1}
             W2: fft(W1)
             W3: fftshift(W2)

             The zero frequency (i.e. DC) value of W2 is the first point.
             The zero frequency of W3 is the 3rd point and appears in the
             middle of the resulting graph.

    Example:
             W1: ravel(1..16, 4)
             W2: fftshift(W1)
             W3: fftshift(W1, 1)
             W4: fftshift(W1, 2)

             W1 == {{1, 5,  9, 13},
                    {2, 6, 10, 14},
                    {3, 7, 11, 15},
                    {4, 8, 12, 16}}

             W2 == {{11, 15,  3,  7},
                    {12, 16,  4,  8},
                    { 9, 13,  1,  5},
                    {10, 14,  2,  6}}

             W3 == {{ 3,  7, 11, 15},
                    { 4,  8, 12, 16},
                    { 1,  5,  9, 13},
                    { 2,  6, 10, 14}}

             W4 == {{ 9, 13,  1,  5},
                    {10, 14,  2,  6},
                    {11, 15,  3,  7},
                    {12, 16,  4,  8}}

             W1 contains a 4x4 array.
             
             W2 shifts the rows and columns of the array. The result is
             the upper left and lower right quadrants are swapped and
             the upper right and lower left quadrands are swapped.

             W3 shifts the rows of the array only.

             W4 shifts the columns of the array only.

    Remarks:
             FFTSHIFT also works on 2D FFT array.

             By default, FFTSHIFT flips the quadrants of an array. Set
             dim to 0 or 1 to shift rows only and not shift the quadrants.

             See IFFTSHIFT to undo the shift of FFTSHIFT.

    See Also:
             Fft
             Fft2
             Ifftshift
#endif


/* shift mid point of series or array */
fftshift(s, dim)
{
        local rmid, rodd, cmid, codd, t, xv;
        local nr, nc, r1, r2, c1, c2, bufsize;

        if (argc < 2)
        {
                if (argc < 1) error(sprintf("%s - input series required", __FUNC__));

                dim = -1;
        }

        /* numrows and numcols */
        (nr, nc) = size(s);

        /* get row midpoints */
        rmid = ceil(nr / 2);
        rodd = nr % 2;

        if (nc > 1 && (dim == 2 || dim < 0) && not(isxy(s)))
        {
                /* array - get col midpoints */
                cmid = ceil(nc / 2 + 1);
                codd = (nc % 2);

                /* array indices */
                if (dim == 2)
                {
                        r1 = 1..rmid;
                        r2 = (rmid + 1)..nr;
                }
                else
                {
                        r1 = (rmid + 1)..nr;
                        r2 = 1..rmid;
                }

                c1 = cmid..nc;
                c2 = 1..(cmid - 1);

                /* flip rows and columns in one operation */
                t = {ravel(s[r1, c1], s[r1, c2]), ravel(s[r2, c1], s[r2, c2])};

                /* fixup offset for even or odd columns */
                setdeltay(t, deltay(s));
                setyoffset(t, yoffset(s) + ((codd - nc) / 2) * deltay(t));
        }
        else if (dim == 1 || dim <= 0)
        {
                /* series - flip each half of all the rows */
                t = concat(extract(s, rmid + 1, -1), extract(s, 1, rmid));

                if (isxy(s))
                {
                        xv = xvals(s);
                        xv = concat(extract(xv, rmid + 1, -1) - deltax(s) * length(s), extract(xv, 1, rmid));
                        t = xy(xv, t);
                }
        }
        else
        {
                t = s;
        }

        /* fixup offset for even or odd rows */
        if (dim <= 2 && not(isxy(s)))
        {
                setdeltax(t, deltax(s));
                setxoffset(t, xoffset(s) + ((rodd - nr) / 2) * deltax(t));
        }

        return(t);
}