!--------( XSLIN : convert AM line data into StdLIN data )--------
  implicit none
  integer :: lwkdir, lrtrim
  real :: hgeoid
  real(8) :: v(30)
  character(80) :: wdr, flog, fnam1, fnam2
  character(150) :: buf;  character(20) :: sdtm, sgmdl, sgmdln
  character(24) :: ca, can;  character(6) :: cu;  character(3) :: cg, cgn
  integer :: ktyp, kg, ka, k1, k2, k3, k4, kmx, ku, kgn, kan
  integer :: ldr, i, lc, l, IER
  real :: calt, valtc
  real(8) :: vlat, vlon, valt, tlat, tlon, talt, wlat, wlon, walt, vmag

  call premsg('XSLIN : convert AM line data into StdLIN 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 premsg('Select AM Line data type :')
  call premsg('  DPAM Line data (1), AMDB-GSJ Line data (2),')
  call premsg('  AMDB-NEDO Line data (3), StdLIN data (4),')
  call gparmi('   or general Random data (0) ==> ', ktyp)
  if ((ktyp < 0) .or. (ktyp > 4)) call abendm('invalid value')
  call gparma('Enter Input AM line data filename ==> ', 81-ldr, fnam1(ldr:80))
  open(1,file=fnam1,status='old')
  call gparmi(' Geodetic ?  in WGS(1) or TokyoDatum(2) ==> ', kg)
  if (kg == 1) then
    cg = 'WGS'
    call gparmi(' Alt. ?  above Ellipsoid(1) or Geoid(2) ==> ', ka)
    if (ka == 1) then
      ca = 'above ITRF-GRS'
    else if (ka == 2) then
      sgmdl = '                    '
      call gparma('  Geoid model ?  0:NGA, or others [Enter NAME]  ==> ',   &
                  20, sgmdl)
      if (sgmdl == '0') sgmdl = 'NGA'
      call sgeoid(sgmdl)
      l = lrtrim(sgmdl)
      if (l <= 12) then
        ca = 'above Geoid ' // sgmdl(1:l)
      else
        ca = 'above Geoid ' // sgmdl(1:9) // '...'
      endif
    else
      call abendm('Alt. reference unknown')
    endif
  else if (kg == 2) then
    cg = 'TYO'
    call premsg(' Input Alt. data is treated as "above Geoid NGA"')
    ca = 'above Geoid NGA';  sgmdl = 'NGA';  call sgeoid(sgmdl)
  else
    call abendm('Geodetic reference unknown')
  endif
  if (ktyp == 1) then
    k1 = 5;  k2 = 6;  k3 = 7;  k4 = 9;  kmx = 9
    ku = 2;  cu = ' Deg. '
  else if (ktyp == 2) then
    k1 = 2;  k2 = 3;  k3 = 0;  k4 = 4;  kmx = 4
    ku = 1;  cu = ' Min. '
    read(1,'(a)') buf
    if ((buf(1:2) /= '##') .or. (buf(27:28) /= 'ft'))   &
      call abendm('seems non AMDB-GSJ data')
    read(buf,'(20x,f6.0)') calt
    valtc = calt * 0.3048
  else if (ktyp == 3) then
    k1 = 3;  k2 = 4;  k3 = 10;  k4 = 8;  kmx = 10
    ku = 1;  cu = ' Min. '
    read(1,'(a)') buf
    if (buf(1:8) /= '##  NEDO') call abendm('seems non AMDB-NEDO data')
  else if (ktyp == 4) then
    k1 = 1;  k2 = 2;  k3 = 3;  k4 = 4;  kmx = 4
    ku = 1;  cu = ' Min. '
  else
!------ 1st data-line display (to answer the following query)
    do
      read(1,'(a)') buf
      call premsg(buf(1:lrtrim(buf)))
      if ((buf(1:1)/='&') .and. (buf(1:1)/='%') .and. (buf(1:1)/='#')) exit
    enddo
    rewind(1)
!------
    call prompt('LatLon unit ?  ')
    call gparmi(' in Min.(1), Deg.(2) or D:Min.(3) ==> ', ku)
    if (ku == 1) then
      cu = ' Min. '
    else if (ku == 2) then
      cu = ' Deg. '
    else if (ku == 3) then
      cu = 'D:Min.'
    else
      call abendm('invalid value')
    endif
    call gparmi2('  data positions of  Lat., Lon. : ', k1,k2)
    if ((k1 <= 0) .or. (k2 <= 0) .or. (k2 == k1))   &
      call abendm('invalid value')
    call gparmi2('  data positions of  Alt., Mag. : ', k3,k4)
    if ((k3 <= 0) .or. (k3 == k1) .or. (k3 == k2) .or.   &
        (k4 <= 0) .or. (k4 == k1) .or. (k4 == k2) .or. (k4 == k3))   &
      call abendm('invalid value')
    if (ku == 3) then
      if (((k1+1) == k2) .or. ((k1+1) == k3) .or. ((k1+1) == k4) .or.   &
          ((k2+1) == k1) .or. ((k2+1) == k3) .or. ((k2+1) == k4))   &
        call abendm('invalid value')
      kmx = max(k1+1, k2+1, k3, k4)
    else
      kmx = max(k1, k2, k3, k4)
    endif
    if (kmx > 30) call abendm('data position unexpected')
  endif

  call gparma('Enter Output StdLIN data filename ==> ', 81-ldr, fnam2(ldr:80))
  open(2,file=fnam2,status='new')
  call gparmi(' Geodetic ?  in WGS(1) or TokyoDatum(2) ==> ', kgn)
  if (kgn == 1) then
    cgn = 'WGS'
    call prompt(' Output Alt. ?  above Ellipsoid(1) or Geoid(2)')
    call gparmi(' ==> ', kan)
    if (kan.eq.1) then
      can = 'above ITRF-GRS  '
    else if (kan.eq.2) then
      sgmdln = '                    '
      call gparma('  Geoid model ?  0:NGA, or others [Enter NAME]  ==> ',   &
                  20, sgmdln)
      if (sgmdln.eq.'0') sgmdln = 'NGA'
      call sgeoid(sgmdln)
      l = lrtrim(sgmdln)
      if (l.le.12) then
        can = 'above Geoid ' // sgmdln(1:l)
      else
        can = 'above Geoid ' // sgmdln(1:9) // '...'
      endif
    else
      call abendm('Alt. reference unknown')
    endif
  else if (kgn == 2) then
    cgn = 'TYO'
    call premsg(' Output Alt. data is "above Geoid NGA"')
    can = 'above Geoid NGA';  sgmdln = 'NGA';  call sgeoid(sgmdln)
  else
    call abendm('Geodetic reference unknown')
  endif
  write(6,'(a)') '----------------------------------------'
  write(6,'(a)') 'XSLIN : convert AM line data into StdLIN data'
  write(6,'(a,a)') ' Input AM line data filename : ', fnam1(1:lrtrim(fnam1))
  write(6,'(8x,i2,5a,i2,2a)')   &
     k1, ': ', cg, '-Lat.(', cu, ')   ', k3, ': Alt.(m) ', ca,   &
     k2, ': ', cg, '-Lon.(', cu, ')   ', k4, ': Mag.(nT)', ''
  write(6,'(a,a)') ' Output StdLIN data filename : ', fnam2(1:lrtrim(fnam2))
  write(6,'(8x,4a)') 'Lat./Lon. in ', cgn, ';  Alt. ', can
  write(2,'(a)')  '# output from XSLIN'
  write(2,'(4a)') '#   Lat./Lon. in ', cgn, ';  Alt. ', can
  write(2,'(a,a)') '# org. line data filename : ', fnam1(1:lrtrim(fnam1))

  do
    read(1,'(a)',iostat=IER) buf
    if (IER < 0) exit
    if ((buf(1:1)=='&') .or. (buf(1:1)=='%') .or. (buf(1:1)=='#')) then
      call premsg(buf(1:lrtrim(buf)))
      write(2,'(a)') buf(1:lrtrim(buf))
      cycle
    endif
    if (ktyp == 4) then
      i = index(buf,'N');  if (i /= 0) buf(i:i) = ' '
      i = index(buf,'E');  if (i /= 0) buf(i:i) = ' '
      i = index(buf,'m');  if (i /= 0) buf(i:i) = ' '
      i = index(buf,'nT');  if (i /= 0) buf(i:i+1) = '  '
    else if (ktyp == 0) then
      do
        lc = index(buf,':');  if (lc == 0) exit
        buf(lc:lc) = ' '
      enddo
    endif
    read(buf,*) (v(i),i=1,kmx)
    vmag = v(k4)
    if (ku == 1) then
      vlat = v(k1)
      vlon = v(k2)
    else if (ku == 2) then
      vlat = v(k1) * 60.
      vlon = v(k2) * 60.
    else
      vlat = v(k1)*60. + v(k1+1)
      vlon = v(k2)*60. + v(k2+1)
    endif
    if (ktyp == 2) valt = valtc
    if (ktyp /= 2) valt = v(k3)
!           [Let vlat/vlon/valt be in Output side spec]
    if (kg == 1) then
      if (kgn == 1) then
!               (vlat/vlon are in WGS, common IN/OUT)
        if (ka == 1) then
          if (kan == 1) then
!               [valt is above WGS Ellips, common IN/OUT]
          else
            valt = valt - hgeoid(sngl(vlat), sngl(vlon))
!               [now valt is above WGS Geoid (conv. from WGS Ellips)]
          endif
        else
          if (kan == 1) then
            valt = valt + hgeoid(sngl(vlat), sngl(vlon))
!               [now valt is above WGS Ellips (conv. from WGS geoid)]
          else
            if (sgmdl.ne.sgmdln) then
              call sgeoid(sgmdl)
              valt = valt + hgeoid(sngl(vlat), sngl(vlon))
              call sgeoid(sgmdln)
              valt = valt - hgeoid(sngl(vlat), sngl(vlon))
!               [now valt is above WGS Geoid (conv. from another model)]
            endif
          endif
        endif
      else
!               (Convert vlat/vlon from WGS to TYO)
        if (ka == 1) then
          valt = valt - hgeoid(sngl(vlat), sngl(vlon))
!               [now valt is above Geoid [TYO] (conv. from WGS Ellips)]
        else
          if (sgmdl.ne.sgmdln) then
            call sgeoid(sgmdl)
            valt = valt + hgeoid(sngl(vlat), sngl(vlon))
            call sgeoid(sgmdln)
            valt = valt - hgeoid(sngl(vlat), sngl(vlon))
!               [now valt is above Geoid [TYO] (conv. from another model)]
          endif
        endif
        call xw84td(vlat, vlon, valt, tlat, tlon, talt)
        vlat = tlat;  vlon = tlon
!               [now vlat/vlon are in TYO (conv. from WGS)]
      endif
    else
      if (kgn == 1) then
!           (Convert vlat/vlon from TYO to WGS)
        call xtw84d(vlat, vlon, valt, wlat, wlon, walt)
        vlat = wlat;  vlon = wlon
!               [now vlat/vlon are in WGS (conv. from TYO)]
        if (kan == 1) then
          valt = valt + hgeoid(sngl(vlat), sngl(vlon))
!               [now valt is above WGS Ellips (conv. from Geoid [TYO])
        else
          if (sgmdl.ne.sgmdln) then
            call sgeoid(sgmdl)
            valt = valt + hgeoid(sngl(vlat), sngl(vlon))
            call sgeoid(sgmdln)
            valt = valt - hgeoid(sngl(vlat), sngl(vlon))
!               [now valt is above WGS Geoid (conv. from Geoid [TYO])]
          endif
        endif
      else
!               (vlat/vlon are in TYO, common IN/OUT)
!               (valt is above Geoid [TYO}, common IN/OUT)
      endif
    endif
    write(2,'(1x,f11.5,a2,f11.5,a2,f8.2,a2,f8.2,a2)')   &
      vlat, 'N ', vlon, 'E ', valt, 'm ', vmag, 'nT'
  enddo
  close(2)
  close(1)

  if (flog(ldr:ldr) /= ' ') then
    open(99,file=flog,access='append')
    write(99,'(//a)') 'XSLIN : convert AM line data into StdLIN data'
    write(99,'(a,a)') ' Input AM line data filename : ', fnam1(1:lrtrim(fnam1))
    write(99,'(8x,i2,5a,i2,2a)')   &
       k1, ': ', cg, '-Lat.(', cu, ')   ', k3, ': Alt.(m) ', ca,   &
       k2, ': ', cg, '-Lon.(', cu, ')   ', k4, ': Mag.(nT)', ''
    write(99,'(a,a)') ' Output StdLIN data filename : ', fnam2(1:lrtrim(fnam2))
    write(99,'(8x,4a)') 'Geodetic in ', cgn, ';  Alt. ', can
    call strdtm(sdtm)
    write(99,'(a,a,a)') '=============== ', sdtm, '==============='
    close(99)
  endif
  stop
  end