subroutine nbo_load1(ztype,zname,zlabel,zunits,zdata,id,zvtot)
 
  ! load a single nubeam output profile into xplasma
 
  use nbi_com
  use xplasma_obj_instance
 
  implicit none
!============
  character*(*),intent(in) :: ztype  ! profile type (see below)
  character*(*),intent(in) :: zname  ! profile xplasma name
  character*(*),intent(in) :: zlabel ! profile label
  character*(*),intent(in) :: zunits ! profile physical units

  real(fp), intent(in) :: zdata(mj)      ! data to write
 
  integer, intent(out) :: id         ! xplasma id
  real(fp), intent(out) :: zvtot       ! total profile contents
  !  (if an integrated profile) or else pt of maximum modulus.

  !  note id=zvtot=0 is returned if zdata=0.0 everywhere
 
  !  profile types (ztype):
  !    density -- store an integrated profile -- spline or piecewise linear
  !       (volume integrated)
  !    curdens -- store an integrated profile -- spline or piecewise linear
  !       (area integrated, for driven currents)
  !   parameter -- store as a profile -- spline or piecewise linear
  !   linparam -- store as a profile -- unsmoothed -- piecewise linear
  !   flow -- store as a boundary oriented profile -- unsmoothed
  !           -- piecewise linear
 
  !---------------------------
  !  local...
 
  integer :: id_rho=0,iorder,i,j,ibc1,ibc2,ierr,iparam
  real(fp) :: zdatwk(mj),zbc1,zbc2,ztmp

  character*32 znamb

#ifdef __DEBUG

  integer :: inxx,inext,iprev
  real(fp), dimension(:), allocatable :: zxx,zf1,zf2

