anvarinquire.f90

Path: src/anvarinquire.f90
Last Update: Wed Jul 20 18:22:22 JST 2005

Copyright (C) GFD Dennou Club, 2000. All rights reserved.

Methods

Included Modules

an_types an_file an_vartable an_generic dc_trace netcdf_f77

Public Instance methods

var :type(an_variable), intent(in)
ndims :integer, intent(out), optional
: 変数の次元数
dimlen :integer, intent(out), optional
: 変数が1次元である場合、次元長
growable :logical, intent(out), optional
: 変数が成長可能次元を持つか
name :character(*), intent(out), optional
: 文字型引数が短いと値の切り詰めが起こりうる ’?’ のあとの変数名
url :character(*), intent(out), optional
: 変数名、少なくともファイル名を含む、なるべく長い名前
xtype :character(*), intent(out), optional
: 変数の型名

問い合わせは型ごとに手続をわけた。

[Source]


subroutine ANVarInquire(var, ndims, dimlen, growable, name, url, xtype)

implicit none
    type(an_variable), intent(in):: var
    ! 変数の次元数
    integer, intent(out), optional:: ndims
    ! 変数が1次元である場合、次元長
    integer, intent(out), optional:: dimlen
    ! 変数が成長可能次元を持つか
    logical, intent(out), optional:: growable
    ! 文字型引数が短いと値の切り詰めが起こりうる
    ! '?' のあとの変数名
    character(*), intent(out), optional:: name
    ! 変数名、少なくともファイル名を含む、なるべく長い名前
    character(*), intent(out), optional:: url
    ! 変数の型名
    character(*), intent(out), optional:: xtype
    ! 内部変数
    type(an_variable_entry):: ent
    integer:: stat, length, i, i_xtype, idim_growable
    character(len = *), parameter:: subname = 'anvarinqurie'
    character(len = nf_max_name):: buffer
    character(len = nf_max_name):: fbuffer
continue
    call beginsub(subname, 'var.id=%d', i=(/var%id/))

    ! フェイルセーフ用にエラー値をまず入れる
    if (present(ndims)) ndims = -1
    if (present(dimlen)) dimlen = -1

    ! 変数実体の探索
    stat = vtable_lookup(var, ent)
    if (stat /= nf_noerr) then
        call endsub(subname, 'var not found')
        return
    endif

    ! 各引数が与えられている場合について値を取得する動作を

    if (present(ndims)) then
        if (associated(ent%dimids)) then
            ndims = size(ent%dimids)
        else
            ndims = 0
        endif
    endif

    if (present(dimlen)) then
        dimlen = 1
        if (ent%dimid > 0) then
            ! 実体に次元としての問い合わせが可能な場合
            stat = nf_inq_dimlen(ent%fileid, ent%dimid, dimlen)
            if (stat /= nf_noerr) then
                dimlen = -1
                call endsub(subname, 'dimlen err')
                return
            endif
        else
            ! 実体が変数として問い合わせるしかない場合
            if (associated(ent%dimids)) then
                do, i = 1, size(ent%dimids)
                    stat = nf_inq_dimlen(ent%fileid, ent%dimids(i), length)
                    if (stat /= nf_noerr) then
                        dimlen = -1
                        exit
                    endif
                    dimlen = dimlen * length
                enddo
            endif
        endif
    endif

    if (present(xtype)) then
        stat = nf_inq_vartype(ent%fileid, ent%varid, xtype=i_xtype)
        if (stat /= NF_NOERR) i_xtype = 0
        call ANXTypeName(i_xtype, xtype)
    endif

    if (present(name)) then
        call local_getname(ent, buffer)
        name = buffer
    endif

    if (present(url)) then
        call local_getname(ent, buffer)
        call DbgMessage('ent%%fileid=%d', i=(/ent%fileid/))
        call inquire(ent%fileid, name=fbuffer)
        url = trim(fbuffer) // '?' // buffer
    endif

    if (present(growable)) then
        growable = .false.
        stat = vtable_lookup(var, ent)
        if (stat /= nf_noerr) return
        stat = nf_inq_unlimdim(ent%fileid, idim_growable)
        if (stat /= nf_noerr) return

        if (ent%varid > 0) then
            if (.not. associated(ent%dimids)) return
            do, i = 1, size(ent%dimids)
                if (ent%dimids(i) == idim_growable) growable = .true.
            enddo
        else
            growable = (ent%dimid == idim_growable)
        endif
    endif

    ! 安全に終った
    call endsub(subname, 'ok')
    return

contains

    subroutine local_getname(ent, varname)

        type(an_variable_entry), intent(in):: ent
        character(len = *), intent(out):: varname
        if (ent%dimid > 0) then
            stat = nf_inq_dimname(ent%fileid, ent%dimid, varname)
        else
            stat = nf_inq_varname(ent%fileid, ent%varid, varname)
        endif
        if (stat /= NF_NOERR) varname = ""
    end subroutine

[Validate]