#include "swancpp.h"
      MODULE WAVES_COUPLER_MOD

#if defined COAWST_COUPLING && defined MCT_LIB
!
!svn $Id: waves_coupler.F 756 2008-09-14 20:18:28Z jcwarner $
!==================================================== John C. Warner ===
!  Copyright (c) 2002-2008 The ROMS/TOMS Group      Hernan G. Arango   !
!   Licensed under a MIT/X style license                               !
!   See License_ROMS.txt                                               !
!=======================================================================
!                                                                      !
!  This module is used to communicate and exchange data between SWAN   !
!  other coupled model(s) using the Model Coupling Toolkit (MCT).      !
!                                                                      !
!=======================================================================
!
!  Componenet model registry.
!
      USE m_MCTWorld, ONLY : MCTWorld_init => init
      USE m_MCTWorld, ONLY : MCTWorld_clean => clean
!
!  Domain decompositin descriptor datatype and assocoiated methods.
!
      USE m_GlobalSegMap, ONLY : GlobalSegMap
      USE m_GlobalSegMap, ONLY : GlobalSegMap_init => init
      USE m_GlobalSegMap, ONLY : GlobalSegMap_lsize => lsize
      USE m_GlobalSegMap, ONLY : GlobalSegMap_clean => clean
      USE m_GlobalSegMap, ONLY : GlobalSegMap_Ordpnts => OrderedPoints
!
!  Field storage data types and associated methods.
!
      USE m_AttrVect, ONLY : AttrVect
      USE m_AttrVect, ONLY : AttrVect_init => init
      USE m_AttrVect, ONLY : AttrVect_zero => zero
      USE m_AttrVect, ONLY : AttrVect_clean => clean
      USE m_AttrVect, ONLY : AttrVect_indxR => indexRA
      USE m_AttrVect, ONLY : AttrVect_importRAttr => importRAttr
      USE m_AttrVect, ONLY : AttrVect_exportRAttr => exportRAttr
!
!  Intercomponent communitcations scheduler.
!
      USE m_Router, ONLY : Router
      USE m_Router, ONLY : Router_init => init
      USE m_Router, ONLY : Router_clean => clean
!
!  Intercomponent transfer.
!
      USE m_Transfer, ONLY : MCT_isend => isend
      USE m_Transfer, ONLY : MCT_irecv => irecv
      USE m_Transfer, ONLY : MCT_waitr => waitrecv
      USE m_Transfer, ONLY : MCT_waits => waitsend
!
# if defined MCT_INTERP_WV2AT
!
!  Sparse Matrix DataType and associated methods.
!
      USE m_SparseMatrix, ONLY : SparseMatrix
      USE m_SparseMatrix, ONLY : SparseMatrix_init => init
      USE m_SparseMatrix, ONLY : SparseMatrix_importGRowInd =>          &
     &                           importGlobalRowIndices
      USE m_SparseMatrix, ONLY : SparseMatrix_importGColInd =>          &
     &                           importGlobalColumnIndices
      USE m_SparseMatrix, ONLY : SparseMatrix_importMatrixElts =>       &
     &                           importMatrixElements
      USE m_SparseMatrix, only : SparseMatrix_lsize => lsize
      USE m_SparseMatrix, only : SparseMatrix_clean => clean
      USE m_SparseMatrixPlus, ONLY : SparseMatrixPlus
      USE m_SparseMatrixPlus, ONLY : SparseMatrixPlus_init => init
      USE m_SparseMatrixPlus, ONLY : SparseMatrixPlus_clean => clean
!
!  Decompose matrix by row.
!
      USE m_SparseMatrixPlus, ONLY : Xonly
!
!  Matrix-Vector multiply methods.
!
      USE m_MatAttrVectMul, ONLY : MCT_MatVecMul => sMatAvMult
# endif

      implicit none
!
      PRIVATE

      PUBLIC :: initialize_wav_coupling
      PUBLIC :: initialize_wav_routers
# ifdef ROMS_COUPLING
      PUBLIC :: wav2ocn_coupling
      PUBLIC :: wavfocn_coupling
# endif
# ifdef WRF_COUPLING
      PUBLIC :: wav2atm_coupling
      PUBLIC :: wavfatm_coupling
# endif
      PUBLIC :: finalize_wav_coupling
!
!  Declarations.
!
      TYPE T_GlobalSegMap_G
        TYPE(GlobalSegMap) :: GSMapSWAN         ! GloabalSegMap variables
      END TYPE T_GlobalSegMap_G
      TYPE (T_GlobalSegMap_G), ALLOCATABLE :: GlobalSegMap_G(:)

      TYPE T_AttrVect_G
# ifdef ROMS_COUPLING
        TYPE(AttrVect) :: wav2ocn_AV            ! AttrVect variables
        TYPE(AttrVect) :: ocn2wav_AV
# endif
# ifdef WRF_COUPLING
        TYPE(AttrVect) :: atm2wav_AV            ! AttrVec variables
        TYPE(AttrVect) :: wav2atm_AV            ! AttrVec variables
# endif
      END TYPE T_AttrVect_G
      TYPE (T_AttrVect_G), ALLOCATABLE :: AttrVect_G(:)

# ifdef ROMS_COUPLING
      TYPE T_Router_O
        type(Router)   :: SWANtoROMS            ! Router variables
      END TYPE T_Router_O
      TYPE (T_Router_O), ALLOCATABLE :: Router_O(:,:)
# endif

# ifdef WRF_COUPLING
      TYPE T_GSMapInterp_A
        TYPE(GlobalSegMap) :: GSMapWRF        ! GloabalSegMap variables
      END TYPE T_GSMapInterp_A
      TYPE (T_GSMapInterp_A), ALLOCATABLE :: GSMapInterp_A(:,:)

      TYPE T_Router_A
        type(Router)   :: SWANtoWRF           ! Router variables
      END TYPE T_Router_A
      TYPE (T_Router_A), ALLOCATABLE :: Router_A(:,:)
# endif

# ifdef MCT_INTERP_WV2AT
      TYPE T_AV2_A
        TYPE(AttrVect) :: wav2atm_AV2         ! AttrVect variables
        TYPE(AttrVect) :: atm2wav_AV2 
      END TYPE T_AV2_A
      TYPE (T_AV2_A), ALLOCATABLE :: AV2_A(:,:)

      TYPE(SparseMatrix) :: sMatW             ! Sparse matrix elements
      TYPE(SparseMatrix) :: sMatA             ! Sparse matrix elements
      TYPE T_SMPlus_G
        TYPE(SparseMatrixPlus) :: A2WMatPlus    ! Sparse matrix plus elements
        TYPE(SparseMatrixPlus) :: W2AMatPlus    ! Sparse matrix plus elements
      END TYPE T_SMPlus_G
      TYPE (T_SMPlus_G), ALLOCATABLE :: SMPlus_G(:,:)
# endif

      CONTAINS

      SUBROUTINE INITIALIZE_WAV_COUPLING (ng)
!
!=======================================================================
!                                                                      !
!  Initialize waves and ocean models coupling stream.  This is the     !
!  training phase use to constuct  MCT  parallel interpolators and     !
!  stablish communication patterns.                                    !
!                                                                      !
!=======================================================================
!
      USE OCPCOMM4
      USE SWCOMM3
      USE M_GENARR
      USE M_PARALL
      USE swan_iounits
      USE mct_coupler_params
# ifdef MCT_INTERP_WV2AT
      USE mod_coupler_iounits
# endif
!
      include 'mpif.h'
      integer, intent(in) :: ng
!
!  Local variable declarations.
!
      integer :: MyError, MyRank
      integer :: gsmsize, nprocs
      integer :: i, j, io, ia, Isize, Jsize, Asize
      integer :: nRows, nCols, num_sparse_elems
      integer :: cid, cad
      character (len=70)  :: nc_name
      character (len=20)  :: to_add
      character (len=120) :: wostring
      character (len=120) :: owstring

      real :: cff

!      integer, dimension(2) :: src_grid_dims, dst_grid_dims
      integer, allocatable :: start(:), length(:)
!
!-----------------------------------------------------------------------
!  Begin initialization phase.
!-----------------------------------------------------------------------
!
!  Get communicator local rank and size.
!
      CALL mpi_comm_rank (WAV_COMM_WORLD, MyRank, MyError)
      CALL mpi_comm_size (WAV_COMM_WORLD, nprocs, MyError)
!
!  Initialize MCT coupled model registry.
!
      IF (ng.eq.1) THEN
        ALLOCATE(GlobalSegMap_G(Nwav_grids))
        ALLOCATE(AttrVect_G(Nwav_grids))
# ifdef MCT_INTERP_WV2AT
        ALLOCATE(SMPlus_G(Nwav_grids,Natm_grids))
        ALLOCATE(AV2_A(Nwav_grids,Natm_grids))
        ALLOCATE(GSMapInterp_A(Nwav_grids,Natm_grids))
# endif
      END IF
!
      WAVid=wavids(ng)
      IF (Nwav_grids.gt.1) THEN
        CALL MCTWorld_init (N_mctmodels, MPI_COMM_WORLD,                &
     &                      WAV_COMM_WORLD,myids=wavids)
      ELSE
        CALL MCTWorld_init (N_mctmodels, MPI_COMM_WORLD,                &
     &                      WAV_COMM_WORLD,WAVid)
      END IF
!
!  Initialize a Global Segment Map for non-haloed transfer of data for
!  SWAN. Determine non-haloed start and length arrays for this
!  processor.
!
      IF (nprocs.eq.1) THEN
        Isize=MXCGL
        Jsize=MYCGL
      ELSE
        IF (MXCGL.gt.MYCGL) THEN
          Isize=MXC-IHALOX*IBLKAD(1)
          Jsize=MYC
        ELSE
          Isize=MXC
          Jsize=MYC-IHALOY*IBLKAD(1)
        END IF
      END IF
!
      allocate ( start(Jsize) )
      allocate ( length(Jsize) )
!
      DO j=1,Jsize
        length(j)=Isize
        IF (MXCGL.gt.MYCGL) THEN
          IF (MyRank.eq.0) THEN
            start(j)=MXF+(j-1)*MXCGL
          ELSE
            start(j)=MXF+(j-1)*MXCGL+IHALOX
          END IF
        ELSE
          IF (MyRank.eq.0) THEN
            start(j)=MYF+(j-1)*MXCGL
          ELSE
            start(j)=(MYF+IHALOY-1)*MXCGL+1+(j-1)*MXCGL
          END IF
        END IF
      END DO
      gsmsize=Isize*Jsize
!
      CALL GlobalSegMap_init (GlobalSegMap_G(ng)%GSMapSWAN, start,      &
     &                        length, 0, WAV_COMM_WORLD, WAVid)
      deallocate (start)
      deallocate (length)

# ifdef MCT_INTERP_WV2AT
!
!  If wave grid and atm grids are different sizes, then
!  develop sparse matrices for interpolation.
!
  35  FORMAT(a3,i1,a7,i1,a11)
      DO ia=1,Natm_grids
