subr. evtrp - j.f.chandler - 1979 nov
derived from subr. sctrp - j.f.chandler, 1977 may
this is a control routine for interpolating a body from its
integration (or n-body or s-body) tape, using everett method
parameters:

jd      integer value of the coord. time of the desired pt
fract   fractional part of the coordinate time of the point
(ignored if jd=0 and nvel=-2)
nvel =0 positions determined
1 velocities determined
-1 velocities determined from position y-vectors
-2 accelerations determined from available y-vectors
lcntl=0 use results of last call if possible
1 do not use results (must recalculate y vectors)
-1 same as 0, except that nb1 is forced to 2
(this option allows for a subsequent call for the same
object for the same observation falling in the
next-earlier tabular interval)
lcntl should, in principle, be set to 1 for the first call,
but this isn't really necessary because ntb1s is initialized
for each series
icall   indicates which is the calling routine
1 - emtrp,  2 - mntrp,  3 - pltrp,  4 - sbtrp,
5 - sctrp,  6 -(prtrp), 7 -(ertrp), 8 - sotrp,
9 - ertsnt, 10 - sotrp for mercury using n-body,
11-19 - sotrp for partials
yv      array for interpolation y-vectors
nb11    pointer to left-hand tabular point in yv
(for objects that may be evaluated at two different
times in the same observation, the later time must
be evaluated first with lcntl=-1, and then the
earlier time with lcntl=0, so that nb11 will be
left over pointing to the correct tabular interval
for the earlier time)
NOTE: this procedure fails if the ephemeris file is
backwards in time and the two needed epochs of
evaluation straddle the epoch of any record
x       output coordinate array
body    coordinates from integration tape
ndim    dimension of 'body' - max. # of partials + 1
jdtb    tabular julian dates
frtb    tabular fractions of a day
p       p-q array to be used

the storage of tabular points in 'body' is in the same time order as
the integration (indicated by 'idir') and is assumed to contain at
least 5 points per record for the three records tagged by jdtb(1-3),
frtb(1-3).  thus, in the case of 5/record the points would be
numbered in order of increasing time as
x              x              x
1, 2, 3, 4, 5, 6, 7, 8, 9,10,11,12,13,14,15  (idir=1)
or
x              x              x
15,14,13,12,11,10, 9, 8, 7, 6, 5, 4, 3, 2, 1  (idir=-1)

where the x's mark the points corresponding to the time tags.  the
interpolation must have at least 5 tabular points bracketing the time
of interest on each side, and the appropriate subroutine is called
if necessary to read the correct records into storage.
bookkeeping is done via common 'trpcom', which is initialized
by subroutine 'trpnit' at the start of each observation series.
bint    directed tabular interval
bintx   absolute value of tabular interval
nb1     same as nb11 above
idirb   (or idir) time direction of integration
nvels=0 y-vectors based on positions
1 y-vectors based on positions and velocities
>1 y-vectors not set up
ntb1s   saved value of ntab1, or 9999 if y-vectors not set up
nbtrp   pointer to left bracketing tabular point (relative
to 2nd or 3rd record, depending on idir)

npto2= half the number of points needed by the interpolator (5 or 7)
nprec= number of tabular points/record (5 or 8 or 10)
nov5 = nprec - npto2
lim1,lim2 = lower and upper limits for left-hand tabular
point for interpolation (relative to time tag of
2nd record if forward or 3rd if backward)

some y-vectors may be retained, shift storage as needed
note: nmo represents time offset; y-vectors
are to be shifted by -nmo.  if that would shift all
present storage off the end, then start fresh.