## Name

CHEBGRID

## Author

Craig B. Markwardt, NASA/GSFC Code 662, Greenbelt, MD 20770

craigm@lheamail.gsfc.nasa.gov

UPDATED VERSIONs can be found on my WEB PAGE:

http://cow.physics.wisc.edu/~craigm/idl/idl.html

## Purpose

Estimate Chebyshev polynomial coefficients of a function on a grid

## Major Topics

Curve and Surface Fitting

## Calling Sequence

p = CHEBGRID(T, X, [ DXDT, NPOINTS=, NPOLY=, NGRANULE= , $

RMS=, DRMS=, RESIDUALS=, DRESIDUALS= , $

XMATRIX=, DXMATRIX=, RESET=,

DERIV_WEIGHT= ] )

## Description

CHEBGRID estimates the coefficients for a finite sum of Chebyshev

polynomials approximating a continuous tabulated function over an

interval. The function (and optionally its derivative) must be

tabulated on a regularly sampled grid. The implementation of this

function is taken from a method described by X. X. Newhall, used

in estimating coefficients for ephemerides in the solar system.

The tabulated function is assumed to be continuous over the entire

interval. A Chebyshev series is fitted to the function over small

segments, called granules. The size of each granule, the number

of points in each granule, and the number of Chebyshev polynomials

are all configurable.

Users may specify either the function alone, or the function and

its first derivative. By also giving the tabulated derivative, a

more accurate Chebyshev polynomial can be developed. Aside from

the constraints mentioned in the next paragraph, the polynomial

that is returned is the best-fit polynomial in a least-squares

sense.

Here is a definition of terms:

GRANULE - a single continuous fitted segment. The length of the

granule, NGRANULE, is specified in units of the tabulated

grid size. Because of the continuity requirements developed

below, granules will always overlap at their endpoints.

Thus, then length of a granule should be a factor of

N_ELEMENTS(X)-1. For simple functions over short intervals,

the granule size can be equal to N_ELEMENTS(X)-1

NUMBER OF POINTS the number of points, NPOINTS, within a

granule to be fitted to the polynomial, not necessarily

equal to the granule size. The greater the number of

points, the more computation time and storage is required.

This number *must* be a factor of NGRANULE. Typically

NPOINTS is a number between 8 and 12. Because of the

single-point overlap between granules (see below), the

actual number of points per fit is NPOINTS+1.

NUMBER OF POLYNOMIALS the number of Chebyshev polynomial terms,

NPOLYNOMIAL, to be fitted per granule. The greater the

number of polynomial terms, the more computation time and

storage is required, but also the greater the approximating

precision of the fit.

The particular set of Chebyshev polynomial coefficients developed

by this function have some special properties. If both the

function and its derivative are specified, then the value and

derivative of the interpolating polynomial at the granule

endpoints will be exactly equal to the tabulated endpoint values.

This feature allows many approximations to be strung together

piecewise, and the function value and first derivative will be

continuous across granule boundaries.

If only the function value is specified, then only the function

value will be continuous at the granule endpoints, and not the

derivative.

An extensive set of statistics are computed to assess the quality

of the Chebyshev polynomial fit. The keywords RESIDUALS and

DRESIDUALS return the residuals of the fit after subtracting the

interpolation. The RMS and DRMS keywords return the root mean

squared deviations between data and model.

If the user does not know how many granules, points, or polynomial

coefficients to use, then he or she should try several

combinations and see which minimizes the r.m.s. value with the

fewest number of coefficients.

If the XMATRIX and DXMATRIX keywords are passed, then CHEBGRID

attempts to avoid recomputing several of the matrices it uses in

estimating the coefficients. If multiple calls to CHEBGRID are to

be made, some compution time savings can be made. In the first

call CHEBGRID the required matrices are computed and returned. In

subsequent calls, CHEBGRID detects the XMATRIX and DXMATRIX

keyword values and uses those values if it can.

The user can also estimate their own coefficients. The matrices

