SCALE-RM
Data Types | Functions/Subroutines
scale_history Module Reference

module HISTORY More...

Functions/Subroutines

subroutine, public hist_setup
 Setup. More...
 
subroutine, public hist_switch (switch)
 set switch More...
 
subroutine, public hist_reg (itemid, zinterp, item, desc, unit, ndim, xdim, ydim, zdim)
 Register/Append variable to history file. More...
 
subroutine, public hist_query (itemid, answer)
 Check time to putting data. More...
 
subroutine hist_put_0d (itemid, var)
 Put 1D data to history buffer. More...
 
subroutine hist_put_1d (itemid, var)
 Put 1D data to history buffer. More...
 
subroutine hist_put_2d (itemid, var, nohalo)
 Put 2D data to history buffer. More...
 
subroutine hist_put_3d (itemid, var, zinterp, xdim, ydim, zdim, nohalo)
 Put 3D data to history buffer. More...
 
subroutine hist_in_0d (var, item, desc, unit)
 Wrapper routine of HIST_reg+HIST_put 0D. More...
 
subroutine hist_get_2d (var, basename, varname, step, allow_missing)
 Get 2D data from file. More...
 
subroutine hist_get_3d (var, basename, varname, step, allow_missing)
 Get 3D data from file. More...
 
subroutine, public hist_write
 Flush history buffer to file. More...
 

Detailed Description

module HISTORY

Description
History output module
Author
Team SCALE
History
  • 2011-12-05 (H.Yashiro) [new]
  • 2012-03-23 (H.Yashiro) [mod] Explicit index parameter inclusion
  • 2012-06-11 (S.Nishizawa) [mod] use gtool_history

Function/Subroutine Documentation

◆ hist_setup()

subroutine, public scale_history::hist_setup ( )

Setup.

Definition at line 87 of file scale_history.F90.

References scale_stdio::h_institute, scale_stdio::h_source, gtool_history::historyinit(), scale_grid_index::ie, scale_grid_index::ieb, scale_grid_index::imax, scale_grid_index::imaxb, scale_stdio::io_fid_conf, scale_stdio::io_fid_log, scale_stdio::io_l, scale_stdio::io_lnml, scale_grid_index::is, scale_grid_index::isb, scale_grid_index::je, scale_grid_index::jeb, scale_grid_index::jmax, scale_grid_index::jmaxb, scale_grid_index::js, scale_grid_index::jsb, scale_grid_index::kmax, scale_rm_process::prc_2drank, scale_process::prc_masterrank, scale_process::prc_mpistop(), scale_process::prc_myrank, scale_prof::prof_rapend(), scale_prof::prof_rapstart(), scale_time::time_dtsec, scale_time::time_nowdate, scale_time::time_nowms, scale_time::time_offset_year, and scale_time::time_startdaysec.

Referenced by mod_rm_driver::scalerm(), and mod_rm_prep::scalerm_prep().

