Bug in INTERPOLATE(/CUBIC) [message #5073] |
Sun, 24 September 1995 00:00 |
landsman
Messages: 93 Registered: August 1991
|
Member |
|
|
I believe that there is a bug in the output of the intrinsic INTERPOLATE
function when used with the /CUBIC keyword. (I have found this problem
in IDL V3.6 and V4.0.1 under both Alpha VMS and SunOS.)
What bothered me was that the use of the /CUBIC keyword never gave as accurate
interpolation as either sinc interpolation or my program quadterp.pro
(available in ftp://idlastro.gsfc.nasa.gov/pub/pro/math), which also does cubic
interpolation using 4 neighboring points. In fact, the use of the /CUBIC
keyword gave results no better than linear interpolation.
I then found that if I *averaged* the results of INTERPOLATE(P,X) and
INTERPOLATE(P,X,/CUBIC), then I would *exactly* match the output of
quadterp.pro. So my guess is that a step is missing from the internal
algorithm for INTERPOLATE(/CUBIC).
I suspect the problem also exists for 2-d interpolation, although I have
done only a quick test of this.
As example of the problem, below I interpolate the function y = exp(x),
first at the grid points x = findgen(20), and then at the point x= 6.5.
IDL> x = findgen(20)
IDL> y = exp(x)
IDL> print,exp(6.5) ;Print the true value at X = 6.5
; 665.142
IDL> print,interpolate(y,6.5) ;Linear interpolation gives a value too large
; 750.031
IDL> print,interpolate(y,6.5,/cubic) ;Cubic keyword gives a value too small
; 546.367
; But the average of the linear and cubic gives a much closer value
IDL> print,(interpolate(y,6.5,/cubic) + interpolate(y,6.5) )/2.
;648.199
;which happens to also be *exactly* equal to the output of quadterp.pro
IDL> quadterp,x,y,6.5,yy & print,yy
;648.199
--Wayne Landsman landsman@stars.gsfc.nasa.gov
|
|
|