anvargetnum.f90

Path: anvargetnum.f90
Last Update: Thu Sep 07 18:21:40 JST 2006

Get AN_VARIABLES

Authors:Yasuhiro MORIKAWA, Eizi TOYODA
Version:$Id: anvargetnum.f90,v 1.5 2006/09/07 09:21:40 morikawa Exp $
Tag Name:$Name: gt4f90io-20070719 $
Copyright:Copyright (C) GFD Dennou Club, 2000-2005. All rights reserved.
License:See COPYRIGHT

以下のサブルーチン、関数は an_generic から gtdata_generic#Get として提供されます。

Methods

Included Modules

an_types an_vartable netcdf_f77 dc_types dc_trace

Public Instance methods

Subroutine :
var :type(AN_VARIABLE), intent(in)
start(:) :integer, intent(in)
cnt(:) :integer, intent(in)
stride(:) :integer, intent(in)
imap(:) :integer, intent(in)
: NetCDF変数と内部データ配列のメモリ内構 造間のマッピングを指定する整数ベクトル. 詳しくは NetCDF マニュアル (NF_PUT_VARM_type 等 を参照のこと)
siz :integer, intent(in)
value(siz) :real(DP), intent(out)
iostat :integer, intent(out)

[Source]

subroutine ANVarGetDouble(var, start, cnt, stride, imap, siz, value, iostat)
  use an_types,    only: AN_VARIABLE
  use an_vartable, only: AN_VARIABLE_ENTRY, vtable_lookup
  use netcdf_f77,  only: NF_NOERR, NF_EINVAL, NF_EINDEFINE, nf_get_varm_Double, nf_get_var1_Double, nf_redef, nf_enddef
  use dc_types,    only: DP
  use dc_trace,    only: BeginSub, EndSub, DbgMessage
  implicit none
  type(AN_VARIABLE), intent(in):: var
  integer, intent(in):: start(:)
  integer, intent(in):: cnt(:)
  integer, intent(in):: stride(:)
  integer, intent(in):: imap(:)
                              ! NetCDF変数と内部データ配列のメモリ内構
                              ! 造間のマッピングを指定する整数ベクトル.
                              ! 詳しくは NetCDF マニュアル
                              ! (NF_PUT_VARM_type 等 を参照のこと)
  integer, intent(in):: siz
  real(DP), intent(out):: value(siz)
  integer, intent(out):: iostat
  integer:: nd ! var が保持する変数が依存する次元変数の数
  integer:: ipos, i
  integer:: defstat
  type(AN_VARIABLE_ENTRY):: ent
  integer, allocatable:: istart(:), istride(:), iimap(:)
continue
  call BeginSub('ANVarGetDouble', fmt='varmap=%d, start=%*d, cnt=%*d, stride=%*d, imap=%*d siz=%d', i=(/var%id, start(:), cnt(:), stride(:), imap(:), siz/), n=(/size(start), size(cnt), size(stride), size(imap)/))
  iostat = vtable_lookup(var, ent)
  if (iostat /= NF_NOERR) goto 999
  ! --- nd check ---
  nd = 0
  if (associated(ent%dimids)) nd = size(ent%dimids)
  if (min(size(start), size(cnt), size(stride), size(imap)) < nd) then
    iostat = NF_EINVAL
    goto 999
  endif
  if (nd == 0) then
    iostat = nf_get_var1_Double(ent%fileid, ent%varid, start, value(1))
    goto 999
  endif
  ! --- stride kakikae buffer ---
  allocate(istart(nd), istride(nd), iimap(nd))
  istart(1:nd) = start(1:nd)
  istride(1:nd) = stride(1:nd)
  iimap(1:nd) = imap(1:nd)
  ipos = 1
  ! --- do read ---
  if (ent%varid <= 0 .or. count(cnt(1:nd) == 0) > 0) then
    call BeginSub('fake_map_get')
    call fake_map_get
    call EndSub('fake_map_get', 'iostat=%d', i=(/iostat/))
  elseif (count(stride(1:nd) < 0) <= 0) then
    iostat = nf_get_varm_Double(ent%fileid, ent%varid, istart, cnt, istride, iimap, value)
    if (iostat == NF_EINDEFINE) then
      defstat = nf_enddef(ent%fileid)
      iostat = nf_get_varm_Double(ent%fileid, ent%varid, istart, cnt, istride, iimap, value)
      defstat = nf_redef(ent%fileid)
    end if
  else
    ! negative stride is not allowed for netcdf
    do, i = 1, nd
      if (stride(i) > 0) cycle
      ipos = ipos + (cnt(i) - 1) * imap(i)
      istart(i) = start(i) + (cnt(i) - 1) * stride(i)
      istride(i) = -stride(i)
      iimap(i) = -imap(i)
      call DbgMessage('dim %d negate: stride->%d start->%d map->%d', i=(/i, istride(i), istart(i), iimap(i)/))
    enddo
    iostat = nf_get_varm_Double(ent%fileid, ent%varid, istart, cnt, istride, iimap, value(ipos))
  endif
  deallocate(istart, istride, iimap)
