* PACKAGE DSPHE  !" 力学 球面調和関数変換
*
**********************************************************************
      SUBROUTINE W2G     !" 球面調和関数変換 (スペクトル→格子)
     M         ( GDATA ,
     I           WDATA ,
     I           HGRAD , HFUNC , KMAXD )
*
*   [PARAM] 
#ifdef SYS_IBMS
      INCLUDE   (ZCDIM)                      !" 格子点数，波数
#else
#include        "zcdim.F"                    !" 格子点数，波数
#endif
*
      INTEGER    KMAXD
*
*   [MODIFY] 
      REAL       GDATA ( IDIM*JDIM, KMAXD )  !" 格子点データ
*
*   [INPUT] 
      REAL       WDATA ( NMDIM, KMAXD    )   !" スペクトルデータ
      CHARACTER  HGRAD*4                     !" 微分フラグ
      CHARACTER  HFUNC*4                     !" 符号フラグ
*
*   [INTERNAL ONCE] 
      REAL       PNM   ( JMXHF*NMDIM )       !" Ｐnm ルジャンドル
      REAL       DPNM  ( JMXHF*NMDIM )       !" Ｐnm μ微分
      REAL       TRIGS ( IDIM*2 )            !" 三角関数表
      INTEGER    IFAX  ( 10 )                !" IMAX の因数分解
      INTEGER    NMOX  ( 2, 0:MMAX, 0:LMAX ) !" スペクトルの添字順番
      REAL       GWX   ( JDIM )              !" ガウス荷重
      REAL       GWDEL ( JDIM )              !" 微分用のガウス荷重
      LOGICAL    OSET                        !" 定数セットフラグ
      DATA       OSET / .FALSE. /
      SAVE       PNM   , DPNM  , TRIGS , IFAX
      SAVE       NMOX  , GWX   , GWDEL , OSET
*
*   [INTERNAL WORK] 
      REAL       ZDATA ( IDIM*JDIM, KMAX )   !" 東西スペクトル
      REAL       WORK  ( IDIM*JDIM, KMAX )   !" ワーク
*
      REAL * 8   QSINLA( JDIM )              !" sin(緯度) (倍精度
      REAL * 8   QGW   ( JDIM )              !" ガウス荷重 ：倍精度
      REAL       ER                          !" 地球半径
      INTEGER    NMO   ( 2, 0:MMAX, 0:LMAX ) !" スペクトルの添字順番
*
      REAL * 8   QPNM  ( 0:NMAX+1, 0:MMAX )  !" Ｐnm ルジャンドル
      REAL * 8   QDPNM ( 0:NMAX+1, 0:MMAX )  !" Ｐnm μ微分
*
      INTEGER    J, L, M, NMH
*
      IF ( .NOT. OSET ) THEN
         CALL MSGDMP( 'E','W2G','CONSTANTS ARE NOT SET.')
      ENDIF
*
      IF ( HGRAD(1:1) .EQ. 'Y' ) THEN
         CALL    SPW2G
     M         ( GDATA ,
     I           WDATA ,
     C           DPNM  , NMOX  , TRIGS , IFAX ,
     F           HGRAD , HFUNC ,
     D           IMAX  , JMAX  , KMAXD , IDIM  , JDIM  ,
     D           LMAX  , MMAX  , NMAX  , MINT  , NMDIM , JMXHF ,
     W           ZDATA , WORK                                   )
      ELSE
         CALL    SPW2G
     M         ( GDATA ,
     I           WDATA ,
     C           PNM   , NMOX  , TRIGS , IFAX ,
     F           HGRAD , HFUNC ,
     D           IMAX  , JMAX  , KMAXD , IDIM  , JDIM  ,
     D           LMAX  , MMAX  , NMAX  , MINT  , NMDIM , JMXHF ,
     W           ZDATA , WORK                                   )
      ENDIF
*
      RETURN
*======================================================================
      ENTRY      G2W      !"  球面調和関数変換 (格子→スペクトル) 
     M         ( WDATA ,
     I           GDATA ,
     I           HGRAD , HFUNC , KMAXD )
*
      IF ( .NOT. OSET ) THEN
         CALL MSGDMP( 'E','G2W','CONSTANTS ARE NOT SET.')
      ENDIF
*
      IF      ( HGRAD(1:1) .EQ. 'Y' ) THEN
         CALL    SPG2W
     M         ( WDATA ,
     I           GDATA ,
     C           DPNM  , NMOX  , TRIGS , IFAX  , GWDEL ,
     F           HGRAD , HFUNC ,
     D           IMAX  , JMAX  , KMAXD , IDIM  , JDIM  ,
     D           LMAX  , MMAX  , NMAX  , MINT  , NMDIM , JMXHF ,
     W           ZDATA , WORK                                   )
      ELSE IF ( HGRAD(1:1) .EQ. 'X' ) THEN
         CALL    SPG2W
     M         ( WDATA ,
     I           GDATA ,
     C           PNM   , NMOX  , TRIGS , IFAX  , GWDEL ,
     F           HGRAD , HFUNC ,
     D           IMAX  , JMAX  , KMAXD , IDIM  , JDIM  ,
     D           LMAX  , MMAX  , NMAX  , MINT  , NMDIM , JMXHF ,
     W           ZDATA , WORK                                   )
      ELSE
         CALL    SPG2W
     M         ( WDATA ,
     I           GDATA ,
     C           PNM   , NMOX  , TRIGS , IFAX  , GWX   ,
     F           HGRAD , HFUNC ,
     D           IMAX  , JMAX  , KMAXD , IDIM  , JDIM  ,
     D           LMAX  , MMAX  , NMAX  , MINT  , NMDIM , JMXHF ,
     W           ZDATA , WORK                                   )
      ENDIF
*
      RETURN
*=====================================================================
      ENTRY      SPSTUP               !" 球面調和関数準備
     I         ( QSINLA, QGW   , ER    ,
     O           NMO                     )
*
*"         < 1. FFTの準備 >
*
      IF ( IMAX .GE. 2 ) THEN
         CALL RFFTIM
     I         ( IMAX  ,
     O           TRIGS , IFAX )
      ENDIF
*
*"         < 2. ルジャンドル関数の計算 >
*
      IF ( JMAX .GE. 2 ) THEN
         CALL SETPNM
     O         ( PNM   , DPNM  , NMO   ,
     I           QSINLA,
     D           JMAX  , MMAX  , LMAX  , NMAX  , MINT  ,
     D           NMDIM , JMXHF ,
     W           QPNM  , QDPNM                          )
      ELSE
         NMH = 0
         DO 2200 M = 0, MMAX, MINT
            NMH = NMH + 1
            NMO ( 1, M, 0 ) = 2* NMH - 1
            NMO ( 2, M, 0 ) = 2* NMH
 2200    CONTINUE
      ENDIF
*
*
*"         < 3. ガウス荷重 >
*
      DO 3100 J = 1, JMAX
         GWX  ( J ) = QGW( J )
         GWDEL( J ) = QGW( J ) / ER / ( 1. - QSINLA(J)**2 )
 3100 CONTINUE
*
*"         < 4. NMO を コピー >
*
      DO 4100 M = 0, MMAX
         DO 4100 L = 0, LMAX
           NMOX ( 1,M,L ) = NMO ( 1,M,L )
           NMOX ( 2,M,L ) = NMO ( 2,M,L )
 4100 CONTINUE
*
      OSET = .TRUE.
*
      RETURN
      END
