/*Translated by FOR_C, v3.4.2 (-), on 07/09/115 at 08:31:57 */
/*FOR_C Options SET: ftn=u io=c no=p op=aimnv s=dbov str=l x=f - prototypes */
#include <math.h>
#include "fcrt.h"
#include "dsval.h"
#include <stdlib.h>
		/* PARAMETER translations */
#define	KMAX	20
#define	ZERO	0.0e0
		/* end of PARAMETER translations */
 
double /*FUNCTION*/ dsval(
long k,
long nc,
double t[],
double bcoef[],
double x,
long ideriv)
{
	long int iderp1, ihi, ilo, imk, ip1, j, jj, km1, kmider, kmj,
	 mode;
	double aj[KMAX], dm[KMAX], dp[KMAX], dsval_v, fkmj;
	static long lefti = 1;
		/* OFFSET Vectors w/subscript range: 1 to dimension */
	double *const Aj = &aj[0] - 1;
	double *const Bcoef = &bcoef[0] - 1;
	double *const Dm = &dm[0] - 1;
	double *const Dp = &dp[0] - 1;
	double *const T = &t[0] - 1;
		/* end of OFFSET VECTORS */
 
	/* Copyright (c) 1996 California Institute of Technology, Pasadena, CA.
	 * ALL RIGHTS RESERVED.
	 * Based on Government Sponsored Research NAS7-03001.
	 *>> 1994-10-20 DSVAL Krogh  Changes to use M77CON
	 *>> 1994-09-26 DSVAL CLL Moved "DSVAL = ZERO" to be 1st executable stmt.
	 *>> 1992-11-12 C. L. Lawson, JPL  Saving LEFTI.
	 *>> 1992-10-27 C. L. Lawson, JPL
	 *>> 1988-03-15 C. L. Lawson, JPL
	 *     Calculates the value at X of the derivative of order IDERIV
	 *     of the spline function represented in B-spline form by
	 *     K, NC, T(), and BCOEF().
	 *
	 *     Based on subroutine BVALUE on pp. 144-145 of A PRACTICAL GUIDE TO
	 *     SPLINES by Carl De Boor, Springer-Verlag, 1978.
	 *     Current version by C. L. Lawson, JPL, March 1988.
	 *     ------------------------------------------------------------------
	 *  K     [in]  Order of the spline functions.  Note that the polynomial
	 *        degree of the segments of the spline is one less than the
	 *        order.  Example:  Cubic splines have order K = 4.
	 *  NC     [in]  Number of B-spline coefficients.  Require NC .ge. 1.
	 *  T()   [in]  Knot sequence, indexed from 1 to NT, with
	 *        NT = NC + K.  Knot values must be nonincreasing.
	 *        Repetition of values is permitted and has the effect of
	 *        decreasing the order of contimuity at the repeated knot.
	 *        Proper function representation by splines of order K is
	 *        supported on the interval from T(K) to T(NC+1).
	 *        Extrapolation can be done outside this interval.
	 *  BCOEF()  [in]  Coefficients of B-spline basis functions, indexed from
	 *        1 to NC.
	 *  X     [in]  Abcissa at which the spline function or one of its
	 *        derivatives is to be evaluated.  The evaluation will use one
	 *        of the polynomial pieces of the spline as follows:
	 *        If X .lt. T(K+1) use the piece associated with [T(K),T(K+1))
	 *        If T(L) .le. X .lt. T(L+1) for some L in [K+1, NC-1] use the
	 *        piece associated with [T(L),T(L+1)).
	 *        If T(NC) .le. X use the piece associated with
	 *        [T(NC),T(NC+1)].
	 *  IDERIV   [in]  Order of derivative to be evaluated.  IDERIV = 0
	 *        selects evaluation of the spline function.  Require
	 *        IDERIV .ge. 0.  All derivatives of orders .ge. K will be zero.
	 *  DSVAL [out] Returned value of the spline function or its requested
	 *        derivative.
	 *     ------------------------------------------------------------------
	 *--D replaces "?": ?SVAL, ?SFIND
	 *     ------------------------------------------------------------------ */
	/*     ------------------------------------------------------------------ */
	dsval_v = ZERO;
	if (k > KMAX)
	{
		ierm1( "DSVAL", 1, 2, "Require KORDER .le. KMAX.", "KORDER"
		 , k, ',' );
		ierv1( "KMAX", KMAX, '.' );
		return( dsval_v );
	}
	kmider = k - ideriv;
	if (kmider <= 0)
		goto L_99;
 
	/*     If T(K) .le. X .lt. NC+1, DSFIND will return LEFTI such that
	 *     T(LEFTI) .le. X .lt. T(LEFTI+1).
	 *     Otherwise if X .lt. T(K), sets LEFTI := K
	 *            or if X .ge. T(NC+1), sets LEFTI := NC.
	 * */
	km1 = k - 1;
	dsfind( t, k, nc + 1, x, &lefti, &mode );
 
	/*        MODE =  -1 if X < T(K)
	 *                  0 if T(1) .le. X .le. T(NC+1)
	 *                 +1 if T(NC+1) .lt. X
	 *
	 *               *** DIFFERENCE THE COEFFICIENTS *IDERIV* TIMES . */
	imk = lefti - k;
	for (j = 1; j <= k; j++)
	{
		Aj[j] = Bcoef[imk + j];
	}
	for (j = 1; j <= ideriv; j++)
	{
		kmj = k - j;
		fkmj = (double)( kmj );
		for (jj = 1; jj <= kmj; jj++)
		{
			ihi = lefti + jj;
			Aj[jj] = (Aj[jj + 1] - Aj[jj])/(T[ihi] - T[ihi - kmj])*
			 fkmj;
		}
	}
 
	/*  *** COMPUTE VALUE AT *X* IN (T(I),T(I+1)) OF IDERIV-TH DERIVATIVE,
	 *      GIVEN ITS RELEVANT B-SPLINE COEFF. IN AJ(1),...,AJ(K-IDERIV).
	 * */
	if (ideriv != km1)
	{
		ip1 = lefti + 1;
		for (j = 1; j <= kmider; j++)
		{
			Dp[j] = T[lefti + j] - x;
			Dm[j] = x - T[ip1 - j];
		}
		iderp1 = ideriv + 1;
		for (j = iderp1; j <= km1; j++)
		{
			kmj = k - j;
			ilo = kmj;
			for (jj = 1; jj <= kmj; jj++)
			{
				Aj[jj] = (Aj[jj + 1]*Dm[ilo] + Aj[jj]*Dp[jj])/(Dm[ilo] +
				 Dp[jj]);
				ilo -= 1;
			}
		}
	}
	dsval_v = Aj[1];
 
L_99:
	;
	return( dsval_v );
} /* end of function */