Class axis_z_mod
In: axis/axis_z.f90

    Copyright (C) GFD Dennou Club, 2005. All rights reserved.

                                                                 !=begin

Module axis_z_mod

  * Developers: Morikawa Yasuhiro
  * Version: $Id: axis_z.f90,v 1.9 2005/01/19 08:52:24 morikawa Exp $
  * Tag Name: $Name:  $
  * Change History:

Overview

This module set axis Z or axis Altitude. Z 軸または高度軸を設定する。

Error Handling

Known Bugs

  • netCDF データから Z 軸を入力する axis_z_netcdf にて、
 および axis_z_half_netcdf にて、元データの units によらず
 データを入力しているが、本来は判定すべき。
 (特に圧力座標を扱う場合には必須となるはず)。
  • 現在は半整数レベルもこのモジュールで読み込んでいるが、
 別モジュールにすべきかも知れない。(「バグ情報」と違うが)。

Note

Future Plans

                                                                 !=end

Methods

Included Modules

type_mod axis_type_mod type_mod nmlfile_mod grid_3d_mod axis_type_mod constants_mod gt4_history dc_types dc_url dc_string dc_trace dc_message axis_type_mod dc_trace axis_type_mod gt4_history dc_url dc_trace axis_type_mod dc_trace type_mod axis_type_mod gt4_history dc_url dc_trace axis_type_mod dc_trace type_mod axis_type_mod gt4_history dc_url dc_trace type_mod constants_mod dc_string dc_trace dc_trace

Public Instance methods

begin

Return Half Level Data of axis Z from NAMELIST

NAMELIST から代入された半整数レベルのデータを半整数レベルの Z 軸データとして返す。

((< axis_z_init >)) の NAMELIST axis_z_half_nml の decision 変数で (({ ‘manual’ })) 以外が与えられた場合は値を代入しないで返す。

[Source]

subroutine axis_z_half_manual (DimHalf)

    use axis_type_mod, only: axis_type_copy
    use dc_trace,   only: DbgMessage, BeginSub, EndSub
                                                                 !=end
    implicit none
                                                                 !=begin
    !==== In/Out
    !
    type(AXISINFO), intent(inout) :: DimHalf        ! 次元情報を包括する変数
                                                                 !=end
    character(STRING), parameter:: subname = "axis_z_half_manual"
  continue
    !-----------------------------------------------------------------
    !   Check Initialization
    !-----------------------------------------------------------------
    call BeginSub( subname )
    if (.not. axis_z_initialized) then
       call EndSub( subname, 'Call axis_z_init before call %c', &
            &       c1=trim(subname) )
       return
    endif

    !----------------------------------------------------------------
    !   decision が 'manual' でない場合は停止して返す。
    !----------------------------------------------------------------
    if (.not. axis_z_half_data_from_manual) then
       call EndSub( subname, &
            & 'In axis_z_half_nml, configurated to not generate from manual')
       return
    endif

    !----------------------------------------------------------------
    !   intent(out) の引数にデータを渡す。
    !----------------------------------------------------------------
    call axis_type_copy(z_DimHalf, DimHalf)
    call EndSub(subname)

end subroutine

begin

Return Half Level Data of axis Z from netCDF file

netCDF データから取得した半整数レベルデータを半整数レベルの Z 軸のデータとして返す。

((< axis_z_init >)) の NAMELIST axis_z_half_nml の decision 変数で gtool4 変数以外が与えられた場合は値を代入しないで返す。

[Source]

subroutine axis_z_half_netcdf (DimHalf)

  !==== Dependency
    use type_mod,      only: STRING
    use axis_type_mod, only: axis_type_copy
    use gt4_history,only: HistoryGet
    use dc_url  ,   only: UrlSplit
    use dc_trace,   only: DbgMessage, BeginSub, EndSub
                                                                 !=end
    implicit none
                                                                 !=begin
    !==== In/Out
    !
    type(AXISINFO), intent(inout) :: DimHalf        ! 次元情報を包括する変数
                                                                 !=end

    character(STRING)           :: file
    character(STRING)           :: varname

    character(STRING), parameter:: subname = "axis_z_half_netcdf"
  continue
    !-----------------------------------------------------------------
    !   Check Initialization
    !-----------------------------------------------------------------
    call BeginSub( subname )
    if (.not. axis_z_initialized) then
       call EndSub( subname, 'Call axis_z_init before call %c', &
            &       c1=trim(subname) )
       return
    endif

    !----------------------------------------------------------------
    !   decision が '***@lon' などでない場合は停止して返す。
    !----------------------------------------------------------------
    if (axis_z_half_data_from_netcdf == '') then
       call EndSub( subname, &
            & 'In axis_z_half_nml, configurated to not generate from netcdf data')
       return
    endif

    !----------------------------------------------------------------
    !   HistoryGet によるデータ入力
    !----------------------------------------------------------------
    call UrlSplit(axis_z_half_data_from_netcdf, file=file, var=varname)
    call HistoryGet(file, varname, z_DimHalf%a_Dim)
    z_DimHalf%stored   = .true.

    !----------------------------------------------------------------
    !   intent(out) の引数にデータを渡す。
    !----------------------------------------------------------------
    call axis_type_copy(z_DimHalf, DimHalf)
    call EndSub(subname)

