#!/usr/bin/env ruby
# -*- f90 -*-
# vi: set sw=4 ts=8:
require("intrinsic_types")
require("optparse")
#
# "historyget.f90" Generator with Ruby.
#
opt = OptionParser.new
opt.on('--histget_dim=VAL') {|v| $histget_dim = v.to_i}
opt.parse!(ARGV)
$histget_dim = 7 unless $histget_dim
print <<"__EndOfFortran90Code__"
!--
#{rb2f90_header_comment}!
!++
!
!= gtool4 netCDF データの入力
!= Input gtool4 netCDF data
!
! Authors:: Yasuhiro MORIKAWA
! Version:: $Id: historyget.rb2f90,v 1.12 2007-09-13 15:38:48 morikawa Exp $
! Tag Name:: $Name: gt4f90io-20080603 $
! Copyright:: Copyright (C) GFD Dennou Club, 2006. All rights reserved.
! License:: See COPYRIGHT[link:../../COPYRIGHT]
!
! 以下のサブルーチン, 関数は gt4_history から gt4_history#HistoryGet
! もしくは gt4_history#HistoryGetPointer として提供されます.
!
! Following subroutines and functions are provided as
! gt4_history#HistoryGet or gt4_history#HistoryGetPointer from
! gt4_history.
!
__EndOfFortran90Code__
types = ["Double", "Real", "Int"]
types.each{ |type|
for num in 0..$histget_dim
print <<"__EndOfFortran90Code__"
subroutine HistoryGet#{type}#{num}(file, varname, array, range, quiet, err)
!
! Note that Japanese and English are described in parallel.
!
#{ifelse(type, "Double", %Q{
#{ifelse(num, 0, %Q{
!
! netCDF データを入力します. gtool4 netCDF 規約に基づくデータを想定
! していますが, 大抵の netCDF データの入力は可能であると期待されます.
!
! デフォルトでは, ファイル内の最新データ, すなわちデータを
! 時刻最大で切り出したものを入力します. 別の時刻または
! 別の次元で切り出したデータを
! 入力したい場合には, 下記の *time* オプションもしくは
! *range* オプションを利用してください.
!
! *file* にファイル名を, *varname* に変数名を与えます.
! *array* にはファイルから入力されたデータが返ります.
! ポインタの配列へデータを入力を行う場合は,
! HistoryGetPointer を利用してください.
!
! ある時刻のデータを明示したい場合には, その時刻を
! *time* に与えます.
! 整数型, 単精度実数型, 倍精度実数型の数値を与えることが可能です.
!
! *range* には gtool4 のコンマ記法
! ({gtool4 netCDF 規約}[link:../xref.htm#label-6] の「5.4 コンマ記法」参照)
! を与えることで, 任意の次元で入力データを切り出すことが可能です.
! *range* に空文字を与えた場合は切り出しを行いません.
!
! *HistoryGet* は複数のサブルーチンの総称名です. *array* には
! 0 〜 #{$histget_dim} 次元の整数型, 単精度実数型,
! 倍精度実数型の配列を与えることが可能です.
!
! デフォルトでは, データの入力時にどのファイルのどの変数が
! どの次元で切り出されて入力されたのかを表示します.
! メッセージ出力が不要な場合は *quiet* に .true. を与えてください.
!
! Input netCDF data. NetCDF data with gtool4 conventions is assumed,
! but most netCDF data is expected to be input.
!
! Give filename to *file*, variable name to *varname*.
! Input data is returned to *array*.
! If you want to use pointer array, use "HistoryGetPointer".
!
! By default, recent data, that is to say, data clipped with
! maximum time is input. In order to input data clipped with
! other time or other dimension,
! use *time* option or *range* option as follows.
!
! In order to get data at certain time, specify the time to *time*.
! Type is integer or single precision real or double precision.
!
! Give gtool4 comma-graphy
! (See "5.4 gtool4 comma-graphy" in {gtool4 netCDF Convention}[link:../xref.htm#label-6])
! to *range*, then input data can be clipped by an arbitrary dimension.
! If null character is given to *range*, data is not clipped.
!
! *HistoryGet* is generic name of multiple subroutines.
! Integer, single precision real, and double precision real
! 0 -- #{$histget_dim} rank array can be given to *array*.
!
! By default, when data is input, filename and variable name and
! clipping information is printed.
! The message is suppressed when .true. is given to *quiet*
!
}, %Q{
!
! 使用方法に関しては HistoryGet を参照してください.
!
! See "HistoryGet", for basic usage.
!
})}
}, %Q{
!
! 使用方法に関しては HistoryGet を参照してください.
!
! See "HistoryGet", for basic usage.
!
})}
!
use gtdata_types, only: GT_VARIABLE
use gtdata_generic, only: Open, Inquire, Close, Get
use dc_string, only: toChar, Split, JoinChar, StoA
use dc_url, only: GT_ATMARK, GT_COMMA, GT_EQUAL, UrlSplit, UrlMerge
use dc_present, only: present_select, present_and_not_empty, present_and_true
use regex, only: match
use dc_types, only: DP, STRING
use dc_message, only: MessageNotify
use dc_trace, only: Beginsub, Endsub, DbgMessage
use dc_error, only: StoreError, DC_NOERR, GT_ENOTURL, GT_ERANKMISMATCH, &
& GT_EARGSIZEMISMATCH
character(*), intent(in):: file
#{ifelse(type, "Double", %Q{
#{ifelse(num, 0, %Q{
! netCDF ファイル名.
! NetCDF filename.
})}
})}
character(*), intent(in):: varname
#{ifelse(type, "Double", %Q{
#{ifelse(num, 0, %Q{
! 変数名.
! Variable name
})}
})}
character(*), intent(in), optional:: range
#{ifelse(type, "Double", %Q{
#{ifelse(num, 0, %Q{
! 切り出し用オプション.
! gtool4 変数のコンマ記法で記述
! {(例: time=100.0,x=10:20,y=^1:^5)}
!
! 詳しくは
! {gtool4 netCDF 規約}[link:../xref.htm#label-6]
! の「5.4 コンマ記法」を参照して
! ください.
!
! Option for clipping.
! Give gtool4 comma-graphy
! {(ex. time=100.0,x=10:20,y=^1:^5)}
!
! For details, see "5.4 gtool4
! comma-graphy" in
! {gtool4 netCDF Convention}[link:../xref.htm#label-6]
})}
})}
logical, intent(in), optional:: quiet
#{$type_intent_out[type]}, intent(out) :: array#{array_colon("#{num}")}
#{ifelse(type, "Double", %Q{
#{ifelse(num, 0, %Q{
! 取得するデータを格納する配列
!
! 型は整数型, 単精度実数型, 倍精度実数型
! かのいづれかです. 取得するデータの空
! 間次元のサイズと配列のサイズとが一致し
! ている必要があります. 入力するデータ
! の型と *array* の型が異なる場合は, 自
! 動的に *array* の型に変換されます.
!
! Array in which input data is store
!
! Type is integer or single precision
! real or double precision. Size of
! array must be identical to input data
! size. When type of input data is
! different from type of *array*, data
! is converted to type of *array*
! automatically.
!
})}
})}
logical, intent(out), optional:: err
#{ifelse(type, "Double", %Q{
#{ifelse(num, 0, %Q{
! 例外処理用フラグ.
! デフォルトでは, この手続き内でエラーが
! 生じた場合, プログラムは強制終了します.
! 引数 *err* が与えられる場合,
! プログラムは強制終了せず, 代わりに
! *err* に .true. が代入されます.
!
! Exception handling flag.
! By default, when error occur in
! this procedure, the program aborts.
! If this *err* argument is given,
! .true. is substituted to *err* and
! the program does not abort.
})}
})}
#{ifelse(num, 0, %Q{
#{$type_internal[type]}:: array_tmp(1)
}, %Q{
#{$type_internal[type]}, allocatable:: array_tmp(:)
integer:: array_allsize
integer:: array_shape(#{num}), data_shape(#{num})
integer:: allcount
})}
type(GT_VARIABLE) :: var
character(STRING) :: url, actual_url
integer:: rank
integer:: domain ! 変数の入出力領域の大きさ
! (= 変数が依存する各次元サイズの積)
integer:: stat
character(STRING):: cause_c
character(*), parameter :: subname = "HistoryGet#{type}#{num}"
interface
subroutine lookup_growable_url(file, varname, url, range, err)
character(*), intent(in) :: file
character(*), intent(in) :: varname
character(*), intent(out) :: url
character(*), intent(in), optional:: range
logical, intent(out), optional :: err
end subroutine lookup_growable_url
end interface
interface
subroutine actual_iorange_dump(url, actual_url, err)
character(*), intent(in) :: url ! 変数 URL
character(*), intent(out), optional :: actual_url
! 正確な入出力範囲指定
logical, intent(out), optional :: err ! エラーのフラグ
end subroutine actual_iorange_dump
end interface
continue
call BeginSub(subname, 'file=%c varname=%c range=%c', &
& c1 = trim(file), c2 = trim(varname), &
& c3 = trim(present_select('', 'no-range', range)))
cause_c = ''
stat = DC_NOERR
#{ifelse(num, 0, %Q{
}, %Q{
!-------------------------------------
! 配列形状の取得
! Get array shape
array_shape = shape( array )
array_allsize = size( array )
})}
!-------------------------------------
! 最新時刻の URL 取得
! Get URL of latest time
call lookup_growable_url(file, varname, url, range, err = err)
if ( present_and_true(err) ) then
stat = GT_ENOTURL
cause_c = url
goto 999
end if
!-------------------------------------
! ファイルオープン
! File open
call Open( var, url, err = err )
if ( present_and_true(err) ) then
stat = GT_ENOTURL
cause_c = url
goto 999
end if
!-------------------------------------
! 配列形状のチェック
! Check array shape
call Inquire( var = var, & ! (in)
& rank = rank ) ! (out)
if ( .not. #{num} == rank ) then
if ( .not. present_and_true(quiet) ) then
call MessageNotify('W', subname, &
& 'Rank of data (%c) is "%d", rank of argument is "%d"', &
& i = (/rank, #{num}/), c1 = trim(url) )
end if
stat = GT_ERANKMISMATCH
cause_c = 'array'
goto 999
end if
#{ifelse(num, 0, %Q{
}, %Q{
#{forloop("\\$dimnum\\$", 1, num, %Q{
call Inquire( var = var , dimord = $dimnum$, & ! (in)
& allcount = allcount ) ! (out)
data_shape($dimnum$) = allcount
})}
if ( .not. all( array_shape == data_shape ) ) then
if ( .not. present_and_true(quiet) ) then
call MessageNotify('W', subname, &
& 'Shape of data (%c) is (%c), shape of argument is (%c)', &
& c1 = trim( url ), &
& c2 = trim( toChar( data_shape ) ), &
& c3 = trim( toChar( array_shape ) ) )
end if
stat = GT_EARGSIZEMISMATCH
cause_c = 'array'
goto 999
end if
})}
!-------------------------------------
! データ取得
! Get data
call Inquire( var = var, & ! (in)
& size = domain ) ! (out)
#{ifelse(num, 0, %Q{
call Get( var = var, & ! (inout)
& nvalue = domain, & ! (in)
& value = array_tmp) ! (out)
array = array_tmp(1)
}, %Q{
if ( allocated( array_tmp ) ) deallocate( array_tmp )
allocate( array_tmp(array_allsize) )
call Get( var, array_tmp, domain )
array = reshape( array_tmp, array_shape )
deallocate( array_tmp )
})}
call Close( var )
!-------------------------------------
! データファイル名と切り出し範囲の印字
! Print data filename and clipping range
call actual_iorange_dump(url, actual_url, err = err)
if ( .not. present_and_true(quiet) ) then
call MessageNotify('M', subname, 'Input %c', c1=trim(actual_url))
end if
999 continue
call StoreError(stat, subname, err, cause_c)
call EndSub(subname)
end subroutine HistoryGet#{type}#{num}
__EndOfFortran90Code__
end
}
types = ["Double", "Real", "Int"]
types.each{ |type|
for num in 0..$histget_dim
print <<"__EndOfFortran90Code__"
subroutine HistoryGet#{type}#{num}Pointer(file, varname, array, range, quiet, err)
!
! Note that Japanese and English are described in parallel.
!
#{ifelse(type, "Double", %Q{
#{ifelse(num, 0, %Q{
! netCDF データを入力します. gtool4 netCDF 規約に基づくデータを想定
! していますが, 大抵の netCDF データの入力は可能であると期待されます.
!
! 基本的な使い方に関しては HistoryGet を参照してください.
! HistoryGet との違いは, *array* にポインタ配列を与えることです.
! *array* には必ず空状態の配列を与えてください.
! すなわち与える配列に対し, 初期値 =>null() を設定するか
! nullify を用いてください.
! 既に割り付けられている場合, もしくは不定状態の場合には
! エラーを返します.
!
! Input netCDF data. NetCDF data with gtool4 conventions is assumed,
! but most netCDF data is expected to be input.
!
! See "HistoryGet", for basic usage.
! Difference from "HistoryGet" is that *array* is pointer array.
! Give null array to *array*.
! More specifically, use '=>null()' as initial value or
! use 'nullify' to the array.
! If the array is allocated already, or undefined,
! error is occurred.
!
}, %Q{
! 使用方法に関しては HistoryGetPointer を参照してください.
!
! See "HistoryGetPointer", for basic usage.
!
})}
}, %Q{
! 使用方法に関しては HistoryGetPointer を参照してください.
!
! See "HistoryGetPointer", for basic usage.
!
})}
!
use gtdata_types, only: GT_VARIABLE
use gtdata_generic, only: Open, Inquire, Close, Get
use dc_string, only: toChar
use dc_present,only: present_select, present_and_true
use dc_types, only: DP, STRING
use dc_message, only: MessageNotify
use dc_trace, only: Beginsub, Endsub, DbgMessage
character(*), intent(in):: file
character(*), intent(in):: varname
character(*), intent(in), optional:: range
logical, intent(in), optional:: quiet
#{$type_intent_out[type]}, pointer :: array#{array_colon("#{num}")} ! (out)
#{ifelse(type, "Double", %Q{
#{ifelse(num, 0, %Q{
!
! 取得するデータを格納する
! ポインタ配列.
!
! 必ず空状態の配列を与えてください.
! すなわち与える配列に対し, 初期値 =>null()
! を設定するか nullify を用いて
! ください. 既に割り付けられている場合,
! もしくは不定状態の場合にはエラーを返し
! ます.
!
! 型は整数型, 単精度実数型, 倍精度実数型
! かのいづれかです. 入力するデータ
! の型と *array* の型が異なる場合は, 自
! 動的に *array* の型に変換されます.
!
! Array in which input data is store.
!
! Give null array to *array*. More
! specifically, use '=>null()' as
! initial value or use 'nullify' to the
! array. If the array is allocated
! already, or undefined, error is
! occurred.
!
! Type is integer or single precision
! real or double precision real. When type of
! input data is different from type of
! *array*, data is converted to type of
! *array* automatically.
!
})}
})}
logical, intent(out), optional:: err
#{ifelse(num, 0, %Q{
#{$type_internal[type]}, target :: array_tmp(1)
integer:: domain ! 変数の入出力領域の大きさ
! (= 変数が依存する各次元サイズの積)
})}
type(GT_VARIABLE) :: var
character(STRING) :: url, actual_url
character(*), parameter :: subname = "HistoryGet#{type}#{num}Pointer"
interface
subroutine lookup_growable_url(file, varname, url, range, err)
character(*), intent(in) :: file ! ファイル名
character(*), intent(in) :: varname ! 変数名
character(*), intent(out) :: url ! gtool変数化した文字列
character(*), intent(in), optional:: range ! 範囲限定や一点切り出し指定
logical, intent(out), optional :: err ! エラーのフラグ
end subroutine lookup_growable_url
end interface
interface
subroutine actual_iorange_dump(url, actual_url, err)
character(*), intent(in) :: url ! 変数 URL
character(*), intent(out), optional :: actual_url
! 正確な入出力範囲指定
logical, intent(out), optional :: err ! エラーのフラグ
end subroutine actual_iorange_dump
end interface
continue
call BeginSub(subname, 'file=%c varname=%c range=%c', &
& c1=trim(file), c2=trim(varname), &
& c3=trim(present_select('', 'no-range', range)))
! 必要な情報を gtool 変数化
call lookup_growable_url(file, varname, url, range, err)
#{ifelse(num, 0, %Q{
allocate(array)
})}
call DbgMessage('@ url =%c', c1=trim(url))
! いよいよデータ取得
call Open(var, url, err)
#{ifelse(num, 0, %Q{
call Inquire(var=var, size=domain)
call Get(var, array_tmp, domain, err)
array = array_tmp(1)
}, %Q{
call Get(var, array, err)
})}
call Close(var, err)
call actual_iorange_dump(url, actual_url, err)
if ( .not. present_and_true(quiet) ) then
call MessageNotify('M', subname, 'Input %c', c1=trim(actual_url))
end if
call EndSub(subname)
end subroutine HistoryGet#{type}#{num}Pointer
__EndOfFortran90Code__
end
}
timetypes = ["Double", "Real", "Int"]
types = ["Double", "Real", "Int"]
fixorptrs = ["", "Pointer"]
def timekind(timetype)
return "D" if timetype == "Double"
return "R" if timetype == "Real"
return "I" if timetype == "Int"
return ""
end
fixorptrs.each{ |fixorptr|
timetypes.each{ |timetype|
types.each{ |type|
for num in 0..$histget_dim
print <<"__EndOfFortran90Code__"
subroutine HistoryGet#{type}#{num}#{fixorptr}Time#{timekind(timetype)}(file, varname, array, time, quiet, err)
!
#{ifelse(fixorptr, "Pointer", %Q{
! Note that Japanese and English are described in parallel.
!
! 数値型引数 *time* で時刻指定可能な HistoryGet です.
! 使用方法に関しては HistoryGet を参照してください.
!
! This is "HistoryGet" with numerical argument *time* for
! specification of time.
! See "HistoryGet", for basic usage.
}, %Q{
! Note that Japanese and English are described in parallel.
!
! 数値型引数 *time* で時刻指定可能な HistoryGetPointer です.
! 使用方法に関しては HistoryGetPointer を参照してください.
!
! This is "HistoryGetPointer" with numerical argument *time* for
! specification of time.
! See "HistoryGetPointer", for basic usage.
})}
!
use dc_string, only: toChar, Split
use dc_types, only: DP, STRING
use dc_trace, only: Beginsub, Endsub, DbgMessage
use dc_url, only: Url_Chop_IOrange, GT_EQUAL
character(*), intent(in):: file, varname
#{$type_intent_in[timetype]}, intent(in):: time
logical, intent(in), optional:: quiet
#{ifelse(fixorptr, "Pointer", %Q{
#{$type_intent_out[type]}, pointer:: array#{array_colon("#{num}")} ! (out)
}, %Q{
#{$type_intent_out[type]}, intent(out):: array#{array_colon("#{num}")}
})}
logical, intent(out), optional:: err
character(STRING):: url, iorange, remainder, timevar_name, range
character(STRING), pointer:: carray (:)
character(*), parameter :: subname = "HistoryGet#{type}#{num}#{fixorptr}Time#{timekind(timetype)}"
interface
subroutine HistoryGet#{type}#{num}#{fixorptr}(file, varname, array, range, quiet, err)
use dc_types, only: DP
character(*), intent(in):: file
character(*), intent(in):: varname
character(*), intent(in), optional:: range
logical, intent(in), optional:: quiet
logical, intent(out), optional:: err
#{ifelse(fixorptr, "Pointer", %Q{
#{$type_intent_out[type]}, pointer:: array#{array_colon("#{num}")} ! (out)
}, %Q{
#{$type_intent_out[type]}, intent(out):: array#{array_colon("#{num}")}
})}
end subroutine HistoryGet#{type}#{num}#{fixorptr}
end interface
interface
subroutine lookup_growable_url(file, varname, url, range, err)
character(*), intent(in) :: file ! ファイル名
character(*), intent(in) :: varname ! 変数名
character(*), intent(out) :: url ! gtool変数化した文字列
character(*), intent(in), optional:: range ! 範囲限定や一点切り出し指定
logical, intent(out), optional :: err ! エラーのフラグ
end subroutine lookup_growable_url
end interface
continue
call BeginSub(subname, 'file=%c varname=%c time=%c', &
& c1=trim(file), c2=trim(varname), c3=toChar(time))
call lookup_growable_url(file = file, varname = varname, &
& url = url, err = err)
call Url_Chop_IOrange( &
& fullname = url, iorange = iorange, remainder = remainder )
call Split( str = iorange, carray = carray, sep = GT_EQUAL )
timevar_name = carray(1)
deallocate( carray )
range = trim(timevar_name) // GT_EQUAL // trim(toChar(time))
call HistoryGet#{type}#{num}#{fixorptr}( file = file, &
& varname = varname, array = array, &
& range = range, quiet = quiet, err = err )
call EndSub(subname)
end subroutine HistoryGet#{type}#{num}#{fixorptr}Time#{timekind(timetype)}
__EndOfFortran90Code__
end
}
}
}
undef timekind
print <<"__EndOfFortran90Code__"
subroutine lookup_growable_url(file, varname, url, range, err)
!
! file の変数 varname が依存する次元の内, 時間の次元
! (growable == .TRUE. のもの, つまり無制限次元) の変数名,
! およびその最後の値を取得し, gtool 変数化
! ("file@varname,time=10.5" みたいな) して返す.
!
! * もしも varname が次元変数である場合は「time=」を付けずに返す.
! * range を与えた場合, 以下のチェックを行った後, それを gtool4
! 変数の iorange 部分に付加する.
! * range に空文字が与えられた場合, range が与えられない場合と
! 同じ動作をする.
! * range 内に時間次元が設定されていない場合は, 自動的に
! 時間次元に関する iorange ("time=0.5") が指定される.
! * 数値のみの文字列 (例: "20", "10.354") が与えられる場合,
! エラーを生じる.
!
use gtdata_types, only: GT_VARIABLE
use gtdata_generic, only: Open, Close, Inquire
use dc_present,only: present_select, present_and_not_empty, present_and_true
use dc_string, only: toChar
use dc_error, only: StoreError, DC_NOERR, GT_ENOUNLIMITDIM, &
& NF_EINVAL, GT_ENOTVAR, GT_EBADGT4COMMAGRAPHY
use dc_url, only: GT_CIRCUMFLEX, GT_COMMA, GT_EQUAL
use dc_url, only: UrlSplit, UrlMerge, UrlSearchIORange
use regex, only: match
use dc_types, only: DP, STRING
use dc_trace, only: Beginsub, Endsub, DbgMessage
character(*), intent(in) :: file ! ファイル名
character(*), intent(in) :: varname ! 変数名
character(*), intent(out) :: url ! gtool変数化した文字列
character(*), intent(in), optional:: range ! 範囲限定や一点切り出し指定
logical, intent(out), optional :: err ! エラーのフラグ
!
type(GT_VARIABLE) :: var
type(GT_VARIABLE), allocatable :: dimvar(:)
character(STRING) :: time_url, time_name, time_iorange
character(STRING) :: iorange, cause_c
logical:: growable, nounlimited
integer:: allcount, timecount, nd, i, stat
integer:: regex_stat, regex_len
character(*), parameter :: subname = "lookup_growable_url"
continue
call BeginSub(subname, '', &
& c1=trim(file), c2=trim(varname), &
& c3=trim(present_select('', 'no-range', range)))
stat = DC_NOERR
cause_c = ""
url = ""
! 引数の正当性をチェック
if (.not. present_and_not_empty(file)) then
stat = NF_EINVAL
cause_c = '"file" is not specified'
goto 999
elseif (.not. present_and_not_empty(varname)) then
stat = NF_EINVAL
cause_c = '"varname" is not specified'
goto 999
end if
! 時刻次元の変数名, およびその最終時刻の
! 探査のために file@varname を open (まだデータを取得しない)
call Open(var, UrlMerge(file, varname), err = err)
if ( present_and_true(err) ) then
stat = GT_ENOTVAR
goto 999
end if
! 次元の数を取得
call Inquire(var=var, alldims=nd)
call DbgMessage('@ alldims = %d', i=(/nd/))
if (allocated(dimvar)) then
deallocate(dimvar)
end if
allocate(dimvar(nd))
!
! 変数が無制限変数を持たない場合には, それに関する iorange を
! 付けないで返すよう, フラグを立てる.
! 無制限次元があれば, .false. にする.
nounlimited = .true.
!
! 各次元毎に情報を取得し, growable == .TRUE. のもの (つまりは時間)
! の変数名 (time_name) を取得する.
call DbgMessage('[%c: growable-dim-search]', c1=trim(subname))
time_name = ''
do, i = 1, nd
call Open(var = dimvar(i), source_var = var, dimord = i, &
& count_compact = .TRUE., err = err)
! まずは変数入り gtool4 変数を time_url に取得
call Inquire(var = dimvar(i), growable = growable, &
& allcount = allcount, url = time_url)
call DbgMessage(' [dim=<%d>: growable=<%y>: url=<%c>]', &
& i = (/i/), L = (/growable/), c1 = trim(time_url))
! 総数 = 最後の数, なので...
if (growable) then
! 変数部分だけ分離
call UrlSplit(fullname = time_url, var = time_name)
timecount = allcount
nounlimited = .false.
endif
call Close(dimvar(i))
end do
! 探査を終了したので閉じる
call Close(var)
if (stat /= DC_NOERR) then
goto 999
end if
! 時刻部分の iorange を作成しておく.
! 格子点情報で取得されているので, 頭に "^" を付加する.
if (nounlimited) then
time_iorange = ''
else
time_iorange = trim(time_name) // GT_EQUAL // &
& GT_CIRCUMFLEX // adjustl(toChar(timecount))
end if
! iorange を指定する.
! 時刻に関しては, range が存在しない場合には
! 自動取得した最後の時刻を付加する.
! range が存在する場合, "=" が含まれなければ, gtool4 のコンマ記法
! として不適切としてエラーを生じる.
! "=" が含まれる場合, iorange としてそのまま iorange になる.
! ただし, その iorange に時刻次元が含まれない場合,
! やはり先ほど自動取得した値が付加される.
! 当然, 時刻次元が存在しない場合には付加しない.
if (.not. present_and_not_empty(range)) then
iorange = time_iorange
else
! range がコンマ記法になっているか, "=" があるかどうかで調べる
call match(GT_EQUAL, range, regex_len, regex_stat)
! コンマ記法になってない場合は無制限次元の値と判定
if (regex_stat < 0) then
cause_c = range
stat = GT_EBADGT4COMMAGRAPHY
goto 999
!!$ iorange = trim(time_name) // GT_EQUAL // adjustl(range)
else
! コンマ記法になっている場合, まずその中に無制限次元が
! 存在しているか調べ, 存在してない場合のみ time_iorange を
! 付加する.
if (trim(UrlSearchIORange(range, time_name)) /= "") then
iorange = range
else
if (trim(time_iorange) /= "") then
iorange = range // GT_COMMA // time_iorange
else
iorange = range
end if
end if
end if
endif
call DbgMessage('@ iorange=%c', c1=trim(iorange))
! file, varname, iorange を gtool変数化
! (「file@varname,time=10.5」のように)
url = UrlMerge(file, varname, '', iorange)
999 continue
call StoreError(stat, subname, err, cause_c)
call EndSub(subname, '', c1=trim(url))
end subroutine lookup_growable_url
subroutine actual_iorange_dump(url, actual_url, err)
!
! 変数 URL *url* に対応するファイル, 変数からデータを取り出す際,
! 入出力範囲指定によって切り出される値の本当の位置を
! 標準出力に出力する. *actual_url* が与えられる場合には
! その引数に値を返し, 標準出力には出力しない.
!
! HistoryGet, HistoryGetPointer が下層で呼び出している
! gtdata_generic#Get は, 入出力範囲が次元データに正確に一致しない
! 場合, 最も近い値を自動的に選択して切り出す. しかしその結果,
! 「本当はどこのデータを入力したか」がわからない場合があるため,
! このサブルーチンによって正確な位置をユーザに知らせる.
!
use dc_types, only: DP, STRING
use dc_string, only: Split, JoinChar, toChar
use dc_url, only: UrlSearchIORange, UrlMerge, UrlSplit
use dc_url, only: GT_COMMA, GT_EQUAL, GT_COLON
use dc_message, only: MessageNotify
use dc_trace, only: Beginsub, Endsub, DbgMessage
use regex, only: match
use gtdata_types, only: GT_VARIABLE
use gtdata_generic, only: Open, Close, Get
use dc_error, only: StoreError, DC_NOERR
character(*), intent(in) :: url ! 変数 URL
character(*), intent(out), optional :: actual_url
! 正確な入出力範囲指定に修正
! された変数 URL
logical, intent(out), optional :: err ! エラーのフラグ
character(STRING), pointer :: iorange_each(:) =>null()
character(STRING), pointer :: range_values(:) =>null()
character(STRING), pointer :: new_iorange_each(:) =>null()
character(STRING), pointer :: new_range_values(:) =>null()
character(STRING) :: new_url, new_iorange, url_tmp, dimname
character(STRING) :: file, varname, range, cause_c
type(GT_VARIABLE) :: var
real :: iorange_value(1)
integer :: i, j, regex_len, regex_stat, stat
character(*), parameter :: subname = "actual_iorange_dump"
continue
call BeginSub(subname, '', c1=trim(url))
new_iorange = ''
cause_c = ''
stat = DC_NOERR
call UrlSplit(url, file, varname, iorange=range)
call Split(range, iorange_each, GT_COMMA)
allocate(new_iorange_each(size(iorange_each)))
do i = 1, size(iorange_each)
call match(GT_EQUAL, iorange_each(i), regex_len, regex_stat)
if (regex_stat < 0 .or. regex_len < 2) then
new_iorange_each(i) = trim(iorange_each(i))
else
dimname = iorange_each(i)(:regex_len-1)
call Split(iorange_each(i)(regex_len+1:), range_values, GT_COLON)
allocate(new_range_values(size(range_values)))
do j = 1, size(range_values)
url_tmp = UrlMerge(file, dimname, '', &
& iorange=trim(dimname) // GT_EQUAL // trim(range_values(j)))
call Open(var, url_tmp)
call Get(var, iorange_value, 1)
call Close(var)
new_range_values(j) = toChar(iorange_value)
end do
new_iorange_each(i) = &
& trim(dimname) // GT_EQUAL // JoinChar(new_range_values, GT_COLON)
deallocate(new_range_values)
deallocate(range_values)
end if
end do
new_iorange = JoinChar(new_iorange_each, GT_COMMA)
deallocate(new_iorange_each)
deallocate(iorange_each)
new_url = UrlMerge(file, varname, '', new_iorange)
if (present(actual_url)) then
actual_url = new_url
else
call MessageNotify('M', subname, 'Input %c', c1=trim(new_url))
end if
999 continue
call StoreError(stat, subname, err, cause_c)
call EndSub(subname, '', c1=trim(new_url))
end subroutine actual_iorange_dump
__EndOfFortran90Code__
print <<"__EndOfFooter__"
!--
! vi:set readonly sw=4 ts=8:
!
#{rb2f90_emacs_readonly}!
!++
__EndOfFooter__