!=begin
!=module gt4_historyauto
!
!  module gt4_history のアプリケーション. 
!  変数毎に時・空間に自由にサンプリングを設定できる.
!  長くなりそうな出力の時分割や並列化に対応. 
!
!  その他の特徴
!  * 一つの変数を切り方や時間間隔を変えて複数出すこともできる。
!  * 並列化対応は、単にプロセス固有のサフィックスをファイル名に付けられる
!    だけという単純なものである。MPIなどを使う際に、各ノードが自分の
!    受け持ちの領域をそれぞれ独立なファイルに出すことを想定している。
!  * 時間積分ループ内でデータを出力するためのコーディングを楽にする工夫
!    * 出力する変数に関する諸々の情報はモジュール内部に取っておき、
!      外部からは名前で指定する。生の gt4f90io の場合、出力初期化でできた
!      GT_HISTORY 構造体を実際の出力場所に渡さなければない。そうすると
!      時間積分のループの内と外をどうやって繋ぐかという悩みが発生するが、
!      ここではその悩みから解放される。
!    * 出力するタイミングも内部で管理されるため時間積分ループ内で毎回出力
!      命令を呼んで構わない。なお、出力のためだけに特別に計算する物理量を
!      無駄に計算しないための工夫もある(関数((<HistoryAutoWhetherPutNow>)))
!  * 地図投影座標のように、多次元の座標（補助）変数を出力したい場合も考慮さ
!    れている(((<HistoryAutoCreate>))のオプション引数 ancilcrdvars)。
!  * 任意の属性が簡単に追加できる(((<HistoryAutoCreate>))のオプション引数 attrs)。
!
!==履歴
!   2004/11/??-19  堀之内 武  作成
!   2005/02/17     堀之内 武  ドキュメントちょこっと修正
!
!==ToDo
!  * 長さ (*) の配列は本当に必要か？ 実は (:) にしたいんだけど、
!    長さチェック実装を面倒がって当座 (*) にしたんでは？？
!    もしそれだけの問題なら、長さチェックをするようにして、まずければ
!    例外が発生するようにすべき。
!  * 現在は内部データの保持は単なる linked list を使っているが、
!    登録変数が多くても高速に検索できるようにするため、2分探索出来る
!    ように変えるべき。
!  * gt4f90io に上位モジュールとして取り込んで貰う
!  * gt4f90io に整合的なドキュメンテーション （あるいは dcpam的に？）
!
!==関数の説明
!
!---subroutine HistoryAutoCreate( name, aryshape, &
!      & longname, units, slfst, sllst, slstp, &
!      & file, title, source, institution, &
!      & dims, axlongnames, axunits, axxtypes, &
!      & time_to_start, put_interval, dt,  &
!      & conventions, gt_version, &
!      & proc, newfile_interval, spcoordvars, &
!      & ancilcrdvars, attrs )
!
!      ヒストリファイル初期化情報の設定。実際のファイル初期化は
!      必要に応じて HistoryAutoPut が行う（時分割するときは適宜
!      クローズと初期化を繰り返さないとならないので、そういう
!      構造になる）。なお、一つのファイルへの出力に対して
!      このサブルーチンを２回以上呼んではならない。複数の
!      変数を一つのファイルに出したい場合は、HistoryAutoCopyCreate
!      を利用せよ。
!      時・空間に自由にサンプリングを設定できる。
!      但し、いずれも等間隔。長い時間積分によって、ファイルが
!      大きくなり過ぎることに対応するため、一定の時間間隔で
!      分割することが可能。また、並列化を念頭に各ノードを特定する
!      文字列を挿入することができる。
!
!---subroutine HistoryAutoCopyCreate( name, longname, units [, file] )
!
!      直前の HistoryAutoCreate を使って、格子及び出力の空間・時間
!      サンプリングが同じ出力を定義する。fileを省略すれば
!      同じファイルを使う。
!
!---subroutine HistoryAutoPut(name, vals, time)
!
!      変数の出力を行う。タイミングは内部で制御するので、全タイム
!      ステップで呼べば良い。なお、下記の HistoryAutoWhetherPutNow
!      を使って呼ぶタイミングを制御しても良い。
!
!---logical function HistoryAutoWhetherPutNow( name, time )
!
!      name の名を持つ出力項目に関し、現在がファイルに出力するタイミ
!      ングかどうかを返す。同名で複数の出力をする場合、どれか一つでも
!      出力するタイミングなら .true. を返す。出力のために特別に計算を
!      を要するようなケースに使うと良い。（ほとんどのステップで無駄に
!      なる計算をするのを避けられる）
!
!==公開データ型とコンストラクター
!
!---GT4_ATTRIBUTE
!     属性を名前と値の組で入れる
!---GT4_ATTRIBUTE function init_gt4_attribute(name,rval,ival,cval)
!     ATTRIBUTEのコンストラクター. 
!     名前 ＆ （実数配列 or 整数配列 or 文字列） を与える
!     
!---GT4_REAL1D
!     配列の配列をつくるための型（実数） 
!---GT4_REAL1D function init_gt4_real1d(ary)
!     REAL1Dのコンストラクター. 
!
!---GT4_NAMED_REALARY
!     名前、次元名、longname, units を持つ実数配列. 配列データは1次元で保持
!---GT4_NAMED_REALARY function init_gt4_named_realary(name,rank,dims,length,ary,longname,units)
!     GT4_NAMED_REALARYのコンストラクター. 
!=end

