orb_interp.f

Go to the documentation of this file.

        subroutine orb_interp(orbit,nlines,timref,time,posi,veli,orbfl)
c
c  prolog(orbit,timref,time,posi,veli)
c
c  purpose: interpolate orbit position and velocity vectors to scan line times
c  
c
c  calling arguments:
c
c  name         Type    i/o     description
c  --------     ----    ---     -----------
c  orbit        struct   i      input structure of filtered orbit data
c  nlines       i*4      i      number of scan line times
c  timref       r*8      i      Size 3 reference time at start line 
c                                of data:  year, day, sec
c  time         r*8      i      array of time in seconds relative to 
c                                timref for every scan line
c  posi         r*4      o      Size 3 by nlines interpolated position
c  veli         r*4      o      Size 3 by nlines interpolated velocity
c  orbfl        i*4      o      Size nlines orbit vector flags(0=good)
c
c
c  by: frederick s. patt, gsc, august 10, 1993
c
c  notes:  method uses cubic polynomial to match positions and velocities
c   at input data points.
c
c  modification history:
c
c  added check for array index exceeded number of scan lines in check for
c  scan lines preceding first scan line.  f.s. patt, saic gsc, may 23, 2000.
c
      implicit none
#include "nav_cnst.fin"
#include "orbit_s.fin"
      type(orbit_struct) :: orbit
c
      integer*4 ind, i, j, i1, nlines, iy, id, jd
      real*8 timref(3), time(nlines)
      real*8 a0(3),a1(3),a2(3),a3(3), t, dt, x, x2, x3, ddifs
      real*4 posi(3,nlines), veli(3,nlines)
      integer*4 orbfl(nlines)
      logical first       
 
      ind = 1
      i1 = 1
      first = .true.
c  make sure that orbit vector precedes first scan line
      iy = timref(1)
      id = timref(2)
      ddifs = (jd(iy,1,id) - jd(orbit%iyr,1,orbit%iday))*864.d2
      t = ddifs + timref(3)+time(i1)
      dowhile((t.lt.orbit%torb(ind)).and.(i1.le.nlines))
        do j=1,3
          posi(j,i1) = 0.0
          veli(j,i1) = 0.0
        end do
        orbfl(i1) = 1
        print *, 'Scan line times before available orbit data'
        i1 = i1 + 1
        t = ddifs + timref(3)+time(i1)
      end do
      
      do i=i1,nlines
        t = ddifs + timref(3)+time(i)
        if ((t.gt.orbit%torb(ind+1)).or.first) then
          dowhile(t.gt.orbit%torb(ind+1))
            ind = ind + 1
            if (ind.gt.orbit%nvec) then
              i1 = 1
              go to 990
            end if
          end do
          first = .false.
c  set up cubic interpolation
          dt = orbit%torb(ind+1) - orbit%torb(ind)
          do j=1,3
            a0(j) = orbit%pos(j,ind)
            a1(j) = orbit%vel(j,ind)*dt
            a2(j) = 3.d0*orbit%pos(j,ind+1) - 3.d0*orbit%pos(j,ind) 
     *              - 2.d0*orbit%vel(j,ind)*dt - orbit%vel(j,ind+1)*dt
            a3(j) = 2.d0*orbit%pos(j,ind) - 2.d0*orbit%pos(j,ind+1) 
     *              + orbit%vel(j,ind)*dt + orbit%vel(j,ind+1)*dt
          end do
        end if
 
c  interpolate orbit position and velocity components to scan line time
        x = (t - orbit%torb(ind))/dt
        x2 = x*x
        x3 = x2*x
        do j=1,3
          posi(j,i) = a0(j) + a1(j)*x + a2(j)*x2 + a3(j)*x3
          veli(j,i) = (a1(j) + 2.*a2(j)*x + 3.*a3(j)*x2)/dt
        end do
        orbfl(i) = 0
 
      end do
      return
 
990   print *, 'Scan line times after available orbit data'
      do i=i1,nlines
        do j=1,3
          posi(j,i) = 0.0
          veli(j,i) = 0.0
        end do
      orbfl(i) = 1
      end do
      return
 
      end