end subroutine

begin

Return Data of axis Z from NAMELIST

NAMELIST から代入されたデータを Z 軸データとして返す。

((< axis_z_init >)) の NAMELIST axis_z_nml の decision 変数で (({ ‘manual’ })) 以外が与えられた場合は値を代入しないで返す。

[Source]

subroutine axis_z_manual (Dim)

    use axis_type_mod, only: axis_type_copy
    use dc_trace,   only: DbgMessage, BeginSub, EndSub
                                                                 !=end
    implicit none
                                                                 !=begin
    !==== In/Out
    !
    type(AXISINFO), intent(inout) :: Dim        ! 次元情報を包括する変数
                                                                 !=end
    character(STRING), parameter:: subname = "axis_z_manual"
  continue
    !-----------------------------------------------------------------
    !   Check Initialization
    !-----------------------------------------------------------------
    call BeginSub( subname )
    if (.not. axis_z_initialized) then
       call EndSub( subname, 'Call axis_z_init before call %c', &
            &       c1=trim(subname) )
       return
    endif

    !----------------------------------------------------------------
    !   decision が 'manual' でない場合は停止して返す。
    !----------------------------------------------------------------
    if (.not. axis_z_data_from_manual) then
       call EndSub( subname, &
            & 'In axis_z_nml, configurated to not generate from manual')
       return
    endif

    !----------------------------------------------------------------
    !   intent(out) の引数にデータを渡す。
    !----------------------------------------------------------------
    call axis_type_copy(z_Dim, Dim)
    call EndSub(subname)

end subroutine

begin

Return Data of axis Z from netCDF file

netCDF データから取得したデータを Z 軸のデータとして返す。

((< axis_z_init >)) の NAMELIST axis_z_nml の decision 変数で gtool4 変数以外が与えられた場合は値を代入しないで返す。

[Source]

subroutine axis_z_netcdf (Dim)

  !==== Dependency
    use type_mod,      only: STRING
    use axis_type_mod, only: axis_type_copy
    use gt4_history,only: HistoryGet
    use dc_url  ,   only: UrlSplit
    use dc_trace,   only: DbgMessage, BeginSub, EndSub
                                                                 !=end
    implicit none
                                                                 !=begin
    !==== In/Out
    !
    type(AXISINFO), intent(inout) :: Dim        ! 次元情報を包括する変数
                                                                 !=end

    character(STRING)           :: file
    character(STRING)           :: varname

    character(STRING), parameter:: subname = "axis_z_netcdf"
  continue
    !-----------------------------------------------------------------
    !   Check Initialization
    !-----------------------------------------------------------------
    call BeginSub( subname )
    if (.not. axis_z_initialized) then
       call EndSub( subname, 'Call axis_z_init before call %c', &
            &       c1=trim(subname) )
       return
    endif

    !----------------------------------------------------------------
    !   decision が '***@lon' などでない場合は停止して返す。
    !----------------------------------------------------------------
    if (axis_z_data_from_netcdf == '') then
       call EndSub( subname, &
            & 'In axis_z_nml, configurated to not generate from netcdf data')
       return
    endif

    !----------------------------------------------------------------
    !   HistoryGet によるデータ入力
    !----------------------------------------------------------------
    call UrlSplit(axis_z_data_from_netcdf, file=file, var=varname)
    call HistoryGet(file, varname, z_Dim%a_Dim)
    z_Dim%stored   = .true.

    !----------------------------------------------------------------
    !   intent(out) の引数にデータを渡す。
    !----------------------------------------------------------------
    call axis_type_copy(z_Dim, Dim)
    call EndSub(subname)

