城市表层土壤重金属污染分析

%建立函数myfun求解参数的初始值 function F = myfun(x) clc

x1=[2883 3692 7 41.7] x2=[4742 7293 9 27.8] x3=[4948 7293 6 34.8] x4=[5006 8846 6 30.3] F =

[x(4)/(4*pi*((0.16*(x1(1)-x(1))*(1+0.0004*(x1(1)-x(1)))^(-0.5)))*(0.14*(x1(1)-x(1))*(1+0.003*(x1(1)-x(1)))^(-0.5)))...

*exp(-(x1(2)-x(2))^2/(2*((0.16*(x1(1)-x(1))*(1+0.0004*(x1(1)-x(1)))^(-0.5)))^2)-(x1(3)-x(3))^2/(2*(0.14*(x1(1)-x(1))*(1+0.003*(x1(1)-x(1)))^(-0.5))^2))-x1(4),

x(4)./(4*pi*((0.16*(x2(1)-x(1)).*(1+0.0004*(x2(1)-x(1))).^(-0.5))).*(0.14*(x2(1)-x(1)).*(1+0.003*(x2(1)-x(1))).^(-0.5)))...

.*exp(-(x2(2)-x(2)).^2./(2*((0.16*(x2(1)-x(1)).*(1+0.0004*(x2(1)-x(1))).^(-0.5)))^2)-(x2(3)-x(3)).^2./(2*(0.14*(x2(1)-x(1)).*(1+0.003*(x2(1)-x(1))).^(-0.5))^2))-x2(4),

x(4)./(4*pi*((0.16*(x3(1)-x(1)).*(1+0.0004*(x3(1)-x(1))).^(-0.5))).*(0.14*(x3(1)-x(1)).*(1+0.003*(x3(1)-x(1))).^(-0.5)))...

.*exp(-(0-x(2)).^2./(2*((0.16*(x3(1)-x(1)).*(1+0.0004*(x3(1)-x(1))).^(-0.5)))^2)-(16-x(3)).^2./(2*(0.14*(x3(1)-x(1)).*(1+0.003*(x3(1)-x(1))).^(-0.5))^2))-x3(4),

x(4)./(4*pi*((0.16*(x4(1)-x(1)).*(1+0.0004*(x4(1)-x(1))).^(-0.5))).*(0.14*(x4(1)-x(1)).*(1+0.003*(x4(1)-x(1))).^(-0.5)))...

.*exp(-(x4(2)-x(2)).^2./(2*((0.16*(x4(1)-x(1)).*(1+0.0004*(x4(1)-x(1))).^(-0.5)))^2)-(x4(3)-x(3)).^2./(2*(0.14*(x4(1)-x(1)).*(1+0.003*(x4(1)-x(1))).^(-0.5))^2))-x4(4)]

%x0=[-5 -5 -5 -5]

%options=optimset('Display','iter'); %[x,fval] = fsolve(@myfun,x0,options)

%对数据回归分析

function f=curvefun1(x,T) x1=T(:,1);y=T(:,2);z=T(:,3);

sigmay=0.16*(x1-x(1)).*(1+0.0004*(x1-x(1))).^(-0.5); sigmaz=0.14*(x1-x(1)).*(1+0.0003*(x1-x(1))).^(-0.5);

f=x(4)./(4*pi*sigmay.*sigmaz).*exp(-(y-x(2)).^2./(2*sigmay)-(z-x(3)).^2./(2*sigmaz))

%[beta，r，J]=nlinfit（T,y,’curvefun1’,x0）