scale_topography.F90

Go to the documentation of this file.

!-------------------------------------------------------------------------------
!-------------------------------------------------------------------------------
#include "scalelib.h"
module scale_topography
  !-----------------------------------------------------------------------------
  !
  !++ used modules
  !
  use scale_precision
  use scale_io
  use scale_prof
  use scale_atmos_grid_cartesc_index
  !-----------------------------------------------------------------------------
  implicit none
  private
  !-----------------------------------------------------------------------------
  !
  !++ Public procedure
  !
  public :: topography_setup
  public :: topography_fillhalo
  public :: topography_write
  public :: topography_calc_tan_slope
 
  !-----------------------------------------------------------------------------
  !
  !++ Public parameters & variables
  !
  logical, public :: topography_exist = .false. 
 
  real(rp), public, allocatable :: topography_zsfc   (:,:)
  real(rp), public, allocatable :: topography_tansl_x(:,:)
  real(rp), public, allocatable :: topography_tansl_y(:,:)
 
  !-----------------------------------------------------------------------------
  !
  !++ Private procedure
  !
  private :: topography_read
 
  !-----------------------------------------------------------------------------
  !
  !++ Private parameters & variables
  !
  character(len=H_LONG),  private :: topography_in_basename = ''
  character(len=H_LONG),  private :: topography_in_varname  = 'topo'
  logical,                private :: topography_in_aggregate
  logical,                private :: topography_in_check_coordinates = .false.
  character(len=H_LONG),  private :: topography_out_basename = ''
  logical,                private :: topography_out_aggregate
  character(len=H_MID),   private :: topography_out_title    = 'SCALE-RM TOPOGRAPHY'
  character(len=H_SHORT), private :: topography_out_dtype    = 'DEFAULT'
 
  !-----------------------------------------------------------------------------
