pgd_grid.F90

Go to the documentation of this file.

!SFX_LIC Copyright 1994-2014 CNRS, Meteo-France and Universite Paul Sabatier
!SFX_LIC This is part of the SURFEX software governed by the CeCILL-C licence
!SFX_LIC version 1. See LICENSE, CeCILL-C_V1-en.txt and CeCILL-C_V1-fr.txt  
!SFX_LIC for details. version 1.
!     #########
      SUBROUTINE pgd_grid (UG, U, GCP, HPROGRAM,HFILE,HFILETYPE,OGRID,HDIR)
!     ##########################################################
!!
!!    PURPOSE
!!    -------
!!   Reads in namelist the grid type and parameters.
!!
!!    METHOD
!!    ------
!!   
!!    EXTERNAL
!!    --------
!!
!!
!!    IMPLICIT ARGUMENTS
!!    ------------------
!!
!!
!!    REFERENCE
!!    ---------
!!
!!    AUTHOR
!!    ------
!!
!!    V. Masson                   Meteo-France
!!
!!    MODIFICATION
!!    ------------
!!
!!    Original     01/2004
!!    E. Martin    10/2007 IGN grid
!!    P. Samuelsson  12/2012 Rotated lonlat
!!    M. Moge      05/02/2015 parallelization (using local sizes, GET_MEAN_OF_COORD_SQRT_ll, SET_NAM_GRID_CONF_PROJ_LOCAL) + MPPDB_CHECK
!!    M. Moge      01/03/2015 call SPLIT_GRID if CPROGRAM == 'PGD   ' + remove SET_NAM_GRID_CONF_PROJ_LOCAL
!!    M. Moge      01/03/2015 change in the input arguments of PGD_GRID_IO_INIT : passing IDXRATIO, IDYRATIO
!----------------------------------------------------------------------------
!
!*    0.     DECLARATION
!            -----------
!
USE modd_surf_atm_grid_n, ONLY : surf_atm_grid_t
USE modd_surf_atm_n, ONLY : surf_atm_t
USE modd_grid_conf_proj_n, ONLY : grid_conf_proj_t 
!
USE modd_surfex_mpi,     ONLY : nsize, npio, nrank, ncomm, nproc
!
USE modd_pgd_grid,       ONLY : nl, xgrid_par, ngrid_par, xmeshlength
USE modn_pgd_grid
USE modd_csts,           ONLY : xpi, xradius
!
USE modi_default_grid
USE modi_grid_from_file
USE modi_open_namelist
USE modi_test_nam_var_surf
USE modi_close_namelist
USE modi_get_luout
USE modi_read_nam_gridtype
USE modi_latlon_grid
USE modi_abor1_sfx
USE modi_pgd_grid_io_init
USE modi_gather_and_write_mpi
!
USE mode_pos_surf
!
#ifdef MNH_PARALLEL
USE modd_conf, ONLY : cprogram
USE mode_tools_ll, ONLY : get_mean_of_coord_sqrt_ll
!
USE modi_get_size_full_n
USE modi_split_grid
#endif
!
USE yomhook   ,ONLY : lhook,   dr_hook
USE parkind1  ,ONLY : jprb
!
IMPLICIT NONE
!
#ifdef SFX_MPI
include "mpif.h"
#endif
!
!*    0.1    Declaration of dummy arguments
!            ------------------------------
!
TYPE(surf_atm_grid_t), INTENT(INOUT) :: UG
TYPE(surf_atm_t), INTENT(INOUT) :: U
TYPE(grid_conf_proj_t),INTENT(INOUT) :: GCP
!
 CHARACTER(LEN=6),  INTENT(IN)   :: HPROGRAM   ! program calling the surface
 CHARACTER(LEN=28), INTENT(IN)   :: HFILE      ! atmospheric file name
 CHARACTER(LEN=6),  INTENT(IN)   :: HFILETYPE  ! atmospheric file type
LOGICAL,           INTENT(IN)   :: OGRID      ! .true. if grid is imposed by atm. model
 CHARACTER(LEN=1), INTENT(IN), OPTIONAL :: HDIR
!
!
!*    0.2    Declaration of local variables
!            ------------------------------
!
REAL, DIMENSION(:), ALLOCATABLE :: ZMESH_SIZE
 CHARACTER(LEN=1) :: YDIR
