!-------------------------------------- LICENCE BEGIN ------------------------------------
!Environment Canada - Atmospheric Science and Technology License/Disclaimer,
!                     version 3; Last Modified: May 7, 2008.
!This is free but copyrighted software; you can use/redistribute/modify it under the terms
!of the Environment Canada - Atmospheric Science and Technology License/Disclaimer
!version 3 or (at your option) any later version that should be found at:
!http://collaboration.cmc.ec.gc.ca/science/rpn.comm/license.html
!
!This software is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY;
!without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.
!See the above mentioned License/Disclaimer for more details.
!You should have received a copy of the License/Disclaimer along with this software;
!if not, you can write to: EC-RPN COMM Group, 2121 TransCanada, suite 500, Dorval (Quebec),
!CANADA, H9P 1J3; or send e-mail to service.rpn@ec.gc.ca
!-------------------------------------- LICENCE END --------------------------------------
!

      SUBROUTINE SPGD 38,6
#if defined (DOC)
*
**s/r SPGD  - Inverse spectral transform
*
* Purpose: builds the grid point model state from the spectral
*          state.
*
*Author  : P.Gauthier/P.Koclas  - *ARMA and CMC  Nov 1997
*     .    Complete rewriting of the transform for a parallel treatment
*     .    based on having the outer loop on the zonal wavenumber
*
*revision: Luc Fillion/Christian Page - Summer 2003
*               - Include NCAR 2D-FFT for Limited-Area analysis
*revision: Luc Fillion - 07 Jul 2004 - Update to NCAR's new fftpack5.
*revision: Luc Fillion - oct 2004 - Add RPN's 2D FFT option.
*revision: Luc Fillion - ARMA/EC - 14 Aug 2007 - Update LAM4D to v_10_0_3.
*
*
*
#endif
      use modstag, only: level2_staggrid
      IMPLICIT NONE
*implicits
#include "pardim.cdk"
#include "comcva.cdk"
#include "comdim.cdk"
#include "comgrd_param.cdk"
#include "comgem.cdk"
#include "comgd0.cdk"
#include "comsp.cdk"
#include "comfftla.cdk"
C
C
      INTEGER  ji,jj,jla,JGL, JK, JLON, ILONMAX
      real*8 zfmla(ni,nksdim,nj), zsp(nla,2)
      real*8 zgdxy(ni,nj)
C
C
C*    1. Inverse Legendre transform
C        --------------------------
 200  CONTINUE
      if (grd_typ.eq.'GU') then                          !
        CALL SPPAR(level2_staggrid)
      endif
                                !
C
C     2.1 Reset to zero the modes that are not part of the truncation
C         -----------------------------------------------------------
C
      if (grd_typ.eq.'GU') then                          !
        ILONMAX = NIEND
        DO JGL = 1, NJ
           DO JLON = 2*(NTRUNC+1)+1, ILONMAX
              DO  JK = 1, NKGDIM
                 GD(JLON,JK,JGL) = 0.
              END DO
           END DO
        END DO
C
        DO JGL =1,NJ
           DO JK = 1, NKGDIM
              GD(0,JK,JGL) = 0.
           END DO
        END DO
      endif
C
C     2.2 Apply the FFT 
C         -------------
C
      if (grd_typ.eq.'GU') then                          !
        CALL FFT3DVAR(GD,NIBEG,NIEND,NKGDIM,NJBEG,NJEND,NKSDIM,+1)
      else if (grd_typ.eq.'LU') then
        do jk = 1, nksdim
          do jla = 1, nla
            zsp(jla,1) = sp(jla,1,jk)
            zsp(jla,2) = sp(jla,2,jk)
          enddo
          if(lrpnfft) then
            call zero(ni*nj,zgdxy)
            call idft2dr(zgdxy,zsp)
            do jj=1,nj
              do ji=1,ni
                gd(ji,jk,jj) = zgdxy(ji,jj)
              enddo
            enddo
          else
!            call idft2d(zfmla,zsp,sdft2d,nindxy,jk,
!     &                   mlen2d,nla,ni,nj,nksdim)
            do jj=1,nj
              do ji=1,ni
                gd(ji,jk,jj) = zfmla(ji,jk,jj)
              enddo
            enddo
          endif
        enddo
        if(nanalvar.eq.4) then
          call sp2rela(tb0,sptb,nflev)
        endif
      endif
C
      END SUBROUTINE SPGD