latlontoxy1d.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 latlontoxy1d(PLAT0,PLON0,PRPK,PBETA,PLATOR,PLONOR, &
                                 pxhatm,pyhatm,plat,plon,kn,pradius)  
!      ################################################
!
!!****  *SM_LATLONTOXY1D * - Routine to compute conformal coordinates
!!
!!
!!     PURPOSE
!!     -------
!        This routine computes the cartesian conformal coordinates 
!      of an array given in latitude-longitude coordinates
!        Three map projections are available: 
!      - polar-stereographic (PRPK=1),
!      - lambert conformal  (0<PRPK<1),
!      - mercator (PRPK=0).
!
!
!!**   METHOD
!!     ------
!!       Spherical earth approximation is used. Longitude origin is 
!!     set in Greenwich, and is positive eastwards. An anticlockwise 
!!     rotation of XBETA degrees is applied to the conformal frame 
!!     with respect to the geographical directions.
!!
!!       WARNING: ALL INPUT AND OUTPUT ANGLES ARE IN DEGREES...
!!
!!     EXTERNAL
!!     --------
!!       None
!!
!!     EXPLICIT ARGUMENTS
!!     ------------------
!!       PXHAT,PYHAT(:)  : 1D arrays of the "velocity" gridpoints
!!                         cartesian conformal coordinates (meters,input).
!!       PLATOR   : Latitude of the (1,1) point of the "mass" grid
!!                      (degrees,input);
!!       PLONOR   : Longitude of the (1,1) point of the "mass" grid
!!                      (degrees,input);
!!       PXHATM   : conformal coordinate x  (meters, mass-grid, input)
!!       PYHATM   : conformal coordinate y  (meters, mass-grid, input)
!!       PLAT    : latitude                (degrees, mass-grid, output)
!!       PLON    : longitude               (degrees, mass-grid, output)
!!
!!     IMPLICIT ARGUMENTS
!!     ------------------
!!       Module MODD_CST         : contains Physical constants
!!          XPI         : Pi;    
!!          XRADIUS     : Earth radius (meters);
!!
!!
!!     REFERENCE
!!     ---------
!!      Asencio N. et al., 1994, "Le projet de modele non-hydrostatique
!!            commun CNRM-LA, specifications techniques", 
!!            Note CNRM/GMME, 26, 139p, (Chapter 2).
!!      Ducrocq V., 1994, "Generation de la grille dans le modele",
!!            Note interne MNH, 5 mai, 3p.
!!      Joly A., 1992, "Geographic parameters for ARPEGE/ALADIN",
!!            Internal note ARPEGE/ALADIN, february 27,28p.
!!      Levallois J., 1970, "Geodesie generale", Tome 2, Collection
!!             de l'IGN, Eyrolles, Paris, 408p.
!!       
!!     AUTHOR
!!     ------
!!      P.M.       *LA*
!!
!!     MODIFICATION
!!     ------------
!!       Original PM  24/05/94
!!       Updated  PM  27/07/94
!!       Updated  VD  23/08/94
!!       Updated  VM  24/10/95 projection from north pole (PRPK<0) and 
!!                             longitudes set between PLON0-180. and PLON0+180.
!!
!-------------------------------------------------------------------------------
!
!*     0.     DECLARATIONS
!             ------------
!
USE modd_csts
!
!
USE yomhook   ,ONLY : lhook,   dr_hook
USE parkind1  ,ONLY : jprb
!
IMPLICIT NONE
!
!*     0.1    Declarations of arguments and results
!
INTEGER,             INTENT(IN) :: kn
REAL,                INTENT(IN) :: plat0
REAL,                INTENT(IN) :: plon0
REAL,                INTENT(IN) :: prpk
REAL,                INTENT(IN) :: pbeta
REAL,                INTENT(IN) :: plator ! Latitude of the origine point
REAL,                INTENT(IN) :: plonor ! Longitude of the origine point
REAL, DIMENSION(KN), INTENT(IN) :: plat,plon
REAL, DIMENSION(KN), INTENT(OUT):: pxhatm,pyhatm
REAL, OPTIONAL,      INTENT(IN) :: pradius
!
!*     0.2    Declarations of local variables
! 
REAL,DIMENSION(KN) :: zlat,zlon
REAL :: zrpk,zlat0,zlon0,zlator,zlonor
REAL :: zrdsdg,zclat0,zslat0,zclator,zslator
REAL :: zro0,zga0,zbeta,zcgam,zsgam
REAL :: zxp,zyp,zraclat0,zxe,zye
REAL :: zradius
!
REAL,DIMENSION(KN) :: zclat,zslat,zro,zga,zxpr,zypr
REAL(KIND=JPRB) :: zhook_handle
!
!
!-------------------------------------------------------------------------------
!
!*     1.     PRELIMINARY CALCULATION FOR ALL PROJECTIONS
!             -------------------------------------------
!
IF (lhook) CALL dr_hook('LATLONTOXY1D',0,zhook_handle)
zrdsdg = xpi/180.         ! Degree to radian conversion factor
!
! By definition, (PLONOR,PLATOR) are the geographical 
! coordinates, and (ZXBM0,ZYBM0) the conformal cartesian 
! coordinates of the (1,1) point of the "mass" grid.
!
!
zlon(:)=plon(:)
zlon(:)=zlon(:)+nint((plon0-zlon(:))/360.)*360.
!
zlonor=plonor
zlonor=zlonor+nint((plon0-zlonor)/360.)*360.
!
zradius = xradius
IF (present(pradius)) zradius = pradius
!
!------------------------------------------------------------------------------
!
!*     2.     POLAR SEREOGRAPHIC AND LAMBERT CONFORMAL CASES
!             ----------------------------------------------
!                   (PRPK=1 P-stereo, 0<PRPK<1 Lambert)
!
IF(prpk  /=  0.) THEN
!
  IF (prpk<0.) THEN     ! projection from north pole
    zrpk=-prpk
    zbeta=-pbeta
    zlat0=-plat0
    zlon0=plon0+180.
    zlator=-plator
    zlonor=zlonor+180.
    zlat(:)=-plat(:)
    zlon(:)=zlon(:)+180.
  ELSE                  ! projection from south pole
    zrpk=prpk
    zbeta=pbeta
    zlat0=plat0
    zlon0=plon0
    zlator=plator
    zlonor=zlonor
    zlat(:)=plat(:)
    zlon(:)=zlon(:)
  ENDIF    
