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 84 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().

84  use gtool_history, only: &
86  use scale_process, only: &
87  prc_mpistop, &
90  use scale_rm_process, only: &
92  use scale_time, only: &
93  time_nowdate, &
94  time_nowms, &
95  time_dtsec, &
98  implicit none
99 
100  character(len=H_MID) :: history_h_title = 'SCALE-RM HISTORY OUTPUT'
101  character(len=H_SHORT) :: history_t_units = 'seconds'
102  character(len=H_MID) :: history_t_since = ''
103 
104  logical :: hist_bnd = .true.
105 
106  namelist / param_hist / &
107  hist_bnd
108 
109  real(DP) :: start_daysec
110  integer :: rankidx(2)
111  integer :: ierr
112  !---------------------------------------------------------------------------
113 
114  if( io_l ) write(io_fid_log,*)
115  if( io_l ) write(io_fid_log,*) '++++++ Module[HISTORY] / Categ[ATMOS-RM IO] / Origin[SCALElib]'
116 
117  !--- read namelist
118  rewind(io_fid_conf)
119  read(io_fid_conf,nml=param_hist,iostat=ierr)
120  if( ierr < 0 ) then !--- missing
121  if( io_l ) write(io_fid_log,*) '*** Not found namelist. Default used.'
122  elseif( ierr > 0 ) then !--- fatal error
123  write(*,*) 'xxx Not appropriate names in namelist PARAM_HIST. Check!'
124  call prc_mpistop
125  endif
126  if( io_lnml ) write(io_fid_log,nml=param_hist)
127 
128 
129  call prof_rapstart('FILE_O_NetCDF', 2)
130 
131  rankidx(1) = prc_2drank(prc_myrank, 1)
132  rankidx(2) = prc_2drank(prc_myrank, 2)
133 
134  start_daysec = time_startdaysec
135 ! if ( TIME_OFFSET_YEAR > 0 ) then
136 ! HISTORY_T_SINCE = '0000-01-01 00:00:00'
137 ! write(HISTORY_T_SINCE(1:4),'(i4)') TIME_OFFSET_YEAR
138 ! end if
139  if ( time_nowdate(1) > 0 ) then
140  write(history_t_since, '(i4.4,a1,i2.2,a1,i2.2,a1,i2.2,a1,i2.2,a1,i2.2)') &
141  time_nowdate(1),'-',time_nowdate(2),'-',time_nowdate(3),' ', &
142  time_nowdate(4),':',time_nowdate(5),':',time_nowdate(6)
143  start_daysec = time_nowms
144  end if
145 
146  if ( hist_bnd ) then
147  im = imaxb
148  jm = jmaxb
149  ims = isb
150  ime = ieb
151  jms = jsb
152  jme = jeb
153  else
154  im = imax
155  jm = jmax
156  ims = is
157  ime = ie
158  jms = js
159  jme = je
160  end if
161 
162  call historyinit( history_h_title, &
163  h_source, &
164  h_institute, &
165  im*jm*kmax, &
166  prc_masterrank, &
167  prc_myrank, &
168  rankidx, &
169  start_daysec, &
170  time_dtsec, &
171  time_units = history_t_units, &
172  time_since = history_t_since, &
173  namelist_fid = io_fid_conf )
174 
175  call hist_put_axes
176 
177  call prof_rapend ('FILE_O_NetCDF', 2)
178 
179  enabled = .true.
180 
181  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
subroutine, public prof_rapstart(rapname_base, level)
Start raptime.
Definition: scale_prof.F90:132
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
subroutine, public prof_rapend(rapname_base, level)
Save raptime.
Definition: scale_prof.F90:178
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 496 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().

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

604  use gtool_history, only: &
606  use scale_time, only: &
608  implicit none
609 
610  integer, intent(in) :: itemid
611  logical, intent(out) :: answer
612 
613  !---------------------------------------------------------------------------
614 
615  answer = .false.
616 
617  if ( .NOT. enabled ) return
618 
619  if ( itemid < 0 ) return
620 
621  call prof_rapstart('FILE_O_NetCDF', 2)
622 
623  call historyquery(itemid, time_nowstep, answer)
624 
625  call prof_rapend ('FILE_O_NetCDF', 2)
626 
627  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 prof_rapstart(rapname_base, level)
Start raptime.
Definition: scale_prof.F90:132
subroutine, public prof_rapend(rapname_base, level)
Save raptime.
Definition: scale_prof.F90:178
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 635 of file scale_history.F90.

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

635  use gtool_history, only: &
636  historyput
637  use scale_time, only: &
639  implicit none
640 
641  integer, intent(in) :: itemid
642  real(RP), intent(in) :: var
643 
644  !---------------------------------------------------------------------------
645 
646  if ( .NOT. enabled ) return
647 
648  if ( itemid < 0 ) return
649 
650  call prof_rapstart('FILE_O_NetCDF', 2)
651 
652  call historyput(itemid, time_nowstep, var)
653 
654  call prof_rapend ('FILE_O_NetCDF', 2)
655 
656  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 prof_rapstart(rapname_base, level)
Start raptime.
Definition: scale_prof.F90:132
subroutine, public prof_rapend(rapname_base, level)
Save raptime.
Definition: scale_prof.F90:178
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 664 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.

