interp Subroutine

private subroutine interp(me, buf, t, ncf, ncm, na, ifl, pv)

this subroutine differentiates and interpolates a set of chebyshev coefficients to give position and velocity.

Type Bound

jpl_ephemeris

Arguments

Type IntentOptional Attributes Name
class(jpl_ephemeris), intent(inout) :: me
real(kind=wp), dimension(ncf,ncm,*) :: buf

1st location of array of d.p. chebyshev coefficients of position
real(kind=wp), intent(in), dimension(2) :: t
integer, intent(in) :: ncf

of coefficients per component
integer, intent(in) :: ncm

of components per set of coefficients
integer, intent(in) :: na

of sets of coefficients in full array

(i.e., # of sub-intervals in full interval)
integer, intent(in) :: ifl

integer flag = 1 for positions only = 2 for pos and vel
real(kind=wp), dimension(ncm,*) :: pv

interpolated quantities requested. dimension expected is pv(ncm,ifl), dp.

Called by

proc~~interp~~CalledByGraph proc~interp jpl_ephemeris_module::jpl_ephemeris%interp proc~state jpl_ephemeris_module::jpl_ephemeris%state proc~state->proc~interp proc~get_state jpl_ephemeris_module::jpl_ephemeris%get_state proc~get_state->proc~state proc~ephemeris_test jpl_ephemeris_module::ephemeris_test proc~ephemeris_test->proc~get_state proc~get_rv_from_jpl_ephemeris jpl_ephemeris_module::jpl_ephemeris%get_rv_from_jpl_ephemeris proc~get_rv_from_jpl_ephemeris->proc~get_state

Source Code

    subroutine interp(me,buf,t,ncf,ncm,na,ifl,pv)

    implicit none

    class(jpl_ephemeris),intent(inout) :: me
    integer,intent(in)                 :: ncf  !! # of coefficients per component
    integer,intent(in)                 :: ncm  !! # of components per set of coefficients
    real(wp),dimension(ncf,ncm,*)      :: buf  !! 1st location of array of d.p. chebyshev coefficients of position
    real(wp),dimension(2),intent(in)   :: t
    integer,intent(in)                 :: na   !! # of sets of coefficients in full array
                                               !! (i.e., # of sub-intervals in full interval)
    integer,intent(in)                 :: ifl  !! integer flag
                                               !! = 1 for positions only
                                               !! = 2 for pos and vel
    real(wp),dimension(ncm,*)          :: pv   !! interpolated quantities requested.  dimension
                                               !! expected is pv(ncm,ifl), dp.

    real(wp) :: dna,dt1,temp,vfac,tc
    integer :: l,i,j

    ! entry point. get correct sub-interval number for this set
    ! of coefficients and then get normalized chebyshev time
    ! within that subinterval.

    dna  = dble(na)
    dt1  = int(t(1))
    temp = dna*t(1)
    l    = int(temp-dt1)+1

    ! tc is the normalized chebyshev time (-1 <= tc <= 1)

    tc=2.0_wp*(mod(temp,1.0_wp)+dt1)-1.0_wp

    ! check to see whether chebyshev time has changed,
    ! and compute new polynomial values if it has.
    ! (the element pc(2) is the value of t1(tc) and hence
    ! contains the value of tc on the previous call.)

    if (tc/=me%pc(2)) then
        me%np    = 2
        me%nv    = 3
        me%pc(2) = tc
        me%twot  = tc+tc
    endif

    ! be sure that at least 'ncf' polynomials have been evaluated
    ! and are stored in the array 'pc'.

    if (me%np<ncf) then
        do i=me%np+1,ncf
            me%pc(i)=me%twot*me%pc(i-1)-me%pc(i-2)
        end do
        me%np=ncf
    endif

    ! interpolate to get position for each component

    do i=1,ncm
        pv(i,1)=0.0_wp
        do j=ncf,1,-1
            pv(i,1)=pv(i,1)+me%pc(j)*buf(j,i,l)
        end do
    end do
    if (ifl<=1) return

    ! if velocity interpolation is wanted, be sure enough
    ! derivative polynomials have been generated and stored.

    vfac=(dna+dna)/t(2)
    me%vc(3)=me%twot+me%twot
    if (me%nv<ncf) then
        do i=me%nv+1,ncf
            me%vc(i)=me%twot*me%vc(i-1)+me%pc(i-1)+me%pc(i-1)-me%vc(i-2)
        end do
        me%nv=ncf
    endif

    ! interpolate to get velocity for each component

    do i=1,ncm
        pv(i,2)=0.0_wp
        do j=ncf,2,-1
            pv(i,2)=pv(i,2)+me%vc(j)*buf(j,i,l)
        end do
        pv(i,2)=pv(i,2)*vfac
    end do

    end subroutine interp