Fitting a Surface to a Function of Two Independent Variables.

This help article provides an example for fitting a surface to a function of two independent variables.

This specific example fits a surface to the function:

f(x,y) = -SIN(2.0 * x) + COS(0.5 * y) .

but can be generalized to work with any function of two variables.

Example:

;++++++++Begin:  The Function to be Fit++++++++
; Please Note:  For better understanding of this
; program, you may want to read the comments to
; the Curve2DFit procedure (below) first.
PRO SurfaceModel, Independent, Parameters, Dependent, $
   Partials
 
; The "Independent" vector is arranged with the numbers
; of the "x" and "y" values in vector at the [0] ("x")
; and [1] ("y") elements.  The remaining part of the
; vector contains "x" then "y".
nX = Independent[0] ;Number of x elements.
nY = Independent[1] ;Number of y elements.
 
; Forming "x" (Xdata) and "y" (Ydata) arrays.
; The (i, j) value of each array coincides with the
; (i, j) value of the "Dependent" ("z" or the surface
; data) variable.
Xdata = Independent[2:nX+1] # REPLICATE(1.0, nY)
Ydata = REPLICATE(1.0, nX) # Independent[nX+2:*]
 
; Calculating the surface data (Dependent) from the
; Xdata and Ydata variables.
Dependent = -SIN(Parameters[0] * Xdata) + $
   COS(Parameters[1] * Ydata)
 
; REFORMing the surface data from an two-dimensional
; array to a vector to be used in CURVEFIT.
Dependent = REFORM(Dependent, N_ELEMENTS(Dependent), $
   /OVERWRITE)
 
; This condition is used to speed up CURVEFIT.  CURVEFIT
; sometimes calls this user-defined procedure to
; obtain data points for both the function and it's
; derivatives; sometimes CURVEFIT just needs the data
; of the function.
if N_PARAMS() ge 4 then begin
   Partials = FLTARR((nX * nY), 2)
   Partials[*,0] = -Xdata * COS(Parameters[0] * Xdata)
   Partials[*,1] = -Ydata * SIN(Parameters[1] * Ydata)
endif
 
END
;+++++++++End:  The Function to be Fit+++++++++
 
;++++++++Begin:  Example Procedure++++++++
PRO Curve2DFit ; Procedure for using CURVEFIT to fit
; two-dimensional data to a curve.
 
; Creating "independent" data.
nX = 31
x = FINDGEN(nX) / 10
nY = 61
y = FINDGEN(nY) / 10
 
; Creating initial guess of parameters.
Parameters = [2.0, 0.5] + 0.5
 
; Calculating and displaying "original" data.
SurfaceModel, [nX, nY, x, y], [2.0, 0.5], z
WINDOW, 0, XSIZE = 350, YSIZE = 350
SURFACE, REFORM(z, nX, nY), x, y, /XSTYLE, /YSTYLE, $
   title = 'Original Data'
 
; Calling CURVEFIT.  All the "dependent" values are
; equally weighted.
Weights = REPLICATE(1.0, N_ELEMENTS(z))
zFit = CURVEFIT([nX, nY, x, y], z, Weights, $
   Parameters, FUNCTION_NAME = 'SurfaceModel')
 
; Calculating and outputting error and the resulting
; parameters.
PRINT, 'Error:', TOTAL( (zFit-z)^2 ) ;Squared error.
PRINT, 'Parameters:', Parameters
 
; REFORMing resulting curve fit into a two-dimensional
; array and displaying it.
zFit = REFORM(zFit, nX, nY, /OVERWRITE)
WINDOW, 1, XSIZE = 350, YSIZE = 350
SURFACE, zFit, x, y, /XSTYLE, /YSTYLE, $
   title = 'Fitted Data'
 
END
;+++++++++End:  Example Procedure+++++++++