c***********************************************************************
c
c transpose subsets of 2D arrays
c
c***********************************************************************

      program f_transp

      implicit none

#include <fdds.h>

c For DDS functions. Arguments and return values
      character prog*(8), Title*(80)
      integer   ier, fd_in, fd_out

      integer*8 ier8,izero8
      parameter (izero8=0)

c Data attributes
      character axes(RANK_MAX)*(AXISNAME_MAX)
      integer n1,n2,n3, rank

c User parameters
      integer i1s,i1e, i2s,i2e

c Internal dynamic array
      real    record_in(1), record_out(1)
      pointer (ptr_record_in,record_in)
      pointer (ptr_record_out,record_out)

c Internal loop counters and other scalar values
      integer i3, irank, n1_out,n2_out
      logical lverb

c-------------------------------------------------------------------
c
c Initialize the Title
c
      prog='f_transp'
      Title=prog//': transpose subsets of 2D arrays'

c
c Open the print file and dump command-line help is requested
c
      ier=fdds_openpr(prog,'')
      if (ier.gt.0) call help(prog,Title)

c
c Open input file
c
      fd_in=fddx_in('in','stdin:',Title)

c
c Retrieve data characteristics
c
      ier=fdds_scanf('size.axis(1)','%d\0',n1)
      ier=fdds_scanf('size.axis(2)','%d\0',n2)
      n3=fdds_axis_prod(3)

      rank=fdds_scank('axis','')
      do irank=1,rank
         ier=fdds_scanf('','%s\0',axes(irank))
      enddo

c
c Retrieve command line and parfile parameters
c
      ier=fdds_dict('par:','scan')

      i1s=1
      i1e=n1
      i2s=1
      i2e=n2
      ier=fdds_scanf('1s','%d\0',i1s)
      ier=fdds_scanf('1e','%d\0',i1e)
      ier=fdds_scanf('2s','%d\0',i2s)
      ier=fdds_scanf('2e','%d\0',i2e)

      lverb=.false.
      if (fdds_scank('verbose',DDS_TRUE_KEY).gt.0) lverb=.true.

c
c Set up the output dictionary
c
      fd_out=fddx_out('out','stdout:',Title,fd_in)

c
c Print modified axis info to output dictionary
c
      ier=fdds_printf('axis',' %s\0',axes(2))
      ier=fdds_printf('',' %s\0',axes(1))
      do irank=3,rank
         ier=fdds_printf('',' %s\0',axes(irank))
      enddo
      ier=fdds_printf('','\n\0')

c
c Twiddle with the dictionary to get axis names to register
c
      ier=fddx_dict(fd_out,'scan')
      ier=fddx_dict(fd_out,'print')

c
c Print modified size info to output dictionary
c
      n1_out=i2e-i2s+1
      n2_out=i1e-i1s+1

      ier=fdds_printf('size.axis(1)','%d\n\0',n1_out)
      ier=fdds_printf('size.axis(2)','%d\n\0',n2_out)

c
c Now actually open the output with a seek
c
      ier8=fdds_lseek8(fd_out,0,izero8,SEEK_SET)

c
c Allocate data arrays
c
      ptr_record_in=fdds_malloc8(dble(n1*n2*SIZEOF_REAL))
      ptr_record_out=fdds_malloc8(dble(n1_out*n2_out*SIZEOF_REAL))

c
c Loop over n3, read-transpose-write
c
      do i3 = 1,n3

         if (lverb .and. mod(i3,10).eq.0) then
            ier=fdds_prtcon('record %5d of %5d\n\0',i3,n3)
         endif

         if (fddx_read(fd_in,record_in,n2).ne.n2) then
            ier=fdds_prterr('Reading record %d\n\0',i3)
            ier=fdds_close(fd_in)
            ier=fdds_close(fd_out)
            ier=fdds_closepr()
            stop
         endif

         call transp2d(record_in, n1,n2, record_out, n1_out,n2_out,
     :                 i1s,i1e, i2s,i2e)

         if (fddx_write(fd_out,record_out,n2_out).ne.n2_out) then
            ier=fdds_prterr('Writing record %d\n\0',i3)
            ier=fdds_close(fd_in)
            ier=fdds_close(fd_out)
            ier=fdds_closepr()
            stop
         endif

      enddo

c
c Clean up
c
      ier=fdds_free(ptr_record_in)
      ier=fdds_free(ptr_record_out)
      ier=fdds_close(fd_in)
      ier=fdds_close(fd_out)
      ier=fdds_closepr()

c
c Stop
c
      stop
      end




c***********************************************************************
c
c subroutine to print help when requested on the command line
c
c***********************************************************************

      subroutine help(prog,Title)

      implicit none

      character prog*(*), Title*(*)
      integer lnblnk

      write(0,*)''
      write(0,*)Title(1:lnblnk(Title))
      write(0,*)''
      write(0,*)' Parameter....Meaning..........................Default'
      write(0,*)''
      write(0,*)'  in=         input dictionary                  stdin:'
      write(0,*)'  out=        output dictionary                stdout:'
      write(0,*)''
      write(0,*)'  1s=         First sample of input to transpose     1'
      write(0,*)'  1e=         Last    "    "    "   "      "        n1'
      write(0,*)'  2s=         First trace  "    "   "      "         1'
      write(0,*)'  2e=         Last    "    "    "   "      "        n2'
      write(0,*)''
      write(0,*)'  verbose=    verbose output                        no'
      write(0,*)''
      write(0,*)'  USAGE: '//prog(1:lnblnk(prog))
      write(0,*)'           in=[in_dict] out=[out_dict]'
      write(0,*)'           [ 1s=[i1s] 1e=[i1e] 2s=[i2s] 2e=[i2e]'
      write(0,*)'             verbose=[y|n] ]'
      write(0,*)''

      stop

      return
      end




c***********************************************************************
c
c transpose a protion of a 2D array
c
c***********************************************************************

      subroutine transp2d
     :           (record_in, n1_in,n2_in, record_out,n1_out,n2_out,
     :            i1s,i1e, i2s,i2e)

      implicit none

c Arguments
      integer n1_in,n2_in, n1_out, n2_out
      integer i1s,i1e, i2s,i2e
      real    record_in(n1_in,n2_in)
      real    record_out(n1_out,n2_out)

c Locals
      integer i1_in,i2_in, i1_out,i2_out

c Go transpose the data
      i1_out=0
      do i2_in = i2s,i2e
         i1_out=i1_out+1
         i2_out=0
         do i1_in = i1s,i1e
            i2_out=i2_out+1
            record_out(i1_out,i2_out) = record_in(i1_in,i2_in)
         enddo
      enddo

      return
      end