📄 spcfun1.m
字号:
function f = spcfun1(params)
global data
global irf
global weight
global fit
global amps
global numTao
global numTaoFixed
global taoFixed
global t
global chisq
global textBox
global spcLeftEnd
global spcRightEnd
global tpch
%This is the function referred to in the call to "leastsq".
%It assumes the form:
%y(x)=a11*exp(-k1*t)+a12*exp(-k2*t)+. . .+a1n*exp(-kn*t)
% a21*exp(-k1*t)+a22*exp(-k2*t)+. . .+a2n*exp(-kn*t)
% .
% .
% .
% am1*exp(-k1*t)+am2*exp(-k2*t)+. . .+amn*exp(-kn*t)
%and makes available the vector of k values and the matrix
%of a values to the script file "fit.m".
dataTemp = data(spcLeftEnd:spcRightEnd,:);
irfTemp = irf(spcLeftEnd:spcRightEnd,:);
numPoints = size(dataTemp,1);
numCurves = size(dataTemp,2);
curve = zeros(numPoints,numTao + numTaoFixed+1);
fit = zeros(size(dataTemp));
a = zeros(numTao + numTaoFixed+1,numCurves);
taoFit = zeros(1,numTao + numTaoFixed);
tao = params(1:numTao);
shift = params(numTao + 1:numTao + numCurves);
shift = shift./tpch;
taoFit(1,1:numTaoFixed) = taoFixed;
taoFit(1,numTaoFixed+1:numTao + numTaoFixed) = tao;
%taoFit = taoFit./tpch;
bottom = floor(shift);
top = ceil(shift);
R = exp(-tpch./taoFit);
for i = 1:numCurves,
if (top(i) < 1)
top(i) = 1;
end
curve(1,:) = irfshift(irfTemp(top(i),i),irfTemp(top(i)+1,i),shift(i));
for count = 1:(numPoints-1),
top(i) = bottom(i) + count;
if (top(i) < 1) top(i) = 1; end
if (top(i) > (numPoints-2)) top(i) = numPoints-2; end
curve(count+1,1:numTao+numTaoFixed) = R.*(0.5*irfshift(irfTemp(top(i),i),irfTemp(top(i)+1,i),shift(i))+...
curve(count,1:numTao+numTaoFixed)) + 0.5 * irfshift(irfTemp(top(i)+1,i),irfTemp(top(i)+2,i),shift(i));
curve(count+1,numTao+numTaoFixed+1)= irfshift(irfTemp(top(i),i),irfTemp(top(i)+1,i),shift(i));
end
%space = curve;
a(:,i) = curve\dataTemp(:,i);
z(:,i) = curve*a(:,i);
end
f = (z - dataTemp)./weight(spcLeftEnd:spcRightEnd,:);
amps = a;
fit=z;
chisq = trace(f'*f)/(numPoints*numCurves-(numTao + numTaoFixed)-(numTao + numTaoFixed)*numCurves-numCurves);
set(textBox,'string',['chisq = ' num2str(chisq)])
drawnow;
function G = irfshift(f1,f2,x0)
G = f1 + (x0 - floor(x0)) * (f2 - f1);
⌨️ 快捷键说明
复制代码
Ctrl + C
搜索代码
Ctrl + F
全屏模式
F11
切换主题
Ctrl + Shift + D
显示快捷键
?
增大字号
Ctrl + =
减小字号
Ctrl + -