subroutine rfxqlo_wavefield_dims(nx,nticrfp1,ierr)

  !  return expected wavefield array dimensions...
  !  this is a convenience, but caller still has to know meaning of radial
  !  and poloidal dimensions...

  use rfxqlo_private
  implicit NONE

  integer,intent(out) :: nx        ! radial dimension
  integer,intent(out) :: nticrfp1  ! poloidal dimension (+1: duplicate endpts)
  integer,intent(out) :: ierr      ! status code (0=OK)

  if(.not.ready) then
     write(lunmsg,*) ' ?rfxqlo_wavefield_dims: RFXQLO module not ready.'
     write(lunmsg,*) '  specify all inputs; call rfxqlo_input_verify first.'
     ierr=1
  else
     ierr=0
     nx=nx_qlo
     nticrfp1 = rfx_grid%nchi_in + 1
  endif

end subroutine rfxqlo_wavefield_dims

subroutine rfxqlo_wavefield_nx(nx,ierr)
  ! return wavefield radial (x) grid size -- Number of Zones
  ! number of Boundaries is nx+1

  use rfxqlo_private
  implicit NONE

  integer,intent(out) :: nx        ! radial dimension
  integer,intent(out) :: ierr      ! status code (0=OK)

  if(.not.ready) then
     write(lunmsg,*) ' ?rfxqlo_wavefield_nx: RFXQLO module not ready.'
     write(lunmsg,*) '  specify all inputs; call rfxqlo_input_verify first.'
     ierr=1
  else
     ierr=0
     nx=nx_qlo
  endif

end subroutine rfxqlo_wavefield_nx

subroutine rfxqlo_wavefield_x(x,ierr)
  ! return wavefield radial grid (x) -- passed array must be large enough.
  !   sqrt(Phi/Philim)

  use rfxqlo_private
  implicit NONE

  real*8,dimension(:) :: x         ! radial grid (x)
  integer,intent(out) :: ierr      ! status code (0=OK)

  if(.not.ready) then
     write(lunmsg,*) ' ?rfxqlo_wavefield_x: RFXQLO module not ready.'
     write(lunmsg,*) '  specify all inputs; call rfxqlo_input_verify first.'
     ierr=1
  else if(size(x).lt.nx_qlo) then
     write(lunmsg,*) ' ?rfxqlo_wavefield_x: passed array not large enough'
     write(lunmsg,*) '  x grid (need=',nx_qlo,' got=',size(x),')'
     ierr=1
  else
     x(1:nx_qlo)=x_qlo(1:nx_qlo)
     ierr=0
  endif

end subroutine rfxqlo_wavefield_x

subroutine rfxqlo_wavefield_xb(xb,ierr)
  ! return wavefield radial grid boundaries (xb) -- passed array must be 
  ! large enough.
  !   sqrt(Phi/Philim)

  use rfxqlo_private
  implicit NONE

  real*8,dimension(:) :: xb        ! radial grid (xb)
  integer,intent(out) :: ierr      ! status code (0=OK)

  if(.not.ready) then
     write(lunmsg,*) ' ?rfxqlo_wavefield_x: RFXQLO module not ready.'
     write(lunmsg,*) '  specify all inputs; call rfxqlo_input_verify first.'
     ierr=1
  else if(size(xb).lt.nx_qlo+1) then
     write(lunmsg,*) ' ?rfxqlo_wavefield_xb: passed array not large enough'
     write(lunmsg,*) '  xb grid (need=',nx_qlo+1,' got=',size(xb),')'
     ierr=1
  else
     xb(1:nx_qlo+1)=xb_qlo(0:nx_qlo)
     ierr=0
  endif

end subroutine rfxqlo_wavefield_xb