INTEGER :: ILUOUT ! output listing logical unit
INTEGER :: ILUNAM ! namelist file  logical unit
INTEGER :: IIMAX_ll, IJMAX_ll ! global size of son model
INTEGER :: ISIZE_FULL
LOGICAL           :: GFOUND ! Flag true if namelist is present
LOGICAL :: GRECT
INTEGER :: INFOMPI
!
REAL(KIND=JPRB) :: ZHOOK_HANDLE
!
!*    0.3    Declaration of namelists
!            ------------------------
!
 INTEGER :: IXOR = 1            ! position of modified bottom left point
 INTEGER :: IYOR = 1            ! according to initial grid
 INTEGER :: IXSIZE = -999       ! number of grid meshes in initial grid to be
 INTEGER :: IYSIZE = -999       ! covered by the modified grid
 INTEGER :: IDXRATIO = 1        ! resolution ratio between modified grid
 INTEGER :: IDYRATIO = 1        ! and initial grid
NAMELIST/nam_inifile_conf_proj/ixor,iyor,ixsize,iysize,idxratio,idyratio
!
!------------------------------------------------------------------------------
!
!*    1.      Defaults
!             --------
!
IF (lhook) CALL dr_hook('PGD_GRID',0,zhook_handle)
!
ydir = 'H'
IF (PRESENT(hdir)) ydir = hdir
!
IF (ASSOCIATED(xgrid_par)) DEALLOCATE(xgrid_par)
!
 CALL default_grid(hprogram,cgrid)
!
yinifile  = '                         '
yinifiletype = '      '
!
IF (ogrid) THEN
  yinifile  = hfile
  yinifiletype = hfiletype
END IF
!
 CALL get_luout(hprogram,iluout)
!------------------------------------------------------------------------------
!
!*    2.      Open namelist
!             -------------
!
IF (.NOT. ogrid) THEN
!
  CALL open_namelist(hprogram,ilunam)
!
!------------------------------------------------------------------------------
!
!*    3.      Read grid type
!             --------------
!
  CALL posnam(ilunam,'NAM_PGD_GRID',gfound,iluout)
  IF (gfound) READ(unit=ilunam,nml=nam_pgd_grid)
!
!------------------------------------------------------------------------------
!
!*    5.      Close namelist
!             --------------
!
  CALL close_namelist(hprogram,ilunam)
!
END IF
!-------------------------------------------------------------------------------
!
!*       4.    check of grid and input file types
!              ----------------------------------
! 
 CALL test_nam_var_surf(iluout,'CGRID',cgrid,'CONF PROJ ','NONE      ','LONLAT REG','CARTESIAN ','GAUSS     ',&
          'IGN       ','LONLATVAL ','LONLAT ROT')  
 CALL test_nam_var_surf(iluout,'YINIFILETYPE',yinifiletype,'      ','MESONH','LFI   ','ASCII ','NC    ')
!
!
!------------------------------------------------------------------------------
!
!*    5.      Initializes grid characteristics
!             --------------------------------
!
!*    5.1     From another file
!             -----------------
!
IF (len_trim(yinifiletype)>0 .AND. len_trim(yinifile)>0 ) THEN
  !
  IF (yinifiletype=='MESONH' .OR. yinifiletype=='LFI   ' .OR. yinifiletype=='ASCII ' .OR. yinifiletype=='NC    ') THEN
    !
    CALL grid_from_file(u,gcp,ug%XGRID_FULL_PAR,hprogram,yinifile,yinifiletype,&
                        ogrid,cgrid,ngrid_par,xgrid_par,nl,ydir)
    !
    ug%G%CGRID  = cgrid
    !
    IF ( cgrid == "IGN       " .OR. cgrid == "GAUSS     " .OR. cgrid == "NONE      " ) THEN
      grect = .false.
    ELSE
      grect = .true.
    ENDIF
    !
    ! on lit la taille globale du modele fils dans la namelist
    CALL open_namelist(hprogram,ilunam)
    CALL posnam(ilunam,'NAM_INIFILE_CONF_PROJ',gfound,iluout)
    IF (gfound) THEN 
      READ(unit=ilunam,nml=nam_inifile_conf_proj)
      iimax_ll = ixsize*idxratio
      ijmax_ll = iysize*idyratio
    ENDIF
    CALL close_namelist(hprogram,ilunam)
    !
    !*    3.      Additional actions for I/O
    !
    IF (gfound) THEN 
#ifdef MNH_PARALLEL
      CALL pgd_grid_io_init(hprogram,ug,ngrid_par,xgrid_par,cgrid,grect,iimax_ll,ijmax_ll,idxratio,idyratio)
#else
      CALL pgd_grid_io_init(hprogram,ug)
#endif
      IF (ydir/='H') u%NDIM_FULL = nl
    ELSE
