src/SURFEX/pgd_grid.F90

   1 !SURFEX_LIC Copyright 1994-2014 Meteo-France
   2 !SURFEX_LIC This is part of the SURFEX software governed by the CeCILL-C  licence
   3 !SURFEX_LIC version 1. See LICENSE, CeCILL-C_V1-en.txt and CeCILL-C_V1-fr.txt
   4 !SURFEX_LIC for details. version 1.
   5 !     #########
   6       SUBROUTINE PGD_GRID(HPROGRAM,HFILE,HFILETYPE,OGRID,HGRID,KGRID_PAR,PGRID_PAR)
   7 !     ##########################################################
   8 !!
   9 !!    PURPOSE
  10 !!    -------
  11 !!   Reads in namelist the grid type and parameters.
  12 !!
  13 !!    METHOD
  14 !!    ------
  15 !!
  16 !!    EXTERNAL
  17 !!    --------
  18 !!
  19 !!
  20 !!    IMPLICIT ARGUMENTS
  21 !!    ------------------
  22 !!
  23 !!
  24 !!    REFERENCE
  25 !!    ---------
  26 !!
  27 !!    AUTHOR
  28 !!    ------
  29 !!
  30 !!    V. Masson                   Meteo-France
  31 !!
  32 !!    MODIFICATION
  33 !!    ------------
  34 !!
  35 !!    Original     01/2004
  36 !!    E. Martin    10/2007 IGN grid
  37 !!    M. Moge      05/02/2015 parallelization (using local sizes, GET_MEAN_OF_COORD_SQRT_ll, SET_NAM_GRID_CONF_PROJ_LOCAL) + MPPDB_CHECK
  38 !!    M. Moge      01/03/2015 call SPLIT_GRID if CPROGRAM == 'PGD   ' + remove SET_NAM_GRID_CONF_PROJ_LOCAL
  39 !!    M. Moge      01/03/2015 change in the input arguments of PGD_GRID_IO_INIT : passing IDXRATIO, IDYRATIO
  40 !----------------------------------------------------------------------------
  41 !
  42 !*    0.     DECLARATION
  43 !            -----------
  44 !
  45 USE MODD_SURFEX_MPI,     ONLY : NSIZE, NINDEX, NPIO, NRANK
  46 USE MODD_SURFEX_OMP,     ONLY : NINDX2, NWORK, XWORK, XWORK2, XWORK3, &
  47                                 NWORK_FULL, XWORK_FULL, XWORK2_FULL
  48 !
  49 USE MODD_PGD_GRID,       ONLY : NL, XGRID_PAR, NGRID_PAR, XMESHLENGTH
  50 USE MODN_PGD_GRID
  51 USE MODD_SURF_ATM_GRID_n, ONLY : XLAT, XLON, XMESH_SIZE, XJPDIR
  52 USE MODD_SURF_ATM_n,     ONLY : NDIM_FULL, NSIZE_FULL
  53 USE MODD_CSTS,           ONLY : XPI, XRADIUS
  54 !
  55 USE MODI_DEFAULT_GRID
  56 USE MODI_GRID_FROM_FILE
  57 USE MODI_OPEN_NAMELIST
  58 USE MODI_TEST_NAM_VAR_SURF
  59 USE MODI_CLOSE_NAMELIST
  60 USE MODI_GET_LUOUT
  61 USE MODI_READ_NAM_GRIDTYPE
  62 USE MODI_LATLON_GRID
  63 !
  64 USE MODE_POS_SURF
  65 !
  66 !
  67 USE YOMHOOK   ,ONLY : LHOOK,   DR_HOOK
  68 USE PARKIND1  ,ONLY : JPRB
  69 !
  70 USE MODI_ABOR1_SFX
  71 !
  72 USE MODI_PGD_GRID_IO_INIT
  73 #ifdef MNH_PARALLEL
  74 USE MODE_TOOLS_ll, ONLY : GET_MEAN_OF_COORD_SQRT_ll
  75 !
  76 USE MODI_GET_SIZE_FULL_n
  77 USE MODI_SPLIT_GRID
  78 USE MODD_CONF, ONLY : CPROGRAM
  79 #endif
  80 !
  81 IMPLICIT NONE
  82 !
  83 !*    0.1    Declaration of dummy arguments
  84 !            ------------------------------
  85 !
  86  CHARACTER(LEN=6),  INTENT(IN)   :: HPROGRAM   ! program calling the surface
  87  CHARACTER(LEN=28), INTENT(IN)   :: HFILE      ! atmospheric file name
  88  CHARACTER(LEN=6),  INTENT(IN)   :: HFILETYPE  ! atmospheric file type
  89 LOGICAL,           INTENT(IN)   :: OGRID      ! .true. if grid is imposed by atm. model
  90  CHARACTER(LEN=10), INTENT(OUT)  :: HGRID      ! grid type
  91 INTEGER,           INTENT(OUT)  :: KGRID_PAR  ! size of PGRID_PAR
  92 REAL, DIMENSION(:), POINTER     :: PGRID_PAR  ! parameters defining this grid
  93 !
  94 !
  95 !*    0.2    Declaration of local variables
  96 !            ------------------------------
  97 !
  98 INTEGER           :: ILUOUT ! output listing logical unit
  99 INTEGER           :: ILUNAM ! namelist file  logical unit
 100 LOGICAL           :: GFOUND ! Flag true if namelist is present
 101 REAL(KIND=JPRB) :: ZHOOK_HANDLE
 102 INTEGER :: IIMAX_ll, IJMAX_ll ! global size of son model
 103 LOGICAL :: GRECT
 104 !
 105  INTEGER :: IXOR = 1            ! position of modified bottom left point
 106  INTEGER :: IYOR = 1            ! according to initial grid
 107  INTEGER :: IXSIZE = -999       ! number of grid meshes in initial grid to be
 108  INTEGER :: IYSIZE = -999       ! covered by the modified grid
 109  INTEGER :: IDXRATIO = 1        ! resolution ratio between modified grid
 110  INTEGER :: IDYRATIO = 1        ! and initial grid
 111 NAMELIST/NAM_INIFILE_CONF_PROJ/IXOR,IYOR,IXSIZE,IYSIZE,IDXRATIO,IDYRATIO
 112 !
 113 !*    0.3    Declaration of namelists
 114 !            ------------------------
 115 !
 116 !------------------------------------------------------------------------------
 117 !
 118 !*    1.      Defaults
 119 !             --------
 120 !
 121 IF (LHOOK) CALL DR_HOOK('PGD_GRID',0,ZHOOK_HANDLE)
 122  CALL DEFAULT_GRID(HPROGRAM,CGRID)
 123 !
 124 YINIFILE  = '                         '
 125 YFILETYPE = '      '
 126 !
 127 IF (OGRID) THEN
 128   YINIFILE  = HFILE
 129   YFILETYPE = HFILETYPE
 130 END IF
 131 !
 132  CALL GET_LUOUT(HPROGRAM,ILUOUT)
 133 !------------------------------------------------------------------------------
 134 !
 135 !*    2.      Open namelist
 136 !             -------------
 137 !
 138 IF (.NOT. OGRID) THEN
 139   CALL OPEN_NAMELIST(HPROGRAM,ILUNAM)
 140 !
 141 !------------------------------------------------------------------------------
 142 !
 143 !*    3.      Read grid type
 144 !             --------------
 145 !
 146   CALL POSNAM(ILUNAM,'NAM_PGD_GRID',GFOUND,ILUOUT)
 147   IF (GFOUND) READ(UNIT=ILUNAM,NML=NAM_PGD_GRID)
 148 !
 149 !------------------------------------------------------------------------------
 150 !
 151 !*    5.      Close namelist
 152 !             --------------
 153 !
 154   CALL CLOSE_NAMELIST(HPROGRAM,ILUNAM)
 155 !
 156 END IF
 157 !-------------------------------------------------------------------------------
 158 !
 159 !*       4.    check of grid and input file types
 160 !              ----------------------------------
 161 !
 162  CALL TEST_NAM_VAR_SURF(ILUOUT,'CGRID',CGRID,'CONF PROJ ','NONE      ','LONLAT REG','CARTESIAN ','GAUSS     ',&
 163           'IGN       ','LONLATVAL ')
 164  CALL TEST_NAM_VAR_SURF(ILUOUT,'YFILETYPE',YFILETYPE,'      ','MESONH','LFI   ','ASCII ')
 165 !
 166 !
 167 !------------------------------------------------------------------------------
 168 !
 169 !*    5.      Initializes grid characteristics
 170 !             --------------------------------
 171 !
 172 !*    5.1     From another file
 173 !             -----------------
 174 !
 175 IF (LEN_TRIM(YFILETYPE)>0 .AND. LEN_TRIM(YINIFILE)>0 ) THEN
 176   IF (YFILETYPE=='MESONH' .OR. YFILETYPE=='LFI   ' .OR. YFILETYPE=='ASCII ') THEN
 177     CALL GRID_FROM_FILE(HPROGRAM,YINIFILE,YFILETYPE,OGRID,CGRID,NGRID_PAR,XGRID_PAR,NL)
 178     HGRID     = CGRID
 179     IF ( HGRID == "IGN       " .OR. HGRID == "GAUSS     " .OR. HGRID == "NONE      " ) THEN
 180       GRECT = .FALSE.
 181     ELSE
 182       GRECT = .TRUE.
 183     ENDIF
 184     ! on lit la taille globale du modele fils dans la namelist
 185     CALL OPEN_NAMELIST(HPROGRAM,ILUNAM)
 186     CALL POSNAM(ILUNAM,'NAM_INIFILE_CONF_PROJ',GFOUND,ILUOUT)
 187     IF (GFOUND) THEN
 188       READ(UNIT=ILUNAM,NML=NAM_INIFILE_CONF_PROJ)
 189       IIMAX_ll = IXSIZE*IDXRATIO
 190       IJMAX_ll = IYSIZE*IDYRATIO
 191     ENDIF
 192     CALL CLOSE_NAMELIST(HPROGRAM,ILUNAM)
 193
 194     !*    3.      Additional actions for I/O
 195     !
 196     IF (GFOUND) THEN
 197 #ifdef MNH_PARALLEL
 198       CALL PGD_GRID_IO_INIT(HPROGRAM,NGRID_PAR,XGRID_PAR, HGRID, GRECT, IIMAX_ll, IJMAX_ll, IDXRATIO, IDYRATIO)
 199 #else
 200       CALL PGD_GRID_IO_INIT(HPROGRAM)
 201 #endif
 202       NDIM_FULL = NL
 203     ELSE
 204 #ifdef MNH_PARALLEL
 205       CALL PGD_GRID_IO_INIT(HPROGRAM,NGRID_PAR,XGRID_PAR, HGRID, GRECT)
 206 #else
 207       CALL PGD_GRID_IO_INIT(HPROGRAM)
 208 #endif
 209     ENDIF
 210     NSIZE     = NDIM_FULL
 211 #ifdef MNH_PARALLEL
 212     CALL GET_SIZE_FULL_n(HPROGRAM,NDIM_FULL,NSIZE_FULL)
 213     NL = NSIZE_FULL
 214 #else
 215     NSIZE_FULL = NL
 216 #endif
 217   ELSE
 218     CALL ABOR1_SFX('PGD_GRID: FILE TYPE NOT SUPPORTED '//HFILETYPE//' FOR FILE '//HFILE)
 219   END IF
 220   !we don't need to call SPLIT_GRID, the grid has been splitted in GRID_FROM_FILE
 221 !
 222 ELSE
 223 !
 224 !*    5.2     Grid not initialized
 225 !             --------------------
 226 !
 227   IF (CGRID=='NONE      ' .OR. CGRID=='          ') THEN
 228     CALL ABOR1_SFX('PGD_GRID: GRID TYPE NOT INITIALIZED, CGRID='//CGRID)
 229
 230 !
 231 !*    5.3     Grid initialized
 232 !             ----------------
 233 !
 234   ELSE
 235 !
 236     CALL READ_NAM_GRIDTYPE(HPROGRAM,CGRID,NGRID_PAR,XGRID_PAR,NL)
 237     HGRID     = CGRID
 238     !*    3.      Additional actions for I/O
 239     !
 240 #ifdef MNH_PARALLEL
 241     CALL PGD_GRID_IO_INIT(HPROGRAM,NGRID_PAR,XGRID_PAR)
 242 #else
 243     CALL PGD_GRID_IO_INIT(HPROGRAM)
 244 #endif
 245     NDIM_FULL = NL
 246     NSIZE     = NDIM_FULL
 247 #ifdef MNH_PARALLEL
 248     CALL GET_SIZE_FULL_n(HPROGRAM,NDIM_FULL,NSIZE_FULL)
 249     NL = NSIZE_FULL
 250 #else
 251     NSIZE_FULL = NL
 252 #endif
 253   END IF
 254 #ifdef MNH_PARALLEL
 255   ! IF we are in PREP_PGD, we need to split the grid. Otherwise, the grid was read in parallel and is already splitted
 256   IF ( CPROGRAM == 'PGD   ') THEN
 257     CALL SPLIT_GRID('MESONH',NGRID_PAR,XGRID_PAR)
 258   ENDIF
 259 #endif
 260
 261 END IF
 262 !
 263 !
 264 IF (.NOT.ALLOCATED(NINDEX)) THEN
 265   ALLOCATE(NINDEX(NDIM_FULL))
 266   NINDEX(:) = 0
 267 ENDIF
 268 NINDX2    = NDIM_FULL
 269 ALLOCATE(NWORK(NDIM_FULL))
 270 ALLOCATE(XWORK(NDIM_FULL))
 271 ALLOCATE(XWORK2(NDIM_FULL,10))
 272 ALLOCATE(XWORK3(NDIM_FULL,10,10))
 273 IF (NRANK==NPIO) THEN
 274   ALLOCATE(NWORK_FULL(NDIM_FULL))
 275   ALLOCATE(XWORK_FULL(NDIM_FULL))
 276   ALLOCATE(XWORK2_FULL(NDIM_FULL,10))
 277 ELSE
 278   ALLOCATE(NWORK_FULL(0))
 279   ALLOCATE(XWORK_FULL(0))
 280   ALLOCATE(XWORK2_FULL(0,0))
 281 ENDIF
 282 !
 283 KGRID_PAR = NGRID_PAR
 284 ALLOCATE(PGRID_PAR(KGRID_PAR))
 285 PGRID_PAR = XGRID_PAR
 286 !
 287 !------------------------------------------------------------------------------
 288 !
 289 !*    6.      Latitude and longitude
 290 !             ----------------------
 291 !
 292 ALLOCATE(XLAT       (NSIZE_FULL))
 293 ALLOCATE(XLON       (NSIZE_FULL))
 294 ALLOCATE(XMESH_SIZE (NSIZE_FULL))
 295 ALLOCATE(XJPDIR     (NSIZE_FULL))
 296  CALL LATLON_GRID(CGRID,NGRID_PAR,NSIZE_FULL,ILUOUT,XGRID_PAR,XLAT,XLON,XMESH_SIZE,XJPDIR)
 297 !
 298 !------------------------------------------------------------------------------
 299 !
 300 !*    7.      Average grid length (in degrees)
 301 !             --------------------------------
 302 !
 303 !* in meters
 304 #ifdef MNH_PARALLEL
 305 CALL GET_MEAN_OF_COORD_SQRT_ll(XMESH_SIZE,NSIZE_FULL,NDIM_FULL,XMESHLENGTH)
 306 #else
 307 XMESHLENGTH = SUM ( SQRT(XMESH_SIZE) ) / NL
 308 #endif
 309 !
 310 !* in degrees (of latitude)
 311 XMESHLENGTH = XMESHLENGTH *180. / XPI / XRADIUS
 312 IF (LHOOK) CALL DR_HOOK('PGD_GRID',1,ZHOOK_HANDLE)
 313 !
 314 !-------------------------------------------------------------------------------
 315 !
 316 END SUBROUTINE PGD_GRID