!!!!!!!!!!!!!!!!!!!!!!
! First work on atm to wave.
!!!!!!!!!!!!!!!!!!!!!!
!
        IF (INODE.eq.MASTER) THEN
          IF (scrip_opt.eq.1) THEN
            write(nc_name,35) 'atm',ia,'_to_wav',ng,'_weights.nc'
          ELSE
            nc_name=A2Wname(ia,ng)
          END IF
          call get_sparse_matrix (ng, nc_name, num_sparse_elems,        &
     &                            src_grid_dims, dst_grid_dims)
!
! Init the sparse matrix.
!
          nRows=dst_grid_dims(1)*dst_grid_dims(2)
          nCols=src_grid_dims(1)*src_grid_dims(2)
!
! Create sparse matrix.
!
!         Sparse rows is the dst address. Multiply the interp weights
!         by the dst masking.
!
          DO i=1,num_sparse_elems
            j=sparse_rows(i)
            cff=REAL(dst_grid_imask(j),m8)
            sparse_weights(i)=sparse_weights(i)*cff
          END DO

          call SparseMatrix_init(sMatA,nRows,nCols,num_sparse_elems)
          call SparseMatrix_importGRowInd(sMatA, sparse_rows,           &
     &                                    size(sparse_rows))
          call SparseMatrix_importGColInd(sMatA, sparse_cols,           &
     &                                    size(sparse_cols))
          call SparseMatrix_importMatrixElts(sMatA, sparse_weights,     &
     &                                       size(sparse_weights))
!
! Deallocate arrays.
!
          deallocate ( sparse_rows )
          deallocate ( sparse_cols )
          deallocate ( sparse_weights )
          deallocate ( dst_grid_imask )

!!!!!!!!!!!!!!!!!!!!!!
! Second work on waves to atm.
!!!!!!!!!!!!!!!!!!!!!!
!
          IF (scrip_opt.eq.1) THEN
            write(nc_name,35) 'wav',ng,'_to_atm',ia,'_weights.nc'
          ELSE
            nc_name=W2Aname(ng,ia)
          END IF
          call get_sparse_matrix (ng, nc_name, num_sparse_elems,        &
     &                            src_grid_dims, dst_grid_dims)
!
! Init the sparse matrix.
!
          nRows=dst_grid_dims(1)*dst_grid_dims(2)
          nCols=src_grid_dims(1)*src_grid_dims(2)
!
! Create sparse matrix.
!
          DO i=1,num_sparse_elems
            j=sparse_rows(i)
            cff=REAL(dst_grid_imask(j),m8)
            sparse_weights(i)=sparse_weights(i)*cff
          END DO
!
! Load the dst grid as a coupling mask.
!
          allocate(W2A_CPLMASK(ng,ia)%dst_mask(nRows))
          DO i=1,nRows
            W2A_CPLMASK(ng,ia)%dst_mask(i)=dst_grid_imask(i)
          END DO
!
          call SparseMatrix_init(sMatW,nRows,nCols,num_sparse_elems)
          call SparseMatrix_importGRowInd(sMatW, sparse_rows,           &
     &                                    size(sparse_rows))
          call SparseMatrix_importGColInd(sMatW, sparse_cols,           &
     &                                    size(sparse_cols))
          call SparseMatrix_importMatrixElts(sMatW, sparse_weights,     &
     &                                       size(sparse_weights))
!
! Deallocate arrays.
!
          deallocate ( sparse_rows )
          deallocate ( sparse_cols )
          deallocate ( sparse_weights )
          deallocate ( dst_grid_imask )
        END IF
        CALL mpi_bcast(dst_grid_dims, 2, MPI_INTEGER, 0,                &
     &                 WAV_COMM_WORLD, MyError)
!
! scatter dst_grid_imask to be used as cpl_mask
!
        IF (INODE.ne.MASTER) THEN
          nRows=dst_grid_dims(1)*dst_grid_dims(2)
          allocate(W2A_CPLMASK(ng,ia)%dst_mask(nRows))
        END IF
        CALL mpi_bcast(W2A_CPLMASK(ng,ia)%dst_mask, nRows,              &
     &                 MPI_INTEGER, 0,                                  &
     &                 WAV_COMM_WORLD, MyError)
!
!  Initialize a Global Segment Map for non-haloed transfer of data
!  for the atmosphere model.
!
        Isize=INT(dst_grid_dims(1)/nprocs)
        IF (MyRank.eq.nprocs-1) THEN
          Isize=dst_grid_dims(1)-Isize*(nprocs-1)
        ENDIF
        IF (.not.allocated(start)) THEN
          allocate ( start(1) )
        END IF
        IF (.not.allocated(length)) THEN
          allocate ( length(1) )
        END IF
        start=(MyRank*INT(dst_grid_dims(1)/nprocs))*dst_grid_dims(2)+1
        length=Isize*dst_grid_dims(2)
!
        CALL GlobalSegMap_init (GSMapInterp_A(ng,ia)%GSMapWRF,          &
     &                          start, length, 0, WAV_COMM_WORLD, WAVid)
        deallocate (start)
        deallocate (length)
        call mpi_barrier(WAV_COMM_WORLD, MyError)
!
! Create ATM sparse matrix plus for interpolation.
! Specify matrix decomposition to be by row.
!
        call SparseMatrixPlus_init(SMPlus_G(ng,ia)%A2WMatPlus, sMatA,   &
     &                             GSMapInterp_A(ng,ia)%GSMapWRF,       &
     &                             GlobalSegMap_G(ng)%GSMapSWAN,        &
     &                             Xonly,0,WAV_COMM_WORLD, WAVid)
        call SparseMatrix_clean(sMatA)
!
! Create Wave sparse matrix plus for interpolation.
! Specify matrix decomposition to be by row.
!
         call SparseMatrixPlus_init(SMPlus_G(ng,ia)%W2AMatPlus, sMatW,  &
     &                              GlobalSegMap_G(ng)%GSMapSWAN,       &
     &                              GSMapInterp_A(ng,ia)%GSMapWRF,      &
     &                              Xonly,0,WAV_COMM_WORLD, WAVid)
        call SparseMatrix_clean(sMatW)
      END DO
# endif
!
!  Initialize attribute vector holding the export data code strings of
!  the wave model.
!
      cad=LEN(wostring)
      DO i=1,cad
        wostring(i:i)=''
      END DO
      cid=1
!
# ifdef ROMS_COUPLING
      to_add='DISBOT'
      cad=LEN_TRIM(to_add)
      write(wostring(cid:cid+cad-1),'(a)') to_add(1:cad)
      cid=cid+cad
!
      to_add=':DISSURF'
      cad=LEN_TRIM(to_add)
      write(wostring(cid:cid+cad-1),'(a)') to_add(1:cad)
      cid=cid+cad
!
      to_add=':DISWCAP'
      cad=LEN_TRIM(to_add)
      write(wostring(cid:cid+cad-1),'(a)') to_add(1:cad)
      cid=cid+cad
!
      to_add=':HSIGN'
      cad=LEN_TRIM(to_add)
      write(wostring(cid:cid+cad-1),'(a)') to_add(1:cad)
      cid=cid+cad
!
      to_add=':RTP'
      cad=LEN_TRIM(to_add)
      write(wostring(cid:cid+cad-1),'(a)') to_add(1:cad)
      cid=cid+cad
!
      to_add=':TMBOT'
      cad=LEN_TRIM(to_add)
      write(wostring(cid:cid+cad-1),'(a)') to_add(1:cad)
      cid=cid+cad
!
      to_add=':UBOT'
      cad=LEN_TRIM(to_add)
      write(wostring(cid:cid+cad-1),'(a)') to_add(1:cad)
      cid=cid+cad
!
      to_add=':DIRE'
      cad=LEN_TRIM(to_add)
      write(wostring(cid:cid+cad-1),'(a)') to_add(1:cad)
      cid=cid+cad

      to_add=':DIRN'
      cad=LEN_TRIM(to_add)
      write(wostring(cid:cid+cad-1),'(a)') to_add(1:cad)
      cid=cid+cad
!
      to_add=':DIREP'
      cad=LEN_TRIM(to_add)
      write(wostring(cid:cid+cad-1),'(a)') to_add(1:cad)
      cid=cid+cad

      to_add=':DIRNP'
      cad=LEN_TRIM(to_add)
      write(wostring(cid:cid+cad-1),'(a)') to_add(1:cad)
      cid=cid+cad
!
      to_add=':WLEN'
      cad=LEN_TRIM(to_add)
      write(wostring(cid:cid+cad-1),'(a)') to_add(1:cad)
      cid=cid+cad
!
      to_add=':WLENP'
      cad=LEN_TRIM(to_add)
      write(wostring(cid:cid+cad-1),'(a)') to_add(1:cad)
      cid=cid+cad
!
      to_add=':QB'
      cad=LEN_TRIM(to_add)
      write(wostring(cid:cid+cad-1),'(a)') to_add(1:cad)
      cid=cid+cad
!
      to_add=':WDSPR'
      cad=LEN_TRIM(to_add)
      write(wostring(cid:cid+cad-1),'(a)') to_add(1:cad)
      cid=cid+cad
!
      to_add=':WQP'
      cad=LEN_TRIM(to_add)
      write(wostring(cid:cid+cad-1),'(a)') to_add(1:cad)
      cid=cid+cad
!
#if defined VEGETATION && defined VEG_SWAN_COUPLING \
      && defined VEG_STREAMING
      to_add=':DISVEG'
      cad=LEN_TRIM(to_add)
      write(wostring(cid:cid+cad-1),'(a)') to_add(1:cad)
      cid=cid+cad
#endif  
!
!  Finalize and remove trailing spaces from the wostring
!  for the rlist.
!
      cad=LEN_TRIM(wostring)
      wostring=wostring(1:cad)
!
      CALL AttrVect_init(AttrVect_G(ng)%wav2ocn_AV,                     &
     &                   rList=TRIM(wostring),lsize=gsmsize)
      CALL AttrVect_zero(AttrVect_G(ng)%wav2ocn_AV)
!
!  Initialize attribute vector holding the export data code string of
!  the ocean model.
!
      cad=LEN(owstring)
      DO i=1,cad
        owstring(i:i)=''
      END DO
      cid=1
!
      to_add='DEPTH'
      cad=LEN_TRIM(to_add)
      write(owstring(cid:cid+cad-1),'(a)') to_add(1:cad)
      cid=cid+cad
!
      to_add=':WLEV'
      cad=LEN_TRIM(to_add)
      write(owstring(cid:cid+cad-1),'(a)') to_add(1:cad)
      cid=cid+cad
!
      to_add=':VELX'
      cad=LEN_TRIM(to_add)
      write(owstring(cid:cid+cad-1),'(a)') to_add(1:cad)
      cid=cid+cad
!
      to_add=':VELY'
      cad=LEN_TRIM(to_add)
      write(owstring(cid:cid+cad-1),'(a)') to_add(1:cad)
      cid=cid+cad
!
      to_add=':ZO'
      cad=LEN_TRIM(to_add)
      write(owstring(cid:cid+cad-1),'(a)') to_add(1:cad)
      cid=cid+cad
#if defined VEGETATION && defined VEG_SWAN_COUPLING
!
      to_add=':VEGDENS'
      cad=LEN_TRIM(to_add)
      write(owstring(cid:cid+cad-1),'(a)') to_add(1:cad)
      cid=cid+cad
#endif
#  ifdef ICE_MODEL
!
      to_add=':SEAICE'
      cad=LEN_TRIM(to_add)
      write(owstring(cid:cid+cad-1),'(a)') to_add(1:cad)
      cid=cid+cad
#  endif
!
!  Finalize and remove trailing spaces from the owstring
!  for the rlist.
!
      cad=LEN_TRIM(owstring)
      owstring=owstring(1:cad)
!
      CALL AttrVect_init (AttrVect_G(ng)%ocn2wav_AV,                    &
     &                    rList=TRIM(owstring),lsize=gsmsize)
      CALL AttrVect_zero (AttrVect_G(ng)%ocn2wav_AV)
!
# endif
# if defined WRF_COUPLING
!
!  Initialize attribute vector holding the export data code string of
!  the atmosphere model.
!
      Asize=GlobalSegMap_lsize(GlobalSegMap_G(ng)%GSMapSWAN,            &
     &                         WAV_COMM_WORLD)
      CALL AttrVect_init (AttrVect_G(ng)%atm2wav_AV, rlist="U10:V10",   &
     &                    lsize=Asize)
      CALL AttrVect_zero (AttrVect_G(ng)%atm2wav_AV)
!
! Initialize atribute vector holding wave data to atm.
!
      CALL AttrVect_init(AttrVect_G(ng)%wav2atm_AV,                     &
     &                   rList="HSIGN:WLENP:RTP",                       &
     &                   lsize=Asize)
      CALL AttrVect_zero(AttrVect_G(ng)%wav2atm_AV)
!
#  if defined MCT_INTERP_WV2AT
      DO ia=1,Natm_grids
!  Initialize attribute vector holding the export data code strings of
!  the atm model. The Asize is the number of grid point on this
!  processor.
!
        Asize=GlobalSegMap_lsize(GSMapInterp_A(ng,ia)%GSMapWRF,         &
     &                           WAV_COMM_WORLD)
        CALL AttrVect_init (AV2_A(ng,ia)%atm2wav_AV2, rlist="U10:V10",  &
     &                      lsize=Asize)
        CALL AttrVect_zero (AV2_A(ng,ia)%atm2wav_AV2)
!
!  Initialize attribute vector holding the export data code string of
!  the wave model.
!
        CALL AttrVect_init(AV2_A(ng,ia)%wav2atm_AV2,                    &
     &                     rList="HSIGN:WLENP:RTP:CPL_MASK",            &
     &                     lsize=Asize)
        CALL AttrVect_zero(AV2_A(ng,ia)%wav2atm_AV2)
      END DO
#  endif
# endif

      RETURN
      END SUBROUTINE INITIALIZE_WAV_COUPLING

      SUBROUTINE INITIALIZE_WAV_ROUTERS
!
!=======================================================================
!                                                                      !
!  Initialize waves routers for wave model.                            !
!                                                                      !
!=======================================================================
!
      USE mct_coupler_params
      USE M_PARALL
!
!      include 'mpif.h'
!
!  Local variable declarations.
!
      integer :: MyError, MyRank
      integer :: ng, iw, ia
!
!  Initialize MCT Routers.
!
# ifdef ROMS_COUPLING
      ALLOCATE(Router_O(Nwav_grids,Nocn_grids))
!
!  Initialize a router to the ocean model component.
!
      DO ng=1,Nocn_grids
        DO iw=1,Nwav_grids
          OCNid=ocnids(ng)
          CALL Router_init (OCNid, GlobalSegMap_G(iw)%GSMapSWAN,        &
     &                      WAV_COMM_WORLD, Router_O(iw,ng)%SWANtoROMS)
        END DO
      END DO
# endif
# ifdef WRF_COUPLING
      ALLOCATE(Router_A(Nwav_grids,Natm_grids))
!
!  Initialize a router to the ocean model component.
!
      DO ia=1,Natm_grids
        DO iw=1,Nwav_grids
          ATMid=atmids(ia)
# if defined MCT_INTERP_WV2AT
          CALL Router_init (ATMid, GSMapInterp_A(iw,ia)%GSMapWRF,       &
     &                      WAV_COMM_WORLD, Router_A(iw,ia)%SWANtoWRF)
# else
          CALL Router_init (ATMid, GlobalSegMap_G(iw)%GSMapSWAN,        &
     &                      WAV_COMM_WORLD, Router_A(iw,ia)%SWANtoWRF)
# endif
        END DO
      END DO
# endif

      RETURN
      END SUBROUTINE INITIALIZE_WAV_ROUTERS

# ifdef ROMS_COUPLING
      SUBROUTINE WAV2OCN_COUPLING (MIP, NVOQP, VOQR, VOQ, IRQ,          &
     &                             IVTYPE, COMPDA, Numcouple, ng, io)
!
!=======================================================================
!                                                                      !
!  This subroutine reads and writes the coupling data streams between  !
!  ocean and wave models. Currently, the following data streams are    !
!  processed:                                                          !
!                                                                      !
!  Fields exported to the OCEAN model:                                 !
!                                                                      !
!     * Wave direction mean (degrees)                                  !
!     * Wave direction peak (degrees)                                  !
!     * Significant wave height (m)                                    !
!     * Average wave length (m)                                        !
!     * Surface wave relative peak period (s)                          !
!     * Bottom wave period (s)                                         !
!     * Percent of breakig waves (nondimensional)                      !
!     * Wave energy dissipation (W/m2)                                 !
!     * Wave bottom orbital velocity (m/s)                             !
!                                                                      !
!  Fields imported from the OCEAN Model:                               !
!                                                                      !
!     * Bathymetry, bottom elevation (m)                               !
!     * Free-surface, water surface elevation (m)                      !
!     * Depth integrated u-momentum (m/s)                              !
!     * Depth integrated v-momentum (m/s)                              !
!                                                                      !
!=======================================================================
!
      USE SWCOMM1
      USE SWCOMM3
      USE SWCOMM4
      USE OUTP_DATA
      USE M_PARALL
      USE M_GENARR
      USE M_MPI
      USE OCPCOMM4
!     USE mod_coupler_kinds
# ifdef COAWST_COUPLING
      USE mct_coupler_params
# endif
!
      implicit none
!
!  Imported variable declarations.
!
      integer :: MIP, IRQ, nvoqp, Numcouple, ng, io
      integer :: VOQR(NMOVAR), IVTYPE, IP, IX, IY

      real :: COMPDA(MCGRD,MCMVAR)
      real :: VOQ(MIP,NVOQP)
!
!  Local variable declarations.
!
      integer :: MyStatus, MyError, MySize, MyRank
      integer :: i, id, j, gsmsize, ierr, indx, Tag
      integer :: Istr, Iend, Jstr, Jend
      integer :: Isize, Jsize, INDXG, NPROCS, OFFSET
      integer, pointer :: points(:)

      real(m8), pointer :: avdata(:)
      real(m8), pointer :: DIRE(:)
      real(m8), pointer :: DIRN(:)
      real(m8), pointer :: DIREP(:)
      real(m8), pointer :: DIRNP(:)
      character (len=40) :: code
!
!-----------------------------------------------------------------------
!  Send wave fields to ROMS.
!-----------------------------------------------------------------------
!
      CALL MPI_COMM_RANK (WAV_COMM_WORLD, MyRank, MyError)
      CALL MPI_COMM_SIZE (WAV_COMM_WORLD, NPROCS, MyError)
!
!  Get the number of grid point on this processor.
!
      gsmsize=GlobalSegMap_lsize(GlobalSegMap_G(ng)%GSMapSWAN,          &
     &                           WAV_COMM_WORLD)
!
!  Allocate attribute vector array used to export/import data.
!
      allocate ( avdata(gsmsize),stat=ierr )
      allocate ( DIRE(gsmsize),stat=ierr )
      allocate ( DIRN(gsmsize),stat=ierr )
      allocate ( DIREP(gsmsize),stat=ierr )
      allocate ( DIRNP(gsmsize),stat=ierr )
      avdata=0.0_m8
      DIRE=0.0_m8
      DIRN=0.0_m8
      DIREP=0.0_m8
      DIRNP=0.0_m8
!
!  Ask for points in this tile.
!
      CALL GlobalSegMap_Ordpnts (GlobalSegMap_G(ng)%GSMapSWAN,          &
     &                           MyRank, points)
!
!  Load SWAN exporting data into MCT storage buffers.  Since this
!  routine is called several times from main, only load field
!  according to the IVTYPE flag.  The data is exported using ROMS
!  definition for real kind m8.
!
      DO IP=1,gsmsize
        avdata(IP)=REAL( VOQ(points(IP),VOQR(IVTYPE)),m8 )
      END DO
!
      IF (IVTYPE.eq.54) THEN
        CALL AttrVect_importRAttr (AttrVect_G(ng)%wav2ocn_AV,           &
     &                             "DISBOT",avdata)
      ELSE IF (IVTYPE.eq.55) THEN
        CALL AttrVect_importRAttr (AttrVect_G(ng)%wav2ocn_AV,           &
     &                             "DISSURF",avdata)
      ELSE IF (IVTYPE.eq.56) THEN
        CALL AttrVect_importRAttr (AttrVect_G(ng)%wav2ocn_AV,           &
     &                             "DISWCAP",avdata)
      ELSE IF (IVTYPE.eq.10) THEN
        CALL AttrVect_importRAttr (AttrVect_G(ng)%wav2ocn_AV,           &
     &                             "HSIGN",avdata)
      ELSE IF (IVTYPE.eq.12) THEN
        CALL AttrVect_importRAttr (AttrVect_G(ng)%wav2ocn_AV,           &
     &                             "RTP",avdata)
      ELSE IF (IVTYPE.eq.50) THEN
        CALL AttrVect_importRAttr (AttrVect_G(ng)%wav2ocn_AV,           &
     &                             "TMBOT",avdata)
      ELSE IF (IVTYPE.eq.6) THEN
        CALL AttrVect_importRAttr (AttrVect_G(ng)%wav2ocn_AV,           &
     &                             "UBOT",avdata)
      ELSE IF (IVTYPE.eq.13) THEN
        DO IP=1,gsmsize
          DIRE(IP)=1.0*SIN(avdata(IP)*PI/180.0)
          DIRN(IP)=1.0*COS(avdata(IP)*PI/180.0)
        END DO
        CALL AttrVect_importRAttr (AttrVect_G(ng)%wav2ocn_AV,           &
     &                             "DIRE",DIRE)
        CALL AttrVect_importRAttr (AttrVect_G(ng)%wav2ocn_AV,           &
     &                             "DIRN",DIRN)
      ELSE IF (IVTYPE.eq.14) THEN
        DO IP=1,gsmsize
          DIREP(IP)=1.0*SIN(avdata(IP)*PI/180.0)
          DIRNP(IP)=1.0*COS(avdata(IP)*PI/180.0)
        END DO
        CALL AttrVect_importRAttr (AttrVect_G(ng)%wav2ocn_AV,           &
     &                             "DIREP",DIREP)
        CALL AttrVect_importRAttr (AttrVect_G(ng)%wav2ocn_AV,           &
     &                             "DIRNP",DIRNP)
      ELSE IF (IVTYPE.eq.17) THEN
        CALL AttrVect_importRAttr (AttrVect_G(ng)%wav2ocn_AV,           &
     &                             "WLEN",avdata)
      ELSE IF (IVTYPE.eq.71) THEN
        CALL AttrVect_importRAttr (AttrVect_G(ng)%wav2ocn_AV,           &
     &                             "WLENP",avdata)
      ELSE IF (IVTYPE.eq.8) THEN
        CALL AttrVect_importRAttr (AttrVect_G(ng)%wav2ocn_AV,           &
     &                             "QB",avdata)
      ELSE IF (IVTYPE.eq.16) THEN
        CALL AttrVect_importRAttr (AttrVect_G(ng)%wav2ocn_AV,           &
     &                             "WDSPR",avdata)
      ELSE IF (IVTYPE.eq.58) THEN
        CALL AttrVect_importRAttr (AttrVect_G(ng)%wav2ocn_AV,           &
     &                             "WQP",avdata)
!
#if defined VEGETATION && defined VEG_SWAN_COUPLING \
      && defined VEG_STREAMING      
      ELSE IF (IVTYPE.eq.57) THEN
        CALL AttrVect_importRAttr (AttrVect_G(ng)%wav2ocn_AV,           &
     &                             "DISVEG",avdata)
#endif   
      END IF
!
      IF (IRQ.eq.Numcouple) THEN
!
!-----------------------------------------------------------------------
!  Create a restart file.
!-----------------------------------------------------------------------
!
!jcw    CALL BACKUP (AC2, SPCSIG, SPCDIR, KGRPNT, XCGRID, YCGRID, ng)
!
!-----------------------------------------------------------------------
!  Send wave fields bundle to ocean model, ROMS.
!-----------------------------------------------------------------------
!
# ifdef WAVES_OCEAN
          Tag=io*100+0*10+ng
          CALL MCT_isend (AttrVect_G(ng)%wav2ocn_AV,                    &
     &                   Router_O(ng,io)%SWANtoROMS, Tag)
          CALL MCT_waits (Router_O(ng,io)%SWANtoROMS)
          IF (MyRank.EQ.0) THEN
            WRITE (SCREEN,36)' == SWAN grid ',ng,                        &
     &                       ' sent wave data to ROMS grid ', io
 36         FORMAT (a14,i2,a29,i2)
          END IF
          IF (MyError.ne.0) THEN
            WRITE (*,*)'coupling send fail swancplr, Error= ', MyError
            CALL FINALIZE_WAV_COUPLING(ng)
          END IF
# endif
      END IF
      deallocate (avdata, points, DIRE, DIRN, DIREP, DIRNP)
!
      RETURN
      END SUBROUTINE WAV2OCN_COUPLING
# endif
# ifdef WRF_COUPLING
      SUBROUTINE WAV2ATM_COUPLING (MIP, NVOQP, VOQR, VOQ, IRQ,          &
     &                             IVTYPE, COMPDA, Numcouple, ng, ia)
!
!=======================================================================
!                                                                      !
!  This subroutine reads and writes the coupling data streams between  !
!  atm and wave models. Currently, the following data streams are      !
!  processed:                                                          !
!                                                                      !
!  Fields exported to the ATM model:                                   !
!                                                                      !
!     * Significant wave height (m)                                    !
!     * Surface wave relative peak period (s)                          !
!     * Surface wave length (m)                                        !
!                                                                      !
!  Fields imported from the ATM Model:                                 !
!                                                                      !
!     * Wind Speed (m/s)                                               !
!                                                                      !
!=======================================================================
!
      USE SWCOMM1
      USE SWCOMM3
      USE SWCOMM4
      USE OUTP_DATA
      USE M_PARALL
      USE M_GENARR
      USE M_MPI
      USE OCPCOMM4
!     USE mod_coupler_kinds
# ifdef COAWST_COUPLING
      USE mct_coupler_params
# endif
!
      implicit none
!
!  Imported variable declarations.
!
      integer :: MIP, IRQ, nvoqp, Numcouple, ng, ia
      integer :: VOQR(NMOVAR), IVTYPE, IP, IX, IY

      real :: COMPDA(MCGRD,MCMVAR)
      real :: VOQ(MIP,NVOQP)
!
!  Local variable declarations.
!
      integer :: MyStatus, MyError, MySize, MyRank
      integer :: i, id, j, gsmsize, ierr, indx, Tag
      integer :: Istr, Iend, Jstr, Jend, Asize
      integer :: Isize, Jsize, INDXG, NPROCS, OFFSET

      integer, pointer :: points(:)
      real(m8), pointer :: avdata(:)
      character (len=40) :: code
#ifdef MCT_INTERP_WV2AT
      integer, pointer :: indices(:)
      real(m8), pointer :: Amask(:)
#endif
!
!-----------------------------------------------------------------------
!  Send wave fields to ROMS.
!-----------------------------------------------------------------------
!
      CALL MPI_COMM_RANK (WAV_COMM_WORLD, MyRank, MyError)
      CALL MPI_COMM_SIZE (WAV_COMM_WORLD, NPROCS, MyError)
!
!  Get the number of grid point on this processor.
!
      gsmsize=GlobalSegMap_lsize(GlobalSegMap_G(ng)%GSMapSWAN,          &
     &                           WAV_COMM_WORLD)
!
!  Allocate attribute vector array used to export/import data.
!
      allocate ( avdata(gsmsize),stat=ierr )
      avdata=0.0_m8
!
!  Ask for points in this tile.
!
      CALL GlobalSegMap_Ordpnts (GlobalSegMap_G(ng)%GSMapSWAN,          &
     &                           MyRank, points)
!
!  Load SWAN exporting data into MCT storage buffers.  Since this
!  routine is called several times from main, only load field
!  according to the IVTYPE flag.  The data is exported using ROMS
!  definition for real kind m8.
!
      DO IP=1,gsmsize
        avdata(IP)=REAL( VOQ(points(IP),VOQR(IVTYPE)),m8 )
      END DO
      deallocate (points)
!
      IF (IVTYPE.eq.10) THEN
        CALL AttrVect_importRAttr (AttrVect_G(ng)%wav2atm_AV,           &
     &                             "HSIGN",avdata)
      ELSE IF (IVTYPE.eq.12) THEN
        CALL AttrVect_importRAttr (AttrVect_G(ng)%wav2atm_AV,           &
     &                             "RTP",avdata)
      ELSE IF (IVTYPE.eq.71) THEN
        CALL AttrVect_importRAttr (AttrVect_G(ng)%wav2atm_AV,           &
     &                             "WLENP",avdata)
      END IF
!
      IF (IRQ.eq.Numcouple) THEN
!
!-----------------------------------------------------------------------
!  Create a restart file.
!-----------------------------------------------------------------------
!
!jcw    CALL BACKUP (AC2, SPCSIG, SPCDIR, KGRPNT, XCGRID, YCGRID, ng)
!
!-----------------------------------------------------------------------
!  Send wave fields bundle to atm model.
!-----------------------------------------------------------------------
!
 35     FORMAT (a14,i2,a23,i5)
!
!  Send fields to atmosphere model.
!
        Tag=ia*10+ng
#  ifdef MCT_INTERP_WV2AT
        CALL MCT_MatVecMul(AttrVect_G(ng)%wav2atm_AV,                   &
     &                     SMPlus_G(ng,ia)%W2AMatPlus,                  &
     &                     AV2_A(ng,ia)%wav2atm_AV2)
!
!  Now add in the CPL_MASK before we send it over to wrf.
!  Get the number of grid points on this processor.
!
        Asize=GlobalSegMap_lsize (GSMapInterp_A(ng,ia)%GSMapWRF,        &
     &                            WAV_COMM_WORLD)
        allocate (Amask(Asize))
        Amask=0.0
!
!  Ask for points in this tile.
!
        CALL GlobalSegMap_Ordpnts (GSMapInterp_A(ng,ia)%GSMapWRF,       &
     &                             MyRank, points)
!
!  Load the dst grid cpl mask into the attr vect.
!
        DO i=1,Asize
          Amask(i)=REAL(W2A_CPLMASK(ng,ia)%dst_mask(points(i)),m8)
        END DO
        deallocate (points)
        CALL AttrVect_importRAttr (AV2_A(ng,ia)%wav2atm_AV2, "CPL_MASK",&
     &                             Amask, Asize)
        CALL MCT_isend (AV2_A(ng,ia)%wav2atm_AV2,                       &
     &                  Router_A(ng,ia)%SWANtoWRF, Tag)
#  else
        CALL MCT_isend (AttrVect_G(ng)%wav2atm_AV,                      &
     &                  Router_A(ng,ia)%SWANtoWRF, Tag)
#  endif
        CALL MCT_waits (Router_A(ng,ia)%SWANtoWRF)
        IF (MyRank.EQ.0) THEN
          WRITE (SCREEN,35) '== SWAN grid ',ng,' sent data to WRF grid '&
     &                      ,ia
        END IF
        IF (MyError.ne.0) THEN
          WRITE (*,*) 'coupl fail swancplr, MyStatus= ', MyError
          CALL FINALIZE_WAV_COUPLING(ng)
        END IF
#ifdef MCT_INTERP_WV2AT
        deallocate (Amask)
        if (associated (indices)) then
          deallocate (indices)
        endif
#endif
      END IF
      deallocate (avdata)
!
      RETURN
      END SUBROUTINE WAV2ATM_COUPLING
# endif
# ifdef ROMS_COUPLING
      SUBROUTINE WAVFOCN_COUPLING (COMPDA, ng, io)
!
!=======================================================================
!                                                                      !
!  This subroutine reads and writes the coupling data streams between  !
!  ocean and wave models. Currently, the following data streams are    !
!  processed:                                                          !
!                                                                      !
!  Fields exported to the OCEAN model:                                 !
!                                                                      !
!     * Mean wave direction (degrees)                                  !
!     * Peak wave direction (degrees)                                  !
!     * Significant wave height (m)                                    !
!     * Average wave length (m)                                        !
!     * Surface wave relative peak period (s)                          !
!     * Bottom wave period (s)                                         !
!     * Percent of breakig waves (nondimensional)                      !
!     * Wave energy dissipation (W/m2)                                 !
!     * Wave bottom orbital velocity (m/s)                             !
!                                                                      !
!  Fields imported from the OCEAN Model:                               !
!                                                                      !
!     * Bathymetry, bottom elevation (m)                               !
!     * Free-surface, water surface elevation (m)                      !
!     * Depth integrated u-momentum (m/s)                              !
!     * Depth integrated v-momentum (m/s)                              !
!                                                                      !
!=======================================================================
!
      USE SWCOMM1
      USE SWCOMM3
      USE SWCOMM4
      USE OUTP_DATA
      USE M_PARALL
      USE M_GENARR
      USE M_MPI
      USE OCPCOMM4
!     USE mod_coupler_kinds
# ifdef COAWST_COUPLING
      USE mct_coupler_params