end subroutine

begin

Return Sigma Full Level and Sigma Half Level (from NAMELIST) as Data of axis Z

NAMELIST から取得した半整数σレベルのデータと、 その半整数σレベルから生成した整数σレベルデータを返す。

((< axis_z_init >)) の NAMELIST axis_z_nml の decision 変数で (({ ‘sigmahalf’ })) が与えられ、且つ NAMELIST axis_z_half_nml の decision 変数で (({ ‘manual’ }))が与えられている場合以外は、値に何も代入せず返す。

[Source]

subroutine axis_z_sigmahalf_manual (Dim, DimHalf)

  !==== Dependency
    use axis_type_mod, only: axis_type_copy
    use dc_trace,   only: DbgMessage, BeginSub, EndSub
                                                                 !=end
    implicit none
                                                                 !=begin
    !==== In/Out
    !
    type(AXISINFO), intent(inout) :: Dim, DimHalf  ! 次元情報を包括する変数
                                                                 !=end
    character(STRING), parameter:: subname = "axis_z_sigmahalf_manual"
  continue
    !-----------------------------------------------------------------
    !   Check Initialization
    !-----------------------------------------------------------------
    call BeginSub( subname )
    if (.not. axis_z_initialized) then
       call EndSub( subname, 'Call axis_z_init before call %c', &
            &       c1=trim(subname) )
       return
    endif

    !----------------------------------------------------------------
    !   axis_z_nml の decision が 'sigmahalf' でなく、
    !   axis_z_sigmahalf_nml の decision が 'manual' でない場合は停止
    !----------------------------------------------------------------
    if (   .not. axis_z_data_from_sigmahalf .or. &
         & .not. axis_z_half_data_from_manual ) then
       call EndSub( subname, &
            & 'In axis_z_nml and in axis_z_sigmahalf_nml, configurated to not generate from sigmahalf (which is generated from NAMELIST)')
       return
    endif

    !----------------------------------------------------------------
    !   DimHalf の値から Dim を作成。
    !----------------------------------------------------------------
    call gen_sigma_from_half(z_DimHalf, z_Dim)

    !----------------------------------------------------------------
    !   intent(out) の引数にデータを渡す。
    !----------------------------------------------------------------
    call axis_type_copy(z_Dim, Dim)
    call axis_type_copy(z_DimHalf, DimHalf)
    call EndSub(subname)

end subroutine

begin

Return Sigma Full Level and Sigma Half Level (from netCDF Data) as Data of axis Z

netCDF データから取得した半整数σレベルのデータと、 その半整数σレベルから生成した整数σレベルデータを返す。

((< axis_z_init >)) の NAMELIST axis_z_nml の decision 変数で (({ ‘sigmahalf’ })) が与えられ、且つ NAMELIST axis_z_half_nml の decision 変数で gtool4変数が与えられている場合以外は、値に何も代入せず返す。

[Source]

subroutine axis_z_sigmahalf_netcdf (Dim, DimHalf)

  !==== Dependency
    use axis_type_mod, only: axis_type_copy
    use gt4_history,only: HistoryGet
    use dc_url  ,   only: UrlSplit
    use dc_trace,   only: DbgMessage, BeginSub, EndSub
                                                                 !=end
    implicit none
                                                                 !=begin
    !==== In/Out
    !
    type(AXISINFO), intent(inout) :: Dim, DimHalf  ! 次元情報を包括する変数
                                                                 !=end
    character(STRING)           :: file
    character(STRING)           :: varname

    character(STRING), parameter:: subname = "axis_z_sigmahalf_netcdf"
  continue
    !-----------------------------------------------------------------
    !   Check Initialization
    !-----------------------------------------------------------------
    call BeginSub( subname )
    if (.not. axis_z_initialized) then
       call EndSub( subname, 'Call axis_z_init before call %c', &
            &       c1=trim(subname) )
       return
    endif

    !----------------------------------------------------------------
    !   axis_z_nml の decision が 'sigmahalf' でなく、
    !   axis_z_sigmahalf_nml の decision が '***@lon' でない時は停止
    !----------------------------------------------------------------
    if (   .not. axis_z_data_from_sigmahalf .or. &
         & .not. axis_z_half_data_from_netcdf /= '') then
       call EndSub( subname, &
            & 'In axis_z_nml and in axis_z_sigmahalf_nml, configurated to not generate from sigmahalf (which is generated from netcdf data)')
       return
    endif

    !----------------------------------------------------------------
    !   HistoryGet によるデータ入力
    !----------------------------------------------------------------
    call UrlSplit(axis_z_half_data_from_netcdf, file=file, var=varname)
    call HistoryGet(file, varname, z_DimHalf%a_Dim)
    z_Dim%stored   = .true.

    !----------------------------------------------------------------
    !   DimHalf の値から Dim を作成。
    !----------------------------------------------------------------
    call gen_sigma_from_half(z_DimHalf, z_Dim)

    !----------------------------------------------------------------
    !   intent(out) の引数にデータを渡す。
    !----------------------------------------------------------------
    call axis_type_copy(z_Dim, Dim)
    call axis_type_copy(z_DimHalf, DimHalf)
    call EndSub(subname)

