c******************************************************
c   XLDPN : Extract Line data  from alog2asc output   *
c              with PostNav Pos. data interporation   *
c                  and calculation of IGRF-Residual   *
c******************************************************
      dimension fid(300), fx(300), fy(300), fz(300)
      real*8 tmag(300), wtm, tms, tme, tmsr, wtmr
      character*80 wdr, flog, finf, fnam, fpos, fout
      character s*44, ss*60, line*150, str*80, name*8, sdtm*20
      call premsg('XLDPN : Extract Line data  from alog2asc output')
      call premsg('           with PostNav Pos. data interporation')
      call premsg('               and calculation of IGRF-Residual')
      call opnpin()
      ldr = lwkdir(50,wdr) + 1
      flog = wdr
      finf = wdr
      fnam = wdr
      fpos = wdr
      fout = wdr
      call gparma('Enter LOG filename ==> ', 81-ldr, flog(ldr:80))
      call gparma('Enter Lines.Inf. filename ==> ',
     *            81-ldr, finf(ldr:80))
      open(10,file=finf,status='old')
      call gparma('Enter Output filename ==> ', 81-ldr, fout(ldr:80))
      open(2,file=fout,status='new')
      call gparmf(' Year of survey  : ', year)
      call gparmi(' Generation number of IGRF : ', ngn)
      call gparmi(' Thin-out Ratio N (crop 1 out of N) ==> ', intv)
      if (intv.le.0) intv = 1
      call clspin()
      write(6,'(a)') '----------------------------------------'
      write(6,'(a)') 'XLDPN : Extract Line data  from alog2asc output'
      write(6,'(8x,a)') '   with PostNav Pos. data interporation'
      write(6,'(8x,a)') '       and calculation of IGRF-Residual'
      write(6,'(2a)') ' Lines.Inf. file : ', finf(1:lrtrim(finf))
      write(6,'(2a)') ' Output filename : ', fout(1:lrtrim(fout))
      write(6,'(a,f7.2)') '   Year of survey  : ', year
      write(6,'(a,i2)')   '   IGRF generation : ', ngn
      write(6,'(a,i3)')   '  Data thin-out ratio :', intv
      call gigrf(ngn, year)
      open(9,status='scratch')
c
      read(10,'(a)') str
      if (str(1:1).ne.'=') call abendm('input file not selected')
   40 is = 1
   41 is = is + 1
      if (is.gt.80) call abendm('Lines.Inf. data illegal')
      if (str(is:is).eq.' ') goto 41
      it = is
   42 it = it + 1
      if (it.gt.80) call abendm('Lines.Inf. data illegal')
      if (str(it:it).ne.' ') goto 42
      fnam(ldr:80) = str(is:it-1)
      is = it
   43 is = is + 1
      if (is.gt.80) call abendm('Lines.Inf. data illegal')
      if (str(is:is).eq.' ') goto 43
      it = is
   44 it = it + 1
      if (it.gt.80) call abendm('Lines.Inf. data illegal')
      if (str(it:it).ne.' ') goto 44
      fpos(ldr:80) = str(is:it-1)
      write(6,'(2a)') ' Input AM filename  : ', fnam(1:lrtrim(fnam))
      open(1,file=fnam,status='old')
      write(9,'(2a)') ' Input AM filename  : ', fnam(1:lrtrim(fnam))
      write(6,'(2a)') ' PNAV data filename : ', fpos(1:lrtrim(fpos))
      open(20,file=fpos,status='old')
      write(9,'(2a)') ' PNAV data filename : ', fpos(1:lrtrim(fpos))
    1 read(1,'(a)') line
      if (line(1:10).ne.'/DateTime:') goto 1
      read(line,'(10x,i9)') iymd
    2 read(1,'(a)') line
      if (line(1:1).eq.'/') goto 2
      read(20,'(a)') ss
      k = index(ss, ':')
      if (k.eq.0 .or.
     *    ss(k+3:k+3).ne.':' .or.
     *    ss(k+6:k+6).ne.'.') call abendm('PNAV data illegal format')
      ss(k:k) = ' '
      ss(k+3:k+3) = ' '
      read(ss,*) ih, im, sec, fidoa, fkdoa, falta
      if (ih.lt.0.or.ih.ge.24 .or.
     *    im.lt.0.or.im.ge.60 .or.
     *    sec.lt.0..or.sec.ge.60.)
     *     call abendm('PNAV data illegal format')
      utca = float(ih*3600 + im*60) + sec
      if (utca.ge.54000.) utca = utca - 86400.
      read(20,'(a)') ss
      k = index(ss, ':')
      if (k.eq.0 .or.
     *    ss(k+3:k+3).ne.':' .or.
     *    ss(k+6:k+6).ne.'.') call abendm('PNAV data illegal format')
      ss(k:k) = ' '
      ss(k+3:k+3) = ' '
      read(ss,*) ih, im, sec, fidob, fkdob, faltb
      utcb = float(ih*3600 + im*60) + sec
      if (utcb.ge.54000.) utcb = utcb - 86400.
      if (utcb.le.utca) call abendm('PNAV data UTC illegal sequence')
      kam = 0
