! dmc apr 2004 -- all xja routines use bmnmx instead of bmidp -- code
!  now in terms of min and max b on flux surfaces, not midplane values
!-----------------------------------------------------------------------
!  xjasbx  extrapolation of orbit intercept function (xjadep subsystem)
!
real(fp) function xjasbx(zx,isav,ze,zb,zrk,zpsik,zvpllk,zbk)
  !
  use nbi_com
  !
  use xjasub_mod
  !
  !--------------------------------
  !
  !
  !  input  zx=pt to extrapolate to
  !       isav= index to starting pt (must be inner or outer plasma
  !             bdy)
  !         ze= "energy" or v**2 of orbit at bdy
  !         zb= b field at bdy
  !
  !    also see xjasubs.blk common
  !
  !  output  zrk= radius extrapolated to
  !        zpsik= psi extrapolated to
  !          zbk= est. b field there
  !       zvpllk= est. vpll there
  !
  !       xjasbx= est. intercept fcn vpll**2+b*mub-v**2
  !               zeroes indicate orbit midplane intercept locations
  !
  !---------------------------------------------------------------------
  !
  !============
  implicit none
  integer :: isav,is1
  !============
  real(fp) :: ze,zb,zrk,zpsik,zvpllk,zbk,zx,zxtrap,zvphi,zf1,zf2,zf
  real(fp) :: zr0
  !============
  zr0=rmnmx(inzp1)
  !
  if(isav.eq.1) then
     is1=isav+1
  else
     is1=isav-1
  end if
  !
  zxtrap=(zx-1.0_fp)/dxi(ledge,2)
  !
  zrk=rmnmx(isav)+(rmnmx(isav)-rmnmx(is1))*zxtrap
  zrk=max(0.001_fp*rmnmx(isav),zrk)
  !
  zpsik=plflxg(lep1)+(plflxg(lep1)-plflxg(ledge))*zxtrap
  !
  zvphi=1.e7_fp*xe*(xzbeami*zpsik-zpsito)/(zmass*zrk)
  zvpllk=zvphi/btbmnmx(isav)
  zbk=zb*rmnmx(isav)/zrk
  !
  zf1=zvpllk**2
  zf2=ze-zmub*zbk
  !
  zf=zf1-zf2
  !
  xjasbx=zf
  !
  return
end function xjasbx
