module strong_fast_global_mod
!$$$   module documentation block
!                .      .    .                                       .
! module:    strong_fast_global_mod
!
! abstract: Contains all routines for fast strong balance constraint
!

! program history log:
!   2008-04-04  safford - add moddule doc block and missing subroutine doc blocks
!
! subroutines included:
!    init_strongvars_2        --
!    initialize_strong_fast_global
!    strong_bal_correction_fast_global
!    strong_bal_correction_fast_global_ad -- adjoint of strong_bal_correction
!    gather_rmstends0         -- get BAL diagnostics
!    gather_rmstends          -- get BAL diagnostics
!    inmi_coupler_sd2ew0      --
!    inmi_coupler_sd2ew1      --
!    inmi_coupler_sd2ew       --
!    inmi_coupler_ew2sd1      --
!    inmi_coupler_ew2sd       --
!    inmi_ew_trans            --
!    inmi_ew_invtrans_ad      --
!    inmi_ew_invtrans         --
!    inmi_ew_trans_ad         --
!    inmi_coupler_ew2ns0      --
!    inmi_coupler_ew2ns1      --
!    inmi_coupler_ew2ns       --
!    inmi_coupler_ns2ew       --
!    inmi_nsuvm2zdm           --
!    spdz2uv_ns               -- compute winds from div and vort for 1 zonal wave number
!    spuv2dz_ns               -- compute div,vort from winds for one zonal wave number
!    spanaly_ns               -- spanaly modified for one zonal wave number
!    spsynth_ns               -- spsynth modified for one zonal wave number
!
! variable definitions:
!
! attributes:
!   language:  f90
!   machine:   ibm RS/6000 SP
!$$$ end documentation block

  use kinds, only: i_kind
  implicit none

! set default to private
  private
! set subroutines to public
  public :: init_strongvars_2
  public :: initialize_strong_fast_global
  public :: strong_bal_correction_fast_global
  public :: strong_bal_correction_fast_global_ad
  public :: gather_rmstends0
  public :: gather_rmstends
  public :: inmi_coupler_sd2ew0
  public :: inmi_coupler_sd2ew1
  public :: inmi_coupler_sd2ew
  public :: inmi_coupler_ew2sd1
  public :: inmi_coupler_ew2sd
  public :: inmi_ew_trans
  public :: inmi_ew_invtrans_ad
  public :: inmi_ew_invtrans
  public :: inmi_ew_trans_ad
  public :: inmi_coupler_ew2ns0
  public :: inmi_coupler_ew2ns1
  public :: inmi_coupler_ew2ns
  public :: inmi_coupler_ns2ew
  public :: inmi_nsuvm2zdm
  public :: spdz2uv_ns
  public :: spuv2dz_ns
  public :: spanaly_ns
  public :: spsynth_ns

  integer(i_kind),allocatable:: mode_list(:,:)
                                                   ! mode_list(1,j) = lat index for ew strip j
                                                   !  mode_list(2,j) = vert mode number for ew strip j
                                                   !  mode_list(3,j) = pe of this lat/vert mode strip 
  integer(i_kind),allocatable:: mmode_list(:,:)
                                                  !  mmode_list(1,j) = m1 (zonal wave number 1) for ns strip
                                                  !  mmode_list(2,j) = m2 (zonal wave number 2) for ns strip
                                                  !  mmode_list(3,j) = vert mode number 1 for ns strip j
                                                  !  mmode_list(4,j) = vert mode number 2 for ns strip j
                                                  !  mmode_list(5,j) = pe for ns strip j
  integer(i_kind) nlatm_0,nlatm_1,m_0,m_1
  integer(i_kind) mthis
  integer(i_kind),allocatable:: mthis0(:),ndisp(:),indexglob(:)
  integer(i_kind),allocatable,dimension(:)::nsend_sd2ew,nrecv_sd2ew
  integer(i_kind),allocatable,dimension(:)::ndsend_sd2ew,ndrecv_sd2ew
  integer(i_kind),allocatable::info_send_sd2ew(:,:),info_recv_sd2ew(:,:)
  integer(i_kind) nallsend_sd2ew,nallrecv_sd2ew
  integer(i_kind),allocatable,dimension(:)::nsend_ew2sd,nrecv_ew2sd
  integer(i_kind),allocatable,dimension(:)::ndsend_ew2sd,ndrecv_ew2sd
  integer(i_kind),allocatable::info_send_ew2sd(:,:),info_recv_ew2sd(:,:)
  integer(i_kind) nallsend_ew2sd,nallrecv_ew2sd
  integer(i_kind) nallsend,nallrecv
  integer(i_kind),allocatable::info_send(:,:),info_recv(:,:)
  integer(i_kind),allocatable::nsend(:),nrecv(:),ndsend(:),ndrecv(:)

contains

  subroutine init_strongvars_2(mype)

!$$$  subprogram documentation block
!                .      .    .
! subprogram:    init_strongvars_2
!
!   prgrmmr:
!
! abstract:
!
! program history log:
!   2008-04-04  safford -- add subprogram doc block
!
!   input argument list:
!     mype     - mpi task id
!
!   output argument list:
!
! attributes:
!   language:  f90
!   machine:   ibm RS/6000 SP
!
!$$$
    implicit none

    integer(i_kind),intent(in   ) :: mype

    call initialize_strong_fast_global(mype)

  end subroutine init_strongvars_2

subroutine initialize_strong_fast_global(mype)
!$$$  subprogram documentation block
!                .      .    .
! subprogram:    initialize_strong_fast_global
!
!   prgrmmr:
!
! abstract:
!
! program history log:
!   2008-04-04  safford -- add subprogram doc block
!   2010-03-31  treadon - replace specmod jcap with sp_a structure
!
!   input argument list:
!     mype     - mpi task id
!
!   output argument list:
!
! attributes:
!   language:  f90
!   machine:   ibm RS/6000 SP
!
!$$$

  use gridmod, only: nlat,sp_a
  use mod_vtrans, only: nvmodes_keep
  implicit none

  integer(i_kind),intent(in   ) :: mype

  allocate(mode_list(3,nlat*nvmodes_keep),mmode_list(5,(sp_a%jcap+1)*nvmodes_keep))
  call inmi_coupler_sd2ew0(mype)
  call inmi_coupler_ew2ns0(mype)
  call gather_rmstends0
  call inmi_coupler_sd2ew1(mype)
  call inmi_coupler_ew2ns1(mype)
  call inmi_coupler_ew2sd1(mype)

end subroutine initialize_strong_fast_global


subroutine strong_bal_correction_fast_global(u_t,v_t,t_t,ps_t,mype,psi,chi,t,ps, &
                                              bal_diagnostic,fullfield,update)

!$$$  subprogram documentation block
!                .      .    .                                       .
! subprogram:    strong_bal_correction  strong balance correction
!   prgmmr: parrish          org: np23                date: 2006-07-15
!
! abstract: given input perturbation tendencies of u,v,t,ps from TLM on gaussian grid,
!           and input perturbation u,v,t,ps, compute balance adjustment to u,v,t,ps
!           which zeroes out input gravity component of perturbation tendencies.
!           also output, for later use, input tendencies projected onto gravity modes.
!
! program history log:
!   2006-07-15  parrish
!   2007-04-16  kleist  - modified for full field or incremental diagnostics
!   2008-10-08  parrish/derber - modify to output streamfunction and vel. pot. and not update time derivatives
!   2009-11-27  parrish - add uv_hyb_ens.  if uv_hyb_ens=true, then
!                          input/output variables psi=u, chi=v.
!   2010-03-31  treadon - replace specmod jcap with sp_a structure
!
!   input argument list:
!     u_t      - input perturbation u tendency on gaussian grid (subdomains)
!     v_t      - input perturbation v tendency on gaussian grid (subdomains)
!     t_t      - input perturbation T tendency on gaussian grid (subdomains)
!     ps_t     - input perturbation surface pressure tendency on gaussian grid (subdomains)
!     mype     - current processor
!     psi      - input perturbation psi on gaussian grid (subdomains)
!     chi      - input perturbation chi on gaussian grid (subdomains)
!     t        - input perturbation T on gaussian grid (subdomains)
!     ps       - input perturbation surface pressure on gaussian grid (subdomains)
!     bal_diagnostic - if true, then compute BAL diagnostic, a measure of amplitude
!                      of balanced gravity mode tendencies
!     fullfield - if true, full field diagnostics
!                 if false, incremental 
!     update   - if false, then do not update u,v,t,ps with balance increment
!
!   output argument list:
!     psi      - output balanced perturbation u on gaussian grid (subdomains)
!     chi      - output balanced perturbation v on gaussian grid (subdomains)
!     t        - output balanced perturbation T on gaussian grid (subdomains)
!     ps       - output balanced perturbation surface pressure on gaussian grid (subdomains)
!
! attributes:
!   language: f90
!   machine:  ibm RS/6000 SP
!
!$$$

  use kinds, only: r_kind
  use mod_vtrans, only: depths,nvmodes_keep,vtrans,vtrans_inv
  use mod_inmi, only: m,gspeed,mmax,dinmi,gproj
  use gridmod, only: nlat,nlon,lat2,lon2,nsig,sp_a
  use constants, only: zero,one,rearth
  use hybrid_ensemble_parameters, only: uv_hyb_ens
  implicit none

  integer(i_kind)                       ,intent(in   ) :: mype
  logical                               ,intent(in   ) :: bal_diagnostic,update,fullfield
  real(r_kind),dimension(lat2,lon2,nsig),intent(in   ) :: u_t,v_t,t_t
  real(r_kind),dimension(lat2,lon2)     ,intent(in   ) :: ps_t
  real(r_kind),dimension(lat2,lon2,nsig),intent(inout) :: psi,chi,t
  real(r_kind),dimension(lat2,lon2)     ,intent(inout) :: ps

  real(r_kind),dimension(nvmodes_keep)::rmstend_uf,rmstend_g_uf
  real(r_kind),dimension(nvmodes_keep)::rmstend_f,rmstend_g_f

  real(r_kind),dimension(lat2,lon2,nsig)::delu,delv,delt
  real(r_kind),dimension(lat2,lon2)::delps
  real(r_kind),dimension(lat2,lon2,nvmodes_keep)::utilde_t,vtilde_t,mtilde_t
  real(r_kind),dimension(lat2,lon2,nvmodes_keep)::delutilde,delvtilde,delmtilde
  real(r_kind),dimension(lat2,lon2,nvmodes_keep)::utilde_t_g,vtilde_t_g,mtilde_t_g
  real(r_kind),allocatable,dimension(:,:)::rmstend_loc_uf,rmstend_g_loc_uf
  real(r_kind),allocatable,dimension(:,:)::rmstend_loc_f,rmstend_g_loc_f
  real(r_kind),dimension(2,0:sp_a%jcap)::divhat,vorthat,mhat,deldivhat,delvorthat,delmhat
  real(r_kind) rmstend_all_uf,rmstend_all_g_uf,rmstend_all_f,rmstend_all_g_f
  real(r_kind),allocatable,dimension(:,:,:,:)::uvm_ew
  real(r_kind),allocatable,dimension(:,:,:,:,:)::uvm_ewtrans,uvm_ns,zdm_hat
  real(r_kind) del2inv

  integer(i_kind) i,ipair,j,k,kk,mode,n
  logical filtered

  filtered=.true.

  mmax=sp_a%jcap

!   1.  u,v,t,ps   -->    utilde,vtilde,mtilde  (vertical mode transform)
!       (subdomains)         (subdomains)

  call vtrans(u_t,v_t,t_t,ps_t,utilde_t,vtilde_t,mtilde_t)

  allocate(uvm_ew(2,3,nlon,nlatm_0:nlatm_1),uvm_ewtrans(2,3,2,0:sp_a%jcap,nlatm_0:nlatm_1))
  allocate(uvm_ns(3,2,nlat,2,m_0:m_1),zdm_hat(3,2,nlat,2,m_0:m_1))
  allocate(rmstend_loc_uf(2,m_0:m_1))
  allocate(rmstend_g_loc_uf(2,m_0:m_1))
  allocate(rmstend_loc_f(2,m_0:m_1))
  allocate(rmstend_g_loc_f(2,m_0:m_1))
  call inmi_coupler_sd2ew(utilde_t,vtilde_t,mtilde_t,utilde_t,vtilde_t,mtilde_t, &
                          uvm_ew,mype)
  call inmi_ew_trans(uvm_ew,uvm_ewtrans)
  call inmi_coupler_ew2ns(uvm_ewtrans,uvm_ns)
  call inmi_nsuvm2zdm(uvm_ns,zdm_hat)
  rmstend_loc_uf=zero
  rmstend_g_loc_uf=zero
  rmstend_loc_f=zero
  rmstend_g_loc_f=zero
  do kk=m_0,m_1
     do ipair=1,2
        m=mmode_list(ipair,kk)
        mode=mmode_list(ipair+2,kk)
        gspeed=sqrt(depths(abs(mode)))
        i=0
        do n=m,sp_a%jcap
           i=i+1
           vorthat(1,n)=zdm_hat(1,1,i,ipair,kk)
           vorthat(2,n)=zdm_hat(1,2,i,ipair,kk)
           divhat( 1,n)=zdm_hat(2,1,i,ipair,kk)
           divhat( 2,n)=zdm_hat(2,2,i,ipair,kk)
           mhat(   1,n)=zdm_hat(3,1,i,ipair,kk)
           mhat(   2,n)=zdm_hat(3,2,i,ipair,kk)
        end do

!   4.  divhat,vorthat,mhat --> deldivhat,delvorthat,delmhat   (inmi correction)
!          (slabs)                        (slabs)

        if(mode >  0) then
!               here, delvorthat, etc contain field corrections necessary to zero out gravity component
!                                         of tendencies
           call dinmi(vorthat(1,m),divhat(1,m),mhat(1,m),delvorthat(1,m),deldivhat(1,m),delmhat(1,m))
        else
!               here, delvorthat, etc contain gravity component of tendencies
           if(bal_diagnostic) &
              call gproj(vorthat(1,m),divhat(1,m),mhat(1,m),delvorthat(1,m),deldivhat(1,m),delmhat(1,m), &
                        rmstend_loc_uf(ipair,kk),rmstend_g_loc_uf(ipair,kk),.not.filtered)
              call gproj(vorthat(1,m),divhat(1,m),mhat(1,m),delvorthat(1,m),deldivhat(1,m),delmhat(1,m), &
                        rmstend_loc_f(ipair,kk),rmstend_g_loc_f(ipair,kk),filtered)
        end if

        if(uv_hyb_ens) then
           i=0
           do n=m,sp_a%jcap
              i=i+1
              zdm_hat(1,1,i,ipair,kk)=delvorthat(1,n)
              zdm_hat(1,2,i,ipair,kk)=delvorthat(2,n)
              zdm_hat(2,1,i,ipair,kk)=deldivhat(1,n)
              zdm_hat(2,2,i,ipair,kk)=deldivhat(2,n)
              zdm_hat(3,1,i,ipair,kk)=delmhat(1,n)
              zdm_hat(3,2,i,ipair,kk)=delmhat(2,n)
           end do
        else
           i=0
           do n=m,sp_a%jcap
              del2inv=zero
              if(n >  0) del2inv=rearth**2/(n*(n+one))
              i=i+1
              zdm_hat(1,1,i,ipair,kk)=-delvorthat(1,n)*del2inv
              zdm_hat(1,2,i,ipair,kk)=-delvorthat(2,n)*del2inv
              zdm_hat(2,1,i,ipair,kk)=-deldivhat(1,n)*del2inv
              zdm_hat(2,2,i,ipair,kk)=-deldivhat(2,n)*del2inv
              zdm_hat(3,1,i,ipair,kk)=delmhat(1,n)
              zdm_hat(3,2,i,ipair,kk)=delmhat(2,n)
           end do
        end if

     end do
  end do

  if(uv_hyb_ens) then
     call inmi_nszdm2uvm(uvm_ns,zdm_hat)
  else
     call inmi_nspcm_hat2pcm(uvm_ns,zdm_hat)
  end if
  call inmi_coupler_ns2ew(uvm_ewtrans,uvm_ns)
  call inmi_ew_invtrans(uvm_ew,uvm_ewtrans)
  call inmi_coupler_ew2sd(delutilde,delvtilde,delmtilde,utilde_t_g,vtilde_t_g,mtilde_t_g,uvm_ew,mype)
  deallocate(uvm_ew,uvm_ewtrans,uvm_ns,zdm_hat)

  if(bal_diagnostic) then
     call gather_rmstends(rmstend_loc_uf,  rmstend_uf)
     call gather_rmstends(rmstend_g_loc_uf,rmstend_g_uf)
     call gather_rmstends(rmstend_loc_f,   rmstend_f)
     call gather_rmstends(rmstend_g_loc_f, rmstend_g_f)
     if(mype == 0) then
        rmstend_all_uf=zero
        rmstend_all_g_uf=zero

        if (fullfield) then
           write(6,*) 'STRONG_FAST_GLOBAL:   FULL FIELD BALANCE DIAGNOSTICS --  '
        else
           write(6,*) 'STRONG_FAST_GLOBAL:   INCREMENTAL BALANCE DIAGNOSTICS --  '
        end if

        do i=1,nvmodes_keep
           rmstend_all_uf=rmstend_all_uf+rmstend_uf(i)
           rmstend_all_g_uf=rmstend_all_g_uf+rmstend_g_uf(i)
           write(6,'(" mode,rmstend_uf,rmstend_g_uf,rat = ",i5,2e14.4,2f10.4)') &
                              i,rmstend_uf(i),rmstend_g_uf(i),&
                              rmstend_g_uf(i)/(rmstend_uf(i)-rmstend_g_uf(i)), &
                              rmstend_g_uf(i)/(rmstend_uf(1)-rmstend_g_uf(1))
        end do
        rmstend_all_f=zero
        rmstend_all_g_f=zero
        do i=1,nvmodes_keep
           rmstend_all_f=rmstend_all_f+rmstend_f(i)
           rmstend_all_g_f=rmstend_all_g_f+rmstend_g_f(i)
           write(6,'(" mode,rmstend_f,rmstend_g_f,rat = ",i5,2e14.4,2f10.4)') &
                              i,rmstend_f(i),rmstend_g_f(i),&
                              rmstend_g_f(i)/(rmstend_f(i)-rmstend_g_f(i)), &
                              rmstend_g_f(i)/(rmstend_f(1)-rmstend_g_f(1))
        end do
        write(6,'(" rmstend_all_uf,g_uf,rat = ",2e14.4,f10.4)') rmstend_all_uf,rmstend_all_g_uf, &
                                                 rmstend_all_g_uf/(rmstend_all_uf-rmstend_all_g_uf) 
        write(6,'(" rmstend_all_f,g_f,rat = ",2e14.4,f10.4)') rmstend_all_f,rmstend_all_g_f, &
                                                 rmstend_all_g_f/(rmstend_all_f-rmstend_all_g_f) 
     end if
  end if
  deallocate(rmstend_loc_uf,rmstend_g_loc_uf,rmstend_loc_f,rmstend_g_loc_f)

