!= 速度から渦度発散を出力するプログラム ! != A program that generates vorticity and divergence from velocity ! ! Authors:: Yasuhiro MORIKAWA ! Version:: $Id: convert_uv2vordiv.f90,v 1.1 2008-04-21 10:28:49 morikawa Exp $ ! Tag Name:: $Name: dcpam4-20080609-1 $ ! Copyright:: Copyright (C) GFD Dennou Club, 2007. All rights reserved. ! License:: See COPYRIGHT[link:../../../COPYRIGHT] ! ! Note that Japanese and English are described in parallel. ! ! 入力ファイルの "U" と "V" という変数の情報を参照し, ! "Vor" と "Div" に変換して出力ファイルへと出力します. ! 入力ファイルには軸データ "lon", "lon_weight", "lat", "lat_weight", ! "sig", "sigm" が格納されている ! 事を想定します. ! ! Variables "U" and "V" in an input file are converted into ! "Vor" and "Div" . The velocities are output to an output file. ! It is expected that axes data "lon", "lon_weight", "lat", "lat_weight", ! "sig", "sigm" are stored in input file. ! program convert_uv2vordiv use constants, only: CONST, Create, Get use dc_types, only: DP, STRING use dc_message, only: MessageNotify use dyn_spectral_as83, only: DYNSPAS83, DynSpAsCreate, DynSpAsClose, & & DynSpAsPutLine, DynSpAsInitialized, UV2VorDiv use dc_string, only: StoA use dc_args, only: ARGS, Open, HelpMsg, Option, Debug, Help, Strict, Close use gt4_history, only: HistoryGet, HistoryGetPointer use gt4_history, only: GT_HISTORY, & & HistoryCreate, HistoryAddVariable, HistoryPut, HistoryClose, & & HistoryAddAttr implicit none !--------------------------------------------------------- ! 実験の表題, モデルの名称, 所属機関名 ! Title of a experiment, name of model, sub-organ !--------------------------------------------------------- character(*), parameter:: title = & & 'convert_uv2vordiv $Name: dcpam4-20080609-1 $ :: ' // & & 'Test program of "dyn_spectral_as83" module' character(*), parameter:: source = & & 'dcpam4 ' // & & '(See http://www.gfd-dennou.org/library/dcpam)' character(*), parameter:: institution = & & 'GFD Dennou Club (See http://www.gfd-dennou.org)' !------------------------------------------------------------------- ! 格子点数・最大全波数 ! Grid points and maximum truncated wavenumber !------------------------------------------------------------------- integer:: nmax ! 最大全波数. ! Maximum truncated wavenumber integer:: imax ! 経度格子点数. ! Number of grid points in longitude integer:: jmax ! 緯度格子点数. ! Number of grid points in latitude integer:: kmax ! 鉛直層数. ! Number of vertical level !----------------------------------------------------------------- ! 物理定数等の設定 ! Configure physical constants etc. !----------------------------------------------------------------- type(CONST):: const_earth real(DP):: PI ! $ \pi $ . 円周率. Circular constant real(DP):: RPlanet ! $ a $ . 惑星半径. Radius of planet real(DP):: Omega ! $ \Omega $ . 回転角速度. Angular velocity !!$ real(DP):: Grav ! $ g $ . 重力加速度. Gravitational acceleration real(DP):: Cp ! $ C_p $ . 大気定圧比熱. Specific heat of air at constant pressure real(DP):: RAir ! $ R $ . 大気気体定数. Gas constant of air real(DP):: EpsV ! $ \epsilon_v $ . 水蒸気分子量比. Molecular weight ratio of water vapor integer:: VisOrder ! 超粘性の次数. Order of hyper-viscosity real(DP):: EFoldTime ! 最大波数に対する e-folding time. E-folding time for maximum wavenumber real(DP):: DelTime !--------------------------------------------------------- ! 軸データ ! Axes data !--------------------------------------------------------- real(DP), pointer:: x_Lon (:) ! 経度. Longitude real(DP), pointer:: x_Lon_Weight(:) ! 経度座標重み. ! Weight of longitude real(DP), pointer:: y_Lat (:) ! 緯度. Latitude real(DP), pointer:: y_Lat_Weight(:) ! 緯度座標重み. ! Weight of latitude real(DP), pointer:: z_Sigma (:) ! $ \sigma $ レベル (整数). ! Full $ \sigma $ level real(DP), pointer:: r_Sigma (:) ! $ \sigma $ レベル (半整数). ! Half $ \sigma $ level !--------------------------------------------------------- ! 物理量 ! Physical values !--------------------------------------------------------- real(DP), allocatable:: xyz_Vor (:,:,:) ! $ \zeta (t-\Delta t) $ . 渦度. Vorticity real(DP), allocatable:: xyz_Div (:,:,:) ! $ D (t-\Delta t) $ . 発散. Divergence real(DP), allocatable:: xyz_U (:,:,:) ! $ U $ . 東西風速. Zonal wind real(DP), allocatable:: xyz_V (:,:,:) ! $ V $ . 南北風速. Meridional wind !--------------------------------------------------------- ! 入出力データ ! Input/Output data !--------------------------------------------------------- type(GT_HISTORY):: gthist character(STRING):: input_data character(STRING):: output_data !--------------------------------------------------------- ! 作業変数 ! Work variables !--------------------------------------------------------- type(ARGS):: arg ! コマンドライン引数. ! Command line arguments !!$ logical:: OPT_namelist ! -N, --namelist オプションの有無. !!$ ! Existence of '-N', '--namelist' option !!$ character(STRING):: VAL_namelist !!$ ! -N, --namelist オプションの値. !!$ ! Value of '-N', '--namelist' option logical:: OPT_input ! -I, --input オプションの有無. ! Existence of '-I', '--input' option character(STRING):: VAL_input ! -I, --input オプションの値. ! Value of '-I', '--input' option logical:: OPT_output ! -O, --output オプションの有無. ! Existence of '-O', '--output' option character(STRING):: VAL_output ! -O, --output オプションの値. ! Value of '-O', '--output' option type(DYNSPAS83):: dyn_sp_as00 logical:: err character(*), parameter:: subname = 'uv2vordiv' continue !--------------------------------------------------------- ! コマンドライン引数の処理 ! Command line arguments handling !--------------------------------------------------------- call Open( arg ) call HelpMsg( arg, 'Title', title ) call HelpMsg( arg, 'Usage', & & './convert_uv2vordiv [Options] --input= --output=' ) call HelpMsg( arg, 'Source', source ) call HelpMsg( arg, 'Institution', institution ) !!$ call Option( arg, StoA('-N', '--namelist'), & !!$ & OPT_namelist, VAL_namelist, help = "NAMELIST filename" ) call Option( arg, StoA('-I', '--input'), & & OPT_input, VAL_input, help = "input data file" ) call Option( arg, StoA('-O', '--output'), & & OPT_output, VAL_output, help = "output data file" ) call Debug( arg ) ; call Help( arg ) ; call Strict( arg, severe = .true. ) call Close( arg ) !--------------------------------------------------------- ! 入出力データ ! Input/Output data !--------------------------------------------------------- input_data = VAL_input output_data = VAL_output if ( trim(input_data) == '' ) then call MessageNotify('E', subname, 'Specify input data like "--input=VorDiv.nc"') end if if ( trim(output_data) == '' ) then call MessageNotify('E', subname, 'Specify output data like "--output=UV.nc"') end if !--------------------------------------------------------- ! 物理定数の準備 ! Prepare physical constants !--------------------------------------------------------- call Create( const_earth ) ! (inout) DelTime = 480.0_DP call Get( constant = const_earth, & ! (inout) & PI = PI, RPlanet = RPlanet, Omega = Omega, & ! (out) & Cp = Cp, RAir = RAir, EpsV = EpsV, & ! (out) & VisOrder = VisOrder, EFoldTime = EFoldTime ) ! (out) !--------------------------------------------------------- ! 軸データ取得 ! Get axes data !--------------------------------------------------------- call HistoryGetPointer ( & & file = input_data, varname = 'lon', & ! (in) & array = x_Lon ) ! (out) call HistoryGetPointer ( & & file = input_data, varname = 'lon_weight', & ! (in) & array = x_Lon_Weight ) ! (out) call HistoryGetPointer ( & & file = input_data, varname = 'lat', & ! (in) & array = y_Lat ) ! (out) call HistoryGetPointer ( & & file = input_data, varname = 'lat_weight', & ! (in) & array = y_Lat_Weight ) ! (out) call HistoryGetPointer ( & & file = input_data, varname = 'sig', & ! (in) & array = z_Sigma ) ! (out) call HistoryGetPointer ( & & file = input_data, varname = 'sigm', & ! (in) & array = r_Sigma ) ! (out) !------------------------------------------------------------------- ! 格子点数・最大全波数の設定 ! Configure grid points and maximum truncated wavenumber !------------------------------------------------------------------- imax = size ( x_Lon ) jmax = size ( y_Lat ) kmax = size ( z_Sigma ) !------------------------------------- ! いい加減な最大全波数自動設定 ! Irresponsible auto-setting of maximum truncated wavenumber if ( imax < 32 .and. jmax < 16 ) then call MessageNotify('E', subname, '"imax=%d" and "jmax=%d" is too small', & & i = (/imax, jmax/) ) elseif ( imax < 64 .and. jmax < 32 ) then nmax = 10 elseif ( imax < 128 .and. jmax < 64 ) then nmax = 21 elseif ( imax < 192 .and. jmax < 96 ) then nmax = 42 else nmax = 63 end if call MessageNotify('M', subname, & & 'Resolution is nmax=%d imax=%d, jmax=%d, kmax=%d', & & i = (/nmax, imax, jmax, kmax/) ) !--------------------------------------------------------- ! 初期設定 ! Initialization !--------------------------------------------------------- call DynSpAsCreate( dyn_sp_as = dyn_sp_as00, & ! (inout) & nmax = nmax, & ! (in) & imax = imax, jmax = jmax, kmax = kmax, & ! (in) & PI = PI, RPlanet = RPlanet, & ! (in) & Omega = Omega, Cp = Cp, RAir = RAir, & ! (in) & EpsV = EpsV, & ! (in) & VisOrder = VisOrder, EFoldTime = EFoldTime, & ! (in) & DelTIme = DelTime ) ! (in) call DynSpAsPutLine( dyn_sp_as00 ) ! (in) !--------------------------------------------------------- ! 渦度発散データ読み込み ! Load vorticity and divergence data !--------------------------------------------------------- allocate( xyz_U(0:imax-1, 0:jmax-1, 0:kmax-1) ) allocate( xyz_V(0:imax-1, 0:jmax-1, 0:kmax-1) ) call HistoryGet ( & & file = input_data, varname = 'U', & ! (in) & array = xyz_U ) ! (out) call HistoryGet ( & & file = input_data, varname = 'V', & ! (in) & array = xyz_V ) ! (out) !--------------------------------------------------------- ! 速度データへ変換 ! Convert into velocities !--------------------------------------------------------- allocate( xyz_Vor(0:imax-1, 0:jmax-1, 0:kmax-1) ) allocate( xyz_Div(0:imax-1, 0:jmax-1, 0:kmax-1) ) call UV2VorDiv( dyn_sp_as = dyn_sp_as00, & ! (inout) & xyz_U = xyz_U, xyz_V = xyz_V, & ! (in) & xyz_Vor = xyz_Vor, xyz_Div = xyz_Div ) ! (out) !--------------------------------------------------------- ! データ出力 ! Output data !--------------------------------------------------------- call HistoryCreate( & & history = gthist, & ! (inout) & file = output_data, & ! (in) & title = 'Vor and Div data converted from "' // & & trim(input_data) // '"', & ! (in) & source = source, institution = institution, & ! (in) & dims = StoA('lon', 'lat', 'sig', 'sigm'), & ! (in) & dimsizes = (/imax, jmax, kmax, kmax + 1/), & ! (in) & longnames = & & StoA('longitude', 'latitude', & & 'sigma at layer midpoints', & & 'sigma at layer end-points (half level))'), & ! (in) & units = StoA('degree_east', 'degree_north', & & '1', '1') ) ! (in) call HistoryPut( & & history = gthist, & ! (inout) & varname = 'lon', array = x_Lon / PI * 180.0_DP ) ! (in) call HistoryPut( & & history = gthist, & ! (inout) & varname = 'lat', array = y_Lat / PI * 180.0_DP ) ! (in) call HistoryPut( & & history = gthist, & ! (inout) & varname = 'sig', array = z_Sigma ) ! (in) call HistoryPut( & & history = gthist, & ! (inout) & varname = 'sigm', array = r_Sigma ) ! (in) call HistoryAddAttr( & & history = gthist, & ! (inout) & varname = 'lon', attrname = 'standard_name', & ! (in) & value = 'longitude' ) ! (in) call HistoryAddAttr( & & history = gthist, & ! (inout) & varname = 'lat', attrname = 'standard_name', & ! (in) & value = 'latitude' ) ! (in) call HistoryAddAttr( & & history = gthist, & ! (inout) & varname = 'sig', attrname = 'standard_name', & ! (in) & value = 'atmosphere_sigma_coordinate' ) ! (in) call HistoryAddAttr( & & history = gthist, & ! (inout) & varname = 'sigm', attrname = 'standard_name', & ! (in) & value = 'atmosphere_sigma_coordinate' ) ! (in) call HistoryAddAttr( & & history = gthist, & ! (inout) & varname = 'sig', attrname = 'positive', & ! (in) & value = 'down' ) ! (in) call HistoryAddAttr( & & history = gthist, & ! (inout) & varname = 'sigm', attrname = 'positive', & ! (in) & value = 'down' ) ! (in) call HistoryAddVariable( & & history = gthist, & ! (inout) & varname = 'lon_weight', & ! (in) & dims = StoA('lon'), & ! (in) & longname = 'weight for integration in longitude', & ! (in) & units = 'radian', xtype = 'double' ) ! (in) call HistoryAddAttr( & & history = gthist, & ! (inout) & varname = 'lon', attrname = 'gt_calc_weight', & ! (in) & value = 'lon_weight' ) ! (in) call HistoryPut( & & history = gthist, & ! (inout) & varname = 'lon_weight', array = x_Lon_Weight ) ! (in) call HistoryAddVariable( & & history = gthist, & ! (inout) & varname = 'lat_weight', & ! (in) & dims = StoA('lat'), & ! (in) & longname = 'weight for integration in latitude', & ! (in) & units = 'radian', xtype = 'double' ) ! (in) call HistoryAddAttr( & & history = gthist, & ! (inout) & varname = 'lat', attrname = 'gt_calc_weight', & ! (in) & value = 'lat_weight' ) ! (in) call HistoryPut( & & history = gthist, & ! (inout) & varname = 'lat_weight', array = y_Lat_Weight ) ! (in) call HistoryAddVariable( & & history = gthist, & ! (inout) & varname = 'Vor', & ! (in) & dims = StoA('lon', 'lat', 'sig'), & ! (in) & longname = 'vorticity', & ! (in) & units = 's-1', xtype = 'double' ) ! (in) call HistoryAddAttr( & & history = gthist, & ! (inout) & varname = 'Vor', attrname = 'standard_name', & ! (in) & value = 'atmosphere_relative_vorticity' ) ! (in) call HistoryAddVariable( & & history = gthist, & ! (inout) & varname = 'Div', & ! (in) & dims = StoA('lon', 'lat', 'sig'), & ! (in) & longname = 'divergence', & ! (in) & units = 's-1', xtype = 'double' ) ! (in) call HistoryAddAttr( & & history = gthist, & ! (inout) & varname = 'Div', attrname = 'standard_name', & ! (in) & value = 'divergence_of_wind' ) ! (in) call HistoryPut( & & history = gthist, & ! (inout) & varname = 'Vor', array = xyz_Vor ) ! (in) call HistoryPut( & & history = gthist, & ! (inout) & varname = 'Div', array = xyz_Div ) ! (in) call HistoryClose( history = gthist ) ! (inout) !--------------------------------------------------------- ! 終了処理 ! Termination !--------------------------------------------------------- call DynSpAsClose( dyn_sp_as00 ) ! (inout) call DynSpAsPutLine( dyn_sp_as00 ) end program convert_uv2vordiv