664  use gtool_history, only: &
665  historyput
666  use scale_time, only: &
668  implicit none
669 
670  integer, intent(in) :: itemid
671  real(RP), intent(in) :: var(:)
672 
673  real(RP) :: var2(kmax)
674  integer :: k
675  !---------------------------------------------------------------------------
676 
677  if ( .NOT. enabled ) return
678 
679  if ( itemid < 0 ) return
680 
681  call prof_rapstart('FILE_O_NetCDF', 2)
682 
683  do k = 1, kmax
684  var2(k) = var(ks+k-1)
685  enddo
686 
687  call historyput(itemid, time_nowstep, var2)
688 
689  call prof_rapend ('FILE_O_NetCDF', 2)
690 
691  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
subroutine, public prof_rapstart(rapname_base, level)
Start raptime.
Definition: scale_prof.F90:132
subroutine, public prof_rapend(rapname_base, level)
Save raptime.
Definition: scale_prof.F90:178
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 700 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.

700  use gtool_file, only: &
701  rmiss
702  use gtool_history, only: &
703  historyput
704  use scale_time, only: &
706  implicit none
707 
708  integer, intent(in) :: itemid
709  real(RP), intent(in) :: var(:,:)
710 
711  logical, intent(in), optional :: nohalo
712 
713  real(RP) :: var2(im*jm)
714  integer :: i, j
715  logical :: nohalo_
716  !---------------------------------------------------------------------------
717 
718  if ( .NOT. enabled ) return
719 
720  if ( itemid < 0 ) return
721 
722  call prof_rapstart('FILE_O_NetCDF', 2)
723 
724  do j = 1, jm
725  do i = 1, im
726  var2(i + (j-1)*im) = var(ims+i-1,jms+j-1)
727  enddo
728  enddo
729 
730  nohalo_ = .false.
731  if ( present(nohalo) ) nohalo_ = nohalo
732 
733  if ( nohalo_ ) then
734  ! W halo
735  do j = 1, jm
736  do i = 1, is-ims
737  var2(i + (j-1)*im) = rmiss
738  enddo
739  enddo
740  ! E halo
741  do j = 1, jm
742  do i = ie-ims+2, ime-ims+1
743  var2(i + (j-1)*im) = rmiss
744  enddo
745  enddo
746  ! S halo
747  do j = 1, js-jms
748  do i = 1, im
749  var2(i + (j-1)*im) = rmiss
750  enddo
751  enddo
752  ! N halo
753  do j = je-jms+2, jme-jms+1
754  do i = 1, im
755  var2(i + (j-1)*im) = rmiss
756  enddo
757  enddo
758  end if
759 
760  call historyput(itemid, time_nowstep, var2)
761 
762  call prof_rapend ('FILE_O_NetCDF', 2)
763 
764  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:110
module Gtool_History
integer, public js
start point of inner domain: y, local
module TIME
Definition: scale_time.F90:15
subroutine, public prof_rapstart(rapname_base, level)
Start raptime.
Definition: scale_prof.F90:132
integer, public ie
end point of inner domain: x, local
subroutine, public prof_rapend(rapname_base, level)
Save raptime.
Definition: scale_prof.F90:178
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 777 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.

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

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

979  implicit none
980 
981  real(RP), intent(in) :: var
982  character(len=*), intent(in) :: item
983  character(len=*), intent(in) :: desc
984  character(len=*), intent(in) :: unit
985 
986  integer :: itemid
987  logical :: zinterp
988  logical :: do_put
989  !---------------------------------------------------------------------------
990 
991  if ( .NOT. enabled ) return
992 
993 
994  call hist_reg( itemid, & ! [OUT]
995  zinterp, & ! [OUT]
996  item, desc, unit, 0 ) ! [IN]
997 
998 
999  call hist_query( itemid, & ! [IN]
1000  do_put ) ! [OUT]
1001 
1002  if ( do_put ) then
1003  call hist_put( itemid, var ) ! [IN]
1004  endif
1005 
1006  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 1202 of file scale_history.F90.

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

1202  use gtool_history, only: &
1203  historyget
1204  implicit none
1205 
1206  real(RP), intent(out) :: var(:,:)
1207  character(len=*), intent(in) :: basename
1208  character(len=*), intent(in) :: varname
1209  integer, intent(in) :: step
1210 
1211  logical, intent(in), optional :: allow_missing
1212 
1213  logical :: am
1214  !---------------------------------------------------------------------------
1215 
1216  call prof_rapstart('FILE_I_NetCDF', 2)
1217 
1218  am = .false.
1219  if( present(allow_missing) ) am = allow_missing
1220 
1221  call historyget( var(:,:), & ! [OUT]
1222  basename, & ! [IN]
1223  varname, & ! [IN]
1224  step, & ! [IN]
1225  allow_missing=am ) ! [IN]
1226 
1227  call prof_rapend ('FILE_I_NetCDF', 2)
1228 
1229  return
module Gtool_History
subroutine, public prof_rapstart(rapname_base, level)
Start raptime.
Definition: scale_prof.F90:132
subroutine, public prof_rapend(rapname_base, level)
Save raptime.
Definition: scale_prof.F90:178
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 1240 of file scale_history.F90.

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

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

◆ hist_write()

subroutine, public scale_history::hist_write ( )

Flush history buffer to file.

Definition at line 1273 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().

1273  use gtool_history, only: &
1275  use scale_time, only: &
1276  time_nowstep
1277  implicit none
1278  !---------------------------------------------------------------------------
1279 
1280  call prof_rapstart('FILE_O_NetCDF', 2)
1281 
1282  call historywriteall( time_nowstep ) ![IN]
1283 
1284  call prof_rapend ('FILE_O_NetCDF', 2)
1285 
1286  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 prof_rapstart(rapname_base, level)
Start raptime.
Definition: scale_prof.F90:132
subroutine, public historywriteall(step_now)
subroutine, public prof_rapend(rapname_base, level)
Save raptime.
Definition: scale_prof.F90:178
Here is the call graph for this function:
Here is the caller graph for this function: