!
! Public domain.
!

!>
!! Module to retrieve a variable from a WRF I/O Internal format file
!! using MPI IO.
!!
!! This module uses an index generated by module_io_int_idx.
!!
!! If the MPI implementations actually supported the external32
!! view on a file we could simply call:
!!
!! MPI_File_set_view(..., "external32", ...)
!!
!! However most do not. Instead we have do a native view and then
!! do a byte order rearrangement.
!!
!! Usage consists of:
!! 1) Creating an index of records in the file
!! 2) Opening the file with MPI
!! 3) Retrieving the data
!! 4) Closing the file with MPI
!!
!! Example:
!!
!! USE module_io_int_idx, only: io_int_index, r_info
!! USE mpiio_rw, only: read_data, write_data
!!
!! type(r_info), pointer :: r(:) => NULL()  ! Define the array of records
!!
!! call io_int_index(filename, r, ierr)     ! Get an index of reocrds
!!
!! call mpi_file_open(mpi_comm_world, trim(filename), &
!!                    mpi_mode_rdonly, mpi_info_null, &
!!                    iunit, ierr)          ! MPI Open the file
!!
!! call read_data(iunit, r, 'TRUELAT1', garb, ierr)
!!                                          ! Get an record (into garb)
!!
!! call mpi_file_close(iunit, ierr)         ! Close the file
!!

module mpiio_rw

    use module_io_int_idx,           only: io_int_loc, r_info
    use, intrinsic :: iso_c_binding, only: c_int32_t

    implicit none

    private
    public :: read_data, write_data

    interface
        integer(c_int32_t)           &
        pure function ntohl(netlong) &
            bind(c, name="ntohl")
            import :: c_int32_t
            integer(c_int32_t), value, intent(in) :: netlong
        end function
        integer(c_int32_t)            &
        pure function htonl(hostlong) &
            bind(c, name="htonl")
            import :: c_int32_t
            integer(c_int32_t), value, intent(in) :: hostlong
        end function
    end interface
    interface read_data
        module procedure  read_i0, read_i1, read_i2, read_i3, &
                          read_r0, read_r1, read_r2, read_r3, &
                          read_c1
    end interface read_data
    interface write_data
        module procedure  write_i0, write_i1, write_i2, write_i3, &
                          write_r0, write_r1, write_r2, write_r3, &
                          write_c1
    end interface write_data

    include "mpif.h"

contains

    !>
    !! Read a single integer
    !
    subroutine read_i0(ifd, records, varname, pos, n, dst, ierr)

    integer,                                 intent(in)  :: ifd
    type(r_info), pointer,                   intent(in)  :: records(:)
    character(len=*),                        intent(in)  :: varname
    integer(kind=mpi_offset_kind), optional, intent(in)  :: pos
    integer, optional,                       intent(in)  :: n
    integer,                                 intent(out) :: dst
    integer,                                 intent(out) :: ierr

    integer(kind=mpi_offset_kind)       :: offset
    integer                             :: count
    integer                             :: tmp

    call io_int_loc(varname, records, offset, count, ierr)
    if (ierr .ne. 0) then
        return
    end if

    if (present(pos)) then
        offset = pos
    end if

    call mpi_file_read_at(ifd, offset, tmp, 1, &
                          mpi_integer4, mpi_status_ignore, ierr)
    if (ierr .ne. 0) then
        write(0,*) 'MPI IO: Unable to read ', varname
        return
    end if
    dst = ntohl(tmp)

    end subroutine read_i0

    !>
    !! Read a 1D integer array.
    !
    subroutine read_i1(ifd, records, varname, pos, n, dst, ierr)

    integer,                                 intent(in)    :: ifd
    type(r_info), pointer,                   intent(in)    :: records(:)
    character(len=*),                        intent(in)    :: varname
    integer(kind=mpi_offset_kind), optional, intent(in)    :: pos
    integer, optional,                       intent(in)    :: n
    integer,                                 intent(inout) :: dst(:)
    integer,                                 intent(out)   :: ierr

    integer(kind=mpi_offset_kind)          :: offset
    integer                                :: count
    integer                                :: num
    integer                                :: i
    integer                                :: its, ite
    integer, allocatable, dimension(:)     :: tmp

    call io_int_loc(varname, records, offset, count, ierr)
    if (ierr .ne. 0) then
        return
    end if

    if (present(pos)) then
        offset = pos
    end if
    if (present(n)) then
        num = n
        its = 1
        ite = num
    else
        its = lbound(dst,1)
        ite = ubound(dst,1)
        num = ite - its + 1
    end if

    allocate(tmp(its:ite), stat=ierr)
    if (ierr .ne. 0) then
        write(0,*) 'Unable to allocate a temporary array'
        return
    end if
    tmp = 0

    call mpi_file_read_at(ifd, offset, tmp, num, &
                          mpi_integer4, mpi_status_ignore, ierr)
    if (ierr .ne. 0) then
        write(0,*) 'MPI IO: Unable to read ', varname
        return
    end if

    forall(i=its:ite)
        dst(i) = transfer(ntohl(tmp(i)), 1.0)
    end forall

    deallocate(tmp)

    end subroutine read_i1

    !>
    !! Read a 2D integer array.
    !
    subroutine read_i2(ifd, records, varname, pos, n, dst, ierr)

    integer,                                 intent(in)    :: ifd
    type(r_info), pointer,                   intent(in)    :: records(:)
    character(len=*),                        intent(in)    :: varname
    integer(kind=mpi_offset_kind), optional, intent(in)    :: pos
    integer, optional,                       intent(in)    :: n
    integer,                                 intent(inout) :: dst(:,:)
    integer,                                 intent(out)   :: ierr

    integer(kind=mpi_offset_kind)          :: offset
    integer                                :: count
    integer                                :: num
    integer                                :: i, j
    integer                                :: its, ite, jts, jte
    integer, allocatable, dimension(:,:)   :: tmp

    call io_int_loc(varname, records, offset, count, ierr)
    if (ierr .ne. 0) then
        return
    end if

    if (present(pos)) then
        offset = pos
    end if

    its = lbound(dst,1)
    ite = ubound(dst,1)
    jts = lbound(dst,2)
    jte = ubound(dst,2)
    num = (ite - its + 1) * (jte - jts +1)

    allocate(tmp(its:ite, jts:jte), stat=ierr)
    if (ierr .ne. 0) then
        write(0,*) 'Unable to allocate a temporary array'
        return
    end if
    tmp = 0

    if (present(n)) then
        num = n
    end if

    call mpi_file_read_at(ifd, offset, tmp, num, &
                          mpi_integer4, mpi_status_ignore, ierr)
    if (ierr .ne. 0) then
        write(0,*) 'MPI IO: Unable to read ', varname
        return
    end if

    forall(i=its:ite, j=jts:jte)
        dst(i,j) = transfer(ntohl(tmp(i,j)), 1.0)
    end forall

    deallocate(tmp)

    end subroutine read_i2

    !>
    !! Read a 3D integer array
    !
    subroutine read_i3(ifd, records, varname, pos, n, dst, ierr)

    integer,                                 intent(in)    :: ifd
    type(r_info), pointer,                   intent(in)    :: records(:)
    character(len=*),                        intent(in)    :: varname
    integer(kind=mpi_offset_kind), optional, intent(in)    :: pos
    integer, optional,                       intent(in)    :: n
    integer,                                 intent(inout) :: dst(:,:,:)
    integer,                                 intent(out)   :: ierr

    integer(kind=mpi_offset_kind)          :: offset
    integer                                :: count
    integer                                :: num
    integer                                :: i, j, k
    integer                                :: its, ite, jts, jte, kts, kte
    integer, allocatable, dimension(:,:,:) :: tmp

    call io_int_loc(varname, records, offset, count, ierr)
    if (ierr .ne. 0) then
        return
    end if

    if (present(pos)) then
        offset = pos
    end if

    its = lbound(dst,1)
    ite = ubound(dst,1)
    jts = lbound(dst,2)
    jte = ubound(dst,2)
    kts = lbound(dst,3)
    kte = ubound(dst,3)
    num = (ite - its + 1) * (jte - jts + 1) * (kte - kts + 1)

    allocate(tmp(its:ite, jts:jte, kts:kte), stat=ierr)
    if (ierr .ne. 0) then
        write(0,*) 'Unable to allocate a temporary array'
        return
    end if
    tmp = 0

    if (present(n)) then
        num = n
    end if

    call mpi_file_read_at(ifd, offset, tmp, num, &
                          mpi_integer4, mpi_status_ignore, ierr)

    if (ierr .ne. 0) then
        write(0,*) 'MPI IO: Unable to read ', varname
        return
    end if

    forall(i=its:ite, j=jts:jte, k=kts:kte)
        dst(i,j,k) = transfer(ntohl(tmp(i,j,k)), 1.0)
    end forall

    deallocate(tmp)

    end subroutine read_i3

    !>
    !! Read a single real.
    !
    subroutine read_r0(ifd, records, varname, pos, n, dst, ierr)

    integer,                                 intent(in)  :: ifd
    type(r_info), pointer,                   intent(in)  :: records(:)
    character(len=*),                        intent(in)  :: varname
    integer(kind=mpi_offset_kind), optional, intent(in)  :: pos
    integer, optional,                       intent(in)  :: n
    real,                                    intent(out) :: dst
    integer,                                 intent(out) :: ierr

    integer(kind=mpi_offset_kind)        :: offset
    integer                              :: count
    integer                              :: tmp

    call io_int_loc(varname, records, offset, count, ierr)
    if (ierr .ne. 0) then
        return
    end if
    if (present(pos)) then
        offset = pos
    end if


    call mpi_file_read_at(ifd, offset, tmp, 1, &
                          mpi_integer4,  mpi_status_ignore, ierr)
    if (ierr .ne. 0) then
        write(0,*) 'MPI IO: Unable to read ', varname
        return
    end if
    dst = transfer(ntohl(tmp), 1.0)

    end subroutine read_r0

    !>
    !! Read a 1D real array.
    !
    subroutine read_r1(ifd, records, varname, pos, n, dst, ierr)

    integer,                                 intent(in)    :: ifd
    type(r_info), pointer,                   intent(in)    :: records(:)
    character(len=*),                        intent(in)    :: varname
    integer(kind=mpi_offset_kind), optional, intent(in)    :: pos
    integer, optional,                       intent(in)    :: n
    real,                                    intent(inout) :: dst(:)
    integer,                                 intent(out)   :: ierr

    integer(kind=mpi_offset_kind)          :: offset
    integer                                :: count
    integer                                :: num
    integer                                :: i
    integer                                :: its, ite
    integer, allocatable, dimension(:)     :: tmp

    call io_int_loc(varname, records, offset, count, ierr)
    if (ierr .ne. 0) then
        return
    end if

    if (present(pos)) then
        offset = pos
    end if
    if (present(n)) then
        num = n
        its = 1
        ite = num
    else
        its = lbound(dst,1)
        ite = ubound(dst,1)
        num = ite - its + 1
    end if

    allocate(tmp(its:ite), stat=ierr)
    if (ierr .ne. 0) then
        write(0,*) 'Unable to allocate a temporary array'
        return
    end if
    tmp = 0

    call mpi_file_read_at(ifd, offset, tmp, num, &
                          mpi_integer4, mpi_status_ignore, ierr)
    if (ierr .ne. 0) then
        write(0,*) 'MPI IO: Unable to read ', varname
        return
    end if

    forall(i=its:ite)
        dst(i) = transfer(ntohl(tmp(i)), 1.0)
    end forall

    deallocate(tmp)

    end subroutine read_r1

    !>
    !! Read a 2D real array.
    !
    subroutine read_r2(ifd, records, varname, pos, n, dst, ierr)

    integer,                                 intent(in)    :: ifd
    type(r_info), pointer,                   intent(in)    :: records(:)
    character(len=*),                        intent(in)    :: varname
    integer(kind=mpi_offset_kind), optional, intent(in)    :: pos
    integer, optional,                       intent(in)    :: n
    real,                                    intent(inout) :: dst(:,:)
    integer,                                 intent(out)   :: ierr

    integer(kind=mpi_offset_kind)          :: offset
    integer                                :: count
    integer                                :: num
    integer                                :: i, j
    integer                                :: its, ite, jts, jte
    integer, allocatable, dimension(:,:)   :: tmp

    call io_int_loc(varname, records, offset, count, ierr)
    if (ierr .ne. 0) then
        return
    end if

    if (present(pos)) then
        offset = pos
    end if

    its = lbound(dst,1)
    ite = ubound(dst,1)
    jts = lbound(dst,2)
    jte = ubound(dst,2)
    num = (ite - its + 1) * (jte - jts + 1)

    allocate(tmp(its:ite, jts:jte), stat=ierr)
    if (ierr .ne. 0) then
        write(0,*) 'Unable to allocate a temporary array'
        return
    end if
    tmp = 0

    if (present(n)) then
        if ( n .le. num) then
            num = n
        end if
    end if
    call mpi_file_read_at(ifd, offset, tmp, num, &
                          mpi_integer4, mpi_status_ignore, ierr)
    if (ierr .ne. 0) then
        write(0,*) 'MPI IO: Unable to read ', varname
        return
    end if

    forall(i=its:ite, j=jts:jte)
        dst(i,j) = transfer(ntohl(tmp(i,j)), 1.0)
    end forall

    deallocate(tmp)

    end subroutine read_r2

    !>
    !! Read a 3D real array
    !
    subroutine read_r3(ifd, records, varname, pos, n, dst, ierr)

    integer,                                 intent(in)    :: ifd
    type(r_info), pointer,                   intent(in)    :: records(:)
    character(len=*),                        intent(in)    :: varname
    integer(kind=mpi_offset_kind), optional, intent(in)    :: pos
    integer, optional,                       intent(in)    :: n
    real,                                    intent(inout) :: dst(:,:,:)
    integer,                                 intent(out)   :: ierr

    integer(kind=mpi_offset_kind)          :: offset
    integer                                :: count
    integer                                :: num
    integer                                :: i, j, k
    integer                                :: its, ite, jts, jte, kts, kte
    integer, allocatable, dimension(:,:,:) :: tmp

    call io_int_loc(varname, records, offset, count, ierr)
    if (ierr .ne. 0) then
        return
    end if

    if (present(pos)) then
        offset = pos
    end if

    its = lbound(dst,1)
    ite = ubound(dst,1)
    jts = lbound(dst,2)
    jte = ubound(dst,2)
    kts = lbound(dst,3)
    kte = ubound(dst,3)
    num = (ite - its + 1) * (jte - jts + 1) * (kte - kts + 1)

    allocate(tmp(its:ite, jts:jte, kts:kte), stat=ierr)
    if (ierr .ne. 0) then
        write(0,*) 'Unable to allocate a temporary array'
        return
    end if
    tmp = 0

    if (present(n)) then
        if ( n .le. num) then
            num = n
        end if
    end if
    call mpi_file_read_at(ifd, offset, tmp, num, &
                          mpi_integer4, mpi_status_ignore, ierr)

    if (ierr .ne. 0) then
        write(0,*) 'MPI IO: Unable to read ', varname
        return
    end if

    forall(i=its:ite, j=jts:jte, k=kts:kte)
        dst(i,j,k) = transfer(ntohl(tmp(i,j,k)), 1.0)
    end forall

    deallocate(tmp)

    end subroutine read_r3

    !>
    !! Read a 1D character array.
    !
    subroutine read_c1(ifd, records, varname, pos, n, dst, ierr)

    integer,                                 intent(in)    :: ifd
    type(r_info), pointer,                   intent(in)    :: records(:)
    character(len=*),                        intent(in)    :: varname
    integer(kind=mpi_offset_kind), optional, intent(in)    :: pos
    integer, optional,                       intent(in)    :: n
    character(len=*),                        intent(inout) :: dst
    integer,                                 intent(out)   :: ierr

    integer(kind=mpi_offset_kind)          :: offset
    integer                                :: count
    integer                                :: num
    integer                                :: i
    integer, allocatable, dimension(:)     :: tmp


    call io_int_loc(varname, records, offset, count, ierr)
    if (ierr .ne. 0) then
        return
    end if

    if (present(pos)) then
        offset = pos
    end if

    num = count

    allocate(tmp(num), stat=ierr)
    if (ierr .ne. 0) then
        write(0,*) 'Unable to allocate a temporary array'
        return
    end if
    tmp = 0

    if (present(n)) then
        if ( n .le. num) then
            num = n
        end if
    end if
    call mpi_file_read_at(ifd, offset, tmp, num, &
                          mpi_integer4, mpi_status_ignore, ierr)
    if (ierr .ne. 0) then
        write(0,*) 'MPI IO: Unable to read ', varname
        return
    end if

! PGI and GNU fortran compilers do not like forall loops over
! character arrays.
! PGI: TPR#19016
! GNU: http://gcc.gnu.org/bugzilla/show_bug.cgi?id=50069
!
! XLF does not like transfer() with the mold being an array,
! so we use achar() instead.
!
!    forall(i=1:num)
!        dst(i:i) = transfer(ntohl(tmp(i)), dst)
!    end forall

    do i=1,num
        dst(i:i) = achar(ntohl(tmp(i)))
    end do

    deallocate(tmp)

    end subroutine read_c1

    !>
    !! write a single integer
    !
    subroutine write_i0(ifd, records, varname, pos, n, src, ierr)

    integer,                                 intent(in)  :: ifd
    type(r_info), pointer,                   intent(in)  :: records(:)
    character(len=*),                        intent(in)  :: varname
    integer(kind=mpi_offset_kind), optional, intent(in)  :: pos
    integer, optional,                       intent(in)  :: n
    integer,                                 intent(in)  :: src
    integer,                                 intent(out) :: ierr

    integer(kind=mpi_offset_kind)       :: offset
    integer                             :: count
    integer                             :: tmp

    call io_int_loc(varname, records, offset, count, ierr)
    if (ierr .ne. 0) then
        return
    end if

    if (present(pos)) then
        offset = pos
    end if

    tmp = htonl(src)

    call mpi_file_write_at(ifd, offset, tmp, 1, &
                          mpi_integer4, mpi_status_ignore, ierr)
    if (ierr .ne. 0) then
        write(0,*) 'MPI IO: Unable to write ', varname
        return
    end if

    end subroutine write_i0

    !>
    !! write a 1D integer array.
    !
    subroutine write_i1(ifd, records, varname, pos, n, src, ierr)

    integer,                                 intent(in)    :: ifd
    type(r_info), pointer,                   intent(in)    :: records(:)
    character(len=*),                        intent(in)    :: varname
    integer(kind=mpi_offset_kind), optional, intent(in)    :: pos
    integer, optional,                       intent(in)    :: n
    integer,                                 intent(in)    :: src(:)
    integer,                                 intent(out)   :: ierr

    integer(kind=mpi_offset_kind)          :: offset
    integer                                :: count
    integer                                :: num
    integer                                :: i
    integer                                :: its, ite
    integer, allocatable, dimension(:)     :: tmp

    call io_int_loc(varname, records, offset, count, ierr)
    if (ierr .ne. 0) then
        return
    end if

    if (present(pos)) then
        offset = pos
    end if
    if (present(n)) then
        num = n
        its = 1
        ite = num
    else
        its = lbound(src,1)
        ite = ubound(src,1)
        num = ite - its + 1
    end if

    allocate(tmp(its:ite), stat=ierr)
    if (ierr .ne. 0) then
        write(0,*) 'Unable to allocate a temporary array'
        return
    end if
    forall(i=its:ite)
        tmp(i) = htonl(src(i))
    end forall

    call mpi_file_write_at(ifd, offset, tmp, num, &
                          mpi_integer4, mpi_status_ignore, ierr)
    if (ierr .ne. 0) then
        write(0,*) 'MPI IO: Unable to write ', varname
        return
    end if

    deallocate(tmp)

    end subroutine write_i1

    !>
    !! write a 2D integer array.
    !
    subroutine write_i2(ifd, records, varname, pos, n, src, ierr)

    integer,                                 intent(in)    :: ifd
    type(r_info), pointer,                   intent(in)    :: records(:)
    character(len=*),                        intent(in)    :: varname
    integer(kind=mpi_offset_kind), optional, intent(in)    :: pos
    integer, optional,                       intent(in)    :: n
    integer,                                 intent(in)    :: src(:,:)
    integer,                                 intent(out)   :: ierr

    integer(kind=mpi_offset_kind)          :: offset
    integer                                :: count
    integer                                :: num
    integer                                :: i, j
    integer                                :: its, ite, jts, jte
    integer, allocatable, dimension(:,:)   :: tmp

    call io_int_loc(varname, records, offset, count, ierr)
    if (ierr .ne. 0) then
        return
    end if

    if (present(pos)) then
        offset = pos
    end if

    its = lbound(src,1)
    ite = ubound(src,1)
    jts = lbound(src,2)
    jte = ubound(src,2)
    num = (ite - its + 1) * (jte - jts +1)

    allocate(tmp(its:ite, jts:jte), stat=ierr)
    if (ierr .ne. 0) then
        write(0,*) 'Unable to allocate a temporary array'
        return
    end if
    forall(i=its:ite, j=jts:jte)
        tmp(i,j) = htonl(src(i,j))
    end forall

    if (present(n)) then
        num = n
    end if

    call mpi_file_write_at(ifd, offset, tmp, num, &
                          mpi_integer4, mpi_status_ignore, ierr)
    if (ierr .ne. 0) then
        write(0,*) 'MPI IO: Unable to write ', varname
        return
    end if

    deallocate(tmp)

    end subroutine write_i2

    !>
    !! write a 3D integer array
    !
    subroutine write_i3(ifd, records, varname, pos, n, src, ierr)

    integer,                                 intent(in)    :: ifd
    type(r_info), pointer,                   intent(in)    :: records(:)
    character(len=*),                        intent(in)    :: varname
    integer(kind=mpi_offset_kind), optional, intent(in)    :: pos
    integer, optional,                       intent(in)    :: n
    integer,                                 intent(in)    :: src(:,:,:)
    integer,                                 intent(out)   :: ierr

    integer(kind=mpi_offset_kind)          :: offset
    integer                                :: count
    integer                                :: num
    integer                                :: i, j, k
    integer                                :: its, ite, jts, jte, kts, kte
    integer, allocatable, dimension(:,:,:) :: tmp

    call io_int_loc(varname, records, offset, count, ierr)
    if (ierr .ne. 0) then
        return
    end if

    if (present(pos)) then
        offset = pos
    end if

    its = lbound(src,1)
    ite = ubound(src,1)
    jts = lbound(src,2)
    jte = ubound(src,2)
    kts = lbound(src,3)
    kte = ubound(src,3)
    num = (ite - its + 1) * (jte - jts + 1) * (kte - kts + 1)

    allocate(tmp(its:ite, jts:jte, kts:kte), stat=ierr)
    if (ierr .ne. 0) then
        write(0,*) 'Unable to allocate a temporary array'
        return
    end if

    forall(i=its:ite, j=jts:jte, k=kts:kte)
        tmp(i,j,k) = htonl(src(i,j,k))
    end forall

    if (present(n)) then
        num = n
    end if

    call mpi_file_write_at(ifd, offset, tmp, num, &
                          mpi_integer4, mpi_status_ignore, ierr)

    if (ierr .ne. 0) then
        write(0,*) 'MPI IO: Unable to write ', varname
        return
    end if

    deallocate(tmp)

    end subroutine write_i3

    !>
    !! write a single real.
    !
    subroutine write_r0(ifd, records, varname, pos, n, src, ierr)

    integer,                                 intent(in)  :: ifd
    type(r_info), pointer,                   intent(in)  :: records(:)
    character(len=*),                        intent(in)  :: varname
    integer(kind=mpi_offset_kind), optional, intent(in)  :: pos
    integer, optional,                       intent(in)  :: n
    real,                                    intent(in)  :: src
    integer,                                 intent(out) :: ierr

    integer(kind=mpi_offset_kind)        :: offset
    integer                              :: count
    integer                              :: tmp

    call io_int_loc(varname, records, offset, count, ierr)
    if (ierr .ne. 0) then
        return
    end if
    if (present(pos)) then
        offset = pos
    end if

    tmp = htonl(transfer(src, 1))

    call mpi_file_write_at(ifd, offset, tmp, 1, &
                          mpi_integer4,  mpi_status_ignore, ierr)
    if (ierr .ne. 0) then
        write(0,*) 'MPI IO: Unable to write ', varname
        return
    end if

    end subroutine write_r0

    !>
    !! write a 1D real array.
    !
    subroutine write_r1(ifd, records, varname, pos, n, src, ierr)

    integer,                                 intent(in)    :: ifd
    type(r_info), pointer,                   intent(in)    :: records(:)
    character(len=*),                        intent(in)    :: varname
    integer(kind=mpi_offset_kind), optional, intent(in)    :: pos
    integer, optional,                       intent(in)    :: n
    real,                                    intent(inout) :: src(:)
    integer,                                 intent(out)   :: ierr

    integer(kind=mpi_offset_kind)          :: offset
    integer                                :: count
    integer                                :: num
    integer                                :: i
    integer                                :: its, ite
    integer, allocatable, dimension(:)     :: tmp

    call io_int_loc(varname, records, offset, count, ierr)
    if (ierr .ne. 0) then
        return
    end if

    if (present(pos)) then
        offset = pos
    end if
    if (present(n)) then
        num = n
        its = 1
        ite = num
    else
        its = lbound(src,1)
        ite = ubound(src,1)
        num = ite - its + 1
    end if

    allocate(tmp(its:ite), stat=ierr)
    if (ierr .ne. 0) then
        write(0,*) 'Unable to allocate a temporary array'
        return
    end if

    forall(i=its:ite)
        tmp(i) = htonl(transfer(src(i), 1))
    end forall

    call mpi_file_write_at(ifd, offset, tmp, num, &
                          mpi_integer4, mpi_status_ignore, ierr)
    if (ierr .ne. 0) then
        write(0,*) 'MPI IO: Unable to write ', varname
        return
    end if

    deallocate(tmp)

    end subroutine write_r1

    !>
    !! write a 2D real array.
    !
    subroutine write_r2(ifd, records, varname, pos, n, src, ierr)

    integer,                                 intent(in)    :: ifd
    type(r_info), pointer,                   intent(in)    :: records(:)
    character(len=*),                        intent(in)    :: varname
    integer(kind=mpi_offset_kind), optional, intent(in)    :: pos
    integer, optional,                       intent(in)    :: n
    real,                                    intent(inout) :: src(:,:)
    integer,                                 intent(out)   :: ierr

    integer(kind=mpi_offset_kind)          :: offset
    integer                                :: count
    integer                                :: num
    integer                                :: i, j
    integer                                :: its, ite, jts, jte
    integer, allocatable, dimension(:,:)   :: tmp

    call io_int_loc(varname, records, offset, count, ierr)
    if (ierr .ne. 0) then
        return
    end if

    if (present(pos)) then
        offset = pos
    end if

    its = lbound(src,1)
    ite = ubound(src,1)
    jts = lbound(src,2)
    jte = ubound(src,2)
    num = (ite - its + 1) * (jte - jts + 1)

    allocate(tmp(its:ite, jts:jte), stat=ierr)
    if (ierr .ne. 0) then
        write(0,*) 'Unable to allocate a temporary array'
        return
    end if

    forall(i=its:ite, j=jts:jte)
        tmp(i,j) = htonl(transfer(src(i,j), 1))
    end forall

    if (present(n)) then
        if ( n .le. num) then
            num = n
        end if
    end if
    call mpi_file_write_at(ifd, offset, tmp, num, &
                          mpi_integer4, mpi_status_ignore, ierr)
    if (ierr .ne. 0) then
        write(0,*) 'MPI IO: Unable to write ', varname
        return
    end if

    deallocate(tmp)

    end subroutine write_r2

    !>
    !! write a 3D real array
    !
    subroutine write_r3(ifd, records, varname, pos, n, src, ierr)

    integer,                                 intent(in)    :: ifd
    type(r_info), pointer,                   intent(in)    :: records(:)
    character(len=*),                        intent(in)    :: varname
    integer(kind=mpi_offset_kind), optional, intent(in)    :: pos
    integer, optional,                       intent(in)    :: n
    real,                                    intent(in)    :: src(:,:,:)
    integer,                                 intent(out)   :: ierr

    integer(kind=mpi_offset_kind)          :: offset
    integer                                :: count
    integer                                :: num
    integer                                :: i, j, k
    integer                                :: its, ite, jts, jte, kts, kte
    integer, allocatable, dimension(:,:,:) :: tmp

    call io_int_loc(varname, records, offset, count, ierr)
    if (ierr .ne. 0) then
        return
    end if

    if (present(pos)) then
        offset = pos
    end if

    its = lbound(src,1)
    ite = ubound(src,1)
    jts = lbound(src,2)
    jte = ubound(src,2)
    kts = lbound(src,3)
    kte = ubound(src,3)
    num = (ite - its + 1) * (jte - jts + 1) * (kte - kts + 1)

    allocate(tmp(its:ite, jts:jte, kts:kte), stat=ierr)
    if (ierr .ne. 0) then
        write(0,*) 'Unable to allocate a temporary array'
        return
    end if

    forall(i=its:ite, j=jts:jte, k=kts:kte)
        tmp(i,j,k) = htonl(transfer(src(i,j,k), 1))
    end forall

    if (present(n)) then
        if ( n .le. num) then
            num = n
        end if
    end if
    call mpi_file_write_at(ifd, offset, tmp, num, &
                          mpi_integer4, mpi_status_ignore, ierr)

    if (ierr .ne. 0) then
        write(0,*) 'MPI IO: Unable to write ', varname
        return
    end if

    deallocate(tmp)

    end subroutine write_r3

    !>
    !! write a 1D character array.
    !
    subroutine write_c1(ifd, records, varname, pos, n, src, ierr)

    integer,                                 intent(in)    :: ifd
    type(r_info), pointer,                   intent(in)    :: records(:)
    character(len=*),                        intent(in)    :: varname
    integer(kind=mpi_offset_kind), optional, intent(in)    :: pos
    integer, optional,                       intent(in)    :: n
    character(len=*),                        intent(in)    :: src
    integer,                                 intent(out)   :: ierr

    integer(kind=mpi_offset_kind)          :: offset
    integer                                :: count
    integer                                :: num
    integer                                :: i
    integer, allocatable, dimension(:)     :: tmp


    call io_int_loc(varname, records, offset, count, ierr)
    if (ierr .ne. 0) then
        return
    end if

    if (present(pos)) then
        offset = pos
    end if

    num = count

    if (present(n)) then
        if ( n .le. num) then
            num = n
        end if
    end if

    allocate(tmp(num), stat=ierr)
    if (ierr .ne. 0) then
        write(0,*) 'Unable to allocate a temporary array'
        return
    end if

    do i=1,num
        tmp(i) = htonl(ichar(src(i:i)))
    end do
    call mpi_file_write_at(ifd, offset, tmp, num, &
                          mpi_integer4, mpi_status_ignore, ierr)
    if (ierr .ne. 0) then
        write(0,*) 'MPI IO: Unable to write ', varname
        return
    end if

    deallocate(tmp)

    end subroutine write_c1

end module mpiio_rw