get_mesh_index_lonlat_rot.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 get_mesh_index_lonlat_rot(KGRID_PAR,KL,PGRID_PAR,PLAT,PLON,KINDEX,KSSO,KISSOX,KISSOY)
!     ###############################################################
!
!!**** *GET_MESH_INDEX_LONLAT_ROT* get the grid mesh where point (lat,lon) is located
!!
!!    Note that this subroutine operates on the rotated grid. Thus incoming
!!    PLON, PLAT on regular grid is rotated.
!!
!!    PURPOSE
!!    -------
!!
!!    AUTHOR
!!    ------
!!
!!    P. Samuelsson  SMHI
!!
!!    MODIFICATION
!!    ------------
!!
!!    Original    12/2012
!!
!----------------------------------------------------------------------------
!
!*    0.     DECLARATION
!            -----------
!
USE modd_get_mesh_index_lonlat_rot, ONLY : xlonlim, xlatlim, nlat, nlon, xlon0, xpolon, xpolat
USE mode_gridtype_lonlat_rot
!
!
USE yomhook   ,ONLY : lhook,   dr_hook
USE parkind1  ,ONLY : jprb
!RJ: missing modi
USE modi_regrot_lonlat_rot
!
IMPLICIT NONE
!
!*    0.1    Declaration of arguments
!            ------------------------
!
INTEGER,                       INTENT(IN)    :: kgrid_par ! size of PGRID_PAR
INTEGER,                       INTENT(IN)    :: kl        ! number of points
REAL,    DIMENSION(KGRID_PAR), INTENT(IN)    :: pgrid_par ! grid parameters
REAL,    DIMENSION(KL),        INTENT(IN)    :: plat      ! latitude of the point
REAL,    DIMENSION(KL),        INTENT(IN)    :: plon      ! longitude of the point
INTEGER, DIMENSION(KL),        INTENT(OUT)   :: kindex    ! index of the grid mesh where the point is
INTEGER,                       INTENT(IN)    :: ksso      ! number of subgrid mesh in each direction
INTEGER, DIMENSION(KL),        INTENT(OUT)   :: kissox    ! X index of the subgrid mesh
INTEGER, DIMENSION(KL),        INTENT(OUT)   :: kissoy    ! Y index of the subgrid mesh
!
!*    0.2    Declaration of other local variables
!            ------------------------------------
!
INTEGER                           :: ji       ! loop counter in x
INTEGER                           :: jj       ! loop counter in y
INTEGER                           :: jl       ! loop counter on input points
!
REAL    :: zwest   ! West longitude in rotated grid (degrees)
REAL    :: zsouth  ! South latitude in rotated grid  (degrees)
REAL    :: zdlon   ! Longitudal grid spacing  (degrees)
REAL    :: zdlat   ! Latitudal grid spacing  (degrees)
!
REAL, DIMENSION(SIZE(PLON)) :: zlon
REAL, DIMENSION(SIZE(PLAT)) :: zlat
REAL(KIND=JPRB) :: zhook_handle
!
!----------------------------------------------------------------------------
!
IF (lhook) CALL dr_hook('GET_MESH_INDEX_LONLAT_ROT',0,zhook_handle)
IF (.NOT. ALLOCATED(xlatlim)) THEN
!
!*    1.     Uncode parameters of the grid
!            -----------------------------
!
  CALL get_gridtype_lonlat_rot(pgrid_par,                                 &
                                 zwest,zsouth,zdlon,zdlat,xpolon,xpolat,  &
                                 nlon,nlat                                )  
!
!----------------------------------------------------------------------------
!
!*    2.     Limits of grid meshes
!            ---------------------
!
  ALLOCATE(xlonlim(nlon+1))
  DO ji=1,nlon+1
    xlonlim(ji) = zwest + float(ji-1)*zdlon
  END DO
  xlonlim = xlonlim - zdlon/2.

  ALLOCATE(xlatlim(nlat+1))
  DO ji=1,nlat+1
    xlatlim(ji) = zsouth + float(ji-1)*zdlat
  END DO
  xlatlim = xlatlim - zdlat/2.
!
  xlon0 = 0.5*(minval(xlonlim)+maxval(xlonlim))
!
END IF
!----------------------------------------------------------------------------
!
!*    3.     Reshifts the longitudes with respect to projection reference point
!            ------------------------------------------------------------------
!
!
 CALL regrot_lonlat_rot(plon,plat,zlon,zlat,    &
                             kl,1,kl,1,        &
                             xpolon,xpolat,1   )  
!
WHERE (zlon(:)>180.) zlon(:)=zlon(:)-360.
WHERE (zlon(:)<-180.) zlon(:)=zlon(:)+360.
!
zlon(:) = zlon(:)+nint((xlon0-zlon(:))/360.)*360.
!
!----------------------------------------------------------------------------
!
!*    4.     Localisation of the data points on (x,y) grid
!            ---------------------------------------------
!
IF (kl/=nlon*nlat) THEN
  kindex = 0
  kissox = 0
  kissoy = 0
END IF
!
!
DO jl=1,SIZE(zlat)
  IF (     zlon(jl)<xlonlim(1) .OR. zlon(jl)>=xlonlim(nlon+1) &
        .OR. zlat(jl)<xlatlim(1) .OR. zlat(jl)>=xlatlim(nlat+1) ) THEN  
    kindex(jl) = 0
    IF (ksso/=0) THEN
      kissox(jl) = 0
      kissoy(jl) = 0
    END IF
    cycle
  END IF
  ji = count(zlon(jl)>=xlonlim(:))
  jj = count(zlat(jl)>=xlatlim(:))
  kindex(jl) = (jj-1) * nlon + ji
!
!
!*    6.     Localisation of the data points on in the subgrid of this mesh
!            --------------------------------------------------------------
!
  IF (ksso/=0) THEN
    kissox(jl) = 1 + int( float(ksso) * (zlon(jl)-xlonlim(ji))/(xlonlim(ji+1)-xlonlim(ji)) )
    kissoy(jl) = 1 + int( float(ksso) * (zlat(jl)-xlatlim(jj))/(xlatlim(jj+1)-xlatlim(jj)) )
  END IF
END DO
IF (lhook) CALL dr_hook('GET_MESH_INDEX_LONLAT_ROT',1,zhook_handle)
!
!-------------------------------------------------------------------------------
!
END SUBROUTINE get_mesh_index_lonlat_rot