#!/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}!
!++
!
!== Input gtool4 netCDF data
!
! Authors:: Yasuhiro MORIKAWA
! Version:: $Id: historyget.rb2f90,v 1.4 2007/06/19 17:49:22 morikawa Exp $
! Tag Name:: $Name: gt4f90io-20070729 $
! Copyright:: Copyright (C) GFD Dennou Club, 2006. All rights reserved.
! License:: See COPYRIGHT[link:../../COPYRIGHT]
!
! 以下のサブルーチン、関数は gt4_history から gt4_history#HistoryGet
! もしくは gt4_history#HistoryGetPointer として提供されます。
__EndOfFortran90Code__
types = ["Double", "Real"]
types.each{ |type|
print <<"__EndOfFortran90Code__"
subroutine HistoryGet#{type}Ex(file, varname, array, range, quiet)
!
! 入力用内部サブルーチン (#{type} の固定長配列用).
! *file* にファイル名を, *varname* に変数名を与えます.
! 外部には提供されず, その他の外部に提供される入力用
! サブルーチンが内部的に呼び出すことを想定しています.
! *array* にはファイルから入力されたデータが返ります.
!
! *range* には gtool4 のコンマ記法を与えることで,
! 入力先のデータの一部を切り出して入力を行うことが
! 可能です. なお, range に数値のみが代入される場合,
! それは時刻の次元 (正確には netCDF の無制限次元) の値として
! 受け取られます.
!
! デフォルトでは, データの入力時にどのファイルのどの変数が
! どの次元で切り出されて入力されたのかを表示します.
! メッセージ出力が不要な場合は *quiet* に .true. を与えてください.
!
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
character(*), intent(in):: file ! ファイル名
character(*), intent(in):: varname ! 変数名
character(*), intent(in), optional:: range
! スライス用オプション.
! gtool4 変数のコンマ記法で記述
! {(例: time=100.0,x=10:20,y=30.5)}
logical, intent(in), optional:: quiet
#{$type_intent_out[type]}, intent(out):: array(*) ! データ
type(GT_VARIABLE) :: var
character(STRING) :: url, actual_url
integer:: domain ! 変数の入出力領域の大きさ
! (= 変数が依存する各次元サイズの積)
character(*), parameter :: subname = "HistoryGet#{type}Ex"
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)))
call lookup_growable_url(file, varname, url, range)
! いよいよデータ取得
call Open(var, url)
call Inquire(var=var, size=domain)
call Get(var, array, domain)
call Close(var)
call actual_iorange_dump(url, actual_url)
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}Ex
__EndOfFortran90Code__
}
types = ["Double", "Real"]
types.each{ |type|
print <<"__EndOfFortran90Code__"
subroutine HistoryGet#{type}0(file, varname, array, range, quiet)
!
#{ifelse(type, "Double", %Q{
!
!== データ入力 (固定長配列用)
!
! 既存の gtool4 データ内の変数のデータを取得します.
!
! ポインタの配列へデータを入力を行う場合は,
! HistoryGetPointer を利用してください.
!
! デフォルトでは, gtool4 データ内の変数データのうち,
! 最も新しい時刻のデータを入力します. その他の時刻データを
! 入力したい場合には, 下記の *range* オプションを利用してください.
!
! *range* には gtool4 のコンマ記法
! ({gtool4 netCDF 規約}[link:../xref.htm#label-6] 参照)
! を与えることで,
! 入力先のデータの一部を切り出して入力を行うことが可能です.
!
! *HistoryGet* は複数のサブルーチンの総称名です. *array* には
! 0 〜 #{$histget_dim} 次元のデータを与えることが可能です.
! 下記のサブルーチンを参照ください.
!
! デフォルトでは, データの入力時にどのファイルのどの変数が
! どの次元で切り出されて入力されたのかを表示します.
! メッセージ出力が不要な場合は *quiet* に .true. を与えてください.
!
})}
!
use dc_string, only: toChar
use dc_present, only: present_select, present_and_not_empty
use dc_types, only: DP, STRING
use dc_trace, only: Beginsub, Endsub, DbgMessage
character(*), intent(in):: file
#{ifelse(type, "Double", %Q{
! 取得するデータが格納された
! netCDF ファイル名
!
})}
character(*), intent(in):: varname
#{ifelse(type, "Double", %Q{
! 取得しようとするデータの変数名
!
! ここで指定する変数名は、 *file*
! で指定されるファイル内に含
! まれていなければなりません。
!
})}
character(*), intent(in), optional:: range
#{ifelse(type, "Double", %Q{
! gtool4 のコンマ記法による
! データの切り出し指定
!
! 詳しくは
! {gtool4 netCDF 規約}[link:../xref.htm#label-6]
! の「5.4 コンマ記法」を参照して
! ください。
!
! 数値のみが与えられる場合
! (※ 引数の型は常に文字型です。
! すなわち
! "10.0"
! のような文字列が与えられた場合
! を指します)、それは時刻の値
! として受け取られます。
! この場合、"^" を
! 先頭に付加することで、
! 時間のステップ数で指定する
! ことも可能です。 (例えば
! "^10"ならば
! 10 ステップ目のデータを取得します)。
!
})}
logical, intent(in), optional:: quiet
#{$type_intent_out[type]}, intent(out) :: array
#{ifelse(type, "Double", %Q{
! 取得するデータを格納する配列
!
! 型は単精度実数型か倍精度実数型
! かのどちらかです。
! 取得するデータの空間次元の
! サイズと配列のサイズとが一
! 致している必要があります。
! *varname* のデータ型と *array*
! の型が異なる場合は、自動的
! に *array* の型に変換されます。
!
})}
#{$type_internal[type]} :: array_tmp(1)
character(*), parameter :: subname = "HistoryGet#{type}0"
interface
subroutine HistoryGet#{type}Ex(file, varname, array, range, quiet)
use dc_types, only: DP
character(*), intent(in):: file ! ファイル名
character(*), intent(in):: varname ! 変数名
character(*), intent(in), optional:: range
logical, intent(in), optional:: quiet
#{$type_intent_out[type]}, intent(out):: array(*) ! データ
end subroutine HistoryGet#{type}Ex
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)))
if (present_and_not_empty(range)) then
call HistoryGet#{type}Ex(file=file, varname=varname, &
& array=array_tmp, range=range, quiet=quiet)
else
call HistoryGet#{type}Ex(file=file, varname=varname, &
& array=array_tmp, quiet=quiet)
end if
array = array_tmp(1)
call EndSub(subname)
end subroutine HistoryGet#{type}0
__EndOfFortran90Code__
}
types = ["Double", "Real"]
types.each{ |type|
for num in 1..$histget_dim
print <<"__EndOfFortran90Code__"
subroutine HistoryGet#{type}#{num}(file, varname, array, range, quiet)
!
! #{type} 型固定配列用 #{num} 次元配列用データ入力サブルーチン
!
!--
! file にファイル名を, varname に変数名を与える.
! array には #{num} 次元配列を与える.
! range には gtool4 のコンマ記法を与えることで,
! 入力先のデータの一部を切り出して入力を行うことが
! 可能である. なお, range に数値のみが代入される場合,
! それは時刻の次元 (正確には netCDF の無制限次元) の値として
! 受け取られる.
!++
use dc_string, only: toChar
use dc_present,only: present_select, present_and_not_empty
use dc_types, only: DP, STRING
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]}, intent(out) :: array#{array_colon("#{num}")}
character(*), parameter :: subname = "HistoryGet#{type}#{num}"
interface
subroutine HistoryGet#{type}Ex(file, varname, array, range, quiet)
use dc_types, only: DP
character(*), intent(in):: file ! ファイル名
character(*), intent(in):: varname ! 変数名
character(*), intent(in), optional:: range
logical, intent(in), optional:: quiet
#{$type_intent_out[type]}, intent(out):: array(*) ! データ
end subroutine HistoryGet#{type}Ex
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)))
if (present_and_not_empty(range)) then
call HistoryGet#{type}Ex(file=file, varname=varname, &
& array=array, range=range, quiet=quiet)
else
call HistoryGet#{type}Ex(file=file, varname=varname, &
& array=array, quiet=quiet)
endif
call EndSub(subname)
end subroutine HistoryGet#{type}#{num}
__EndOfFortran90Code__
end
}
types = ["Double", "Real"]
types.each{ |type|
print <<"__EndOfFortran90Code__"
subroutine HistoryGet#{type}0TimeR(file, varname, array, time, quiet)
!
! 将来廃止予定です。*time* は *range* に置き換え、文字型で与えて下さい。
!
use dc_string, only: toChar
use dc_types, only: DP, STRING
use dc_trace, only: Beginsub, Endsub, DbgMessage
character(*), intent(in):: file, varname
real, intent(in) :: time
logical, intent(in), optional:: quiet
#{$type_intent_out[type]}, intent(out) :: array
#{$type_internal[type]} :: array_tmp(1)
character(*), parameter :: subname = "HistoryGet#{type}0TimeR"
interface
subroutine HistoryGet#{type}Ex(file, varname, array, range, quiet)
use dc_types, only: DP
character(*), intent(in):: file ! ファイル名
character(*), intent(in):: varname ! 変数名
character(*), intent(in), optional:: range
logical, intent(in), optional:: quiet
#{$type_intent_out[type]}, intent(out):: array(*) ! データ
end subroutine HistoryGet#{type}Ex
end interface
continue
call BeginSub(subname, 'file=%c varname=%c time=%c', &
& c1=trim(file), c2=trim(varname), c3=toChar(time))
call HistoryGet#{type}Ex(file=file, varname=varname, &
& array=array_tmp, range=toChar(time), quiet=quiet)
call EndSub(subname)
array = array_tmp(1)
end subroutine
__EndOfFortran90Code__
}
types = ["Double", "Real"]
types.each{ |type|
print <<"__EndOfFortran90Code__"
subroutine HistoryGet#{type}0TimeD(file, varname, array, time, quiet)
!
! 将来廃止予定です。*time* は *range* に置き換え、文字型で与えて下さい。
!
use dc_string, only: toChar
use dc_types, only: DP
use dc_trace, only: Beginsub, Endsub, DbgMessage
character(*), intent(in):: file, varname
real(DP), intent(in) :: time
logical, intent(in), optional:: quiet
#{$type_intent_out[type]}, intent(out) :: array
#{$type_internal[type]} :: array_tmp(1)
character(*), parameter :: subname = "HistoryGet#{type}0TimeD"
interface
subroutine HistoryGet#{type}Ex(file, varname, array, range, quiet)
use dc_types, only: DP
character(*), intent(in):: file ! ファイル名
character(*), intent(in):: varname ! 変数名
character(*), intent(in), optional:: range
logical, intent(in), optional:: quiet
#{$type_intent_out[type]}, intent(out):: array(*) ! データ
end subroutine HistoryGet#{type}Ex
end interface
continue
call BeginSub(subname, 'file=%c varname=%c time=%c', &
& c1=trim(file), c2=trim(varname), c3=toChar(time))
call HistoryGet#{type}Ex(file=file, varname=varname, &
& array=array_tmp, range=toChar(time), quiet=quiet)
call EndSub(subname)
array = array_tmp(1)
end subroutine
__EndOfFortran90Code__
}
types = ["Double", "Real"]
types.each{ |type|
for num in 1..$histget_dim
print <<"__EndOfFortran90Code__"
subroutine HistoryGet#{type}#{num}TimeR(file, varname, array, time, quiet)
!
! 将来廃止予定です。*time* は *range* に置き換え、文字型で与えて下さい。
!
use dc_string, only: toChar
use dc_types, only: DP, STRING
use dc_trace, only: Beginsub, Endsub, DbgMessage
character(*), intent(in):: file, varname
real, intent(in) :: time
logical, intent(in), optional:: quiet
#{$type_intent_out[type]}, intent(out) :: array#{array_colon("#{num}")}
character(*), parameter :: subname = "HistoryGet#{type}#{num}TimeR"
interface
subroutine HistoryGet#{type}Ex(file, varname, array, range, quiet)
use dc_types, only: DP
character(*), intent(in):: file ! ファイル名
character(*), intent(in):: varname ! 変数名
character(*), intent(in), optional:: range
logical, intent(in), optional:: quiet
#{$type_intent_out[type]}, intent(out):: array(*) ! データ
end subroutine HistoryGet#{type}Ex
end interface
continue
call BeginSub(subname, 'file=%c varname=%c time=%c', &
& c1=trim(file), c2=trim(varname), c3=toChar(time))
call HistoryGet#{type}Ex(file=file, varname=varname, &
& array=array, range=toChar(time), quiet=quiet)
call EndSub(subname)
end subroutine
__EndOfFortran90Code__
end
}
types = ["Double", "Real"]
types.each{ |type|
for num in 1..$histget_dim
print <<"__EndOfFortran90Code__"
subroutine HistoryGet#{type}#{num}TimeD(file, varname, array, time, quiet)
!
! 将来廃止予定です。*time* は *range* に置き換え、文字型で与えて下さい。
!
use dc_string, only: toChar
use dc_types, only: DP, STRING
use dc_trace, only: Beginsub, Endsub, DbgMessage
character(*), intent(in):: file, varname
real(DP), intent(in) :: time
logical, intent(in), optional:: quiet
#{$type_intent_out[type]}, intent(out) :: array#{array_colon("#{num}")}
character(*), parameter :: subname = "HistoryGet#{type}#{num}TimeD"
interface
subroutine HistoryGet#{type}Ex(file, varname, array, range, quiet)
use dc_types, only: DP
character(*), intent(in):: file ! ファイル名
character(*), intent(in):: varname ! 変数名
character(*), intent(in), optional:: range
logical, intent(in), optional:: quiet
#{$type_intent_out[type]}, intent(out):: array(*) ! データ
end subroutine HistoryGet#{type}Ex
end interface
continue
call BeginSub(subname, 'file=%c varname=%c time=%c', &
& c1=trim(file), c2=trim(varname), c3=toChar(time))
call HistoryGet#{type}Ex(file=file, varname=varname, &
& array=array, range=toChar(time), quiet=quiet)
call EndSub(subname)
end subroutine
__EndOfFortran90Code__
end
}
types = ["Double", "Real"]
types.each{ |type|
print <<"__EndOfFortran90Code__"
subroutine HistoryGet#{type}0Pointer(file, varname, array, range, quiet)
!
#{ifelse(type, "Double", %Q{
!
!== データ入力 (ポインタ配列用)
!
! 既存の gtool4 データ内の変数のデータをポインタ配列へ取得します.
!
! 固定長の配列へデータを入力を行う場合は,
! HistoryGet を利用してください.
!
! *array* は必ず空状態で与えてください.
! すなわち初期値 =>null() を設定するか nullify を用いてください.
! 既に割り付けられている場合, もしくは不定状態の場合には
! エラーを返します.
!
! デフォルトでは, gtool4 データ内の変数データのうち,
! 最も新しい時刻のデータを入力します. その他の時刻データを
! 入力したい場合には, 下記の *range* オプションを利用してください.
!
! *range* には gtool4 のコンマ記法
! ({gtool4 netCDF 規約}[link:../xref.htm#label-6] 参照)
! を与えることで,
! 入力先のデータの一部を切り出して入力を行うことが可能です.
!
! *HistoryGetPointer* は複数のサブルーチンの総称名です. *array* には
! 0 〜 #{$histget_dim} 次元のデータを与えることが可能です.
! 下記のサブルーチンを参照ください.
!
! デフォルトでは, データの入力時にどのファイルのどの変数が
! どの次元で切り出されて入力されたのかを表示します.
! メッセージ出力が不要な場合は *quiet* に .true. を与えてください.
!
})}
!
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
#{ifelse(type, "Double", %Q{
! 取得するデータが格納された
! netCDF ファイル名
!
})}
character(*), intent(in):: varname
#{ifelse(type, "Double", %Q{
! 取得しようとするデータの変数名
!
! ここで指定する変数名は、 *file*
! で指定されるファイル内に含
! まれていなければなりません。
!
})}
character(*), intent(in), optional:: range
#{ifelse(type, "Double", %Q{
! gtool4 のコンマ記法による
! データの切り出し指定
!
! 詳しくは
! {gtool4 netCDF 規約}[link:../xref.htm#label-6]
! の「5.4 コンマ記法」を参照して
! ください。
!
! 数値のみが与えられる場合
! (※ 引数の型は常に文字型です。
! すなわち
! "10.0"
! のような文字列が与えられた場合
! を指します)、それは時刻の値
! として受け取られます。
! この場合、"^" を
! 先頭に付加することで、
! 時間のステップ数で指定する
! ことも可能です。 (例えば
! "^10"ならば
! 10 ステップ目のデータを取得します)。
!
})}
logical, intent(in), optional:: quiet
#{$type_intent_out[type]}, pointer :: array
#{ifelse(type, "Double", %Q{
! (out)
! 取得するデータを格納する
! ポインタの配列
!
! 必ず空状態で与えてください。
! すなわち初期値 =>null() を
! 設定するか nullify を用い
! ください。
!
! 型は単精度実数型か倍精度実数型
! かのどちらかです。
! 取得するデータの空間次元の
! サイズに合わせ、配列のサイズは
! 自動的に設定されます。
! *varname* のデータ型と *array*
! の型が異なる場合は、自動的
! に *array* の型に変換されます。
!
})}
#{$type_internal[type]}, target :: array_tmp(1)
type(GT_VARIABLE) :: var
character(STRING) :: url, actual_url
integer:: domain ! 変数の入出力領域の大きさ
! (= 変数が依存する各次元サイズの積)
character(*), parameter :: subname = "HistoryGet#{type}0Pointer"
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)
allocate(array)
! いよいよデータ取得
call Open(var, url)
call Inquire(var=var, size=domain)
call Get(var, array_tmp, domain)
call Close(var)
array = array_tmp(1)
call actual_iorange_dump(url, actual_url)
if ( .not. present_and_true(quiet) ) then
call MessageNotify('M', subname, 'Input %c', c1=trim(actual_url))
end if
call EndSub(subname)
end subroutine
__EndOfFortran90Code__
}
types = ["Double", "Real"]
types.each{ |type|
print <<"__EndOfFortran90Code__"
subroutine HistoryGet#{type}0PointerTimeR(file, varname, array, time, quiet)
!
! 将来廃止予定です。*time* は *range* に置き換え、文字型で与えて下さい。
!
use gt4_history, only: HistoryGetPointer
use dc_string, only: toChar
use dc_types, only: DP, STRING
use dc_trace, only: Beginsub, Endsub, DbgMessage
character(*), intent(in) :: file, varname
real, intent(in) :: time
logical, intent(in), optional:: quiet
#{$type_intent_out[type]}, pointer :: array ! (out)
character(*), parameter :: subname = "HistoryGet#{type}0PointerTimeR"
continue
call BeginSub(subname, 'file=%c varname=%c time=%c', &
& c1=trim(file), c2=trim(varname), c3=toChar(time))
call HistoryGetPointer(file=file, &
& varname=varname, array=array, range=toChar(time), quiet=quiet)
call EndSub(subname)
end subroutine
__EndOfFortran90Code__
}
types = ["Double", "Real"]
types.each{ |type|
print <<"__EndOfFortran90Code__"
subroutine HistoryGet#{type}0PointerTimeD(file, varname, array, time, quiet)
!
! 将来廃止予定です。*time* は *range* に置き換え、文字型で与えて下さい。
!
use gt4_history, only: HistoryGetPointer
use dc_string, only: toChar
use dc_types, only: DP, STRING
use dc_trace, only: Beginsub, Endsub, DbgMessage
character(*), intent(in) :: file, varname
real(DP), intent(in) :: time
logical, intent(in), optional:: quiet
#{$type_intent_out[type]}, pointer :: array ! (out)
character(*), parameter :: subname = "HistoryGet#{type}0PointerTimeD"
continue
call BeginSub(subname, 'file=%c varname=%c time=%c', &
& c1=trim(file), c2=trim(varname), c3=toChar(time))
call HistoryGetPointer(file=file, &
& varname=varname, array=array, range=toChar(time), quiet=quiet)
call EndSub(subname)
end subroutine
__EndOfFortran90Code__
}
types = ["Double", "Real"]
types.each{ |type|
for num in 1..$histget_dim
print <<"__EndOfFortran90Code__"
subroutine HistoryGet#{type}#{num}Pointer(file, varname, array, range, quiet)
!
! #{type} 型ポインタ配列用 #{num} 次元配列用データ入力サブルーチン
!
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, varname
character(*), intent(in), optional:: range
logical, intent(in), optional:: quiet
#{$type_intent_out[type]}, pointer :: array#{array_colon("#{num}")} ! (out)
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)
call DbgMessage('@ url =%c', c1=trim(url))
! いよいよデータ取得
call Open(var, url)
call Get(var, array)
call Close(var)
call actual_iorange_dump(url, actual_url)
if ( .not. present_and_true(quiet) ) then
call MessageNotify('M', subname, 'Input %c', c1=trim(actual_url))
end if
call EndSub(subname)
end subroutine
__EndOfFortran90Code__
end
}
types = ["Double", "Real"]
types.each{ |type|
for num in 1..$histget_dim
print <<"__EndOfFortran90Code__"
subroutine HistoryGet#{type}#{num}PointerTimeR(file, varname, array, time, quiet)
!
! 将来廃止予定です。*time* は *range* に置き換え、文字型で与えて下さい。
!
use gt4_history, only: HistoryGetPointer
use dc_string, only: toChar
use dc_types, only: DP, STRING
use dc_trace, only: Beginsub, Endsub, DbgMessage
character(*), intent(in):: file, varname
real, intent(in) :: time
logical, intent(in), optional:: quiet
#{$type_intent_out[type]}, pointer :: array#{array_colon("#{num}")} ! (out)
character(*), parameter :: subname = "HistoryGet#{type}#{num}PointerTimeR"
continue
call BeginSub(subname, 'file=%c varname=%c time=%c', &
& c1=trim(file), c2=trim(varname), c3=toChar(time))
call HistoryGetPointer( &
& file=file, varname=varname, array=array, range=toChar(time), quiet=quiet)
call EndSub(subname)
end subroutine
__EndOfFortran90Code__
end
}
types = ["Double", "Real"]
types.each{ |type|
for num in 1..$histget_dim
print <<"__EndOfFortran90Code__"
subroutine HistoryGet#{type}#{num}PointerTimeD(file, varname, array, time, quiet)
!
! 将来廃止予定です。*time* は *range* に置き換え、文字型で与えて下さい。
!
use gt4_history, only: HistoryGetPointer
use dc_string, only: toChar
use dc_types, only: DP, STRING
use dc_trace, only: Beginsub, Endsub, DbgMessage
character(*), intent(in):: file, varname
real(DP), intent(in) :: time
logical, intent(in), optional:: quiet
#{$type_intent_out[type]}, pointer :: array#{array_colon("#{num}")} ! (out)
character(*), parameter :: subname = "HistoryGet#{type}#{num}PointerTimeD"
continue
call BeginSub(subname, 'file=%c varname=%c time=%c', &
& c1=trim(file), c2=trim(varname), c3=toChar(time))
call HistoryGetPointer( &
& file=file, varname=varname, array=array, range=toChar(time), quiet=quiet)
call EndSub(subname)
end subroutine
__EndOfFortran90Code__
end
}
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 部分に付加する.
! * 数値のみが与えられる場合, 時間次元の値として認識され,
! その値を用いた url が返る.
! * range 内に時間次元が設定されていない場合は, 自動的に
! 時間次元に関する iorange ("time=0.5") が指定される.
!
use gtdata_types, only: GT_VARIABLE
use gtdata_generic, only: Open, Close, Inquire
use dc_present,only: present_select, present_and_not_empty
use dc_string, only: toChar
use dc_error, only: StoreError, DC_NOERR, GT_ENOUNLIMITDIM, NF_EINVAL
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 = ""
! 引数の正当性をチェック
if (.not. present_and_not_empty(file)) then
stat = NF_EINVAL
cause_c = '"file" is not specified'
url = ""
goto 999
elseif (.not. present_and_not_empty(varname)) then
stat = NF_EINVAL
cause_c = '"varname" is not specified'
url = ""
goto 999
end if
! 時刻次元の変数名, およびその最終時刻の
! 探査のために file@varname を open (まだデータを取得しない)
call Open(var, UrlMerge(file, varname))
! 次元の数を取得
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 が存在する場合, "=" が含まれなければ単に時刻の
! 値として取得.
! "=" が含まれる場合, 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
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__