Class | axis_z_mod |
In: |
axis/axis_z.f90
|
Copyright (C) GFD Dennou Club, 2005. All rights reserved.
!=begin
* Developers: Morikawa Yasuhiro * Version: $Id: axis_z.f90,v 1.9 2005/01/19 08:52:24 morikawa Exp $ * Tag Name: $Name: $ * Change History:
This module set axis Z or axis Altitude. Z 軸または高度軸を設定する。
および axis_z_half_netcdf にて、元データの units によらず データを入力しているが、本来は判定すべき。 (特に圧力座標を扱う場合には必須となるはず)。
別モジュールにすべきかも知れない。(「バグ情報」と違うが)。
!=end
NAMELIST から代入された半整数レベルのデータを半整数レベルの Z 軸データとして返す。
((< axis_z_init >)) の NAMELIST axis_z_half_nml の decision 変数で (({ ‘manual’ })) 以外が与えられた場合は値を代入しないで返す。
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
netCDF データから取得した半整数レベルデータを半整数レベルの Z 軸のデータとして返す。
((< axis_z_init >)) の NAMELIST axis_z_half_nml の decision 変数で gtool4 変数以外が与えられた場合は値を代入しないで返す。
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
NAMELIST から代入されたデータを Z 軸データとして返す。
((< axis_z_init >)) の NAMELIST axis_z_nml の decision 変数で (({ ‘manual’ })) 以外が与えられた場合は値を代入しないで返す。
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
netCDF データから取得したデータを Z 軸のデータとして返す。
((< axis_z_init >)) の NAMELIST axis_z_nml の decision 変数で gtool4 変数以外が与えられた場合は値を代入しないで返す。
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
NAMELIST から取得した半整数σレベルのデータと、 その半整数σレベルから生成した整数σレベルデータを返す。
((< axis_z_init >)) の NAMELIST axis_z_nml の decision 変数で (({ ‘sigmahalf’ })) が与えられ、且つ NAMELIST axis_z_half_nml の decision 変数で (({ ‘manual’ }))が与えられている場合以外は、値に何も代入せず返す。
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
netCDF データから取得した半整数σレベルのデータと、 その半整数σレベルから生成した整数σレベルデータを返す。
((< axis_z_init >)) の NAMELIST axis_z_nml の decision 変数で (({ ‘sigmahalf’ })) が与えられ、且つ NAMELIST axis_z_half_nml の decision 変数で gtool4変数が与えられている場合以外は、値に何も代入せず返す。
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
axis_z_mod の内部サブルーチン。
引数 DimHalf を半整数σレベルとして、その値から整数σレベルを生成し、 引数 Dim として返す。
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