87  use gtool_history, only: &
89  use scale_process, only: &
90  prc_mpistop, &
93  use scale_rm_process, only: &
95  use scale_time, only: &
96  time_nowdate, &
97  time_nowms, &
98  time_dtsec, &
101  implicit none
102 
103  character(len=H_MID) :: history_h_title = 'SCALE-RM HISTORY OUTPUT'
104  character(len=H_SHORT) :: history_t_units = 'seconds'
105  character(len=H_MID) :: history_t_since = ''
106 
107  logical :: hist_bnd = .false.
108 
109  namelist / param_hist / &
110  hist_bnd
111 
112  real(DP) :: start_daysec
113  integer :: rankidx(2)
114  integer :: ierr
115  !---------------------------------------------------------------------------
116 
117  if( io_l ) write(io_fid_log,*)
118  if( io_l ) write(io_fid_log,*) '++++++ Module[HISTORY] / Categ[ATMOS-RM IO] / Origin[SCALElib]'
119 
120  !--- read namelist
121  rewind(io_fid_conf)
122  read(io_fid_conf,nml=param_hist,iostat=ierr)
123  if( ierr < 0 ) then !--- missing
124  if( io_l ) write(io_fid_log,*) '*** Not found namelist. Default used.'
125  elseif( ierr > 0 ) then !--- fatal error
126  write(*,*) 'xxx Not appropriate names in namelist PARAM_HIST. Check!'
127  call prc_mpistop
128  endif
129  if( io_lnml ) write(io_fid_log,nml=param_hist)
130 
131  call prof_rapstart('FILE_O_NetCDF', 2)
132 
133  rankidx(1) = prc_2drank(prc_myrank, 1)
134  rankidx(2) = prc_2drank(prc_myrank, 2)
135 
136  start_daysec = time_startdaysec
137  if ( time_nowdate(1) > 0 ) then
138  write(history_t_since, '(I4.4,5(A1,I2.2))') &
139  time_nowdate(1),'-',time_nowdate(2),'-',time_nowdate(3),' ', &
140  time_nowdate(4),':',time_nowdate(5),':',time_nowdate(6)
141  start_daysec = time_nowms
142  endif
143 
144  if ( hist_bnd ) then
145  im = imaxb
146  jm = jmaxb
147  ims = isb
148  ime = ieb
149  jms = jsb
150  jme = jeb
151  else
152  im = imax
153  jm = jmax
154  ims = is
155  ime = ie
156  jms = js
157  jme = je
158  endif
159 
160  call historyinit( history_h_title, & ! [IN]
161  h_source, & ! [IN]
162  h_institute, & ! [IN]
163  im*jm*kmax, & ! [IN]
164  prc_masterrank, & ! [IN]
165  prc_myrank, & ! [IN]
166  rankidx, & ! [IN]
167  start_daysec, & ! [IN]
168  time_dtsec, & ! [IN]
169  time_units = history_t_units, & ! [IN]
170  time_since = history_t_since, & ! [IN]
171  namelist_fid = io_fid_conf ) ! [IN]
172 
173  call hist_put_axes
174 
175  call prof_rapend ('FILE_O_NetCDF', 2)
176 
177  enabled = .true.
178 
179  return
integer, public imax
of computational cells: x
integer, public is
start point of inner domain: x, local
integer, public je
end point of inner domain: y, local
subroutine, public prc_mpistop
Abort MPI.
integer, public jeb
real(dp), public time_nowms
subsecond part of current time [millisec]
Definition: scale_time.F90:66
real(dp), public time_startdaysec
second of start time [sec]
Definition: scale_time.F90:74
integer, public imaxb
integer, public jmaxb
module Gtool_History
integer, public ieb
real(dp), public time_dtsec
time interval of model [sec]
Definition: scale_time.F90:36
integer, public time_offset_year
time offset [year]
Definition: scale_time.F90:73
integer, public kmax
of computational cells: z
subroutine, public historyinit(title, source, institution, array_size, master, myrank, rankidx, time_start, time_interval, time_units, time_since, default_basename, default_tinterval, default_tunit, default_taverage, default_zinterp, default_datatype, namelist_filename, namelist_fid)
integer, public js
start point of inner domain: y, local
module TIME
Definition: scale_time.F90:15
module PROCESS
integer, parameter, public prc_masterrank
master process in each communicator
integer, public prc_myrank
process num in local communicator
module RM PROCESS
integer, public ie
end point of inner domain: x, local
integer, dimension(:,:), allocatable, public prc_2drank
node index in 2D topology
integer, public isb
integer, dimension(6), public time_nowdate
current time [YYYY MM DD HH MM SS]
Definition: scale_time.F90:65
integer, public io_fid_conf
Config file ID.
Definition: scale_stdio.F90:55
integer, public jsb
integer, public jmax
of computational cells: y
Here is the call graph for this function:
Here is the caller graph for this function:

◆ hist_switch()

subroutine, public scale_history::hist_switch ( logical, intent(in)  switch)

◆ hist_reg()