!   7.  delutilde,delvtilde,delmtilde  -->  delu,delv,delt,delps   (vertical mode inverse transform)
!       (subdomains)                      (subdomains)

  delu=zero
  delv=zero
  delt=zero
  delps=zero
  call vtrans_inv(delutilde,delvtilde,delmtilde,delu,delv,delt,delps)
!  call vtrans_inv(utilde_t_g,vtilde_t_g,mtilde_t_g,u_t_g,v_t_g,t_t_g,ps_t_g)


!  update u,v,t,ps


  if(update) then
     do k=1,nsig
        do j=1,lon2
           do i=1,lat2
              psi(i,j,k)=psi(i,j,k)+delu(i,j,k)
              chi(i,j,k)=chi(i,j,k)+delv(i,j,k)
              t(i,j,k)=t(i,j,k)+delt(i,j,k)
           end do
        end do
     end do
     do j=1,lon2
        do i=1,lat2
           ps(i,j)=ps(i,j)+delps(i,j)
        end do
     end do
  end if

end subroutine strong_bal_correction_fast_global

subroutine strong_bal_correction_fast_global_ad(u_t,v_t,t_t,ps_t,mype,psi,chi,t,ps)

!$$$  subprogram documentation block
!                .      .    .                                       .
! subprogram:    strong_bal_correction_ad adjoint of strong_bal_correction
!   prgmmr: parrish          org: np23                date: 2006-07-15
!
! abstract: adjoint of strong_bal_correction
!
! program history log:
!   2006-07-15  parrish
!   2008-10-08  parrish/derber - modify to output streamfunction and vel. pot. and not update time derivatives
!   2009-11-27  parrish - add uv_hyb_ens.  if present and true, then
!                          input/output variables psi=u, chi=v.
!   2010-03-31  treadon - replace specmod jcap with sp_a structure
!
!   input argument list:
!     u_t      - input perturbation u tendency on gaussian grid (subdomains)
!     v_t      - input perturbation v tendency on gaussian grid (subdomains)
!     t_t      - input perturbation T tendency on gaussian grid (subdomains)
!     ps_t     - input perturbation surface pressure tendency on gaussian grid (subdomains)
!     mype     - current processor
!     psi      - input perturbation psi on gaussian grid (subdomains)
!     chi      - input perturbation chi on gaussian grid (subdomains)
!     t        - input perturbation T on gaussian grid (subdomains)
!     ps       - input perturbation surface pressure on gaussian grid (subdomains)
!
!   output argument list:
!     u_t      - output perturbation u tendency on gaussian grid (subdomains)
!     v_t      - output perturbation v tendency on gaussian grid (subdomains)
!     t_t      - output perturbation T tendency on gaussian grid (subdomains)
!     ps_t     - output perturbation surface pressure tendency on gaussian grid (subdomains)
!
! attributes:
!   language: f90
!   machine:  ibm RS/6000 SP
!
!$$$

  use kinds, only: r_kind
  use mod_vtrans, only: depths,nvmodes_keep,vtrans_ad,vtrans_inv_ad
  use mod_inmi, only: m,gspeed,mmax,dinmi_ad,gproj_ad
  use gridmod, only: nlat,nlon,lat2,lon2,nsig,sp_a
  use constants, only: zero,one,rearth
  use hybrid_ensemble_parameters, only: uv_hyb_ens
  implicit none

  integer(i_kind)                       ,intent(in   ) :: mype
  real(r_kind),dimension(lat2,lon2,nsig),intent(inout) :: u_t,v_t,t_t
  real(r_kind),dimension(lat2,lon2)     ,intent(inout) :: ps_t
  real(r_kind),dimension(lat2,lon2,nsig),intent(in   ) :: psi,chi,t
  real(r_kind),dimension(lat2,lon2)     ,intent(in   ) :: ps

  real(r_kind),dimension(lat2,lon2,nsig)::delu,delv,delt
  real(r_kind),dimension(lat2,lon2)::delps
  real(r_kind),dimension(lat2,lon2,nvmodes_keep)::utilde_t,vtilde_t,mtilde_t
  real(r_kind),dimension(lat2,lon2,nvmodes_keep)::utilde_t2,vtilde_t2,mtilde_t2
  real(r_kind),dimension(lat2,lon2,nvmodes_keep)::delutilde,delvtilde,delmtilde
  real(r_kind),dimension(lat2,lon2,nvmodes_keep)::utilde_t_g,vtilde_t_g,mtilde_t_g
  real(r_kind),dimension(2,0:sp_a%jcap)::divhat,vorthat,mhat,deldivhat,delvorthat,delmhat
  real(r_kind),allocatable,dimension(:,:,:,:)::uvm_ew
  real(r_kind),allocatable,dimension(:,:,:,:,:)::uvm_ewtrans,uvm_ns,zdm_hat
  real(r_kind) del2inv

  integer(i_kind) i,ipair,j,k,kk,mode,n

  mmax=sp_a%jcap

!  adjoint of update u,v,t,ps

  do j=1,lon2
     do i=1,lat2
        delps(i,j)=ps(i,j)
     end do
  end do

  do k=1,nsig
     do j=1,lon2
        do i=1,lat2
           delu(i,j,k)=psi(i,j,k)
           delv(i,j,k)=chi(i,j,k)
           delt(i,j,k)=t(i,j,k)
        end do
     end do
  end do

!   7.  adjoint of delutilde,delvtilde,delmtilde  -->  delu,delv,delt,delps  (vert mode inverse transform)
!       (subdomains)                      (subdomains)

  delutilde=zero ; delvtilde=zero ; delmtilde=zero
  utilde_t_g=zero ; vtilde_t_g=zero ; mtilde_t_g=zero
  call vtrans_inv_ad(delutilde,delvtilde,delmtilde,delu,delv,delt,delps)

  allocate(uvm_ew(2,3,nlon,nlatm_0:nlatm_1),uvm_ewtrans(2,3,2,0:sp_a%jcap,nlatm_0:nlatm_1))
  allocate(uvm_ns(3,2,nlat,2,m_0:m_1),zdm_hat(3,2,nlat,2,m_0:m_1))
  call inmi_coupler_sd2ew(delutilde,delvtilde,delmtilde,utilde_t_g,vtilde_t_g,mtilde_t_g,uvm_ew,mype)
  call inmi_ew_invtrans_ad(uvm_ew,uvm_ewtrans)
  call inmi_coupler_ew2ns(uvm_ewtrans,uvm_ns)
  if(uv_hyb_ens) then
     call inmi_nszdm2uvm_ad(uvm_ns,zdm_hat)
  else
     call inmi_nspcm_hat2pcm_ad(uvm_ns,zdm_hat)
  end if
  do kk=m_0,m_1
     do ipair=1,2
        m=mmode_list(ipair,kk)
        mode=mmode_list(ipair+2,kk)
        gspeed=sqrt(depths(abs(mode)))
        do n=m,sp_a%jcap
           vorthat(1,n)=zero
           vorthat(2,n)=zero
           divhat( 1,n)=zero
           divhat( 2,n)=zero
           mhat(   1,n)=zero
           mhat(   2,n)=zero
        end do
        if(uv_hyb_ens) then
           i=0
           do n=m,sp_a%jcap
              i=i+1
              delvorthat(1,n)=zdm_hat(1,1,i,ipair,kk)
              delvorthat(2,n)=zdm_hat(1,2,i,ipair,kk)
              deldivhat(1,n)=zdm_hat(2,1,i,ipair,kk)
              deldivhat(2,n)=zdm_hat(2,2,i,ipair,kk)
              delmhat(1,n)=zdm_hat(3,1,i,ipair,kk)
              delmhat(2,n)=zdm_hat(3,2,i,ipair,kk)
           end do
        else
           i=0
           do n=m,sp_a%jcap
              del2inv=zero
              if(n >  0) del2inv=rearth**2/(n*(n+one))
              i=i+1
              delvorthat(1,n)=-zdm_hat(1,1,i,ipair,kk)*del2inv
              delvorthat(2,n)=-zdm_hat(1,2,i,ipair,kk)*del2inv
              deldivhat(1,n)=-zdm_hat(2,1,i,ipair,kk)*del2inv
              deldivhat(2,n)=-zdm_hat(2,2,i,ipair,kk)*del2inv
              delmhat(1,n)=zdm_hat(3,1,i,ipair,kk)
              delmhat(2,n)=zdm_hat(3,2,i,ipair,kk)
           end do
        end if
        if(mode >  0) then
           call dinmi_ad(vorthat(1,m),divhat(1,m),mhat(1,m),&
                       delvorthat(1,m)   ,   deldivhat(1,m),   delmhat(1,m))
        else
           call gproj_ad(vorthat(1,m),divhat(1,m),mhat(1,m),delvorthat(1,m),deldivhat(1,m),delmhat(1,m))
        end if

        i=0
        do n=m,sp_a%jcap
           i=i+1
           zdm_hat(1,1,i,ipair,kk)=vorthat(1,n)
           zdm_hat(1,2,i,ipair,kk)=vorthat(2,n)
           zdm_hat(2,1,i,ipair,kk)=divhat(1,n)
           zdm_hat(2,2,i,ipair,kk)=divhat(2,n)
           zdm_hat(3,1,i,ipair,kk)=mhat(1,n)
           zdm_hat(3,2,i,ipair,kk)=mhat(2,n)
        end do

     end do
  end do

  call inmi_nsuvm2zdm_ad(uvm_ns,zdm_hat)
  call inmi_coupler_ns2ew(uvm_ewtrans,uvm_ns)
  call inmi_ew_trans_ad(uvm_ew,uvm_ewtrans)
  call inmi_coupler_ew2sd(utilde_t,vtilde_t,mtilde_t,utilde_t2,vtilde_t2,mtilde_t2,uvm_ew,mype)

  utilde_t=utilde_t+utilde_t2
  vtilde_t=vtilde_t+vtilde_t2
  mtilde_t=mtilde_t+mtilde_t2
!
!    1.  adjoint of u,v,t,ps   -->    utilde,vtilde,mtilde  (vertical mode transform)
!                      (subdomains)         (subdomains)

  call vtrans_ad(u_t,v_t,t_t,ps_t,utilde_t,vtilde_t,mtilde_t)
  deallocate(uvm_ew,uvm_ewtrans,uvm_ns,zdm_hat)

end subroutine strong_bal_correction_fast_global_ad

subroutine gather_rmstends0

!$$$  subprogram documentation block
!                .      .    .                                       .
! subprogram:    gather_rmstends0  get BAL diagnostics
!   prgmmr: parrish          org: np23                date: 2006-08-03
!
! abstract: compute BAL diagnostics which give amplitude of 
!           gravity projection of energy norm of tendencies
!
! program history log:
!   2006-08-03  parrish
!   2008-04-04  safford  - rm unused uses
!   2010-03-31  treadon - replace specmod jcap with sp_a structure
!
!   input argument list:
!
!   output argument list:
!     rmstend  -  all vertical modes of rmstend assembled across all processors
!
! attributes:
!   language: f90
!   machine:  ibm RS/6000 SP
!
!$$$

  use mpimod, only: ierror,mpi_comm_world,mpi_integer4,npe
  use mod_vtrans, only: nvmodes_keep
  use gridmod, only: sp_a
  implicit none

  integer(i_kind) indexloc(m_0:m_1)
  integer(i_kind) i
  
  allocate(mthis0(npe),ndisp(npe+1),indexglob((sp_a%jcap+1)*nvmodes_keep))
  do i=m_0,m_1
     indexloc(i)=i
  end do
  mthis=m_1-m_0+1
  call mpi_allgather(mthis,1,mpi_integer4,mthis0,1,mpi_integer4,mpi_comm_world,ierror)
  ndisp(1)=0
  do i=2,npe+1
     ndisp(i)=ndisp(i-1)+mthis0(i-1)
  end do
  call mpi_allgatherv(indexloc,mthis,mpi_integer4, &
                      indexglob,mthis0,ndisp,mpi_integer4,mpi_comm_world,ierror)

end subroutine gather_rmstends0

subroutine gather_rmstends(rmstend_loc,rmstend)

!$$$  subprogram documentation block
!                .      .    .                                       .
! subprogram:    gather_rmstends  get BAL diagnostics
!   prgmmr: parrish          org: np23                date: 2006-08-03
!
! abstract: compute BAL diagnostics which give amplitude of 
!           gravity projection of energy norm of tendencies
!
! program history log:
!   2006-08-03  parrish
!   2008-04-04  safford  - rm unused uses
!
!   input argument list:
!     rmstend_loc - previously computed energy norms of vertical modes
!                   as distributed on local processors
!
!   output argument list:
!     rmstend  -  all vertical modes of rmstend assembled across all processors
!
! attributes:
!   language: f90
!   machine:  ibm RS/6000 SP
!
!$$$

  use kinds, only: r_kind
  use mpimod, only: ierror,mpi_comm_world,mpi_rtype,npe
  use mod_vtrans, only: nvmodes_keep
  use constants, only: zero
  use gridmod, only: sp_a
  implicit none

  real(r_kind),intent(in   ) :: rmstend_loc(2,m_0:m_1)
  real(r_kind),intent(  out) :: rmstend(nvmodes_keep)

  integer(i_kind) i,ii,mode,mpi_string1
  real(r_kind) work(2,(sp_a%jcap+1)*nvmodes_keep)
  
  call mpi_type_contiguous(2,mpi_rtype,mpi_string1,ierror)
  call mpi_type_commit(mpi_string1,ierror)
  call mpi_allgatherv(rmstend_loc,mthis,mpi_string1, &
                     work,mthis0,ndisp,mpi_string1,mpi_comm_world,ierror)
  call mpi_type_free(mpi_string1,ierror)
  rmstend=zero
  do i=1,ndisp(npe+1)
     ii=indexglob(i)
     mode=mmode_list(3,ii)
     if(mode <  0) rmstend(-mode)=work(1,ii)+rmstend(-mode)
     mode=mmode_list(4,ii)
     if(mode <  0) rmstend(-mode)=work(2,ii)+rmstend(-mode)
  end do

end subroutine gather_rmstends