#ifdef MNH_PARALLEL
      CALL pgd_grid_io_init(hprogram,ug,ngrid_par,xgrid_par,cgrid,grect)
#else
      CALL pgd_grid_io_init(hprogram,ug)
#endif
    ENDIF
    !
#ifdef MNH_PARALLEL
    CALL get_size_full_n(hprogram,u%NDIM_FULL,u%NSIZE_FULL,isize_full)
    u%NSIZE_FULL = isize_full
    nl = u%NSIZE_FULL
#else
    u%NSIZE_FULL = nl    
#endif
  ELSE
    CALL abor1_sfx('PGD_GRID: FILE TYPE NOT SUPPORTED '//hfiletype//' FOR FILE '//hfile)
  END IF
!
ELSE
!
!*    5.2     Grid not initialized
!             --------------------
!
  IF (cgrid=='NONE      ' .OR. cgrid=='          ') THEN
    CALL abor1_sfx('PGD_GRID: GRID TYPE NOT INITIALIZED, CGRID='//cgrid)

!
!*    5.3     Grid initialized
!             ----------------
!
  ELSE
!
    CALL read_nam_gridtype(gcp,ug%XGRID_FULL_PAR,u%NDIM_FULL,hprogram,cgrid,ngrid_par,xgrid_par,nl,ydir)
!
    !*    3.      Additional actions for I/O
    !
    ug%G%CGRID  = cgrid
    !
#ifdef MNH_PARALLEL
    CALL pgd_grid_io_init(hprogram,ug,ngrid_par,xgrid_par)
#else
    CALL pgd_grid_io_init(hprogram,ug)
#endif
    IF (ydir/='H') u%NDIM_FULL = nl
#ifdef MNH_PARALLEL
    CALL get_size_full_n(hprogram,nl,u%NSIZE_FULL,isize_full)
    u%NSIZE_FULL = isize_full
    nl = u%NSIZE_FULL
#endif
    !
  END IF  

#ifdef MNH_PARALLEL
  ! IF we are in PREP_PGD, we need to split the grid. Otherwise, the grid was read in parallel and is already splitted
  IF ( cprogram == 'PGD   ') THEN
    CALL split_grid(ug,u,'MESONH',ngrid_par,xgrid_par)
  ENDIF
#endif
    
END IF
!
IF (hdir=='A') THEN
  ug%NGRID_FULL_PAR = ngrid_par
  ALLOCATE(ug%XGRID_FULL_PAR(ug%NGRID_FULL_PAR))
  ug%XGRID_FULL_PAR = xgrid_par        
ELSE
  ug%G%NGRID_PAR = ngrid_par
  ALLOCATE(ug%G%XGRID_PAR(ug%G%NGRID_PAR))
  ug%G%XGRID_PAR = xgrid_par
ENDIF
!
IF (ydir=='A') THEN  
  CALL dr_hook('PGD_GRID',1,zhook_handle)
  RETURN
ENDIF
!
!------------------------------------------------------------------------------
!
!*    6.      Latitude and longitude
!             ----------------------
!
ALLOCATE(ug%G%XLAT       (nl))
ALLOCATE(ug%G%XLON       (nl))
ALLOCATE(ug%G%XMESH_SIZE (nl))
ALLOCATE(ug%XJPDIR       (nl))
 CALL latlon_grid(ug%G,nl,ug%XJPDIR)
!
!------------------------------------------------------------------------------
!
!*    7.      Average grid length (in degrees)
!             --------------------------------
!
!* in meters
#ifdef MNH_PARALLEL
 CALL get_mean_of_coord_sqrt_ll(ug%G%XMESH_SIZE,u%NSIZE_FULL,u%NDIM_FULL,xmeshlength)
#else
IF (nrank==npio) ALLOCATE(zmesh_size(u%NDIM_FULL))
 CALL gather_and_write_mpi(ug%G%XMESH_SIZE,zmesh_size)
IF (nrank==npio) xmeshlength = sum( sqrt(zmesh_size) ) / max(u%NDIM_FULL,1)
IF (nproc>1) THEN
#ifdef SFX_MPI  
  CALL mpi_bcast(xmeshlength,kind(xmeshlength)/4,mpi_real,npio,ncomm,infompi)
#endif
ENDIF
IF (nrank==npio) DEALLOCATE(zmesh_size)
#endif
!
!* in degrees (of latitude)
xmeshlength = xmeshlength *180. / xpi / xradius
IF (lhook) CALL dr_hook('PGD_GRID',1,zhook_handle)
!
!-------------------------------------------------------------------------------
!
END SUBROUTINE pgd_grid