! != netCDF ファイルへ次元変数作成 ! ! Authors:: Eizi TOYODA, Yasuhiro MORIKAWA ! Version:: $Id: anvarcreated.f90,v 1.3 2006/06/07 16:35:22 morikawa Exp $ ! Tag Name:: $Name: gt4f90io-20070116 $ ! Copyright:: Copyright (C) GFD Dennou Club, 2000-2005. All rights reserved. ! License:: See COPYRIGHT[link:../../COPYRIGHT] ! ! 以下のサブルーチン、関数は an_generic から an_generic#Create ! として提供されます。 subroutine ANVarCreateD(var, url, xtype, length, overwrite, err) ! !== 次元変数作成 ! ! 変数 URL *url* に次元変数を作成します. ! 次元変数の長さを *length* に与えます. ! 返される引数 *var* には変数 ID などの情報が格納されます. ! ! *overwrite* に .true. を設定すると上書き可能なモードになります. ! デフォルトは上書き不可です. ! *err* を与える場合, 次元変数生成時にエラーが生じても ! プログラムを終了せず, *err* に .false. が返ります. ! use an_types, only: an_variable, an_variable_search use an_vartable, only: vtable_add use dc_string, only: strieq use dc_types, only: string use dc_url, only: UrlSplit use dc_trace, only: BeginSub, EndSub, DbgMessage use netcdf_f77, only: nf_noerr, nf_real, nf_int, nf_double, & & nf_def_var, nf_def_dim use an_file, only: anfileopen, anfiledefinemode use dc_error, only: StoreError, gt_enomem implicit none type(an_variable), intent(out):: var character(len = *), intent(in):: url character(len = *), intent(in):: xtype integer, intent(in):: length logical, intent(in), optional:: overwrite logical, intent(out), optional:: err type(an_variable_search):: ent character(len = string):: filename, varname, cause_c integer:: stat integer:: nc_xtype character(len = *), parameter:: subname = "ANVarCreateD" continue call BeginSub(subname, 'url=<%c>, xtype=<%c>, length=<%d>', & & c1=trim(url), c2=trim(xtype), i=(/length/)) cause_c = trim(url) ! ! --- ファイルを用意 --- call UrlSplit(url, file=filename, var=varname) call ANFileOpen(ent%fileid, filename, stat=stat, writable=.TRUE., & & overwrite=overwrite) if (stat /= NF_NOERR) goto 999 stat = ANFileDefineMode(ent%fileid) if (stat /= NF_NOERR) goto 999 ! ! --- 型の決定 --- nc_xtype = NF_REAL if (strieq(xtype, "double") .or. & & strieq(xtype, "DOUBLEPRECISION")) then nc_xtype = NF_DOUBLE endif if (strieq(xtype, "int") .or. strieq(xtype, "INTEGER")) then nc_xtype = NF_INT endif ! ! --- 次元変数の作成 --- stat = nf_def_dim(ent%fileid, trim(varname), len=length, & & dimid=ent%dimid) if (stat /= NF_NOERR) goto 999 stat = nf_def_var(ent%fileid, trim(varname), & & xtype=nc_xtype, ndims=1, dimids=(/ent%dimid/), varid=ent%varid) if (stat /= NF_NOERR) goto 999 ! stat = vtable_add(var, ent) if (stat /= NF_NOERR) goto 999 999 continue call StoreError(stat, subname, err, cause_c=cause_c) if (stat /= NF_NOERR) var = an_variable(-1) call EndSub(subname, 'stat=%d', i=(/stat/)) end subroutine ANVarCreateD