#endif
 
  !---------------------------
 
  ierr = 0
 
  znamb='nbi_'//trim(zname)
  oprefix='nbi_'
  zdatwk=0.0_fp 

  !---------------------------
  !  (a) initiate x axes if necessary
 
  call nbo_load1_xizb
 
  call nbo_load1_xizc
 
  !---------------------------
  !  (b) store profile.  buf if identically zero, skip...
 
  iparam = 0

  if(maxval(abs(zdata(lcentr:lep1))).eq.0.0_fp) then
 
     id=0
     zvtot=0
 
  else
 
     ! there is data... process according to type...
 
     if((ztype.eq.'curdens').or.(ztype.eq.'density').or. &
          (ztype.eq.'parameter')) then
        if(dxbsmoo.le.0.0_fp) then
           iorder=0  ! piecewise linear (replacement allowed)
        else
           iorder=2  ! spline (replacement allowed)
        endif
     else
        iorder=0     ! always piecewise linear (replacement allowed)
     endif
 
     ibc1 = 0        ! default bc's for splines: not a knot
     ibc2 = 0
     zbc1 = 0
     zbc2 = 0
 
     if((ztype.eq.'density').or.(ztype.eq.'curdens')) then
 
        !  "density" profile -- spline the integral...
        !  (note that nubeam internal data is not volume-normalized)
        !  (so units are #, #/sec, watts, amps, etc;
        !   not #/cm3, #/cm3/sec, amps/cm2, etc.)
 
        zdatwk(1)=0.0_fp
        i=1
        do j=lcentr,ledge
           zdatwk(i+1)=zdatwk(i)+zdata(j)
           i=i+1
        enddo
 
        if(iorder.gt.0) then
           ibc1=1  ! f(rho)*v(rho)=0 => slope of integral =0 @ axis
        endif
 
        call xplasma_create_1dprof(s,znamb,id_rhozb,zdatwk(1:nzones+1),id, &
             ierr, &
             ispline=iorder,ibca=ibc1,zbca=zbc1,ibcb=ibc2,zbcb=zbc2, &
             label=zlabel,units=zunits)
 
        zvtot = zdatwk(nzones+1)
        if(zvtot.eq.0.0_fp) zvtot = 1.0e-100_fp  ! do not allow zero

     else if((ztype.eq.'parameter').or.(ztype.eq.'linparam')) then
 
        iparam = 1

        !  parameter profile -- supply the data itself; extend to edges
 
        zdatwk(2:nzones+1) = zdata(lcentr:ledge)
        zdatwk(1)=zdatwk(2)
        zdatwk(nzones+2)=zdatwk(nzones+1)
 
        call xplasma_create_1dprof(s,znamb,id_rhozc,zdatwk(1:nzones+2),id, &
             ierr, &
             ispline=iorder,ibca=ibc1,zbca=zbc1,ibcb=ibc2,zbcb=zbc2, &
             label=zlabel,units=zunits)
 
        zvtot = 0
        do i=1,nzones+2
           ztmp=abs(zdatwk(i))
           if(ztmp.gt.abs(zvtot)) zvtot=zdatwk(i)
        enddo
        if(zvtot.eq.0.0_fp) zvtot = 1.0e-100_fp  ! do not allow zero

     else if(ztype.eq.'flow') then
 
        !  flow profile -- assume zero on axis
 
        zdatwk(1:nzones+1)=zdata(lcentr:lep1)
 
        call xplasma_create_1dprof(s,znamb,id_rhozb,zdatwk(1:nzones+1),id, &
             ierr, &
             ispline=iorder,ibca=ibc1,zbca=zbc1,ibcb=ibc2,zbcb=zbc2, &
             label=zlabel,units=zunits)
 
        zvtot = zdatwk(nzones+1)
        if(zvtot.eq.0.0_fp) zvtot = 1.0e-100_fp  ! do not allow zero
 
     else
        ierr=88888
        write(nonlin,*) ' ?nbo_load1:  unexpected ztype: "'//ztype//'"'
        write(nonlin,*) '  failed to load: ',trim(zname)
     endif
  endif
 
  if((ierr.ne.0).and.(ierr.ne.88888)) then
     write(nonlin,*) ' ??nbo_load1:  unexpected error.'
     write(nonlin,*) &
          ' ?nbo_load1: xplasma error on "'//ztype//'" data: "'//zname,'".'
     write(nonlin,*) '  failed to load: ',trim(zname)
     call xplasma_error(s,ierr,nonlin)
  endif
 
#ifdef __DEBUG
  if((ierr.eq.0).and.(id.gt.0)) then
     allocate(zxx(10*nzones+1),zf1(10*nzones+1),zf2(10*nzones+1))

     zxx(1)=xizb(lcentr)
     zf1(1)=zdatwk(1)
     inext=1
     do i=1,nzones
        iprev=inext
        inext=iprev+10
        zxx(inext)=xizb(lcentr+i)
        zf1(inext)=zdatwk(1+i)
        do j=1,9
           zxx(iprev+j)=zxx(iprev)+j*(xizb(lcentr+i)-xizb(lcentr+i-1))/10
           if(iparam.eq.1) then
              zf1(iprev+j)=zdatwk(1+i)
           else
              zf1(iprev+j)=zdatwk(i)+j*(zdatwk(i+1)-zdatwk(i))/10
           endif
        enddo
     enddo

     call xplasma_eval_prof(s,id,zxx,zf2,ierr)

!     call r8_grafx2(zxx,zf1,zf2,10*nzones+1,'x',zunits, &
!          trim(runid)//' nbo_load1 debug plot',zlabel,zname)

     deallocate(zxx,zf1,zf2)
  endif
#endif
end subroutine nbo_load1
 
subroutine nbo_load1_xizb
 
  !  define boundary oriented radial grid for nubeam
  !  only do this if the axis is not already defined.
 
  use nbi_com
  use xplasma_obj_instance
 
  implicit none
  integer :: id_rho=0,ierr
  real(fp) :: zxwk(mj),ztol
 
  !-----------------------------------
 
  call xplasma_gridid(s,'rhozb_nubeam',id_rhozb)
  if(id_rhozb.ne.0) then
     call xplasma_grid_size(s,id_rhozb,num_rhozb,ierr)
     if(num_rhozb.eq.nzones+1) then
        return  ! ok -- axis already defined
     endif
  endif
 
  call xplasma_gridid(s,'__rho',id_rho)
 
  ! setup grid for boundary oriented quantities
 
  num_rhozb=nzones+1
  zxwk(1)=0.0_fp; zxwk(num_rhozb)=1.0_fp      ! end pts
  zxwk(2:nzones)=xizb_r8(lcp1:ledge)          ! interior pts
 
  call xplasma_author_set(s,'nubeam',ierr)
  if(ierr.ne.0) then
     write(nonlin,*) &
          ' ?nbo_load1_xizb: unexpected error: xplasma_author_set:'
     call xplasma_error(s,ierr,nonlin)
     return
  endif

  ztol=1.0e-4_fp
  call xplasma_create_grid(s,'rhozb_nubeam',xplasma_rho_coord, &
       zxwk(1:num_rhozb), &
       id_rhozb,ierr)
 
  if(ierr.ne.0) then
     write(nonlin,*) &
          ' ?nbo_load1_xizb: unexpected error: xplasma_create_grid(rhozb...)'
     call xplasma_error(s,ierr,nonlin)
  endif
 
  call xplasma_author_clear(s,'nubeam',ierr)
  if(ierr.ne.0) then
     write(nonlin,*) &
          ' ?nbo_load1_xizb: unexpected error: xplasma_author_set:'
     call xplasma_error(s,ierr,nonlin)
  endif

end subroutine nbo_load1_xizb
 
subroutine nbo_load1_xizc
 
  !  define zone center oriented radial grid for nubeam
  !  only do this if the axis is not already defined.
 
  use nbi_com
  use xplasma_obj_instance
 
  implicit none
  integer :: id_rho=0,ierr
  real(fp) zxwk(mj),ztol
 
  !-----------------------------------
 
  call xplasma_gridid(s,'rhozc_nubeam',id_rhozc)
  if(id_rhozc.ne.0) then
     call xplasma_grid_size(s,id_rhozc,num_rhozc,ierr)
     if(num_rhozc.eq.nzones+2) then
        return  ! ok -- axis already defined
     endif
  endif
 
  call xplasma_gridid(s,'__rho',id_rho)
 
  ! setup grid for zone oriented quantities
 
  num_rhozc=nzones+2
  zxwk(1)=0.0_fp; zxwk(num_rhozc)=1.0_fp      ! end pts
  zxwk(2:nzones+1)=xizc_r8(lcentr:ledge)          ! interior pts
 
  call xplasma_author_set(s,'nubeam',ierr)
  if(ierr.ne.0) then
     write(nonlin,*) &
          ' ?nbo_load1_xizc: unexpected error: xplasma_author_set:'
     call xplasma_error(s,ierr,nonlin)
     return
  endif

  ztol=1.0e-4_fp
  call xplasma_create_grid(s,'rhozc_nubeam',xplasma_rho_coord, &
       zxwk(1:num_rhozc), &
       id_rhozc,ierr)
 
  if(ierr.ne.0) then
     write(nonlin,*) &
          ' ?nbo_load1_xizc: unexpected error: xplasma_create_grid(rhozc...)'
     call xplasma_error(s,ierr,nonlin)
  endif
 
  call xplasma_author_clear(s,'nubeam',ierr)
  if(ierr.ne.0) then
     write(nonlin,*) &
          ' ?nbo_load1_xizc: unexpected error: xplasma_author_clear:'
     call xplasma_error(s,ierr,nonlin)
  endif
 
end subroutine nbo_load1_xizc