module gt4_historyauto

   use dc_types, only: STRING, TOKEN
   use dc_trace, only: BeginSub, EndSub, DbgMessage
   !!use gtdata_types, only: GT_VARIABLE
   use gt4_history, only: GT_HISTORY, HistoryCreate, HistoryPut, &
	& HistoryClose, HistoryAddVariable, HistorySetTime, &
	& HistoryAddAttr, HistoryInquire, Inquire, GT_HISTORY_VARINFO
   implicit none
   private
   public :: GT4_ATTRIBUTE, init_gt4_attribute
   public :: GT4_REAL1D, init_gt4_real1d
   public :: GT4_NAMED_REALARY, init_gt4_named_realary
   public :: HistoryAutoPut
   public :: HistoryAutoCreate
   public :: HistoryAutoCopyCreate
   public :: HistoryAutoWhetherPutNow
   public :: HistoryAutoSetRunInfo
   public :: HistoryAutoSetGrid

    interface HistoryAutoCreate
        module procedure HistoryAutoCreate1
        module procedure HistoryAutoCreate2
    end interface

   type GT4_ATTRIBUTE
      character(len=TOKEN)      :: name
      real,pointer              :: rval(:)
      integer,pointer           :: ival(:)
      character(len=STRING)     :: cval
   end type GT4_ATTRIBUTE

   type GT4_REAL1D    ! to make an array of 1D arrays
      real,pointer :: ary(:)
   end type GT4_REAL1D

   type GT4_NAMED_REALARY
      integer :: rank
      character(len=TOKEN) :: name
      character(len=TOKEN) :: dims(3)    ! to support up to 3D
      character(len=STRING) :: longname
      character(len=STRING) :: units
      real,pointer :: ary(:)
   end type GT4_NAMED_REALARY

   type GTHP
      type(GT_HISTORY),pointer:: hs
   end type GTHP
   

   type HIST_EACHVAR
      character(len=TOKEN):: name
      type(GTHP), pointer    :: h
      character(len=STRING)  :: longname
      character(len=STRING)  :: units
      integer              :: size
      integer              :: aryshape(3)
      integer              :: slfst(3)
      integer              :: sllst(3)
      integer              :: slstp(3)
      ! needed for HistoryCreate:
      character(len=STRING)  :: file
      character(len=TOKEN)   :: proc
      real                   :: newfile_interval    ! negative-> no renew
      character(len=STRING)  :: title
      character(len=STRING)  :: source
      character(len=STRING)  :: institution
      integer                :: sprank
      character(len=TOKEN)   :: dims(4)
      integer                :: dimsizes(4)
      character(len=STRING)  :: axlongnames(4)
      character(len=STRING)  :: axunits(4)
      character(len=TOKEN)   :: axxtypes(4)
      real                   :: time_last
      real                   :: time_to_start
      real                   :: put_interval    ! output time interval
      real                   :: dt        ! model time step (to quantize upon)
      character(len=STRING)  :: conventions
      character(len=TOKEN)   :: gt_version
      ! spastal coordinate variables:
      type(GT4_REAL1D)           :: spcoordvars(3) 
      type(GT4_NAMED_REALARY),pointer  :: ancilcrdvars(:)
      ! attributes of the data
      type(GT4_ATTRIBUTE),pointer      :: attrs(:)
   end type HIST_EACHVAR

   ! initial values of some components of HIST_EACHVAR:
!   real,parameter         :: time_last_inival = -1e35

   type HIST_LINK
      character(len=TOKEN):: name
      type(HIST_EACHVAR) :: hist
      type(HIST_LINK),pointer :: next
   end type HIST_LINK

   type(HIST_LINK),save,pointer:: HISTPOOL => NULL()  ! Fortran95 feature
   type(HIST_LINK),save,pointer:: HISTGRIDPOOL => NULL()  ! Fortran95 feature

   integer,parameter           :: max_char_len = 200
   character(len=max_char_len) :: com_proc=''
   character(len=max_char_len) :: com_title='', com_source='', com_institution=''
   character(len=max_char_len) :: com_conventions='', com_gt_version='4.2'

