!--------( PLTRK : plot trackline path from AM line data )--------
  implicit none
  integer :: lwkdir, lrtrim
  external cviken
  real :: v(30)
  character(80) :: wdr, flog, fnin, fnam
  character(150) :: buf;  character(20) :: sdtm
  character(8) :: area, aval;  character(6) :: cu
  character(4) :: anc;  character(3) :: cg, spl;  character(1) :: yn
  integer :: ldr, nc, ncc, iorg, korg, ispa, ispb, IER, i1, i2
  integer :: ktyp, kg, ku, k1, k2, kmx, mc, kpl, kscl, klll, kshl, k
  integer :: km, iinc, kinc, kp, lc, i, isf, itf, ksf, ktf, idf, kdf
  real :: xs, ys, rx, ry, fscl, high, wide, ccm, sclr, sclu, dscl
  real :: vlat, vlon, tlat, tlon, wlat, wlon, ye, xn, px, py, alt, x, y
  real :: fi1, fi2, fi3, fi4, fk1, fk2, fk3, fk4

  call premsg('PLTRK : plot trackline path from AM line data')
  call opnpin()
  ldr = lwkdir(50,wdr) + 1
  flog = wdr;  fnin = wdr;  fnam = 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 line/random data filename ==> ',   &
              81-ldr, fnin(ldr:80))
  open(1,file=fnin,status='old')
  call gparmi(' Lat./Lon. ?   in WGS(1) or TokyoDatum(2) ==> ', kg)
  if (kg == 1) then
    cg = 'WGS'
  else if (kg == 2) then
    cg = 'TYO'
  else
    call abendm('invalid value')
  endif
  ku = 1;  cu = ' Min. '
  if (ktyp == 1) then
    k1 = 5;  k2 = 6;  kmx = 6;  ku = 2;  cu = ' Deg. '
  else if (ktyp == 2) then
    k1 = 2;  k2 = 3;  kmx = 3
    read(1,'(a)') buf
    if ((buf(1:2) /= '##') .or. (buf(27:28) /= 'ft'))   &
      call abendm(' seems non AMDB-GSJ data')
  else if (ktyp == 3) then
    k1 = 3;  k2 = 4;  kmx = 4
    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;  kmx = 2
  else
    do
      read(1,'(a)') buf
      call premsg(buf(1:lrtrim(buf)))
      if ((buf(1:1) == '&') .or. (buf(1:1) == '%')) cycle
      if (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')
    kmx = max(k1, k2)
    if (ku == 3) then
      if (((k1+1) == k2) .or. (k1 == (k2+1))) call abendm('invalid value')
      kmx = max(k1+1, k2+1)
    endif
    if (kmx > 30) call abendm('too large data position')
  endif

  call gparma('Enter Output PS filename  ==> ', 81-ldr, fnam(ldr:80))
  if (fnam(ldr:ldr) == ' ') fnam(1:1) = char(0)
  if (fnam(80:80) /= ' ') call abendm('too long filename')
  call gparma('  Areaname (8 chrs max) : ', 8, area)
  call gparmi('Select Coordinate Number (+200 for WGS) ==> ', ncc)
  if (ncc >= 800) then
    mc = 0;  nc = ncc - 800
  else
    mc = 1;  nc = ncc
  endif
  ispa = 0;  ispb = 0
  if ((nc >= 1) .and. (nc <= 62)) then
    iorg = 0
    korg = (nc-30)*360 - 180
    if (nc > 60) korg = 0
    if (nc == 61) iorg = +5400
    if (nc == 62) iorg = -5400
  else if (nc == 65) then
    iorg = 0
    call gparmi2('  Central Longitude ?  (Deg. Min.) : ', i1,i2)
    korg = i1*60 + i2
  else
    if (nc == 72) call abendm('nc=72 not implemented')
    call gparmi2('  Origin Latitude ?   (Deg. Min.) : ', i1,i2)
    iorg = i1*60 + i2
    call gparmi2('  Origin Longitude ?  (Deg. Min.) : ', i1,i2)
    korg = i1*60 + i2
  endif
  call premsg('  Southwest corner Coord. in km ?')
  call gparmf('    Northing : ', xs)
  call gparmf('    Easting  : ', ys)
  call premsg('  Area size in km ?')
  call gparmf('    N-S range : ', rx)
  if (rx <= 0.) call abendm('invalid value')
  call gparmf('    E-W range : ', ry)
  if (ry <= 0.) call abendm('invalid value')
  call gparmi('  Select  1 (Portrait) or 2 (Landscape) ==> ', kpl)
  if (kpl == 1) then
    spl = 'A4P'
  else if (kpl == 2) then
    spl = 'A4L'
  else
    call abendm('invalid value')
  endif
  call psopn(fnam, spl)
  if (kpl == 1) then
    call gparmf('  Width of Drawing (in cm) ?  ==> ', wide)
    fscl = ry*1000. / wide
    high = rx*1000. / fscl
  else
    call gparmf('  Height of Drawing (in cm) ?  ==> ', high)
    fscl = rx*1000. / high
    wide = ry*1000. / fscl
  endif
  fscl = fscl * 100.
  ccm = fscl/100000.
  kscl = 0;  sclr = 0.;  sclu = 0.;  klll = 0;  kshl = 0
  km = nint(fscl)/25000
  if (km < 2) then;  km = 1
  else if (km < 5) then;  km = 2
  else if (km < 10) then;  km = 5
  else if (km < 20) then;  km = 10
  else if (km < 50) then;  km = 20
  else;  km = 50
  endif
  iinc = km;  kinc = km + km/2
  if (km == 50) iinc = 60
  if ((km == 1) .or. (km == 5) .or. (km == 50)) kinc = iinc
  dscl = float(km)/ccm
  call prompt('  Draw  Scale Bar  ?  ')
  call gparma('  "y" for Yes, othewise No ==> ', 1, yn)
  if ((yn == 'y') .or. (yn == 'Y')) then
    kscl = 1
    write(anc,'(f4.1)') dscl
    call premsg('      size :' // anc // 'cm Wide, 0.8cm High')
    call prompt('    Lower-Left Pos. ? (right and up)')
    call gparmf2(' in cm ==> ', sclr, sclu)
  endif
  call prompt('  Draw Lat.Lon.Lines ?')
  call gparma('  "y" for Yes, othewise No ==> ', 1, yn)
  if ((yn == 'y') .or. (yn == 'Y')) klll = 1
  call prompt('  Draw Coastlines ?   ')
  call gparma('  "y" for Yes, othewise No ==> ', 1, yn)
  if ((yn == 'y') .or. (yn == 'Y')) then
    kshl = 2
    call prompt('       Rivers  ?      ')
    call gparma('  "y" for Yes, othewise No ==> ', 1, yn)
    if ((yn == 'y') .or. (yn == 'Y')) kshl = 12
    call prompt('       Pref.Boundary ?')
    call gparma('  "y" for Yes, othewise No ==> ', 1, yn)
    if ((yn == 'y') .or. (yn == 'Y')) then
      kshl = kshl + 200
      if (kshl == 202) kshl = 102
    endif
  endif
  call clspin()
  write(6,'(a)') '----------------------------------------'
  write(6,'(a)') 'PLTRK : plot trackline path from AM line data'
  write(6,'(a,a)') ' Input AM line data filename : ', fnin(1:lrtrim(fnin))
  write(6,'(8x,i2,5a,3x,i2,5a)')   &
    k1, ': ', cg, '-Lat.(', cu, ')', k2, ': ', cg, '-Lon.(', cu, ')'
  write(6,'(a,a)') ' Output PS filename  : ', fnam(1:lrtrim(fnam))
  write(6,'(a12,a8,6x,a6,i4,4x,4i8)')   &
    'Areaname :  ', area, 'Coord.', ncc, iorg,korg, ispa,ispb
  write(6,'(12x,a,f8.2,a,f8.2)') 'SW corner :      x =', xs, '   y =', ys
  write(6,'(12x,a,f8.2,a,f8.2)') '  area size :    NS=', rx, '   EW=', ry
  write(6,'(6x,a,f8.0,4x,a3,f8.1,a1,f8.1,a1)')   &
    'scale  1 /', fscl, spl, high, 'H', wide, 'W'
  if ((kscl /= 0) .or. (klll /= 0) .or. (kshl /= 0)) then
    write(6,'(a,$)') '      with '
    if (kshl /= 0) write(6,'(a,$)') '  Coast'
    if (mod(kshl/10,10) /= 0) write(6,'(a,$)') '+River'
    if ((kshl/100) /= 0) write(6,'(a,$)') '+PrefB'
    if (klll /= 0) write(6,'(a,$)') '  LatLonLines'
    if (kscl /= 0) then
      write(6,'(a,f5.1,$)') '  ScaleBar [', sclr
      write(6,'(a,f5.1,a,$)') ',', sclu, ' ]'
    endif
    write(6,'(a)')
  endif

  call cvinit(ncc,float(iorg),float(korg),float(ispa),float(ispb))
  call plots(1.5, 1.5)
  call wrect(0., 0., wide, high)
  call lstyle('T', 0.4, 0., 0, -255)
  write(aval,'(f8.1)') xs
  call ptext(aval, 8, -0.1, 0., 1)
  write(aval,'(f8.1)') xs+rx
  call ptext(aval, 8, -0.1, high-0.28, 1)
  write(aval,'(f8.1)') ys
  call ptext(aval, 8, 0.7, -0.4, 1)
  write(aval,'(f8.1)') ys+ry
  call ptext(aval, 8, wide, -0.4, 1)
  kp = 3
  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)))
      kp = 3;  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) = ' '
      read(buf,*) v(1), v(2)
    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)
    endif
    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 ((kg == 1) .and. (ncc >= 800)) then
      call xw84t(vlat, vlon, 0., tlat, tlon, alt)
      call cviken(tlat, tlon, ye, xn)
    else if ((kg == 2) .and. (ncc < 800)) then
      call xtw84(vlat, vlon, 0., wlat, wlon, alt)
      call cviken(wlat, wlon, ye, xn)
    else
      call cviken(vlat, vlon, ye, xn)
    endif
    px = (ye-ys) / ccm;  py = (xn-xs) / ccm
    call plot(px, py, kp)
    kp = 2
  enddo
  close(1)
  if (kscl /= 0) then
    call plot(sclr,      sclu+0.1, 3)
    call plot(sclr+dscl, sclu+0.1, 2)
    call newpen(2)
    call plot(sclr+dscl, sclu+0.3, 3)
    call plot(sclr+dscl, sclu,     2)
    call plot(sclr,      sclu,     2)
    call plot(sclr,      sclu+0.3, 2)
    write(anc,'(i2,a2)') km, 'km'
    call lstyle('T', 0.5, 0., 0, -255)
    call ptext(anc, 4, sclr+dscl/2., sclu+0.5, 2)
    call newpen(1)
  endif
  call scisor(0., 0., wide, high)
  if ((klll /= 0) .or. (kshl /= 0)) then
    call cvinit(ncc,float(iorg),float(korg),float(ispa),float(ispb))
    call cvenik(ys,    xs,    fi1, fk1)
    call cvenik(ys+ry, xs,    fi2, fk2)
    call cvenik(ys,    xs+rx, fi3, fk3)
    call cvenik(ys+ry, xs+rx, fi4, fk4)
    isf = nint(amin1(fi1,fi2) - 0.5)
    itf = nint(amax1(fi3,fi4) + 0.5)
    ksf = nint(amin1(fk1,fk3) - 0.5)
    ktf = nint(amax1(fk2,fk4) + 0.5)
  endif
  if (klll /= 0) then
    do i=isf/iinc,itf/iinc
      idf = i * iinc
      if ((idf < isf) .or. (idf > itf)) cycle
      call cviken(float(idf), float(ksf), y, x)
      call plot((y-ys)/ccm, (x-xs)/ccm, 3)
      do kdf=ksf+1,ktf
        call cviken(float(idf), float(kdf), y, x)
        call plot((y-ys)/ccm, (x-xs)/ccm, 2)
      enddo
    enddo
    do k=ksf/kinc,ktf/kinc
      kdf = k * kinc
      if ((kdf < ksf) .or. (kdf > ktf)) cycle
      call cviken(float(isf), float(kdf), y, x)
      call plot((y-ys)/ccm, (x-xs)/ccm, 3)
      do idf=isf+1,itf
        call cviken(float(idf), float(kdf), y, x)
        call plot((y-ys)/ccm, (x-xs)/ccm, 2)
      enddo
    enddo
  endif
  if (kshl /= 0) then
    call rshore(mc, isf/60, (itf+59)/60, ksf/60, (ktf+59)/60)
    call pshore(-ys/ccm, -xs/ccm, kshl, ccm, cviken)
  endif
  call scisor(0., 0., 0., 0.)
  call lstyle('T', 0.4, 0., 0, -255)
  call ptext(area, 8, -1., -1.5, 0)
  call ptext(fnin, 50, 2., -1.5, 0)
  call plote()
  call pscls()

  if (flog(ldr:ldr) /= ' ') then
    open(99,file=flog,access='append')
    write(99,'(//a,$)') 'PLTRK :'
    write(99,'(a)')   ' plot trackline path from AM line data'
    write(99,'(a,a)') ' Input AM line data filename : ', fnin(1:lrtrim(fnin))
    write(99,'(8x,i2,5a,3x,i2,5a)')   &
      k1, ': ', cg, '-Lat.(', cu, ')', k2, ': ', cg, '-Lon.(', cu, ')'
    write(99,'(a,a)') ' Output PS filename  : ', fnam(1:lrtrim(fnam))
    write(99,'(a12,a8,6x,a6,i4,4x,4i8)')   &
      'Areaname :  ', area, 'Coord.', ncc, iorg,korg,ispa,ispb
    write(99,'(12x,a,f8.2,a,f8.2)') 'SW corner :      x =', xs, '   y =', ys
    write(99,'(12x,a,f8.2,a,f8.2)') '  area size :    NS=', rx, '   EW=', ry
    write(99,'(6x,a,f8.0,4x,a3,f8.1,a1,f8.1,a1)')   &
      'scale  1 /', fscl, spl, high, 'H', wide, 'W'
    if ((kscl /= 0) .or. (klll /= 0) .or. (kshl /= 0)) then
      write(99,'(a,$)') '      with '
      if (kshl /= 0) write(99,'(a,$)') '  Coast'
      if (mod(kshl/10,10) /= 0) write(99,'(a,$)') '+River'
      if ((kshl/100) /= 0) write(99,'(a,$)') '+PrefB'
      if (klll /= 0) write(99,'(a,$)') '  LatLonLines'
      if (kscl /= 0) then
        write(99,'(a,f5.1,$)') '  ScaleBar [', sclr
        write(99,'(a,f5.1,a,$)') ',', sclu, ' ]'
      endif
      write(99,'(a)')
    endif
    call strdtm(sdtm)
    write(99,'(a,a,a)') '=============== ', sdtm, '==============='
    close(99)
  endif
  stop
  end