c
c     =====================================================
      subroutine mapc2p(xc,yc,xp,yp)
c     =====================================================
c
c     # on input,  (xc,yc) is a computational grid point
c     # on output, (xp,yp) is corresponding point in physical space
c
      implicit double precision (a-h,o-z)
      common /cgrid/ igrid
c
c
      go to (100,200,300) igrid

  100 continue
c     # Cartesian (no mapping):
      xp = xc
      yp = yc
      go to 900

  200 continue
c     # Polar coordinates, xc = r,  yc = theta
      xp = xc * dcos(yc) 
      yp = xc * dsin(yc)
      go to 900

  300 continue
c     # radial projection mapping
      r1 = 1.d0
      d = dmax1(dabs(xc), dabs(yc))
      r = dmax1(dsqrt(xc**2 + yc**2), 1.d-10)
      xp = r1*d*xc/r
      yp = r1*d*yc/r

c     # interpolate between above mapping and Cartesian grid:
      w = dmax1(dabs(xc),dabs(yc))
      w = dmin1(w, 1.d0)      !# for ghost cells
      xp = w*xp + (1.d0-w)*xc
      yp = w*yp + (1.d0-w)*yc
      go to 900
c
  900 continue
      return
      end