contains

   logical function HistoryHasVariable(history, varname) result(result)
      implicit none
      type(GT_HISTORY), intent(in):: history
      character(len = *):: varname
      type(GT_HISTORY_VARINFO), pointer :: varinfo(:)
      integer:: i
      logical :: err
      character(STRING) :: name

      result = .false.
      call Inquire(history, err = err, varinfo = varinfo)
      do i = 1, size(varinfo)
	 call Inquire(varinfo(i), name=name)
	 if (name == varname) then
	    result = .true.
	    return
	 endif
      end do
      return
   end function HistoryHasVariable

   function init_gt4_attribute(name,rval,ival,cval) result(result)
      implicit none
      type(GT4_ATTRIBUTE)                       :: result
      character(len=*),intent(in)           :: name
      real,intent(in),optional              :: rval(:)
      integer,intent(in),optional           :: ival(:)
      character(len=*),intent(in),optional  :: cval
      result%name = name
      if(present(rval)) then
	 allocate(result%rval(size(rval)))
	 result%rval = rval
	 nullify(result%ival)
      else if (present(ival)) then
	 allocate(result%ival(size(ival)))
	 result%ival = ival
	 nullify(result%rval)
      else if (present(cval)) then
	 nullify(result%rval)
	 nullify(result%ival)
	 result%cval = cval
      endif
   end function init_gt4_attribute

   function init_gt4_real1d(ary) result(result)
      implicit none
      type(GT4_REAL1D)      :: result
      real,intent(in)   :: ary(:)
      if(associated(result%ary)) deallocate(result%ary)
      allocate(result%ary(size(ary)))
      result%ary = ary
   end function init_gt4_real1d

   function init_gt4_named_realary(name,rank,dims,length,ary,longname,units) &
	& result(result)
      use dc_error, only: USR_ECHAR, StoreError
      implicit none
      type(GT4_NAMED_REALARY)            :: result
      !
      character(len=*),intent(in)    :: name
      integer,intent(in)             :: rank
      character(len = *),intent(in)  :: dims(rank)
      integer,intent(in)             :: length
      real,intent(in)                :: ary(length)
      character(len=*),intent(in)    :: longname
      character(len=*),intent(in)    :: units
      !
      character(len = *), parameter:: subname = 'init_gt4_named_realary'
      !
      call BeginSub(subname)
      if(rank>3 .or. rank<0) call StoreError(USR_ECHAR, subname, &
	      &   cause_c='rank must be <= 3 and >=1')
      result%rank = rank
      result%name = name
      result%dims(1:rank) = dims(1:rank)
      allocate(result%ary(length))     ! always new allocation
      result%ary(1:length) = ary(1:length)
      result%longname = longname
      result%units = units
      call EndSub(subname)
   end function init_gt4_named_realary

   function make_slice(vals, rank, aryshape, slfst, sllst, slstp) &
	& result(result)
      implicit none
      real,pointer,dimension(:)          :: result
      real,intent(in)                    :: vals(:)
      integer,intent(in)                 :: rank
      integer,intent(in)                 :: aryshape(*)
      integer,intent(in)                 :: slfst(*)
      integer,intent(in)                 :: sllst(*)
      integer,intent(in)                 :: slstp(*)
      !
      integer  :: i,slsize
      logical  :: slicing_needed
      real,pointer :: v1(:),v2(:,:),v3(:,:,:)
      character(len = *), parameter:: subname = 'make_slice'
      !
      call BeginSub(subname)
      !
      slicing_needed = .false.
      do i=1,rank
	 if (slfst(i)/=1)           slicing_needed = .true.
	 if (sllst(i)/=aryshape(i)) slicing_needed = .true.
	 if (slstp(i)/=1)           slicing_needed = .true.
      enddo

      slsize = 1
      do i=1,rank
	 slsize = slsize * ( (sllst(i)-slfst(i))/slstp(i) + 1 )
      enddo


      if(.not.slicing_needed) then
	 nullify(result)
      else
	 if(associated(result)) deallocate(result)
	 allocate(result(slsize))
	 select case(rank)
	 case (1)
	    if(associated(v1)) deallocate(v1)
	    allocate(v1(aryshape(1)))
	    v1 = reshape(vals,aryshape(1:1))
	    result = v1(slfst(1):sllst(1):slstp(1))
	 case (2)
	    if(associated(v2)) deallocate(v2)
	    allocate(v2(aryshape(1),aryshape(2)))
	    v2 = reshape(vals,aryshape(1:2))
	    result = reshape( &
                 & v2(slfst(1):sllst(1):slstp(1), slfst(2):sllst(2):slstp(2)),&
		 & (/slsize/) )
	 case (3)
	    if(associated(v3)) deallocate(v3)
	    allocate(v3(aryshape(1),aryshape(2),aryshape(3)))
	    v3 = reshape(vals,aryshape(1:3))
	    result = reshape( &
 	         & v3(slfst(1):sllst(1):slstp(1), slfst(2):sllst(2):slstp(2),&
		 &    slfst(3):sllst(3):slstp(3)), (/slsize/) )
	 end select
      endif
      call EndSub(subname)
   end function make_slice

   subroutine add_gt4_attribute(hst, attr)
      implicit none
      type(HIST_EACHVAR),intent(inout)   :: hst
      type(GT4_ATTRIBUTE),intent(in)         :: attr
      !
      if( hst%name /= "" ) then
	 if(associated(attr%rval)) then
	    call HistoryAddAttr(hst%name, attr%name, attr%rval, hst%h%hs)
	 else if(associated(attr%rval)) then
	    call HistoryAddAttr(hst%name, attr%name, attr%ival, hst%h%hs)
	 else
	    call HistoryAddAttr(hst%name, attr%name, trim(attr%cval), hst%h%hs)
	 endif
      endif
   end subroutine add_gt4_attribute

   subroutine add_ancilcrdvar(hst, var)
      implicit none
      type(HIST_EACHVAR),intent(inout)   :: hst
      type(GT4_NAMED_REALARY),intent(in)     :: var
      !
      integer     :: rank
      !
      rank = var%rank
      call HistoryAddVariable(var%name, var%dims(1:rank), &
	   & trim(var%longname), trim(var%units), &
	   & history=hst%h%hs)
   end subroutine add_ancilcrdvar

   subroutine put_ancilcrdvar(hst, var)
      use dc_error, only: GT_ENOMATCHDIM, StoreError
      implicit none
      type(HIST_EACHVAR),intent(inout)   :: hst
      type(GT4_NAMED_REALARY),intent(in)     :: var
      !
      integer                      :: rank,i,j
      integer,allocatable          :: idx(:)
      real,pointer                 :: subset(:)
      character(len=*), parameter  :: subname = 'put_ancilcrdvar'
      !
      call BeginSub(subname)
      rank = var%rank
      allocate(idx(rank))
      loopi: do i=1,rank
	 do j=1,hst%sprank
	    if ( var%dims(i) == hst%dims(j) ) then
	       idx(i)=j
	       cycle loopi
	    endif
	 enddo
	 call StoreError(GT_ENOMATCHDIM, subname)
      enddo loopi

      subset => make_slice(var%ary, rank, hst%aryshape(idx), &
	                 & hst%slfst(idx), hst%sllst(idx), hst%slstp(idx))
      if (associated(subset)) then
	 call HistoryPut(var%name, subset, hst%h%hs)
      else
	 call HistoryPut(var%name, var%ary, hst%h%hs)
      endif
      deallocate(idx)
      call EndSub(subname)
   end subroutine put_ancilcrdvar

   subroutine HistoryAutoPut(name, vals, time)
      implicit none
      character(len=*), intent(in)     :: name
      real                             :: vals(*)
      real                             :: time
      !
      type(HIST_EACHVAR),pointer      :: hst
      integer                         :: ith, j, rank
      character(len=STRING)           :: file_actual
      real                            :: eps=3e-7, newest
      type(GT_HISTORY),pointer        :: hist
      logical                         :: put_now
      integer                         :: arysize
      real,pointer                    :: subset(:)
      character(len = *), parameter:: subname = 'HistoryAutoPut'
      !
      call BeginSub(subname)
      ith = 1
      do while( histpl_find(HISTPOOL, name, ith, hst) )
	 put_now = whether_to_put_now( time, hst%time_last, &
	                   & hst%time_to_start, hst%put_interval, hst%dt )
	 if ( put_now ) then
	    if ( associated(hst%h%hs) ) then
	       call HistoryInquire(hst%h%hs, newest=newest)
	       if ( hst%newfile_interval > 0 .and. &
		    & time >= hst%time_to_start+hst%newfile_interval*(1.0-eps)&
		    & .and. newest < time) then
                  ! to make a new file
		  hst%time_to_start = hst%time_to_start + hst%newfile_interval
		  call HistoryClose(hst%h%hs)
		  nullify(hst%h%hs)
	       endif
	    endif
	    if (.not.associated(hst%h%hs)) then
	       if (hst%newfile_interval > 0) then
		  file_actual = merge_file_proc_time(hst%file,hst%proc, &
		       & hst%time_to_start)
	       else
		  file_actual = merge_file_proc_time(hst%file,hst%proc)
	       endif
	       rank = hst%sprank + 1
	       allocate(hist)   ! always new allocataion
	       call HistoryCreate( file_actual, trim(hst%title), &
		    & trim(hst%source), trim(hst%institution), &
		    & hst%dims(1:rank), hst%dimsizes(1:rank), &
		    & hst%axlongnames(1:rank), hst%axunits(1:rank), &
		    & hst%time_to_start, hst%put_interval, &
		    & hst%axxtypes(1:rank), hist, &
		    & trim(hst%conventions),  trim(hst%gt_version))
	       hst%h%hs => hist
	       call HistoryAddVariable(name, hst%dims(1:rank), &
                           & trim(hst%longname), trim(hst%units), &
			   & history=hst%h%hs)
	       if (associated(hst%attrs)) then
		  do j=1,size(hst%attrs)
		     call add_gt4_attribute(hst, hst%attrs(j))
		  enddo
	       endif
	       if (associated(hst%ancilcrdvars)) then
		  do j=1,size(hst%ancilcrdvars)
		     call add_ancilcrdvar(hst, hst%ancilcrdvars(j))
		  enddo
	       endif
	       do j=1,hst%sprank
		  subset => make_slice(hst%spcoordvars(j)%ary, &
		       & 1, hst%aryshape(j), &
		       & hst%slfst(j), hst%sllst(j), hst%slstp(j) )
		  if (associated(subset)) then
		     call HistoryPut(hst%dims(j), &
			  & subset, hst%h%hs)
		  else
		     call HistoryPut(hst%dims(j), &
			  & hst%spcoordvars(j)%ary, hst%h%hs)
		  endif
	       enddo
	       if (associated(hst%ancilcrdvars)) then
		  do j=1,size(hst%ancilcrdvars)
		     call put_ancilcrdvar(hst, hst%ancilcrdvars(j))
		  enddo
	       endif

	       call HistorySetTime(time, hst%h%hs)
	    else
	       if ( .not. HistoryHasVariable(hst%h%hs, name) ) then
		  call HistoryAddVariable(name, hst%dims(1:rank), &
                           & trim(hst%longname), trim(hst%units), &
			   & history=hst%h%hs)
		  !" ここで HistorySetTime すると問題が起きるので前回に従う
		  if (associated(hst%attrs)) then
		     do j=1,size(hst%attrs)
			call add_gt4_attribute(hst, hst%attrs(j))
		     enddo
		  endif
	       else
		  call HistorySetTime(time, hst%h%hs)
	       endif
	    endif

	    arysize = product(hst%aryshape(1:hst%sprank))
	    subset => make_slice(vals(1:arysize), hst%sprank, hst%aryshape, &
		           & hst%slfst, hst%sllst, hst%slstp)
	    if (associated(subset)) then
	       call HistoryPut(name, subset, hst%h%hs)
	    else
	       call HistoryPut(name, vals(1:hst%size), hst%h%hs)
	    endif
	    hst%time_last = time

	 endif
      enddo
      call EndSub(subname)
   end subroutine HistoryAutoPut

   function HistoryAutoWhetherPutNow( name, time ) result(result)
      implicit none
      logical                         :: result
      character(len=*), intent(in)    :: name
      real, intent(in)                :: time
      !
      integer                      :: ith
      type(HIST_EACHVAR),pointer   :: hst
      character(len = *), parameter:: subname = 'HistoryAutoWhetherPutNow'
      logical                      :: put_now
      !
      call BeginSub(subname)
      result = .false.
      ith = 1
      do while( histpl_find(HISTPOOL, name, ith, hst) )
	 put_now = whether_to_put_now( time, hst%time_last, &
	                   & hst%time_to_start, hst%put_interval, hst%dt )
	 if (put_now) then
	    result = .true.
	    exit
	 endif
      enddo
      call EndSub(subname)
   end function HistoryAutoWhetherPutNow

   function whether_to_put_now( time_now, time_last, &
	time_to_start, put_interval, dt ) result(result)
      implicit none
      logical           :: result
      real, intent(in)  :: time_now
      real, intent(in)  :: time_last
      real, intent(in)  :: time_to_start
      real, intent(in)  :: put_interval
      real, intent(in)  :: dt
      !
      real     :: next_put_time
      real     :: eps
      character(len = *), parameter:: subname = 'whether_to_put_now'

      call BeginSub(subname)
      eps = dt * 1e-3         ! allowable error in time in float

      if (time_now < time_to_start - eps) then
	 result = .false.
	 return
      end if

      next_put_time = time_last + put_interval  ! initially very small because
                                                ! of the init val of time_last

      if ( time_now >= (next_put_time - eps) ) then
	 result = .true.
      else
	 result = .false.
      endif
      call EndSub(subname)
   end function whether_to_put_now

   subroutine HistoryAutoCopyCreate( name, longname, units, file )
      ! use the result of the latest call of HistoryAutoCreate
      character(len=*), intent(in)          :: name
      character(len=*), intent(in)          :: longname
      character(len=*), intent(in)          :: units
      character(len=*), intent(in),optional :: file
      !
      type(HIST_EACHVAR)              :: hist
      type(HIST_EACHVAR),pointer      :: histpt
      character(len = *), parameter:: subname = 'HistoryAutoCopyCreate'
      !
      call BeginSub(subname)
      histpt => histpl_last(HISTPOOL)
      hist = histpt         ! copy the contents
      if(present_and_not_empty(file)) then
	 hist%file = file
	 allocate(hist%h)       ! always new allocation
	 nullify(hist%h%hs)
      else
	 hist%h => histpt%h
      endif
      hist%name = name
      hist%longname = longname
      hist%units = units
      call histpl_push(HISTPOOL, hist)
      call EndSub(subname)
   end subroutine HistoryAutoCopyCreate

   subroutine HistoryAutoSetRunInfo( &
	& title, source, institution, proc, &
        & conventions, gt_version  )
      use dc_error, only: USR_EINT, StoreError
      implicit none
      character(len=*), intent(in), optional     :: title, source, institution
      character(len=*), intent(in), optional     :: proc
      character(len=*), intent(in), optional     :: conventions, gt_version
      !
      character(len = *),parameter  :: subname = "HistoryAutoSetRunInfo"
   continue
      call BeginSub(subname)

      if (present(title)) com_title = title
      if (present(source)) com_source = source
      if (present(institution)) com_institution = institution
      if (present(proc)) com_proc = proc
      if (present(conventions)) com_conventions = conventions
      if (present(gt_version)) com_gt_version = gt_version

      call EndSub(subname)
   end subroutine HistoryAutoSetRunInfo

   subroutine HistoryAutoSetGrid( grid_label, aryshape, &
        & dims, axlongnames, axunits, axxtypes, &
	& coord1, coord2, coord3, ancilcrdvars )
      use dc_error, only: USR_EINT, StoreError
      implicit none
      character(len=*), intent(in)     :: grid_label
      integer, intent(in)              :: aryshape(:) ! size <= 3 (--> sprank)
      character(len=*), intent(in)     :: dims(:)           !size == sprank+1
      character(len=*), intent(in)     :: axlongnames(:)    !size == sprank+1
      character(len=*), intent(in)     :: axunits(:)        !size == sprank+1
      character(len=*), intent(in)     :: axxtypes(:)
      real, intent(in),optional        :: coord1(:) ! must present if sprank>=1
      real, intent(in),optional        :: coord2(:) ! must present if sprank>=2
      real, intent(in),optional        :: coord3(:) ! must present if sprank>=3
      type(GT4_NAMED_REALARY),intent(in),optional  :: ancilcrdvars(:)
      !
      type(HIST_EACHVAR)            :: hist
      integer                       :: sprank
      character(len = *),parameter  :: subname = "HistoryAutoSetGrid"
   continue
      call BeginSub(subname)

      sprank = min( size(aryshape), 3 )

      hist%sprank = sprank
      hist%name = grid_label
      hist%aryshape(1:sprank) = aryshape(1:sprank)
      hist%dims(1:sprank+1) = dims(1:sprank+1)
      hist%axlongnames(1:sprank+1) = axlongnames(1:sprank+1)
      hist%axunits(1:sprank+1) = axunits(1:sprank+1)
      hist%axxtypes(1:sprank+1) = axxtypes(1:sprank+1)

      if (sprank >= 1) hist%spcoordvars(1) = init_gt4_real1d( coord1 )
      if (sprank >= 2) hist%spcoordvars(2) = init_gt4_real1d( coord2 )
      if (sprank >= 3) hist%spcoordvars(3) = init_gt4_real1d( coord3 )

      if(.not. present(ancilcrdvars)) then
 	 nullify(hist%ancilcrdvars)
      else if ( size(ancilcrdvars)==0 )then
 	 nullify(hist%ancilcrdvars)
      else
	 allocate(hist%ancilcrdvars(size(ancilcrdvars))) ! always new alloc
	 hist%ancilcrdvars = ancilcrdvars
      endif

      call histpl_push(HISTGRIDPOOL, hist)
      call EndSub(subname)
   end subroutine HistoryAutoSetGrid

   ! Create a history using a grid data set with a previous call of
   ! HistoryAutoSetGrid. 
   subroutine HistoryAutoCreate2( name, longname, units, file, &
	& slfst, sllst, slstp, &
	& time_to_start, put_interval, dt, newfile_interval, &
	& attrs, &
        & grid_label, &
	& title, source, institution, conventions, gt_version, proc )
      use dc_error, only: USR_ECHAR, USR_EINT, StoreError
      implicit none
      character(len=*), intent(in)     :: name
      character(len=*), intent(in)     :: longname
      character(len=*), intent(in)     :: units
      character(len=*), intent(in)     :: file
      integer, intent(in)              :: slfst(*)     ! size == sprank
      integer, intent(in)              :: sllst(*)     ! size == sprank
      integer, intent(in)              :: slstp(*)     ! size == sprank
      real,             intent(in)     :: time_to_start, put_interval, dt
      real, intent(in)                 :: newfile_interval
      type(GT4_ATTRIBUTE),intent(in),optional  :: attrs(:)
      character(len=*), intent(in)     :: grid_label  ! <-- HistoryAutoSetGrid
      character(len=*), intent(in),optional  :: proc
      character(len=*), intent(in),optional  :: title, source, institution
      character(len=*), intent(in),optional  :: conventions, gt_version
      !
      type(HIST_EACHVAR),pointer    :: hist
      integer                       :: ith
      character(len = *),parameter  :: subname = "HistoryAutoCreate2"

      call BeginSub(subname)
      ith = 1
      if (.not.histpl_find(HISTGRIDPOOL, grid_label, ith, hist)) then
	 call StoreError(USR_ECHAR, subname, &
	      & cause_c='grid '//trim(subname)//' not found')
      endif
      
      call HistoryAutoCreate1( name, longname, units, file, &
	& slfst, sllst, slstp, &
	& time_to_start, put_interval, dt, newfile_interval, &
	& attrs, &
        & hist%aryshape(1:hist%sprank), &
	& hist%dims(1:hist%sprank+1), hist%axlongnames(1:hist%sprank+1), &
	& hist%axunits(1:hist%sprank+1), hist%axxtypes(1:hist%sprank+1), &
	& hist%spcoordvars, hist%ancilcrdvars, &
	& title, source, institution, conventions, gt_version, proc )

      call EndSub(subname)
   end subroutine HistoryAutoCreate2


   subroutine HistoryAutoCreate1( name, longname, units, file, &
	& slfst, sllst, slstp, &
	& time_to_start, put_interval, dt, newfile_interval, &
	& attrs, &
        & aryshape, dims, axlongnames, axunits, axxtypes, &
	& spcoordvars, ancilcrdvars, &
	& title, source, institution, conventions, gt_version, proc )
      use dc_error, only: USR_EINT, StoreError
      implicit none
      character(len=*), intent(in)     :: name
      character(len=*), intent(in)     :: longname
      character(len=*), intent(in)     :: units
      character(len=*), intent(in)     :: file
      integer, intent(in)              :: slfst(*)     ! size == sprank
      integer, intent(in)              :: sllst(*)     ! size == sprank
      integer, intent(in)              :: slstp(*)     ! size == sprank
      real,             intent(in)     :: time_to_start, put_interval, dt
      real, intent(in)                 :: newfile_interval
      type(GT4_ATTRIBUTE),intent(in),optional  :: attrs(:)
      integer, intent(in)              :: aryshape(:) ! size <= 3 (--> sprank)
      character(len=*), intent(in)     :: dims(*)           !size == sprank+1
      character(len=*), intent(in)     :: axlongnames(*)    !size == sprank+1
      character(len=*), intent(in)     :: axunits(*)        !size == sprank+1
      character(len=*), intent(in)     :: axxtypes(*)
      type(GT4_REAL1D), intent(in)         :: spcoordvars(*)    ! size == sprank
      type(GT4_NAMED_REALARY),intent(in),optional  :: ancilcrdvars(:)
      character(len=*), intent(in),optional  :: proc
      character(len=*), intent(in),optional  :: title, source, institution
      character(len=*), intent(in),optional  :: conventions, gt_version
      !
      type(HIST_EACHVAR)            :: hist
      integer                       :: sprank,i
      character(len = *),parameter  :: subname = "HistoryAutoCreate1"

      call BeginSub(subname)

      !< initialize hist except hist%h -- actual creation is deferred >
 
      hist%time_last = -1e35 ! time_last_inival
      allocate(hist%h)       ! always new allocation
      nullify(hist%h%hs)

      hist%name = name
      hist%longname = longname
      hist%units = units

      sprank = min( size(aryshape), 3 )
      hist%sprank = sprank

      hist%size = 1
      do i=1,sprank
	 hist%aryshape(i) = aryshape(i)
	 if(slfst(i) > 0) then
	    hist%slfst(i) = slfst(i)
	 else
	    hist%slfst(i) = slfst(i) + aryshape(i)
	 endif
	 if(sllst(i) > 0) then
	    hist%sllst(i) = sllst(i)
	 else
	    hist%sllst(i) = sllst(i) + aryshape(i)
	 endif
	 if(slstp(i) > 0) then
	    hist%slstp(i) = slstp(i)
	 else
	    hist%slstp(i) = 1
	 endif
	 hist%dimsizes(i) = (hist%sllst(i)-hist%slfst(i))/hist%slstp(i) + 1
	 if (hist%slfst(i)<=0 .or. hist%slfst(i)>aryshape(i)) &
	      & call StoreError(USR_EINT, subname, &
	      &   cause_c='str not within the index range for dim:', cause_i=i)
	 if (hist%sllst(i)<=0 .or. hist%sllst(i)>aryshape(i)) &
	      & call StoreError(USR_EINT, subname, &
	      &   cause_c='end not within the index range for dim:', cause_i=i)
	 if (hist%slstp(i)<=0) &
	      & call StoreError(USR_EINT, subname, &
	      &   cause_c='step not positive for dim:', cause_i=i)
	 if (hist%dimsizes(i)<=0) &
	      & call StoreError(USR_EINT, subname, &
	      &   cause_c='negative dimsize for dim:', cause_i=i)
	 hist%size = hist%size * hist%dimsizes(i)
      enddo
      hist%dimsizes(sprank+1) = 0   ! unlimited dimension

      hist%file = file
      hist%newfile_interval = newfile_interval
      hist%dims(1:sprank+1) = dims(1:sprank+1)
      hist%axlongnames(1:sprank+1) = axlongnames(1:sprank+1)
      hist%axunits(1:sprank+1) = axunits(1:sprank+1)
      hist%time_to_start = time_to_start
      hist%put_interval = put_interval
      hist%dt = dt
      hist%axxtypes(1:sprank+1) = axxtypes(1:sprank+1)
      if(present(title)) then
	 hist%title = title
      else
	 hist%title = com_title
      endif
      if(present(source)) then
	 hist%source = source
      else
	 hist%source = com_source
      endif
      if(present(institution)) then
	 hist%institution = institution
      else
	 hist%institution = com_institution
      endif
      if(present(conventions)) then
	 hist%conventions = conventions
      else
	 hist%conventions = com_conventions
      endif
      if(present(gt_version)) then
	 hist%gt_version = gt_version
      else
	 hist%gt_version = com_gt_version
      endif
      if(present(proc)) then
	 hist%proc = proc
      else
	 hist%proc = com_proc
      endif

      hist%spcoordvars(1:sprank) = spcoordvars(1:sprank)
      if(.not. present(ancilcrdvars)) then
 	 nullify(hist%ancilcrdvars)
      else if ( size(ancilcrdvars)==0 )then
 	 nullify(hist%ancilcrdvars)
      else
	 allocate(hist%ancilcrdvars(size(ancilcrdvars))) ! always new alloc
	 hist%ancilcrdvars = ancilcrdvars
      endif

      if(.not. present(attrs)) then
 	 nullify(hist%attrs)
      else if ( size(attrs)==0 )then
 	 nullify(hist%attrs)
      else
	 allocate(hist%attrs(size(attrs))) ! always new alloc
	 hist%attrs = attrs
      endif

      call histpl_push(HISTPOOL, hist)

      call EndSub(subname)
   end subroutine HistoryAutoCreate1

   function histpl_find(histpl, name,ith,hist) result(result)
      implicit none
      logical                        :: result   ! .true. if found
      type(HIST_LINK),pointer        :: histpl ! intent(in)
      character(len=*), intent(in)   :: name
      integer, intent(inout)         :: ith   ! ith+=1 when return (to iterate)
      type(HIST_EACHVAR),pointer     :: hist
      !
      type(HIST_LINK),pointer,save :: hp
      integer,save :: cnt=1
      character(len=TOKEN),save   :: name_save = ''
      !
      if(name/=name_save .or. ith<cnt) then
	 cnt = 1
	 hp => histpl
      endif
      do while (associated(hp))
	 !!print *,trim(name),ith,cnt,trim(hp%name)
	 if (hp%name == name) then
	    if(cnt==ith) then
	       !!print *,' ...found'
	       hist => hp%hist
	       result = .true.      ! found
	       name_save = name     !    save the name found
	       ith = ith+1          !    stepped forward for the next search
	       cnt = cnt + 1        !    stepped forward for the next search
	       hp => hp%next        !    stepped forward for the next search
	       return
	    endif
	    cnt = cnt + 1
	 endif
	 hp => hp%next
      end do
      !!print *,'    ...not found'
      result = .false.     ! not found
      name_save = ''       !     initialize
      cnt = 1              !     initialize
   end function histpl_find

   function histpl_to_the_end(histpl) result(result)
      type(HIST_LINK),pointer            :: result
      type(HIST_LINK),pointer        :: histpl ! intent(in)
      result => histpl
      do while (associated(result))
	 if (associated(result%next)) then
	    result => result%next
	 else
	    exit
	 endif
      end do
   end function histpl_to_the_end

   function histpl_last(histpl) result(result)
      implicit none
      type(HIST_EACHVAR),pointer         :: result
      type(HIST_LINK),pointer        :: histpl ! intent(in)
      !
      type(HIST_LINK),pointer            :: hp
      hp => histpl_to_the_end(histpl)
      result => hp%hist
   end function histpl_last

   subroutine histpl_push(histpl, hist)
      implicit none
      type(HIST_LINK),pointer        :: histpl ! intent(in)
      type(HIST_EACHVAR),intent(in)        :: hist
      !
      type(HIST_LINK),pointer               :: hp, nxt
      hp => histpl_to_the_end(histpl)
      if ( .not. associated(hp) ) then
         ! must be the first time
	 allocate(hp)   ! always new allocation
	 histpl => hp
      else
	 allocate(nxt)
	 hp%next => nxt
	 hp => nxt
      endif
      hp%hist = hist
      hp%name = hist%name
   end subroutine histpl_push

   function merge_file_proc_time(file,proc,time) result(result)
      implicit none
      character(len=STRING)         :: result
      character(len=*), intent(in)  :: file
      character(len=*), intent(in)  :: proc
      real, intent(in), optional    :: time
      !
      integer                       :: idx
      character(len=TOKEN)          :: ctime
      character(len=10)             :: fmt
      !
      if(.not.present(time)) then
	 ctime = ""
      else
	 if (aint(time) == time) then
	    fmt = "(I)"
	    write(ctime,fmt=fmt) nint(time)
	    idx = index(ctime, '.')
	    if (idx>0) ctime = ctime(1:idx-1)
	 else
	    write(ctime,*) time
	 endif
	 ctime = '_t'//adjustl(ctime)//'-'
      endif
      !
      if (proc == "") then
	 result = file
      else
	 idx = index(file, '.nc', .true.)  ! tru -> search the right-most match
	 if (idx == 0) then
	    result = trim(file) // trim(adjustl(proc))
	 else if (idx /= 1) then
	    result = file(1:idx-1) // trim(adjustl(proc)) // '.nc'
	 else
	    result = trim(adjustl(proc)) // '.nc'
	 endif
      endif
      if (ctime == "") then
	 ! do nothing
      else
	 idx = index(result, '.nc', .true.)  !tru-> search the right-most match
	 if (idx == 0) then
	    result = trim(result) // trim(adjustl(ctime))
	 else if (idx /= 1) then
	    result = result(1:idx-1) // trim(adjustl(ctime)) // '.nc'
	 else
	    result = trim(adjustl(ctime)) // '.nc'
	 endif
      endif
   end function merge_file_proc_time

   function present_and_not_empty(arg) result(result)
      logical :: result
      character(len=*),intent(in),optional :: arg
      if(present(arg)) then
         if(arg/="") then
            result=.true.
         else
            result=.false.
         endif
      else
         result=.false.
      endif
   end function present_and_not_empty

end module gt4_historyauto