999 continue
  call EndSub('ANVarGetDouble', 'iostat=%d', i=(/iostat/))
  return
contains

  subroutine fake_map_get
    !
    ! ent%varid <= 0 の場合, または整数型配列 cnt に 0 が含まれる
    ! 場合に呼ばれるサブルーチン.
    !
    ! 本来は nf_get_varm_Double で行う 1 次元データの入力を, 
    ! nf_get_var1_Double で 1 要素ごとに入力し, value に格納する.
    !
    ! iostat にはデータ入出力時のステータスが返る.
    !
    use netcdf_f77, only: NF_EINDEFINE, nf_enddef, nf_redef, nf_get_var1_Double
    implicit none
    integer:: ofs(nd), here(nd)
    integer:: j, defstat
  continue
    iostat = NF_NOERR
    ofs(1:nd) = 0
    do
      j = ipos + dot_product(ofs(1:nd), imap(1:nd))
      here(1:nd) = istart(1:nd) + ofs(1:nd) * istride(1:nd)
      if (j < 1 .or. j > siz) then
        iostat = NF_EINVAL
        call DbgMessage('nf_get_var1_Double(ncid=%d, varid=%d, indx=[%*d], out-ofs=%d)', i=(/ent%fileid, ent%varid, here(1:nd), j/), n=(/nd/))
        return
      endif
      if (ent%varid == 0) then
        value(j) = j
        iostat = NF_NOERR
      else
        iostat = nf_get_var1_Double(ent%fileid, ent%varid, here(1), value(j))
        if (iostat == NF_EINDEFINE) then
          defstat = nf_enddef(ent%fileid)
          iostat = nf_get_var1_Double(ent%fileid, ent%varid, here(1), value(j))
          defstat = nf_redef(ent%fileid)
        end if
      endif
      if (iostat /= NF_NOERR) return
      ofs(1) = ofs(1) + 1
      do, j = 1, nd - 1
        if (ofs(j) < cnt(j)) exit
        ofs(j) = 0
        ofs(j + 1) = ofs(j + 1) + 1
      enddo
      if (ofs(nd) >= cnt(nd)) exit
    enddo
  end subroutine fake_map_get

end subroutine ANVarGetDouble
Subroutine :
var :type(AN_VARIABLE), intent(in)
start(:) :integer, intent(in)
cnt(:) :integer, intent(in)
stride(:) :integer, intent(in)
imap(:) :integer, intent(in)
: NetCDF変数と内部データ配列のメモリ内構 造間のマッピングを指定する整数ベクトル. 詳しくは NetCDF マニュアル (NF_PUT_VARM_type 等 を参照のこと)
siz :integer, intent(in)
value(siz) :real, intent(out)
iostat :integer, intent(out)

[Source]

subroutine ANVarGetReal(var, start, cnt, stride, imap, siz, value, iostat)
  use an_types,    only: AN_VARIABLE
  use an_vartable, only: AN_VARIABLE_ENTRY, vtable_lookup
  use netcdf_f77,  only: NF_NOERR, NF_EINVAL, NF_EINDEFINE, nf_get_varm_Real, nf_get_var1_Real, nf_redef, nf_enddef
  use dc_types,    only: DP
  use dc_trace,    only: BeginSub, EndSub, DbgMessage
  implicit none
  type(AN_VARIABLE), intent(in):: var
  integer, intent(in):: start(:)
  integer, intent(in):: cnt(:)
  integer, intent(in):: stride(:)
  integer, intent(in):: imap(:)
                              ! NetCDF変数と内部データ配列のメモリ内構
                              ! 造間のマッピングを指定する整数ベクトル.
                              ! 詳しくは NetCDF マニュアル
                              ! (NF_PUT_VARM_type 等 を参照のこと)
  integer, intent(in):: siz
  real, intent(out):: value(siz)
  integer, intent(out):: iostat
  integer:: nd ! var が保持する変数が依存する次元変数の数
  integer:: ipos, i
  integer:: defstat
  type(AN_VARIABLE_ENTRY):: ent
  integer, allocatable:: istart(:), istride(:), iimap(:)