end subroutine

begin

Generate Sigma Level from Sigma Half Level.

axis_z_mod の内部サブルーチン。

引数 DimHalf を半整数σレベルとして、その値から整数σレベルを生成し、 引数 Dim として返す。

[Source]

subroutine gen_sigma_from_half (DimHalf, Dim)

  !==== Dependency
    use type_mod,   only: DBKIND, STRING, REKIND
    use constants_mod,only: RAir, Cp
    use dc_string,  only: toChar
    use dc_trace,   only: DbgMessage, BeginSub, EndSub
                                                                 !=end
    implicit none
                                                                 !=begin
    !==== Input
    !
    type(AXISINFO), intent(in)   :: DimHalf ! 半整数σレベル

    !==== In/Out
    !
    type(AXISINFO), intent(inout):: Dim     ! 整数σレベル
                                                                 !=end
    !----- 作業用内部変数 -----
    real(REKIND), allocatable :: z_DelSigma(:) ! Δσ(整数レベル)
    real(DBKIND)              :: Kappa  ! κ=R/Cp (気体定数/定圧比熱)

    ! do ループ用作業変数
    integer(INTKIND)   :: k

    character(STRING), parameter:: subname = "gen_sigma_from_half"
  continue
    !-----------------------------------------------------------------
    !   Check Initialization
    !-----------------------------------------------------------------
    call BeginSub( subname )
    if (.not. axis_z_initialized) then
       call EndSub( subname, 'Call axis_z_init before call %c', &
            &       c1=trim(subname) )
       return
    endif

    !-----------------------------------------------------------------
    !   z_DelSigma (Δσ) の生成
    !-----------------------------------------------------------------
    allocate( z_DelSigma(Dim%axisinfo%length) )

    do k = 1, Dim%axisinfo%length
       z_DelSigma(k) = z_DimHalf%a_Dim(k) - z_DimHalf%a_Dim(k+1)
    enddo

    call DbgMessage( 'DimHalf-a_Dim(:)=<%c>', c1=trim( toChar(DimHalf%a_Dim(:)) ) )
    call DbgMessage( 'Dim-length=<%c>', c1=trim( toChar(Dim%axisinfo%length) )  )
    call DbgMessage( 'z_DelSigma(:)=<%c>', c1=trim( toChar(z_DelSigma(:)) ) )

    !-----------------------------------------------------------------
    !   κ=R/Cp (気体定数/定圧比熱) の生成
    !-----------------------------------------------------------------
    Kappa = RAir / Cp

    call DbgMessage( 'Kappa=<%c>', c1=trim( toChar(Kappa) ) )

    !-----------------------------------------------------------------
    !   整数レベルを作成
    !-----------------------------------------------------------------
    if (.not. allocated(Dim%a_Dim) ) then
       allocate( Dim%a_Dim(Dim%axisinfo%length) )
    endif

    do k = 1, Dim%axisinfo%length
       Dim%a_Dim(k)= (                                             &
            &           (  z_DimHalf%a_Dim ( k )**( 1. + Kappa )   &
            &            - z_DimHalf%a_Dim (k+1)**( 1. + Kappa ) ) &
            &                                                      &
            &           / (   z_DelSigma ( k ) *  ( 1. + Kappa ) ) &
            &       )                                              &
            &        **( 1. /Kappa )
    enddo

    call DbgMessage( 'z_Sigma(:)=<%c>', c1=trim( toChar(Dim%a_Dim(:)) ) )
    call EndSub(subname)

end subroutine

[Validate]