subroutine, public scale_history::hist_reg ( integer, intent(out)  itemid,
logical, intent(out)  zinterp,
character(len=*), intent(in)  item,
character(len=*), intent(in)  desc,
character(len=*), intent(in)  unit,
integer, intent(in)  ndim,
character(len=*), intent(in), optional  xdim,
character(len=*), intent(in), optional  ydim,
character(len=*), intent(in), optional  zdim 
)

Register/Append variable to history file.

Parameters
[out]itemidindex number of the item
[out]zinterpz* -> z flag of the item
[in]itemname of the item
[in]descdescription of the item
[in]unitunit of the item
[in]ndimdimension of the item

Definition at line 494 of file scale_history.F90.

References gtool_history::historyaddvariable(), scale_prof::prof_rapend(), scale_prof::prof_rapstart(), and scale_time::time_nowstep.

Referenced by scale_atmos_phy_mp_tomita08::atmos_phy_mp_tomita08_setup(), mod_atmos_vars::atmos_vars_setup(), and hist_in_0d().

494  use gtool_history, only: &
496  use scale_time, only: &
497  nowstep => time_nowstep
498  implicit none
499 
500  integer, intent(out) :: itemid
501  logical, intent(out) :: zinterp
502  character(len=*), intent(in) :: item
503  character(len=*), intent(in) :: desc
504  character(len=*), intent(in) :: unit
505  integer, intent(in) :: ndim
506 
507  character(len=*), intent(in), optional :: xdim
508  character(len=*), intent(in), optional :: ydim
509  character(len=*), intent(in), optional :: zdim
510 
511  logical :: xh
512  logical :: yh
513 
514  logical :: existed
515 
516  character(len=H_SHORT) :: dims(3)
517  !---------------------------------------------------------------------------
518 
519  call prof_rapstart('FILE_O_NetCDF', 2)
520 
521  if ( ndim == 1 ) then
522 
523  dims(1) = "z"
524  if ( present(zdim) ) then
525  if ( zdim=='half' ) then
526  dims(1) = "zh"
527  endif
528  endif
529 
530  else
531 
532  xh = .false.
533  if ( present(xdim) ) then
534  if ( xdim=='half' ) xh = .true.
535  endif
536  yh = .false.
537  if ( present(ydim) ) then
538  if ( ydim=='half' ) yh = .true.
539  endif
540 
541  if ( xh .AND. yh ) then
542  dims(1) = 'lon_uv'
543  dims(2) = 'lat_uv'
544  dims(3) = 'height_uvz'
545  else if ( xh ) then
546  dims(1) = 'lon_uy'
547  dims(2) = 'lat_uy'
548  dims(3) = 'height_uyz'
549  else if ( yh ) then
550  dims(1) = 'lon_xv'
551  dims(2) = 'lat_xv'
552  dims(3) = 'height_xvz'
553  else
554  dims(1) = 'lon'
555  dims(2) = 'lat'
556  dims(3) = 'height'
557  endif
558 
559  if ( present(zdim) ) then
560  if ( zdim=='land' ) then
561  dims(3) = 'lz'
562  else if ( zdim=='landhalf' ) then
563  dims(3) = 'lzh'
564  else if ( zdim=='urban' ) then
565  dims(3) = 'uz'
566  else if ( zdim=='urbanhalf' ) then
567  dims(3) = 'uzh'
568  else if ( zdim=='half' ) then
569  if ( xh .AND. yh ) then
570  dims(3) = 'height_uvw'
571  else if ( xh ) then
572  dims(3) = 'height_uyw'
573  else if ( yh ) then
574  dims(3) = 'height_xvw'
575  else
576  dims(3) = 'height_xyw'
577  endif
578  endif
579  endif
580 
581  endif
582 
583  call historyaddvariable( item, & ! [IN]
584  dims(1:ndim), & ! [IN]
585  desc, & ! [IN]
586  unit, & ! [IN]
587  nowstep, & ! [IN]
588  itemid, & ! [OUT]
589  zinterp, & ! [OUT]
590  existed ) ! [OUT]
591 
592  call prof_rapend ('FILE_O_NetCDF', 2)
593 
594  return
integer, public time_nowstep
current step [number]
Definition: scale_time.F90:70
module Gtool_History
subroutine, public historyaddvariable(varname, dims, desc, units, now_step, id, zinterp, existed, options)
module TIME
Definition: scale_time.F90:15
Here is the call graph for this function:
Here is the caller graph for this function:

◆ hist_query()

subroutine, public scale_history::hist_query ( integer, intent(in)  itemid,
logical, intent(out)  answer 
)

Check time to putting data.

Parameters
[in]itemidindex number of the item
[out]answeris it time to store?

Definition at line 602 of file scale_history.F90.

References gtool_history::historyquery(), scale_prof::prof_rapend(), scale_prof::prof_rapstart(), and scale_time::time_nowstep.

Referenced by scale_atmos_phy_mp_tomita08::atmos_phy_mp_tomita08(), and hist_in_0d().

602  use gtool_history, only: &
604  use scale_time, only: &
606  implicit none
607 
608  integer, intent(in) :: itemid
609  logical, intent(out) :: answer
610 
611  !---------------------------------------------------------------------------
612 
613  answer = .false.
614 
615  if ( .NOT. enabled ) return
616 
617  if ( itemid < 0 ) return
618 
619  call prof_rapstart('FILE_O_NetCDF', 2)
620 
621  call historyquery(itemid, time_nowstep, answer)
622 
623  call prof_rapend ('FILE_O_NetCDF', 2)
624 
625  return
integer, public time_nowstep
current step [number]
Definition: scale_time.F90:70
module Gtool_History
module TIME
Definition: scale_time.F90:15
subroutine, public historyquery(itemid, step_next, answer)
Here is the call graph for this function:
Here is the caller graph for this function:

◆ hist_put_0d()

subroutine scale_history::hist_put_0d ( integer, intent(in)  itemid,
real(rp), intent(in)  var 
)

Put 1D data to history buffer.

Parameters
[in]itemidindex number of the item
[in]varvalue

Definition at line 633 of file scale_history.F90.

References scale_prof::prof_rapend(), scale_prof::prof_rapstart(), and scale_time::time_nowstep.

633  use gtool_history, only: &
634  historyput
635  use scale_time, only: &
637  implicit none
638 
639  integer, intent(in) :: itemid
640  real(RP), intent(in) :: var
641 
642  !---------------------------------------------------------------------------
643 
644  if ( .NOT. enabled ) return
645 
646  if ( itemid < 0 ) return
647 
648  call prof_rapstart('FILE_O_NetCDF', 2)
649 
650  call historyput(itemid, time_nowstep, var)
651 
652  call prof_rapend ('FILE_O_NetCDF', 2)
653 
654  return
integer, public time_nowstep
current step [number]
Definition: scale_time.F90:70
module Gtool_History
module TIME
Definition: scale_time.F90:15
Here is the call graph for this function:

◆ hist_put_1d()

subroutine scale_history::hist_put_1d ( integer, intent(in)  itemid,
real(rp), dimension(:), intent(in)  var 
)

Put 1D data to history buffer.

Parameters
[in]itemidindex number of the item
[in]varvalue

Definition at line 662 of file scale_history.F90.

References scale_grid_index::kmax, scale_grid_index::ks, scale_prof::prof_rapend(), scale_prof::prof_rapstart(), and scale_time::time_nowstep.

662  use gtool_history, only: &
663  historyput
664  use scale_time, only: &
666  implicit none
667 
668  integer, intent(in) :: itemid
669  real(RP), intent(in) :: var(:)
670 
671  real(RP) :: var2(kmax)
672  integer :: k
673  !---------------------------------------------------------------------------
674 
675  if ( .NOT. enabled ) return
676 
677  if ( itemid < 0 ) return
678 
679  call prof_rapstart('FILE_O_NetCDF', 2)
680 
681  do k = 1, kmax
682  var2(k) = var(ks+k-1)
683  enddo
684 
685  call historyput(itemid, time_nowstep, var2)
686 
687  call prof_rapend ('FILE_O_NetCDF', 2)
688 
689  return
integer, public time_nowstep
current step [number]
Definition: scale_time.F90:70
module Gtool_History
integer, public kmax
of computational cells: z
module TIME
Definition: scale_time.F90:15
integer, public ks
start point of inner domain: z, local
Here is the call graph for this function:

◆ hist_put_2d()

subroutine scale_history::hist_put_2d ( integer, intent(in)  itemid,
real(rp), dimension(:,:), intent(in)  var,
logical, intent(in), optional  nohalo 
)

Put 2D data to history buffer.

Parameters
[in]itemidindex number of the item
[in]varvalue

Definition at line 698 of file scale_history.F90.

References scale_grid_index::ie, scale_grid_index::is, scale_grid_index::je, scale_grid_index::js, scale_prof::prof_rapend(), scale_prof::prof_rapstart(), gtool_file::rmiss, and scale_time::time_nowstep.

698  use gtool_file, only: &
699  rmiss
700  use gtool_history, only: &
701  historyput
702  use scale_time, only: &
704  implicit none
705 
706  integer, intent(in) :: itemid
707  real(RP), intent(in) :: var(:,:)
708 
709  logical, intent(in), optional :: nohalo
710 
711  real(RP) :: var2(im*jm)
712  integer :: i, j
713  logical :: nohalo_
714  !---------------------------------------------------------------------------
715 
716  if ( .NOT. enabled ) return
717 
718  if ( itemid < 0 ) return
719 
720  call prof_rapstart('FILE_O_NetCDF', 2)
721 
722  do j = 1, jm
723  do i = 1, im
724  var2(i + (j-1)*im) = var(ims+i-1,jms+j-1)
725  enddo
726  enddo
727 
728  nohalo_ = .false.
729  if ( present(nohalo) ) nohalo_ = nohalo
730 
731  if ( nohalo_ ) then
732  ! W halo
733  do j = 1, jm
734  do i = 1, is-ims
735  var2(i + (j-1)*im) = rmiss
736  enddo
737  enddo
738  ! E halo
739  do j = 1, jm
740  do i = ie-ims+2, ime-ims+1
741  var2(i + (j-1)*im) = rmiss
742  enddo
743  enddo
744  ! S halo
745  do j = 1, js-jms
746  do i = 1, im
747  var2(i + (j-1)*im) = rmiss
748  enddo
749  enddo
750  ! N halo
751  do j = je-jms+2, jme-jms+1
752  do i = 1, im
753  var2(i + (j-1)*im) = rmiss
754  enddo
755  enddo
756  endif
757 
758  call historyput(itemid, time_nowstep, var2)
759 
760  call prof_rapend ('FILE_O_NetCDF', 2)
761 
762  return
integer, public time_nowstep
current step [number]
Definition: scale_time.F90:70
integer, public is
start point of inner domain: x, local
module GTOOL_FILE
Definition: gtool_file.f90:17
integer, public je
end point of inner domain: y, local
real(dp), parameter, public rmiss
Definition: gtool_file.f90:141
module Gtool_History
integer, public js
start point of inner domain: y, local
module TIME
Definition: scale_time.F90:15
integer, public ie
end point of inner domain: x, local
Here is the call graph for this function:

◆ hist_put_3d()

subroutine scale_history::hist_put_3d ( integer, intent(in)  itemid,
real(rp), dimension(:,:,:), intent(in)  var,
logical, intent(in)  zinterp,
character(len=*), intent(in), optional  xdim,
character(len=*), intent(in), optional  ydim,
character(len=*), intent(in), optional  zdim,
logical, intent(in), optional  nohalo 
)

Put 3D data to history buffer.

Parameters
[in]itemidindex number of the item
[in]varvalue
[in]zinterpvertical interpolation?

Definition at line 775 of file scale_history.F90.

References scale_grid_index::ia, scale_grid_index::ie, scale_interpolation::interp_available, scale_interpolation::interp_vertical_xi2z(), scale_grid_index::is, scale_grid_index::ja, scale_grid_index::je, scale_grid_index::js, scale_grid_index::ka, scale_grid_index::kmax, scale_grid_index::ks, scale_land_grid_index::lke, scale_land_grid_index::lkmax, scale_land_grid_index::lks, scale_prof::prof_rapend(), scale_prof::prof_rapstart(), gtool_file::rmiss, scale_time::time_nowstep, scale_urban_grid_index::uke, scale_urban_grid_index::ukmax, and scale_urban_grid_index::uks.