continue
  call BeginSub('ANVarGetReal', fmt='varmap=%d, start=%*d, cnt=%*d, stride=%*d, imap=%*d siz=%d', i=(/var%id, start(:), cnt(:), stride(:), imap(:), siz/), n=(/size(start), size(cnt), size(stride), size(imap)/))
  iostat = vtable_lookup(var, ent)
  if (iostat /= NF_NOERR) goto 999
  ! --- nd check ---
  nd = 0
  if (associated(ent%dimids)) nd = size(ent%dimids)
  if (min(size(start), size(cnt), size(stride), size(imap)) < nd) then
    iostat = NF_EINVAL
    goto 999
  endif
  if (nd == 0) then
    iostat = nf_get_var1_Real(ent%fileid, ent%varid, start, value(1))
    goto 999
  endif
  ! --- stride kakikae buffer ---
  allocate(istart(nd), istride(nd), iimap(nd))
  istart(1:nd) = start(1:nd)
  istride(1:nd) = stride(1:nd)
  iimap(1:nd) = imap(1:nd)
  ipos = 1
  ! --- do read ---
  if (ent%varid <= 0 .or. count(cnt(1:nd) == 0) > 0) then
    call BeginSub('fake_map_get')
    call fake_map_get
    call EndSub('fake_map_get', 'iostat=%d', i=(/iostat/))
  elseif (count(stride(1:nd) < 0) <= 0) then
    iostat = nf_get_varm_Real(ent%fileid, ent%varid, istart, cnt, istride, iimap, value)
    if (iostat == NF_EINDEFINE) then
      defstat = nf_enddef(ent%fileid)
      iostat = nf_get_varm_Real(ent%fileid, ent%varid, istart, cnt, istride, iimap, value)
      defstat = nf_redef(ent%fileid)
    end if
  else
    ! negative stride is not allowed for netcdf
    do, i = 1, nd
      if (stride(i) > 0) cycle
      ipos = ipos + (cnt(i) - 1) * imap(i)
      istart(i) = start(i) + (cnt(i) - 1) * stride(i)
      istride(i) = -stride(i)
      iimap(i) = -imap(i)
      call DbgMessage('dim %d negate: stride->%d start->%d map->%d', i=(/i, istride(i), istart(i), iimap(i)/))
    enddo
    iostat = nf_get_varm_Real(ent%fileid, ent%varid, istart, cnt, istride, iimap, value(ipos))
  endif
  deallocate(istart, istride, iimap)
999 continue
  call EndSub('ANVarGetReal', 'iostat=%d', i=(/iostat/))
  return
contains

  subroutine fake_map_get
    !
    ! ent%varid <= 0 の場合, または整数型配列 cnt に 0 が含まれる
    ! 場合に呼ばれるサブルーチン.
    !
    ! 本来は nf_get_varm_Real で行う 1 次元データの入力を, 
    ! nf_get_var1_Real で 1 要素ごとに入力し, value に格納する.
    !
    ! iostat にはデータ入出力時のステータスが返る.
    !
    use netcdf_f77, only: NF_EINDEFINE, nf_enddef, nf_redef, nf_get_var1_Real
    implicit none
    integer:: ofs(nd), here(nd)
    integer:: j, defstat
  continue
    iostat = NF_NOERR
    ofs(1:nd) = 0
    do
      j = ipos + dot_product(ofs(1:nd), imap(1:nd))
      here(1:nd) = istart(1:nd) + ofs(1:nd) * istride(1:nd)
      if (j < 1 .or. j > siz) then
        iostat = NF_EINVAL
        call DbgMessage('nf_get_var1_Real(ncid=%d, varid=%d, indx=[%*d], out-ofs=%d)', i=(/ent%fileid, ent%varid, here(1:nd), j/), n=(/nd/))
        return
      endif
      if (ent%varid == 0) then
        value(j) = j
        iostat = NF_NOERR
      else
        iostat = nf_get_var1_Real(ent%fileid, ent%varid, here(1), value(j))
        if (iostat == NF_EINDEFINE) then
          defstat = nf_enddef(ent%fileid)
          iostat = nf_get_var1_Real(ent%fileid, ent%varid, here(1), value(j))
          defstat = nf_redef(ent%fileid)
        end if
      endif
      if (iostat /= NF_NOERR) return
      ofs(1) = ofs(1) + 1
      do, j = 1, nd - 1
        if (ofs(j) < cnt(j)) exit
        ofs(j) = 0
        ofs(j + 1) = ofs(j + 1) + 1
      enddo
      if (ofs(nd) >= cnt(nd)) exit
    enddo
  end subroutine fake_map_get

end subroutine ANVarGetReal

Private Instance methods

Subroutine :

ent%varid <= 0 の場合, または整数型配列 cnt に 0 が含まれる 場合に呼ばれるサブルーチン.

