!--------( RAERX : extract a rear set data        )--------
!--------(           from multiple sets data file )--------
  implicit none
  integer :: lwkdir, lrtrim
  real, allocatable :: f(:,:)
  character(80) :: wdr, flog, fnam1, fnam2, buf
  character(72) :: out;  character(20) :: sdtm;  character(10) :: anul
  character(8) :: area, salt
  integer :: ldr, i, j, ns, n
  integer :: nc, ncc, iorg, korg, ispa, ispb
  integer :: ixs, iys, mszx, mszy, mx, my
  real :: vnul, alt, valt, xs, ys, dmy

  call premsg('REARX : extract a RAER set data')
  call opnpin()
  ldr = lwkdir(50,wdr) + 1
  flog = wdr;  fnam1 = wdr;  fnam2 = wdr
  call gparma('Enter LOG filename ==> ', 81-ldr, flog(ldr:80))
  call gparma('Enter input Source data filename ==> ',   &
              81-ldr, fnam1(ldr:80))
  open(1,file=fnam1,status='old')
  call prompt('  Select Seq.No. of the SET to extract ')
  call gparmi('(1: Topmost set) ==> ', ns)
  if (ns <= 0) call abendm('illegal set number')
  do n=1,ns-1
    read(1,'(a)') buf
    do while (buf(1:1) == '#')
      read(1,'(a)') buf
    enddo
    read(1,'(a)') buf
    read(buf,*) ixs, iys, mszx, mszy, mx, my, vnul, alt
    read(1,*) ((dmy,i=1,mx),j=1,my)
  enddo
  read(1,'(a)') buf
  do while (buf(1:1) == '#')
    read(1,'(a)') buf
  enddo
  if (buf(13:16) == '    ') then        !-- v2005
    read(buf,'(a8,i4,4x,4i8)') area, ncc, iorg, korg, ispa, ispb
    if ((ncc<0) .or. (ncc>=400)) call abendm('invalid coord.no.')
    nc = ncc
    if (ncc < 200) then
      nc = ncc; ncc = nc + 800
    else
      ncc = ncc - 200; nc = ncc
    endif
  else                                  !-- v2018
    read(buf,'(a8,4x,i4,4i8)') area, ncc, iorg, korg, ispa, ispb
    nc = ncc
    if ((ncc >= 800) .and. (ncc < 1000)) nc = ncc - 800
  endif
  if ((nc >= 1) .and. (nc <= 62)) then
    iorg = 0;  ispa = 0;  ispb = 0
    korg = (nc-30)*360 - 180
    if (nc > 60) korg = 0
    if (nc == 61) iorg = +5400
    if (nc == 62) iorg = -5400
  endif
  read(1,'(a)') buf
  read(buf,*) ixs, iys, mszx, mszy, mx, my, vnul, alt
  write(salt,'(1x,f7.0)') alt
  valt = alt
  if (alt <= 0.) then
    valt = -1.
    if (alt == 0.) salt = '    var.'
    if (alt < 0.) salt = '   undef'
  endif
  xs = float(ixs)/1000.
  ys = float(iys)/1000.

  allocate(f(mx,my))

  write(anul,'(f9.2)') vnul
  if ((anul(9:9) /= '0') .or. (anul(1:1) /= ' '))   &
    call abendm('unable to handle vnul value')
  write(out,'(a14,2a10,2x,2a8,2a6,a8)')   &
     '---   (mesh)  ', 'xs', 'ys', 'msz-x', 'msz-y', 'mx', 'my', 'alt '
  call premsg(out)
  write(out,'(a16,2f10.2,2i8,2i6,a8)')   &
     '---   source    ', xs, ys, mszx, mszy, mx, my, salt
  call premsg(out)
  read(1,*) ((f(i,j),i=1,mx),j=1,my)
  close(1)

  call premsg('')
  call gparma('Enter output new data filename ==> ', 81-ldr, fnam2(ldr:80))
  call clspin()
  open(2,file=fnam2,status='new')
  write(6,'(a)') '----------------------------------------'
  write(6,'(a)') 'REARX : extract a REAR set data'
  write(6,'(a,a)')  ' Source infile   : ', fnam1(1:lrtrim(fnam1))
  write(6,'(a,i3)') '  Seq.No. of the set to extract :', ns
  write(6,'(a,a)')  '    areaname : ', area
  write(6,'(a,a)')  ' Output new file : ', fnam2(1:lrtrim(fnam2))
  write(6,'(a12,a8,6x,a6,i4,4x,4i8)')   &
    'Areaname :  ', area, 'Coord.', ncc, iorg,korg, ispa,ispb
  write(6,'(4x,a)')   &
    '  (mesh)          xs        ys     msz-x   msz-y    mx    my    alt '
  write(6,'(4x,a12,2f10.2,2i8,2i6,a8)')   &
    'grid out=in ', xs, ys, mszx, mszy, mx, my, salt
  if (alt == 0.) then
    write(6,'(35x,a)') 'var. alt inf. is lost (set to undef)'
  endif

  write(2,'(a)')    '# prog.REARX applied'
  write(2,'(a,a)')  '# Source: ', fnam1(1:lrtrim(fnam1))
  write(2,'(a,i3)') '#  Seq.No. of the set to extract :', ns
  write(2,'(a8,4x,i4,4i8)') area, ncc, iorg,korg, ispa,ispb
  write(2,'(2i12,4i6,1x,f7.1,1x,f7.0)') ixs,iys, mszx,mszy, mx,my, vnul, valt
  do j=1,my
    write(2,'((f7.1,9(1x,f7.1)))') (f(i,j),i=1,mx)
  enddo
  close(2)

  if (flog(ldr:ldr) /= ' ') then
    open(99,file=flog,access='append')
    write(99,'(//a)')  'REARX : extract a REAR set data'
    write(99,'(a,a)')  ' Source infile  : ', fnam1(1:lrtrim(fnam1))
    write(99,'(a,i3)') '  Seq.No. of the set to extract :', ns
    write(99,'(a,a)')  '    areaname : ', area
    write(99,'(a,a)')  ' Output new file: ', fnam2(1:lrtrim(fnam2))
    write(99,'(a12,a8,6x,a6,i4,4x,4i8)')   &
      'Areaname :  ', area, 'Coord.', ncc, iorg,korg, ispa,ispb
    write(99,'(4x,a)')   &
      '  (mesh)          xs        ys     msz-x   msz-y    mx    my    alt '
    write(99,'(4x,a12,2f10.2,2i8,2i6,a8)')   &
      'grid out=in ', xs, ys, mszx, mszy, mx, my, salt
    if (alt == 0.) then
      write(99,'(35x,a)') 'var. alt inf. is lost (set to undef)'
    endif
    call strdtm(sdtm)
    write(99,'(a,a,a)') '=============== ', sdtm, '==============='
    close(99)
  endif
  stop
  end