subroutine inmi_coupler_sd2ew0(mype)
!$$$  subprogram documentation block
!                .      .    .                                       .
! subprogram:    inmi_coupler_sd2ew0
!   prgmmr:
!
! abstract:  create ew (lat strips) subdivision for use with inmi
!
! program history log:
!   2008-04-04  safford  - add doc block
!
!   input argument list:
!     mype        - current processor number
!
!   output argument list:
!
! attributes:
!   language: f90
!   machine:  ibm RS/6000 SP
!
!$$$

  use mod_vtrans, only: nvmodes_keep
  use gridmod, only: nlat
  use mpimod, only: npe
  implicit none

  integer(i_kind),intent(in   ) :: mype

  integer(i_kind) i,k,kchk,kk,n,nn,nlatm_this

  nlatm_this=nlat*nvmodes_keep/npe
  if(mod(nlat*nvmodes_keep,npe)/=0) nlatm_this=nlatm_this+1
  if(mod(nlat*nvmodes_keep,npe)==0) then
     kchk=npe
  else
     kchk=mod(nlat*nvmodes_keep,npe)
  end if

  nn=0
  do k=1,nvmodes_keep
     do i=1,nlat
        nn=nn+1
        mode_list(1,nn)=i
        mode_list(2,nn)=k
        mode_list(3,nn)=-1
     end do
  end do
  
  nlatm_0=-1
  nlatm_1=-2
  nn=0
  do n=1,npe
     if(n <= kchk) then
        kk=nlatm_this
     else
        kk=nlatm_this-1
     end if
     if(kk >  0) then
        if(mype+1 == n) then
           nlatm_0=nn+1
           nlatm_1=nn+kk
        end if
        do k=1,kk
           nn=nn+1
           mode_list(3,nn)=n
        end do
     end if
  end do

end subroutine inmi_coupler_sd2ew0

subroutine inmi_coupler_sd2ew1(mype)
!$$$  subprogram documentation block
!                .      .    .
! subprogram:    inmi_coupler_sd2ew1
!
!   prgrmmr:
!
! abstract: use mpi_alltoallv to move u_sd,v_sd,m_sd (subdomains) to uvm_ew (lat strips)
!
! program history log:
!   2008-04-04  safford -- add subprogram doc block
!
!   input argument list:
!     mype     - mpi task id
!
!   output argument list:
!
! attributes:
!   language:  f90
!   machine:   ibm RS/6000 SP
!
!$$$

  use mod_vtrans, only: nvmodes_keep
  use gridmod, only: nlat,lon2,lat2,jstart,istart
  use mpimod, only: npe,mpi_comm_world,ierror,mpi_integer4
  implicit none


  integer(i_kind),intent(in   ) :: mype

  integer(i_kind) mode2_list(nlat,nvmodes_keep)
  integer(i_kind) i,ii,ii0,ilat,imode,j,mm1,nn,nlonloc,ipe,ilatm,ilon,mpi_string1

  allocate(nsend_sd2ew(npe),nrecv_sd2ew(npe))
  allocate(ndsend_sd2ew(npe+1),ndrecv_sd2ew(npe+1))
  mm1=mype+1

  nn=0
  mode2_list=0
  do j=1,nlat*nvmodes_keep
     ilat=mode_list(1,j)
     imode=mode_list(2,j)
     if(mode2_list(ilat,imode) /= 0) then
        if(mype == 0) write(6,*)' problem in inmi_coupler_sd2ew'
        call mpi_finalize(i)
        stop
     end if
     mode2_list(ilat,imode)=j
  end do
  do imode=1,nvmodes_keep
     if(imode == 0) cycle
     do ilat=1,nlat
        if(mode2_list(ilat,imode) == 0) then
           if(mype == 0) write(6,*)' problem in inmi_coupler_sd2ew'
           call mpi_finalize(i)
           stop
        end if
     end do
  end do

!  obtain counts of points to send to each pe from this pe

  nsend_sd2ew=0
  nlonloc=lon2-2
  do imode=1,nvmodes_keep
     if(imode == 0) cycle
     do i=2,lat2-1
        ilat=i+istart(mm1)-2
        j=mode2_list(ilat,imode)
        ipe=mode_list(3,j)
        nsend_sd2ew(ipe)=nsend_sd2ew(ipe)+nlonloc
     end do
  end do

  ndsend_sd2ew(1)=0
  do i=2,npe+1
     ndsend_sd2ew(i)=ndsend_sd2ew(i-1)+nsend_sd2ew(i-1)
  end do
  nallsend_sd2ew=ndsend_sd2ew(npe+1)
  allocate(info_send_sd2ew(2,nallsend_sd2ew))
  nsend_sd2ew=0
  do imode=1,nvmodes_keep
     if(imode == 0) cycle
     do i=2,lat2-1
        ilat=i+istart(mm1)-2
        ilatm=mode2_list(ilat,imode)
        ipe=mode_list(3,ilatm)
        do ii=2,lon2-1
           ilon=ii+jstart(mm1)-2
           nsend_sd2ew(ipe)=nsend_sd2ew(ipe)+1
           ii0=ndsend_sd2ew(ipe)+nsend_sd2ew(ipe)
           info_send_sd2ew(1,ii0)=ilon
           info_send_sd2ew(2,ii0)=ilatm
        end do
     end do
  end do

  call mpi_alltoall(nsend_sd2ew,1,mpi_integer4,nrecv_sd2ew,1,mpi_integer4,mpi_comm_world,ierror)
  ndrecv_sd2ew(1)=0
  do i=2,npe+1
     ndrecv_sd2ew(i)=ndrecv_sd2ew(i-1)+nrecv_sd2ew(i-1)
  end do
  nallrecv_sd2ew=ndrecv_sd2ew(npe+1)
  allocate(info_recv_sd2ew(2,nallrecv_sd2ew))
  call mpi_type_contiguous(2,mpi_integer4,mpi_string1,ierror)
  call mpi_type_commit(mpi_string1,ierror)
  call mpi_alltoallv(info_send_sd2ew,nsend_sd2ew,ndsend_sd2ew,mpi_string1, &
                     info_recv_sd2ew,nrecv_sd2ew,ndrecv_sd2ew,mpi_string1,mpi_comm_world,ierror)
  call mpi_type_free(mpi_string1,ierror)
    
end subroutine inmi_coupler_sd2ew1

subroutine inmi_coupler_sd2ew(u_sd1,v_sd1,m_sd1,u_sd2,v_sd2,m_sd2,uvm_ew,mype)

!$$$  subprogram documentation block
!                .      .    .
! subprogram:    init_coupler_sd2ew
!
!   prgrmmr:
!
! abstract: use mpi_alltoallv to move u_sd,v_sd,m_sd (subdomains) to uvm_ew (lat strips)
!
! program history log:
!   2008-04-04  safford -- add subprogram doc block, rm unused vars
!
!   input argument list:
!     mype     - mpi task id
!     u_sd1    -
!     v_sd1    -
!     m_sd1    -
!     u_sd2    -
!     v_sd2    -
!     m_sd2    -
!
!   output argument list:
!     uvm_ew   -
!
! attributes:
!   language:  f90
!   machine:   ibm RS/6000 SP
!
!$$$

  use kinds, only: r_kind
  use mod_vtrans, only: nvmodes_keep
  use gridmod, only: nlon,lon2,lat2,jstart,istart
  use mpimod, only: mpi_comm_world,ierror,mpi_rtype
  implicit none


  integer(i_kind)                                 ,intent(in   ) :: mype
  real(r_kind),dimension(lat2,lon2,nvmodes_keep)  ,intent(in   ) :: u_sd1,v_sd1,m_sd1
  real(r_kind),dimension(lat2,lon2,nvmodes_keep)  ,intent(in   ) :: u_sd2,v_sd2,m_sd2
  real(r_kind),dimension(2,3,nlon,nlatm_0:nlatm_1),intent(  out) :: uvm_ew

  real(r_kind),allocatable::sendbuf(:,:,:),recvbuf(:,:,:)
  integer(i_kind) ilat,imode,j,mm1,ilatm,ilon,mpi_string1

  mm1=mype+1

  allocate(sendbuf(2,3,nallsend_sd2ew))
  do j=1,nallsend_sd2ew
     ilon=info_send_sd2ew(1,j)
     ilatm=info_send_sd2ew(2,j)
     ilat=mode_list(1,ilatm)
     imode=mode_list(2,ilatm)
     sendbuf(1,1,j)=u_sd1(ilat-istart(mm1)+2,ilon-jstart(mm1)+2,imode)
     sendbuf(1,2,j)=v_sd1(ilat-istart(mm1)+2,ilon-jstart(mm1)+2,imode)
     sendbuf(1,3,j)=m_sd1(ilat-istart(mm1)+2,ilon-jstart(mm1)+2,imode)
     sendbuf(2,1,j)=u_sd2(ilat-istart(mm1)+2,ilon-jstart(mm1)+2,imode)
     sendbuf(2,2,j)=v_sd2(ilat-istart(mm1)+2,ilon-jstart(mm1)+2,imode)
     sendbuf(2,3,j)=m_sd2(ilat-istart(mm1)+2,ilon-jstart(mm1)+2,imode)
  end do
  allocate(recvbuf(2,3,nallrecv_sd2ew))
  call mpi_type_contiguous(6,mpi_rtype,mpi_string1,ierror)
  call mpi_type_commit(mpi_string1,ierror)
  call mpi_alltoallv(sendbuf,nsend_sd2ew,ndsend_sd2ew,mpi_string1, &
                     recvbuf,nrecv_sd2ew,ndrecv_sd2ew,mpi_string1,mpi_comm_world,ierror)
  call mpi_type_free(mpi_string1,ierror)
  deallocate(sendbuf)
    
  do j=1,nallrecv_sd2ew
     ilon=info_recv_sd2ew(1,j)
     ilatm=info_recv_sd2ew(2,j)
     uvm_ew(1,1,ilon,ilatm)=recvbuf(1,1,j)
     uvm_ew(1,2,ilon,ilatm)=recvbuf(1,2,j)
     uvm_ew(1,3,ilon,ilatm)=recvbuf(1,3,j)
     uvm_ew(2,1,ilon,ilatm)=recvbuf(2,1,j)
     uvm_ew(2,2,ilon,ilatm)=recvbuf(2,2,j)
     uvm_ew(2,3,ilon,ilatm)=recvbuf(2,3,j)
  end do
  deallocate(recvbuf)
    
end subroutine inmi_coupler_sd2ew

subroutine inmi_coupler_ew2sd1(mype)
!$$$  subprogram documentation block
!                .      .    .
! subprogram:    inmi_coupler_ew2sd1
!
!   prgrmmr:
!
! abstract: use mpi_alltoallv to move uvm_ew (lat strips) to u_sd,v_sd,m_sd (subdomains)
!
! program history log:
!   2008-04-04  safford -- add subprogram doc block, rm unused vars and uses
!
!   input argument list:
!     mype     - mpi task id
!
!   output argument list:
!
! attributes:
!   language:  f90
!   machine:   ibm RS/6000 SP
!
!$$$

  use gridmod, only: nlon,jstart,istart,ilat1,jlon1
  use mpimod, only: npe,mpi_comm_world,ierror,mpi_integer4
  implicit none


  integer(i_kind),intent(in   ) :: mype

  integer(i_kind) i,ilat,imode,j,k,mm1,ipe,ilon,mpi_string1,nn

  allocate(nsend_ew2sd(npe),nrecv_ew2sd(npe))
  allocate(ndsend_ew2sd(npe+1),ndrecv_ew2sd(npe+1))
  mm1=mype+1

!      1.  for each pe, gather up list of points from this set of lat strips destined
!             for subdomain of pe
  do ipe=1,npe
     nn=0
     do k=nlatm_0,nlatm_1
        ilat=mode_list(1,k)
        imode=mode_list(2,k)
        i=ilat-istart(ipe)+2
        if(i <  1.or.i >  ilat1(ipe)+2) cycle
        do j=1,jlon1(ipe)+2
           ilon=j+jstart(ipe)-2
           if(ilon <  1) ilon=ilon+nlon
           if(ilon >  nlon) ilon=ilon-nlon
           nn=nn+1
        end do
     end do
     nsend_ew2sd(ipe)=nn
  end do

  ndsend_ew2sd(1)=0
  do i=2,npe+1
     ndsend_ew2sd(i)=ndsend_ew2sd(i-1)+nsend_ew2sd(i-1)
  end do
  nallsend_ew2sd=ndsend_ew2sd(npe+1)
  allocate(info_send_ew2sd(3,nallsend_ew2sd))
  nn=0
  do ipe=1,npe
     do k=nlatm_0,nlatm_1
        ilat=mode_list(1,k)
        imode=mode_list(2,k)
        i=ilat-istart(ipe)+2
        if(i <  1.or.i >  ilat1(ipe)+2) cycle
        do j=1,jlon1(ipe)+2
           ilon=j+jstart(ipe)-2
           if(ilon <  1) ilon=ilon+nlon
           if(ilon >  nlon) ilon=ilon-nlon
           nn=nn+1
           info_send_ew2sd(1,nn)=ilon
           info_send_ew2sd(2,nn)=j
           info_send_ew2sd(3,nn)=k
        end do
     end do
  end do

  call mpi_alltoall(nsend_ew2sd,1,mpi_integer4,nrecv_ew2sd,1,mpi_integer4,mpi_comm_world,ierror)
  ndrecv_ew2sd(1)=0
  do i=2,npe+1
     ndrecv_ew2sd(i)=ndrecv_ew2sd(i-1)+nrecv_ew2sd(i-1)
  end do
  nallrecv_ew2sd=ndrecv_ew2sd(npe+1)
  allocate(info_recv_ew2sd(3,nallrecv_ew2sd))
  call mpi_type_contiguous(3,mpi_integer4,mpi_string1,ierror)
  call mpi_type_commit(mpi_string1,ierror)
  call mpi_alltoallv(info_send_ew2sd,nsend_ew2sd,ndsend_ew2sd,mpi_string1, &
                     info_recv_ew2sd,nrecv_ew2sd,ndrecv_ew2sd,mpi_string1,mpi_comm_world,ierror)
  call mpi_type_free(mpi_string1,ierror)

end subroutine inmi_coupler_ew2sd1

subroutine inmi_coupler_ew2sd(u_sd1,v_sd1,m_sd1,u_sd2,v_sd2,m_sd2,uvm_ew,mype)
!$$$  subprogram documentation block
!                .      .    .
! subprogram:    inmi_coupler_ew2sd
!
!   prgrmmr:
!
! abstract: use mpi_alltoallv to move uvm_ew (lat strips) to u_sd,v_sd,m_sd (subdomains)
!
! program history log:
!   2008-04-04  safford -- add subprogram doc block, rm unused vars
!
!   input argument list:
!     mype     - mpi task id
!     uvm_ew   -
!
!   output argument list:
!     u_sd1    -
!     v_sd1    -
!     m_sd1    -
!     u_sd2    -
!     v_sd2    -
!     m_sd2    -
!
! attributes:
!   language:  f90
!   machine:   ibm RS/6000 SP
!
!$$$

  use kinds, only: r_kind
  use mod_vtrans, only: nvmodes_keep
  use gridmod, only: nlat,nlon,lon2,lat2,istart
  use mpimod, only: mpi_comm_world,ierror,mpi_rtype
  implicit none


  integer(i_kind)                                 ,intent(in   ) :: mype
  real(r_kind),dimension(lat2,lon2,nvmodes_keep)  ,intent(  out) :: u_sd1,v_sd1,m_sd1
  real(r_kind),dimension(lat2,lon2,nvmodes_keep)  ,intent(  out) :: u_sd2,v_sd2,m_sd2
  real(r_kind),dimension(2,3,nlon,nlatm_0:nlatm_1),intent(in   ) :: uvm_ew

  real(r_kind),allocatable::sendbuf(:,:,:),recvbuf(:,:,:)
  integer(i_kind) ilat,imode,j,mm1,ilatm,ilon,mpi_string1,ilonloc

  mm1=mype+1

  allocate(sendbuf(2,3,nallsend_ew2sd))
  do j=1,nallsend_ew2sd
     ilon=info_send_ew2sd(1,j)
     ilatm=info_send_ew2sd(3,j)
     sendbuf(1,1,j)=uvm_ew(1,1,ilon,ilatm)
     sendbuf(1,2,j)=uvm_ew(1,2,ilon,ilatm)
     sendbuf(1,3,j)=uvm_ew(1,3,ilon,ilatm)
     sendbuf(2,1,j)=uvm_ew(2,1,ilon,ilatm)
     sendbuf(2,2,j)=uvm_ew(2,2,ilon,ilatm)
     sendbuf(2,3,j)=uvm_ew(2,3,ilon,ilatm)
  end do
  allocate(recvbuf(2,3,nallrecv_ew2sd))
  call mpi_type_contiguous(6,mpi_rtype,mpi_string1,ierror)
  call mpi_type_commit(mpi_string1,ierror)
  call mpi_alltoallv(sendbuf,nsend_ew2sd,ndsend_ew2sd,mpi_string1, &
                     recvbuf,nrecv_ew2sd,ndrecv_ew2sd,mpi_string1,mpi_comm_world,ierror)
  call mpi_type_free(mpi_string1,ierror)
  deallocate(sendbuf)
  do j=1,nallrecv_ew2sd
     ilonloc=info_recv_ew2sd(2,j)
     ilatm=info_recv_ew2sd(3,j)
     ilat=mode_list(1,ilatm)
     imode=mode_list(2,ilatm)
     u_sd1(ilat-istart(mm1)+2,ilonloc,imode)=recvbuf(1,1,j)
     v_sd1(ilat-istart(mm1)+2,ilonloc,imode)=recvbuf(1,2,j)
     m_sd1(ilat-istart(mm1)+2,ilonloc,imode)=recvbuf(1,3,j)
     u_sd2(ilat-istart(mm1)+2,ilonloc,imode)=recvbuf(2,1,j)
     v_sd2(ilat-istart(mm1)+2,ilonloc,imode)=recvbuf(2,2,j)
     m_sd2(ilat-istart(mm1)+2,ilonloc,imode)=recvbuf(2,3,j)