subroutine rfxqlo_wavefield_xpsi(xpsi,ierr)
  ! return wavefield radial grid (x) -- passed array must be large enough.
  ! at Zone Boundaries! 
  !   sqrt(Psi/Psilim)

  use rfxqlo_private
  implicit NONE

  real*8,dimension(:) :: xpsi      ! radial grid (xpsi)
  integer,intent(out) :: ierr      ! status code (0=OK)

  if(.not.ready) then
     write(lunmsg,*) ' ?rfxqlo_wavefield_x: RFXQLO module not ready.'
     write(lunmsg,*) '  specify all inputs; call rfxqlo_input_verify first.'
     ierr=1
  else if(size(xpsi).lt.nx_qlo+1) then
     write(lunmsg,*) ' ?rfxqlo_wavefield_x: passed array not large enough'
     write(lunmsg,*) '  x grid (need=',nx_qlo,' got=',size(xpsi),')'
     ierr=1
  else
     call eq_rgetf(nx_qlo+1,xb_qlo,rfx_grid%id_xpsi,0,xpsi,ierr)
     ierr=0
  endif

end subroutine rfxqlo_wavefield_xpsi
  
subroutine rfxqlo_wavefield_id_xpsi(id_xpsi,ierr)
  ! just return the Xplasma id of the profile sqrt(Psi/Psilim)

  use rfxqlo_private
  implicit NONE

  integer, intent(out) :: id_xpsi    ! id returned (0 if error occurs)
  integer, intent(out) :: ierr       ! status code (0 if no error)

  if(.not.ready) then
     write(lunmsg,*) ' ?rfxqlo_wavefield_x: RFXQLO module not ready.'
     write(lunmsg,*) '  specify all inputs; call rfxqlo_input_verify first.'
     ierr=1
     id_xpsi=0
  else
     ierr=0
     id_xpsi=rfx_grid%id_xpsi
  endif

end subroutine rfxqlo_wavefield_id_xpsi

subroutine rfxqlo_set_coupling(mode_coupling,ierr)

  !  set the mode coupling (or load factor in units such as MW/kA^2)
  !  if RF code sets for one n_phi mode, it must set for all.  If none
  !  are set, all couplings are at unity and the vacuum power spectrum
  !  is retained.

  use rfxqlo_private
  implicit NONE

  real*8, intent(in) :: mode_coupling
  integer, intent(out) :: ierr    ! status code (0=OK)

  !------------------------
  integer :: iant, iphi
  !------------------------

  if(.not.ready) then
     write(lunmsg,*) ' ?rfxqlo_set_coupling: RFXQLO module not ready.'
     write(lunmsg,*) '  specify all inputs; call rfxqlo_input_verify first.'
     ierr=1

  else
     ierr=0

     iant = selected_antenna
     iphi = selected_nphi

     wavefields(iphi,iant)%mode_coupling = mode_coupling
  endif
end subroutine rfxqlo_set_coupling

subroutine rfxqlo_get_mode_power(mode_power,istate,ierr)

  !  get the mode power associated with currently selected wave component; 
  !  (a single nphi mode on a single antenna); also return state code.

  use rfxqlo_private
  implicit NONE

  real*8, intent(out) :: mode_power  ! power in wavefield source

  integer, intent(out) :: istate     ! =0: power is zero, none was calculated;
  !            =1: vacuum normalization;
  !            =2: global renormalization, accounting for coupling efficiency.

  integer, intent(out) :: ierr    ! status code (0=OK)

  !------------------------
  integer :: iant, iphi
  !------------------------

  if(.not.ready) then
     write(lunmsg,*) ' ?rfxqlo_set_coupling: RFXQLO module not ready.'
     write(lunmsg,*) '  specify all inputs; call rfxqlo_input_verify first.'
     ierr=1
     mode_power = 0.0d0
     istate = 0

  else
     ierr=0

     iant = selected_antenna
     iphi = selected_nphi

     mode_power = wavefields(iphi,iant)%mode_power
     istate = wavefields(iphi,iant)%nstate
  endif

end subroutine rfxqlo_get_mode_power

