c--------( DESPIKE : suppress magnetic spike in line data )--------
      parameter (ph=3.)
      dimension fid(10), fido(10), fkdo(10), falt(10), tres(10)
      dimension fx(10), fy(10), fz(10), utc(10)
      dimension iymd(10), km(10), kf(10)
      real*8 wtm(10), tmag(10), ar, cr, d
      character*80 wdr, flog, fnam, fout
      character buf*150, sdtm*20
      call premsg('DESPIKE : suppress magnetic spike in line data')
      call opnpin()
      ldr = lwkdir(50,wdr) + 1
      flog = wdr
      fnam = wdr
      fout = wdr
      call gparma('Enter LOG filename ==> ', 81-ldr, flog(ldr:80))
      call gparma('Enter Input data filename   ==> ',
     *            81-ldr, fnam(ldr:80))
      open(1,file=fnam,status='old')
      call gparma('Enter Output data filename  ==> ',
     *            81-ldr, fout(ldr:80))
      call clspin()
      open(2,file=fout,status='new')
      write(6,'(a)') '----------------------------------------'
      write(6,'(a)') 'DESPIKE : suppress magnetic spike in line data'
      write(6,'(a,a)') ' Input filename   : ', fnam(1:lrtrim(fnam))
      write(6,'(a,a)') ' Output filename  : ', fout(1:lrtrim(fout))
c
    1 do 100 i=1,9
        read(1,'(a)',end=10) buf
        if (buf(1:1).eq.'&'.or.buf(1:1).eq.'%') goto 2
        read(buf,*) fid(i),iymd(i),wtm(i),km(i),fido(i),fkdo(i),
     *            falt(i),tmag(i),tres(i),fx(i),fy(i),fz(i),utc(i)
        kf(i) = 0
  100 continue
      goto 3
    2 n = i-1
      do 200 i=1,n
        write(2,1000) fid(i),iymd(i),wtm(i),km(i),fido(i),fkdo(i),
     *            falt(i),tmag(i),tres(i),fx(i),fy(i),fz(i),utc(i)
  200 continue
      write(2,'(a)') buf(1:lrtrim(buf))
      goto 1
c
    3 j = 0
      dm = ph
      call estim(tmag, kf, ar, cr)
      do 300 i=1,9
        if (kf(i).ne.0) goto 300
        d = ar*float(i-5) + cr - tmag(i)
        if (d.le.dm .and. d.ge.(-dm)) goto 300
        j = i
        dm = abs(d)
  300 continue
      if (j.eq.0) goto 4
      kf(j) = 1
      goto 3
    4 do 400 i=5,1,-1
        if (kf(i).eq.0) goto 400
        d = ar*float(i-5) + cr - tmag(i)
        tres(i) = tres(i) + d
        tmag(i) = tmag(i) + d
        kf(i) = 0
        call estim(tmag, kf, ar, cr)
  400 continue
      do 500 i=1,5
        write(2,1000) fid(i),iymd(i),wtm(i),km(i),fido(i),fkdo(i),
     *            falt(i),tmag(i),tres(i),fx(i),fy(i),fz(i),utc(i)
  500 continue
    5 read(1,'(a)',end=20) buf
      if (buf(1:1).eq.'&'.or.buf(1:1).eq.'%') goto 6
      read(buf,*) fid(10),iymd(10),wtm(10),km(10),fido(10),fkdo(10),
     *       falt(10),tmag(10),tres(10),fx(10),fy(10),fz(10),utc(10)
      kf(10) = 0
      do 600 i=1,9
        fid(i) = fid(i+1)
        iymd(i) = iymd(i+1)
        wtm(i) = wtm(i+1)
        km(i) = km(i+1)
        fido(i) = fido(i+1)
        fkdo(i) = fkdo(i+1)
        falt(i) = falt(i+1)
        tmag(i) = tmag(i+1)
        tres(i) = tres(i+1)
        fx(i) = fx(i+1)
        fy(i) = fy(i+1)
        fz(i) = fz(i+1)
        utc(i) = utc(i+1)
        kf(i) = kf(i+1)
  600 continue
      call estim(tmag, kf, ar, cr)
      d = cr - tmag(5)
      if (kf(5).ne.0 .or. d.gt.ph .or. d.lt.(-ph)) then
        tres(5) = tres(5) + d
        tmag(5) = tmag(5) + d
        kf(5) = 0
      endif
      d = ar*4. + cr - tmag(9)
      if (d.gt.ph .or. d.lt.(-ph)) kf(9) = 1
      write(2,1000) fid(5),iymd(5),wtm(5),km(5),fido(5),fkdo(5),
     *          falt(5),tmag(5),tres(5),fx(5),fy(5),fz(5),utc(5)
      goto 5
    6 do 700 i=6,9
        if (kf(i).ne.0) then
          d = ar*float(i-5) + cr - tmag(i)
          tres(i) = tres(i) + d
          tmag(i) = tmag(i) + d
          kf(i) = 0
          call estim(tmag, kf, ar, cr)
        endif
        write(2,1000) fid(i),iymd(i),wtm(i),km(i),fido(i),fkdo(i),
     *            falt(i),tmag(i),tres(i),fx(i),fy(i),fz(i),utc(i)
  700 continue
      write(2,'(a)') buf(1:lrtrim(buf))
      goto 1
c
   10 n = i-1
      do 800 i=1,n
        write(2,1000) fid(i),iymd(i),wtm(i),km(i),fido(i),fkdo(i),
     *            falt(i),tmag(i),tres(i),fx(i),fy(i),fz(i),utc(i)
  800 continue
      goto 8
   20 do 900 i=6,9
        if (kf(i).ne.0) then
          d = ar*float(i-5) + cr - tmag(i)
          tres(i) = tres(i) + d
          tmag(i) = tmag(i) + d
          kf(i) = 0
          call estim(tmag, kf, ar, cr)
        endif
        write(2,1000) fid(i),iymd(i),wtm(i),km(i),fido(i),fkdo(i),
     *            falt(i),tmag(i),tres(i),fx(i),fy(i),fz(i),utc(i)
  900 continue
    8 close(2)
c
      if(flog(ldr:ldr).ne.' ') then
        open(99,file=flog,access='append')
        write(99,'(//a)')
     *       'DESPIKE : suppress magnetic spike in line data'
        write(99,'(a,a)') ' Input filename   : ', fnam(1:lrtrim(fnam))
        write(99,'(a,a)') ' Output filename  : ', fout(1:lrtrim(fout))
        call strdtm(sdtm)
        write(99,'(a,a,a)') '=============== ', sdtm, '==============='
        close(99)
      endif
      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, 4(1x,f8.2))
      end
c---
      subroutine estim(tmag, kf, ar, cr)
      real*8 tmag(9), ar, cr
      dimension kf(9)
      sxx = 0.
      sxy = 0.
      sx = 0.
      sy = 0.
      n = 0
      do 80 i=1,9
        if (kf(i).ne.0) goto 80
        x = float(i-5)
        y = tmag(i) - tmag(5)
        sxx = sxx + x*x
        sxy = sxy + x*y
        sx = sx + x
        sy = sy + y
        n = n + 1
   80 continue
      denom = sxx*float(n) - sx*sx
      ar = (sxy*float(n) - sx*sy) / denom
      cr = (sxx*sy - sxy*sx) / denom + tmag(5)
      return
      end