!--------------check for north or south pole
     ilat=-1
     if(mode_list(1,ilatm) == nlat) ilat=nlat+1
     if(mode_list(1,ilatm) == 1) ilat=0
     if(ilat == -1) cycle
!-----------------do repeat rows for north/south pole
     u_sd1(ilat-istart(mm1)+2,ilonloc,imode)=recvbuf(1,1,j)
     v_sd1(ilat-istart(mm1)+2,ilonloc,imode)=recvbuf(1,2,j)
     m_sd1(ilat-istart(mm1)+2,ilonloc,imode)=recvbuf(1,3,j)
     u_sd2(ilat-istart(mm1)+2,ilonloc,imode)=recvbuf(2,1,j)
     v_sd2(ilat-istart(mm1)+2,ilonloc,imode)=recvbuf(2,2,j)
     m_sd2(ilat-istart(mm1)+2,ilonloc,imode)=recvbuf(2,3,j)
  end do
  deallocate(recvbuf)

end subroutine inmi_coupler_ew2sd

subroutine inmi_ew_trans(uvm_ew,uvm_ewtrans)

!$$$  subprogram documentation block
!                .      .    .
! subprogram:    inmi_ew_trans
!
!   prgrmmr:
!
! abstract:
!
! program history log:
!   2008-04-04  safford -- add subprogram doc block, rm unused uses
!   2010-03-31  treadon - replace specmod components with sp_a structure
!
!   input argument list:
!     uvm_ew   -
!
!   output argument list:
!     uvm_ewtrans -
!
! attributes:
!   language:  f90
!   machine:   ibm RS/6000 SP
!
!$$$

  use kinds, only: r_kind
  use gridmod, only: nlon,sp_a
  implicit none

  real(r_kind),dimension(2,3,nlon,nlatm_0:nlatm_1)    ,intent(in   ) :: uvm_ew
  real(r_kind),dimension(2,3,2,0:sp_a%jcap,nlatm_0:nlatm_1),intent(  out) :: uvm_ewtrans

  integer(i_kind) i,j,k
  real(r_kind) grid(nlon,2),halfwave(2,0:nlon/2,2)

  do k=nlatm_0,nlatm_1
     do j=1,3
        do i=1,nlon
           grid(i,1)=uvm_ew(1,j,i,k)
           grid(i,2)=uvm_ew(2,j,i,k)
        end do
        call spffte(nlon,1+nlon/2,nlon,2,halfwave,grid,-1,sp_a%afft)
        do i=0,sp_a%jcap
           uvm_ewtrans(1,j,1,i,k)=halfwave(1,i,1)
           uvm_ewtrans(1,j,2,i,k)=halfwave(2,i,1)
           uvm_ewtrans(2,j,1,i,k)=halfwave(1,i,2)
           uvm_ewtrans(2,j,2,i,k)=halfwave(2,i,2)
        end do
     end do
  end do

end subroutine inmi_ew_trans

subroutine inmi_ew_invtrans_ad(uvm_ew,uvm_ewtrans)
!$$$  subprogram documentation block
!                .      .    .
! subprogram:    inmi_ew_invtrans_ad
!
!   prgrmmr:
!
! abstract:
!
! program history log:
!   2008-04-04  safford -- add subprogram doc block, rm unused uses
!   2010-03-31  treadon - replace specmod components with sp_a structure
!
!   input argument list:
!     uvm_ew   -
!
!   output argument list:
!     uvm_ewtrans -
!
! attributes:
!   language:  f90
!   machine:   ibm RS/6000 SP
!
!$$$

  use kinds, only: r_kind
  use gridmod, only: nlon,sp_a
  use constants, only: two
  implicit none

  real(r_kind),dimension(2,3,nlon,nlatm_0:nlatm_1)    ,intent(in   ) :: uvm_ew
  real(r_kind),dimension(2,3,2,0:sp_a%jcap,nlatm_0:nlatm_1),intent(  out) :: uvm_ewtrans

  integer(i_kind) i,j,k
  real(r_kind) grid(nlon,2),halfwave(2,0:nlon/2,2)

  do k=nlatm_0,nlatm_1
     do j=1,3
        do i=1,nlon
           grid(i,1)=uvm_ew(1,j,i,k)
           grid(i,2)=uvm_ew(2,j,i,k)
        end do
        call spffte(nlon,1+nlon/2,nlon,2,halfwave,grid,-1,sp_a%afft)
        uvm_ewtrans(1,j,1,0,k)=halfwave(1,0,1)*float(nlon)
        uvm_ewtrans(1,j,2,0,k)=halfwave(2,0,1)*float(nlon)
        uvm_ewtrans(2,j,1,0,k)=halfwave(1,0,2)*float(nlon)
        uvm_ewtrans(2,j,2,0,k)=halfwave(2,0,2)*float(nlon)
        do i=1,sp_a%jcap
           uvm_ewtrans(1,j,1,i,k)=halfwave(1,i,1)*two*float(nlon)
           uvm_ewtrans(1,j,2,i,k)=halfwave(2,i,1)*two*float(nlon)
           uvm_ewtrans(2,j,1,i,k)=halfwave(1,i,2)*two*float(nlon)
           uvm_ewtrans(2,j,2,i,k)=halfwave(2,i,2)*two*float(nlon)
        end do
     end do
  end do

end subroutine inmi_ew_invtrans_ad

subroutine inmi_ew_invtrans(uvm_ew,uvm_ewtrans)

!$$$  subprogram documentation block
!                .      .    .
! subprogram:    inmi_ew_invtrans
!
!   prgrmmr:
!
! abstract:
!
! program history log:
!   2008-04-04  safford -- add subprogram doc block, rm unused vars
!   2010-03-31  treadon - replace specmod components with sp_a structure
!
!   input argument list:
!     uvm_ewtrans -
!
!   output argument list:
!     uvm_ew      -
!
! attributes:
!   language:  f90
!   machine:   ibm RS/6000 SP
!
!$$$

  use kinds, only: r_kind
  use gridmod, only: nlon,sp_a
  use constants, only: zero
  implicit none

  real(r_kind),dimension(2,3,nlon,nlatm_0:nlatm_1)    ,intent(  out) :: uvm_ew
  real(r_kind),dimension(2,3,2,0:sp_a%jcap,nlatm_0:nlatm_1),intent(in   ) :: uvm_ewtrans

  integer(i_kind) i,j,k
  real(r_kind) grid(nlon,2),halfwave(2,0:nlon/2,2)

  do k=nlatm_0,nlatm_1
     do j=1,3
        do i=0,sp_a%jcap
           halfwave(1,i,1)=uvm_ewtrans(1,j,1,i,k)
           halfwave(2,i,1)=uvm_ewtrans(1,j,2,i,k)
           halfwave(1,i,2)=uvm_ewtrans(2,j,1,i,k)
           halfwave(2,i,2)=uvm_ewtrans(2,j,2,i,k)
        end do
        do i=sp_a%jcap+1,nlon/2
           halfwave(1,i,1)=zero
           halfwave(2,i,1)=zero
           halfwave(1,i,2)=zero
           halfwave(2,i,2)=zero
        end do
        call spffte(nlon,1+nlon/2,nlon,2,halfwave,grid,1,sp_a%afft)
        do i=1,nlon
           uvm_ew(1,j,i,k)=grid(i,1)
           uvm_ew(2,j,i,k)=grid(i,2)
        end do
     end do
  end do

end subroutine inmi_ew_invtrans

subroutine inmi_ew_trans_ad(uvm_ew,uvm_ewtrans)
!$$$  subprogram documentation block
!                .      .    .
! subprogram:    inmi_ew_trans_ad
!
!   prgrmmr:
!
! abstract:
!
! program history log:
!   2008-04-04  safford -- add subprogram doc block, rm unused vars
!   2010-03-31  treadon - replace specmod components with sp_a structure
!
!   input argument list:
!     uvm_ewtrans -
!
!   output argument list:
!   output argument list:
!     uvm_ew      -
!
! attributes:
!   language:  f90
!   machine:   ibm RS/6000 SP
!
!$$$


  use kinds, only: r_kind
  use gridmod, only: nlon,sp_a
  use constants, only: zero,half
  implicit none

  real(r_kind),dimension(2,3,nlon,nlatm_0:nlatm_1)    ,intent(  out) :: uvm_ew
  real(r_kind),dimension(2,3,2,0:sp_a%jcap,nlatm_0:nlatm_1),intent(in   ) :: uvm_ewtrans

  integer(i_kind) i,j,k
  real(r_kind) grid(nlon,2),halfwave(2,0:nlon/2,2)

  do k=nlatm_0,nlatm_1
     do j=1,3
        halfwave(1,0,1)=uvm_ewtrans(1,j,1,0,k)/float(nlon)
        halfwave(2,0,1)=zero
        halfwave(1,0,2)=uvm_ewtrans(2,j,1,0,k)/float(nlon)
        halfwave(2,0,2)=zero
        do i=1,sp_a%jcap
           halfwave(1,i,1)=half*uvm_ewtrans(1,j,1,i,k)/float(nlon)
           halfwave(2,i,1)=half*uvm_ewtrans(1,j,2,i,k)/float(nlon)
           halfwave(1,i,2)=half*uvm_ewtrans(2,j,1,i,k)/float(nlon)
           halfwave(2,i,2)=half*uvm_ewtrans(2,j,2,i,k)/float(nlon)
        end do
        do i=sp_a%jcap+1,nlon/2
           halfwave(1,i,1)=zero
           halfwave(2,i,1)=zero
           halfwave(1,i,2)=zero
           halfwave(2,i,2)=zero
        end do
        call spffte(nlon,1+nlon/2,nlon,2,halfwave,grid,1,sp_a%afft)
        do i=1,nlon
           uvm_ew(1,j,i,k)=grid(i,1)
           uvm_ew(2,j,i,k)=grid(i,2)
        end do
     end do
  end do

end subroutine inmi_ew_trans_ad

subroutine inmi_coupler_ew2ns0(mype)

!$$$  subprogram documentation block
!                .      .    .
! subprogram:    inmi_coupler_ew2ns0
!
!   prgrmmr:
!
! abstract:
!
! program history log:
!   2008-04-04  safford -- add subprogram doc block, rm unused uses
!   2010-03-31  treadon - replace specmod jcap with sp_a structure
!
!   input argument list:
!     mype     - mpi task id
!
!   output argument list:
!
! attributes:
!   language:  f90
!   machine:   ibm RS/6000 SP
!
!$$$

  use mod_vtrans, only: nvmodes_keep
  use mpimod, only: npe
  use gridmod, only: sp_a
  implicit none

  integer(i_kind),intent(in   ) :: mype

  integer(i_kind) k,kk,m,n,num_per_pe,total_groups,nn,kchk

!   in laying out by zonal wave number/vertical mode, have two types of groupings:

!   1.  jcap odd: 

!     group zonal wave numbers in pairs as follows:

!      (0,(jcap+1)/2)    then    ( m,jcap+1-m)  for m=1,(jcap-1)/2

!   2.  jcap even:

!     0  (+,- mode pair together),  then   (m,jcap+1-m ) for m=1,jcap/2  single modes

  total_groups=(sp_a%jcap+1)*nvmodes_keep
  num_per_pe=total_groups/npe
  if(mod(total_groups,npe)/=0) num_per_pe=num_per_pe+1
  if(mod(total_groups,npe)==0) then
     kchk=npe
  else
     kchk=mod(total_groups,npe)
  end if

  if(mod(sp_a%jcap,2) /= 0) then

!    case  jcap odd:
   
     nn=0
     do k=1,nvmodes_keep
        nn=nn+1
        mmode_list(1,nn)=0
        mmode_list(2,nn)=(sp_a%jcap+1)/2
        mmode_list(3,nn)=k
        mmode_list(4,nn)=k
        mmode_list(5,nn)=-1
        do m=1,(sp_a%jcap-1)/2
           nn=nn+1
           mmode_list(1,nn)=m
           mmode_list(2,nn)=sp_a%jcap+1-m
           mmode_list(3,nn)=k
           mmode_list(4,nn)=k
           mmode_list(5,nn)=-1
        end do
     end do
     do k=1,nvmodes_keep
        nn=nn+1
        mmode_list(1,nn)=0
        mmode_list(2,nn)=(sp_a%jcap+1)/2
        mmode_list(3,nn)=-k
        mmode_list(4,nn)=-k
        mmode_list(5,nn)=-1
        do m=1,(sp_a%jcap-1)/2
           nn=nn+1
           mmode_list(1,nn)=m
           mmode_list(2,nn)=sp_a%jcap+1-m
           mmode_list(3,nn)=-k
           mmode_list(4,nn)=-k
           mmode_list(5,nn)=-1
        end do
     end do

  else

!    case  jcap even:
   
     nn=0
     do k=1,nvmodes_keep
        nn=nn+1
        mmode_list(1,nn)=0
        mmode_list(2,nn)=0
        mmode_list(3,nn)=k
        mmode_list(4,nn)=-k
        mmode_list(5,nn)=-1
        do m=1,sp_a%jcap/2
           nn=nn+1
           mmode_list(1,nn)=m
           mmode_list(2,nn)=sp_a%jcap+1-m
           mmode_list(3,nn)=k
           mmode_list(4,nn)=k
           mmode_list(5,nn)=-1
        end do
     end do
     do k=1,nvmodes_keep
        do m=1,sp_a%jcap/2
           nn=nn+1
           mmode_list(1,nn)=m
           mmode_list(2,nn)=sp_a%jcap+1-m
           mmode_list(3,nn)=-k
           mmode_list(4,nn)=-k
           mmode_list(5,nn)=-1
        end do
     end do

  end if

  m_0=-1
  m_1=-2
  nn=0
  do n=1,npe
     if(n <= kchk) then
        kk=num_per_pe
     else
        kk=num_per_pe-1
     end if
     if(kk >  0) then
        if(mype+1 == n) then
           m_0=nn+1
           m_1=nn+kk
        end if
        do k=1,kk
           nn=nn+1
           mmode_list(5,nn)=n
        end do
     end if
  end do

end subroutine inmi_coupler_ew2ns0

subroutine inmi_coupler_ew2ns1(mype)

!$$$  subprogram documentation block
!                .      .    .
! subprogram:    inmi_coupler_ew2ns1
!
!   prgrmmr:
!
! abstract:
!
! program history log:
!   2008-04-04  safford -- add subprogram doc block
!   2010-03-31  treadon - replace specmod jcap with sp_a structure
!
!   input argument list:
!     mype     - mpi task id
!
!   output argument list:
!
! attributes:
!   language:  f90
!   machine:   ibm RS/6000 SP
!
!$$$

  use mod_vtrans, only: nvmodes_keep
  use mpimod, only: npe,mpi_comm_world,ierror,mpi_integer4,mpi_sum
  use gridmod, only: sp_a
  implicit none

  integer(i_kind),intent(in   )::mype

  integer(i_kind) mmode2_list(0:sp_a%jcap,-nvmodes_keep:nvmodes_keep)
  integer(i_kind) i,ip12,ipe,j,k,m,nn,m1,m2,ilat,imode,imode1,imode2
  integer(i_kind) mpi_string1
  integer(i_kind) ibad,ibad0,loop

  allocate(nsend(npe),nrecv(npe),ndsend(npe+1),ndrecv(npe+1))
  nn=0
  mmode2_list=0
  do j=1,(sp_a%jcap+1)*nvmodes_keep
     m1=mmode_list(1,j)
     m2=mmode_list(2,j)
     imode1=mmode_list(3,j)
     imode2=mmode_list(4,j)
     if(imode1 == imode2) then
        if(mmode2_list(m1,imode1) /= 0.or.mmode2_list(m2,imode1) /= 0) then
           if(mype == 0) write(6,*)' problem in inmi_coupler_ew2ns'
           call mpi_finalize(i)
           stop
        end if
        mmode2_list(m1,imode1)=j
        mmode2_list(m2,imode1)=j
     end if
     if(m1 == m2) then
        if(mmode2_list(m1,imode1) /= 0.or.mmode2_list(m1,imode2) /= 0) then
           if(mype == 0) write(6,*)' problem in inmi_coupler_ew2ns'
           call mpi_finalize(i)
           stop
        end if
        mmode2_list(m1,imode1)=j
        mmode2_list(m1,imode2)=j
     end if
  end do
  do imode=-nvmodes_keep,nvmodes_keep
     if(imode == 0) cycle
     do m=0,sp_a%jcap
        if(mmode2_list(m,imode) == 0) then
           if(mype == 0) write(6,*)' problem in inmi_coupler_ew2ns'
           call mpi_finalize(i)
           stop
        end if
     end do
  end do
                