subroutine rfxqlo_wavefield_set1(skperp,eplus,crfpol,ierr)

  !  "nstate=1" setup of wavefield for current (antenna, nphi)

  use rfxqlo_private
  implicit NONE

  real*8, dimension(:,:), intent(in) :: skperp, eplus
  complex*16, dimension(:,:), intent(in) :: crfpol

  integer, intent(out) :: ierr    ! status code (0=OK)

  !------------------------
  integer :: iant, iphi, nticrf
  !------------------------

  if(.not.ready) then
     write(lunmsg,*) ' ?rfxqlo_wavefield_set1: RFXQLO module not ready.'
     write(lunmsg,*) '  specify all inputs; call rfxqlo_input_verify first.'
     ierr=1
  else
     ierr=0

     iant = selected_antenna
     iphi = selected_nphi

     nticrf=rfx_grid%nchi_in

     call ckdims('skperp',skperp,ierr)
     if(ierr.ne.0) return

     call ckdims('eplus',eplus,ierr)
     if(ierr.ne.0) return

     call cckdims('crfpol',crfpol,ierr)
     if(ierr.ne.0) return

     wavefields(iphi,iant)%nstate=1
     wavefields(iphi,iant)%skperp = skperp
     wavefields(iphi,iant)%eplus = eplus
     wavefields(iphi,iant)%crfpol = crfpol

     wavefields(iphi,iant)%mode_power = antenna_array(iant)%power * &
          antenna_array(iant)%wnphi_vac(iphi)

  endif

  contains

    subroutine ckdims(lbl,data,ierr)

      character*(*) :: lbl   ! name of data
      real*8, dimension(:,:) :: data
      integer, intent(out) :: ierr

      ierr=0

      if(size(data,1).ne.nx_qlo) then
         write(lunmsg,*) ' ?rfxqlo_wavefield_set1(ckdims): ',trim(lbl), &
              ' 1st (radial) dimension mismatch:'
         write(lunmsg,*) '  expected: ',nx_qlo,'  got: ',size(data,1)
         ierr=1
      endif

      if(size(data,2).ne.nticrf+1) then
         write(lunmsg,*) ' ?rfxqlo_wavefield_set1(ckdims): ',trim(lbl), &
              ' 2nd (poloidal) dimension mismatch:'
         write(lunmsg,*) '  expected: ',nticrf+1,'  got: ',size(data,2)
         ierr=1
      endif

    end subroutine ckdims

    subroutine cckdims(lbl,data,ierr)

      character*(*) :: lbl   ! name of data
      complex*16, dimension(:,:) :: data
      integer, intent(out) :: ierr

      ierr=0

      if(size(data,1).ne.nx_qlo) then
         write(lunmsg,*) ' ?rfxqlo_wavefield_set1(cckdims): ',trim(lbl), &
              ' 1st (radial) dimension mismatch:'
         write(lunmsg,*) '  expected: ',nx_qlo,'  got: ',size(data,1)
         ierr=1
      endif

      if(size(data,2).ne.nticrf+1) then
         write(lunmsg,*) ' ?rfxqlo_wavefield_set1(cckdims): ',trim(lbl), &
              ' 2nd (poloidal) dimension mismatch:'
         write(lunmsg,*) '  expected: ',nticrf+1,'  got: ',size(data,2)
         ierr=1
      endif

    end subroutine cckdims
end subroutine rfxqlo_wavefield_set1