# endif
!
      implicit none
!
!  Imported variable declarations.
!
      integer :: ng, io
      integer :: IP, IX, IY

      real :: COMPDA(MCGRD,MCMVAR)
!
!  Local variable declarations.
!
      integer :: MyStatus, MyError, MySize, MyRank
      integer :: i, id, j, gsmsize, ierr, indx, Tag
      integer :: Istr, Iend, Jstr, Jend
      integer :: Isize, Jsize, INDXG, NPROCS, OFFSET
      integer :: NUMTRANSFER, NNEIGH, HALOSIZE, NUMSENT, INB
      integer :: WHICHWAY, GDEST, GSRC, TAGT, TAGB, TAGR, TAGL
      integer :: TREQUEST,BREQUEST,RREQUEST,LREQUEST,MSIZE
      integer :: iddep, idwlv, idvlx, idvly, idveg
      integer :: idruf, idice
      integer, dimension(MPI_STATUS_SIZE,4) :: status

      real :: cff
# if !defined BBL_MODEL
      real :: cff2
# endif
      real, parameter ::  Large = 1.0E+20
      real, dimension(2) :: range
      real(m8), pointer :: avdata(:)

      real, pointer :: TEMPMCT(:,:)
      real, pointer :: GRECVT(:), GRECVB(:), GRECVR(:), GRECVL(:)
      real, pointer :: GSENDT(:), GSENDB(:), GSENDR(:), GSENDL(:)

      character (len=40) :: code
!
!-----------------------------------------------------------------------
!  Send wave fields to ROMS.
!-----------------------------------------------------------------------
!
      CALL MPI_COMM_RANK (WAV_COMM_WORLD, MyRank, MyError)
      CALL MPI_COMM_SIZE (WAV_COMM_WORLD, NPROCS, MyError)
!
!  Get the number of grid point on this processor.
!
      gsmsize=GlobalSegMap_lsize(GlobalSegMap_G(ng)%GSMapSWAN,          &
     &                           WAV_COMM_WORLD)
!
!  Allocate attribute vector array used to export/import data.
!
      allocate ( avdata(gsmsize),stat=ierr )
      avdata=0.0_m8
!
!-----------------------------------------------------------------------
!  Receive from ROMS: Depth, Water Level, VELX, and VELY.
!-----------------------------------------------------------------------
!
!  Schedule receiving field from ocean model.
!
      Tag=io*100+0*10+ng
      CALL MCT_irecv (AttrVect_G(ng)%ocn2wav_AV,                        &
     &                Router_O(ng,io)%SWANtoROMS, Tag)
!
!     Wait to make sure the OCN data has arrived.
!
      CALL MCT_waitr (AttrVect_G(ng)%ocn2wav_AV,                        &
     &                 Router_O(ng,io)%SWANtoROMS)
!
      IF (MyRank.EQ.0) THEN
        WRITE (SCREEN,35) ' == SWAN grid ',ng,                          &
     &                    ' recv data from ROMS grid ', io
      END IF
      IF (MyError.ne.0) THEN
       WRITE (*,*) 'coupling fail swancplr, MyStatus= ', MyError
        CALL FINALIZE_WAV_COUPLING(ng)
      END IF
 35   FORMAT (a14,i2,a26,i2)
!
! Pass the non-halo data from MCT into tempmct array.
!
        NNEIGH = IBLKAD(1)
        IF (NPROCS.eq.1) THEN
          Istr=1
          Iend=MXC
          Jstr=1
          Jend=MYC
        ELSE
          IF (MXCGL.GT.MYCGL) THEN
            IF (MyRank.eq.0) THEN
              Istr=1
            ELSE
              Istr=IHALOX+1
            END IF
            Isize=MXC-IHALOX*IBLKAD(1)
            Iend=Istr+Isize-1
            Jstr=1
            Jend=MYC
            HALOSIZE=IHALOX*MYC
          ELSE
            IF (MyRank.eq.0) THEN
              Jstr=1
            ELSE
              Jstr=IHALOY+1
            END IF
            Jsize=MYC-IHALOY*IBLKAD(1)
            Jend=Jstr+Jsize-1
            Istr=1
            Iend=MXC
            HALOSIZE=IHALOY*MXC
          END IF
        END IF
!
!  Determine the amount of fields and assign id numbers.
!
        IP=0
# ifdef WAVES_OCEAN
        IP=IP+1
        iddep=IP
!
        IP=IP+1
        idwlv=IP
!
        IP=IP+1
        idvlx=IP
!
        IP=IP+1
        idvly=IP
!
        IP=IP+1
        idruf=IP
!
#  if defined VEGETATION && defined VEG_SWAN_COUPLING
        IP=IP+1
        idveg=IP
#  endif
!
#  ifdef ICE_MODEL
        IP=IP+1
        idice=IP
#  endif
# endif
        NUMTRANSFER=IP
!
 40     FORMAT (a36,1x,2(1pe14.6))
        allocate ( TEMPMCT(MXC*MYC,NUMTRANSFER),stat=ierr )
        TEMPMCT=0.0
# ifdef WAVES_OCEAN
!
!  Bottom elevation.
!
        CALL AttrVect_exportRAttr (AttrVect_G(ng)%ocn2wav_AV,           &
     &                             "DEPTH",avdata,gsmsize)
        range(1)= Large
        range(2)=-Large
        IP=0
        DO IY=Jstr,Jend
          DO IX=Istr,Iend
            IP=IP+1
            INDXG=(IY-1)*MXC+IX
            TEMPMCT(INDXG,iddep)=avdata(IP)
            range(1)=MIN(range(1),REAL(avdata(IP)))
            range(2)=MAX(range(2),REAL(avdata(IP)))
          END DO
        END DO
        CALL SWREDUCE(range(1),1,SWREAL,SWMIN)
        CALL SWREDUCE(range(2),1,SWREAL,SWMAX)
        IF (MyRank.eq.0) THEN
          write(SCREEN,40) 'ROMStoSWAN Min/Max DEPTH   (m):     ',      &
     &                      range(1),range(2)
        END IF
!
!  Water surface elevation.
!
        CALL AttrVect_exportRAttr (AttrVect_G(ng)%ocn2wav_AV,           &
     &                             "WLEV",avdata,gsmsize)
        range(1)= Large
        range(2)=-Large
        IP=0
        DO IY=Jstr,Jend
          DO IX=Istr,Iend
            IP=IP+1
            INDXG=(IY-1)*MXC+IX
            TEMPMCT(INDXG,idwlv)=avdata(IP)
            range(1)=MIN(range(1),REAL(avdata(IP)))
            range(2)=MAX(range(2),REAL(avdata(IP)))
          END DO
        END DO
        CALL SWREDUCE(range(1),1,SWREAL,SWMIN)
        CALL SWREDUCE(range(2),1,SWREAL,SWMAX)
        IF (MyRank.eq.0) THEN
          write(SCREEN,40) 'ROMStoSWAN Min/Max WLEV    (m):     ',      &
     &                      range(1),range(2)
        END IF
!
!  Depth-integrated u-velocity.
!
        CALL AttrVect_exportRAttr(AttrVect_G(ng)%ocn2wav_AV,           &
     &                             "VELX",avdata,gsmsize)
        range(1)= Large
        range(2)=-Large
        IP=0
        DO IY=Jstr,Jend
          DO IX=Istr,Iend
            IP=IP+1
            INDXG=(IY-1)*MXC+IX
            TEMPMCT(INDXG,idvlx)=avdata(IP)
            range(1)=MIN(range(1),REAL(avdata(IP)))
            range(2)=MAX(range(2),REAL(avdata(IP)))
          END DO
        END DO
        CALL SWREDUCE(range(1),1,SWREAL,SWMIN)
        CALL SWREDUCE(range(2),1,SWREAL,SWMAX)
        IF (MyRank.eq.0) THEN
          write(SCREEN,40) 'ROMStoSWAN Min/Max VELX    (ms-1):  ',      &
     &                      range(1),range(2)
        END IF
!
!  Depth-integrated v-velocity.
!
        CALL AttrVect_exportRAttr(AttrVect_G(ng)%ocn2wav_AV,            &
     &                            "VELY",avdata,gsmsize)
        range(1)= Large
        range(2)=-Large
        IP=0
        DO IY=Jstr,Jend
          DO IX=Istr,Iend
            IP=IP+1
            INDXG=(IY-1)*MXC+IX
            TEMPMCT(INDXG,idvly)=avdata(IP)
            range(1)=MIN(range(1),REAL(avdata(IP)))
            range(2)=MAX(range(2),REAL(avdata(IP)))
          END DO
        END DO
        CALL SWREDUCE(range(1),1,SWREAL,SWMIN)
        CALL SWREDUCE(range(2),1,SWREAL,SWMAX)
        IF (MyRank.eq.0) THEN
          write(SCREEN,40) 'ROMStoSWAN Min/Max VELY    (ms-1):  ',      &
     &                      range(1),range(2)
        END IF
!
!  Bottom roughness.
!
        CALL AttrVect_exportRAttr(AttrVect_G(ng)%ocn2wav_AV,           &
     &                            "ZO",avdata,gsmsize)
        range(1)= Large
        range(2)=-Large
        IP=0
        DO IY=Jstr,Jend
          DO IX=Istr,Iend
            IP=IP+1
            INDXG=(IY-1)*MXC+IX
            TEMPMCT(INDXG,idruf)=avdata(IP)
            range(1)=MIN(range(1),REAL(avdata(IP)))
            range(2)=MAX(range(2),REAL(avdata(IP)))
          END DO
        END DO
        CALL SWREDUCE(range(1),1,SWREAL,SWMIN)
        CALL SWREDUCE(range(2),1,SWREAL,SWMAX)
        IF (MyRank.eq.0) THEN
          write(SCREEN,40) 'ROMStoSWAN Min/Max ZO      (m):      ',    &
     &                      range(1),range(2)
        END IF
#  if defined VEGETATION && defined VEG_SWAN_COUPLING
!
!  Vegetation density.
!
        CALL AttrVect_exportRAttr(AttrVect_G(ng)%ocn2wav_AV,           &
     &                            "VEGDENS",avdata,gsmsize)
        range(1)= Large
        range(2)=-Large
        IP=0
        DO IY=Jstr,Jend
          DO IX=Istr,Iend
            IP=IP+1
            INDXG=(IY-1)*MXC+IX
            TEMPMCT(INDXG,idveg)=avdata(IP)                     
            range(1)=MIN(range(1),REAL(avdata(IP)))
            range(2)=MAX(range(2),REAL(avdata(IP)))
          END DO
        END DO
        CALL SWREDUCE(range(1),1,SWREAL,SWMIN)
        CALL SWREDUCE(range(2),1,SWREAL,SWMAX)
        IF (MyRank.eq.0) THEN
          write(SCREEN,40) 'ROMStoSWAN Min/Max VEGDENS (indv/m-2): ',   &
     &                      range(1),range(2)
        END IF