c
   30 read(10,'(a)',end=3) str
      if (str(1:1).eq.'=') then
        close(1)
        close(20)
        goto 40
      endif
      is = 0
   31 is = is + 1
      if (is.gt.80) call abendm('Lines.Inf. data illegal')
      if (str(is:is).eq.' ') goto 31
      it = is
   32 it = it + 1
      if (it.gt.80) call abendm('Lines.Inf. data illegal')
      if (str(it:it).ne.' ') goto 32
      name = str(is:it-1)
      read(str(it:80),*) tms, tme
      write(str,'(a8,2x,f9.2,1x,f9.2)') name, tms, tme
      call premsg('   Line: ' // str(1:29))
      write(2,'(a1,a)') '&', str(1:29)
      kst = 0
      kintv = 1
      if (kam.ne.0) goto 15
c
   10 read(1,'(f7.1,22x,f10.3,3f8.3,40x,a44)',end=19)
     *     fid(1),tmag(1),fx(1),fy(1),fz(1),s
      if (s(41:41).ne.'1') goto 10
      read(s,'(f9.5)') utc1
      utc1 = utc1 * 3600.
      if (utc1.ge.54000.) utc1 = utc1 - 86400.
      if (utc1.lt.utca) goto 10
   20 if (utc1.le.utcb) goto 21
      utca = utcb
      fidoa = fidob
      fkdoa = fkdob
      falta = faltb
      read(20,'(a)',end=29) ss
      k = index(ss, ':')
      if (k.eq.0 .or.
     *    ss(k+3:k+3).ne.':' .or.
     *    ss(k+6:k+6).ne.'.') call abendm('PNAV data illegal format')
      ss(k:k) = ' '
      ss(k+3:k+3) = ' '
      read(ss,*) ih, im, sec, fidob, fkdob, faltb
      utcb = float(ih*3600 + im*60) + sec
      if (utcb.ge.54000.) utcb = utcb - 86400.
      if (utcb.le.utca) call abendm('PNAV data UTC illegal sequence')
      goto 20
   21 f = (utc1-utca) / (utcb-utca)
      fido1 = fidoa + (fidob-fidoa)*f
      fkdo1 = fkdoa + (fkdob-fkdoa)*f
      falt1 = falta + (faltb-falta)*f
      m = 0
      goto 12
c
   11 fid(1) = fid(n)
      tmag(1) = tmag(n)
      fx(1) = fx(n)
      fy(1) = fy(n)
      fz(1) = fz(n)
      utc1 = utcn
      fido1 = fidon
      fkdo1 = fkdon
      falt1 = faltn
      m = 1
   12 do 100 i=2,300
        read(1,'(f7.1,22x,f10.3,3f8.3,40x,a44)',end=19)
     *       fid(i),tmag(i),fx(i),fy(i),fz(i),s
        if (s(41:41).eq.'1') goto 13
  100 continue
      write(str(1:8),'(f8.1)') fid(1)
      call premsg('--- GPS data lost after Fid.=' // str(1:8))
      if (kst.ne.0) call abendm('unable to recover position')
      goto 10
   13 n = i
      read(s,'(f9.5)') utcn
      utcn = utcn * 3600.
      if (utcn.ge.54000.) utcn = utcn - 86400.
      if (utcn.le.utc1) call abendm('Input data UTC illegal sequence')
   23 if (utcn.le.utcb) goto 24
      utca = utcb
      fidoa = fidob
      fkdoa = fkdob
      falta = faltb
      read(20,'(a)',end=29) ss
      k = index(ss, ':')
      if (k.eq.0 .or.
     *    ss(k+3:k+3).ne.':' .or.
     *    ss(k+6:k+6).ne.'.') call abendm('PNAV data illegal format')
      ss(k:k) = ' '
      ss(k+3:k+3) = ' '
      read(ss,*) ih, im, sec, fidob, fkdob, faltb
      utcb = float(ih*3600 + im*60) + sec
      if (utcb.ge.54000.) utcb = utcb - 86400.
      if (utcb.le.utca) call abendm('PNAV data UTC illegal sequence')
      goto 23
   24 f = (utcn-utca) / (utcb-utca)
      fidon = fidoa + (fidob-fidoa)*f
      fkdon = fkdoa + (fkdob-fkdoa)*f
      faltn = falta + (faltb-falta)*f
   14 m = m + 1
      if (m.gt.n) goto 11
      f = float(m-1) / float(n-1)
      utc = utc1 + (utcn-utc1)*f
      wtm = utc + 32400.
      nhm = int(wtm/60.)
      nh = nhm / 60
      nm = nhm - nh*60
      wtm = wtm + float(nm*40 + nh*6400)
   15 if (wtm.lt.tms) goto 14
      if (kst.eq.0) then
        tmsr = wtm
        kst = 1
        if (wtm.gt.tme) call abendm('No data Extracted')
      endif
      if (wtm.gt.tme) goto 16
      kintv = kintv - 1
      if (kintv.ge.1) goto 14
      fido = fido1 + (fidon-fido1)*f
      fkdo = fkdo1 + (fkdon-fkdo1)*f
      falt = falt1 + (faltn-falt1)*f
      call igrfc(fido, fkdo, falt, tigrf)
      tres = tmag(m) - tigrf
      write(2,1000) fid(m), iymd, wtm, 3, fido, fkdo, falt,
     *              tmag(m),tres, fx(m),fy(m),fz(m), utc
      wtmr = wtm
      kintv = intv
      goto 14
   16 write(str,'(a17,f9.2,1x,f9.2,a10,i4)')
     *     '    Extracted :  ', tmsr, wtmr, '    intv =', intv
      write(6,'(a)') str(1:50)
      write(9,'(a)') str(1:50)
      kam = 1
      goto 30
c
   19 call abendm('EOF reached on AM input')
   29 call abendm('EOF reached on PNAV input')
c
    3 close(2)
      close(10)
      close(1)
      close(20)
      endfile(9)
      rewind(9)
c
      if(flog(ldr:ldr).ne.' ') then
        open(99,file=flog,access='append')
        write(99,'(//a)')
     *        'XLDPN : Extract Line data  from alog2asc output'
        write(99,'(8x,a)') '   with PostNav Pos. data interporation'
        write(99,'(8x,a)') '       and calculation of IGRF-Residual'
        write(99,'(2a)') ' Lines.Inf. file : ', finf(1:lrtrim(finf))
        write(99,'(2a)') ' Output filename : ', fout(1:lrtrim(fout))
        write(99,'(a,f7.2)') '   Year of survey  : ', year
        write(99,'(a,i2)')   '   IGRF generation : ', ngn
        write(99,'(a,i3)')   '  Data thin-out ratio :', intv
   92   read(9,'(a)',end=93) line
        write(99,'(a)') line(1:lrtrim(line))
        goto 92
   93   call strdtm(sdtm)
        write(99,'(a,a,a)') '=============== ', sdtm, '==============='
        close(99)
      endif
      close(9)
      stop
c
 1000 format(f7.1, 1x,i8, 1x,f9.2, 1x,i2, 1x,f9.5, 1x,f10.5, 1x,f6.1,
     *       2(1x,f8.2), 3(1x,f7.3), 1x,f9.2)
      end