!  obtain counts of points to send to each pe from this pe

  nsend=0
  do k=nlatm_0,nlatm_1
     imode=mode_list(2,k)
     do m=0,sp_a%jcap
        j=mmode2_list(m,imode)
        ipe=mmode_list(5,j)
        nsend(ipe)=nsend(ipe)+1
        j=mmode2_list(m,-imode)
        ipe=mmode_list(5,j)
        nsend(ipe)=nsend(ipe)+1
     end do
  end do

  ndsend(1)=0
  do i=2,npe+1
     ndsend(i)=ndsend(i-1)+nsend(i-1)
  end do
  nallsend=ndsend(npe+1)
  allocate(info_send(6,nallsend))
  nsend=0
  ibad =0
  do k=nlatm_0,nlatm_1
     ilat=mode_list(1,k)
     do loop=1,2
        imode=mode_list(2,k)
        if(loop == 2) imode=-mode_list(2,k)
        do m=0,sp_a%jcap
           j=mmode2_list(m,imode)
           m1=mmode_list(1,j)
           m2=mmode_list(2,j)
           imode1=mmode_list(3,j)
           imode2=mmode_list(4,j)
           ipe=mmode_list(5,j)
           ip12=0
           if(m1 == m2.and.imode == imode1) ip12=1
           if(m1 == m2.and.imode == imode2) ip12=2
           if(imode1 == imode2.and.m == m1) ip12=1
           if(imode1 == imode2.and.m == m2) ip12=2
           if(ip12 == 0) ibad=ibad+1
           ipe=mmode_list(5,j)
           nsend(ipe)=nsend(ipe)+1
           info_send(1,ndsend(ipe)+nsend(ipe))=k
           info_send(2,ndsend(ipe)+nsend(ipe))=ilat
           info_send(3,ndsend(ipe)+nsend(ipe))=imode
           info_send(4,ndsend(ipe)+nsend(ipe))=m
           info_send(5,ndsend(ipe)+nsend(ipe))=j
           info_send(6,ndsend(ipe)+nsend(ipe))=ip12
        end do
     end do
  end do

  call mpi_allreduce(ibad,ibad0,1,mpi_integer4,mpi_sum,mpi_comm_world,ierror)
  if(ibad0 >  0) then
     if(mype == 0) write(0,*)' ibad = ',ibad0,'  inconsistency in inmi_coupler_ew2ns1'
     call mpi_finalize(ierror)
     stop
  end if

  call mpi_alltoall(nsend,1,mpi_integer4,nrecv,1,mpi_integer4,mpi_comm_world,ierror)
  ndrecv(1)=0
  do i=2,npe+1
     ndrecv(i)=ndrecv(i-1)+nrecv(i-1)
  end do
  nallrecv=ndrecv(npe+1)
  allocate(info_recv(6,nallrecv))
  call mpi_type_contiguous(6,mpi_integer4,mpi_string1,ierror)
  call mpi_type_commit(mpi_string1,ierror)
  call mpi_alltoallv(info_send,nsend,ndsend,mpi_string1, &
                     info_recv,nrecv,ndrecv,mpi_string1,mpi_comm_world,ierror)
  call mpi_type_free(mpi_string1,ierror)

end subroutine inmi_coupler_ew2ns1

subroutine inmi_coupler_ew2ns(uvm_ewtrans,uvm_ns)

!$$$  subprogram documentation block
!                .      .    .
! subprogram:    inmi_coupler_ew2ns
!
!   prgrmmr:
!
! abstract:
!
! program history log:
!   2008-04-04  safford -- add subprogram doc block, rm unused vars and uses
!   2010-03-31  treadon - replace specmod components with sp_a structure
!
!   input argument list:
!     uvm_ewtrans -
!
!   output argument list:
!     uvm_ns      -
!
! attributes:
!   language:  f90
!   machine:   ibm RS/6000 SP
!
!$$$

  use kinds, only: r_kind
  use gridmod, only: nlat
  use mpimod, only: mpi_comm_world,ierror,mpi_rtype
  use gridmod, only: sp_a
  implicit none

  real(r_kind),dimension(2,3,2,0:sp_a%jcap,nlatm_0:nlatm_1),intent(in   ) :: uvm_ewtrans
  real(r_kind),dimension(3,2,nlat,2,m_0:m_1)          ,intent(  out) :: uvm_ns

  integer(i_kind) ip12,j,m,mm,ilat,ilatm,imode
  integer(i_kind) mpi_string1
  real(r_kind),allocatable::sendbuf(:,:,:),recvbuf(:,:,:)
  integer(i_kind) loop

  allocate(sendbuf(3,2,nallsend))
  do j=1,nallsend
     ilatm=info_send(1,j)
     imode=info_send(3,j)
     m=info_send(4,j)
     loop=1
     if(imode <  0) loop=2
     sendbuf(1,1,j)=uvm_ewtrans(loop,1,1,m,ilatm)
     sendbuf(2,1,j)=uvm_ewtrans(loop,2,1,m,ilatm)
     sendbuf(3,1,j)=uvm_ewtrans(loop,3,1,m,ilatm)
     sendbuf(1,2,j)=uvm_ewtrans(loop,1,2,m,ilatm)
     sendbuf(2,2,j)=uvm_ewtrans(loop,2,2,m,ilatm)
     sendbuf(3,2,j)=uvm_ewtrans(loop,3,2,m,ilatm)
  end do
  allocate(recvbuf(3,2,nallrecv))
  call mpi_type_contiguous(6,mpi_rtype,mpi_string1,ierror)
  call mpi_type_commit(mpi_string1,ierror)
  call mpi_alltoallv(sendbuf,nsend,ndsend,mpi_string1, &
                     recvbuf,nrecv,ndrecv,mpi_string1,mpi_comm_world,ierror)
  call mpi_type_free(mpi_string1,ierror)
  deallocate(sendbuf)

  do j=1,nallrecv
     ilat=info_recv(2,j)
     mm=info_recv(5,j)
     ip12=info_recv(6,j)
     uvm_ns(1,1,ilat,ip12,mm)=recvbuf(1,1,j)
     uvm_ns(2,1,ilat,ip12,mm)=recvbuf(2,1,j)
     uvm_ns(3,1,ilat,ip12,mm)=recvbuf(3,1,j)
     uvm_ns(1,2,ilat,ip12,mm)=recvbuf(1,2,j)
     uvm_ns(2,2,ilat,ip12,mm)=recvbuf(2,2,j)
     uvm_ns(3,2,ilat,ip12,mm)=recvbuf(3,2,j)
  end do
  deallocate(recvbuf)

end subroutine inmi_coupler_ew2ns

subroutine inmi_coupler_ns2ew(uvm_ewtrans,uvm_ns)

!$$$  subprogram documentation block
!                .      .    .
! subprogram:    inmi_coupler_ns2ew
!
!   prgrmmr:
!
! abstract:
!
! program history log:
!   2008-04-04  safford -- add subprogram doc block, rm unused vars and uses
!   2010-03-31  treadon - replace specmod jcap with sp_a structure
!
!   input argument list:
!     uvm_ns   -
!
!   output argument list:
!     uvm_ewtrans -
!
! attributes:
!   language:  f90
!   machine:   ibm RS/6000 SP
!
!$$$

  use kinds, only: r_kind
  use gridmod, only: nlat
  use mpimod, only: mpi_comm_world,ierror,mpi_rtype
  use gridmod, only: sp_a
  implicit none

  real(r_kind),dimension(2,3,2,0:sp_a%jcap,nlatm_0:nlatm_1),intent(  out) :: uvm_ewtrans
  real(r_kind),dimension(3,2,nlat,2,m_0:m_1)          ,intent(in   ) :: uvm_ns

  integer(i_kind) ip12,j,m,mm,ilat,ilatm,imode
  integer(i_kind) mpi_string1
  real(r_kind),allocatable::sendbuf(:,:,:),recvbuf(:,:,:)
  integer(i_kind) loop


  allocate(recvbuf(3,2,nallrecv))
  do j=1,nallrecv
     ilat=info_recv(2,j)
     mm=info_recv(5,j)
     ip12=info_recv(6,j)
     recvbuf(1,1,j)=uvm_ns(1,1,ilat,ip12,mm)
     recvbuf(2,1,j)=uvm_ns(2,1,ilat,ip12,mm)
     recvbuf(3,1,j)=uvm_ns(3,1,ilat,ip12,mm)
     recvbuf(1,2,j)=uvm_ns(1,2,ilat,ip12,mm)
     recvbuf(2,2,j)=uvm_ns(2,2,ilat,ip12,mm)
     recvbuf(3,2,j)=uvm_ns(3,2,ilat,ip12,mm)
  end do
  allocate(sendbuf(3,2,nallsend))
  call mpi_type_contiguous(6,mpi_rtype,mpi_string1,ierror)
  call mpi_type_commit(mpi_string1,ierror)
  call mpi_alltoallv(recvbuf,nrecv,ndrecv,mpi_string1, &
                     sendbuf,nsend,ndsend,mpi_string1,mpi_comm_world,ierror)
  call mpi_type_free(mpi_string1,ierror)
  deallocate(recvbuf)
  do j=1,nallsend
     ilatm=info_send(1,j)
     imode=info_send(3,j)
     m=info_send(4,j)
     loop=1
     if(imode <  0) loop=2
     uvm_ewtrans(loop,1,1,m,ilatm)=sendbuf(1,1,j)
     uvm_ewtrans(loop,2,1,m,ilatm)=sendbuf(2,1,j)
     uvm_ewtrans(loop,3,1,m,ilatm)=sendbuf(3,1,j)
     uvm_ewtrans(loop,1,2,m,ilatm)=sendbuf(1,2,j)
     uvm_ewtrans(loop,2,2,m,ilatm)=sendbuf(2,2,j)
     uvm_ewtrans(loop,3,2,m,ilatm)=sendbuf(3,2,j)
  end do
  deallocate(sendbuf)

end subroutine inmi_coupler_ns2ew

subroutine inmi_nsuvm2zdm(uvm_ns,zdm_hat)

!$$$  subprogram documentation block
!                .      .    .
! subprogram:    inmi_nsuv2zdm
!
!   prgrmmr:
!
! abstract:
!
! program history log:
!   2008-04-04  safford -- add subprogram doc block, rm unused uses
!   2010-03-31  treadon - replace specmod components with sp_a structure
!
!   input argument list:
!     uvm_ns   -
!
!   output argument list:
!     zdm_had  -
!
! attributes:
!   language:  f90
!   machine:   ibm RS/6000 SP
!   machine:   ibm RS/6000 SP
!
!$$$

  use kinds, only: r_kind
  use gridmod, only: nlat,sp_a
  use constants, only: zero
  implicit none

  real(r_kind),dimension(3,2,nlat,2,m_0:m_1),intent(in   ) :: uvm_ns
  real(r_kind),dimension(3,2,nlat,2,m_0:m_1),intent(  out) :: zdm_hat

  integer(i_kind) i,ics,j,jnorth,jsouth,m,mm,n,ipair
  real(r_kind) spcz(2,0:sp_a%jcap),spcd(2,0:sp_a%jcap),spcp(2,0:sp_a%jcap),spcu(2,0:sp_a%jcap+1),spcv(2,0:sp_a%jcap+1)
  real(r_kind) plnloc(0:sp_a%jcap+1)
  real(r_kind) fu(2,2),fv(2,2),fp(2,2)
  real(r_kind):: enn1_loc(sp_a%ncd2),elonn1_loc(sp_a%ncd2),eon_loc(sp_a%ncd2),eontop_loc(sp_a%jcap+1)

  enn1_loc  = sp_a%enn1
  elonn1_loc= sp_a%elonn1
  eon_loc   = sp_a%eon
  eontop_loc= sp_a%eontop

!$omp parallel do  schedule(dynamic,1) private(mm,ipair,m,ics,i,n) &
!$omp private(spcz,spcd,spcp,spcu,spcv,j,jnorth,jsouth,plnloc,fu,fv,fp)
  do mm=m_0,m_1
     do ipair=1,2
        m=mmode_list(ipair,mm)
        ics=1+m*(2*sp_a%jcap+3-m)/2

        do n=m,sp_a%jcap
           spcp(1,n)=zero
           spcp(2,n)=zero
        end do
        do n=m,sp_a%jcap+1
           spcu(1,n)=zero
           spcu(2,n)=zero
           spcv(1,n)=zero
           spcv(2,n)=zero
        end do

        do j=sp_a%jb,sp_a%je
           jsouth=1+j
           jnorth=nlat-j

           fu(1,1)=uvm_ns(1,1,jnorth,ipair,mm)/sp_a%clat(j)**2
           fu(2,1)=uvm_ns(1,2,jnorth,ipair,mm)/sp_a%clat(j)**2
           fu(1,2)=uvm_ns(1,1,jsouth,ipair,mm)/sp_a%clat(j)**2
           fu(2,2)=uvm_ns(1,2,jsouth,ipair,mm)/sp_a%clat(j)**2
           fv(1,1)=uvm_ns(2,1,jnorth,ipair,mm)/sp_a%clat(j)**2
           fv(2,1)=uvm_ns(2,2,jnorth,ipair,mm)/sp_a%clat(j)**2
           fv(1,2)=uvm_ns(2,1,jsouth,ipair,mm)/sp_a%clat(j)**2
           fv(2,2)=uvm_ns(2,2,jsouth,ipair,mm)/sp_a%clat(j)**2
           fp(1,1)=uvm_ns(3,1,jnorth,ipair,mm)
           fp(2,1)=uvm_ns(3,2,jnorth,ipair,mm)
           fp(1,2)=uvm_ns(3,1,jsouth,ipair,mm)
           fp(2,2)=uvm_ns(3,2,jsouth,ipair,mm)
!           create plnloc

           do n=m,sp_a%jcap
              plnloc(n)=sp_a%pln(ics+n-m,j)
           end do
           plnloc(sp_a%jcap+1)=sp_a%plntop(m+1,j)
 
           call spanaly_ns(sp_a%jcap,m,sp_a%wlat(j),sp_a%clat(j),plnloc(m),1 ,fu,spcu(1,m))
           call spanaly_ns(sp_a%jcap,m,sp_a%wlat(j),sp_a%clat(j),plnloc(m),1 ,fv,spcv(1,m))
           call spanaly_ns(sp_a%jcap,m,sp_a%wlat(j),sp_a%clat(j),plnloc(m),0,fp,spcp(1,m))

        end do

        call spuv2dz_ns(sp_a%jcap,m,enn1_loc(ics),elonn1_loc(ics),eon_loc(ics),eontop_loc(m+1), &
                spcu(1,m),spcv(1,m),spcu(1,sp_a%jcap+1),spcv(1,sp_a%jcap+1),spcd(1,m),spcz(1,m))

        i=0
        do n=m,sp_a%jcap
           i=i+1
           zdm_hat(1,1,i,ipair,mm)=spcz(1,n)
           zdm_hat(1,2,i,ipair,mm)=spcz(2,n)
           zdm_hat(2,1,i,ipair,mm)=spcd(1,n)
           zdm_hat(2,2,i,ipair,mm)=spcd(2,n)
           zdm_hat(3,1,i,ipair,mm)=spcp(1,n)
           zdm_hat(3,2,i,ipair,mm)=spcp(2,n)
        end do

     end do

  end do

end subroutine inmi_nsuvm2zdm

subroutine inmi_nszdm2uvm_ad(uvm_ns,zdm_hat)
!$$$  subprogram documentation block
!                .      .    .
! subprogram:    inmi_nszdm2uvm_ad
!
!   prgrmmr:
!
! abstract:
!
! program history log:
!   2008-04-04  safford -- add subprogram doc block, rm unused uses
!   2010-03-31  treadon - replace specmod components with sp_a structure
!
!   input argument list:
!     uvm_ns   -
!
!   output argument list:
!     zdm_had  -
!
! attributes:
!   language:  f90
!   machine:   ibm RS/6000 SP
!
!$$$

  use kinds, only: r_kind
  use gridmod, only: nlat,sp_a
  use constants, only: zero
  implicit none

  real(r_kind),dimension(3,2,nlat,2,m_0:m_1),intent(in   ) :: uvm_ns
  real(r_kind),dimension(3,2,nlat,2,m_0:m_1),intent(  out) :: zdm_hat

  real(r_kind),dimension(3,2,nlat,2,m_0:m_1)::uvm_ns_temp
  integer(i_kind) i,ics,j,jnorth,jsouth,m,mm,n,ipair
  real(r_kind) spcz(2,0:sp_a%jcap),spcd(2,0:sp_a%jcap),spcp(2,0:sp_a%jcap),spcu(2,0:sp_a%jcap+1),spcv(2,0:sp_a%jcap+1)
  real(r_kind) plnloc(0:sp_a%jcap+1)
  real(r_kind) fu(2,2),fv(2,2),fp(2,2)
  real(r_kind):: enn1_loc(sp_a%ncd2),elonn1_loc(sp_a%ncd2),eon_loc(sp_a%ncd2),eontop_loc(sp_a%jcap+1)

  enn1_loc  = sp_a%enn1
  elonn1_loc= sp_a%elonn1
  eon_loc   = sp_a%eon
  eontop_loc= sp_a%eontop
  uvm_ns_temp=uvm_ns