本来は nf_get_varm_Double で行う 1 次元データの入力を, nf_get_var1_Double で 1 要素ごとに入力し, value に格納する.

iostat にはデータ入出力時のステータスが返る.

[Source]

  subroutine fake_map_get
    !
    ! ent%varid <= 0 の場合, または整数型配列 cnt に 0 が含まれる
    ! 場合に呼ばれるサブルーチン.
    !
    ! 本来は nf_get_varm_Double で行う 1 次元データの入力を, 
    ! nf_get_var1_Double で 1 要素ごとに入力し, value に格納する.
    !
    ! iostat にはデータ入出力時のステータスが返る.
    !
    use netcdf_f77, only: NF_EINDEFINE, nf_enddef, nf_redef, nf_get_var1_Double
    implicit none
    integer:: ofs(nd), here(nd)
    integer:: j, defstat
  continue
    iostat = NF_NOERR
    ofs(1:nd) = 0
    do
      j = ipos + dot_product(ofs(1:nd), imap(1:nd))
      here(1:nd) = istart(1:nd) + ofs(1:nd) * istride(1:nd)
      if (j < 1 .or. j > siz) then
        iostat = NF_EINVAL
        call DbgMessage('nf_get_var1_Double(ncid=%d, varid=%d, indx=[%*d], out-ofs=%d)', i=(/ent%fileid, ent%varid, here(1:nd), j/), n=(/nd/))
        return
      endif
      if (ent%varid == 0) then
        value(j) = j
        iostat = NF_NOERR
      else
        iostat = nf_get_var1_Double(ent%fileid, ent%varid, here(1), value(j))
        if (iostat == NF_EINDEFINE) then
          defstat = nf_enddef(ent%fileid)
          iostat = nf_get_var1_Double(ent%fileid, ent%varid, here(1), value(j))
          defstat = nf_redef(ent%fileid)
        end if
      endif
      if (iostat /= NF_NOERR) return
      ofs(1) = ofs(1) + 1
      do, j = 1, nd - 1
        if (ofs(j) < cnt(j)) exit
        ofs(j) = 0
        ofs(j + 1) = ofs(j + 1) + 1
      enddo
      if (ofs(nd) >= cnt(nd)) exit
    enddo
  end subroutine fake_map_get
Subroutine :

ent%varid <= 0 の場合, または整数型配列 cnt に 0 が含まれる 場合に呼ばれるサブルーチン.

本来は nf_get_varm_Real で行う 1 次元データの入力を, nf_get_var1_Real で 1 要素ごとに入力し, value に格納する.

iostat にはデータ入出力時のステータスが返る.

[Source]

  subroutine fake_map_get
    !
    ! ent%varid <= 0 の場合, または整数型配列 cnt に 0 が含まれる
    ! 場合に呼ばれるサブルーチン.
    !
    ! 本来は nf_get_varm_Real で行う 1 次元データの入力を, 
    ! nf_get_var1_Real で 1 要素ごとに入力し, value に格納する.
    !
    ! iostat にはデータ入出力時のステータスが返る.
    !
    use netcdf_f77, only: NF_EINDEFINE, nf_enddef, nf_redef, nf_get_var1_Real
    implicit none
    integer:: ofs(nd), here(nd)
    integer:: j, defstat
  continue
    iostat = NF_NOERR
    ofs(1:nd) = 0
    do
      j = ipos + dot_product(ofs(1:nd), imap(1:nd))
      here(1:nd) = istart(1:nd) + ofs(1:nd) * istride(1:nd)
      if (j < 1 .or. j > siz) then
        iostat = NF_EINVAL
        call DbgMessage('nf_get_var1_Real(ncid=%d, varid=%d, indx=[%*d], out-ofs=%d)', i=(/ent%fileid, ent%varid, here(1:nd), j/), n=(/nd/))
        return
      endif
      if (ent%varid == 0) then
        value(j) = j
        iostat = NF_NOERR
      else
        iostat = nf_get_var1_Real(ent%fileid, ent%varid, here(1), value(j))
        if (iostat == NF_EINDEFINE) then
          defstat = nf_enddef(ent%fileid)
          iostat = nf_get_var1_Real(ent%fileid, ent%varid, here(1), value(j))
          defstat = nf_redef(ent%fileid)
        end if
      endif
      if (iostat /= NF_NOERR) return
      ofs(1) = ofs(1) + 1
      do, j = 1, nd - 1
        if (ofs(j) < cnt(j)) exit
        ofs(j) = 0
        ofs(j + 1) = ofs(j + 1) + 1
      enddo
      if (ofs(nd) >= cnt(nd)) exit
    enddo
  end subroutine fake_map_get

[Validate]