今天早上接到师姐的紧急传呼，说是将一组数据按照确定的公式进行指数拟合。可是来北京后就没用过matlap,很多基础的东西都忘 了，没事，有事无事找度娘啊！话不多说，说正事。我仔细阅读了要求文档，核心有两个，1：方程的构建，:2：计算流程
先拿公式：
（方程1）
已知数据：每个波段所对应的ad(λ)值和ad(440)（波长为440的值），λ为400-700 nm的波长，间隔1 nm
求解sd
波段数据如下：
数据特点是由很多site,每个site都要分别计算，想了下简单，一个循环的事，难点是如何动态设置已知变量ad(440)，查了一些，自己试验了一下发现也挺简单
。
不多说先放代码，下面是计算每个site拟合方程的代码：
function []=CaculaterSd()
data=xlsread(''); %这是波段的原始数据
[r,c]=size(data);
result=zeros(r+2,c);
for i=2:c
a=data(1,i);
x=data(2:r,1);
y=data(2:r,i);
x0=0; %因为估算了一下，0值比较接近真实的sd值，所以选择0了
%定义拟合函数
f=@(sd,x)a*exp(sd*(440-x));
sd=lsqcurvefit(f,x0,x,y);
y1=a*exp(sd*(440-x));
% -----计算r方方法一------------
r2=norm(y1-mean(y1))^2/norm(y-mean(y))^2;
% -----计算r方 方法二------------
% ss0=0;
% ss1=0;
% for k=1:r-1
% ss0=(y(k)-mean(y))*(y1(k)-mean(y1))+ss0;
% ss1=(y(k)-mean(y))^2*(y1(k)-mean(y1))^2+ss1;
% end
% r2=(ss0/sqrt(ss1))^2;
%-----------------------------------
%存入sd、ad、r2
result(r+1,i)=sd;
result(r+2,i)=r2;
result(2:r,i)=y1;
end
xlswrite('result_sd_a8.xls',result)
end
计算完后，得查看一下拟合效果吧
这个也简单，上代码
function []=plotfun()
data=xlsread(''); %读取原始数据
ExcuteData=xlsread('result_sd_a7.xls'); %读取处理数据
[r,c]=size(data);
x=data(2:r,1);%获取波段号
y=data(2:r,3); %读取site2的波段值
y1=ExcuteData(2:r,3);%读取处理后site2的波段值
plot(x,y,'-or',x,y1,'-*b');legend('实际值','计算值')
end
看了下图，拟合效果很不错啊
