!--------------------------------------------------------------------- ! Copyright (C) GFD Dennou Club, 2005. All rights reserved. !--------------------------------------------------------------------- !=begin != Module axis_x_mod ! ! * Developers: Morikawa Yasuhiro ! * Version: $Id: axis_x.f90,v 1.4 2006/09/29 14:19:30 morikawa Exp $ ! * Tag Name: $Name: $ ! * Change History: ! !== Overview ! !This module set axis X or axis Longitude. ! !X 軸または経度軸を設定する。 ! !== Error Handling ! !== Known Bugs ! !netCDF データから X 軸を入力する axis_x_netcdf にて、 !元データが radians でも degrees でも、そのまま入力されるように !なっている。本来は元データの units から判定すべき。 ! !== Note ! !== Future Plans ! ! !=end module axis_x_mod !=begin !== Dependency use type_mod, only : INTKIND, STRING use axis_type_mod, only : AXISINFO !=end implicit none !=begin !== Public Interface private public :: axis_x_init, axis_x_weight, axis_x_spectral ! subroutines public :: axis_x_manual, axis_x_netcdf, axis_x_end ! subroutines !=end type(AXISINFO) ,save:: x_Dim ! 次元情報を包括する変数 type(AXISINFO) ,save:: x_Dim_Weight ! 次元重み情報を包括する変数 logical ,save:: axis_x_data_from_spectral = .false. logical ,save:: axis_x_data_from_manual = .false. character(STRING),save:: axis_x_data_from_netcdf = '' logical ,save:: axis_x_initialized = .false. character(STRING),parameter:: version = & & '$Id: axis_x.f90,v 1.4 2006/09/29 14:19:30 morikawa Exp $' character(STRING),parameter:: tagname = '$Name: $' contains !=begin !== Procedure Interface ! !=== Initialize module and acquire NAMELIST ! !モジュールを初期化し、NAMELIST から値を取得する。 !NAMELIST から値が取得できないものに関しては上記のデフォルト値が !用いられる。 ! !NAMELIST ファイルは、メインプログラムにて ((< nmlfile_mod >)) の !((< nmlfile_init >)) で指定されることが想定されているが、 !もしもこの初期化ルーチンより以前に指定されていなければ、 !((< nmlfile_init >)) のデフォルトで指定される NAMELIST ファイルを !読む。 ! subroutine axis_x_init !==== Dependency use type_mod, only: STRING, TOKEN, INTKIND, REKIND, DBKIND, NMLARRAY use nmlfile_mod,only: nmlfile_init, nmlfile_open, nmlfile_close use grid_3d_mod,only: im, grid_3d_init use spml_mod, only: spml_init ! use axis_type_mod, only : axis_attrs_copy, axis_attrs_init ! use gt4_history, only: GT_HISTORY_ATTR use dc_types, only: GT_TOKEN => TOKEN, GT_STRING => STRING use dc_url, only: GT_ATMARK, GT_QUESTION use dc_trace, only: DbgMessage, BeginSub, EndSub use dc_message,only: MessageNotify use gt4_history, only: HistoryAxisCreate, HistoryAxisAddAttr use dc_string, only: StrHead !=end implicit none !------------------------------------------------------------------- ! 変数定義 !------------------------------------------------------------------- !=begin ! !==== NAMELIST ! !X 軸の次元変数に関する情報を与える。 !値を与えないものに関しては以下のデフォルトの値が用いられる。 ! !変数 decision には X 軸のデータをどのように与えるかを指定する。 ! !* (({ 'manual' })) ! * Data 配列に格納したデータをそのまま X 軸として与える。 ! !* (({ 'spectral' })) ! * スペクトル法を用いる事を想定し、 length に応じて自動的に ! X 軸のデータが決まる。 ! !* gtool4 変数 (例えば (({ 'foo.nc@lon' })) など) ! * 該当する変数から X 軸のデータを取得する。 ! !変数 length には、 (()) の公開要素 ((< im >)) と同じ !値を与えなければならない。 ! character(TOKEN) :: name = 'lon' ! 次元変数名 integer(INTKIND) :: length = 64 ! 次元長 (配列サイズ) character(STRING) :: longname = 'Longitude' ! 次元変数の記述的名称 character(STRING) :: units = 'degrees_east' ! 次元変数の単位 character(TOKEN) :: xtype = 'float' ! 次元変数の型 character(STRING) :: decision = 'spectral' ! 次元データの取得方法 real(REKIND) :: Data(NMLARRAY) = 0.0 ! 次元データ入力用 namelist /axis_x_nml/ & & name , & ! 次元変数名 & length , & ! 次元長 (配列サイズ) & longname , & ! 次元変数の記述的名称 & units , & ! 次元変数の単位 & xtype , & ! 次元変数の型 & decision , & ! 次元データの取得方法 & Data ! 次元データ !=end !=begin ! !X 軸の次元変数の属性に関する情報を与える。 !NAMELIST に複数の axis_x_attr_nml を用意しておく事で !複数の情報を与える事が可能である。 !与えない場合には属性情報は付加されない。 ! !attrtype には与える属性値の種類を設定する。 !(()) !を参照せよ。なお、arraysize に 1 以上の値を設定すると、 !配列データが優先されて属性値に設定される。 ! character(GT_TOKEN) :: attrname = '' ! 属性名 character(GT_TOKEN) :: attrtype = '' ! 属性値の型 character(GT_STRING) :: cvalue = '' ! 属性の値 (文字) integer(INTKIND) :: ivalue = 0 ! 属性の値 (整数) real(REKIND) :: rvalue = 0.0 ! 属性の値 (単精度実数) real(DBKIND) :: dvalue = 0.0d0 ! 属性の値 (倍精度実数) logical :: lvalue = .false.! 属性の値 (論理) integer(INTKIND) :: arraysize= 0 ! 配列のサイズ integer(INTKIND) :: iarray(NMLARRAY) = 0 ! 属性の値 (整数) real(REKIND) :: rarray(NMLARRAY) = 0.0 ! 属性の値 (単精度実数) real(DBKIND) :: darray(NMLARRAY) = 0.0d0! 属性の値 (倍精度実数) namelist /axis_x_attr_nml/ & & attrname , & ! 属性名 & attrtype , & ! 属性値の型 & cvalue , & ! 属性の値 (文字) & ivalue , & ! 属性の値 (整数) & rvalue , & ! 属性の値 (単精度実数) & dvalue , & ! 属性の値 (倍精度実数) & lvalue , & ! 属性の値 (論理) & arraysize , & ! 配列のサイズ & iarray , & ! 属性の値 (整数) & rarray , & ! 属性の値 (単精度実数) & darray ! 属性の値 (倍精度実数) !=end !----------------------------------------------------------------- ! 変数情報の一時格納変数 !----------------------------------------------------------------- ! type(GT_HISTORY_ATTR), allocatable :: attrs_tmp(:) !----------------------------------------------------------------- ! 汎用変数 !----------------------------------------------------------------- integer(INTKIND) :: i, k integer(INTKIND) :: nmlstat, nmlunit logical :: nmlreadable, next character(TOKEN) :: position character(STRING), parameter:: subname = "axis_x_init" continue !---------------------------------------------------------------- ! Check Initialization !---------------------------------------------------------------- call BeginSub(subname) if (axis_x_initialized) then call EndSub( subname, '%c is already called', c1=trim(subname) ) return else axis_x_initialized = .true. endif !---------------------------------------------------------------- ! Version identifier !---------------------------------------------------------------- call DbgMessage('%c :: %c', c1=trim(version), c2=trim(tagname)) !----------------------------------------------------------------- ! Initialize Dependent Modules !----------------------------------------------------------------- call grid_3d_init call spml_init !---------------------------------------------------------------- ! read axis_x_nml !---------------------------------------------------------------- call nmlfile_init call nmlfile_open(nmlunit, nmlreadable) if (nmlreadable) then read(nmlunit, nml=axis_x_nml, iostat=nmlstat) call DbgMessage('Stat of NAMELIST axis_x_nml Input is <%d>', & & i=(/nmlstat/)) write(0, nml=axis_x_nml) else call DbgMessage('Not Read NAMELIST axis_x_nml') call MessageNotify('W', subname, & & 'Can not Read NAMELIST axis_x_nml. Force Use Default Value.') end if call nmlfile_close x_Dim%stored = .false. ! 次元変数の情報を構造型変数 x_Dim への代入 call HistoryAxisCreate(x_Dim % axisinfo, & & name, length, longname, units, xtype) !!$ x_Dim%axisinfo%name = name !!$ x_Dim%axisinfo%length = length !!$ x_Dim%axisinfo%longname = longname !!$ x_Dim%axisinfo%units = units !!$ x_Dim%axisinfo%xtype = xtype !!$ allocate( x_Dim%a_Dim(x_Dim%axisinfo%length) ) allocate( x_Dim%a_Dim(length) ) ! 次元変数の情報を構造型変数 x_Dim への代入 select case(decision) ! manual: NAMELIST の Data で入力 case('manual') x_Dim%a_Dim(:) = 0 x_Dim%a_Dim(1:length) = Data(1:length) x_Dim%stored = .true. axis_x_data_from_manual = .true. ! spectral: SPMODEL の wa_module により生成 case('spectral') axis_x_data_from_spectral = .true. x_Dim%stored = .false. ! foo.nc@lon: foo.nc ファイルの lon 変数から取得 ! その他 : spectral と同じに case default ! 文字の中に '@' か '?' が含まれる場合は gtool4 変数として ! 認識し、その nc ファイルから変数情報をコピーする。 if ( index(decision, GT_ATMARK) > 0 .or. & & index(decision, GT_QUESTION) > 0) then axis_x_data_from_netcdf = decision x_Dim%stored = .false. ! それ以外は 'spectral' と同じように処理 else axis_x_data_from_spectral = .true. x_Dim%stored = .false. endif end select !---------------------------------------------------------------- ! read axis_x_attr_nml !---------------------------------------------------------------- call nmlfile_init call nmlfile_open(nmlunit, nmlreadable) if (.not. nmlreadable) then call DbgMessage('Not Read NAMELIST axis_x_attr_nml') call MessageNotify('W', subname, & & 'Can not Read NAMELIST axis_x_attr_nml.') else i = 0 next = .false. axis_x_attr_nml_input : do i = i + 1 call DbgMessage('NAMELIST axis_x_attr_nml Input, <%d> time', i=(/i/)) ! 初期化 attrname = '' ! 属性名 attrtype = '' ! 属性値の型 cvalue = '' ! 属性の値 (文字) ivalue = 0 ! 属性の値 (整数) rvalue = 0.0 ! 属性の値 (単精度実数) dvalue = 0.0d0 ! 属性の値 (倍精度実数) lvalue = .false.! 属性の値 (論理) arraysize = 0 ! 配列のサイズ iarray(:) = 0 ! 属性の値 (整数) rarray(:) = 0.0 ! 属性の値 (単精度実数) darray(:) = 0.0d0 ! 属性の値 (倍精度実数) ! read nml read(nmlunit, nml=axis_x_attr_nml, iostat=nmlstat) call DbgMessage('Stat of NAMELIST axis_x_attr_nml Input is <%d>', & & i=(/nmlstat/)) write(0, nml=axis_x_attr_nml) ! Inquire access position inquire(nmlunit, position=position) if ( trim(position) /= 'APPEND' .and. nmlstat == 0 ) then next = .true. else next = .false. endif ! 有効でない値を含むものに関しては無視。 if (attrname == '') then call DbgMessage('attrname is blank. so this axis_x_attr_nml is ignored.') if (next) cycle if (.not. next) exit elseif (attrtype == '') then call DbgMessage('attrtype is blank. so this axis_x_attr_nml is ignored.') if (next) cycle if (.not. next) exit endif !----------------------------------------------------------- ! x_Dim%attrs への格納 !----------------------------------------------------------- ! attrs(:) の拡張 !!$ if ( .not. associated(x_Dim%attrs) ) then !!$ allocate( x_Dim%attrs(1) ) !!$ k = 1 !!$ else !!$ k = size( x_Dim%attrs ) + 1 !!$ ! 配列データの領域確保 !!$ allocate( attrs_tmp(k-1) ) !!$ call axis_attrs_copy(from=x_Dim%attrs(1:k-1), to=attrs_tmp(1:k-1)) !!$ deallocate( x_Dim%attrs ) !!$ allocate( x_Dim%attrs(k) ) !!$ call axis_attrs_copy(from=attrs_tmp(1:k-1), to=x_Dim%attrs(1:k-1)) !!$ deallocate( attrs_tmp ) !!$ endif !!$ !!$ if (arraysize > 0) then !!$ call axis_attrs_init(x_Dim%attrs(k)) !!$ !!$ deallocate( x_Dim%attrs(k)%iarray ) !!$ deallocate( x_Dim%attrs(k)%rarray ) !!$ deallocate( x_Dim%attrs(k)%darray ) !!$ !!$ allocate( x_Dim%attrs(k)%iarray( arraysize ) ) !!$ allocate( x_Dim%attrs(k)%rarray( arraysize ) ) !!$ allocate( x_Dim%attrs(k)%darray( arraysize ) ) !!$ !!$ x_Dim%attrs(k)%array = .true. !!$ !!$ else !!$ call axis_attrs_init(x_Dim%attrs(k)) !!$ endif if (StrHead(attrtype, 'char')) then call HistoryAxisAddAttr(x_Dim % axisinfo, attrname, cvalue) elseif (StrHead(attrtype, 'logical')) then call HistoryAxisAddAttr(x_Dim % axisinfo, attrname, lvalue) elseif (StrHead(attrtype, 'int')) then if (arraysize > 0) then call HistoryAxisAddAttr(x_Dim % axisinfo, attrname, iarray) else call HistoryAxisAddAttr(x_Dim % axisinfo, attrname, ivalue) end if elseif (StrHead(attrtype, 'real')) then if (arraysize > 0) then call HistoryAxisAddAttr(x_Dim % axisinfo, attrname, rarray) else call HistoryAxisAddAttr(x_Dim % axisinfo, attrname, rvalue) end if elseif (StrHead(attrtype, 'double')) then if (arraysize > 0) then call HistoryAxisAddAttr(x_Dim % axisinfo, attrname, darray) else call HistoryAxisAddAttr(x_Dim % axisinfo, attrname, dvalue) end if end if !!$ x_Dim%attrs(k)%attrname = attrname !!$ x_Dim%attrs(k)%attrtype = attrtype !!$ x_Dim%attrs(k)%cvalue = cvalue !!$ x_Dim%attrs(k)%ivalue = ivalue !!$ x_Dim%attrs(k)%rvalue = rvalue !!$ x_Dim%attrs(k)%dvalue = dvalue !!$ x_Dim%attrs(k)%lvalue = lvalue !!$ x_Dim%attrs(k)%iarray(1:max(1,arraysize)) = iarray(1:max(1,arraysize)) !!$ x_Dim%attrs(k)%rarray(1:max(1,arraysize)) = rarray(1:max(1,arraysize)) !!$ x_Dim%attrs(k)%darray(1:max(1,arraysize)) = darray(1:max(1,arraysize)) !!$ call DbgMessage('x_Dim-attrs(%d) [attrname=<%c> ' // & !!$ & 'attrtype=<%c> array=<%b> cvalue=<%c> ' // & !!$ & 'ivalue=<%d> rvalue=<%r> dvalue=<%f> ' // & !!$ & 'iarray(1:%d)=<%d, ...> ' // & !!$ & 'rarray(1:%d)=<%r, ...> darray(1:%d)=<%f, ...>' , & !!$ & c1=trim( x_Dim%attrs(k)%attrname ) , & !!$ & c2=trim( x_Dim%attrs(k)%attrtype ) , & !!$ & c3=trim( x_Dim%attrs(k)%cvalue ) , & !!$ & i=(/ k, x_Dim%attrs(k)%ivalue , & !!$ & size(x_Dim%attrs(k)%iarray) , & !!$ & x_Dim%attrs(k)%iarray , & !!$ & size(x_Dim%attrs(k)%rarray) , & !!$ & size(x_Dim%attrs(k)%darray) & !!$ & /) , & !!$ & r=(/x_Dim%attrs(k)%rvalue, x_Dim%attrs(k)%rarray/) , & !!$ & d=(/x_Dim%attrs(k)%dvalue, x_Dim%attrs(k)%darray/) , & !!$ & l=(/x_Dim%attrs(k)%lvalue/) ) if (.not. next) exit axis_x_attr_nml_input next = .false. ! 次回のための初期化 enddo axis_x_attr_nml_input end if call nmlfile_close !---------------------------------------------------------------- ! grid_3d_mod との整合性をチェック !---------------------------------------------------------------- call grid_3d_init if (length /= im) then call MessageNotify('E', subname, & & message='axis length is inconsistent with im in grid_3d_mod') endif !---------------------------------------------------------------- ! 例外処理 !---------------------------------------------------------------- if (length < 1) then call MessageNotify('E', subname, message='Invalid grid number.') endif call EndSub( subname ) end subroutine axis_x_init !=begin !=== Return Weight of axis X ! !重みデータとその付加情報を返す。 !((< axis_x_init >)) の NAMELIST axis_x_nml の decision 変数で !(({ 'spectral' })) 以外が与えられた場合は値を代入しないで返す。 ! !また、X 軸の次元変数に重みデータに関する付加情報を加える。 ! subroutine axis_x_weight(Dim_Weight) !==== Dependency use constants_mod, only: constants_init, pi use axis_type_mod, only: axis_type_copy ! use axis_type_mod, only: axis_type_copy, axis_attrs_copy, axis_attrs_init use spml_mod, only: wa_module_x_Lon_Weight => x_Lon_Weight use grid_3d_mod,only: im ! use gt4_history,only: GT_HISTORY_ATTR use dc_trace, only: DbgMessage, BeginSub, EndSub use gt4_history,only: HistoryAxisCopy, HistoryAxisInquire !=end implicit none !=begin !==== Output ! type(AXISINFO), intent(out) :: Dim_Weight ! 次元情報を包括する変数 !=end ! type(GT_HISTORY_ATTR), allocatable :: attrs_tmp(:) ! 次元情報を一時格納 character(STRING) :: name, longname character(STRING), parameter:: subname = "axis_x_weight" continue !---------------------------------------------------------------- ! Check Initialization !---------------------------------------------------------------- call BeginSub( subname ) if (.not. axis_x_initialized) then call EndSub( subname, 'Call axis_x_init before call %c', & & c1=trim(subname) ) return endif !---------------------------------------------------------------- ! decision が 'spectral' でない場合は停止して返す。 !---------------------------------------------------------------- if (.not. axis_x_data_from_spectral) then call EndSub( subname, & & 'In axis_x_nml, configurated to not generate from spmodel') return endif !---------------------------------------------------------------- ! Initialize Dependent module !---------------------------------------------------------------- call constants_init() !---------------------------------------------------------------- ! x_Dim の情報から x_Dim_Weight の情報を生成。 !---------------------------------------------------------------- call HistoryAxisInquire(x_Dim % axisinfo, name=name, longname=longname) call HistoryAxisCopy(& & axis_dest=x_Dim_Weight % axisinfo, & & axis_src=x_Dim % axisinfo, & & name=trim(name) // '_weight', longname=trim(longname) // ' weight') !!$ x_Dim_Weight%axisinfo%name = trim(x_Dim%axisinfo%name) // '_weight' !!$ x_Dim_Weight%axisinfo%length = x_Dim%axisinfo%length !!$ x_Dim_Weight%axisinfo%longname = trim(x_Dim%axisinfo%longname) // ' weight' !!$ x_Dim_Weight%axisinfo%units = trim(x_Dim%axisinfo%units) !!$ x_Dim_Weight%axisinfo%xtype = trim(x_Dim%axisinfo%xtype) !!$ if ( associated(x_Dim_Weight%attrs) ) then !!$ deallocate( x_Dim_Weight%attrs ) !!$ endif !!$ !!$ ! x_Dim にも情報付加 !!$ if ( associated(x_Dim%attrs) ) then !!$ allocate( attrs_tmp(size(x_Dim%attrs)) ) !!$ call axis_attrs_copy(from=x_Dim%attrs(:), to=attrs_tmp(:)) !!$ deallocate( x_Dim%attrs ) !!$ allocate( x_Dim%attrs(size(attrs_tmp)+1) ) !!$ call axis_attrs_copy(from=attrs_tmp(:), to=x_Dim%attrs(1:size(attrs_tmp)) ) !!$ !!$ call axis_attrs_init( x_Dim%attrs(size(attrs_tmp)+1) ) !!$ x_Dim%attrs(size(attrs_tmp)+1)%attrname = 'gt_calc_weight' !!$ x_Dim%attrs(size(attrs_tmp)+1)%attrtype = 'char' !!$ x_Dim%attrs(size(attrs_tmp)+1)%cvalue = x_Dim_Weight%axisinfo%name !!$ x_Dim%attrs(size(attrs_tmp)+1)%array = .false. !!$ else !!$ allocate( x_Dim%attrs(1) ) !!$ call axis_attrs_init( x_Dim%attrs(1) ) !!$ x_Dim%attrs(1)%attrname = 'gt_calc_weight' !!$ x_Dim%attrs(1)%attrtype = 'char' !!$ x_Dim%attrs(1)%cvalue = x_Dim_Weight%axisinfo%name !!$ x_Dim%attrs(1)%array = .false. !!$ endif !---------------------------------------------------------------- ! wa_module からスペクトルデータを受け取る。 !---------------------------------------------------------------- allocate( x_Dim_Weight%a_Dim(im) ) ! ラジアンを度数に変換して代入 x_Dim_Weight%a_Dim(:) = wa_module_x_Lon_Weight(:) * 180. / pi x_Dim_Weight%stored = .true. !---------------------------------------------------------------- ! intent(out) の引数にデータを渡す。 !---------------------------------------------------------------- call axis_type_copy(x_Dim_Weight, Dim_Weight) call EndSub( subname ) end subroutine axis_x_weight !=begin !=== Return Data of axis X for Spectral Method. ! !スペクトル法を用いる場合を想定した X 軸のデータを返す。 ! !((< axis_x_init >)) の NAMELIST axis_x_nml の units に !(({ 'radian' })) または (({ 'rad.' })) を与える場合には !単位がラジアンでデータが返される。それ以外では度数でデータが返る。 ! !((< axis_x_init >)) の NAMELIST axis_x_nml の decision 変数で !(({ 'spectral' })) 以外が与えられた場合は値を代入しないで返す。 ! subroutine axis_x_spectral(Dim) !==== Dependency use axis_type_mod, only: axis_type_copy use constants_mod, only: constants_init, pi use spml_mod, only: wa_module_x_Lon => x_Lon use dc_string, only: toChar, StrHead, LChar use dc_trace, only: DbgMessage, BeginSub, EndSub use gt4_history, only: HistoryAxisInquire !=end implicit none !=begin !==== In/Out ! type(AXISINFO), intent(inout) :: Dim ! 次元情報を包括する変数 !=end character(STRING) :: units character(STRING), parameter:: subname = "axis_x_spectral" continue !---------------------------------------------------------------- ! Check Initialization !---------------------------------------------------------------- call BeginSub( subname ) if (.not. axis_x_initialized) then call EndSub( subname, 'Call axis_x_init before call %c', & & c1=trim(subname) ) return endif !---------------------------------------------------------------- ! decision が 'spectral' でない場合は停止して返す。 !---------------------------------------------------------------- if (.not. axis_x_data_from_spectral) then call EndSub( subname, & & 'In axis_x_nml, configurated to not generate from spmodel') return endif !---------------------------------------------------------------- ! Initialize Dependent module !---------------------------------------------------------------- call constants_init() !---------------------------------------------------------------- ! wa_module からスペクトルデータを受け取る。 !---------------------------------------------------------------- call HistoryAxisInquire(x_Dim % axisinfo, units=units) if ( StrHead( 'radians', trim(LChar(units)) ) .or.& & StrHead( 'rad.', trim(LChar(units)) ) ) then ! radian でそのまま代入 x_Dim%a_Dim(:) = wa_module_x_Lon(:) else ! radian を degree に変換して代入 x_Dim%a_Dim(:) = wa_module_x_Lon(:) * 180. / pi end if x_Dim%stored = .true. call DbgMessage('x_Lon=<%c>', c1=trim(toChar(x_Dim%a_Dim)) ) !---------------------------------------------------------------- ! intent(out) の引数にデータを渡す。 !---------------------------------------------------------------- call axis_type_copy(x_Dim, Dim) call EndSub( subname ) end subroutine axis_x_spectral !=begin !=== Return Data of axis X from NAMELIST ! !NAMELIST から代入されたデータを X 軸データとして返す。 ! !((< axis_x_init >)) の NAMELIST axis_x_nml の decision 変数で !(({ 'manual' })) 以外が与えられた場合は値を代入しないで返す。 ! subroutine axis_x_manual(Dim) !==== Dependency use axis_type_mod, only: axis_type_copy use spml_mod, only: wa_module_x_Lon => x_Lon use dc_trace, only: DbgMessage, BeginSub, EndSub !=end implicit none !=begin !==== In/Out ! type(AXISINFO), intent(inout) :: Dim ! 次元情報を包括する変数 !=end character(STRING), parameter:: subname = "axis_x_manual" continue !---------------------------------------------------------------- ! Check Initialization !---------------------------------------------------------------- call BeginSub( subname ) if (.not. axis_x_initialized) then call EndSub( subname, 'Call axis_x_init before call %c', & & c1=trim(subname) ) return endif !---------------------------------------------------------------- ! decision が 'manual' でない場合は停止して返す。 !---------------------------------------------------------------- if (.not. axis_x_data_from_manual) then call EndSub( subname, & & 'In axis_x_nml, configurated to not generate from manual') return endif !---------------------------------------------------------------- ! intent(out) の引数にデータを渡す。 !---------------------------------------------------------------- call axis_type_copy(x_Dim, Dim) call EndSub( subname ) end subroutine axis_x_manual !=begin !=== Return Data of axis X from netCDF file ! !netCDF データから取得したデータを X 軸のデータとして返す。 ! !現在、取得先のデータの単位に関わらず、そのままデータが !入力される。(()) 参照。 ! !((< axis_x_init >)) の NAMELIST axis_x_nml の decision 変数で !gtool4 変数以外が与えられた場合は値を代入しないで返す。 ! subroutine axis_x_netcdf(Dim) !==== 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 ! 次元情報を包括する変数 !=end character(STRING) :: file character(STRING) :: varname character(STRING), parameter:: subname = "axis_x_netcdf" continue !---------------------------------------------------------------- ! Check Initialization !---------------------------------------------------------------- call BeginSub( subname ) if (.not. axis_x_initialized) then call EndSub( subname, 'Call axis_x_init before call %c', & & c1=trim(subname) ) return endif !---------------------------------------------------------------- ! decision が '***@lon' などでない場合は停止して返す。 !---------------------------------------------------------------- if (axis_x_data_from_netcdf == '') then call EndSub( subname, & & 'In axis_x_nml, configurated to not generate from netcdf data') return endif !---------------------------------------------------------------- ! HistoryGet によるデータ入力 !---------------------------------------------------------------- call UrlSplit(axis_x_data_from_netcdf, file=file, var=varname) call HistoryGet(file, varname, x_Dim%a_Dim) x_Dim%stored = .true. ! 入力するデータの units が 'degree' か 'radian' か調べ、 ! 一時、'radian' に変換し、その後、x_Dim%axisinfo%units が ! 'degree' か 'radian' かに応じて 180. / pi を掛けるかどうか ! 判定すべき !---------------------------------------------------------------- ! intent(out) の引数にデータを渡す。 !---------------------------------------------------------------- call axis_type_copy(x_Dim, Dim) call EndSub( subname ) end subroutine axis_x_netcdf !=begin !=== Terminate module ! !(()) で設定された値を破棄し、デフォルトに戻す。 ! subroutine axis_x_end() !==== Dependency use dc_trace, only: DbgMessage, BeginSub, EndSub use gt4_history, only: HistoryAxisClear !=end implicit none character(STRING), parameter:: subname = "axis_x_end" continue !---------------------------------------------------------------- ! Check Initialization !---------------------------------------------------------------- call BeginSub(subname) if ( .not. axis_x_initialized) then call EndSub( subname, 'axis_x_init was not called', & & c1=trim(subname) ) return else axis_x_initialized = .false. axis_x_data_from_spectral = .false. axis_x_data_from_manual = .false. axis_x_data_from_netcdf = '' !-------------------------------------------------------------- ! Initialize x_Dim !-------------------------------------------------------------- x_Dim%stored = .false. call HistoryAxisClear(x_Dim % axisinfo) !!$ x_Dim%axisinfo%name = '' !!$ x_Dim%axisinfo%length = 0 !!$ x_Dim%axisinfo%longname = '' !!$ x_Dim%axisinfo%units = '' !!$ x_Dim%axisinfo%xtype = '' !!$ if ( associated(x_Dim%a_Dim) ) then !!$ deallocate( x_Dim%a_Dim ) !!$ endif !!$ !!$ if ( associated(x_Dim%attrs) ) then !!$ deallocate( x_Dim%attrs ) !!$ end if !-------------------------------------------------------------- ! Initialize x_Dim_Weight !-------------------------------------------------------------- call HistoryAxisClear(x_Dim_Weight % axisinfo) !!$ x_Dim_Weight%axisinfo%name = '' !!$ x_Dim_Weight%axisinfo%length = 0 !!$ x_Dim_Weight%axisinfo%longname = '' !!$ x_Dim_Weight%axisinfo%units = '' !!$ x_Dim_Weight%axisinfo%xtype = '' !!$ !!$ if ( associated(x_Dim%a_Dim) ) then !!$ deallocate( x_Dim%a_Dim ) !!$ end if !!$ !!$ if ( associated(x_Dim%attrs) ) then !!$ deallocate( x_Dim%attrs ) !!$ end if endif call EndSub( subname ) end subroutine axis_x_end end module axis_x_mod