subroutine nbmomry_xs(zxii,zthi,ijac,zr,zz,zjac)
!added aug.24, 2015 
  use iso_c_binding, only: fp => c_double
  use xplasma_obj_instance
 
  !  !use xplasma spline maps to evaluate what was done by the legacy
  !  momry routine
 
  implicit none
  real(fp), intent(in) :: zxii,zthi        ! mag. coordinates in
  integer, intent(in) :: ijac        ! =0: no jacobian; .gt.0: jacobian
  real(fp), intent(out) :: zr,zz         ! r,z coordinates out
  real(fp), intent(out) :: zjac(2,2)     ! 2x2 jacobian out
 
  !
  !  zjac(1,1) = dr/dxi     zjac(1,2) = dr/dth
  !  zjac(2,1) = dz/dxi     zjac(2,2) = dz/dth
  !
 
  real(fp), parameter :: ceps4 = 1.0e-4_fp
  real(fp), parameter :: ceps13= 1.0e-13_fp
  real(fp), parameter :: tol_par=1.e-6_fp
  real(fp) :: ztol,zbdy_tol,zjac11,zjac12,zjac21,zjac22
  logical, save :: ilset=.false.
  integer :: ierr
 
  !---------------------------------------------
  ! normalize range
  if(.not.ilset) then
     ztol=max(ceps13,min(ceps4,tol_par))
     call xplasma_bdytol_set(s,ztol,ierr) 
     call xplasma_global_info(s,ierr,bdytol=zbdy_tol)
     ilset=.true.
  endif
  if(ijac.eq.0)then
     call xplasma_ctrans(s, ierr, rho_in = zxii, theta_in = zthi, tol=tol_par, & 
          &  r_out=zr, z_out=zz) !r,z in m this routine about 10 times slower than nbmomry
     zr=zr*1.0e+2_fp !cm
     zz=zz*1.0e+2_fp !cm

  else
     call xplasma_rzjac1(s, zxii, zthi, ierr,r1=zr,z1=zz,drdrho1=zjac11,&
          &  dzdrho1=zjac21,drdtheta1=zjac12,dzdtheta1=zjac22)
     zjac(1,1)=zjac11*1.0e+2_fp !cm	  
     zjac(2,1)=zjac21*1.0e+2_fp !cm	  
     zjac(1,2)=zjac12*1.0e+2_fp !cm	  
     zjac(2,2)=zjac22*1.0e+2_fp !cm	  
     zr=zr*1.0e+2_fp !cm
     zz=zz*1.0e+2_fp !cm
          
  endif

  return 
end subroutine nbmomry_xs
! 11sep2004 fgtok -s r8_precision.sub nubeam_test_r8.sub "r8con.csh conversion"