#  endif
#  ifdef ICE_MODEL
!
!  sea ice fraction.
!
        CALL AttrVect_exportRAttr(AttrVect_G(ng)%ocn2wav_AV,           &
     &                            "SEAICE",avdata,gsmsize)
        range(1)= Large
        range(2)=-Large
        IP=0
        DO IY=Jstr,Jend
          DO IX=Istr,Iend
            IP=IP+1
            INDXG=(IY-1)*MXC+IX
            TEMPMCT(INDXG,idice)=avdata(IP)
            range(1)=MIN(range(1),REAL(avdata(IP)))
            range(2)=MAX(range(2),REAL(avdata(IP)))
          END DO
        END DO
        CALL SWREDUCE(range(1),1,SWREAL,SWMIN)
        CALL SWREDUCE(range(2),1,SWREAL,SWMAX)
        IF (MyRank.eq.0) THEN
          write(SCREEN,40) 'ROMStoSWAN Min/Max DEPTH   (-):     ',      &
     &                      range(1),range(2)
        END IF
#  endif
# endif
!
!  Pack and send halo regions to be exchanged with adjacent tiles.
!  IBLKAD contains the tile data.
!  WHICHWAY: [top, bot, right, left] = [1 2 3 4]
!
        IF (NPROCS.GT.1) THEN
          MSIZE=HALOSIZE*NUMTRANSFER
          IF (MXCGL.GT.MYCGL) THEN
            allocate ( GSENDR(MSIZE),stat=ierr )
            allocate ( GSENDL(MSIZE),stat=ierr )
            allocate ( GRECVR(MSIZE),stat=ierr )
            allocate ( GRECVL(MSIZE),stat=ierr )
            GSENDR=0.0
            GSENDL=0.0
            GRECVR=0.0
            GRECVL=0.0
          ELSE
            allocate ( GSENDT(MSIZE),stat=ierr )
            allocate ( GSENDB(MSIZE),stat=ierr )
            allocate ( GRECVT(MSIZE),stat=ierr )
            allocate ( GRECVB(MSIZE),stat=ierr )
            GSENDT=0.0
            GSENDB=0.0
            GRECVT=0.0
            GRECVB=0.0
          END IF
          TAGT=1
          TAGB=2
          TAGR=3
          TAGL=4
          DO INB=1,NNEIGH
            OFFSET=0
            WHICHWAY=IBLKAD(3*INB)
            DO NUMSENT=1,NUMTRANSFER
              IP=OFFSET
              IF (WHICHWAY.EQ.1) THEN
                DO IY=MYC-IHALOX-2,MYC-3
                  DO IX=1,MXC
                    IP=IP+1
                    INDXG=(IY-1)*MXC+IX
                    GSENDT(IP)=TEMPMCT(INDXG,NUMSENT)
                  END DO
                END DO
              ELSE IF (WHICHWAY.EQ.2) THEN
                DO IY=IHALOY+1,IHALOY+3
                  DO IX=1,MXC
                    IP=IP+1
                    INDXG=(IY-1)*MXC+IX
                    GSENDB(IP)=TEMPMCT(INDXG,NUMSENT)
                  END DO
                END DO
              ELSE IF (WHICHWAY.EQ.3) THEN
                DO IY=1,MYC
                  DO IX=MXC-IHALOX-2,MXC-3
                    IP=IP+1
                    INDXG=(IY-1)*MXC+IX
                    GSENDR(IP)=TEMPMCT(INDXG,NUMSENT)
                  END DO
                END DO
              ELSE IF (WHICHWAY.EQ.4) THEN
                DO IY=1,MYC
                  DO IX=IHALOX+1,IHALOX+3
                    IP=IP+1
                    INDXG=(IY-1)*MXC+IX
                    GSENDL(IP)=TEMPMCT(INDXG,NUMSENT)
                  END DO
                END DO
              END IF
              OFFSET=OFFSET+HALOSIZE
            END DO
          END DO
          DO INB=1,NNEIGH
            GSRC=IBLKAD(3*INB-1)-1
            WHICHWAY=IBLKAD(3*INB)
            IF (WHICHWAY.EQ.1) THEN
              CALL mpi_irecv (GRECVT,MSIZE,SWREAL,                      &
     &                        GSRC,TAGB,WAV_COMM_WORLD,TREQUEST,MyError)
            ELSE IF (WHICHWAY.EQ.2) THEN
              CALL mpi_irecv (GRECVB,MSIZE,SWREAL,                      &
     &                        GSRC,TAGT,WAV_COMM_WORLD,BREQUEST,MyError)
            ELSE IF (WHICHWAY.EQ.3) THEN
              CALL mpi_irecv (GRECVR,MSIZE,SWREAL,                      &
     &                        GSRC,TAGL,WAV_COMM_WORLD,RREQUEST,MyError)
            ELSE IF (WHICHWAY.EQ.4) THEN
              CALL mpi_irecv (GRECVL,MSIZE,SWREAL,                      &
     &                        GSRC,TAGR,WAV_COMM_WORLD,LREQUEST,MyError)
            END IF
          END DO
          DO INB=1,NNEIGH
            GDEST=IBLKAD(3*INB-1)-1
            WHICHWAY=IBLKAD(3*INB)
            IF (WHICHWAY.EQ.1) THEN
              CALL mpi_send (GSENDT,MSIZE,SWREAL,                       &
     &                       GDEST,TAGT,WAV_COMM_WORLD,MyError)
            ELSE IF (WHICHWAY.EQ.2) THEN
              CALL mpi_send (GSENDB,MSIZE,SWREAL,                       &
     &                       GDEST,TAGB,WAV_COMM_WORLD,MyError)
            ELSE IF (WHICHWAY.EQ.4) THEN
              CALL mpi_send (GSENDL,MSIZE,SWREAL,                       &
     &                       GDEST,TAGL,WAV_COMM_WORLD,MyError)
            ELSE IF (WHICHWAY.EQ.3) THEN
              CALL mpi_send (GSENDR,MSIZE,SWREAL,                       &
     &                       GDEST,TAGR,WAV_COMM_WORLD,MyError)
            END IF
          END DO
!
! Receive and unpack halo regions exchanged with adjacent tiles.
! [top, bot, right, left] = [1 2 3 4]
!
          DO INB=1,NNEIGH
            WHICHWAY=IBLKAD(3*INB)
            IF (WHICHWAY.EQ.1) THEN
              CALL mpi_wait (TREQUEST,status(1,1),MyError)
            ELSE IF (WHICHWAY.EQ.2) THEN
              CALL mpi_wait (BREQUEST,status(1,2),MyError)
            ELSE IF (WHICHWAY.EQ.3) THEN
              CALL mpi_wait (RREQUEST,status(1,3),MyError)
            ELSE IF (WHICHWAY.EQ.4) THEN
              CALL mpi_wait (LREQUEST,status(1,4),MyError)
            END IF
          END DO
!
          DO INB=1,NNEIGH
            OFFSET=0
            WHICHWAY=IBLKAD(3*INB)
            IF (WHICHWAY.EQ.1) THEN
              DO NUMSENT=1,NUMTRANSFER
                IP=OFFSET
                DO IY=MYC-2,MYC
                  DO IX=1,MXC
                    IP=IP+1
                    INDXG=(IY-1)*MXC+IX
                    TEMPMCT(INDXG,NUMSENT)=GRECVT(IP)
                  END DO
                END DO
                OFFSET=OFFSET+HALOSIZE
              END DO
            ELSE IF (WHICHWAY.EQ.2) THEN
              DO NUMSENT=1,NUMTRANSFER
                IP=OFFSET
                DO IY=1,IHALOY
                  DO IX=1,MXC
                    IP=IP+1
                    INDXG=(IY-1)*MXC+IX
                    TEMPMCT(INDXG,NUMSENT)=GRECVB(IP)
                  END DO
                END DO
                OFFSET=OFFSET+HALOSIZE
              END DO
            ELSE IF (WHICHWAY.EQ.3) THEN
              DO NUMSENT=1,NUMTRANSFER
                IP=OFFSET
                DO IY=1,MYC
                  DO IX=MXC-2,MXC
                    IP=IP+1
                    INDXG=(IY-1)*MXC+IX
                    TEMPMCT(INDXG,NUMSENT)=GRECVR(IP)
                  END DO
                END DO
                OFFSET=OFFSET+HALOSIZE
              END DO
            ELSE IF (WHICHWAY.EQ.4) THEN
              DO NUMSENT=1,NUMTRANSFER
                IP=OFFSET
                DO IY=1,MYC
                  DO IX=1,IHALOX
                    IP=IP+1
                    INDXG=(IY-1)*MXC+IX
                    TEMPMCT(INDXG,NUMSENT)=GRECVL(IP)
                  END DO
                END DO
                OFFSET=OFFSET+HALOSIZE
              END DO
            END IF
          END DO
          IF (MXCGL.GT.MYCGL) THEN
            deallocate (GRECVR,GRECVL,GSENDR,GSENDL)
          ELSE
            deallocate (GRECVT,GRECVB,GSENDT,GSENDB)
          END IF
        END IF
!
! Finally insert the full (MXC*MYC) TEMPMCT array into the SWAN
! array for DEPTH and computational array COMPDA. Only insert
! active (wet points) using array KGRPNT.
!
# if defined WAVES_OCEAN && defined SED_MORPH
! Insert depth into SWAN array.
! Here the depth is really the bottom level. 
! When the user requests depth for output, that value is written
! as COMPDA(INDX,JWLV2) which is DEPTH (really botlev) + watlev.
!
        IP=0
        DO IY=1,MYC
          DO IX=1,MXC
            IP=IP+1
            IF (MXCGL.gt.MYCGL) THEN
              OFFSET=(MXF-1)+(IY-1)*MXCGL+IX
            ELSE
              OFFSET=(MYF-1)*MXCGL+(IY-1)*MXCGL+IX
            END IF
            cff=TEMPMCT(IP,1)
            IF (cff.gt.0.) THEN
              DEPTH(OFFSET)=TEMPMCT(IP,1)
            END IF
          END DO
        END DO
# endif
!
! Move values at 'present' time level 2 to 'old' time level 1.
! MCGRD = MXC*MYC+1-#masked cells.
! MXC = # cells x-dir in this tile including halox.
! MYC = # cells y-dir in this tile including haloy.
! COMPDA has only active wet points + 1.
!
# ifdef WAVES_OCEAN
        IF (io.eq.1) THEN
          DO INDX = 1, MCGRD
            COMPDA(INDX,JWLV1)=COMPDA(INDX,JWLV2)
            COMPDA(INDX,JVX1) =COMPDA(INDX,JVX2)
            COMPDA(INDX,JVY1) =COMPDA(INDX,JVY2)
            COMPDA(INDX,JFRC3)=COMPDA(INDX,JFRC2)
#  if defined VEGETATION && defined VEG_SWAN_COUPLING
            COMPDA(INDX,JNPLA3)=COMPDA(INDX,JNPLA2)
#  endif
          END DO
        END IF
