polynom <- function (x, norder=1, ctr=midrange, nderiv) { # This computes values of the polynomials, # P_l(x) = (x-ctr)^l, l=0,...,NORDER-1 # or their derivatives. # The degree of the highest order polynomial is one less than NORDER. # The default is the constant function. # Arguments are as follows: # X ... array of values at which the polynomials are to # evaluated # NORDER ... the polynomial basis object # NDERIV ... highest order derivative. 0 means only function values # are returned. # CTR ... a constant to be subtracted from a value prior to taking # its power # Return is a matrix with length(X) rows and NORDER columns containing # the values of the polynomials # last modified 17 January 2006 x <- as.vector(x) n <- length(x) ndegree <- norder - 1 if (nderiv > ndegree) stop('NDERIV exceeds highest degree of polynomial.') rangex <- range(x) midrange <- mean(rangex) lfac <- 1 if (nderiv > 1) for (l in 2:nderiv) lfac <- lfac*l polyval <- matrix(0,n,norder) polyval[,nderiv+1] <- lfac if (norder > nderiv+1) for (l in (nderiv+2):norder) polyval[,l] <- polyval[,l-1]*(x-ctr)*(l-1)/(l-nderiv-1) return (polyval) }