* PACKAGE XMKINI   !" 初期値ファイルの作成(球面発散順圧モデル用)
*
*" [HIS]  92/06/11 
*         96/11/15 (takepiro)
*******************************************************************
      PROGRAM XMKINI
*
#if   SYS_IBMS
      INCLUDE   (ZCDIM)
      INCLUDE   (ZCCOM)
#else
#include        "zcdim.F"
#include        "zccom.F"
#endif
*
*"     ========== 変数 =========================================
*
*"     ------- 格子点データ ------------------------------------
*
      REAL       GAU   ( IDIM, JDIM, KMAX )  !" 西風  ｕ(t)
      REAL       GBU   ( IDIM, JDIM, KMAX )  !" 西風  ｕ(t-Δt)
      REAL       GAV   ( IDIM, JDIM, KMAX )  !" 南風  ｖ(t)
      REAL       GBV   ( IDIM, JDIM, KMAX )  !" 南風  ｖ(t-Δt)
      REAL       GAPHI ( IDIM, JDIM, KMAX )  !" 層厚  φ(t)
      REAL       GBPHI ( IDIM, JDIM, KMAX )  !" 層厚  φ(t-Δt)
*
*"     ========= 定数 ==========================================
*
*"     ------ 座標格子定数 -------------------------------------
      REAL       ALON  ( IDIM )              !" 経度
      REAL       DLON  ( IDIM )              !" 経度荷重
      REAL       ALAT  ( JDIM )              !" 緯度
      REAL       DLAT  ( JDIM )              !" 緯度荷重
*
*"     ==========  内部変数 ====================================
*
*"     ======= 内部定数(DSETCで確定) ============================
*
*"     ------ 球関数展開定数 -----------------------------------
*
      INTEGER    NMO   ( 2, 0:MMAX, 0:LMAX ) !" スペクトルの添字順番
      REAL       FLAPLA( NMDIM  )            !" ラプラシアンの係数
      REAL       EDEL  ( NMDIM  )            !" ζ，Ｄ→Ｕ，Ｖ
*
*"     ------ 力学定数 -----------------------------------------
*
      REAL       CORIOL ( IDIM, JDIM )       !" コリオリ係数 ｆ
      REAL       UVFACT ( IDIM, JDIM )       !" u→U のファクター
*
*"     ------ 時間積分定数 -------------------------------------
*
*"     ------- 時刻変数 ----------------------------------------
*
      INTEGER    ISTEP
      INTEGER    ITA, ITB
      INTEGER    IDATEA( 6 ) , IDATEB( 6 )
*
*"     ------ 実験管理パラメーター -----------------------------
*
      REAL       BGEP                            !" 平均ジオポテンシャル
      REAL       DGEP                            !" 擾乱ジオポテンシャル
      REAL       ALAT0, ALON0, RAD               !" 擾乱の位置と大きさ
*
      DATA       ISTEP        / 0  /
      DATA       ITA, ITB     / 0 , 0 /
      DATA       IDATEA       / 1990, 1, 1, 0, 0, 0 /
      DATA       IDATEB       / 1990, 1, 1, 0, 0, 0 /
*
      CALL YPREP
*
      WRITE( 6, * ) ' STRENGTH OF MEAN GEO-POT?'
      READ ( 5, * ) BGEP
      WRITE( 6, * ) ' STRENGTH OF DISTURBANCE OF GEO-POT?'
      READ ( 5, * ) DGEP
      WRITE( 6, * ) ' CENTER OF DIST. IN DEGREE (lat,lon) ? ' 
      READ ( 5, * ) ALAT0, ALON0
      WRITE( 6, * ) ' RADIUS OF DIST. IN DEGREE '
      READ ( 5, * ) RAD
*
*"          < 0. 定数のセット >
*
      CALL       PCONST
*
      CALL ACRGET
     O         ( ALON  , DLON  ,
     O           ALAT  , DLAT  )
*
c$$$      CALL DCSPHE     !"  力学パートのスペクトル定数セット
c$$$     O     ( NMO   , FLAPLA , EDEL  , UVFACT, 
c$$$     C       ALAT  , DLAT                     )
*
      CALL INISET
     O         ( GAU   , GAV   , GAPHI  ,
     I           BGEP  , DGEP  , ALAT0 , ALON0 , RAD  ,
     C           ALAT  , ALON                   )
*
      CALL COPY     ( GBU  , GAU  , IDIM*JDIM*KMAX )
      CALL COPY     ( GBV  , GAV  , IDIM*JDIM*KMAX )
      CALL COPY     ( GBPHI, GAPHI, IDIM*JDIM*KMAX )
*
      CALL IHSTRT
*
      CALL IWRSTR
     I         ( GAU   , GAV   , GAPHI ,
     I           GBU   , GBV   , GBPHI ,
     I           ITA   , ITB   , IDATEA, IDATEB, ISTEP )
*
      STOP
      END
**********************************************************************
      SUBROUTINE INISET
     O         ( GAU   , GAV   , GAPHI  ,
     I           BGEP  , DGEP  , ALAT0 , ALON0 , RAD  ,
     C           ALAT  , ALON                   )
*
*"           < 初期擾乱 
*"                 geo potential 強度 = DGEP 
*"                 中心 = (RIC,RJC) ; 半径 = RAD >
*
#if   SYS_IBMS
      INCLUDE   (ZCDIM)
#else
#include        "zcdim.F"
#endif
*
* [OUTPUT]
      REAL       GAU   ( IDIM, JDIM, KMAX )  !" 西風  ｕ
      REAL       GAV   ( IDIM, JDIM, KMAX )  !" 南風  ｖ
      REAL       GAPHI ( IDIM, JDIM, KMAX )  !" 層厚  φ
*
* [INPUT]
      REAL       BGEP                            !" 平均ジオポテンシャル
      REAL       DGEP                            !" 擾乱ジオポテンシャル
      REAL       ALAT0, ALON0, RAD               !" 擾乱の位置と大きさ
*
      REAL       ALON  ( IDIM )                  !" 経度
      REAL       ALAT  ( JDIM )                  !" 緯度
*
*  [INTERNAL WORK]
      INTEGER    I, J, K                         !" ループカウンター
      REAL       PI                              !" 円周率
*
*
      PI = ATAN( 1. )*4
      ALAT0 = ALAT0*PI/180.
      ALON0 = ALON0*PI/180.
      RAD   = RAD*PI/180.
*
      CALL RESET ( GAU   , IDIM*JDIM*KMAX )
      CALL RESET ( GAV   , IDIM*JDIM*KMAX )
      CALL RESET ( GAPHI , IDIM*JDIM*KMAX )
*
      DO 100 J = 1 , JDIM
        DO 100 I = 1, IDIM
           DO 100 K = 1, KMAX
              GAPHI ( I,J,K ) 
     &      =   BGEP 
     &        + DGEP * EXP( -(   ( ALAT(J) - ALAT0 )**2
     &                         + ( ALON(I) - ALON0 )**2  )/RAD**2 )
  100 CONTINUE
*
      RETURN
      END