returned are (NPOINTS+1)x(NPOLYNOMIAL). The coefficients from a

NPOINTS+1 tabulation, X, are found by:

PCHEB = XMATRIX ## X + DXMATRIX ## DXDT

if derivative information is known, or

PCHEB = XMATRIX ## X

if no derivative information is known. [ Note: the matrices are

different, depending on whether derivative information is known or

not. ]

## Inputs

T - array of regularly sampled *independent* variables. The number

of elements in T should be a multiple of NGRANULE, plus one.

X - array of regularly sampled *dependent* variables. The number

of elements in X should be equal to the number of elements in

T.

DXDT - optionally, a tabulated array of first derivatives of X

with respect to T, at the same grid points.

## Keyword Parameters

NGRANULE - size of a "granule", in grid intervals. NGRANULE must

be at least 2, and a factor of N_ELEMENTS(T)-1.

Default: 8

NPOINTS - number of points per granule that are fitted. NPOINTS

must be at least 2, and a factor of NGRANULE.

Default: NGRANULE

NPOLYNOMIAL - number of Chebyshev polynomial terms per fit.

NPOLYNOMIAL must be at least 2 and less than

2*(NPOINTS+1), when derivative information is

specified; or less than NPOINTS+1, when no

derivative information is specified.

Default: 7

RESIDUALS - upon return, an array of size N_ELEMENTS(T), with

residuals of the tabulated function minus the

interpolated function.

DRESIDUALS - same as RESIDUALS, but for the function's first

derivative.

RMS - upon return, the root mean square of the function value

residuals.

DRMS - same as RMS, but for the function's first derivative.

XMATRIX - upon return, the matrix used to compute Chebyshev

polynomial coefficients from the function value.

Upon input, CHEBGRID determines if XMATRIX will apply to

the data, and if so, XMATRIX is reused rather than

computed. If XMATRIX cannot be reused, then it is

computed afresh, and the new value is returned in the

XMATRIX keyword.

The user should not modify the contents of this array.

DXMATRIX - same as XMATRIX, but for the function's first

derivative.

RESET - if set, force a recomputation of XMATRIX and/or DXMATRIX.

DERIV_WEIGHT - amount of weight to give to function derivative,

relative to the function value.

Default: 0.16d

## Returns

An array of coefficient values. The dimensions of the array are

NPOLYNOMIALxNSEGS, where NSEGS is the number of granules in the

entire interval.

## Example

;; Estimate Chebyshev coefficients for the function SIN(X), on the

;; interval [-1,+1].

xx = dindgen(9)/4d - 1d ;; Regular grid from -1 to 1 (9 points)

yy = sin(xx) ;; Function values, sin(x), ...

dy = cos(xx) ;; ... and derivatives

;; Estimate coefficients using CHEBGRID (single granule of 8 intervals)

p = chebgrid(xx, yy, dy, npoints=8, ngranule=8, npoly=10)

xxx = dindgen(1001)/500 - 1d ;; New grid for testing

res = sin(xxx) - chebeval(xxx, p)

plot, xxx, res

;; Same as example above, except extended range to [-1, +15],

using eight granules.

xx2 = dindgen(65)/4d - 1

yy2 = sin(xx2)

dy2 = cos(xx2)

p = chebgrid(xx2, yy2, dy2, ngranule=8, npoint=8, npoly=10)

help, p

P DOUBLE = Array[10, 8]

;; (i.e., 10 polynomial coefficients over 8 granules)

## References

Abramowitz, M. & Stegun, I., 1965, *Handbook of Mathematical

Functions*, 1965, U.S. Government Printing Office, Washington,

D.C. (Applied Mathematical Series 55)

Newhall, X. X. 1989, Celestial Mechanics, 45, p. 305-310

## Modification History

Written, CM, Feb 2002

Documented, CM, 24 Mar 2002

Corrected documentation, CM, 28 Apr 2002

Typo correction, CM, 10 Oct 2002