!$omp parallel do  schedule(dynamic,1) private(mm,ipair,m,ics,i,n) &
!$omp private(spcz,spcd,spcp,spcu,spcv,j,jnorth,jsouth,plnloc,fu,fv,fp)
  do mm=m_0,m_1
     do ipair=1,2
        m=mmode_list(ipair,mm)
        ics=1+m*(2*sp_a%jcap+3-m)/2

        do n=m,sp_a%jcap
           spcp(1,n)=zero
           spcp(2,n)=zero
        end do
        do n=m,sp_a%jcap+1
           spcu(1,n)=zero
           spcu(2,n)=zero
           spcv(1,n)=zero
           spcv(2,n)=zero
        end do

!  adjoint of set pole values
        if(m == 0) then
           uvm_ns_temp(3,1,2,ipair,mm)=uvm_ns_temp(3,1,1,ipair,mm)+uvm_ns_temp(3,1,2,ipair,mm)
           uvm_ns_temp(3,1,nlat-1,ipair,mm)=uvm_ns_temp(3,1,nlat,ipair,mm)+uvm_ns_temp(3,1,nlat-1,ipair,mm)
        else if(m == 1) then
           uvm_ns_temp(1,1,2,ipair,mm)=uvm_ns_temp(1,1,1,ipair,mm)+uvm_ns_temp(1,1,2,ipair,mm)
           uvm_ns_temp(1,2,2,ipair,mm)=uvm_ns_temp(1,2,1,ipair,mm)+uvm_ns_temp(1,2,2,ipair,mm)
           uvm_ns_temp(2,1,2,ipair,mm)=uvm_ns_temp(2,1,1,ipair,mm)+uvm_ns_temp(2,1,2,ipair,mm)
           uvm_ns_temp(2,2,2,ipair,mm)=uvm_ns_temp(2,2,1,ipair,mm)+uvm_ns_temp(2,2,2,ipair,mm)
           uvm_ns_temp(1,1,nlat-1,ipair,mm)=uvm_ns_temp(1,1,nlat,ipair,mm)+uvm_ns_temp(1,1,nlat-1,ipair,mm)
           uvm_ns_temp(1,2,nlat-1,ipair,mm)=uvm_ns_temp(1,2,nlat,ipair,mm)+uvm_ns_temp(1,2,nlat-1,ipair,mm)
           uvm_ns_temp(2,1,nlat-1,ipair,mm)=uvm_ns_temp(2,1,nlat,ipair,mm)+uvm_ns_temp(2,1,nlat-1,ipair,mm)
           uvm_ns_temp(2,2,nlat-1,ipair,mm)=uvm_ns_temp(2,2,nlat,ipair,mm)+uvm_ns_temp(2,2,nlat-1,ipair,mm)
        end if

        do j=sp_a%jb,sp_a%je
           jsouth=1+j
           jnorth=nlat-j

           fu(1,1)=uvm_ns_temp(1,1,jnorth,ipair,mm)/(sp_a%wlat(j)*sp_a%clat(j)**2)
           fu(2,1)=uvm_ns_temp(1,2,jnorth,ipair,mm)/(sp_a%wlat(j)*sp_a%clat(j)**2)
           fu(1,2)=uvm_ns_temp(1,1,jsouth,ipair,mm)/(sp_a%wlat(j)*sp_a%clat(j)**2)
           fu(2,2)=uvm_ns_temp(1,2,jsouth,ipair,mm)/(sp_a%wlat(j)*sp_a%clat(j)**2)
           fv(1,1)=uvm_ns_temp(2,1,jnorth,ipair,mm)/(sp_a%wlat(j)*sp_a%clat(j)**2)
           fv(2,1)=uvm_ns_temp(2,2,jnorth,ipair,mm)/(sp_a%wlat(j)*sp_a%clat(j)**2)
           fv(1,2)=uvm_ns_temp(2,1,jsouth,ipair,mm)/(sp_a%wlat(j)*sp_a%clat(j)**2)
           fv(2,2)=uvm_ns_temp(2,2,jsouth,ipair,mm)/(sp_a%wlat(j)*sp_a%clat(j)**2)
           fp(1,1)=uvm_ns_temp(3,1,jnorth,ipair,mm)/sp_a%wlat(j)
           fp(2,1)=uvm_ns_temp(3,2,jnorth,ipair,mm)/sp_a%wlat(j)
           fp(1,2)=uvm_ns_temp(3,1,jsouth,ipair,mm)/sp_a%wlat(j)
           fp(2,2)=uvm_ns_temp(3,2,jsouth,ipair,mm)/sp_a%wlat(j)
!           create plnloc

           do n=m,sp_a%jcap
              plnloc(n)=sp_a%pln(ics+n-m,j)
           end do
           plnloc(sp_a%jcap+1)=sp_a%plntop(m+1,j)

           call spanaly_ns(sp_a%jcap,m,sp_a%wlat(j),sp_a%clat(j),plnloc(m),1 ,fu,spcu(1,m))
           call spanaly_ns(sp_a%jcap,m,sp_a%wlat(j),sp_a%clat(j),plnloc(m),1 ,fv,spcv(1,m))
           call spanaly_ns(sp_a%jcap,m,sp_a%wlat(j),sp_a%clat(j),plnloc(m),0,fp,spcp(1,m))

        end do

        call spuv2dz_ns(sp_a%jcap,m,enn1_loc(ics),elonn1_loc(ics),eon_loc(ics),eontop_loc(m+1), &
                spcu(1,m),spcv(1,m),spcu(1,sp_a%jcap+1),spcv(1,sp_a%jcap+1),spcd(1,m),spcz(1,m))

        i=0
        if(m == 0) then
           i=i+1
           zdm_hat(1,1,i,ipair,mm)=zero
           zdm_hat(1,2,i,ipair,mm)=zero
           zdm_hat(2,1,i,ipair,mm)=zero
           zdm_hat(2,2,i,ipair,mm)=zero
           zdm_hat(3,1,i,ipair,mm)=spcp(1,0)
           zdm_hat(3,2,i,ipair,mm)=spcp(2,0)
        end if
        do n=max(1,m),sp_a%jcap
           i=i+1
           zdm_hat(1,1,i,ipair,mm)=spcz(1,n)/sp_a%enn1(ics+n-m)
           zdm_hat(1,2,i,ipair,mm)=spcz(2,n)/sp_a%enn1(ics+n-m)
           zdm_hat(2,1,i,ipair,mm)=spcd(1,n)/sp_a%enn1(ics+n-m)
           zdm_hat(2,2,i,ipair,mm)=spcd(2,n)/sp_a%enn1(ics+n-m)
           zdm_hat(3,1,i,ipair,mm)=spcp(1,n)
           zdm_hat(3,2,i,ipair,mm)=spcp(2,n)
        end do

     end do

  end do

end subroutine inmi_nszdm2uvm_ad

subroutine inmi_nszdm2uvm(uvm_ns,zdm_hat)
!$$$  subprogram documentation block
!                .      .    .
! subprogram:    inmi_szdm2uvm
!
!   prgrmmr:
!
! abstract:
!
! program history log:
!   2008-04-04  safford -- add subprogram doc block, rm unused uses
!   2010-03-31  treadon - replace specmod components with sp_a structure
!
!   input argument list:
!     uvm_ns   -
!
!   output argument list:
!     zdm_had  -
!
! attributes:
!   language:  f90
!   machine:   ibm RS/6000 SP
!
!$$$

  use kinds, only: r_kind
  use gridmod, only: nlat,sp_a
  use constants, only: zero
  implicit none

  real(r_kind),dimension(3,2,nlat,2,m_0:m_1),intent(  out) :: uvm_ns
  real(r_kind),dimension(3,2,nlat,2,m_0:m_1),intent(in   ) :: zdm_hat

  integer(i_kind) i,ics,j,jnorth,jsouth,m,mm,n,ipair
  real(r_kind) spcz(2,0:sp_a%jcap),spcd(2,0:sp_a%jcap),spcp(2,0:sp_a%jcap),spcu(2,0:sp_a%jcap+1),spcv(2,0:sp_a%jcap+1)
  real(r_kind) plnloc(0:sp_a%jcap+1)
  real(r_kind) fu(2,2),fv(2,2),fp(2,2)

  real(r_kind):: elonn1_loc(sp_a%ncd2),eon_loc(sp_a%ncd2),eontop_loc(sp_a%jcap+1)

  elonn1_loc= sp_a%elonn1
  eon_loc   = sp_a%eon
  eontop_loc= sp_a%eontop

!$omp parallel do  schedule(dynamic,1) private(mm,ipair,m,ics,i,n) &
!$omp private(spcz,spcd,spcp,spcu,spcv,j,jnorth,jsouth,plnloc,fu,fv,fp)
  do mm=m_0,m_1
     do ipair=1,2
        m=mmode_list(ipair,mm)
        ics=1+m*(2*sp_a%jcap+3-m)/2

!           gather up spcz, spcd, spcp

        i=0
        do n=m,sp_a%jcap
           i=i+1
           spcz(1,n)=zdm_hat(1,1,i,ipair,mm)
           spcz(2,n)=zdm_hat(1,2,i,ipair,mm)
           spcd(1,n)=zdm_hat(2,1,i,ipair,mm)
           spcd(2,n)=zdm_hat(2,2,i,ipair,mm)
           spcp(1,n)=zdm_hat(3,1,i,ipair,mm)
           spcp(2,n)=zdm_hat(3,2,i,ipair,mm)
        end do

!           convert to spcu, spcv

        call spdz2uv_ns(sp_a%jcap,m,elonn1_loc(ics),eon_loc(ics),eontop_loc(m+1), &
                spcd(1,m),spcz(1,m),spcu(1,m),spcv(1,m),spcu(1,sp_a%jcap+1),spcv(1,sp_a%jcap+1))

        do j=sp_a%jb,sp_a%je
           jsouth=1+j
           jnorth=nlat-j

!           create plnloc

           do n=m,sp_a%jcap
              plnloc(n)=sp_a%pln(ics+n-m,j)
           end do
           plnloc(sp_a%jcap+1)=sp_a%plntop(m+1,j)
 
!          obtain f

           call spsynth_ns(sp_a%jcap,m,sp_a%clat(j),plnloc(m),1,spcu(1,m),fu)
           call spsynth_ns(sp_a%jcap,m,sp_a%clat(j),plnloc(m),1,spcv(1,m),fv)
           call spsynth_ns(sp_a%jcap,m,sp_a%clat(j),plnloc(m),0,spcp(1,m),fp)

!          scatter back to output pairs of lats
 
           uvm_ns(1,1,jnorth,ipair,mm)=fu(1,1)
           uvm_ns(1,2,jnorth,ipair,mm)=fu(2,1)
           uvm_ns(1,1,jsouth,ipair,mm)=fu(1,2)
           uvm_ns(1,2,jsouth,ipair,mm)=fu(2,2)
           uvm_ns(2,1,jnorth,ipair,mm)=fv(1,1)
           uvm_ns(2,2,jnorth,ipair,mm)=fv(2,1)
           uvm_ns(2,1,jsouth,ipair,mm)=fv(1,2)
           uvm_ns(2,2,jsouth,ipair,mm)=fv(2,2)
           uvm_ns(3,1,jnorth,ipair,mm)=fp(1,1)
           uvm_ns(3,2,jnorth,ipair,mm)=fp(2,1)
           uvm_ns(3,1,jsouth,ipair,mm)=fp(1,2)
           uvm_ns(3,2,jsouth,ipair,mm)=fp(2,2)

        end do

!  set pole values
        if(m == 0) then
           uvm_ns(1,1,1,ipair,mm)=zero
           uvm_ns(1,2,1,ipair,mm)=zero
           uvm_ns(2,1,1,ipair,mm)=zero
           uvm_ns(2,2,1,ipair,mm)=zero
           uvm_ns(3,1,1,ipair,mm)=uvm_ns(3,1,2,ipair,mm)
           uvm_ns(3,2,1,ipair,mm)=zero
           uvm_ns(1,1,nlat,ipair,mm)=zero
           uvm_ns(1,2,nlat,ipair,mm)=zero
           uvm_ns(2,1,nlat,ipair,mm)=zero
           uvm_ns(2,2,nlat,ipair,mm)=zero
           uvm_ns(3,1,nlat,ipair,mm)=uvm_ns(3,1,nlat-1,ipair,mm)
           uvm_ns(3,2,nlat,ipair,mm)=zero
        else if(m == 1) then
           uvm_ns(1,1,1,ipair,mm)=uvm_ns(1,1,2,ipair,mm)
           uvm_ns(1,2,1,ipair,mm)=uvm_ns(1,2,2,ipair,mm)
           uvm_ns(2,1,1,ipair,mm)=uvm_ns(2,1,2,ipair,mm)
           uvm_ns(2,2,1,ipair,mm)=uvm_ns(2,2,2,ipair,mm)
           uvm_ns(3,1,1,ipair,mm)=zero
           uvm_ns(3,2,1,ipair,mm)=zero
           uvm_ns(1,1,nlat,ipair,mm)=uvm_ns(1,1,nlat-1,ipair,mm)
           uvm_ns(1,2,nlat,ipair,mm)=uvm_ns(1,2,nlat-1,ipair,mm)
           uvm_ns(2,1,nlat,ipair,mm)=uvm_ns(2,1,nlat-1,ipair,mm)
           uvm_ns(2,2,nlat,ipair,mm)=uvm_ns(2,2,nlat-1,ipair,mm)
           uvm_ns(3,1,nlat,ipair,mm)=zero
           uvm_ns(3,2,nlat,ipair,mm)=zero
        else
           uvm_ns(1,1,1,ipair,mm)=zero
           uvm_ns(1,2,1,ipair,mm)=zero
           uvm_ns(2,1,1,ipair,mm)=zero
           uvm_ns(2,2,1,ipair,mm)=zero
           uvm_ns(3,1,1,ipair,mm)=zero
           uvm_ns(3,2,1,ipair,mm)=zero
           uvm_ns(1,1,nlat,ipair,mm)=zero
           uvm_ns(1,2,nlat,ipair,mm)=zero
           uvm_ns(2,1,nlat,ipair,mm)=zero
           uvm_ns(2,2,nlat,ipair,mm)=zero
           uvm_ns(3,1,nlat,ipair,mm)=zero
           uvm_ns(3,2,nlat,ipair,mm)=zero
        end if

     end do

  end do

end subroutine inmi_nszdm2uvm

subroutine inmi_nspcm_hat2pcm(pcm_ns,pcm_hat)

!$$$  subprogram documentation block
!                .      .    .                                       .
! subprogram:    inmi_nspcm_hat2pcm
!   prgmmr:
!
! abstract:
!
! program history log:
!   2009-08-13  lueken - added subprogram doc block
!   2010-03-31  treadon - replace specmod components with sp_a structure
!
!   input argument list:
!    pcm_hat
!
!   output argument list:
!    pcm_ns
!
! attributes:
!   language: f90
!   machine:
!
!$$$ end documentation block

  use kinds, only: r_kind
  use gridmod, only: nlat,sp_a
  use constants, only: zero
  implicit none

  real(r_kind),dimension(3,2,nlat,2,m_0:m_1),intent(  out) :: pcm_ns
  real(r_kind),dimension(3,2,nlat,2,m_0:m_1),intent(in   ) :: pcm_hat

  integer(i_kind) i,ics,j,jnorth,jsouth,m,mm,n,ipair
  real(r_kind) spcp(2,0:sp_a%jcap),spcc(2,0:sp_a%jcap),spcm(2,0:sp_a%jcap)
  real(r_kind) plnloc(0:sp_a%jcap+1)
  real(r_kind) fp(2,2),fc(2,2),fm(2,2)

!$omp parallel do  schedule(dynamic,1) private(mm,ipair,m,ics,i,n) &
!$omp private(spcc,spcm,spcp,j,jnorth,jsouth,plnloc,fc,fm,fp)
  do mm=m_0,m_1
     do ipair=1,2
        m=mmode_list(ipair,mm)
        ics=1+m*(2*sp_a%jcap+3-m)/2

!           gather up spcp, spcc, spcm

        i=0
        do n=m,sp_a%jcap
           i=i+1
           spcp(1,n)=pcm_hat(1,1,i,ipair,mm)
           spcp(2,n)=pcm_hat(1,2,i,ipair,mm)
           spcc(1,n)=pcm_hat(2,1,i,ipair,mm)
           spcc(2,n)=pcm_hat(2,2,i,ipair,mm)
           spcm(1,n)=pcm_hat(3,1,i,ipair,mm)
           spcm(2,n)=pcm_hat(3,2,i,ipair,mm)
        end do

        do j=sp_a%jb,sp_a%je
           jsouth=1+j
           jnorth=nlat-j