subroutine rfxqlo_wavefield_get1(nstate, skperp,eplus,crfpol,ierr)

  !  retrieve wavefield for current (antenna, nphi)

  use rfxqlo_private
  implicit NONE

  integer, intent(out) :: nstate
  real*8, dimension(:,:), intent(out) :: skperp, eplus
  complex*16, dimension(:,:), intent(out) :: crfpol

  integer, intent(out) :: ierr    ! status code (0=OK)

  !------------------------
  integer :: iant, iphi, nticrf
  !------------------------

  if(.not.ready) then
     write(lunmsg,*) ' ?rfxqlo_wavefield_get1: RFXQLO module not ready.'
     write(lunmsg,*) '  specify all inputs; call rfxqlo_input_verify first.'
     ierr=1
  else
     ierr=0

     iant = selected_antenna
     iphi = selected_nphi

     nticrf=rfx_grid%nchi_in

     call ckdims('skperp',skperp,ierr)
     if(ierr.ne.0) return

     call ckdims('eplus',eplus,ierr)
     if(ierr.ne.0) return

     call cckdims('crfpol',crfpol,ierr)
     if(ierr.ne.0) return

     nstate = wavefields(iphi,iant)%nstate
     skperp = wavefields(iphi,iant)%skperp
     eplus = wavefields(iphi,iant)%eplus
     crfpol = wavefields(iphi,iant)%crfpol

  endif

  contains

    subroutine ckdims(lbl,data,ierr)

      character*(*) :: lbl   ! name of data
      real*8, dimension(:,:) :: data
      integer, intent(out) :: ierr

      ierr=0

      if(size(data,1).ne.nx_qlo) then
         write(lunmsg,*) ' ?rfxqlo_wavefield_get1(ckdims): ',trim(lbl), &
              ' 1st (radial) dimension mismatch:'
         write(lunmsg,*) '  expected: ',nx_qlo,'  got: ',size(data,1)
         ierr=1
      endif

      if(size(data,2).ne.nticrf+1) then
         write(lunmsg,*) ' ?rfxqlo_wavefield_get1(ckdims): ',trim(lbl), &
              ' 2nd (poloidal) dimension mismatch:'
         write(lunmsg,*) '  expected: ',nticrf+1,'  got: ',size(data,2)
         ierr=1
      endif

    end subroutine ckdims

    subroutine cckdims(lbl,data,ierr)

      character*(*) :: lbl   ! name of data
      complex*16, dimension(:,:) :: data
      integer, intent(out) :: ierr

      ierr=0

      if(size(data,1).ne.nx_qlo) then
         write(lunmsg,*) ' ?rfxqlo_wavefield_get1(cckdims): ',trim(lbl), &
              ' 1st (radial) dimension mismatch:'
         write(lunmsg,*) '  expected: ',nx_qlo,'  got: ',size(data,1)
         ierr=1
      endif

      if(size(data,2).ne.nticrf+1) then
         write(lunmsg,*) ' ?rfxqlo_wavefield_get1(cckdims): ',trim(lbl), &
              ' 2nd (poloidal) dimension mismatch:'
         write(lunmsg,*) '  expected: ',nticrf+1,'  got: ',size(data,2)
         ierr=1
      endif

    end subroutine cckdims
end subroutine rfxqlo_wavefield_get1

!-----------------------------------------------------------------------
subroutine rfxqlo_wavefields_alloc(ierr)

  !  initialize & allocate wavefield structures and contained arrays...
  !  INTERNAL rfxqlo call, normally...

  use rfxqlo_private
  implicit NONE

  integer, intent(out) :: ierr

  !-----------------
  integer :: iant,iphi,nticrf
  !-----------------

  if(.not.ready) then
     write(lunmsg,*) ' ?rfxqlo_wavefields_init: RFXQLO module not ready.'
     write(lunmsg,*) '  specify all inputs; call rfxqlo_input_verify first.'
     ierr=1
  else
     ierr=0
     nticrf = rfx_grid%nchi_in
     call rfxqlo_wavefields_free
     allocate(wavefields(max_nphi,n_antenna))
     do iant=1,n_antenna
        do iphi=1,max_nphi
           wavefields(iphi,iant)%nstate = 0
           wavefields(iphi,iant)%mode_power = 0.0d0
           wavefields(iphi,iant)%mode_coupling = 1.0d0
           allocate(wavefields(iphi,iant)%skperp(nx_qlo,nticrf+1))
           allocate(wavefields(iphi,iant)%eplus(nx_qlo,nticrf+1))
           allocate(wavefields(iphi,iant)%crfpol(nx_qlo,nticrf+1))
        enddo
     enddo
  endif

end subroutine rfxqlo_wavefields_alloc

subroutine rfxqlo_wavefields_free

  !  conditionally deallocate wavefield structures
  !  INTERNAL rfxqlo call, normally...

  use rfxqlo_private
  implicit NONE

  integer :: idim1,idim2,i1,i2

  if(allocated(wavefields)) then

     idim1 = size(wavefields,1)
     idim2 = size(wavefields,2)

     do i2=1,idim2
        do i1=1,idim1

           deallocate(wavefields(i1,i2)%skperp)
           deallocate(wavefields(i1,i2)%eplus)
           deallocate(wavefields(i1,i2)%crfpol)

        enddo
     enddo

     deallocate(wavefields)
  endif

end subroutine rfxqlo_wavefields_free