contains
  !-----------------------------------------------------------------------------
  subroutine topography_setup
    use scale_file, only: &
       file_aggregate
    use scale_prc, only: &
       prc_abort
    implicit none
 
    namelist / param_topography / &
       topography_in_basename,          &
       topography_in_varname,           &
       topography_in_aggregate,         &
       topography_in_check_coordinates, &
       topography_out_basename,         &
       topography_out_aggregate,        &
       topography_out_dtype
 
    integer :: ierr
    !---------------------------------------------------------------------------
 
    log_newline
    log_info("TOPOGRAPHY_setup",*) 'Setup'
 
    topography_in_aggregate  = file_aggregate
    topography_out_aggregate = file_aggregate
 
    !--- read namelist
    rewind(io_fid_conf)
    read(io_fid_conf,nml=param_topography,iostat=ierr)
    if( ierr < 0 ) then !--- missing
       log_info("TOPOGRAPHY_setup",*) 'Not found namelist. Default used.'
    elseif( ierr > 0 ) then !--- fatal error
       log_error("TOPOGRAPHY_setup",*) 'Not appropriate names in namelist PARAM_TOPOGRAPHY. Check!'
       call prc_abort
    endif
    log_nml(param_topography)
 
    allocate( topography_zsfc(ia,ja) )
    allocate( topography_tansl_x(ia,ja) )
    allocate( topography_tansl_y(ia,ja) )
    topography_zsfc(:,:) = 0.0_rp
    topography_tansl_x(:,:) = 0.0_rp
    topography_tansl_y(:,:) = 0.0_rp
 
    ! read from file
    call topography_read
 
    return
  end subroutine topography_setup
 
  !-----------------------------------------------------------------------------
  subroutine topography_fillhalo( Zsfc, FILL_BND )
    use scale_comm_cartesc, only: &
       comm_vars8, &
       comm_wait
    implicit none
 
    real(rp), intent(inout), optional :: zsfc(ia,ja)
    logical,  intent(in),    optional :: fill_bnd
 
    logical :: fill_bnd_
    !---------------------------------------------------------------------------
 
    fill_bnd_ = .false.
    if ( present(fill_bnd) ) fill_bnd_ = fill_bnd
 
    if ( present(zsfc) ) then
       call comm_vars8( zsfc(:,:), 1 )
       call comm_wait ( zsfc(:,:), 1, fill_bnd_ )
    else
       call comm_vars8( topography_zsfc(:,:), 1 )
       call comm_wait ( topography_zsfc(:,:), 1, fill_bnd_ )
    end if
 
    return
  end subroutine topography_fillhalo
 
  !-----------------------------------------------------------------------------
  subroutine topography_read
    use scale_file_cartesc, only: &
       file_cartesc_open, &
       file_cartesc_read, &
       file_cartesc_flush, &
       file_cartesc_check_coordinates, &
       file_cartesc_close
    use scale_prc, only: &
       prc_abort
    implicit none
 
    integer :: fid
    !---------------------------------------------------------------------------
 
    log_newline
    log_info("TOPOGRAPHY_read",*) 'Input topography file '
 
    if ( topography_in_basename /= '' ) then
 
       call file_cartesc_open( topography_in_basename, fid, aggregate=topography_in_aggregate )
       call file_cartesc_read( fid, topography_in_varname, 'XY', topography_zsfc(:,:) )
 
       call file_cartesc_flush( fid )
 
       if ( topography_in_check_coordinates ) then
          call file_cartesc_check_coordinates( fid )
       end if
 
       call file_cartesc_close( fid )
 
       call topography_fillhalo( fill_bnd=.false. )
 
       topography_exist = .true.
 
    else
       log_info_cont(*) 'topography file is not specified.'
 
       topography_exist = .false.
    endif
 
    return
  end subroutine topography_read
 
  !-----------------------------------------------------------------------------
  subroutine topography_write
    use scale_file_cartesc, only: &
       file_cartesc_create,    &
       file_cartesc_def_var,   &
       file_cartesc_enddef,    &
       file_cartesc_write_var, &
       file_cartesc_close
    implicit none
 
    integer :: fid, vid
    !---------------------------------------------------------------------------
 
    if ( topography_out_basename /= '' .and. topography_out_basename /= topography_in_basename ) then
 
       log_newline
       log_info("TOPOGRAPHY_write",*) 'Output topography file '
 
       call topography_fillhalo( fill_bnd=.false. )
 
       call file_cartesc_create( topography_out_basename, topography_out_title, topography_out_dtype, & ! [IN]
                                 fid,                                                                 & ! [OUT]
                                 haszcoord=.false., aggregate=topography_out_aggregate                ) ! [IN]
 
       call file_cartesc_def_var( fid, 'topo', 'Topography', 'm', 'XY', topography_out_dtype, & ! [IN]
                                  vid,                                                        & ! [OUT]
                                  standard_name="surface_altitude"                            ) ! [IN]
 
       call file_cartesc_enddef( fid )
 
       call file_cartesc_write_var( fid, vid, topography_zsfc(:,:), 'topo', 'XY' ) ! [IN]
 
 
       call file_cartesc_close( fid )
 
    endif
 
    return
  end subroutine topography_write
 
  subroutine topography_calc_tan_slope( &
       IA, IS, IE, JA, JS, JE, &
       RCDX, RCDY, MAPF )
    use scale_prc_cartesc, only: &
       prc_twod
    implicit none
       integer,  intent(in) :: ia, is, ie
       integer,  intent(in) :: ja, js, je
       real(rp), intent(in) :: rcdx(ia), rcdy(ja)
       real(rp), intent(in) :: mapf(ia,ja,2)
 
       integer :: i, j
 
       if ( prc_twod ) then
          !$omp parallel do
          do j = js, je
          do i = is, ie
             topography_tansl_x(i,j) = 0.0_rp
          end do
          end do
       else
          !$omp parallel do
          do j = js, je
          do i = is, ie
             topography_tansl_x(i,j) = ( ( topography_zsfc(i+1,j) + topography_zsfc(i  ,j) ) * 0.5_rp &
                                       - ( topography_zsfc(i  ,j) + topography_zsfc(i-1,j) ) * 0.5_rp ) &
                                       * rcdx(i) * mapf(i,j,1)
          end do
          end do
       end if
       !$omp parallel do
       do j = js, je
       do i = is, ie
          topography_tansl_y(i,j) = ( ( topography_zsfc(i,j+1) + topography_zsfc(i,j  ) ) * 0.5_rp &
                                    - ( topography_zsfc(i,j  ) + topography_zsfc(i,j-1) ) * 0.5_rp ) &
                                  * rcdy(j) * mapf(i,j,2)
       end do
       end do
 
       call topography_fillhalo( topography_tansl_x(:,:), .true. )
       call topography_fillhalo( topography_tansl_y(:,:), .true. )
 
       return
     end subroutine topography_calc_tan_slope
 
end module scale_topography