*-----------------------------------------------------------------------
*     MAP PROJECTION (POLAR STEREO)                    93/02/20 S.SAKAI
*-----------------------------------------------------------------------
*     Copyright (C) 2000-2016 GFD Dennou Club. All rights reserved.
*-----------------------------------------------------------------------
      SUBROUTINE MPFPST(XLON, YLAT, X, Y)

      PARAMETER (EPSL = 1.E-4)

      EXTERNAL  RFPI


      PI = RFPI()
      TH  = PI/2 - YLAT
      IF (TH .GE. PI-EPSL) TH = PI - EPSL
      R = 2*TAN(TH/2)
      X =  R*SIN(XLON)
      Y = -R*COS(XLON)

      RETURN
*-----------------------------------------------------------------------
      ENTRY MPIPST(X, Y, XLON, YLAT)

      PI = RFPI()
      R   = SQRT(X*X + Y*Y)
      YLAT = PI/2 - 2*ATAN(R/2)
      IF (R.EQ.0) THEN
        XLON = 0
      ELSE
        XLON = ATAN2(X, -Y)
      END IF

      END