775  use gtool_file, only: &
776  rmiss
777  use gtool_history, only: &
778  historyput
779  use scale_time, only: &
781  use scale_interpolation, only: &
784  implicit none
785 
786  integer, intent(in) :: itemid
787  real(RP), intent(in) :: var(:,:,:)
788  logical, intent(in) :: zinterp
789 
790  character(len=*), intent(in), optional :: xdim
791  character(len=*), intent(in), optional :: ydim
792  character(len=*), intent(in), optional :: zdim
793  logical, intent(in), optional :: nohalo
794 
795  intrinsic shape
796  integer :: s(3)
797 
798  character(len=H_SHORT) :: xd, yd, zd
799 
800  real(RP), allocatable :: var_z(:,:,:)
801  real(RP), allocatable :: var2 (:)
802 
803  integer :: i, j, k
804  integer :: isize, jsize, ksize
805  integer :: iall, jall, kall
806  integer :: istart, jstart, kstart
807 
808  logical :: nohalo_
809  !---------------------------------------------------------------------------
810 
811  if ( .NOT. enabled ) return
812 
813  if ( itemid < 0 ) return
814 
815  call prof_rapstart('FILE_O_NetCDF', 2)
816 
817  xd = ''
818  yd = ''
819  zd = ''
820  if( present(xdim) ) xd = xdim
821  if( present(ydim) ) yd = ydim
822  if( present(zdim) ) zd = zdim
823 
824  nohalo_ = .false.
825  if ( present(nohalo) ) nohalo_ = nohalo
826 
827 
828  ! select dimension
829  select case ( xd )
830  case default
831  isize = im
832  iall = ia
833  istart = ims
834  end select
835  select case ( yd )
836  case default
837  jsize = jm
838  jall = ja
839  jstart = jms
840  end select
841  select case ( zd )
842  case ('land')
843  ksize = lkmax
844  kall = lke
845  kstart = lks
846  case ('urban')
847  ksize = ukmax
848  kall = uke
849  kstart = uks
850  case default
851  ksize = kmax
852  kall = ka
853  kstart = ks
854  end select
855 
856  allocate( var_z(kall,iall,jall) )
857  allocate( var2(ksize*isize*jsize) )
858 
859  s = shape(var)
860  if ( s(1) == 1 ) then
861 
862  do j = 1, jsize
863  do i = 1, isize
864  var2(i + (j-1)*isize) = var(1,istart+i-1,jstart+j-1)
865  enddo
866  enddo
867 
868  if ( nohalo_ ) then
869  ! W halo
870  do j = 1, jsize
871  do i = 1, is-ims
872  var2(i + (j-1)*isize) = rmiss
873  enddo
874  enddo
875  ! E halo
876  do j = 1, jsize
877  do i = ie-ims+2, ime-ims+1
878  var2(i + (j-1)*isize) = rmiss
879  enddo
880  enddo
881  ! S halo
882  do j = 1, js-jms
883  do i = 1, isize
884  var2(i + (j-1)*isize) = rmiss
885  enddo
886  enddo
887  ! N halo
888  do j = je-jms+2, jme-jms+1
889  do i = 1, isize
890  var2(i + (j-1)*isize) = rmiss
891  enddo
892  enddo
893  endif
894 
895  call historyput(itemid, time_nowstep, var2(1:isize*jsize))
896 
897  else
898  if ( ksize == kmax &
899  .AND. zinterp &
900  .AND. interp_available ) then
901  call prof_rapstart('FILE_O_interp', 2)
902  call interp_vertical_xi2z( var(:,:,:), & ! [IN]
903  var_z(:,:,:) ) ! [OUT]
904  call prof_rapend ('FILE_O_interp', 2)
905 
906  do k = 1, ksize
907  do j = 1, jsize
908  do i = 1, isize
909  var2(i + (j-1)*isize + (k-1)*jsize*isize) = var_z(kstart+k-1,istart+i-1,jstart+j-1)
910  enddo
911  enddo
912  enddo
913  else
914  do k = 1, ksize
915  do j = 1, jsize
916  do i = 1, isize
917  var2(i + (j-1)*isize + (k-1)*jsize*isize) = var(kstart+k-1,istart+i-1,jstart+j-1)
918  enddo
919  enddo
920  enddo
921  endif
922 
923  if ( nohalo_ ) then
924  ! W halo
925  do k = 1, ksize
926  do j = 1, jsize
927  do i = 1, is-ims
928  var2(i + (j-1)*isize + (k-1)*jsize*isize) = rmiss
929  enddo
930  enddo
931  enddo
932  ! E halo
933  do k = 1, ksize
934  do j = 1, jsize
935  do i = ie-ims+2, ime-ims+1
936  var2(i + (j-1)*isize + (k-1)*jsize*isize) = rmiss
937  enddo
938  enddo
939  enddo
940  ! S halo
941  do k = 1, ksize
942  do j = 1, js-jms
943  do i = 1, isize
944  var2(i + (j-1)*isize + (k-1)*jsize*isize) = rmiss
945  enddo
946  enddo
947  enddo
948  ! N halo
949  do k = 1, ksize
950  do j = je-jms+2, jme-jms+1
951  do i = 1, isize
952  var2(i + (j-1)*isize + (k-1)*jsize*isize) = rmiss
953  enddo
954  enddo
955  enddo
956  endif
957 
958  call historyput(itemid, time_nowstep, var2)
959 
960  endif
961 
962  deallocate( var_z )
963  deallocate( var2 )
964 
965  call prof_rapend ('FILE_O_NetCDF', 2)
966 
967  return
integer, public time_nowstep
current step [number]
Definition: scale_time.F90:70
integer, public is
start point of inner domain: x, local
module GTOOL_FILE
Definition: gtool_file.f90:17
integer, public je
end point of inner domain: y, local
real(dp), parameter, public rmiss
Definition: gtool_file.f90:141
subroutine, public interp_vertical_xi2z(var, var_Z)
Reset random seed.
module Gtool_History
integer, public ia
of x whole cells (local, with HALO)
integer, public ka
of z whole cells (local, with HALO)
integer, public kmax
of computational cells: z
integer, public js
start point of inner domain: y, local
module TIME
Definition: scale_time.F90:15
integer, public ks
start point of inner domain: z, local
integer, public ie
end point of inner domain: x, local
logical, public interp_available
vertical interpolation is available? (have meaning?)
module INTERPOLATION
integer, public ja
of y whole cells (local, with HALO)
Here is the call graph for this function:

◆ hist_in_0d()

subroutine scale_history::hist_in_0d ( real(rp), intent(in)  var,
character(len=*), intent(in)  item,
character(len=*), intent(in)  desc,
character(len=*), intent(in)  unit 
)

Wrapper routine of HIST_reg+HIST_put 0D.

Definition at line 977 of file scale_history.F90.

References hist_query(), hist_reg(), scale_prof::prof_rapend(), and scale_prof::prof_rapstart().

977  implicit none
978 
979  real(RP), intent(in) :: var
980  character(len=*), intent(in) :: item
981  character(len=*), intent(in) :: desc
982  character(len=*), intent(in) :: unit
983 
984  integer :: itemid
985  logical :: zinterp
986  logical :: do_put
987  !---------------------------------------------------------------------------
988 
989  if ( .NOT. enabled ) return
990 
991 
992  call hist_reg( itemid, & ! [OUT]
993  zinterp, & ! [OUT]
994  item, desc, unit, 0 ) ! [IN]
995 
996 
997  call hist_query( itemid, & ! [IN]
998  do_put ) ! [OUT]
999 
1000  if ( do_put ) then
1001  call hist_put( itemid, var ) ! [IN]
1002  endif
1003 
1004  return
Here is the call graph for this function:

◆ hist_get_2d()

