*"    ルジャンドル変換ルーチンの初期化（沼口バージョン用）
*"    [HIS]  堀之内  93/09/01
*"    DSETC をもとにしたルーチン DSETCH または DSETC そのものを内部で呼ぶ
***********************************************************************
      SUBROUTINE LTINIT

#ifdef SYS_IBMS
      INCLUDE   (ZCDIM)                      !" 格子点数, 波数
      INCLUDE   (ZHDIM)                      !" 文字列文字数
#else
#include        "zcdim.F"                    !" 格子点数, 波数
#include        "zhdim.F"                    !" 文字列文字数
#endif

      REAL       ALAT  ( JDIM )        !" 緯度
      REAL       DLAT  ( JDIM )        !" 緯度荷重
*
*dsetc      REAL       CORIOL ( IDIM, JDIM )       !" コリオリ係数 ｆ
*dsetc      REAL       UVFACT ( IDIM, JDIM )       !" u→U のファクター
*dsetc      REAL       DIFV  ( NMDIM  )            !" 運動量水平拡散係数
*dsetc      REAL       DIFT  ( NMDIM  )            !" 熱，水水平拡散係数
*
      INTEGER    NMO   ( 2, 0:MMAX, 0:LMAX ) !" スペクトルの添字順番
      REAL       FLAPLA( NMDIM  )            !" ラプラシアンの係数
      REAL       EDEL  ( NMDIM  )            !" ζ，Ｄ→Ｕ，Ｖ
*
*     -- for entry LTGET --{
      CHARACTER HP   *(*)                    !" パラメーターの名前
      INTEGER   IPRDIM                       !" 配列IPARAの全寸法
*"   (INTEGER)  IPARA(*)                     !" パラメーターの内容
      INTEGER   IIPARA( 2, 0:MMAX, 0:LMAX )  !" パラメーターの内容
      REAL      FIPARA(*)                    !" パラメーターの内容
      INTEGER   J,M,L
*     }-- for entry LTGET --

      CALL ALTDLT     !"  モデル座標定数の計算 
     O         ( ALAT  , DLAT )

      CALL    DSETCH     !"  球面調和関数展開パラメーターのセット
     O         ( NMO   , FLAPLA, EDEL  ,
     C           ALAT  , DLAT                   )
*dsetc      CALL    DSETC     !"  球面調和関数展開パラメーターのセット
*dsetc     O         ( NMO   , FLAPLA, EDEL  ,
*dsetc     O           CORIOL, UVFACT, DIFV  , DIFT  ,
*dsetc     C           ALAT  , DLAT                   )

      RETURN
*======================================================================
*" パラメーター（配列）の取得。
      ENTRY LTGET 
     I         ( HP,IPRDIM,
     O           IPARA )
*
      IF ( HP .EQ. 'NMO' ) THEN
        CALL LTGETI
     I         ( HP,IPRDIM,
     O           IPARA )
      ELSE 
        CALL LTGETF
     I         ( HP,IPRDIM,
     O           IPARA )
      ENDIF
      RETURN
*======================================================================
*  when ipara is NMO
      ENTRY  LTGETI
     I         ( HP,
     D           IPRDIM,
     O           IIPARA )
      IF( IPRDIM .LT. 2*(MMAX+1)*(LMAX+1) ) THEN
        CALL MSGDMP( 'E', 'LTGETI', 'ARRAY IS INSUFFICIENT' )
      ENDIF
      DO 100 L = 0, LMAX
      DO 100 M = 0, MMAX
      DO 100 J = 1, 2
        IIPARA( J, M, L ) = NMO( J, M, L )
  100 CONTINUE
      RETURN
*======================================================================
*  when ipara is Real
      ENTRY  LTGETF
     I         ( HP,
     O           FIPARA )
      CALL MSGDMP( 'E', 'LTGETF', 'Mada dekite nainoyo gomen.' )
*-----" ALAT などを出力できるようにするつもり。--------
      RETURN
      END
***********************************************************************
      SUBROUTINE ALTDLT     !"  モデル座標定数（緯度方向）を得る
     O         ( ALAT  , DLAT )
*
#ifdef SYS_IBMS
      INCLUDE   (ZCDIM)                      !" 格子点数, 波数
      INCLUDE   (ZHDIM)                      !" 文字列文字数
#else
#include        "zcdim.F"                    !" 格子点数, 波数
#include        "zhdim.F"                    !" 文字列文字数
#endif
      REAL       ALAT  ( JDIM )              !" 緯度
      REAL       DLAT  ( JDIM )              !" 緯度荷重
      CHARACTER  HALAT  *(NCC)               !" 緯度軸名称
*
*"         < 1. 座標値のセット >
*
      CALL SETLAT
     O         ( ALAT  , DLAT  , HALAT ,
     D           JMAX  , JDIM            )
*
      RETURN
      END
***********************************************************************
**** PACKAGE DSETC  !" 力学 定数設定(水平差分)  --->
* SUBROUTINE DSETCH  
*
*"  [HIS] 90/08/31(numaguti)
*"        92/05/18(masaaki)     SETCOF: EDEL(LMAX=0)
*"        92/07/20(takepiro)    SETDIF              エラーチェック
*"        93/09/01(horinout)    DSETC ---> DSETCH
*
***********************************************************************
      SUBROUTINE DSETCH     !"  球面調和関数展開パラメーターのセット
     O         ( NMO   , FLAPLA, EDEL  ,
     C           ALAT  , DLAT                   )
*
*   [PARAM] 
#ifdef SYS_IBMS
      INCLUDE   (ZCDIM)                      !" 格子点数，波数
****      INCLUDE   (ZCCOM)                      !" 標準物理定数
#else
#include        "zcdim.F"                    !" 格子点数，波数
****#include        "zccom.F"                    !" 標準物理定数
#endif
*
      REAL ER !" 地球半径：本来は zccom.F の COMMON BLOCK にあるんだけど
      DATA ER / 6370.E+3 /
      SAVE ER
*   [OUTPUT] 
      INTEGER    NMO   ( 2, 0:MMAX, 0:LMAX ) !" スペクトルの添字順番
      REAL       FLAPLA( NMDIM  )            !" ラプラシアンの係数
      REAL       EDEL  ( NMDIM  )            !" ζ，Ｄ→Ｕ，Ｖ
*"       : 球関数展開定数
*
      REAL       ALAT  ( JDIM )              !" 緯度
      REAL       DLAT  ( JDIM )              !" 緯度荷重
*        : 座標格子定数
*
*   [INTERNAL WORK] 
*
      REAL * 8   QSINLA( JDIM )              !" sin（緯度）：倍精度
      REAL * 8   QGW   ( JDIM )              !" ガウス荷重 ：倍精度
*
*"         < 2. 球面調和関数展開パラメーター >
*
      IF ( JMAX .GE. 2 ) THEN
         CALL GAUSS
     O         ( QSINLA, QGW   ,
     I           JMAX            )
      ELSE
         QSINLA( 1 ) = SIN ( ALAT( 1 ) )
         QGW   ( 1 ) = DLAT( 1 )
      ENDIF
**          write(6,*) 'alat',alat
**          write(6,*) 'dlat',dlat
**          write(6,*) 'qsinla',qsinla
**          write(6,*) 'qgw',qgw
**          write(6,*) 'er',er
*
      CALL SPSTUP
     I         ( QSINLA, QGW   , ER    ,
     O           NMO                    )
*
      CALL SETCOF
     O         ( FLAPLA, EDEL  ,
     I           ER    , NMO     )
*
      RETURN
      END
