subroutine map_xe(stcprm, x_map, xi, eta,flag)
      parameter (pi=3.14159265358979323846,pi_2=pi/2.)
      character * 1 flag
      dimension stcprm(15),x_map(3)
        if (abs(x_map(3)) .ge. 1.) then
C* Projection pole or its antipodes
            xi = 0.
          if (stcprm(1) * x_map(3) .gt. 0.) then
            eta = 1./stcprm(1)
          else
            if (x_map(3) .gt. 0) then
              eta = 1.0e10
            else
              eta = -1.0e10
            endif
          endif
          return
        else
          if (abs(stcprm(1)) .eq. 1.) then
C* Stereographic Case, away from pole
            fact = 1. / (stcprm(1) + x_map(3) )
            xi = x_map(2)*fact
            eta = 1./stcprm(1) - x_map(1) * fact
          else
            ymerc = .5 * log( ( 1. + x_map(3) ) / (1. - x_map(3) ) )
C* This Projection has a cut.  Select theta according to the rule
C* If cutflag = 'e', -PI/2 <= theta < 3PI/2, if cutflag = 'E',
C* 0 <= theta < 2PI, if cutflag = 'w', -3Pi/2 <= theta < PI/2,
C* cutflag = 'W', -2PI <= theta < 0., else -PI<=theta<PI. 
          if (flag .eq. 'E') then
            theta = atan2(- x_map(2), - x_map(1) ) + pi
          else if (flag .eq. 'W') then
            theta = atan2(- x_map(2), - x_map(1) ) - pi
          else if (flag .eq. 'e') then
            theta = atan2(- x_map(1), x_map(2) )  + pi_2
          else if (flag .eq. 'w') then
            theta = atan2(x_map(1), - x_map(2) )  - pi_2
          else
            theta = atan2(x_map(2), x_map(1) )
          endif
            rhog = xpabva( stcprm(1), - ymerc)
C* rhog = ( exp( - gamma * ymerc ) - 1. ) / gamma */
            xi = (1. + stcprm(1) * rhog) * snabva(stcprm(1), theta)
            eta = stcprm(1) * (1. + stcprm(1) * rhog) *
     A                        csabva(stcprm(1), theta) - rhog
          end if
        end if
        return
      end