subroutine scale_history::hist_get_2d ( real(rp), dimension(:,:), intent(out)  var,
character(len=*), intent(in)  basename,
character(len=*), intent(in)  varname,
integer, intent(in)  step,
logical, intent(in), optional  allow_missing 
)

Get 2D data from file.

Parameters
[out]varvalue
[in]basenamebasename of the file
[in]varnamename of the variable
[in]stepstep number

Definition at line 1201 of file scale_history.F90.

References scale_prof::prof_rapend(), and scale_prof::prof_rapstart().

1201  use gtool_history, only: &
1202  historyget
1203  implicit none
1204 
1205  real(RP), intent(out) :: var(:,:)
1206  character(len=*), intent(in) :: basename
1207  character(len=*), intent(in) :: varname
1208  integer, intent(in) :: step
1209 
1210  logical, intent(in), optional :: allow_missing
1211 
1212  logical :: am
1213  !---------------------------------------------------------------------------
1214 
1215  call prof_rapstart('FILE_I_NetCDF', 2)
1216 
1217  am = .false.
1218  if( present(allow_missing) ) am = allow_missing
1219 
1220  call historyget( var(:,:), & ! [OUT]
1221  basename, & ! [IN]
1222  varname, & ! [IN]
1223  step, & ! [IN]
1224  allow_missing=am ) ! [IN]
1225 
1226  call prof_rapend ('FILE_I_NetCDF', 2)
1227 
1228  return
module Gtool_History
Here is the call graph for this function:

◆ hist_get_3d()

subroutine scale_history::hist_get_3d ( real(rp), dimension(:,:,:), intent(out)  var,
character(len=*), intent(in)  basename,
character(len=*), intent(in)  varname,
integer, intent(in)  step,
logical, intent(in), optional  allow_missing 
)

Get 3D data from file.

Parameters
[out]varvalue
[in]basenamebasename of the file
[in]varnamename of the variable
[in]stepstep number

Definition at line 1239 of file scale_history.F90.

References scale_prof::prof_rapend(), and scale_prof::prof_rapstart().

1239  use gtool_history, only: &
1240  historyget
1241  implicit none
1242 
1243  real(RP), intent(out) :: var(:,:,:)
1244  character(len=*), intent(in) :: basename
1245  character(len=*), intent(in) :: varname
1246  integer, intent(in) :: step
1247 
1248  logical, intent(in), optional :: allow_missing
1249 
1250  logical :: am
1251  !---------------------------------------------------------------------------
1252 
1253  call prof_rapstart('FILE_I_NetCDF', 2)
1254 
1255  am = .false.
1256  if( present(allow_missing) ) am = allow_missing
1257 
1258  call historyget( var(:,:,:), & ! [OUT]
1259  basename, & ! [IN]
1260  varname, & ! [IN]
1261  step, & ! [IN]
1262  allow_missing=am ) ! [IN]
1263 
1264  call prof_rapend ('FILE_I_NetCDF', 2)
1265 
1266  return
module Gtool_History
Here is the call graph for this function:

◆ hist_write()

subroutine, public scale_history::hist_write ( )

Flush history buffer to file.

Definition at line 1272 of file scale_history.F90.

References gtool_history::historywriteall(), scale_prof::prof_rapend(), scale_prof::prof_rapstart(), and scale_time::time_nowstep.

Referenced by mod_rm_driver::scalerm().

1272  use gtool_history, only: &
1274  use scale_time, only: &
1275  time_nowstep
1276  implicit none
1277  !---------------------------------------------------------------------------
1278 
1279  call prof_rapstart('FILE_O_NetCDF', 2)
1280 
1281  call historywriteall( time_nowstep ) ![IN]
1282 
1283  call prof_rapend ('FILE_O_NetCDF', 2)
1284 
1285  return
integer, public time_nowstep
current step [number]
Definition: scale_time.F90:70
module Gtool_History
module TIME
Definition: scale_time.F90:15
subroutine, public historywriteall(step_now)
Here is the call graph for this function:
Here is the caller graph for this function: