第1个回答 2008-06-18
二维数据拟合
试试这个
function p = polyfit2d(x,y,z,n,m)
% P= POLYFIT2D(x,y,z,n,m) finds the coefficients of a
% polynomial function of 2 variables formed from the
% data in vectors x and y of degrees n and m, respectively,
% that fit the data in vector z in a least-squares sense.
%
% The regression problem is formulated in matrix format as:
%
% z = A*P So that if the polynomial is cubic in
% x and linear in y, the problem becomes:
%
% z = [y.*x.^3 y.*x.^2 y.*x y x.^3 x.^2 x ones(length(x),1)]*
% [p31 p21 p11 p01 p30 p20 p10 p00]'
%
% Note that the various xy products are column vectors of length(x).
%
% The coefficents of the output p
% matrix are arranged as shown:
%
% p31 p30
% p21 p20
% p11 p10
% p01 p00
%
% The indices on the elements of p correspond to the
% order of x and y associated with that element.
%
% For a solution to exist, the number of ordered
% triples [x,y,z] must equal or exceed (n+1)*(m+1).
% Note that m or n may be zero.
%
% To evaluate the resulting polynominal function,
% use POLYVAL2D.
% Perry W. Stout June 29, 1995
% 4829 Rockland Way
% Fair Oaks, CA 95628
% (916) 966-0236
% Based on the Matlab function polyfit.
if any((size(x) ~= size(y)) | (size(z) ~= size(y)))
error('X, Y,and Z vectors must be the same size')
end
x = x(:); y = y(:); z= z(:); % Switches vectors to columns--matrices, too
if length(x) < (n+1)*(m+1)
error('Number of points must equal or exceed order of polynomial function.')
end
n = n + 1;
m = m + 1; % Increments n and m to equal row, col numbers of p.
a = zeros(max(size(x)),n*m);
% Construct the extended Vandermonde matrix, containing all xy products.
for i1= 1:m
for j1=1:n
a(:,j1+(i1-1)*n) = (x.^(n-j1)).*(y.^(m-i1));
end
end
p1 = (a\z);
% Reform p as a matrix.
p=[];
for i1=1:m
p=[p, p1((n*(i1-1)+1):(n*i1))];
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function z = polyval2d(c,x,y)
% z = POLYVAL2D(V,sx,sy) Two dimensional polynomial evaluation.
%
% If V is a matrix whose elements are the coefficients of a
% polynomial function of 2 variables, then POLYVAL2D(V,sx,sy)
% is the value of the polynomial evaluated at [sx,sy]. Row
% numbers in V correspond to powers of x, while column numbers
% in V correspond to powers of y. If sx and sy are matrices or
% vectors,the polynomial function is evaluated at all points
% in [sx, sy].
%
% If V is one dimensional, POLYVAL2D returns the same result as
% POLYVAL.
%
% Use POLYFIT2D to generate appropriate polynomial matrices from
% f(x,y) data using a least squares method.
% Perry W. Stout June 28, 1995
% 4829 Rockland Way
% Fair Oaks, CA 95628
% (916) 966-0236
% Based on the Matlab function POLYVAL.
% Polynomial evaluation c(x,y) is implemented using Horner's slick
% method. Note use of the filter function to speed evaluation when
% the ordered pair [sx,sy] is single valued.
if nargin==2
y=ones(x); end
if any(size(x) ~= size(y))
error('x and y must have the same dimensions.')
end
[m,n] = size(x);
[rown,coln]= size(c);
if (m+n) == 2
% Use the built-in filter function when [sx,sy] is single valued
% to implement Hoerner's method.
z= 0;
for i1= 1:coln
ccol= c(:,coln+1-i1);
w = filter(1,[1 -x],ccol);
w = w(rown)*(y^(i1-1));
z= z+w;
end
return
end % Of the scalar computation
% Do general case where X and Y are arrays
z = zeros(m,n);
for i1=1:coln
ccol= c(:,coln+1-i1);
w= zeros(m,n);
for j1=1:rown
w = x.*w + ccol(j1) * ones(m,n);
end
z= z+w.*(y.^(i1-1));
end本回答被网友采纳