aboutsummaryrefslogtreecommitdiff
path: root/math/iminterp/doc
diff options
context:
space:
mode:
Diffstat (limited to 'math/iminterp/doc')
-rw-r--r--math/iminterp/doc/arbpix.hlp57
-rw-r--r--math/iminterp/doc/arider.hlp59
-rw-r--r--math/iminterp/doc/arieval.hlp48
-rw-r--r--math/iminterp/doc/asider.hlp52
-rw-r--r--math/iminterp/doc/asieval.hlp44
-rw-r--r--math/iminterp/doc/asifit.hlp40
-rw-r--r--math/iminterp/doc/asifree.hlp25
-rw-r--r--math/iminterp/doc/asigeti.hlp36
-rw-r--r--math/iminterp/doc/asigetr.hlp36
-rw-r--r--math/iminterp/doc/asigrl.hlp40
-rw-r--r--math/iminterp/doc/asiinit.hlp39
-rw-r--r--math/iminterp/doc/asirestore.hlp36
-rw-r--r--math/iminterp/doc/asisave.hlp39
-rw-r--r--math/iminterp/doc/asisinit.hlp60
-rw-r--r--math/iminterp/doc/asitype.hlp95
-rw-r--r--math/iminterp/doc/asivector.hlp52
-rw-r--r--math/iminterp/doc/im1dinterp.spc525
-rw-r--r--math/iminterp/doc/im2dinterp.spc432
-rw-r--r--math/iminterp/doc/iminterp.hd37
-rw-r--r--math/iminterp/doc/iminterp.hlp234
-rw-r--r--math/iminterp/doc/iminterp.men32
-rw-r--r--math/iminterp/doc/iminterp.spc525
-rw-r--r--math/iminterp/doc/mrider.hlp79
-rw-r--r--math/iminterp/doc/mrieval.hlp57
-rw-r--r--math/iminterp/doc/msider.hlp52
-rw-r--r--math/iminterp/doc/msieval.hlp46
-rw-r--r--math/iminterp/doc/msifit.hlp45
-rw-r--r--math/iminterp/doc/msifree.hlp26
-rw-r--r--math/iminterp/doc/msigeti.hlp35
-rw-r--r--math/iminterp/doc/msigetr.hlp37
-rw-r--r--math/iminterp/doc/msigrid.hlp51
-rw-r--r--math/iminterp/doc/msigrl.hlp43
-rw-r--r--math/iminterp/doc/msiinit.hlp41
-rw-r--r--math/iminterp/doc/msirestore.hlp36
-rw-r--r--math/iminterp/doc/msisave.hlp38
-rw-r--r--math/iminterp/doc/msisinit.hlp61
-rw-r--r--math/iminterp/doc/msisqgrl.hlp38
-rw-r--r--math/iminterp/doc/msitype.hlp95
-rw-r--r--math/iminterp/doc/msivector.hlp54
39 files changed, 3377 insertions, 0 deletions
diff --git a/math/iminterp/doc/arbpix.hlp b/math/iminterp/doc/arbpix.hlp
new file mode 100644
index 00000000..0d7ec9ec
--- /dev/null
+++ b/math/iminterp/doc/arbpix.hlp
@@ -0,0 +1,57 @@
+.help arbpix Dec98 "Image Interpolator Package"
+.ih
+NAME
+arbpix -- replace INDEF valued pixels with interpolated values
+.ih
+SYNOPSIS
+include <math/iminterp.h>
+
+arbpix (datain, dataout, npix, interp_type, boundary_type)
+
+.nf
+ real datain[npix] #I input data
+ real dataout[npix] #O output array, dataout != datain
+ int npix #I number of data points
+ int interp_type #I type of interpolant
+ int boundary_type #I type of boundary condition
+.fi
+.ih
+ARGUMENTS
+.ls datain
+Array of input data containing 0 or more INDEF valued pixels.
+.le
+.ls dataout
+Array of output data with INDEFS replaced by interpolated values.
+The dataout array must be different from the datain array.
+.le
+.ls npix
+Number of data points.
+.le
+.ls interp_type
+Type of interpolant. Options are II_NEAREST, II_LINEAR, II_POLY3, II_POLY5,
+II_SPLINE3, II_SINC / II_LSINC, and II_DRIZZLE. The look-up table sinc
+interpolant is not supported, and defaults to the sinc interpolant.
+The sinc interpolant width is 31 pixels. The drizzle interpolant is not
+supported and defaults to the linear interpolant. The interpolant type
+definitions are stored in the file math/iminterp.h.
+.le
+.ls boundary_type
+Type of boundary extension. The only supported option is II_BOUNDARYEXT.
+Polynomial interpolants of lower order are used if there are not enough
+good pixels to define the requested interpolant. Nearest neighbor boundary
+extension is used if there are not enough good points to define the sinc
+interpolant. The boundary type definitions are stored in the header file
+math/iminterp.h.
+.le
+.ih
+DESCRIPTION
+If there are no good points in datain, ARBPIX returns INDEFS in dataout.
+Points below and above the first and last good point are replaced by the
+first and last good point values respectively.
+.ih
+NOTES
+The sinc function actually evaluates the interpolant by computing the
+average of two interpolations at +-0.05 pixels about the bad pixel since
+the interpolant is undefined exactly at a pixel.
+.ih
+SEE ALSO
diff --git a/math/iminterp/doc/arider.hlp b/math/iminterp/doc/arider.hlp
new file mode 100644
index 00000000..b8631eaa
--- /dev/null
+++ b/math/iminterp/doc/arider.hlp
@@ -0,0 +1,59 @@
+.help arider Dec98 "Image Interpolator Package"
+.ih
+NAME
+arider -- calculate the interpolant derivatives at x
+.ih
+SYNOPSIS
+include <math/iminterp.h>
+
+arider (x, datain, npix, der, nder, interp_type)
+
+.nf
+ real x[2] #I x value, 1 <= x[1-2] <= npts
+ real datain[npix] #I array of data points
+ int npix #I number of data points
+ real der[nder] #O derivatives, der[1] = function value
+ int nder #I number of derivatives, 1 + max order
+ int interp_type #I interpolant type
+.fi
+.ih
+ARGUMENTS
+.ls x
+Single X value, or pair of X values defining a range in the case of the
+drizzle interpolant.
+.le
+.ls datain
+Array of data values.
+.le
+.ls npix
+Number of data points.
+.le
+.ls der
+Array of derivatives. Der[1] contains the function value, der[2] the
+first derivative, and so on.
+.le
+.ls nder
+Number of derivatives. ARIDER checks that the requested number of derivatives
+is sensible. The sinc interpolant returns the function value and the first
+two derivatives. The drizzle interpolant returns the function and the first
+derivative.
+.le
+.ls interp_type
+Interpolant type. The options are II_NEAREST, II_LINEAR, II_POLY3, II_POLY5,
+II_SPLINE3, II_SINC / II_LSINC, and II_DRIZZLE. The look-up table sinc
+is not supported and defaults to sinc. The sinc interpolant width is 31 pixels.
+The drizzle pixel fraction is 1.0. The interpolant type definitions are found
+in the package header file math/iminterp.h.
+.le
+.ih
+DESCRIPTION
+ARIDER permits the evaluation of the interpolant at a few randomly spaced
+points within datain without the storage requirements of the sequential
+version.
+.ih
+NOTES
+Checking for INDEF valued or out of bounds pixels is the responsibility
+of the user.
+.ih
+SEE ALSO
+asider
diff --git a/math/iminterp/doc/arieval.hlp b/math/iminterp/doc/arieval.hlp
new file mode 100644
index 00000000..52c62148
--- /dev/null
+++ b/math/iminterp/doc/arieval.hlp
@@ -0,0 +1,48 @@
+.help arieval Dec98 "Image Interpolator Package"
+.ih
+NAME
+arieval -- evaluate the interpolant at x
+.ih
+SYNOPSIS
+include <math/iminterp.h>
+
+y = arieval (x, datain, npix, interp_type)
+
+.nf
+ real x[2] #I x value, 1 <= x[1-2] <= npix
+ real datain[npix] #I data values
+ int npix #I number of data values
+ int interp_type #I interpolant type
+.fi
+.ih
+ARGUMENTS
+.ls x
+Single X value, or a pair of X values specifying a range in the case
+of the drizzle interpolant.
+.le
+.ls datain
+Array of input data.
+.le
+.ls npix
+Number of data points.
+.le
+.ls interp_type
+Interpolant type. Options are II_NEAREST, II_LINEAR, II_POLY3, II_POLY5,
+II_SPLINE3, II_SINC / II_LSINC, and II_DRIZZLE, for nearest neighbor,
+linear, 3rd and fifth order polynomials, cubic spline, sinc, look-up
+table sinc, and drizzle interpolants respectively. The look-up table sinc
+interpolant is not supported and defaults to the sinc interpolant. The sinc
+width is 31 pixels. The drizzle pixel fraction is 1.0. The interpolant
+type definitions are contained in the package header file math/iminterp.h
+.le
+.ih
+DESCRIPTION
+ARIEVAL allows the evaluation of a few interpolated points without the
+storage required for the sequential interpolant.
+.ih
+NOTES
+Checking for out of bounds and INDEF valued pixels is the responsibility of
+the user.
+.ih
+SEE ALSO
+arider, asieval, asivector
diff --git a/math/iminterp/doc/asider.hlp b/math/iminterp/doc/asider.hlp
new file mode 100644
index 00000000..0c27ffbc
--- /dev/null
+++ b/math/iminterp/doc/asider.hlp
@@ -0,0 +1,52 @@
+.help asider Dec98 "Image Interpolator Package"
+.ih
+NAME
+asider -- evaluate the interpolant derivatives at x
+.ih
+SYNOPSIS
+asider (asi, x, der, nder)
+
+.nf
+ pointer asi #I interpolant descriptor
+ real x[2] #I x value, 1 <= x[1-2] <= npix
+ real der[] #O der[1] = interpolant, der[2] = 1st derivative
+ int nder #I number of derivatives
+.fi
+.ih
+ARGUMENTS
+.ls asi
+Pointer to the sequential interpolant descriptor.
+.le
+.ls x
+Single X value, or pair of X values defining a range in the case of
+the drizzle interpolant.
+.le
+.ls der
+Array containing the derivatives. Der[1] = interpolant at x, der[2] the
+first derivative of the interpolant at x and so on.
+.le
+.ls nder
+Number of derivatives. Nder = 1 + order of the maximum desired derivative.
+ASIDER checks that nder is reasonable. The sinc interpolant returns the
+interpolant value and first two derivatives. The drizzle interpolant returns
+the interpolant value and the first derivative.
+.le
+.ih
+DESCRIPTION
+The polynomial coefficients are evaluated directly from the data points
+for the polynomial interpolants and from the B-spline coefficients
+for the cubic spline interpolant. The derivatives are evaluated from
+the polynomial coefficients using nested multiplication. The sinc
+derivatives are analytic but are defined only for the first two derivatives.
+The drizzle derivative is an approximation defined for the first derivative
+only.
+.ih
+NOTES
+ASIDER checks that the number of derivatives requested is reasonable.
+Checking for out of bounds and INDEF valued pixels is the responsibility
+of the user. ASIINIT or ASISINIT and ASIFIT must be called before ASIDER
+is called.
+.ih
+SEE ALSO
+asieval, asivector, arieval, arider
+.endhelp
diff --git a/math/iminterp/doc/asieval.hlp b/math/iminterp/doc/asieval.hlp
new file mode 100644
index 00000000..20f70abe
--- /dev/null
+++ b/math/iminterp/doc/asieval.hlp
@@ -0,0 +1,44 @@
+.help asieval Dec98 "Image Interpolator Package"
+.ih
+NAME
+asieval -- procedure to evaluate interpolant at x
+.ih
+SYNOPSIS
+y = asieval (asi, x)
+
+.nf
+ pointer asi #I interpolant descriptor
+ real x[2] #I x value, 1 <= x[1-2] <= npts
+.fi
+.ih
+ARGUMENTS
+.ls asi
+Pointer to the sequential interpolant descriptor structure.
+.le
+.ls x
+Single X value, or pair of X values defining a range in the case of the
+drizzle interpolant.
+.le
+.ih
+DESCRIPTION
+The polynomial coefficients are calculated directly from the data points
+for the polynomial interpolants, and from the B-spline coefficients for
+the cubic spline interpolant. The actual calculation is done by adding and
+multiplying terms according to Everett's central difference interpolation
+formula. The boundary extension algorithm is projection.
+
+The sinc interpolant is computed using a range of data points around
+the desired position. Look-up table sinc interpolation is computed
+using the most appropriate entry in a precomputed look-up table.
+The boundary extension algorithm is nearest neighbor.
+
+The drizzle interpolant is computed by summing the data over the user
+supplied X interval.
+.ih
+NOTES
+Checking for out of bounds and INDEF valued pixels is the responsibility of
+the user. ASIINIT or ASISINIT and ASIFIT must be called before using ASIEVAL.
+.ih
+SEE ALSO
+asivector, arieval
+.endhelp
diff --git a/math/iminterp/doc/asifit.hlp b/math/iminterp/doc/asifit.hlp
new file mode 100644
index 00000000..bcd1fdc8
--- /dev/null
+++ b/math/iminterp/doc/asifit.hlp
@@ -0,0 +1,40 @@
+.help asifit Dec98 "Image Interpolator Package"
+.ih
+NAME
+asifit - fit the interpolant to data
+.ih
+SYNOPSIS
+asifit (asi, datain, npix)
+
+.nf
+ pointer asi #I interpolant descriptor
+ real datain[npix] #I input data
+ int npix #I the number of data points
+.fi
+.ih
+ARGUMENTS
+.ls asi
+Pointer to sequential interpolant descriptor structure.
+.le
+.ls datain
+Array of input data.
+.le
+.ls npix
+Number of data points.
+.le
+.ih
+DESCRIPTION
+The datain array is checked for size, memory is allocated for the coefficient
+array, and the end conditions are specified. The interior polynomial, sinc and
+drizzle interpolants are saved as the data points. The polynomial coefficients
+are calculated directly from the data points in the evaluation stage. The
+B-spline coefficients are calculated in ASIFIT as they depend on the entire
+data array.
+.ih
+NOTES
+Checking for INDEF valued and out of bounds pixels is the responsibility
+of the user. ASIINIT or ASISINIT and ASIFIT must be called before using
+ASIEVAL, ASIVECTOR, ASIDER or ASIGRL.
+.ih
+SEE ALSO
+.endhelp
diff --git a/math/iminterp/doc/asifree.hlp b/math/iminterp/doc/asifree.hlp
new file mode 100644
index 00000000..c61f2ce0
--- /dev/null
+++ b/math/iminterp/doc/asifree.hlp
@@ -0,0 +1,25 @@
+.help asifree Dec98 "Image Interpolator Package"
+.ih
+NAME
+asifree - free sequential interpolant descriptor
+.ih
+SYNOPSIS
+asifree (asi)
+
+.nf
+pointer asi #U interpolant descriptor
+.fi
+.ih
+ARGUMENTS
+.ls asi
+Pointer to the sequential interpolant descriptor structure.
+.le
+.ih
+DESCRIPTION
+ASIFREE frees the sequential interpolant descriptor structure allocated by
+ASIINIT or ASISINIT. ASIFREE should be called when interpolation is complete.
+.ih
+NOTES
+.ih
+SEE ALSO
+asiinit, asisinit
diff --git a/math/iminterp/doc/asigeti.hlp b/math/iminterp/doc/asigeti.hlp
new file mode 100644
index 00000000..0e9d1964
--- /dev/null
+++ b/math/iminterp/doc/asigeti.hlp
@@ -0,0 +1,36 @@
+.help asigeti Dec98 asigeti.hlp
+.ih
+NAME
+asigeti -- fetch an asi integer parameter
+.ih
+SYNOPSIS
+include <math/iminterp.h>
+
+value = asigeti (asi, param)
+
+.nf
+ pointer asi #I interpolant descriptor
+ int param #I parameter
+.fi
+.ih
+ARGUMENTS
+.ls asi
+Pointer to the sequential interpolant descriptor structure.
+.le
+.ls param
+The parameter to be fetched. The choices are: II_ASITYPE the interpolant
+type, II_ASINSAVE the length of the saved coefficient array, and
+II_ASINSINC the half-width of the sinc interpolant. The parameter
+definitions are contained in the package header file math/iminterp.h.
+.le
+.ih
+DESCRIPTION
+ASIGETI is used to determine the size of the coefficient array that
+must be allocated to save the sequential interpolant description structure,
+and to fetch selected sequential interpolant parameters.
+.ih
+NOTES
+.ih
+SEE ALSO
+asiinit, asisinit, asigetr
+.endhelp
diff --git a/math/iminterp/doc/asigetr.hlp b/math/iminterp/doc/asigetr.hlp
new file mode 100644
index 00000000..8a31deef
--- /dev/null
+++ b/math/iminterp/doc/asigetr.hlp
@@ -0,0 +1,36 @@
+.help asigetr Dec98 asigetr.hlp
+.ih
+NAME
+asigetr -- fetch an asi integer parameter
+.ih
+SYNOPSIS
+include <math/iminterp.h>
+
+value = asigetr (asi, param)
+
+.nf
+ pointer asi #I interpolant descriptor
+ int param #I parameter
+.fi
+.ih
+ARGUMENTS
+.ls asi
+Pointer to the sequential interpolant descriptor structure.
+.le
+.ls param
+The parameter to be fetched. The choices are: II_ASIBADVAL the undefined
+pixel value for the drizzle interpolant. The parameter definitions are
+contained in the package header file math/iminterp.h.
+.le
+.ih
+DESCRIPTION
+ASIGETR is used to set the value of undefined drizzle interpolant pixels.
+Undefined pixels are those for which the interpolation coordinates do not
+overlap the input coordinates, but are still, within the boundaries of the input
+image, a situation which may occur when the pixel fraction is < 1.0.
+.ih
+NOTES
+.ih
+SEE ALSO
+asiinit, asisinit, asigeti
+.endhelp
diff --git a/math/iminterp/doc/asigrl.hlp b/math/iminterp/doc/asigrl.hlp
new file mode 100644
index 00000000..4c3087bb
--- /dev/null
+++ b/math/iminterp/doc/asigrl.hlp
@@ -0,0 +1,40 @@
+.help asigrl Dec98 "Image Interpolator Package"
+.ih
+NAME
+asigrl -- integrate interpolant from a to b
+.ih
+SYNOPSIS
+integral = asigrl (asi, a, b)
+
+.nf
+ pointer asi #I interpolant descriptor
+ real a #I lower limit for integral, 1 <= a <= npix
+ real b #I upper limit for integral, 1 <= b <= npix
+.fi
+.ih
+ARGUMENTS
+.ls asi
+Pointer to the sequential interpolant descriptor structure.
+.le
+.ls a
+Lower limit to the integral, where 1 <= a <= npix.
+.le
+.ls b
+Upper limit to the integral, where 1 <= b <= npix.
+.le
+.ih
+DESCRIPTION
+The integral is calculated assuming that the interior polynomial, sinc, and
+drizzle interpolants are stored as the data points, and that the spline
+interpolant is stored as an array of B-spline coefficients.
+
+The integral of the sinc interpolant is computed by dividing the integration
+interval into a number of equal size subintervals which are at most one pixel
+wide. The integral of each subinterval is the central value times the interval
+width. The look-up table sinc interpolant is not supported and defaults to
+the sinc interpolant.
+.ih
+NOTES
+ASIINIT or ASISINIT and ASIFIT must be called before using ASIGRL.
+.ih
+SEE ALSO
diff --git a/math/iminterp/doc/asiinit.hlp b/math/iminterp/doc/asiinit.hlp
new file mode 100644
index 00000000..711d7969
--- /dev/null
+++ b/math/iminterp/doc/asiinit.hlp
@@ -0,0 +1,39 @@
+.help asiinit Dec98 "Image Interpolator Package"
+.ih
+NAME
+asiinit -- initialize the sequential interpolant descriptor using default parameters
+.ih
+SYNOPSIS
+include <math/iminterp.h>
+
+asiinit (asi, interp_type)
+
+.nf
+ pointer asi #O interpolant descriptor
+ int interp_type #I interpolant type
+.fi
+
+.ih
+ARGUMENTS
+.ls asi
+Pointer to sequential interpolant descriptor.
+.le
+.ls interp_type
+Interpolant type. The options are II_NEAREST, II_LINEAR, II_POLY3, II_POLY5,
+II_SPLINE3, II_SINC, II_LSINC, and II_DRIZZLE for nearest neighbor, linear,
+3rd order polynomial, 5th order polynomial, cubic spline, sinc, look-up
+table sinc, and drizzle respectively. The interpolant type definitions are
+found in the package header file math/iminterp.h.
+.le
+.ih
+DESCRIPTION
+The interpolant descriptor is allocated and initialized. The pointer asi is
+returned by ASIINIT. The sinc interpolant width defaults to 31 pixels. The
+look-up table sinc interpolant resolution defaults to 20 intervals or
+0.05 pixels. The drizzle pixel fraction defaults to 1.0.
+.ih
+NOTES
+ASIINIT or ASISINIT must be called before using any other ASI routines.
+.ih
+SEE ALSO
+asisinit, asifree
diff --git a/math/iminterp/doc/asirestore.hlp b/math/iminterp/doc/asirestore.hlp
new file mode 100644
index 00000000..6dd70262
--- /dev/null
+++ b/math/iminterp/doc/asirestore.hlp
@@ -0,0 +1,36 @@
+.help asirestore Dec98 "Image Interpolator Package"
+.ih
+NAME
+asirestore -- restore interpolant
+.ih
+SYNOPSIS
+asirestore (asi, interpolant)
+
+.nf
+ pointer asi #O interpolant descriptor
+ real interpolant[] #I array containing interpolant
+.fi
+.ih
+ARGUMENTS
+.ls asi
+Pointer to the interpolant descriptor structure.
+.le
+.ls interpolant
+Array containing the interpolant. The length of interpolant can be determined
+by a call to ASIGETI.
+.le
+
+.nf
+ len_interpolant = asigeti (asi, II_ASINSAVE)
+.fi
+.ih
+DESCRIPTION
+ASIRESTORE allocates space for the interpolant descriptor and restores the
+parameters and coefficients stored in the interpolant array to a structure
+for use with ASIEVAL, ASIVECTOR, ASIDER and ASIGRL.
+.ih
+NOTES
+.ih
+SEE ALSO
+asisave
+.endhelp
diff --git a/math/iminterp/doc/asisave.hlp b/math/iminterp/doc/asisave.hlp
new file mode 100644
index 00000000..7c8ff37a
--- /dev/null
+++ b/math/iminterp/doc/asisave.hlp
@@ -0,0 +1,39 @@
+.help asisave Dec98 "Image Interpolator Package"
+.ih
+NAME
+asisave -- save interpolant
+.ih
+SYNOPSIS
+asisave (asi, interpolant)
+
+.nf
+ pointer asi #I interpolant descriptor
+ real interpolant[] #O array containing the interpolant
+.fi
+.ih
+ARGUMENTS
+.ls asi
+Pointer to the interpolant descriptor structure.
+.le
+.ls interpolant
+Array where the interpolant is stored. The size of interpolant can be
+determined by a call to asigeti.
+.le
+
+.nf
+ len_interpolant = asigeti (asi, II_ASINSAVE)
+.fi
+.ih
+DESCRIPTION
+The interpolant type, number of coefficients and the position of
+the first data point in the coefficient array, along with various
+parameters such as the sinc interpolant width, sinc look-up table
+resolution, and drizzle pixel fraction, are stored in the first
+7 elements of the interpolant array. The remaining elements contain
+the coefficients calculated by ASIFIT.
+.ih
+NOTES
+.ih
+SEE ALSO
+asirestore
+.endhelp
diff --git a/math/iminterp/doc/asisinit.hlp b/math/iminterp/doc/asisinit.hlp
new file mode 100644
index 00000000..1ec9bc4e
--- /dev/null
+++ b/math/iminterp/doc/asisinit.hlp
@@ -0,0 +1,60 @@
+.help asisinit Dec98 "Image Interpolator Package"
+.ih
+NAME
+asisinit -- initialize the sequential interpolant descriptor using user parameters
+.ih
+SYNOPSIS
+include <math/iminterp.h>
+
+asisinit (asi, interp_type, nsinc, nincr, rparam, badval)
+
+.nf
+ pointer asi #O interpolant descriptor
+ int interp_type #I interpolant type
+ int nsinc #I sinc interpolant width in pixels
+ int nincr #I sinc look-up table resolution
+ real pixfrac #I sinc shift or drizzle pixel fraction
+ real badval #I drizzle undefined pixel value
+.fi
+
+.ih
+ARGUMENTS
+.ls asi
+Pointer to sequential interpolant descriptor.
+.le
+.ls interp_type
+Interpolant type. The options are II_NEAREST, II_LINEAR, II_POLY3, II_POLY5,
+II_SPLINE3, II_SINC, II_LSINC, and II_DRIZZLE for the nearest neighbour, linear,
+3rd order polynomial, 5th order polynomial, cubic spline, sinc, look-up
+table sinc, and drizzle interpolants respectively. The interpolant type
+definitions are found in the package header file math/iminterp.h.
+.le
+.ls nsinc
+The sinc and look-up table sinc interpolant width in pixels. Nsinc is
+rounded up internally to the nearest odd number.
+.le
+.ls nincr
+The look-up table sinc resolution in number of entries. Nincr = 10 implies
+a pixel resolution of 0.1 pixels, nincr = 20 a pixel resolution of 0.05
+pixels, etc.
+.le
+.ls pixfrac
+The look-up table sinc fractional pixel shift if nincr = 1 in which case
+-0.5 <= pixfrac <= 0.5 , or the drizzle pixel fraction in which case
+0.0 <= pixfrac <= 1.0. A minimum value of 0.001 is imposed on pixfrac.
+.le
+.ls badval
+The undefined pixel value for the drizzle interpolant. Pixels within
+the boundaries of the input image which do not overlap the input image
+pixels are assigned a value of badval.
+.le
+.ih
+DESCRIPTION
+The interpolant descriptor is allocated and initialized. The pointer asi is
+returned by ASISINIT.
+.ih
+NOTES
+ASIINIT or ASISINIT must be called before using any other ASI routines.
+.ih
+SEE ALSO
+asisinit, asifree
diff --git a/math/iminterp/doc/asitype.hlp b/math/iminterp/doc/asitype.hlp
new file mode 100644
index 00000000..d8c78b44
--- /dev/null
+++ b/math/iminterp/doc/asitype.hlp
@@ -0,0 +1,95 @@
+.help asitype Dec98 "Image Interpolator Package"
+.ih
+NAME
+asitype -- decode an interpolant string
+.ih
+SYNOPSIS
+include <math/iminterp.h>
+
+asitype (interpstr, interp_type, nsinc, nincr, rparam)
+
+.nf
+ char interpstr #I the input interpolant string
+ int interp_type #O interpolant type
+ int nsinc #O sinc interpolant width in pixels
+ int nincr #O sinc look-up table resolution
+ real rparam #O sinc or drizzle pixel fraction
+.fi
+
+.ih
+ARGUMENTS
+.ls interpstr
+The user supplied interpolant string to be decoded. The options are
+.ls nearest
+Nearest neighbor interpolation.
+.le
+.ls linear
+Linear interpolation
+.le
+.ls poly3
+Cubic polynomial interpolation.
+.le
+.ls poly5
+Quintic polynomial interpolation.
+.le
+.ls spline3
+Cubic spline interpolation.
+.le
+.ls sinc
+Sinc interpolation. Users can specify the sinc interpolant width by
+appending a width value to the interpolant string, e.g. sinc51 specifies
+a 51 pixel wide sinc interpolant. The sinc width will be rounded up to
+the nearest odd number. The default sinc width is 31.
+.le
+.ls lsinc
+Look-up table sinc interpolation. Users can specify the look-up table sinc
+interpolant width by appending a width value to the interpolant string, e.g.
+lsinc51 specifies a 51 pixel wide look-up table sinc interpolant. The user
+supplied sinc width will be rounded up to the nearest odd number. The default
+sinc width is 31 pixels. Users can specify the resolution of the sinc lookup
+up table by appending the look-up table size in square brackets to the
+interpolant string, e.g. lsinc51[20] specifies a 20 element sinc look-up
+table interpolant with a pixel resolution of 0.05 pixels. The default
+look-up table size and resolution are 20 and 0.05 pixels respectively.
+The fractional pixel shift for a 1 element look-up table sinc may be
+specified by replacing the number of lookup-table elements with the fractional
+shift, e.g. sinc51[0.5] will precompute a lookup table for a 0.5 pixel shift.
+.le
+.ls drizzle
+Drizzle interpolation. Users can specify the drizzle pixel fraction by
+appending the pixel fraction value to the interpolant string in square
+brackets, e.g. drizzle[0.5] specifies a pixel fraction of 0.5.
+The default pixel fraction is 1.0.
+.le
+.le
+.ls interp_type
+The output interpolant type. The options are II_NEAREST, II_LINEAR, II_POLY3,
+II_POLY5, II_SPLINE3, II_SINC, II_LSINC, and II_DRIZZLE for the nearest
+neighbor, linear, 3rd order polynomial, 5th order polynomial, cubic spline,
+sinc, look-up table sinc, and drizzle interpolants respectively. The
+interpolant type definitions are found in the package header file
+math/iminterp.h.
+.le
+.ls nsinc
+The output sinc and look-up table sinc interpolant width in pixels. The
+default value is 31 pixels.
+.le
+.ls nincr
+The output sinc look-up table size. Nincr = 10 implies a pixel resolution
+of 0.1 pixels, nincr = 20 a pixel resolution of 0.05 pixels, etc. The
+default value of nincr is 20.
+.le
+.ls rparam
+The output look-up table sinc fractional pixel shift if nincr = 1 in which case
+-0.5 <= rparam <= 0.5 , or the drizzle pixel fraction in which case
+0.0 <= rparam <= 1.0.
+.le
+.ih
+DESCRIPTION
+The interpolant string is decoded into values suitable for the ASISINIT
+or ASIINIT routines.
+.ih
+NOTES
+.ih
+SEE ALSO
+asinit, asisinit, asifree
diff --git a/math/iminterp/doc/asivector.hlp b/math/iminterp/doc/asivector.hlp
new file mode 100644
index 00000000..95bac138
--- /dev/null
+++ b/math/iminterp/doc/asivector.hlp
@@ -0,0 +1,52 @@
+.help asivector Dec98 "Image Interpolator Package"
+.ih
+NAME
+asivector -- evaluate the interpolant
+.ih
+SYNOPSIS
+asivector (asi, x, y, npix)
+
+.nf
+ pointer asi #I interpolator descriptor
+ real x[npix/2*npix] #I x array, 1 <= x[i] <= npix
+ real y[npix] #O array of interpolated values
+ int npix #I number of x values
+.fi
+.ih
+ARGUMENTS
+.ls asi
+Pointer to the sequential interpolator descriptor structure.
+.le
+.ls x
+Array of npix x values, or array of npix x ranges if the interpolant is
+drizzle.
+.le
+.ls y
+Array of interpolated values.
+.le
+.ls npix
+The number of x values.
+.le
+.ih
+DESCRIPTION
+The polynomial coefficients are calculated directly from the data points
+for the polynomial interpolants, and from the B-spline coefficients for
+the cubic spline interpolant. The actual calculation is done by adding and
+multiplying terms according to Everett's central difference interpolation
+formula. The boundary extension algorithm is projection.
+
+The sinc interpolant is computed using a range of data points around
+the desired position. Look-up table sinc interpolation is computed
+using the most appropriate entry in a precomputed look-up table.
+The boundary extension algorithm is nearest neighbor.
+
+The drizzle interpolant is computed by summing the data over the user
+supplied X intervals.
+.ih
+NOTES
+Checking for out of bounds and INDEF valued pixels is the responsibility of the
+user. ASIINIT or ASISINIT and ASIFIT must be called before calling ASIVECTOR.
+.ih
+SEE ALSO
+asieval, asider, arieval, arider
+.endhelp
diff --git a/math/iminterp/doc/im1dinterp.spc b/math/iminterp/doc/im1dinterp.spc
new file mode 100644
index 00000000..ce5b8680
--- /dev/null
+++ b/math/iminterp/doc/im1dinterp.spc
@@ -0,0 +1,525 @@
+.help iminterp Jul84 "Math Package"
+
+.ce
+Specifications for the Image Interpolator Package
+.ce
+Lindsey Davis
+.ce
+Vesa Junkkarinen
+.ce
+August 1984
+
+.sh
+1. Introduction
+
+ One of the most common operations in image processing is
+interpolation in a data array. Due to the large amount of data involved,
+efficiency is highly important. The advantage of having locally written
+interpolators, includes the ability to optimize for uniformly spaced data
+and the possibility of adding features that are useful to the final
+application.
+
+.sh
+2. Requirements
+
+.ls (1)
+The package shall take as input a one-dimensional array containing image
+data. The pixels are assumed to be equally spaced along a line.
+The coordinates of a pixel are assumed to be
+the same as the subscript of the pixel in the data array.
+The coordinate of the first pixel in the array and the spacing between pixels
+is assumed to be 1.0. All pixels are assumed to be good.
+Checking for INDEF valued and out of bounds pixels is the responsibility of the
+user. A routine to remove INDEF valued pixels from a data array shall be
+included in the package.
+.le
+.ls (2)
+The package is divided into array sequential interpolators and array
+random interpolators. The sequential interpolators have been optimized
+for returning many values as is the case when an array is shifted, or
+oversampled at many points in order to produce a
+smooth plot.
+The random interpolators allow the evaluation of a few interpolated
+points without the computing time and storage overhead required for
+setting up the sequential version.
+.le
+.ls (3)
+The quality of the interpolant will be set at run time. The options are:
+
+.nf
+ II_NEAREST - nearest neighbour
+ II_LINEAR - linear interpolation
+ II_POLY3 - 3rd order divided differences
+ II_POLY5 - 5th order divided differences
+ II_SPLINE3 - cubic spline
+.fi
+
+The calling sequences shall be invariant to the interpolant selected.
+Routines should be designed so that new interpolants can be added
+with minimal changes to the code and no change to the calling sequences.
+.le
+.ls (4)
+The interpolant parameters and the arrays necessary to store the coefficients
+are stored in a structure referenced by a pointer. The pointer is returned
+to the user program by the initial call to ASIINIT or ASIRESTORE and freed
+by a call to ASIFREE (see section 3.1).
+.le
+.ls (5)
+The package routines shall be able to:
+.ls o
+Calculate the coefficients of the interpolant and store these coefficients in
+the appropriate part of the interpolant descriptor structure.
+.le
+.ls o
+Evaluate the interplant at a given x(s) coordinate(s).
+.le
+.ls o
+Calculate the derivatives of the interpolant at a given value of x.
+.le
+.ls o
+Integrate the interpolant over a specified x interval.
+.le
+.ls o
+Store the interpolant in a user supplied array. Restore the saved interpolant
+to the interpolant descriptor structure for later use by ASIEVAL, ASIVECTOR,
+ASIDER and ASIGRL.
+.le
+
+.sh
+3. Specifications
+
+.sh
+3.1. The Array Sequential Interpolator Routines
+
+ The package prefix is asi and the package routines are:
+
+.nf
+ asiinit (asi, interp_type)
+ asifit (asi, datain, npix)
+ y = asieval (asi, x)
+ asivector (asi, x, yfit, npix)
+ asider (asi, x, der, nder)
+ v = asigrl (asi, a, b)
+ asisave (asi, interpolant)
+ asirestore (asi, interpolant)
+ asifree (asi)
+.fi
+
+.sh
+3.2. The Array Random Interpolator Routines
+
+ The package prefix is ari and the package routines are:
+
+.nf
+ y = arieval (x, datain, npix, interp_type)
+ arider (x, datain, npix, der, nder, interp_type)
+.fi
+
+.sh
+3.3. Miscellaneous
+
+ A routine has been included in the package to remove INDEF valued
+pixels from an array.
+
+.nf
+ arbpix (datain, dataout, npix, interp_type, boundary_type)
+.fi
+
+.sh
+3.4. Algorithms
+
+.sh
+3.4.1. Coefficients
+
+ The coefficient array used by the asi routines is calculated by ASIFIT.
+ASIFIT accepts an array of data, checks that the number
+of data points is appropriate for the interpolant selected, allocates
+space for the interpolant, and calculates the coefficients.
+Boundary coefficient values are calculated
+using boundary projection. With the exception of the cubic spline interpolant,
+the coefficients are stored as the data points.
+The B-spline coefficients are
+calculated using natural end conditions (Prenter 1975).
+After a call to ASIFIT the coefficient array contains the following.
+
+.nf
+ case II_NEAREST:
+
+ # no boundary conditions necessary
+ coeff[1] = datain[1]
+ .
+ .
+ .
+ coeff[npts] = datain[npix]
+
+ case II_LINEAR:
+
+ # coeff[npxix+1] required if x = npix
+ coeff[1] = datain[1]
+ .
+ .
+ .
+ coeff[npix] = datain[npix]
+ coeff[npix+1] = 2. * datain[npix] - datain[npix-1]
+
+ case II_POLY3:
+
+ # coeff[0] required if x = 1
+ # coeff[npix+1], coeff[npix+2] required if x = npix
+ coeff[0] = 2. * datain[1] - datain[2]
+ coeff[1] = datain[1]
+ .
+ .
+ .
+ coeff[npix] = datain[npix]
+ coeff[npix+1] = 2. * datain[npix] - datain[npix-1]
+ coeff[npix+2] = 2. * datain[npix] - datain[npix-2]
+
+ case II_POLY5:
+
+ # coeff[1], coeff[0] reqired if x = 1
+ # coeff[npix+1], coeff[npix+2], coeff[npix=3]
+ # required if x = npix
+
+ coeff[-1] = 2. * datain[1] - datain[3]
+ coeff[0] = 2. * datain[1] - datain[2]
+ coeff[1] = datain[1]
+ .
+ .
+ .
+ coeff[npix] = datain[npix]
+ coeff[npix+1] = 2. * datain[npix] - datain[npix-1]
+ coeff[npix+2] = 2. * datain[npix] - datain[npix-2]
+ coeff[npix+3] = 2. * datain[npix] - datain[npix-3]
+
+ case SPLINE3:
+
+ # coeff[0] = 2nd der at x = 1, coeff[0] = 0.
+ # coeff[npix+1] = 2nd der at x = npts, coeff[npix+1] = 0.
+ # coeff[npix+2] = 0., required if x = npix
+ coeff[0] = b[1]
+ coeff[1] = b[2]
+ .
+ .
+ .
+ coeff[npix] = b[npix+1]
+ coeff[npix+1] = b[npix+2]
+ coeff[npix+2] = 0.
+.fi
+
+.sh
+3.4.2. Evaluation
+
+ The ASIEVAL and ASIVECTOR routines have been optimized to be as efficient
+as possible. The values of the II_NEAREST and II_LINEAR interpolants
+are calculated directly. The II_SPLINE3 interpolant is evaluated using
+polynomial coefficients calculated directly from the B-spline coefficients
+(de Boor 1978). Values of the higher order polynomial interpolants
+are calculated using central differences. The equations for each case are
+listed below.
+
+.nf
+case II_NEAREST:
+
+ y = coeff[int (x + 0.5)]
+
+case II_LINEAR:
+
+ nx = x
+ y = (x - nx) * coeff[nx+1] + (nx + 1 - x) * coeff[nx]
+
+case II_POLY3:
+
+ nx = x
+ s = x - nx
+ t = 1. - s
+
+ # second central differences
+ cd20 = 1./6. * (coeff[nx+1] - 2. * coeff[nx] + coeff[nx-1])
+ cd21 = 1./6. * (coeff[nx+2] - 2. * coeff[nx+1] + coeff[nx])
+
+ y = s * (coeff[nx+1] + (s * s - 1.) * cd21) + t * (coeff[nx] +
+ (t * t - 1.) * cd20)
+
+case II_POLY5:
+
+ nx = x
+ s = x - nx
+ t = 1. - s
+
+ # second central differences
+ cd20 = 1./6. * (coeff[nx+1] - 2. * coeff[nx] + coeff[nx-1])
+ cd21 = 1./6. * (coeff[nx+2] - 2. * coeff[nx+1] + coeff[nx])
+
+ # fourth central diffreences
+ cd40 = 1./120. * (coeff[nx-2] - 4. * coeff[nx-1] + 6. * coeff[nx] - 4. *
+ coeff[nx+1] + a[nx+2])
+ cd41 = 1./120. * (coeff[nx-1] - 4. * coeff[nx] + 6. * coeff[nx+1] - 4. *
+ coeff[nx+2] + coeff[nx+3]
+
+ y = s * (coeff[nx+1] + (s * s - 1.) * (cd21 + (s * s - 4.) * cd41)) +
+ t * (coeff[nx] + (t * t - 1.) * (cd20 + (t * t - 4.) * cd40))
+
+case II_SPLINE3:
+
+ nx = x
+ s = x - nx
+
+ pc[1] = coeff[nx-1] + 4. * coeff[nx] + coeff[nx+1]
+ pc[2] = 3. * (coeff[nx+1] - coeff[nx-1])
+ pc[3] = 3. * (coeff[nx-1] - 2. * coeff[nx] + coeff[nx+1])
+ pc[4] = -coeff[nx-1] + 3. * coeff[nx] - 3. * coeff[nx+1] + coeff[nx+2]
+
+ y = pc[1] + s * (pc[2] + s * (pc[3] + s * pc[4]))
+.fi
+
+
+ The ARIEVAL routine uses the expressions above to evaluate the
+interpolant. However unlike ASIEVAL, ARIEVAL does not use a previously
+calculated coefficient array. Instead ARIEVAL selects the appropriate
+portion of the data array, calculates the coefficients and boundary
+coefficients if necessary, and evaluates the interpolant at the time it
+is called. The cubic spline interpolant uses at most SPLTS (currently 16)
+data points to calculate the B-spline coefficients.
+
+.sh
+3.4.3. Derivatives
+
+ Derivatives of the interpolant are calculated by evaluating the
+derivatives of the interpolating polynomial. For all interpolants der[1]
+equals the value of the interpolant at x.
+For the sake of efficiency the derivatives
+of the II_NEAREST and II_LINEAR interpolants are calculated directly.
+
+.nf
+ case II_NEAREST:
+
+ der[1] = coeff[int (x+0.5)]
+
+ case II_LINEAR:
+
+ der[1] = (x - nx) * coeff [nx+1] + (nx + 1 - x) * coeff[nx]
+ der[2] = coeff[nx+1] - coeff[nx]
+.fi
+
+ In order to calculate the derivatives of the cubic spline and
+polynomial interpolants
+the coefficients of the interpolating polynomial must be calculated.
+The polynomial
+coefficients for the cubic spline interpolant are computed directly from the
+B-spline coefficients (see 3.4.2.). The higher order polynomial
+interpolant coefficients are calculated as follows.
+
+First the appropriate portion of the coefficient array is loaded.
+
+.nf
+ do i = 1, nterms
+ d[i] = coeff[nx - nterms/2 + i]
+.fi
+
+Next the divided differences are calculated (Conte and de Boor 1972).
+
+.nf
+ do k = 1, nterms - 1
+ do i = 1, nterms - k
+ d[i] = (d[i+1] - d[i]) / k
+.fi
+
+The d[i] are the coefficients of an interpolating polynomial of the
+following form. The x[i] are the nterms data points surrounding the
+point of interest.
+
+.nf
+ p(x) = d[1] * (x-x[1]) * ... * (x-x[nterms-1) +
+ d[2] * (x-x[2]) * ... * (x-x[nterms-1]) + ... + d[nterms]
+.fi
+
+Next a new set of polynomial coefficients are calculated
+(Conte and de Boor 1972).
+
+.nf
+ do k = nterms, 2, -1
+ do i = 2, k
+ d[i] = d[i] + d[i-1] * (k - i - nterms/2)
+.fi
+
+The new d[i] are the coefficients of the follwoing polynomial.
+
+.nf
+ nx = x
+ p(x) = d[1] * (x-nx) ** (nterms-1) + d[2] * (x-nx) ** (nterms-2) + ...
+ d[nterms]
+.fi
+
+The d[i] array is flipped. The value and derivatives
+of the interpolant are then calculated using the d[i] array and
+nested multiplication.
+
+.nf
+ s = x - nx
+
+ do k = 1, nder {
+
+ accum = d[nterms-k+1]
+
+ do j = nterms - k, 1, -1
+ accum = d[j] + s * accum
+
+ der[k] = accum
+
+ # differnetiate
+ do j = 1, nterms - k
+ d[j] = j * d[j + 1]
+ }
+.fi
+
+ ARIDER calculates the derivatives of the interpolant using the same
+technique ASIDER. However ARIDER does not use a previously calculated
+coefficient array like ASIDER. Instead ARIDER selects the appropriate portion
+of the data array, calculates the coefficients and boundary coefficients,
+and computes the derivatives at the time it is called.
+
+.sh
+3.4.5. Integration
+
+ ASIGRL calculates the integral of the interpolant between fixed limits
+by integrating the interpolating polynomial. The coefficients of the
+interpolating polynomial are calculated as discussed in section 3.4.4.
+
+.sh
+4. Usage
+
+.sh
+4.1. User Notes
+
+The following series of steps illustrates the use of the package.
+
+.ls 4
+.ls (1)
+Insert an include <iminterp.h> statement in the calling program to make
+the IINTERP definitions available to the user program.
+.le
+.ls (2)
+Remove INDEF valued pixels from the data using ARBPIX.
+.le
+.ls (3)
+Call ASIINIT to initialize the interpolant parameters.
+.le
+.ls (4)
+Call ASIFIT to calculate the coefficients of the interpolant.
+.le
+.ls (5)
+Evaluate the interpolant at a given value of x(s) using ASIEVAL or
+ASIVECTOR.
+.le
+.ls (6)
+Calculate the derivatives and integral or the interpolant using
+ASIDER and ASIGRL.
+.le
+.ls (7)
+Free the interpolator structure by calling ASIFREE.
+.le
+.le
+
+ The interpolant can be saved and restored using ASISAVE and ASIRESTORE.
+If the values and derivatives of only a few points in an array are desired
+ARIEVAL and ARIDER can be called.
+
+.sh
+4.2. Examples
+
+.nf
+Example 1: Shift a data array by a constant amount
+
+ include <iminterp.h>
+ ...
+ call asiinit (asi, II_POLY5)
+ call asifit (asi, inrow, npix)
+
+ do i = 1, npix
+ outrow[i] = asieval (asi, i + shift)
+
+ call asifree (asi)
+ ...
+
+Example 2: Calculate the integral under the data array
+
+ include <iminterp.h>
+ ...
+ call asiinit (asi, II_POLY5)
+ call asifit (asi, datain, npix)
+
+ integral = asigrl (asi, 1. real (npix))
+
+ call asifree (asi)
+ ...
+
+Example 2: Store interpolant for later use by ASIEVAL
+ LEN_INTERP must be at least npix + 8 units long where npix is
+ is defined in the call to ASIFIT.
+
+ include <iminterp.h>
+
+ real interpolant[LEN_INTERP]
+ ...
+ call asiinit (asi, II_POLY3)
+ call asifit (asi, datain, npix)
+ call asisave (asi, interpolant)
+ call asifree (asi)
+ ...
+ call asirestore (asi, interpolant)
+ do i = 1, npts
+ yfit[i] = asieval (asi, x[i])
+ call asifree (asi)
+ ...
+.fi
+.sh
+5. Detailed Design
+
+.sh
+5.1. Interpolator Descriptor
+
+ The interpolant parameters and coefficients are stored in a
+structure listed below.
+
+.nf
+ define LEN_ASISTRUCT 4 # Length in structure units of
+ # interpolant descriptor
+
+ define ASI_TYPE Memi[$1] # Interpolant type
+ define ASI_NCOEFF Memi[$1+1] # No. of coefficients
+ define ASI_OFFSET Memi[$1+2] # First "data" point in
+ # coefficient array
+ define ASI_COEFF Memi[$1+3] # Pointer to coefficient array
+.fi
+
+.sh
+5.2. Storage Requirements
+
+ The interpolant descriptor requires LEN_ASISTRUCT storage units. The
+coefficient array requires ASI_NCOEFF(asi) real storage elements, where
+ASI_NCOEFF(asi) is defined as follows.
+
+.nf
+ ASI_NCOEFF(asi) = npix # II_NEAREST
+ ASI_NCOEFF(asi) = npix+1 # II_LINEAR
+ ASI_NCOEFF(asi) = npix+3 # II_POLY3
+ ASI_NCOEFF(asi) = npix+5 # II_POLY5
+ ASI_NCOEFF(asi) = npix+3 # II_SPLINE3
+.fi
+
+.sh
+6. References
+
+.ls (1)
+Carl de Boor, "A Practical Guide to Splines", 1978, Springer-Verlag New
+York Inc.
+.le
+.ls (2)
+S.D. Conte and C. de Boor, "Elementary Numerical Analysis", 1972, McGraw-Hill,
+Inc.
+.le
+.ls (3)
+P.M. Prenter, "Splines and Variational Methods", 1975, John Wiley and Sons Inc.
+.le
+.endhelp
diff --git a/math/iminterp/doc/im2dinterp.spc b/math/iminterp/doc/im2dinterp.spc
new file mode 100644
index 00000000..e4d88be3
--- /dev/null
+++ b/math/iminterp/doc/im2dinterp.spc
@@ -0,0 +1,432 @@
+.help iinterp Dec84 "Math Package"
+.ce
+Specifications for the 2D-Image Interpolator Package
+.ce
+Lindsey Davis
+.ce
+December 1984
+
+.sh
+1. Introduction
+
+One of the most common operations required in image processing is
+two-dimensional interpolation in a data array. Any image operator which
+physically moves pixels around requires image interpolation to calculate
+the new gray levels of the pixels. The advantage
+of having a locally written image interpolation package includes the ability to
+optimize for uniformly spaced data and the possibility of adding features
+that are useful to the final application.
+
+.sh
+2. Requirements
+
+.ls (1)
+The package shall take as input a 2-D array containing an image or image
+section. The pixels are assumed to be equally spaced on a rectangular grid.
+The coordinates of the pixels
+are assumed to be the same as the subscripts of the pixel in the data array.
+Therefore the coordinates of the first pixel in the array are assumed
+to be (1,1). For operations on image sections the calling program must
+keep track of any transformations between image and section coordinates.
+All pixels are assumed to be good.
+Checking for INDEF and out of bounds pixels is
+the responsibility of the user. A routine to remove INDEF valued pixels
+ARBPIX is available in the 1-D package.
+.le
+.ls (2)
+The package is divided into array sequential interpolants and array random
+interpolants. The sequential interpolants have been optimized for returning
+many values as when an array is shifted or when an array is oversampled
+at many points in order to produce a smooth surface plot.
+The random interpolants
+allow the evaluation of a few interpolated points without the computing
+time and storage overhead required for setting up the sequential version.
+.le
+.ls (3)
+The quality of the interpolant will be set at run time. The options are:
+
+.nf
+ II_BINEAREST # nearest neighbour
+ II_BILINEAR # bilinear
+ II_BIPOLY3 # bicubic interior polynomial
+ II_BIPOLY5 # biquintic interior polynomial
+ II_BISPLINE3 # bicubic spline
+.fi
+
+The calling sequences shall be invariant to the interpolant selected.
+Routines should be designed so that new interpolants can be added with
+minimal changes to the code and no changes to the calling sequences.
+.le
+.ls (4)
+The interpolant parameters and the arrays necessary to store the
+coefficients are stored in a structure referenced by a pointer. The pointer
+is returned to the user program by the initial call to MSIINIT or
+MSIRESTORE and freed by a call to MSIFREE.
+.le
+.ls (5)
+The package routines should be able to:
+.ls o
+Calculate the coefficients of the interpolant and store these coefficients
+in the appropriate part of the interpolant descriptor structure.
+.le
+.ls o
+Evaluate the interpolant at a given x-y(s) coordinate.
+.le
+.ls o
+Evaluate the first nder derivatives at a given value of x-y.
+.le
+.ls o
+Integrate the interpolant over a specified interval in x-y.
+.le
+
+.sh
+3. Specifications
+
+.sh
+3.1. The Matrix Sequential Interpolator Routines
+
+The package prefix is msi and the package routines are:
+
+.nf
+ msiinit (msi, interp_type)
+ msifit (msi, datain, nxpix, nypix, len_datain)
+ y = msieval (msi, x, y)
+ msivector (msi, x, y, zfit, npix)
+ msigrid (msi, x, y, zfit, nx, ny, len_zfit)
+ msider (msi, x, y, der, nxder, nyder, len_der)
+ v = msigrl (msi, x, y, npts)
+ v = msisqgrl (msi, x1, x2, y1, y2)
+ msisave (msi, interpolant)
+ msirestore (msi, interpolant)
+ msifree (msi)
+.fi
+
+.sh
+3.2. The Matrix Random Interpolator Routines
+
+The package prefix is mri and the package routines are:
+
+.nf
+ y = mrieval (x, y, datain, nx, ny, len_datain, interp_type)
+ mrider (x, y, datain, nx, ny, len_datain, der, nxder, nyder,
+ len_der, interp_type)
+.fi
+
+.sh
+3.3. Miscellaneous
+
+A routine ARBPIX will remove INDEF valued pixels from an
+array and is available in the 1-D package.
+
+.nf
+ arbpix (datain, dataout, npix, interp_type, boundary_type)
+.fi
+
+.sh
+3.4. Algorithms
+
+.sh
+3.4.1. Coefficients
+
+The coefficients for the msi routines are calculated by MSIFIT. MSIFIT
+accepts the input data, checks that the number of data pixels is
+appropriate for the interpolant selected, and allocates space for the
+interpolant. Boundary coefficient values are calculated using boundary
+projection. With the exception of the II_BISPLINE3 option, the interpolant
+is stored as the data points. The B-spline coefficients are calculated
+using the "natural" end conditions in two steps.
+First the B-spline coefficients in x at each value of y are calculated.
+The B-x coefficients are then solved for the B-spline coefficients in x-y.
+After a call to MSIFIT the coefficient
+array contains the following.
+
+.nf
+ CASE II_BINEAREST:
+
+ # no boundary extension required, coeff[1:nx,1:ny]
+
+ coeff[i,j] = data[i,j] i = 1,...,nx
+ j = 1,...,ny
+
+ case II_BILINEAR:
+
+ # must extend nx by 1 and ny by 1, coeff[1:nx+1,1:ny+1]
+
+ coeff[i,j] = data[i,j] i = 1,...,nx
+ j = 1,...,ny
+
+ coeff[nx+1,j] = 2 * data[nx] - data[nx-1] j = 1,...,ny
+ coeff[i,ny+1] = 2 * data[ny] - data[ny-1] i = 1,...,nx
+
+ coeff[nx+1,ny+1] = 2 * coeff[nx+1,ny] - data[nx,ny]
+
+ case II_BIPOLY3:
+
+ # must extend nx by -1 and 2 and ny by -1 and 2, coeff[0:nx+2,0:ny+2]
+
+ coeff[i,j] = data[i,j] i = 1,...,nx
+ j = 1,...,ny
+
+ coeff[0,j] = 2 * data[1,j] - data[2,j] j = 1,...,ny
+ coeff[nx+1,j] = 2 * data[nx,j] - data[nx-1,j] j = 1,...,ny
+ coeff[nx+2,j] = 2 * data[nx,j] - data[nx-2,j] j = 1,...,ny
+
+ coeff[i,0] = 2 * data[i,1] - data[i,2] i = 1,...,ny
+ coeff[i,ny+1] = 2 * data[i,ny] - data[i,ny-1] i = 1,...,nx
+ coeff[i,ny+2] = 2 * data[i,ny] - data[i,ny-2] i = 1,...,nx
+
+ # plus remaining points
+
+ case II_BIPOLY5:
+
+ # extend -2 and 3 in nx and -2 and 3 in ny, coeff[-1:nx+3,-1:ny+3]
+
+ coeff[i,j] = data[i,j] i = 1,...,nx
+ j = 1,...,ny
+
+ coeff[-1,j] = 2 * data[1,j] - data[3,j] j = 1,...,ny
+ coeff[0,j] = 2 * data[1,j] - data[2,j] j = 1,...,ny
+ coeff[nx+1,j] = 2 * data[nx,j] - data[nx-1,j] j = 1,...,ny
+ coeff[nx+2,j] = 2 * data[nx,j] - data[nx-2,j] j = 1,...,ny
+ coeff[nx+3,j] = 2 * data[nx,j] - data[nx-3,j] j = 1,...,ny
+
+ coeff[i,-1] = 2 * data[i,1] - data[i,3] i = 1,...,nx
+ coeff[i,0] = 2 * data[i,1] - data[i,2] i = 1,...,nx
+ coeff[i,ny+1] = 2 * data[i,ny] - data[i,ny-1] i = 1,...,nx
+ coeff[i,ny+2] = 2 * data[i,ny] - data[i,ny-2] i = 1,...,nx
+ coeff[i,ny+3] = 2 * data[i,ny] - data[i,ny-3] i = 1,...,nx
+
+ # plus remaining conditions
+
+ case II_BISPLINE3:
+
+ # the natural boundary conditions are used, coeff[0:nx+2,0:ny+2]
+
+ coeff[i,j] = data[i,j] i = 1,...,nx
+ j = 1,...,ny
+
+ coeff[i,0] = 0. i = 0,...,nx+2
+ coeff[i,ny+1] = 0. i = 0,...,nx+2
+ coeff[i,ny+2] = 0. i = 0,...,nx+2
+
+ coeff[0,j] = 0. j = 1,...,ny
+ coeff[nx+1,j] = 0. j = 1,...,ny
+ coeff[nx+2,j] = 0. j = 1,...,ny
+
+ # plus remaining coefficients
+
+.fi
+
+.sh
+3.4.2. Evaluation
+
+The MSIEVAL and MSIVECTOR routines will be optimized to be as efficient
+as possible for evaluating arbitrarily spaced data. A special function
+MSIGRID is included for evaluating closely spaced points on a rectangular grid.
+For the options II_BINEAREST and II_BILINEAR the value of the
+interpolant is calculated directly. The II_BISPLINE3 interpolant is
+evaluated using polynomial coefficients calculated directly from the
+B-spline coefficients. Values of the higher order polynomial interpolants
+are calculated using Everett's central difference formula. The equations
+are listed below.
+
+.nf
+case II_BINEAREST:
+
+ z = coeff[int (x + 0.5), int (y + 0.5))
+
+case II_BILINEAR:
+
+ sx = x - nx sy = y - ny
+ tx = 1. - sx ty = 1. - sy
+
+ z = tx * ty * coeff[nx,ny] + sx * ty * coeff[nx+1,ny] +
+ sy * tx * coeff[nx,ny+1] + sx * sy * coeff[nx+1,ny+1]
+
+
+case II_BIPOLY3:
+
+ nx = x
+ sx = x - nx
+ tx = 1. - sx
+
+ # interpolate in x
+
+ i = 1
+
+ do j = ny-1, ny+2 {
+
+ cd20[i] = 1./6. * (coeff[nx+1,j] - 2. * coeff[nx,j] + coeff[nx-1,j])
+ cd21[i] = 1./6. * (coeff[nx+2,j] - 2. * coeff[nx+1,j] + coeff[nx,j])
+
+ z[i] = sx * (coeff[nx+1,j] + (sx * sx - 1.) * cd21[i]) +
+ tx * (coeff[nx,j] + (tx * tx - 1.) * cd20[i])
+
+ i = i + 1
+ }
+
+ ny = y
+ sy = y - ny
+ ty = 1. - sy
+
+ # interpolate in y
+
+ cd20y = 1./6. * (z[3] - 2. * z[2] + z[1])
+ cd21y = 1./6. * (z[4] - 2. * z[3] + z[2])
+
+ value = sy * (z[3] + (sy * sy - 1.) * cd21y) +
+ ty * (z[2] + (ty * ty - 1.) *cd20y)
+
+
+case II_BIPOLY5:
+
+ nx = x
+ sx = x - nx
+ sx2 = sx * sx
+ tx = 1. - sx
+ tx2 = tx * tx
+
+ # interpolate in x
+
+ i = 1
+
+ do j = ny-2, ny+3 {
+
+ cd20[i] = 1./6. * (coeff[nx+1,j] - 2. * coeff[nx,j] + coeff[nx-1,j])
+ cd21[i] = 1./6. * (coeff[nx+2,j] - 2. * coeff[nx+1,j] + coeff[nx,j])
+ cd40[i] = 1./120. * (coeff[nx-2,j] - 4. * coeff[nx-1,j] +
+ 6. * coeff[nx,j] - 4. * coeff[nx+1,j] + coeff[nx+2,j])
+ cd41[i] = 1./120. * (coeff[nx-1,j] - 4. * coeff[nx,j] +
+ 6. * coeff[nx+1,j] - 4. * coeff[nx+2,j] + coeff[nx+3,j])
+
+ z[i] = sx * (coeff[nx+1,j] + (sx2 - 1.) * (cd21[j] + (sx2 - 4.) *
+ cd41[j])) + tx * (coeff[nx,j] + (tx2 - 1.) *
+ (cd20[j] + (tx2 - 4.) * cd40[j]))
+
+ i = i + 1
+ }
+
+ ny = y
+ sy = y - ny
+ sy2 = sy * sy
+ ty = 1. - sy
+ ty2 = ty * ty
+
+ # interpolate in y
+
+ cd20y = 1./6. * (z[3] - 2. * z[2] + z[1])
+ cd21y = 1./6. * (z[4] - 2. * z[3] + z[2])
+ cd40y = 1./120. * (z[1] - 4. * z[2] + 6. * z[3] - 4. * z[4] + z[5])
+ cd41y = 1./120. * (z[2] - 4. * z[3] + 6. * z[4] - 4. * z[5] + z[6])
+
+ value = sy * (z[4] + (sy2 - 1.) * (cd21y + (sy2 - 4.) * cd41y)) +
+ ty * (z[3] + (ty2 - 1.) * (cd20y + (ty2 - 4.) * cd40y))
+
+case II_BISPLINE3:
+
+ # use the B-spline representation
+
+ value = coeff[i,j] * B[i](x) * B[j](y)
+
+ # the B-splines are the following
+
+ B1(x) = (x - x1) ** 3
+ B2(x) = 1 + 3 * (x - x2) + 3 * (x - x2) ** 2 - 3 * (x - x2) ** 3
+ B3(x) = 1 + 3 * (x3 - x) + 3 * (x3 - x) ** 2 - 3 * (x3 - x) ** 3
+ B4(x) = (x4 - x) ** 3
+
+ # the y B-splines are identical
+
+.fi
+
+.sh
+3.4.3. Derivatives
+
+The derivatives are calculated by evaluating the derivatives of the
+interpolating polynomial. The 0-th derivative is the value of the
+interpolant. The 1-st derivatives are d/dx and d/dy. The second are
+d2/dx2, d2/dy2 and d2/dxdy. The derivatives in the case II_BINEAREST
+and II_BILINEAR are calculated directly.
+
+.nf
+case II_BINEAREST:
+
+ der[1,1] = value # see previous section
+
+case II_BILINEAR:
+
+ der[1] = value # see previous section
+
+ # d/dx
+ der[2,1] = -ty * coeff[nx,ny] + ty * coeff[nx+1,ny] -
+ sy * coeff[nx,ny+1] + sy * coeff[nx+1,ny+1]
+
+ # d/dy
+ der[1,2] = -tx * coeff[nx,ny] + sx * coeff[nx+1,ny] +
+ tx * coeff[nx,ny+1] + sx * coeff[nx+1,ny+1]
+
+ # d2/dxdy
+ der[2,2] = coeff[nx,ny] - coeff[nx+1,ny] - coeff[nx,ny+1] + coeff[nx+1,ny+1]
+
+case II_BIPOLY3, II_BIPOLY5, II_BISPLINE3:
+.fi
+
+For the higher order interpolants the coefficients of the interpolating
+polynomial in x and y are calculated. In the case of II_BIPOLY3 and II_BIPOLY5
+this is the Everett polynomial of the 3-rd and 5-th order respectively.
+In the case of II_BISPLINE3 the pp-representation is calculated for
+the B-spline coefficients. The value of the interpolant and its
+derivatives is calculated using nested multiplication.
+
+.sh
+3.4.5. Integration
+
+Integration is most easily accomplished by integrating the interpolant in
+x for each value of y. The resulting function of y can then be integrated
+in y to derive the 2-D integral. The limits of the integral are assumed
+to be the corners of a polygon. At minimum of three points which define
+a triangular region in x-y are required. The integral will be an approximation.
+A special function for integrating over a rectangular region is also
+provided.
+
+.sh
+4. Detailed Design
+
+.sh
+4.1. Interpolant Descriptor
+
+The interpolant parameters and coefficients will be stored in a structure
+as listed below.
+
+.nf
+ define LEN_MSISTRUCT 10
+
+ struct {
+
+ int msi_type # interpolant type
+ int msi_nxcoeff # size of coefficient array in x
+ int msi_nycoeff # size of coefficient array in y
+ int msi_fstpnt # first datapoint in coefficient array
+ int asi_coeff # pointer to 1-D coefficient array
+
+ } msistruct
+
+.fi
+
+.sh
+4.2. Storage Requirements
+
+The interpolant descriptor requires LEN_MSISTRUCT storage units.
+The coefficient array storage is dynamically allocated and requires
+msi_nxcoeff * msi_nycoeff real storage elements. The requirements for
+each interpolant type are listed below.
+
+.nf
+ II_BINEAREST nxpix * nypix
+ II_BILINEAR (nxpix + 1) * (nypix + 1)
+ II_BIPOLY3 (nxpix + 3) * (nypix + 3)
+ II_BIPOLY5 (nxpix + 5) * (nypix + 5)
+ II_BISPLINE3 (nxpix + 3) * (nypix + 3)
+.fi
+
+.endhelp
diff --git a/math/iminterp/doc/iminterp.hd b/math/iminterp/doc/iminterp.hd
new file mode 100644
index 00000000..de836c3e
--- /dev/null
+++ b/math/iminterp/doc/iminterp.hd
@@ -0,0 +1,37 @@
+# Help directory for the IMINTERP (image interpolator) package.
+
+$iminterp = "math$iminterp/"
+
+arbpix hlp = arbpix.hlp, src = iminterp$arbpix.x
+arider hlp = arider.hlp, src = iminterp$arider.x
+arieval hlp = arieval.hlp, src = iminterp$arieval.x
+asider hlp = asider.hlp, src = iminterp$asider.x
+asieval hlp = asieval.hlp, src = iminterp$asieval.x
+asifit hlp = asifit.hlp, src = iminterp$asifit.x
+asifree hlp = asifree.hlp, src = iminterp$asifree.x
+asigeti hlp = asigeti.hlp, src = iminterp$asigeti.x
+asigetr hlp = asigetr.hlp, src = iminterp$asigetr.x
+asigrl hlp = asigrl.hlp, src = iminterp$asigrl.x
+asiinit hlp = asiinit.hlp, src = iminterp$asiinit.x
+asirestore hlp = asirestore.hlp, src = iminterp$asirestore.x
+asisave hlp = asisave.hlp, src = iminterp$asisave.x
+asisinit hlp = asisinit.hlp, src = iminterp$asisinit.x
+asivector hlp = asivector.hlp, src = iminterp$asivector.x
+asitype hlp = asitype.hlp, src = iminterp$asitype.x
+
+mrider hlp = mrider.hlp, src = iminterp$mrider.x
+mrieval hlp = mrieval.hlp, src = iminterp$mrieval.x
+msider hlp = msider.hlp, src = iminterp$msider.x
+msieval hlp = msieval.hlp, src = iminterp$msieval.x
+msifit hlp = msifit.hlp, src = iminterp$msifit.x
+msifree hlp = msifree.hlp, src = iminterp$msifree.x
+msigeti hlp = msigeti.hlp, src = iminterp$msigeti.x
+msigetr hlp = msigetr.hlp, src = iminterp$msigetr.x
+msigrid hlp = msigrid.hlp, src = iminterp$msigrid.x
+msigrl hlp = msigrl.hlp, src = iminterp$msigrl.x
+msiinit hlp = msiinit.hlp, src = iminterp$msiinit.x
+msirestore hlp = msirestore.hlp, src = iminterp$msirestore.x
+msisave hlp = msisave.hlp, src = iminterp$msisave.x
+msisinit hlp = msisinit.hlp, src = iminterp$msisinit.x
+msitype hlp = msitype.hlp, src = iminterp$msitype.x
+msivector hlp = msivector.hlp, src = iminterp$msivector.x
diff --git a/math/iminterp/doc/iminterp.hlp b/math/iminterp/doc/iminterp.hlp
new file mode 100644
index 00000000..c602b43e
--- /dev/null
+++ b/math/iminterp/doc/iminterp.hlp
@@ -0,0 +1,234 @@
+.help iminterp Dec98 "Math Package"
+.ih
+NAME
+iminterp -- image interpolator package
+.ih
+SYNOPSIS
+
+.nf
+ asitype (interpstr, interp_type, nsinc, nincr, rparam)
+ asiinit (asi, interp_type)
+ asisinit (asi, interp_type, nsinc, nincr, rparam)
+ asifit (asi, datain, npix)
+ivalue = asigeti (asi, param)
+rvalue = asigetr (asi, param)
+ y = asieval (asi, x)
+ asivector (asi, x, yfit, npix)
+ asider (asi, x, der, nder)
+ v = asigrl (asi, a, b)
+ asisave (asi, interpolant)
+ asirestore (asi, interpolant)
+ asifree (asi)
+
+ y = arieval (x, datain, npix, interp_type)
+ arider (x, datain, npix, der, nder, interp_type)
+
+ arbpix (datain, dataout, npix, interp_type, boundary_type)
+.fi
+
+.nf
+ msitype (interpstr, interp_type, nsinc, nincr, rparam)
+ msisinit (msi, interp_type, nsinc, nxincr, nyincr, rparam1, rparam2)
+ msiinit (msi, interp_type)
+ msifit (msi, datain, nxpix, nypix, len_datain)
+ivalue = msigeti (msi, param)
+rvalue = msigetr (msi, param)
+ y = msieval (msi, x, y)
+ msivector (msi, x, y, zfit, npts)
+ msider (msi, x, y, der, nxder, nyder, len_der)
+ v = msigrl (msi, x, y, npts)
+ v = msisqgrl (msi, x1, x2, y1, y2)
+ msisave (msi, interpolant)
+ msirestore (msi, interpolant)
+ msifree (msi)
+
+ y = mrieval (x, y, datain, nxpix, nypix, len_dataina, interp_type)
+ mrider (x, y, datain, nxpix, nypix, len_datain, der, nxder, nyder,
+ len_der, interp_type)
+.fi
+
+.ih
+DESCRIPTION
+The iminterp package provides a set of routines for interpolating uniformly
+spaced data assuming that the spacing between data points is 1.0. The
+package is divided into 1D and 2D array sequential interpolants,
+prefixes asi and msi, and 1D and 2D
+array random interpolants, prefixes ari and mri.
+The sequential interpolants have
+been optimized for returning many values as is the case when an array is
+shifted. The random interpolants allow evaluation of a few interpolated
+points without the computing time and storage overhead required for
+setting up the sequential version.
+.ih
+NOTES
+The interpolant is chosen at run time from the following list.
+
+.nf
+ II_NEAREST # nearest neighbour in x
+ II_LINEAR # linear interpolation in x
+ II_POLY3 # 3rd order interior polynomial in x
+ II_POLY5 # fifth order interior polynomial in x
+ II_SPLINE3 # cubic spline in x
+ II_SINC # sinc interpolation in x
+ II_LSINC # look-up table sinc interpolation in x
+ II_DRIZZLE # drizzle interpolation in x
+
+ II_BINEAREST # nearest neighbour in x and y
+ II_BILINEAR # bilinear interpolation
+ II_BIPOLY3 # 3rd order interior polynomial in x and y
+ II_BIPOLY5 # 5th order interior polynomial in x and y
+ II_BISPLINE3 # bicubic spline
+ II_BISINC # sinc interpolation in x and y
+ II_BILSINC # look-up table sinc interpolation in x and y
+ II_BIDRIZZLE # drizzle interpolation in x and y
+.fi
+
+The routines assume that all x (1D, 2D) and y (2D) values of interest lie in
+the region 1 <= x <= nxpix, 1 <= y <= nypix.
+Checking for out of bounds x and/or y values is the responsibility
+of the calling program. The asi, ari, msi, and mri routines assume that INDEF
+valued pixels have been removed from the data prior to entering the
+package. The routine ARBPIX has been added to the package to facilitate
+INDEF valued pixel removal.
+
+In order to make the package definitions available to the calling program
+an include <math/iminterp.h> statement must appear in the calling program.
+Either ASIINIT, ASISINIT or ASIRESTORE must be called before using the
+asi routines. ASIFREE frees the space used by the asi routines. For the
+msi routines the corresponding examples are MSIINIT, MSISINIT, MSIRESTORE
+and MSIFREE.
+.ih
+EXAMPLES
+.nf
+Example 1: Shift a 1D data array by a constant amount using a 5th order
+polynomial interpolant and the drizzle routine respectively. Note that
+in this example the drizzle interpolant is equivalent to the linear
+interpolant since the default drizzle pixel fraction is 1.0 and there
+is no scale change. Out-of-bounds pixels are set to 0.0
+
+ include <math/iminterp.h>
+ ...
+ call asiinit (asi, II_POLY5)
+ call asifit (asi, inrow, npix)
+
+ do i = 1, npix
+ if (i + xshift < 1.0 || i + xshift > npix)
+ outrow[i] = 0.0
+ else
+ outrow[i] = asieval (asi, i + xshift)
+
+ call asifree (asi)
+ ...
+
+
+ include <math/iminterp.h>
+
+ real tmpx[2]
+ ...
+ call asiinit (asi, II_DRIZZLE)
+ call asifit (asi, inrow, npix)
+
+ do i = 1, npix
+ tmpx[1] = i + xshift - 0.5
+ tmpx[2] = i + xshift + 0.5
+ if (tmpx[1] < 1 || tmpx[2] > npix)
+ outrow[i] = 0.0
+ else
+ outrow[i] = asieval (asi, tmpx)
+
+ call asifree (asi)
+ ...
+
+
+Example 2: Shift a 2D array by a constant amount using a 3rd order polynomial
+interpolant and the drizzle interpolant respectively. Note that
+in this example the drizzle interpolant is equivalent to the linear
+interpolant since the default drizzle pixel fraction is 1.0 and there
+is no scale change. Out-of-bounds pixels are set to 0.0.
+
+ include <math/iminterp.h>
+ ...
+ call msiinit (msi, II_BIPOLY3)
+ call msifit (msi, insection, nxpix, nypix, nxpix)
+
+ do j = 1, nypix
+ if (j + yshift < 1 || j + yshift > nypix)
+ do i = 1, nxpix
+ outsection[i,j] = 0.0
+ else
+ do i = 1, nxpix
+ if (i + xshift < 1 || i + xshift > nxpix)
+ outsection[i,j] = 0.0
+ else
+ outsection[i,j] = msieval (msi, i + xshift, j + yshift)
+
+ call msifree (msi)
+ ...
+
+
+ include <math/iminterp.h>
+ ...
+ real tmpx[4], tmpy[4]
+ ...
+ call msiinit (msi, II_BIDRIZZLE)
+ call msifit (msi, insection, nxpix, nypix, nxpix)
+
+ do j = 1, nypix {
+ tmpy[1] = j + yshift - 0.5
+ tmpy[2] = j + yshift - 0.5
+ tmpy[3] = j + yshift + 0.5
+ tmpy[4] = j + yshift + 0.5
+ if (tmpy[1] < 1 || tmpy[4] > nypix)
+ do i = 1, nxpix
+ outsection[i,j] = 0.0
+ else
+ do i = 1, nxpix
+ tmpx[1] = i + xshift - 0.5
+ tmpx[2] = i + xshift + 0.5
+ tmpx[3] = i + xshift + 0.5
+ tmpx[4] = i + xshift - 0.5
+ if (tmpx[1] < 1 || tmpx[2] > nxpix)
+ outsection[i,j] = 0.0
+ else
+ outsection[i,j] = msieval (msi, tmpx, tmpy)
+ }
+
+ call msifree (msi)
+ ...
+
+
+Example 3: Calculate the integral under a 1D data array
+
+ include <math/iminterp.h>
+ ...
+ call asiinit (asi, II_POLY5)
+ call asifit (asi, datain, npix)
+
+ integral = asigrl (asi, 1. real (npix))
+
+ call asifree (asi)
+ ...
+
+Example 4: Store a 1D interpolant for later use by ASIEVAL
+
+ include <math/iminterp.h>
+
+ ...
+ call asiinit (asi, II_POLY3)
+ call asifit (asi, datain, npix)
+
+ len_interpolant = asigeti (asi, ASINSAVE)
+ call salloc (interpolant, len_interpolant, TY_REAL)
+ call asisave (asi, Memr[interpolant])
+
+ call asifree (asi)
+ ...
+ call asirestore (asi, Memr[interpolant])
+
+ do i = 1, npts
+ yfit[i] = asieval (asi, x[i])
+
+ call asifree (asi)
+ ...
+.fi
+.endhelp
diff --git a/math/iminterp/doc/iminterp.men b/math/iminterp/doc/iminterp.men
new file mode 100644
index 00000000..9daf5883
--- /dev/null
+++ b/math/iminterp/doc/iminterp.men
@@ -0,0 +1,32 @@
+ arbpix - Replace bad pixels in an array
+ arider - Calculate derivatives for a few points in a 1-D array
+ arieval - Evaluate the interpolant at a few points in a 1-D array
+ asider - Evaluate derivatives at x
+ asieval - Evaluate the interpolant at x
+ asifit - Fit the 1-D interpolant
+ asifree - Free space allocated by asiinit or asisinit
+ asigeti - Fetch an integer parameter
+ asigrl - Calculate the integral under an array
+ asiinit - Initialize 1-D interpolant using default parameters
+ asirestore - Restore interpolant parameters and coefficients
+ asisave - Save interpolant parameters and coefficients
+ asisinit - Initialize 1-D interpolant using user parameters
+ asitype - Decode string into interpolant type and parameters
+ asivector - Evaluate interpolant at an array of x
+
+ mrider - Calculate derivatives for a few points in a 2-D array
+ mrieval - Evaluate the interpolant at a few points in a 2-D array
+ msider - Evaluate derivatives at x and y
+ msieval - Evaluate the interpolant at x and y
+ msifit - Fit the 2-D interpolant
+ msifree - Free space allocated by msiinit or msisinit
+ msigeti - Fetch an integer parameter
+ msigetr - Fetch a real parameter
+ msigrl - Calculate the integral inside a polygon
+ msiinit - Initialize the 2-D interpolant using default parameters
+ msirestore - Restore interpolant parameters and coefficients
+ msisave - Save 2-D interpolant parameters and coefficients
+ msisinit - Initialize the 2-D interpolant using user parameters
+ msisqgrl - Calculate the integral inside a rectangle
+ msitype - Decode string into interpolant type and parameters
+ msivector - Evaluate interpolant at an array of x and y
diff --git a/math/iminterp/doc/iminterp.spc b/math/iminterp/doc/iminterp.spc
new file mode 100644
index 00000000..ce5b8680
--- /dev/null
+++ b/math/iminterp/doc/iminterp.spc
@@ -0,0 +1,525 @@
+.help iminterp Jul84 "Math Package"
+
+.ce
+Specifications for the Image Interpolator Package
+.ce
+Lindsey Davis
+.ce
+Vesa Junkkarinen
+.ce
+August 1984
+
+.sh
+1. Introduction
+
+ One of the most common operations in image processing is
+interpolation in a data array. Due to the large amount of data involved,
+efficiency is highly important. The advantage of having locally written
+interpolators, includes the ability to optimize for uniformly spaced data
+and the possibility of adding features that are useful to the final
+application.
+
+.sh
+2. Requirements
+
+.ls (1)
+The package shall take as input a one-dimensional array containing image
+data. The pixels are assumed to be equally spaced along a line.
+The coordinates of a pixel are assumed to be
+the same as the subscript of the pixel in the data array.
+The coordinate of the first pixel in the array and the spacing between pixels
+is assumed to be 1.0. All pixels are assumed to be good.
+Checking for INDEF valued and out of bounds pixels is the responsibility of the
+user. A routine to remove INDEF valued pixels from a data array shall be
+included in the package.
+.le
+.ls (2)
+The package is divided into array sequential interpolators and array
+random interpolators. The sequential interpolators have been optimized
+for returning many values as is the case when an array is shifted, or
+oversampled at many points in order to produce a
+smooth plot.
+The random interpolators allow the evaluation of a few interpolated
+points without the computing time and storage overhead required for
+setting up the sequential version.
+.le
+.ls (3)
+The quality of the interpolant will be set at run time. The options are:
+
+.nf
+ II_NEAREST - nearest neighbour
+ II_LINEAR - linear interpolation
+ II_POLY3 - 3rd order divided differences
+ II_POLY5 - 5th order divided differences
+ II_SPLINE3 - cubic spline
+.fi
+
+The calling sequences shall be invariant to the interpolant selected.
+Routines should be designed so that new interpolants can be added
+with minimal changes to the code and no change to the calling sequences.
+.le
+.ls (4)
+The interpolant parameters and the arrays necessary to store the coefficients
+are stored in a structure referenced by a pointer. The pointer is returned
+to the user program by the initial call to ASIINIT or ASIRESTORE and freed
+by a call to ASIFREE (see section 3.1).
+.le
+.ls (5)
+The package routines shall be able to:
+.ls o
+Calculate the coefficients of the interpolant and store these coefficients in
+the appropriate part of the interpolant descriptor structure.
+.le
+.ls o
+Evaluate the interplant at a given x(s) coordinate(s).
+.le
+.ls o
+Calculate the derivatives of the interpolant at a given value of x.
+.le
+.ls o
+Integrate the interpolant over a specified x interval.
+.le
+.ls o
+Store the interpolant in a user supplied array. Restore the saved interpolant
+to the interpolant descriptor structure for later use by ASIEVAL, ASIVECTOR,
+ASIDER and ASIGRL.
+.le
+
+.sh
+3. Specifications
+
+.sh
+3.1. The Array Sequential Interpolator Routines
+
+ The package prefix is asi and the package routines are:
+
+.nf
+ asiinit (asi, interp_type)
+ asifit (asi, datain, npix)
+ y = asieval (asi, x)
+ asivector (asi, x, yfit, npix)
+ asider (asi, x, der, nder)
+ v = asigrl (asi, a, b)
+ asisave (asi, interpolant)
+ asirestore (asi, interpolant)
+ asifree (asi)
+.fi
+
+.sh
+3.2. The Array Random Interpolator Routines
+
+ The package prefix is ari and the package routines are:
+
+.nf
+ y = arieval (x, datain, npix, interp_type)
+ arider (x, datain, npix, der, nder, interp_type)
+.fi
+
+.sh
+3.3. Miscellaneous
+
+ A routine has been included in the package to remove INDEF valued
+pixels from an array.
+
+.nf
+ arbpix (datain, dataout, npix, interp_type, boundary_type)
+.fi
+
+.sh
+3.4. Algorithms
+
+.sh
+3.4.1. Coefficients
+
+ The coefficient array used by the asi routines is calculated by ASIFIT.
+ASIFIT accepts an array of data, checks that the number
+of data points is appropriate for the interpolant selected, allocates
+space for the interpolant, and calculates the coefficients.
+Boundary coefficient values are calculated
+using boundary projection. With the exception of the cubic spline interpolant,
+the coefficients are stored as the data points.
+The B-spline coefficients are
+calculated using natural end conditions (Prenter 1975).
+After a call to ASIFIT the coefficient array contains the following.
+
+.nf
+ case II_NEAREST:
+
+ # no boundary conditions necessary
+ coeff[1] = datain[1]
+ .
+ .
+ .
+ coeff[npts] = datain[npix]
+
+ case II_LINEAR:
+
+ # coeff[npxix+1] required if x = npix
+ coeff[1] = datain[1]
+ .
+ .
+ .
+ coeff[npix] = datain[npix]
+ coeff[npix+1] = 2. * datain[npix] - datain[npix-1]
+
+ case II_POLY3:
+
+ # coeff[0] required if x = 1
+ # coeff[npix+1], coeff[npix+2] required if x = npix
+ coeff[0] = 2. * datain[1] - datain[2]
+ coeff[1] = datain[1]
+ .
+ .
+ .
+ coeff[npix] = datain[npix]
+ coeff[npix+1] = 2. * datain[npix] - datain[npix-1]
+ coeff[npix+2] = 2. * datain[npix] - datain[npix-2]
+
+ case II_POLY5:
+
+ # coeff[1], coeff[0] reqired if x = 1
+ # coeff[npix+1], coeff[npix+2], coeff[npix=3]
+ # required if x = npix
+
+ coeff[-1] = 2. * datain[1] - datain[3]
+ coeff[0] = 2. * datain[1] - datain[2]
+ coeff[1] = datain[1]
+ .
+ .
+ .
+ coeff[npix] = datain[npix]
+ coeff[npix+1] = 2. * datain[npix] - datain[npix-1]
+ coeff[npix+2] = 2. * datain[npix] - datain[npix-2]
+ coeff[npix+3] = 2. * datain[npix] - datain[npix-3]
+
+ case SPLINE3:
+
+ # coeff[0] = 2nd der at x = 1, coeff[0] = 0.
+ # coeff[npix+1] = 2nd der at x = npts, coeff[npix+1] = 0.
+ # coeff[npix+2] = 0., required if x = npix
+ coeff[0] = b[1]
+ coeff[1] = b[2]
+ .
+ .
+ .
+ coeff[npix] = b[npix+1]
+ coeff[npix+1] = b[npix+2]
+ coeff[npix+2] = 0.
+.fi
+
+.sh
+3.4.2. Evaluation
+
+ The ASIEVAL and ASIVECTOR routines have been optimized to be as efficient
+as possible. The values of the II_NEAREST and II_LINEAR interpolants
+are calculated directly. The II_SPLINE3 interpolant is evaluated using
+polynomial coefficients calculated directly from the B-spline coefficients
+(de Boor 1978). Values of the higher order polynomial interpolants
+are calculated using central differences. The equations for each case are
+listed below.
+
+.nf
+case II_NEAREST:
+
+ y = coeff[int (x + 0.5)]
+
+case II_LINEAR:
+
+ nx = x
+ y = (x - nx) * coeff[nx+1] + (nx + 1 - x) * coeff[nx]
+
+case II_POLY3:
+
+ nx = x
+ s = x - nx
+ t = 1. - s
+
+ # second central differences
+ cd20 = 1./6. * (coeff[nx+1] - 2. * coeff[nx] + coeff[nx-1])
+ cd21 = 1./6. * (coeff[nx+2] - 2. * coeff[nx+1] + coeff[nx])
+
+ y = s * (coeff[nx+1] + (s * s - 1.) * cd21) + t * (coeff[nx] +
+ (t * t - 1.) * cd20)
+
+case II_POLY5:
+
+ nx = x
+ s = x - nx
+ t = 1. - s
+
+ # second central differences
+ cd20 = 1./6. * (coeff[nx+1] - 2. * coeff[nx] + coeff[nx-1])
+ cd21 = 1./6. * (coeff[nx+2] - 2. * coeff[nx+1] + coeff[nx])
+
+ # fourth central diffreences
+ cd40 = 1./120. * (coeff[nx-2] - 4. * coeff[nx-1] + 6. * coeff[nx] - 4. *
+ coeff[nx+1] + a[nx+2])
+ cd41 = 1./120. * (coeff[nx-1] - 4. * coeff[nx] + 6. * coeff[nx+1] - 4. *
+ coeff[nx+2] + coeff[nx+3]
+
+ y = s * (coeff[nx+1] + (s * s - 1.) * (cd21 + (s * s - 4.) * cd41)) +
+ t * (coeff[nx] + (t * t - 1.) * (cd20 + (t * t - 4.) * cd40))
+
+case II_SPLINE3:
+
+ nx = x
+ s = x - nx
+
+ pc[1] = coeff[nx-1] + 4. * coeff[nx] + coeff[nx+1]
+ pc[2] = 3. * (coeff[nx+1] - coeff[nx-1])
+ pc[3] = 3. * (coeff[nx-1] - 2. * coeff[nx] + coeff[nx+1])
+ pc[4] = -coeff[nx-1] + 3. * coeff[nx] - 3. * coeff[nx+1] + coeff[nx+2]
+
+ y = pc[1] + s * (pc[2] + s * (pc[3] + s * pc[4]))
+.fi
+
+
+ The ARIEVAL routine uses the expressions above to evaluate the
+interpolant. However unlike ASIEVAL, ARIEVAL does not use a previously
+calculated coefficient array. Instead ARIEVAL selects the appropriate
+portion of the data array, calculates the coefficients and boundary
+coefficients if necessary, and evaluates the interpolant at the time it
+is called. The cubic spline interpolant uses at most SPLTS (currently 16)
+data points to calculate the B-spline coefficients.
+
+.sh
+3.4.3. Derivatives
+
+ Derivatives of the interpolant are calculated by evaluating the
+derivatives of the interpolating polynomial. For all interpolants der[1]
+equals the value of the interpolant at x.
+For the sake of efficiency the derivatives
+of the II_NEAREST and II_LINEAR interpolants are calculated directly.
+
+.nf
+ case II_NEAREST:
+
+ der[1] = coeff[int (x+0.5)]
+
+ case II_LINEAR:
+
+ der[1] = (x - nx) * coeff [nx+1] + (nx + 1 - x) * coeff[nx]
+ der[2] = coeff[nx+1] - coeff[nx]
+.fi
+
+ In order to calculate the derivatives of the cubic spline and
+polynomial interpolants
+the coefficients of the interpolating polynomial must be calculated.
+The polynomial
+coefficients for the cubic spline interpolant are computed directly from the
+B-spline coefficients (see 3.4.2.). The higher order polynomial
+interpolant coefficients are calculated as follows.
+
+First the appropriate portion of the coefficient array is loaded.
+
+.nf
+ do i = 1, nterms
+ d[i] = coeff[nx - nterms/2 + i]
+.fi
+
+Next the divided differences are calculated (Conte and de Boor 1972).
+
+.nf
+ do k = 1, nterms - 1
+ do i = 1, nterms - k
+ d[i] = (d[i+1] - d[i]) / k
+.fi
+
+The d[i] are the coefficients of an interpolating polynomial of the
+following form. The x[i] are the nterms data points surrounding the
+point of interest.
+
+.nf
+ p(x) = d[1] * (x-x[1]) * ... * (x-x[nterms-1) +
+ d[2] * (x-x[2]) * ... * (x-x[nterms-1]) + ... + d[nterms]
+.fi
+
+Next a new set of polynomial coefficients are calculated
+(Conte and de Boor 1972).
+
+.nf
+ do k = nterms, 2, -1
+ do i = 2, k
+ d[i] = d[i] + d[i-1] * (k - i - nterms/2)
+.fi
+
+The new d[i] are the coefficients of the follwoing polynomial.
+
+.nf
+ nx = x
+ p(x) = d[1] * (x-nx) ** (nterms-1) + d[2] * (x-nx) ** (nterms-2) + ...
+ d[nterms]
+.fi
+
+The d[i] array is flipped. The value and derivatives
+of the interpolant are then calculated using the d[i] array and
+nested multiplication.
+
+.nf
+ s = x - nx
+
+ do k = 1, nder {
+
+ accum = d[nterms-k+1]
+
+ do j = nterms - k, 1, -1
+ accum = d[j] + s * accum
+
+ der[k] = accum
+
+ # differnetiate
+ do j = 1, nterms - k
+ d[j] = j * d[j + 1]
+ }
+.fi
+
+ ARIDER calculates the derivatives of the interpolant using the same
+technique ASIDER. However ARIDER does not use a previously calculated
+coefficient array like ASIDER. Instead ARIDER selects the appropriate portion
+of the data array, calculates the coefficients and boundary coefficients,
+and computes the derivatives at the time it is called.
+
+.sh
+3.4.5. Integration
+
+ ASIGRL calculates the integral of the interpolant between fixed limits
+by integrating the interpolating polynomial. The coefficients of the
+interpolating polynomial are calculated as discussed in section 3.4.4.
+
+.sh
+4. Usage
+
+.sh
+4.1. User Notes
+
+The following series of steps illustrates the use of the package.
+
+.ls 4
+.ls (1)
+Insert an include <iminterp.h> statement in the calling program to make
+the IINTERP definitions available to the user program.
+.le
+.ls (2)
+Remove INDEF valued pixels from the data using ARBPIX.
+.le
+.ls (3)
+Call ASIINIT to initialize the interpolant parameters.
+.le
+.ls (4)
+Call ASIFIT to calculate the coefficients of the interpolant.
+.le
+.ls (5)
+Evaluate the interpolant at a given value of x(s) using ASIEVAL or
+ASIVECTOR.
+.le
+.ls (6)
+Calculate the derivatives and integral or the interpolant using
+ASIDER and ASIGRL.
+.le
+.ls (7)
+Free the interpolator structure by calling ASIFREE.
+.le
+.le
+
+ The interpolant can be saved and restored using ASISAVE and ASIRESTORE.
+If the values and derivatives of only a few points in an array are desired
+ARIEVAL and ARIDER can be called.
+
+.sh
+4.2. Examples
+
+.nf
+Example 1: Shift a data array by a constant amount
+
+ include <iminterp.h>
+ ...
+ call asiinit (asi, II_POLY5)
+ call asifit (asi, inrow, npix)
+
+ do i = 1, npix
+ outrow[i] = asieval (asi, i + shift)
+
+ call asifree (asi)
+ ...
+
+Example 2: Calculate the integral under the data array
+
+ include <iminterp.h>
+ ...
+ call asiinit (asi, II_POLY5)
+ call asifit (asi, datain, npix)
+
+ integral = asigrl (asi, 1. real (npix))
+
+ call asifree (asi)
+ ...
+
+Example 2: Store interpolant for later use by ASIEVAL
+ LEN_INTERP must be at least npix + 8 units long where npix is
+ is defined in the call to ASIFIT.
+
+ include <iminterp.h>
+
+ real interpolant[LEN_INTERP]
+ ...
+ call asiinit (asi, II_POLY3)
+ call asifit (asi, datain, npix)
+ call asisave (asi, interpolant)
+ call asifree (asi)
+ ...
+ call asirestore (asi, interpolant)
+ do i = 1, npts
+ yfit[i] = asieval (asi, x[i])
+ call asifree (asi)
+ ...
+.fi
+.sh
+5. Detailed Design
+
+.sh
+5.1. Interpolator Descriptor
+
+ The interpolant parameters and coefficients are stored in a
+structure listed below.
+
+.nf
+ define LEN_ASISTRUCT 4 # Length in structure units of
+ # interpolant descriptor
+
+ define ASI_TYPE Memi[$1] # Interpolant type
+ define ASI_NCOEFF Memi[$1+1] # No. of coefficients
+ define ASI_OFFSET Memi[$1+2] # First "data" point in
+ # coefficient array
+ define ASI_COEFF Memi[$1+3] # Pointer to coefficient array
+.fi
+
+.sh
+5.2. Storage Requirements
+
+ The interpolant descriptor requires LEN_ASISTRUCT storage units. The
+coefficient array requires ASI_NCOEFF(asi) real storage elements, where
+ASI_NCOEFF(asi) is defined as follows.
+
+.nf
+ ASI_NCOEFF(asi) = npix # II_NEAREST
+ ASI_NCOEFF(asi) = npix+1 # II_LINEAR
+ ASI_NCOEFF(asi) = npix+3 # II_POLY3
+ ASI_NCOEFF(asi) = npix+5 # II_POLY5
+ ASI_NCOEFF(asi) = npix+3 # II_SPLINE3
+.fi
+
+.sh
+6. References
+
+.ls (1)
+Carl de Boor, "A Practical Guide to Splines", 1978, Springer-Verlag New
+York Inc.
+.le
+.ls (2)
+S.D. Conte and C. de Boor, "Elementary Numerical Analysis", 1972, McGraw-Hill,
+Inc.
+.le
+.ls (3)
+P.M. Prenter, "Splines and Variational Methods", 1975, John Wiley and Sons Inc.
+.le
+.endhelp
diff --git a/math/iminterp/doc/mrider.hlp b/math/iminterp/doc/mrider.hlp
new file mode 100644
index 00000000..47515c48
--- /dev/null
+++ b/math/iminterp/doc/mrider.hlp
@@ -0,0 +1,79 @@
+.help mrider Dec98 "Image Interpolation Package"
+.ih
+NAME
+mrider -- calculate the derivatives at x and y
+.ih
+SYNOPSIS
+include <math/iminterp.h>
+
+.nf
+mrider (x, y, datain, nxpix, nypix, len_datain, der, nxder, nyder, len_der,
+ interp_type)
+.fi
+
+.nf
+real x[4] #I x value, 1. <= x[1-4] <= nxpix
+real y[4] #I y value, 1. <= y[1-4] <= nypix
+real datain[len_datain, ARB] #I data array
+int nxpix #I number of data pixels in x
+int nypix #I number of data pixels in y
+int len_datain #I length of datain, len_datain >= nxpix
+real der[len_der, ARB] #O derivative array
+int nxder #I x order of the derivatives
+int nyder #I y order of the derivatives
+int len_der #I row length of der, len_der >= nxder
+int interp_type #I interpolant type
+.fi
+.ih
+ARGUMENTS
+.ls x, y
+The single x and y points or in the case of the drizzle interpolant the
+single quadrilateral at / over which the derivatives are to be evaluated.
+The quadrilateral vertices may be stored in clock-wise or counter-clockwise
+order.
+.le
+.ls datain
+Array of data values.
+.le
+.ls nxpix, nypix
+The number of data values in the x and y directions
+.le
+.ls len_datain
+The row length of the datain array. Len_datain must be >= nxpix.
+.le
+.ls der
+The derivative array. Der[1,1] equals the function value at x and y and
+der[2,1], der[1,2] are the first derivatives with respect to x and y
+respectively.
+.le
+.ls nxder, nyder
+The number of the derivatives in x and y to be returned. MRIDER checks
+that the requested number of derivatives is sensible. The sinc interpolants
+return the interpolant value and all the first and second order derivatives.
+The drizzle interpolant returns the interpolant value and the first
+derivative in x and y.
+.le
+.ls len_der
+The row length of the derivative array. Len_der must be >= nxder.
+.le
+.ls interp_type
+Interpolant type. The options are II_BINEAREST, II_BILINEAR, II_BIPOLY3,
+II_BIPOLY5, II_BISPLINE3, II_SINC / II_LSINC, and II_DRIZZLE. The look-up
+table sinc is not supported and defaults to the sinc interpolant. The
+interpolant width is 31 pixels. The drizzle pixel fraction is 1.0. The
+interpolant type definitions are found in the package header file
+math/iminterp.h.
+.le
+.ih
+DESCRIPTION
+MRIDER is useful for evaluating the function and derivatives at a few
+widely spaced points in a data array without the storage space required
+by the sequential version.
+.ih
+NOTES
+Checking for out of bounds and INDEF valued pixels is the
+responsibility of the user.
+.ih
+SEE ALSO
+msider
+.endhelp
diff --git a/math/iminterp/doc/mrieval.hlp b/math/iminterp/doc/mrieval.hlp
new file mode 100644
index 00000000..35477614
--- /dev/null
+++ b/math/iminterp/doc/mrieval.hlp
@@ -0,0 +1,57 @@
+.help mrieval Dec98 "Image Interpolation Package"
+.ih
+NAME
+mrieval -- evaluate the interpolant at x and y
+.ih
+SYNOPSIS
+include <math/iminterp.h>
+
+y = mrieval (x, y, datain, nxpix, nypix, len_datain, interp_type)
+
+.nf
+ real x[4] #I x value, 1 <= x[1-4] <= nxpix
+ real y[4] #I y value, 1 <= y[1-4] <= nypix
+ real datain[len_datain, ARB] #I data array
+ int nxpix #I number of x values
+ int nypix #I number of y values
+ int len_datain #I length datain, len_datain >= nxpix
+ int interp_type #I interpolant type
+.fi
+.ih
+ARGUMENTS
+.ls x, y
+The single x and y values or in the case of the drizzle interpolant the
+single quadrilateral at / over which the interpolant is to be evaluated.
+The vertices of the quadilateral must be defined in clock-wise or
+counter-clockwise order.
+.le
+.ls datain
+The array of data values.
+.le
+.ls nxpix, nypix
+The number of data pixels in x and y.
+.le
+.ls len_datain
+The row length of datain. Len_datain must be >= nxpix.
+.le
+.ls interp_type
+Interpolant type. The options are II_BINEAREST, II_BILINEAR, II_BIPOLY3,
+II_BIPOLY5, II_BISPLINE3, II_SINC / II_LSINC, and II_DRIZZLE. The look-up
+table sinc interpolant is not supported and defaults to the sinc interpolant.
+The sinc interpolant width is 31 pixels. The drizzle pixel fraction is 1.0.
+The interpolant type definitions are found in the package header file
+math/iminterp.h.
+.le
+.ih
+DESCRIPTION
+MRIEVAL is useful for evaluating the interpolant at a few selected points
+in the datain array without the storage overhead required for the sequential
+version.
+.ih
+NOTES
+Checking for INDEF valued or out of bounds pixels is the
+responsibility of the user.
+.ih
+SEE ALSO
+msieval, msivector, mrider
+.endhelp
diff --git a/math/iminterp/doc/msider.hlp b/math/iminterp/doc/msider.hlp
new file mode 100644
index 00000000..0139c0a0
--- /dev/null
+++ b/math/iminterp/doc/msider.hlp
@@ -0,0 +1,52 @@
+.help msider Dec98 "Image Interpolation Package"
+.ih
+NAME
+msider -- evaluate the interpolant derivatives at x and y
+.ih
+SYNOPSIS
+msider (msi, x, y, der, nxder, nyder, len_der)
+
+.nf
+ pointer msi #I interpolant descriptor
+ real x[4] #I x value, 1 <= x[1-4] <= nxpix
+ real y[4] #I y value, 1 <= y[1-4] <= nypix
+ real der[len_der, ARB] #O derivative array
+ int nxder #I number of x derivatives
+ int nyder #I number of y derivatives
+ int len_der #I row length of der, len_der >= nxder
+.fi
+.ih
+ARGUMENTS
+.ls msi
+Pointer to the 2D sequential interpolant descriptor.
+.le
+.ls x, y
+The single x and y values or in the case of the drizzle interpolant the
+single quadrilateral at / over which the point is to be evaluated.
+.le
+.ls der
+The array containing the derivatives. Der[1,1] contains the value of
+the interpolant at x and y. Der[2,1] and der[1,2] contain the 1st
+derivatives of x and y respectively.
+.le
+.ls nxder, nyder
+The number derivatives in x and y.
+.le
+.ls len_der
+The row length of der. Len_der must be >= nxder.
+.le
+.ih
+DESCRIPTION
+The polynomial and spline interpolants are evaluated using the polynomial
+coefficients and nested multiplication. The polynomial interpolants are
+stored as the data points. The spline interpolant is stored as a set of
+B-spline coefficients.
+.ih
+NOTES
+MRIDER checks that the number of derivatives requested is reasonable.
+Checking for out of bounds and INDEF valued pixels is the responsibility of the
+user. MSIINIT and MSIFIT must be called before using MSIDER.
+.ih
+SEE ALSO
+mrider
+.endhelp
diff --git a/math/iminterp/doc/msieval.hlp b/math/iminterp/doc/msieval.hlp
new file mode 100644
index 00000000..9a77c006
--- /dev/null
+++ b/math/iminterp/doc/msieval.hlp
@@ -0,0 +1,46 @@
+.help msieval Dec98 "Image Interpolation Package"
+.ih
+NAME
+msieval -- procedure to evaluate the interpolant at x and y
+.ih
+SYNOPSIS
+z = msieval (msi, x, y)
+
+.nf
+pointer msi #I interpolant descriptor
+real x[4] #I x value, 1 <= x[1-4] <= nxpix
+real y[4] #I y value, 1 <= y[1-4] <= nypix
+.fi
+.ih
+ARGUMENTS
+.ls msi
+The pointer to the sequential interpolant descriptor structure.
+.le
+.ls x, y
+The single x and y values of or in the case of the drizzle interpolant
+the single quadrilateral over which the point is to be evaluated.
+.le
+.ih
+DESCRIPTION
+The polynomial coefficients are calculated from the data points in the
+case of the polynomial interpolants and the B-spline coefficients in
+the case of the spline interpolant. The polynomial interpolants
+are evaluated using Everett's central difference formula. The boundary
+extension algorithm is projection.
+
+The sinc interpolant is evaluated using an array of data points around
+the desired position. The look-up table sinc interpolant is computed
+using an a pre-computed look--up table entry. The boundary extension
+algorithm is nerest neighbor.
+
+The drizzle interpolant is computed by computing the mean value of the
+data within the user supplied quadrilateral.
+.ih
+NOTES
+Checking for out of bounds and INDEF valued pixels is the responsibility of
+the user. MSIINIT or MSISINIT and MSIFIT must be called before calling
+MSIEVAL.
+.ih
+SEE ALSO
+msivector, mrieval, mrider
+.endhelp
diff --git a/math/iminterp/doc/msifit.hlp b/math/iminterp/doc/msifit.hlp
new file mode 100644
index 00000000..9a63ae2d
--- /dev/null
+++ b/math/iminterp/doc/msifit.hlp
@@ -0,0 +1,45 @@
+.help msifit Dec98 "Image Interpolation Package"
+.ih
+NAME
+msifit - fit the interpolant to the data
+.ih
+SYNOPSIS
+msifit (msi, datain, nxpix, nypix, len_datain)
+
+.nf
+ pointer msi #I interpolant descriptor
+ real datain[len_datain,ARB] #I data array
+ int nxpix #I number of x pixels
+ int nypix #I number of y pixels
+ int len_datain #I length of datain, len_datain >= nxpix
+.fi
+.ih
+ARGUMENTS
+.ls msi
+Pointer to the sequential interpolant descriptor.
+.le
+.ls datain
+Array containing the data.
+.le
+.ls nxpix, nypix
+The number of pixels in x and y.
+.le
+.ls len_datain
+The row length of the datain array. Len_datain must be >= nxpix.
+.le
+.ih
+DESCRIPTION
+The datain array is checked for size, memory is allocated for the coefficient
+array and the end conditions are specified. The interior polynomial, sinc,
+and drizzle interpolants are saved as the data points. The polynomial
+coefficients are calculated from the data points in the evaluation stage.
+The B-spline coefficients are calculated in MSIFIT as they depend on the
+entire data array.
+.ih
+NOTES
+Checking for INDEF valued pixels is the responsibility of the user.
+MSIINIT or MSISINIT must be called before using MSIFIT. MSIFIT must be
+called before using MSIEVAL, MSIVECTOR, MSIDER, MSIGRL or MSISQGRL.
+.ih
+SEE ALSO
+.endhelp
diff --git a/math/iminterp/doc/msifree.hlp b/math/iminterp/doc/msifree.hlp
new file mode 100644
index 00000000..79b1966d
--- /dev/null
+++ b/math/iminterp/doc/msifree.hlp
@@ -0,0 +1,26 @@
+.help msifree Dec98 "Image Interpolation Package"
+.ih
+NAME
+msifree -- free sequential interpolant descriptor
+.ih
+SYNOPSIS
+msifree (msi)
+
+.nf
+ pointer msi #U interpolant descriptor
+.fi
+.ih
+ARGUMENTS
+.ls msi
+Pointer to the sequential interpolant descriptor structure.
+.le
+.ih
+DESCRIPTION
+MSIFREE frees the sequential interpolant descriptor structure.
+MSIFREE should be called when interpolation is complete.
+.ih
+NOTES
+.ih
+SEE ALSO
+msiinit, msisinit
+.endhelp
diff --git a/math/iminterp/doc/msigeti.hlp b/math/iminterp/doc/msigeti.hlp
new file mode 100644
index 00000000..5ab46156
--- /dev/null
+++ b/math/iminterp/doc/msigeti.hlp
@@ -0,0 +1,35 @@
+.help msigeti Dec98 msigeti.hlp
+.ih
+NAME
+msigeti -- fetch an msi integer parameter
+.ih
+SYNOPSIS
+include <math/iminterp.h>
+
+ivalue = msigeti (msi, param)
+
+.nf
+ pointer msi #I interpolant descriptor
+ int param #I parameter
+.fi
+.ih
+ARGUMENTS
+.ls msi
+Pointer to the sequential interpolant descriptor structure.
+.le
+.ls param
+The parameter to be fetched. The choices are: II_MSITYPE, the interpolant
+type, II_MSINSAVE, the length of the saved coefficient array, and
+II_MSINSINC, the half-width of the sinc interpolant.
+.le
+.ih
+DESCRIPTION
+MSIGETI is used to determine the size of the coefficient array that
+must be allocated to save the sequential interpolant description structure,
+and to fetch selected sequential interpolant parameters.
+.ih
+NOTES
+.ih
+SEE ALSO
+msiinit, msisinit, msigetr
+.endhelp
diff --git a/math/iminterp/doc/msigetr.hlp b/math/iminterp/doc/msigetr.hlp
new file mode 100644
index 00000000..dadac5c2
--- /dev/null
+++ b/math/iminterp/doc/msigetr.hlp
@@ -0,0 +1,37 @@
+.help msigetr Dec98 msigetr.hlp
+.ih
+NAME
+msigetr -- fetch an msi real parameter
+.ih
+SYNOPSIS
+include <math/iminterp.h>
+
+rvalue = msigetr (msi, param)
+
+.nf
+ pointer msi #I interpolant descriptor
+ int param #I parameter
+.fi
+.ih
+ARGUMENTS
+.ls msi
+Pointer to the sequential interpolant descriptor structure.
+.le
+.ls param
+The parameter to be fetched. The choices are: II_MSIBADVAL, the undefined
+pixel value for the drizzle interpolant. The parameter definitions are
+contained in the package header file math/iminterp.h.
+.le
+.ih
+DESCRIPTION
+MSIGETR is used to set the value of undefined drizzle interpolant pixels.
+Undefined pixels are those for which the interpolation coordinates do not
+overlap the input coordinates, but are still, within the boundaries of the input
+image, a situation which may occur when the pixel fraction is < 1.0.
+.ih
+.ih
+NOTES
+.ih
+SEE ALSO
+msiinit, msisinit, msigeti
+.endhelp
diff --git a/math/iminterp/doc/msigrid.hlp b/math/iminterp/doc/msigrid.hlp
new file mode 100644
index 00000000..6bc110d5
--- /dev/null
+++ b/math/iminterp/doc/msigrid.hlp
@@ -0,0 +1,51 @@
+.help msigrid Dec98 "Image Interpolation Package"
+.ih
+NAME
+msigrid -- evaluate the interpolant on a grid of points
+.ih
+SYNOPSIS
+msigrid (msi, x, y, zfit, nx, ny, len_zfit)
+
+.nf
+ pointer msi #I interpolant descriptor
+ real x[2*nx] #I x values, 1 <= x[i] <= nx
+ real y[2*ny] #I y values, 1 <= y[i] <= ny
+ real zfit[len_zfit,ARB] #O grid of interpolated values
+ int nx #I number of x points
+ int ny #I number of y points
+ int len_zfit #I length zfit, len_zfit >= nx
+.fi
+.ih
+ARGUMENTS
+.ls msi
+Pointer to the interpolant descriptor structure.
+.le
+.ls x, y
+The x and y arrays of points to be evaluated, or in the case of the drizzle
+interpolant the x and y ranges over which the points are to be evaluated.
+The x and y arrays must be ordered in increasing values of x and y respectively.
+.le
+.ls zfit
+The array of interpolated points.
+.le
+.ls nx, ny
+The number of points in the x and y directions respectively.
+.le
+.ls len_zfit
+The row length of the zfit array. Len_zfit >= nx.
+.le
+.ih
+DESCRIPTION
+MSIGRID evaluates the interpolant at a set of x and y values on a
+rectangular grid or in the case of the drizzle interpolant within
+rectangular regions. It is most efficient for evaluating the interpolant
+at many values which are closely spaced in x and y. For widely spaced
+points MSIVECTOR should be used.
+.ih
+NOTES
+Checking for out of bounds and INDEF valued pixels is the responsibility
+of the user.
+.ih
+SEE ALSO
+msieval, msivector, msider, mrieval, mrider
+.endhelp
diff --git a/math/iminterp/doc/msigrl.hlp b/math/iminterp/doc/msigrl.hlp
new file mode 100644
index 00000000..f6e2326a
--- /dev/null
+++ b/math/iminterp/doc/msigrl.hlp
@@ -0,0 +1,43 @@
+.help msigrl Dec98 "Image Interpolation Package"
+.ih
+NAME
+msigrl -- integrate the interpolant inside a polygon
+.ih
+SYNOPSIS
+y = msigrl (msi, x, y, npts)
+
+.nf
+ pointer msi #I interpolant descriptor
+ real x[npts] #I x values, 1 <= x <= npts, x[1] = x[npts]
+ real y[npts] #I y values, 1 <= y <= npts, y[1] = y[npts]
+ int npts #I number of points
+.fi
+.ih
+ARGUMENTS
+.ls msi
+Pointer to the sequential interpolant descriptor structure.
+.le
+.ls x, y
+An array of x and y values describing a polygon, where x[1] = x[npts] and
+y[1] = y[npts]. X and y describe a closed curve where any horizontal line
+segment intersects the domain of integration at at most one point.
+.le
+.ls npts
+The number of points describing the polygon. Npts must >= 4 (triangle).
+.le
+.ih
+DESCRIPTION
+MSIGRL integrates the interpolant exactly for rectangular domains
+of integration. For more irregular regions of integration MSIGRL
+returns an approximation whose accuracy depends on the size of the
+integration region and the shape of the polygon.
+.ih
+NOTES
+Checking for out of bound integration regimes is the responsibility of
+the user. Non-rectangular partial pixel domains of integration default
+to rectangular regions. MSIINIT or MSISINIT and MSIFIT must be called
+before using MSIGRL.
+.ih
+SEE ALSO
+msisqrgl
+.endhelp
diff --git a/math/iminterp/doc/msiinit.hlp b/math/iminterp/doc/msiinit.hlp
new file mode 100644
index 00000000..68dba684
--- /dev/null
+++ b/math/iminterp/doc/msiinit.hlp
@@ -0,0 +1,41 @@
+.help msiinit Dec98 "Image Interpolation Package"
+.ih
+NAME
+msiinit -- initialize the sequential interpolant descriptor
+.ih
+SYNOPSIS
+include <math/iminterp.h>
+
+msiinit (msi, interp_type)
+
+.nf
+ pointer msi #U interpolant descriptor
+ int interp_type #I interpolant type
+.fi
+.ih
+ARGUMENTS
+.ls msi
+Pointer to the sequential interpolant descriptor.
+.le
+.ls interp_type
+Interpolant type. The options are II_BINEAREST, II_BILINEAR, II_BIPOLY3,
+II_BIPOLY5, II_BISPLINE3, II_BISINC, II_BILSINC, and II_BIDRIZZLE, for
+nearest neighbour, bilinear, 3rd and 5th order interior polynomials, bicubic
+spline, sinc, look-up table sinc, and drizzle respectively. The interpolant
+definitions are found in the package header file math/iminterp.h.
+.le
+.ih
+DESCRIPTION
+The interpolant type is allocated and initialized. The pointer msi is
+returned by MSIINIT. The sinc interpolant width defaults to 31 pixels
+in x and y. The look-up table sinc resolution defaults to 20 resolution
+elements or 0.05 pixels in x and y. The drizzle pixel fraction defaults
+to 1.0.
+.ih
+NOTES
+MSIINIT, MSISINIT or MSIRESTORE must be called before using any other
+MSI routines.
+.ih
+SEE ALSO
+msirestore, msifree
+.endhelp
diff --git a/math/iminterp/doc/msirestore.hlp b/math/iminterp/doc/msirestore.hlp
new file mode 100644
index 00000000..664c9b50
--- /dev/null
+++ b/math/iminterp/doc/msirestore.hlp
@@ -0,0 +1,36 @@
+.help msirestore Dec98 "Image Interpolator Package"
+.ih
+NAME
+msirestore -- restore interpolant
+.ih
+SYNOPSIS
+msirestore (msi, interpolant)
+
+.nf
+ pointer msi #U interpolant descriptor
+ real interpolant[] #I array containing interpolant
+.fi
+.ih
+ARGUMENTS
+.ls msi
+Pointer to the interpolant descriptor structure.
+.le
+.ls interpolant
+Array containing the interpolant. The amount of space required by interpolant
+can be determined by a call to msigeti.
+.le
+
+.nf
+ len_interpolant = msigeti (msi, II_MSINSAVE)
+.fi
+.ih
+DESCRIPTION
+MSIRESTORE allocates space for the interpolant descriptor and restores the
+parameters and coefficients stored in the interpolant array to the
+interpolant structure for use by MSIEVAL, MSIVECTOR, MSIDER and MSIGRL.
+.ih
+NOTES
+.ih
+SEE ALSO
+msisave
+.endhelp
diff --git a/math/iminterp/doc/msisave.hlp b/math/iminterp/doc/msisave.hlp
new file mode 100644
index 00000000..7d5f67ef
--- /dev/null
+++ b/math/iminterp/doc/msisave.hlp
@@ -0,0 +1,38 @@
+.help msisave Dec98 "Image Interpolator Package"
+.ih
+NAME
+msisave -- save interpolant
+.ih
+SYNOPSIS
+msisave (msi, interpolant)
+
+.nf
+ pointer msi #I interpolant descriptor
+ real interpolant[] #O array containing the interpolant
+.fi
+.ih
+ARGUMENTS
+.ls msi
+Pointer to the interpolant descriptor structure.
+.le
+.ls interpolant
+Array where the interpolant is stored. The required interpolant array length
+required can be determined by a call to msigeti.
+.le
+
+.nf
+ len_interpolant = msigeti (msi, II_MSINSAVE)
+.fi
+.ih
+DESCRIPTION
+The interpolant type, number of coefficients in x and y, the position of
+the first data point in the coefficient array, and the sinc and drizzle
+interpolant parameters are stored in the first eleven elements of interpolant.
+The remaining elements contain the coefficients and look-up tables
+calculated by MSIFIT.
+.ih
+NOTES
+.ih
+SEE ALSO
+msirestore
+.endhelp
diff --git a/math/iminterp/doc/msisinit.hlp b/math/iminterp/doc/msisinit.hlp
new file mode 100644
index 00000000..0a05e7ab
--- /dev/null
+++ b/math/iminterp/doc/msisinit.hlp
@@ -0,0 +1,61 @@
+.help msisinit Dec98 "Image Interpolator Package"
+.ih
+NAME
+msisinit -- initialize the sequential interpolant descriptor using user parameters
+.ih
+SYNOPSIS
+include <math/iminterp.h>
+
+msisinit (msi, interp_type, nsinc, nxincr, nyincr, pixfrac1, pixfrac2, badval)
+
+.nf
+ pointer msi #O interpolant descriptor
+ int interp_type #I interpolant type
+ int nsinc #I sinc interpolant width in pixels
+ int nxincr,nyincr #I sinc look-up table resolution
+ real pixfrac1,pixfrac2 #I sinc or drizzle pixel fractions
+ real badval #I drizzle undefined pixel value
+.fi
+
+.ih
+ARGUMENTS
+.ls msi
+Pointer to sequential interpolant descriptor.
+.le
+.ls interp_type
+Interpolant type. The options are II_BINEAREST, II_BILINEAR, II_BIPOLY3,
+II_BIPOLY5, II_BISPLINE3, II_BISINC, II_BILSINC, and II_BIDRIZZLE for the
+nearest neighbour, linear, 3rd order polynomial, 5th order polynomial,
+cubic spline, sinc, look-up table sinc, and drizzle interpolants respectively.
+The interpolant type definitions are found in the package header file
+math/iminterp.h.
+.le
+.ls nsinc
+The sinc and look-up table sinc interpolant width in pixels. Nsinc is
+rounded up internally to the nearest odd number.
+.le
+.ls nxincr, nyincr
+The look-up table sinc resolution in x and y in number of entries. Nxincr = 10
+implies a pixel resolution of 0.1 pixels in x, nxincr = 20 a pixel resolution
+of 0.05 pixels in x, etc. The default value of nxincr and nyincr are 20 and 20
+.le
+.ls pixfrac1, pixfrac2
+The look-up table sinc fractional pixel shifts in x and y if nincr = 1 in
+which case -0.5 <= rparam[1/2] <= 0.5 , or the drizzle pixel fractions in
+which case 0.0 <= rparam[1/2] <= 1.0.
+.le
+.ls badval
+The undefined pixel value for the drizzle interpolant. Pixels within
+the boundaries of the input image which do not overlap the input image
+pixels are assigned a value of badval.
+.le
+.ih
+DESCRIPTION
+The interpolant descriptor is allocated and initialized. The pointer msi is
+returned by MSISINIT.
+.ih
+NOTES
+MSIINIT or MSISINIT must be called before using any other MSI routines.
+.ih
+SEE ALSO
+msisinit, msifree
diff --git a/math/iminterp/doc/msisqgrl.hlp b/math/iminterp/doc/msisqgrl.hlp
new file mode 100644
index 00000000..fc49514d
--- /dev/null
+++ b/math/iminterp/doc/msisqgrl.hlp
@@ -0,0 +1,38 @@
+.help msisqgrl Dec98 "Image Interpolation Package"
+.ih
+NAME
+msisqgrl -- integrate the interpolant over a rectangular region
+.ih
+SYNOPSIS
+y = msisqgrl (msi, x1, x2, y1, y2)
+
+.nf
+ pointer msi #I interpolant descriptor
+ real x1 #I lower x limit, 1 <= x1 <= nxpix
+ real x2 #I upper x limit, 1 <= x2 <= nxpix
+ real y1 #I lower y limit, 1 <= y1 <= nypix
+ real y2 #I upper y limit, 1 <= y2 <= nypix
+.fi
+.ih
+ARGUMENTS
+.ls msi
+Pointer to the sequential interpolant descriptor structure.
+.le
+.ls x1, x2
+The x limits of integration
+.le
+.ls y1, y2
+The y limits of integration.
+.le
+.ih
+DESCRIPTION
+MSISQGRL integrates the interpolant exactly for rectangular domains
+of integration, including partial pixel regions.
+.ih
+NOTES
+Checking for out of bound integration regimes is the responsibility of
+the user. MSIFIT must be called before using MSISQGRL.
+.ih
+SEE ALSO
+msigrl
+.endhelp
diff --git a/math/iminterp/doc/msitype.hlp b/math/iminterp/doc/msitype.hlp
new file mode 100644
index 00000000..5da6dc6d
--- /dev/null
+++ b/math/iminterp/doc/msitype.hlp
@@ -0,0 +1,95 @@
+.help msitype Dec98 "Image Interpolator Package"
+.ih
+NAME
+msitype -- decode an interpolant string
+.ih
+SYNOPSIS
+include <math/iminterp.h>
+
+msitype (interpstr, interp_type, nsinc, nincr, pixfrac)
+
+.nf
+ char interpstr #I the input interpolant string
+ int interp_type #O interpolant type
+ int nsinc #O sinc interpolant width in pixels
+ int nincr #O sinc look-up table resolution
+ real pixfrac #O sinc or drizzle pixel fraction
+.fi
+
+.ih
+ARGUMENTS
+.ls interpstr
+The user supplied interpolant string to be decoded. The options are
+.ls nearest
+Nearest neighbor interpolation.
+.le
+.ls linear
+Bilinear interpolation
+.le
+.ls poly3
+Bicubic polynomial interpolation.
+.le
+.ls poly5
+Biquintic polynomial interpolation.
+.le
+.ls spline3
+Bicubic spline interpolation.
+.le
+.ls sinc
+2D sinc interpolation. Users can specify the sinc interpolant width by
+appending a width value to the interpolant string, e.g. sinc51 specifies
+a 51 by 51 pixel wide sinc interpolant. The sinc width will be rounded up to
+the nearest odd number. The default sinc width is 31 by 31.
+.le
+.ls lsinc
+Look-up table sinc interpolation. Users can specify the look-up table sinc
+interpolant width by appending a width value to the interpolant string, e.g.
+lsinc51 specifies a 51 by 51 pixel wide look-up table sinc interpolant. The user
+supplied sinc width will be rounded up to the nearest odd number. The default
+sinc width is 31 by 31 pixels. Users can specify the resolution of the lookup
+table sinc by appending the look-up table size in square brackets to the
+interpolant string, e.g. lsinc51[20] specifies a 20 by 20 element sinc
+look-up table interpolant with a pixel resolution of 0.05 pixels in x and y.
+The default look-up table size and resolution are 20 by 20 and 0.05 pixels
+in x and y respectively.
+.le
+.ls drizzle
+Drizzle interpolation. Users can specify the drizzle pixel fraction by
+appending the pixel fraction value to the interpolant string in square
+brackets, e.g. drizzle[0.5] specifies a pixel fraction of 0.5 in x and y.
+The default pixel fraction is 1.0. If either of the x or y pixel
+fractions are 0.0, then both are set to 0.0. A minimum value of 0.001
+is imposed on the actual value of pixfrac.
+.le
+.le
+.ls interp_type
+The output interpolant type. The options are II_BINEAREST, II_BILINEAR,
+II_BIPOLY3, II_BIPOLY5, II_BISPLINE3, II_BISINC, II_BILSINC, and II_BIDRIZZLE
+for the nearest neighbor, linear, 3rd order polynomial, 5th order polynomial,
+cubic spline, sinc, look-up table sinc, and drizzle interpolants respectively.
+The interpolant type definitions are found in the package header file
+math/iminterp.h.
+.le
+.ls nsinc
+The output sinc and look-up table sinc interpolant width in pixels. The
+default value is 31 pixels in x and y.
+.le
+.ls nincr
+The output sinc look-up table size. Nincr = 10 implies a pixel resolution
+of 0.1 pixels in x, nincr = 20 a pixel resolution of 0.05 pixels in y, etc. The
+default value of nincr is 20.
+.le
+.ls pixfrac
+The output look-up table sinc fractional pixel shift if nincr = 1
+in which case -0.5 <= pixfrac <= 0.5 , or the drizzle pixel
+fraction in which case 0.0 <= pixfrac <= 1.0.
+.le
+.ih
+DESCRIPTION
+The interpolant string is decoded into values suitable for the MSISINIT
+or MSIINIT routines.
+.ih
+NOTES
+.ih
+SEE ALSO
+msinit, msisinit, msifree
diff --git a/math/iminterp/doc/msivector.hlp b/math/iminterp/doc/msivector.hlp
new file mode 100644
index 00000000..d6561be7
--- /dev/null
+++ b/math/iminterp/doc/msivector.hlp
@@ -0,0 +1,54 @@
+.help msivector Dec98 "Image Interpolation Package"
+.ih
+NAME
+msivector -- evaluate the interpolant at an array of x and y points
+.ih
+SYNOPSIS
+msivector (msi, x, y, zfit, npts)
+
+.nf
+ pointer msi #I interpolant descriptor
+ real x[npts/4*npts] #I x values, 1 <= x <= nxpix
+ real y[npts/4*npts] #I y values, 1 <= y <= nypix
+ real zfit[npts] #O interpolated values
+ int npts #I number of points
+.fi
+.ih
+ARGUMENTS
+.ls msi
+The pointer to the sequential interpolant descriptor
+.le
+.ls x, y
+The array of x and y values at or in the case tof the drizzle interpolant the
+array of quadrilaterals over which to evaluate the interpolant.
+.le
+.ls zfit
+The interpolated values.
+.le
+.ls npts
+The number of points.
+.le
+.ih
+DESCRIPTION
+The polynomial coefficients are calculated directly from the data points,
+The polynomial interpolants are evaluated using Everett's central difference
+formula. The spline interpolant uses the B-spline coefficients
+calculated using the MSIFIT routine. The boundary extension algorithm is
+projection.
+
+The sinc interpolant is evaluated using a array of data points around
+the point in question. The look-up table since is computed by convolving
+the data with a pre-computed look-up table entry. The boundary extension
+algorithm is nearest neighbor.
+
+The drizzle interpolant is evaluated by summing the data over the
+list of user supplied quadrilaterals.
+.ih
+NOTES
+Checking for out of bounds and INDEF valued pixels is the responsibility
+of the user. MSIINIT or MSISINIT and MSIFIT must be called before using
+MSIVECTOR.
+.ih
+SEE ALSO
+msieval, mrieval
+.endhelp