subroutine ptclrgen(zxi,zth,zmu,zpphi,imech_flag,zvion,zksid,zvperp,zvpll,zmo)
  !
  !  from conserved quantities compute velocities for upcoming
  !  timestep ... see related routine ptclsave.f90
  !  dmc nov 16 2001
  use nbi_com
  
  implicit none
  real(fp), intent(in) :: zxi,zth         ! flux coordinate position of ptcl
  real(fp), intent(in) :: zmu             ! mu (cgs)
  real(fp), intent(in) :: zpphi          ! canonical or mechanical momentum (cgs)
 
  logical, intent(in) :: imech_flag     ! .true. if zpphi is mechanical
  ! momentum and mechanical, not canonical, momentum is to be preserved

  real(fp), intent(out) :: zvion,zksid    ! v,vpll/v (cm/sec & dimensionless)
  real(fp), intent(out) :: zvperp,zvpll   ! vperp,vpll (cm/sec)
  real(fp), intent(out) :: zmo            ! mechanical momentum (cgs)
 
  !--------------------------------------
 
  real(fp) zpsi_fac
 
  real(fp) zbi    ! 1/b, 1/gauss
 
  !--------------------------------------
 
  integer :: jx       ! radial zone index
 
  real(fp)    :: zx0      ! xi at radial zone for spline interpolation
  real(fp)    :: zh       ! offset into radial zone
 
  real(fp)    :: zkonst   ! toroidal flux / pi with normalization
 
  real(fp)    :: zg       ! (r*bt) (from spline)
  real(fp)    :: zpsi     ! poloidal flux (from spline integration)
 
  real(fp)    :: zb       ! z of beam ion
 
  real(fp) :: zdum1,zdum2,zdum3,zdum4,zdum5
  real(fp) :: zf1d,zpsint
  !--------------------------------------

  if(imech_flag) then
     zpsi_fac = 0.0_fp
  else
     zpsi_fac = 1.0_fp
  endif

  !--------------------------------------
  ! (1/b)
 
  call nbsplb(zxi,zth,0,zbi,zdum1,zdum2,zdum3,zdum4,zdum5)
 
  ! vperp from mu
 
  zvperp = sqrt(zmu/(0.5_fp*zmp*abeam*zbi))
 
  ! get g (r*bt) & psi for pphi
 
  jx  = int(zxi/dxi(lcentr,2)) + 1
  zx0 = real((jx-1),fp)*dxi(lcentr,2)
  zh  = zxi - zx0
 
  zkonst = 1.e8_fp*tflux/xpi
 
  zg = zf1d(zh,splg(1,jx),splg(2,jx),splg(3,jx),splg(4,jx))
  zpsi = zpsint(plflxi(jx),zkonst,zx0,zh,spliot(1,jx),  &
       spliot(2,jx),spliot(3,jx),spliot(4,jx))
 
  zpsi = zpsi_fac*zpsi

  ! momentum (input)
  !!! zpphi = abeam*zmp*zvpll*zg*zbi - xzbeam*zel*zpsi/zc
 
  zb = xzbeams(lsbeam)
  zvpll = (zpphi + zb*zel*zpsi/zc)/(abeam*zmp*zg*zbi)
 
  zmo = abeam*zmp*zvpll*zg*zbi
 
  zvion = sqrt(zvperp*zvperp+zvpll*zvpll)
  zksid = zvpll/zvion
 
  return
end subroutine ptclrgen
! 11sep2004 fgtok -s r8_precision.sub nubeam_test_r8.sub "r8con.csh conversion"
! 11sep2004 fgtok