# endif
!
! Insert bot level, water level, velx, vely, fric, and winds 
! into SWAN arrays.
!
! If not using a BBL model, then determine the non-spatially 
! varying friction coeff to enter into the JFRC2 array.
!
# if defined WAVES_OCEAN && !defined BBL_MODEL
        cff2=0.05_m8                  ! default
        IF (IBOT.EQ.1) THEN           ! Jonswap
          cff2=PBOT(3)
        ELSE IF (IBOT.EQ.2) THEN      ! Collins
          cff2=PBOT(2)
        ELSE IF (IBOT.EQ.3) THEN      ! Madsen
          cff2=PBOT(5)
        END IF
# endif
        IP=0
        DO IY=1,MYC
          DO IX=1,MXC
            IP=IP+1
            INDX = KGRPNT(IX,IY)
            IF (INDX.GT.1) THEN
# ifdef WAVES_OCEAN
#  ifdef BBL_MODEL
              IF (io.eq.1) THEN
                COMPDA(INDX,JFRC2)=MAX(REAL(TEMPMCT(IP,idruf)), 0.0001)
              ELSE
                COMPDA(INDX,JFRC2)=COMPDA(INDX,JFRC2)+                  &
     &                             MAX(REAL(TEMPMCT(IP,idruf)), 0.0001)
              END IF
#  else
              COMPDA(INDX,JFRC2)=cff2
#  endif
!             cff=TEMPMCT(IP,iddep)+TEMPMCT(IP,idwlv)
!             IF (cff.gt.0.) THEN
                IF (io.eq.1) THEN
#  ifdef SED_MORPH
                  COMPDA(INDX,JBOTLV)=TEMPMCT(IP,iddep)
#  endif
                  COMPDA(INDX,JWLV2)=TEMPMCT(IP,idwlv)
                  COMPDA(INDX,JVX2)=TEMPMCT(IP,idvlx)
                  COMPDA(INDX,JVY2)=TEMPMCT(IP,idvly)
#  if defined VEGETATION && defined VEG_SWAN_COUPLING 
                  COMPDA(INDX,JNPLA2)=TEMPMCT(IP,idveg)
#  endif
                ELSE
#  ifdef SED_MORPH
                  COMPDA(INDX,JBOTLV)=COMPDA(INDX,JBOTLV)+              &
     &                                TEMPMCT(IP,iddep)
#  endif
                  COMPDA(INDX,JWLV2)=COMPDA(INDX,JWLV2)+                &
     &                               TEMPMCT(IP,idwlv)
                  COMPDA(INDX,JVX2)=COMPDA(INDX,JVX2)+TEMPMCT(IP,idvlx)
                  COMPDA(INDX,JVY2)=COMPDA(INDX,JVY2)+TEMPMCT(IP,idvly)
#  if defined VEGETATION && defined VEG_SWAN_COUPLING
                  COMPDA(INDX,JNPLA2)=COMPDA(INDX,JNPLA2)+              &
     &                                TEMPMCT(IP,idveg)
#  endif

                END IF
!             END IF
# endif
            END IF
          END DO
        END DO
!
      deallocate (TEMPMCT)
      deallocate (avdata)
!
      RETURN
      END SUBROUTINE WAVFOCN_COUPLING
# endif
# ifdef WRF_COUPLING
      SUBROUTINE WAVFATM_COUPLING (COMPDA, ng, ia)
!
!=======================================================================
!                                                                      !
!  This subroutine reads and writes the coupling data streams between  !
!  ocean and wave models. Currently, the following data streams are    !
!  processed:                                                          !
!                                                                      !
!                                                                      !
!  Fields imported from the ATM Model:                                 !
!                                                                      !
!     * Wind E and N (m/s)                                             !
!                                                                      !
!=======================================================================
!
      USE SWCOMM1
      USE SWCOMM3
      USE SWCOMM4
      USE OUTP_DATA
      USE M_PARALL
      USE M_GENARR
      USE M_MPI
      USE OCPCOMM4
!     USE mod_coupler_kinds
# ifdef COAWST_COUPLING
      USE mct_coupler_params
# endif
!
      implicit none
!
!  Imported variable declarations.
!
      integer :: ng, ia
      integer :: IP, IX, IY

      real :: COMPDA(MCGRD,MCMVAR)
!
!  Local variable declarations.
!
      integer :: MyStatus, MyError, MySize, MyRank
      integer :: i, id, j, gsmsize, ierr, indx, Tag
      integer :: Istr, Iend, Jstr, Jend
      integer :: Isize, Jsize, INDXG, NPROCS, OFFSET
      integer :: NUMTRANSFER, NNEIGH, HALOSIZE, NUMSENT, INB
      integer :: WHICHWAY, GDEST, GSRC, TAGT, TAGB, TAGR, TAGL
      integer :: TREQUEST,BREQUEST,RREQUEST,LREQUEST,MSIZE
      integer :: iddep, idwlv, idvlx, idvly
      integer :: idu10, idv10
      integer, dimension(MPI_STATUS_SIZE,4) :: status
# ifdef MCT_INTERP_WV2AT
      integer, pointer :: indices(:)
# endif

      real :: cff
# if !defined BBL_MODEL
      real :: cff2
# endif
      real, pointer :: TEMPMCT(:,:)
      real, pointer :: GRECVT(:), GRECVB(:), GRECVR(:), GRECVL(:)
      real, pointer :: GSENDT(:), GSENDB(:), GSENDR(:), GSENDL(:)
      real(m8), pointer :: avdata(:)
      real, parameter ::  Large = 1.0E+20
      real, dimension(2) :: range

      character (len=40) :: code
!
!-----------------------------------------------------------------------
!  Send wave fields to ROMS.
!-----------------------------------------------------------------------
!
      CALL MPI_COMM_RANK (WAV_COMM_WORLD, MyRank, MyError)
      CALL MPI_COMM_SIZE (WAV_COMM_WORLD, NPROCS, MyError)
!
!  Get the number of grid point on this processor.
!
      gsmsize=GlobalSegMap_lsize(GlobalSegMap_G(ng)%GSMapSWAN,          &
     &                           WAV_COMM_WORLD)
!
!  Allocate attribute vector array used to export/import data.
!
      allocate ( avdata(gsmsize),stat=ierr )
      avdata=0.0_m8
      MyError=0
!
!-----------------------------------------------------------------------
!  Send wave fields bundle to ocean model, ROMS.
!-----------------------------------------------------------------------
!
 35   FORMAT (a14,i2,a24,i2)
!
!  Receive fields from atmosphere model.
!
      Tag=0*100+ia*10+ng
# ifdef MCT_INTERP_WV2AT
      CALL MCT_irecv (AV2_A(ng,ia)%atm2wav_AV2,                         &
     &                Router_A(ng,ia)%SWANtoWRF, Tag)
!     Wait to make sure the WRF data has arrived.
      CALL MCT_waitr (AV2_A(ng,ia)%atm2wav_AV2,                         &
     &                Router_A(ng,ia)%SWANtoWRF)
      CALL MCT_MatVecMul(AV2_A(ng,ia)%atm2wav_AV2,                      &
     &                   SMPlus_G(ng,ia)%A2WMatPlus,                    &
     &                   AttrVect_G(ng)%atm2wav_AV)
# else
      CALL MCT_irecv (AttrVect_G(ng)%atm2wav_AV,                        &
     &                Router_A(ng,ia)%SWANtoWRF, Tag)
!     Wait to make sure the WRF data has arrived.
      CALL MCT_waitr (AttrVect_G(ng)%atm2wav_AV,                        &
     &                Router_A(ng,ia)%SWANtoWRF)
# endif
        IF (MyRank.EQ.0) THEN
          WRITE (SCREEN,35)'== SWAN grid ',ng,' recv data from WRF grid'&
     &                     ,ia
        END IF
        IF (MyError.ne.0) THEN
          WRITE (*,*) 'coupling fail swancplr, MyStatus= ', MyError
          CALL FINALIZE_WAV_COUPLING(ng)
        END IF
!
! Pass the non-halo data from MCT into tempmct array.
!
        NNEIGH = IBLKAD(1)
        IF (NPROCS.eq.1) THEN
          Istr=1
          Iend=MXC
          Jstr=1
          Jend=MYC
        ELSE
          IF (MXCGL.GT.MYCGL) THEN
            IF (MyRank.eq.0) THEN
              Istr=1
            ELSE
              Istr=IHALOX+1
            END IF
            Isize=MXC-IHALOX*IBLKAD(1)
            Iend=Istr+Isize-1
            Jstr=1
            Jend=MYC
            HALOSIZE=IHALOX*MYC
          ELSE
            IF (MyRank.eq.0) THEN
              Jstr=1
            ELSE
              Jstr=IHALOY+1
            END IF
            Jsize=MYC-IHALOY*IBLKAD(1)
            Jend=Jstr+Jsize-1
            Istr=1
            Iend=MXC
            HALOSIZE=IHALOY*MXC
          END IF
        END IF
!
!  Determine the amount of fields and assign id numbers.
!
        IP=0
        IP=IP+1
        idu10=IP
!
        IP=IP+1
        idv10=IP
        NUMTRANSFER=IP
!
 40     FORMAT (a36,1x,2(1pe14.6))
        allocate ( TEMPMCT(MXC*MYC,NUMTRANSFER),stat=ierr )
        TEMPMCT=0.0
# ifdef AIR_WAVES
!
!  U-wind.
!
        CALL AttrVect_exportRAttr(AttrVect_G(ng)%atm2wav_AV,"U10",      &
     &                            avdata,gsmsize)
        range(1)= Large
        range(2)=-Large
        IP=0
        DO IY=Jstr,Jend
          DO IX=Istr,Iend
            IP=IP+1
            INDXG=(IY-1)*MXC+IX
            TEMPMCT(INDXG,idu10)=avdata(IP)
            range(1)=MIN(range(1),REAL(avdata(IP)))
            range(2)=MAX(range(2),REAL(avdata(IP)))
          END DO
        END DO
        CALL SWREDUCE(range(1),1,SWREAL,SWMIN)
        CALL SWREDUCE(range(2),1,SWREAL,SWMAX)
        IF (MyRank.eq.0) THEN
          write(SCREEN,40) 'WRFtoSWAN Min/Max  U10     (ms-1):  ',      &
     &                      range(1),range(2)
        END IF
!
!  V-wind.
!
        CALL AttrVect_exportRAttr(AttrVect_G(ng)%atm2wav_AV,"V10",      &
     &                            avdata,gsmsize)
        range(1)= Large
        range(2)=-Large
        IP=0
        DO IY=Jstr,Jend
          DO IX=Istr,Iend
            IP=IP+1
            INDXG=(IY-1)*MXC+IX
            TEMPMCT(INDXG,idv10)=avdata(IP)
            range(1)=MIN(range(1),REAL(avdata(IP)))
            range(2)=MAX(range(2),REAL(avdata(IP)))
          END DO
        END DO
        CALL SWREDUCE(range(1),1,SWREAL,SWMIN)
        CALL SWREDUCE(range(2),1,SWREAL,SWMAX)
        IF (MyRank.eq.0) THEN
          write(SCREEN,40) 'WRFtoSWAN Min/Max  V10     (ms-1):  ',      &
     &                      range(1),range(2)
        END IF