!           create plnloc

           do n=m,sp_a%jcap
              plnloc(n)=sp_a%pln(ics+n-m,j)
           end do
           plnloc(sp_a%jcap+1)=sp_a%plntop(m+1,j)

!          obtain f

           call spsynth_ns(sp_a%jcap,m,sp_a%clat(j),plnloc(m),0,spcp(1,m),fp)
           call spsynth_ns(sp_a%jcap,m,sp_a%clat(j),plnloc(m),0,spcc(1,m),fc)
           call spsynth_ns(sp_a%jcap,m,sp_a%clat(j),plnloc(m),0,spcm(1,m),fm)

!          scatter back to output pairs of lats

           pcm_ns(1,1,jnorth,ipair,mm)=fp(1,1)
           pcm_ns(1,2,jnorth,ipair,mm)=fp(2,1)
           pcm_ns(1,1,jsouth,ipair,mm)=fp(1,2)
           pcm_ns(1,2,jsouth,ipair,mm)=fp(2,2)
           pcm_ns(2,1,jnorth,ipair,mm)=fc(1,1)
           pcm_ns(2,2,jnorth,ipair,mm)=fc(2,1)
           pcm_ns(2,1,jsouth,ipair,mm)=fc(1,2)
           pcm_ns(2,2,jsouth,ipair,mm)=fc(2,2)
           pcm_ns(3,1,jnorth,ipair,mm)=fm(1,1)
           pcm_ns(3,2,jnorth,ipair,mm)=fm(2,1)
           pcm_ns(3,1,jsouth,ipair,mm)=fm(1,2)
           pcm_ns(3,2,jsouth,ipair,mm)=fm(2,2)
 
        end do

!  set pole values
        if(m == 0) then
           pcm_ns(1,1,1,ipair,mm)=pcm_ns(1,1,2,ipair,mm)
           pcm_ns(1,2,1,ipair,mm)=zero
           pcm_ns(1,1,nlat,ipair,mm)=pcm_ns(1,1,nlat-1,ipair,mm)
           pcm_ns(1,2,nlat,ipair,mm)=zero
           pcm_ns(2,1,1,ipair,mm)=pcm_ns(2,1,2,ipair,mm)
           pcm_ns(2,2,1,ipair,mm)=zero
           pcm_ns(2,1,nlat,ipair,mm)=pcm_ns(2,1,nlat-1,ipair,mm)
           pcm_ns(2,2,nlat,ipair,mm)=zero
           pcm_ns(3,1,1,ipair,mm)=pcm_ns(3,1,2,ipair,mm)
           pcm_ns(3,2,1,ipair,mm)=zero
           pcm_ns(3,1,nlat,ipair,mm)=pcm_ns(3,1,nlat-1,ipair,mm)
           pcm_ns(3,2,nlat,ipair,mm)=zero
        else
           pcm_ns(1,1,1,ipair,mm)=zero
           pcm_ns(1,2,1,ipair,mm)=zero
           pcm_ns(2,1,1,ipair,mm)=zero
           pcm_ns(2,2,1,ipair,mm)=zero
           pcm_ns(3,1,1,ipair,mm)=zero
           pcm_ns(3,2,1,ipair,mm)=zero
           pcm_ns(1,1,nlat,ipair,mm)=zero
           pcm_ns(1,2,nlat,ipair,mm)=zero
           pcm_ns(2,1,nlat,ipair,mm)=zero
           pcm_ns(2,2,nlat,ipair,mm)=zero
           pcm_ns(3,1,nlat,ipair,mm)=zero
           pcm_ns(3,2,nlat,ipair,mm)=zero
        end if

     end do

  end do

end subroutine inmi_nspcm_hat2pcm

subroutine inmi_nspcm_hat2pcm_ad(pcm_ns,pcm_hat)

!$$$  subprogram documentation block
!                .      .    .                                       .
! subprogram:    inmi_nspcm_hat2pcm_ad
!   prgmmr:
!
! abstract:
!
! program history log:
!   2009-08-13  lueken - added subprogram doc block
!   2010-03-31  treadon - replace specmod components with sp_a structure
!
!   input argument list:
!    pcm_ns
!
!   output argument list:
!    pcm_hat
!
! attributes:
!   language: f90
!   machine:
!
!$$$ end documentation block
  use kinds, only: r_kind
  use gridmod, only: nlat,sp_a
  use constants, only: zero
  implicit none

  real(r_kind),dimension(3,2,nlat,2,m_0:m_1),intent(in   ) :: pcm_ns
  real(r_kind),dimension(3,2,nlat,2,m_0:m_1),intent(  out) :: pcm_hat

  real(r_kind),dimension(3,2,nlat,2,m_0:m_1)::pcm_ns_temp
  integer(i_kind) i,ics,j,jnorth,jsouth,m,mm,n,ipair
  real(r_kind) spcp(2,0:sp_a%jcap),spcc(2,0:sp_a%jcap),spcm(2,0:sp_a%jcap)
  real(r_kind) plnloc(0:sp_a%jcap+1)
  real(r_kind) fp(2,2),fc(2,2),fm(2,2)

  pcm_ns_temp=pcm_ns
  pcm_hat=zero
!$omp parallel do  schedule(dynamic,1) private(mm,ipair,m,ics,i,n) &
!$omp private(spcc,spcm,spcp,j,jnorth,jsouth,plnloc,fc,fm,fp)
  do mm=m_0,m_1
     do ipair=1,2
        m=mmode_list(ipair,mm)
        ics=1+m*(2*sp_a%jcap+3-m)/2
 
        do n=m,sp_a%jcap
           spcp(1,n)=zero
           spcp(2,n)=zero
           spcc(1,n)=zero
           spcc(2,n)=zero
           spcm(1,n)=zero
           spcm(2,n)=zero
        end do

!  adjoint of set pole values

        if(m == 0) then
           pcm_ns_temp(1,1,        2,ipair,mm)=pcm_ns_temp(1,1,   1,ipair,mm)+pcm_ns_temp(1,1,        2,ipair,mm)
           pcm_ns_temp(1,1,nlat-1,ipair,mm)=pcm_ns_temp(1,1,nlat,ipair,mm)+pcm_ns_temp(1,1,nlat-1,ipair,mm)
           pcm_ns_temp(2,1,        2,ipair,mm)=pcm_ns_temp(2,1,   1,ipair,mm)+pcm_ns_temp(2,1,        2,ipair,mm)
           pcm_ns_temp(2,1,nlat-1,ipair,mm)=pcm_ns_temp(2,1,nlat,ipair,mm)+pcm_ns_temp(2,1,nlat-1,ipair,mm)
           pcm_ns_temp(3,1,        2,ipair,mm)=pcm_ns_temp(3,1,   1,ipair,mm)+pcm_ns_temp(3,1,        2,ipair,mm)
           pcm_ns_temp(3,1,nlat-1,ipair,mm)=pcm_ns_temp(3,1,nlat,ipair,mm)+pcm_ns_temp(3,1,nlat-1,ipair,mm)
        end if

        do j=sp_a%jb,sp_a%je
           jsouth=1+j
           jnorth=nlat-j
 
!          adjoint of scatter back to output pairs of lats

           fp(1,1)=pcm_ns_temp(1,1,jnorth,ipair,mm)/sp_a%wlat(j)
           fp(2,1)=pcm_ns_temp(1,2,jnorth,ipair,mm)/sp_a%wlat(j)
           fp(1,2)=pcm_ns_temp(1,1,jsouth,ipair,mm)/sp_a%wlat(j)
           fp(2,2)=pcm_ns_temp(1,2,jsouth,ipair,mm)/sp_a%wlat(j)
           fc(1,1)=pcm_ns_temp(2,1,jnorth,ipair,mm)/sp_a%wlat(j)
           fc(2,1)=pcm_ns_temp(2,2,jnorth,ipair,mm)/sp_a%wlat(j)
           fc(1,2)=pcm_ns_temp(2,1,jsouth,ipair,mm)/sp_a%wlat(j)
           fc(2,2)=pcm_ns_temp(2,2,jsouth,ipair,mm)/sp_a%wlat(j)
           fm(1,1)=pcm_ns_temp(3,1,jnorth,ipair,mm)/sp_a%wlat(j)
           fm(2,1)=pcm_ns_temp(3,2,jnorth,ipair,mm)/sp_a%wlat(j)
           fm(1,2)=pcm_ns_temp(3,1,jsouth,ipair,mm)/sp_a%wlat(j)
           fm(2,2)=pcm_ns_temp(3,2,jsouth,ipair,mm)/sp_a%wlat(j)
 
!           create plnloc

           do n=m,sp_a%jcap
              plnloc(n)=sp_a%pln(ics+n-m,j)
           end do
           plnloc(sp_a%jcap+1)=sp_a%plntop(m+1,j)

!          adjoint of obtain f

           call spanaly_ns(sp_a%jcap,m,sp_a%wlat(j),sp_a%clat(j),plnloc(m),0,fp,spcp(1,m))
           call spanaly_ns(sp_a%jcap,m,sp_a%wlat(j),sp_a%clat(j),plnloc(m),0,fc,spcc(1,m))
           call spanaly_ns(sp_a%jcap,m,sp_a%wlat(j),sp_a%clat(j),plnloc(m),0,fm,spcm(1,m))
 
        end do

!       adjoint of gather up spcp, spcc, spcm

        i=0
        do n=m,sp_a%jcap
           i=i+1
           pcm_hat(1,1,i,ipair,mm)=spcp(1,n)
           pcm_hat(1,2,i,ipair,mm)=spcp(2,n)
           pcm_hat(2,1,i,ipair,mm)=spcc(1,n)
           pcm_hat(2,2,i,ipair,mm)=spcc(2,n)
           pcm_hat(3,1,i,ipair,mm)=spcm(1,n)
           pcm_hat(3,2,i,ipair,mm)=spcm(2,n)
        end do

     end do

  end do

end subroutine inmi_nspcm_hat2pcm_ad

subroutine inmi_nsuvm2zdm_ad(uvm_ns,zdm_hat)

!$$$  subprogram documentation block
!                .      .    .                                       .
! subprogram:    inmi_nsuvm2zdm_ad
!   prgmmr:
!
! abstract:
!
! program history log:
!   2009-08-13  lueken - added subprogram doc block
!   2010-03-31  treadon - replace specmod components with sp_a structure
!
!   input argument list:
!    uvm_ns
!    zdm_hat
!
!   output argument list:
!    uvm_ns
!
! attributes:
!   language: f90
!   machine:
!
!$$$ end documentation block

  use kinds, only: r_kind
  use gridmod, only: nlat,sp_a
  implicit none

  real(r_kind),dimension(3,2,nlat,2,m_0:m_1),intent(inout) :: uvm_ns
  real(r_kind),dimension(3,2,nlat,2,m_0:m_1),intent(in   ) :: zdm_hat

  integer(i_kind) i,ics,j,jnorth,jsouth,m,mm,n,ipair
  real(r_kind) spcz(2,0:sp_a%jcap),spcd(2,0:sp_a%jcap),spcp(2,0:sp_a%jcap),spcu(2,0:sp_a%jcap+1),spcv(2,0:sp_a%jcap+1)
  real(r_kind) plnloc(0:sp_a%jcap+1)
  real(r_kind) fu(2,2),fv(2,2),fp(2,2)

  real(r_kind):: elonn1_loc(sp_a%ncd2),eon_loc(sp_a%ncd2),eontop_loc(sp_a%jcap+1)

  elonn1_loc= sp_a%elonn1
  eon_loc   = sp_a%eon
  eontop_loc= sp_a%eontop

!$omp parallel do  schedule(dynamic,1) private(mm,ipair,m,ics,i,n) &
!$omp private(spcz,spcd,spcp,spcu,spcv,j,jnorth,jsouth,plnloc,fu,fv,fp)
  do mm=m_0,m_1
     do ipair=1,2
        m=mmode_list(ipair,mm)
        ics=1+m*(2*sp_a%jcap+3-m)/2

!           gather up spcz, spcd, spcp

        i=0
        do n=m,sp_a%jcap
           i=i+1
           spcz(1,n)=zdm_hat(1,1,i,ipair,mm)*sp_a%enn1(ics+n-m)
           spcz(2,n)=zdm_hat(1,2,i,ipair,mm)*sp_a%enn1(ics+n-m)
           spcd(1,n)=zdm_hat(2,1,i,ipair,mm)*sp_a%enn1(ics+n-m)
           spcd(2,n)=zdm_hat(2,2,i,ipair,mm)*sp_a%enn1(ics+n-m)
           spcp(1,n)=zdm_hat(3,1,i,ipair,mm)
           spcp(2,n)=zdm_hat(3,2,i,ipair,mm)
        end do

!           convert to spcu, spcv

        call spdz2uv_ns(sp_a%jcap,m,elonn1_loc(ics),eon_loc(ics),eontop_loc(m+1), &
                spcd(1,m),spcz(1,m),spcu(1,m),spcv(1,m),spcu(1,sp_a%jcap+1),spcv(1,sp_a%jcap+1))

        do j=sp_a%jb,sp_a%je
           jsouth=1+j
           jnorth=nlat-j
 
!           create plnloc

           do n=m,sp_a%jcap
              plnloc(n)=sp_a%pln(ics+n-m,j)
           end do
           plnloc(sp_a%jcap+1)=sp_a%plntop(m+1,j)

!          obtain f

           call spsynth_ns(sp_a%jcap,m,sp_a%clat(j),plnloc(m),1,spcu(1,m),fu)
           call spsynth_ns(sp_a%jcap,m,sp_a%clat(j),plnloc(m),1,spcv(1,m),fv)
           call spsynth_ns(sp_a%jcap,m,sp_a%clat(j),plnloc(m),0,spcp(1,m),fp)

!          scatter back to output pairs of lats

           uvm_ns(1,1,jnorth,ipair,mm)=fu(1,1)*sp_a%wlat(j)
           uvm_ns(1,2,jnorth,ipair,mm)=fu(2,1)*sp_a%wlat(j)
           uvm_ns(1,1,jsouth,ipair,mm)=fu(1,2)*sp_a%wlat(j)
           uvm_ns(1,2,jsouth,ipair,mm)=fu(2,2)*sp_a%wlat(j)
           uvm_ns(2,1,jnorth,ipair,mm)=fv(1,1)*sp_a%wlat(j)
           uvm_ns(2,2,jnorth,ipair,mm)=fv(2,1)*sp_a%wlat(j)
           uvm_ns(2,1,jsouth,ipair,mm)=fv(1,2)*sp_a%wlat(j)
           uvm_ns(2,2,jsouth,ipair,mm)=fv(2,2)*sp_a%wlat(j)
           uvm_ns(3,1,jnorth,ipair,mm)=fp(1,1)*sp_a%wlat(j)
           uvm_ns(3,2,jnorth,ipair,mm)=fp(2,1)*sp_a%wlat(j)
           uvm_ns(3,1,jsouth,ipair,mm)=fp(1,2)*sp_a%wlat(j)
           uvm_ns(3,2,jsouth,ipair,mm)=fp(2,2)*sp_a%wlat(j)

        end do

     end do

  end do

end subroutine inmi_nsuvm2zdm_ad
      subroutine spdz2uv_ns(M,L,ELONN1,EON,EONTOP,D,Z,U,V,UTOP,VTOP)
