subroutine pair_difference( NSets, Epoch, Incl, Node, E, Omega,
     &   M, N, bStar, j1, j2, dT, dR, dV, dRdT, dVdT )
!
!     Compute the position and velocity errors propagating from tle(j1)
!     to tle(j2)
!
!     B. Knapp, 1999-08-19
!
      implicit none
!
!     Input
      integer NSets, j1, j2
      real*8 Epoch(NSets), Incl(NSets), Node(NSets), E(NSets),
     &   Omega(NSets), M(NSets), N(NSets), bStar(NSets)
!
!     Output
      real*8 dT, dR, dV, dRdT, dVdT
!
!     Local
      integer i, IFlag, Status
      real*8 PVj1(6), PVj2(6)
!
!     Subroutines
      external SGP4
!
      dT = Epoch(j2)-Epoch(j1)
      IFlag = 1
      call SGP4( IFlag, dT*1440.d0, M(j1), Node(j1), Omega(j1),
     &   E(j1), Incl(j1), N(j1), bStar(j1), PVj1, Status )
      IFlag = 1
      call SGP4( IFlag, 0.0d0,      M(j2), Node(j2), Omega(j2),
     &   E(j2), Incl(j2), N(j2), bStar(j2), PVj2, Status )
      dR = 0.0d0
      dV = 0.0d0
      do i=1,3
         dR = dR + (PVj1(i)-PVj2(i))**2
         dV = dV + (PVj1(i+3)-PVj2(i+3))**2
      enddo
      dR = sqrt(dR)
      dV = sqrt(dV)
      dRdT = dR/dT
      dVdT = dV/dT
!
      return
      end