# endif
!
!  Pack and send halo regions to be exchanged with adjacent tiles.
!  IBLKAD contains the tile data.
!  WHICHWAY: [top, bot, right, left] = [1 2 3 4]
!
        IF (NPROCS.GT.1) THEN
          MSIZE=HALOSIZE*NUMTRANSFER
          IF (MXCGL.GT.MYCGL) THEN
            allocate ( GSENDR(MSIZE),stat=ierr )
            allocate ( GSENDL(MSIZE),stat=ierr )
            allocate ( GRECVR(MSIZE),stat=ierr )
            allocate ( GRECVL(MSIZE),stat=ierr )
            GSENDR=0.0
            GSENDL=0.0
            GRECVR=0.0
            GRECVL=0.0
          ELSE
            allocate ( GSENDT(MSIZE),stat=ierr )
            allocate ( GSENDB(MSIZE),stat=ierr )
            allocate ( GRECVT(MSIZE),stat=ierr )
            allocate ( GRECVB(MSIZE),stat=ierr )
            GSENDT=0.0
            GSENDB=0.0
            GRECVT=0.0
            GRECVB=0.0
          END IF
          TAGT=1
          TAGB=2
          TAGR=3
          TAGL=4
          DO INB=1,NNEIGH
            OFFSET=0
            WHICHWAY=IBLKAD(3*INB)
            DO NUMSENT=1,NUMTRANSFER
              IP=OFFSET
              IF (WHICHWAY.EQ.1) THEN
                DO IY=MYC-IHALOX-2,MYC-3
                  DO IX=1,MXC
                    IP=IP+1
                    INDXG=(IY-1)*MXC+IX
                    GSENDT(IP)=TEMPMCT(INDXG,NUMSENT)
                  END DO
                END DO
              ELSE IF (WHICHWAY.EQ.2) THEN
                DO IY=IHALOY+1,IHALOY+3
                  DO IX=1,MXC
                    IP=IP+1
                    INDXG=(IY-1)*MXC+IX
                    GSENDB(IP)=TEMPMCT(INDXG,NUMSENT)
                  END DO
                END DO
              ELSE IF (WHICHWAY.EQ.3) THEN
                DO IY=1,MYC
                  DO IX=MXC-IHALOX-2,MXC-3
                    IP=IP+1
                    INDXG=(IY-1)*MXC+IX
                    GSENDR(IP)=TEMPMCT(INDXG,NUMSENT)
                  END DO
                END DO
              ELSE IF (WHICHWAY.EQ.4) THEN
                DO IY=1,MYC
                  DO IX=IHALOX+1,IHALOX+3
                    IP=IP+1
                    INDXG=(IY-1)*MXC+IX
                    GSENDL(IP)=TEMPMCT(INDXG,NUMSENT)
                  END DO
                END DO
              END IF
              OFFSET=OFFSET+HALOSIZE
            END DO
          END DO
          DO INB=1,NNEIGH
            GSRC=IBLKAD(3*INB-1)-1
            WHICHWAY=IBLKAD(3*INB)
            IF (WHICHWAY.EQ.1) THEN
              CALL mpi_irecv (GRECVT,MSIZE,SWREAL,                      &
     &                        GSRC,TAGB,WAV_COMM_WORLD,TREQUEST,MyError)
            ELSE IF (WHICHWAY.EQ.2) THEN
              CALL mpi_irecv (GRECVB,MSIZE,SWREAL,                      &
     &                        GSRC,TAGT,WAV_COMM_WORLD,BREQUEST,MyError)
            ELSE IF (WHICHWAY.EQ.3) THEN
              CALL mpi_irecv (GRECVR,MSIZE,SWREAL,                      &
     &                        GSRC,TAGL,WAV_COMM_WORLD,RREQUEST,MyError)
            ELSE IF (WHICHWAY.EQ.4) THEN
              CALL mpi_irecv (GRECVL,MSIZE,SWREAL,                      &
     &                        GSRC,TAGR,WAV_COMM_WORLD,LREQUEST,MyError)
            END IF
          END DO
          DO INB=1,NNEIGH
            GDEST=IBLKAD(3*INB-1)-1
            WHICHWAY=IBLKAD(3*INB)
            IF (WHICHWAY.EQ.1) THEN
              CALL mpi_send (GSENDT,MSIZE,SWREAL,                       &
     &                       GDEST,TAGT,WAV_COMM_WORLD,MyError)
            ELSE IF (WHICHWAY.EQ.2) THEN
              CALL mpi_send (GSENDB,MSIZE,SWREAL,                       &
     &                       GDEST,TAGB,WAV_COMM_WORLD,MyError)
            ELSE IF (WHICHWAY.EQ.4) THEN
              CALL mpi_send (GSENDL,MSIZE,SWREAL,                       &
     &                       GDEST,TAGL,WAV_COMM_WORLD,MyError)
            ELSE IF (WHICHWAY.EQ.3) THEN
              CALL mpi_send (GSENDR,MSIZE,SWREAL,                       &
     &                       GDEST,TAGR,WAV_COMM_WORLD,MyError)
            END IF
          END DO
!
! Receive and unpack halo regions exchanged with adjacent tiles.
! [top, bot, right, left] = [1 2 3 4]
!
          DO INB=1,NNEIGH
            WHICHWAY=IBLKAD(3*INB)
            IF (WHICHWAY.EQ.1) THEN
              CALL mpi_wait (TREQUEST,status(1,1),MyError)
            ELSE IF (WHICHWAY.EQ.2) THEN
              CALL mpi_wait (BREQUEST,status(1,2),MyError)
            ELSE IF (WHICHWAY.EQ.3) THEN
              CALL mpi_wait (RREQUEST,status(1,3),MyError)
            ELSE IF (WHICHWAY.EQ.4) THEN
              CALL mpi_wait (LREQUEST,status(1,4),MyError)
            END IF
          END DO
!
          DO INB=1,NNEIGH
            OFFSET=0
            WHICHWAY=IBLKAD(3*INB)
            IF (WHICHWAY.EQ.1) THEN
              DO NUMSENT=1,NUMTRANSFER
                IP=OFFSET
                DO IY=MYC-2,MYC
                  DO IX=1,MXC
                    IP=IP+1
                    INDXG=(IY-1)*MXC+IX
                    TEMPMCT(INDXG,NUMSENT)=GRECVT(IP)
                  END DO
                END DO
                OFFSET=OFFSET+HALOSIZE
              END DO
            ELSE IF (WHICHWAY.EQ.2) THEN
              DO NUMSENT=1,NUMTRANSFER
                IP=OFFSET
                DO IY=1,IHALOY
                  DO IX=1,MXC
                    IP=IP+1
                    INDXG=(IY-1)*MXC+IX
                    TEMPMCT(INDXG,NUMSENT)=GRECVB(IP)
                  END DO
                END DO
                OFFSET=OFFSET+HALOSIZE
              END DO
            ELSE IF (WHICHWAY.EQ.3) THEN
              DO NUMSENT=1,NUMTRANSFER
                IP=OFFSET
                DO IY=1,MYC
                  DO IX=MXC-2,MXC
                    IP=IP+1
                    INDXG=(IY-1)*MXC+IX
                    TEMPMCT(INDXG,NUMSENT)=GRECVR(IP)
                  END DO
                END DO
                OFFSET=OFFSET+HALOSIZE
              END DO
            ELSE IF (WHICHWAY.EQ.4) THEN
              DO NUMSENT=1,NUMTRANSFER
                IP=OFFSET
                DO IY=1,MYC
                  DO IX=1,IHALOX
                    IP=IP+1
                    INDXG=(IY-1)*MXC+IX
                    TEMPMCT(INDXG,NUMSENT)=GRECVL(IP)
                  END DO
                END DO
                OFFSET=OFFSET+HALOSIZE
              END DO
            END IF
          END DO
          IF (MXCGL.GT.MYCGL) THEN
            deallocate (GRECVR,GRECVL,GSENDR,GSENDL)
          ELSE
            deallocate (GRECVT,GRECVB,GSENDT,GSENDB)
          END IF
        END IF
!
! Finally insert the full (MXC*MYC) TEMPMCT array into the SWAN
! array for DEPTH and computational array COMPDA. Only insert
! active (wet points) using array KGRPNT.
!
!
# ifdef AIR_WAVES
! Move values at 'present' time level 2 to 'old' time level 1.
! MCGRD = MXC*MYC+1-#masked cells.
! MXC = # cells x-dir in this tile including halox.
! MYC = # cells y-dir in this tile including haloy.
! COMPDA has only active wet points + 1.
!
        IF (ia.eq.1) THEN
          DO INDX = 1, MCGRD
            COMPDA(INDX,JWX3) =COMPDA(INDX,JWX2)
            COMPDA(INDX,JWY3) =COMPDA(INDX,JWY2)
          END DO
        END IF
# endif
!
! Insert bot level, water level, velx, vely, fric, and winds 
! into SWAN arrays.
!
! If not using a BBL model, then determine the non-spatially 
! varying friction coeff to enter into the JFRC2 array.
!
        IP=0
        DO IY=1,MYC
          DO IX=1,MXC
            IP=IP+1
            INDX = KGRPNT(IX,IY)
            IF (INDX.GT.1) THEN
              IF (ia.eq.1) THEN
                COMPDA(INDX,JWX2)=TEMPMCT(IP,idu10)
                COMPDA(INDX,JWY2)=TEMPMCT(IP,idv10)
              ELSE
                COMPDA(INDX,JWX2)=COMPDA(INDX,JWX2)+TEMPMCT(IP,idu10)
                COMPDA(INDX,JWY2)=COMPDA(INDX,JWY2)+TEMPMCT(IP,idv10)
              END IF
            END IF
          END DO
        END DO
!
      deallocate (TEMPMCT)
      deallocate (avdata)
# ifdef MCT_INTERP_WV2AT
      if (associated (indices)) then
        deallocate (indices)
      endif
# endif
!
      RETURN
      END SUBROUTINE WAVFATM_COUPLING
# endif

      SUBROUTINE FINALIZE_WAV_COUPLING(ng)
!
!=======================================================================
!                                                                    ===
!  This routines terminates execution during coupling error.         ===
!                                                                    ===
!=======================================================================
      USE mct_coupler_params
!
!  Local variable declarations.
!
      integer :: ng, io, ia, MyError
!
!-----------------------------------------------------------------------
!  Deallocate MCT environment.
!-----------------------------------------------------------------------
!
# ifdef ROMS_COUPLING
      DO io=1,Nocn_grids
        CALL Router_clean (Router_O(ng,io)%SWANtoROMS, MyError)
      END DO
      CALL AttrVect_clean (AttrVect_G(ng)%wav2ocn_AV, MyError)
# endif
# ifdef WRF_COUPLING
      DO ia=1,Natm_grids
        CALL Router_clean (Router_A(ng,ia)%SWANtoWRF, MyError)
      END DO
      CALL AttrVect_clean (AttrVect_G(ng)%atm2wav_AV, MyError)
# endif
      CALL GlobalSegMap_clean (GlobalSegMap_G(ng)%GSMapSWAN, MyError)

      END SUBROUTINE FINALIZE_WAV_COUPLING
#endif
      END MODULE WAVES_COUPLER_MOD