!$$$  SUBPROGRAM DOCUMENTATION BLOCK
!
! SUBPROGRAM:    SPDZ2UV_ns  COMPUTE WINDS FROM div and vort for one zonal wave number
!   PRGMMR: IREDELL          ORG: W/NMC23     DATE: 92-10-31
!
! ABSTRACT: COMPUTES THE WIND COMPONENTS FROM DIVERGENCE AND VORTICITY
!           IN SPECTRAL SPACE.
!           SUBPROGRAM SPEPS SHOULD BE CALLED ALREADY.
!           IF L IS THE ZONAL WAVENUMBER, N IS THE TOTAL WAVENUMBER,
!           EPS(L,N)=SQRT((N**2-L**2)/(4*N**2-1)) AND A IS EARTH RADIUS,
!           THEN THE ZONAL WIND COMPONENT U IS COMPUTED AS
!             U(L,N)=-I*L/(N*(N+1))*A*D(L,N)
!                    +EPS(L,N+1)/(N+1)*A*Z(L,N+1)-EPS(L,N)/N*A*Z(L,N-1)
!           AND THE MERIDIONAL WIND COMPONENT V IS COMPUTED AS
!             V(L,N)=-I*L/(N*(N+1))*A*Z(L,N)
!                    -EPS(L,N+1)/(N+1)*A*D(L,N+1)+EPS(L,N)/N*A*D(L,N-1)
!           WHERE D IS DIVERGENCE AND Z IS VORTICITY.
!           U AND V ARE WEIGHTED BY THE COSINE OF LATITUDE.
!           EXTRA TERMS ARE COMPUTED OVER TOP OF THE SPECTRAL DOMAIN.
!           ADVANTAGE IS TAKEN OF THE FACT THAT EPS(L,L)=0
!           IN ORDER TO VECTORIZE OVER THE ENTIRE SPECTRAL DOMAIN.
!           Triangular truncation only
!
! PROGRAM HISTORY LOG:
!   91-10-31  MARK IREDELL
!   2006-09-05 parrish -- modify to do one zonal wave number only for parallel
!                         computation across processors by zonal wave number.
!   2010-05-14 derber  - make triangular truncation only
!
! USAGE:    CALL SPDZ2UV_ns(M,L,ELONN1,EON,EONTOP,D,Z,U,V,UTOP,VTOP)
!
!   INPUT ARGUMENT LIST:
!     M        - INTEGER SPECTRAL TRUNCATION
!     L        - zonal wave number
!     ELONN1   - REAL (L:M) L/(N*(N+1))*A
!     EON      - REAL (L:M) EPSILON/N*A
!     EONTOP   - REAL       EPSILON/N*A OVER TOP
!     D        - REAL (2,L:M) DIVERGENCE for zonal wave number L
!     Z        - REAL (2,L:M) VORTICITY for zonal wave number L
!
!   OUTPUT ARGUMENT LIST:
!     U        - REAL (2,L:M) ZONAL WIND (TIMES COSLAT) for zonal wave number L
!     V        - REAL (2,L:M) MERID WIND (TIMES COSLAT) for zonal wave number L
!     UTOP     - REAL (2) ZONAL WIND (TIMES COSLAT) OVER TOP for zonal wave number L
!     VTOP     - REAL (2) MERID WIND (TIMES COSLAT) OVER TOP for zonal wave number L
!
! ATTRIBUTES:
!   LANGUAGE: CRAY FORTRAN
!
!$$$
      use kinds, only: r_kind
      implicit none

      integer(i_kind),intent(in   ) :: m,l
      REAL(r_kind)   ,intent(in   ) :: ELONN1(L:M)
      REAL(r_kind)   ,intent(in   ) :: EON(L:M),EONTOP
      REAL(r_kind)   ,intent(in   ) :: D(2,L:M),Z(2,L:M)
      REAL(r_kind)   ,intent(  out) :: U(2,L:M),V(2,L:M)
      REAL(r_kind)   ,intent(  out) :: UTOP(2),VTOP(2)

      integer(i_kind) n
! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
!  COMPUTE WINDS IN THE SPECTRAL DOMAIN

      do n=l,m
         U(1,n)= elonn1(n)*d(2,n)
         U(2,n)=-elonn1(n)*d(1,n)
         V(1,n)= elonn1(n)*z(2,n)
         V(2,n)=-elonn1(n)*z(1,n)
      end do
      do n=l,m-1
         U(1,n)=u(1,n)+EON(n+1)*Z(1,n+1)
         U(2,n)=u(2,n)+EON(n+1)*Z(2,n+1)
         V(1,n)=v(1,n)-EON(n+1)*D(1,n+1)
         V(2,n)=v(2,n)-EON(n+1)*D(2,n+1)
      end do
      do n=l+1,m
         U(1,n)=u(1,n)-EON(n)*Z(1,n-1)
         U(2,n)=u(2,n)-EON(n)*Z(2,n-1)
         V(1,n)=v(1,n)+EON(n)*D(1,n-1)
         V(2,n)=v(2,n)+EON(n)*D(2,n-1)
      end do
! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
!  COMPUTE WINDS OVER TOP OF THE SPECTRAL DOMAIN
      UTOP(1)=-EONTOP*Z(1,m)
      UTOP(2)=-EONTOP*Z(2,m)
      VTOP(1)= EONTOP*D(1,m)
      VTOP(2)= EONTOP*D(2,m)
! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
      RETURN
      END subroutine spdz2uv_ns
!-----------------------------------------------------------------------
      subroutine spuv2dz_ns(M,L,ENN1,ELONN1,EON,EONTOP,U,V,UTOP,VTOP,D,Z)
!$$$  SUBPROGRAM DOCUMENTATION BLOCK
!
! SUBPROGRAM:    spuv2dz_ns  COMPUTE DIV,VORT FROM WINDS for one zonal wave number
!   PRGMMR: IREDELL          ORG: W/NMC23     DATE: 92-10-31
!
! ABSTRACT: COMPUTES THE DIVERGENCE AND VORTICITY FROM WIND COMPONENTS
!           IN SPECTRAL SPACE.
!           SUBPROGRAM SPEPS SHOULD BE CALLED ALREADY.
!           IF L IS THE ZONAL WAVENUMBER, N IS THE TOTAL WAVENUMBER,
!           EPS(L,N)=SQRT((N**2-L**2)/(4*N**2-1)) AND A IS EARTH RADIUS,
!           THEN THE DIVERGENCE D IS COMPUTED AS
!             D(L,N)=I*L*A*U(L,N)
!                    +EPS(L,N+1)*N*A*V(L,N+1)-EPS(L,N)*(N+1)*A*V(L,N-1)
!           AND THE VORTICITY Z IS COMPUTED AS
!             Z(L,N)=I*L*A*V(L,N)
!                    -EPS(L,N+1)*N*A*U(L,N+1)+EPS(L,N)*(N+1)*A*U(L,N-1)
!           WHERE U IS THE ZONAL WIND AND V IS THE MERIDIONAL WIND.
!           U AND V ARE WEIGHTED BY THE SECANT OF LATITUDE.
!           EXTRA TERMS ARE USED OVER TOP OF THE SPECTRAL DOMAIN.
!           ADVANTAGE IS TAKEN OF THE FACT THAT EPS(L,L)=0
!           IN ORDER TO VECTORIZE OVER THE ENTIRE SPECTRAL DOMAIN.
!           Triangular Truncation only
!
! PROGRAM HISTORY LOG:
!   91-10-31  MARK IREDELL
!   2006-09-05 parrish -- modify to do one zonal wave number only for parallel
!                         computation across processors by zonal wave number.
!   2010-05-14 derber  - Triangular truncation only
!
! USAGE:    CALL spuv2dz_ns(M,L,ENN1,ELONN1,EON,EONTOP,U,V,UTOP,VTOP,D,Z)
!
!   INPUT ARGUMENT LIST:
!     L        - zonal wave number
!     M        - INTEGER SPECTRAL TRUNCATION
!     ENN1     - REAL (L:M) N*(N+1)/A**2
!     ELONN1   - REAL (L:M) L/(N*(N+1))*A
!     EON      - REAL (L:M) EPSILON/N*A
!     EONTOP   - REAL    EPSILON/N*A OVER TOP
!     U        - REAL (2,L:M) ZONAL WIND (OVER COSLAT)
!     V        - REAL (2,L:M) MERID WIND (OVER COSLAT)
!     UTOP     - REAL (2) ZONAL WIND (OVER COSLAT) OVER TOP
!     VTOP     - REAL (2) MERID WIND (OVER COSLAT) OVER TOP
!
!   OUTPUT ARGUMENT LIST:
!     D        - REAL (2,L:M) DIVERGENCE
!     Z        - REAL (2,L:M) VORTICITY
!
! ATTRIBUTES:
!   LANGUAGE: CRAY FORTRAN
!
!$$$
      use kinds, only: r_kind
      implicit none

      integer(i_kind),intent(in   ) :: m,l
      REAL(r_kind)   ,intent(in   ) :: ENN1(L:M),ELONN1(L:M)
      REAL(r_kind)   ,intent(in   ) :: EON(L:M),EONTOP
      REAL(r_kind)   ,intent(in   ) :: U(2,L:M),V(2,L:M)
      REAL(r_kind)   ,intent(in   ) :: UTOP(2),VTOP(2)
      REAL(r_kind)   ,intent(  out) :: D(2,L:M),Z(2,L:M)

      integer(i_kind) n
! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
!  COMPUTE TERMS FROM THE SPECTRAL DOMAIN

      do n=l,m
         D(1,n)=-ELONN1(n)*U(2,n)
         D(2,n)= ELONN1(n)*U(1,n)
         Z(1,n)=-ELONN1(n)*V(2,n)
         Z(2,n)= ELONN1(n)*V(1,n)
      end do
      do n=l,m-1
         D(1,n)=d(1,n)+EON(n+1)*V(1,n+1)
         D(2,n)=d(2,n)+EON(n+1)*V(2,n+1)
         Z(1,n)=z(1,n)-EON(n+1)*U(1,n+1)
         Z(2,n)=z(2,n)-EON(n+1)*U(2,n+1)
      end do
      do n=l+1,m
         D(1,n)=d(1,n)-EON(n)*V(1,n-1)
         D(2,n)=d(2,n)-EON(n)*V(2,n-1)
         Z(1,n)=z(1,n)+EON(n)*U(1,n-1)
         Z(2,n)=z(2,n)+EON(n)*U(2,n-1)
      end do

!  COMPUTE TERMS FROM OVER TOP OF THE SPECTRAL DOMAIN
      n=m
      d(1,n)=d(1,n)+eontop*vtop(1)
      d(2,n)=d(2,n)+eontop*vtop(2)
      z(1,n)=z(1,n)-eontop*utop(1)
      z(2,n)=z(2,n)-eontop*utop(2)

!  MULTIPLY BY LAPLACIAN TERM
      do n=l,m
         d(1,n)=d(1,n)*enn1(n)
         d(2,n)=d(2,n)*enn1(n)
         z(1,n)=z(1,n)*enn1(n)
         z(2,n)=z(2,n)*enn1(n)
      end do

      RETURN
      end subroutine spuv2dz_ns

!-----------------------------------------------------------------------
      subroutine spanaly_ns(M,L,WGT,CLAT,PLN,MP,F,SPC)
!$$$  SUBPROGRAM DOCUMENTATION BLOCK
!
! SUBPROGRAM:    spanaly     ANALYZE SPECTRAL FROM FOURIER
!   PRGMMR: IREDELL          ORG: W/NMC23     DATE: 92-10-31
!
! ABSTRACT: ANALYZES SPECTRAL COEFFICIENTS FROM FOURIER COEFFICIENTS
!           FOR A LATITUDE PAIR (NORTHERN AND SOUTHERN HEMISPHERES).
!           VECTOR COMPONENTS ARE MULTIPLIED BY COSINE OF LATITUDE.
!           Triangular truncation only
!
! PROGRAM HISTORY LOG:
!   91-10-31  MARK IREDELL
!   94-08-01  MARK IREDELL   MOVED ZONAL WAVENUMBER LOOP INSIDE
! 1998-12-15  IREDELL  OPENMP DIRECTIVES INSERTED
! 2006-09-10  parrish -- modify to do one zonal wave number only for parallel
!                        computation across processors by zonal wave number
! 2010-05-14  derber  - triangular truncation only
!
! USAGE:    CALL spanaly_ns(M,L,WGT,CLAT,PLN,MP,F,SPC)
!
!   INPUT ARGUMENT LIST:
!     M        - INTEGER SPECTRAL TRUNCATION
!     L        - zonal wave number to process for this call
!     WGT      - REAL GAUSSIAN WEIGHT
!     CLAT     - REAL COSINE OF LATITUDE
!     PLN      - REAL (L:M+MP) LEGENDRE POLYNOMIALS
!     MP       - INTEGER  IDENTIFIERS (0 FOR SCALAR, 1 FOR VECTOR)
!     F        - REAL (2,2) input zonal wave number coefficients for this lat pair
!
!   OUTPUT ARGUMENT LIST:
!     SPC      - REAL (2,L:M+MP) SPECTRAL COEFFICIENTS
!
! ATTRIBUTES:
!   LANGUAGE: CRAY FORTRAN
!
!$$$
      use kinds, only: r_kind
      implicit none

      integer(i_kind),intent(in   ) :: m,l
      integer(i_kind),intent(in   ) :: mp
      real(r_kind)   ,intent(in   ) :: wgt,clat
      real(r_kind)   ,intent(in   ) :: pln(l:m+mp)
      real(r_kind)   ,intent(in   ) :: f(2,2)
      real(r_kind)   ,intent(inout) :: spc(2,l:m+mp)

      real(r_kind) f11,f21,f12,f22
      real(r_kind) wgtloc
      integer(i_kind) n
! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
      if(mp == 0) then
         wgtloc=wgt
      else
         wgtloc=wgt*clat
      end if
      f11=wgtloc*(f(1,1)+f(1,2))
      f21=wgtloc*(f(2,1)+f(2,2))
      f12=wgtloc*(f(1,1)-f(1,2))
      f22=wgtloc*(f(2,1)-f(2,2))
      do n=l,m+mp,2
         spc(1,n)=spc(1,n)+pln(n)*f11
         spc(2,n)=spc(2,n)+pln(n)*f21
      end do
      do n=l+1,m+mp,2
         spc(1,n)=spc(1,n)+pln(n)*f12
         spc(2,n)=spc(2,n)+pln(n)*f22
      end do
! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
      return
      end subroutine spanaly_ns
!-----------------------------------------------------------------------
      subroutine spsynth_ns(m,l,clat,pln,mp,spc,f)
!$$$  SUBPROGRAM DOCUMENTATION BLOCK
!
! SUBPROGRAM:    spsynth_ns  spsynth modified for one zonal wave number triangular trunkation only
!   PRGMMR: IREDELL          ORG: W/NMC23     DATE: 92-10-31
!
! ABSTRACT: SYNTHESIZES FOURIER COEFFICIENTS FROM SPECTRAL COEFFICIENTS
!           FOR A LATITUDE PAIR (NORTHERN AND SOUTHERN HEMISPHERES).
!           VECTOR COMPONENTS ARE DIVIDED BY COSINE OF LATITUDE.
!
! PROGRAM HISTORY LOG:
!   91-10-31  MARK IREDELL
! 1998-12-18  MARK IREDELL  INCLUDE SCALAR AND GRADIENT OPTION
! 2006-09-06  parrish -- modify to do one zonal wave number only for parallel
!                        computation across processors by zonal wave number
! 2010-05-14  derber - modify of triangular truncation only
!
! USAGE:    CALL spsynth_ns(M,L,CLAT,PLN,MP,SPC,F)
!
!   INPUT ARGUMENT LIST:
!     M        - INTEGER SPECTRAL TRUNCATION
!     L        - zonal wave number to process for this call
!     CLAT     - REAL COSINE OF LATITUDE
!     PLN      - REAL (L:M+MP) LEGENDRE POLYNOMIAL
!     SPC      - REAL (2,L:M+MP) SPECTRAL COEFFICIENTS
!     MP       - INTEGER  IDENTIFIERS (0 FOR SCALAR, 1 FOR VECTOR)
!
!   OUTPUT ARGUMENT LIST:
!     F        - REAL (2,2) zonal wave number for this LATITUDE PAIR
!
! ATTRIBUTES:
!   LANGUAGE: CRAY FORTRAN
!
!$$$
      use kinds, only: r_kind
      use constants, only: zero
      implicit none
 
      integer(i_kind),intent(in   ) :: m,l
      integer(i_kind),intent(in   ) :: mp
      real(r_kind)   ,intent(in   ) :: clat
      real(r_kind)   ,intent(in   ) :: pln(l:m+mp)
      real(r_kind)   ,intent(in   ) :: spc(2,l:m+mp)
      real(r_kind)   ,intent(  out) :: f(2,2)

      integer(i_kind) n
      real(r_kind) f1r,f1i,f2r,f2i
! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
!  ZERO OUT FOURIER COEFFICIENTS.
      f1r=zero
      f1i=zero
      f2r=zero
      f2i=zero
! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
!  SYNTHESIS OVER FINITE LATITUDE.
!  FOR EACH ZONAL WAVENUMBER, SYNTHESIZE TERMS OVER TOTAL WAVENUMBER.
!  SYNTHESIZE EVEN AND ODD POLYNOMIALS SEPARATELY.
      do n=l,m+mp,2
         f1r=f1r+pln(n)*spc(1,n)
         f1i=f1i+pln(n)*spc(2,n)
      enddo
      do n=l+1,m+mp,2
         f2r=f2r+pln(n)*spc(1,n)
         f2i=f2i+pln(n)*spc(2,n)
      enddo
!  SEPARATE FOURIER COEFFICIENTS FROM EACH HEMISPHERE.
!  ODD POLYNOMIALS CONTRIBUTE NEGATIVELY TO THE SOUTHERN HEMISPHERE.
!  DIVIDE VECTOR COMPONENTS BY COSINE LATITUDE.
      if(mp == 1) then
         f(1,1)=(f1r+f2r)/clat
         f(2,1)=(f1i+f2i)/clat
         f(1,2)=(f1r-f2r)/clat
         f(2,2)=(f1i-f2i)/clat
      else
         f(1,1)=f1r+f2r
         f(2,1)=f1i+f2i
         f(1,2)=f1r-f2r
         f(2,2)=f1i-f2i
      endif
! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
      return
      end subroutine spsynth_ns

end module strong_fast_global_mod