!
!*     2.1    Preliminary calculations
!
  zclat0  = cos(zrdsdg*zlat0)
  zslat0  = sin(zrdsdg*zlat0)
  zclator = cos(zrdsdg*zlator)
  zslator = sin(zrdsdg*zlator)
  zro0    = (zradius/zrpk)*(abs(zclat0))**(1.-zrpk)    &
            * ((1.+zslat0)*abs(zclator)/(1.+zslator))**zrpk  
  zga0    = (zrpk*(zlonor-zlon0)-zbeta)*zrdsdg
  zxp     = -zro0*sin(zga0)
  zyp     = zro0*cos(zga0)
!
!*    2.2    Conformal coordinates in meters
!
  zclat(:)  = cos(zrdsdg*zlat(:))
  zslat(:)  = sin(zrdsdg*zlat(:))
  zro(:)    = (zradius/zrpk)*(abs(zclat0))**(1.-zrpk)    &
               * ((1.+zslat0)*abs(zclat(:))/(1.+zslat(:)))**zrpk  
  zga(:)    = (zrpk*(zlon(:)-zlon0)-zbeta)*zrdsdg
!
  pxhatm(:) = zxp+zro(:)*sin(zga(:))
  pyhatm(:) = zyp-zro(:)*cos(zga(:))
!
  IF (prpk<0.) THEN     ! projection from north pole
    pyhatm(:)=-pyhatm(:)
  ENDIF
!
!-------------------------------------------------------------------------------
!
!*  3.        MERCATOR PROJECTION WITH ROTATION
!             ---------------------------------
!                       (PRPK=0)
!
ELSE
!
!*  3.1       Preliminary calculations
!
  zcgam    = cos(-zrdsdg*pbeta)
  zsgam    = sin(-zrdsdg*pbeta)
  zraclat0 = zradius*cos(zrdsdg*plat0)
  zxe      = - zraclat0*(plonor-plon0)*zrdsdg 
  zye      = - zraclat0*log(tan(xpi/4.+plator*zrdsdg/2.))
!
!*  3.2       Conformal coordinates
!
  zxpr(:)   = zraclat0*(zlon(:)-plon0)*zrdsdg+zxe
  zypr(:)   = zraclat0*log(tan(xpi/4.+plat(:)*zrdsdg/2.))+zye
  !
  pxhatm = zxpr(:)*zcgam-zypr(:)*zsgam
  pyhatm = zxpr(:)*zsgam+zypr(:)*zcgam
!
!-------------------------------------------------------------------------------
!
!*  4.        EXIT
!             ----
!
END IF
IF (lhook) CALL dr_hook('LATLONTOXY1D',1,zhook_handle)
!-------------------------------------------------------------------------------
END SUBROUTINE latlontoxy1d