Fit multiple systems from FOR loop
Posted: Wed Apr 04, 2018 1:13 am
I am working on the problem previously described in viewtopic.php?f=3&t=282. I am now implementing an automated fitting of EPR spectra based on a set of subspectra automatically generated by a FOR loop, based on a theoretical model.
My main function is:
My custom function is:
The problem is that each subspectrum seems to be almost identical, despite the variations of the Jij values, which should greatly influence the gperp components. I suspect that there is a problem in the invocation of
I have double-checked with a model with two manually created subsystems and it works as predicted. Does anyone see an error that eludes me?
Thanks!
My main function is:
Code: Select all
clear all; close all; format long;
[B,spc,params] = eprload('AKB1004A_4.3K_X_powder_01.DTA');
%--------------------------------------------------------------------------
Si = 3/2; % The spin of the individual ions
Jav=-10.9; % Isotropic exchange in cm-1 according to the -2JSiSj formalism
r=0.05; % Radial distortion coefficient
gz=1.971; %g|| common for all ions
gxy=1.971; %gperpendicular common for all ions
iter=3; %Number of subspectra
wG=0; %Gaussian linewidth
wL=2; %Lorentzian linewidth
Gx = 0; Gy = 0; Gz = 0.2; % Antisymmetric exchange pseudovector components
SysA.Si = Si;
SysA.Jav = Jav;
SysA.r = r;
SysA.gz = gz;
SysA.gxy = gxy;
SysA.iter = iter;
SysA.wG = wG;
SysA.wL = wL;
SysA.Gx = Gx;
SysA.Gy = Gy;
SysA.Gz = Gz;
VaryA.Jav = 0;
VaryA.r = SysA.r;
VaryA.Gz = Gz;
VaryA.gz = 0.0;
VaryA.gxy = 0.0;
VaryA.wG = 0;
VaryA.wL = 0;
Exp.Temperature = 4.3; Exp.mwFreq=params.MWFQ/1e9; Exp.CenterSweep=[params.A1CT*1e3 params.A1SW*1e3]; Exp.nPoints=params.XPTS;
FitOpt.Method = 'simplex fcn'; % simplex algorithm, data as is
FitOpt.Scaling = 'lsq0';
esfit('constrainJg_equilateral_dynamic_EPR_fit',spc,SysA,VaryA,Exp,[],FitOpt);Code: Select all
function y = constrainJg_equilateral_dynamic_EPR_fit(SysA,Exp,Opt)
cm=100*clight/1e6; % Conversion constant from cm-1 to MHz
%Define the theta angles (y-O-Mj) for each metal center
theta1=0;
theta2=-2*pi/3;
theta3=2*pi/3;
GG = [0 SysA.Gz -SysA.Gy; -SysA.Gz 0 SysA.Gx; SysA.Gy -SysA.Gx 0];
%Variation of phi angles of the central oxide
for i=1:SysA.iter
phi(i)=((i-1)/SysA.iter)*(2*pi/3);
%Construct multiplicative factors for each spin pair and for each iteration
a12(i)=SysA.r*(cos(phi(i)-theta1)+cos(phi(i)-theta2));
a13(i)=SysA.r*(cos(phi(i)-theta1)+cos(phi(i)-theta3));
a23(i)=SysA.r*(cos(phi(i)-theta2)+cos(phi(i)-theta3));
% Individual isotropic exchange constants
J12(i)=SysA.Jav*(1+a12(i));
J13(i)=SysA.Jav*(1+a13(i));
J23(i)=SysA.Jav*(1+a23(i));
%The individual systems
%--------------------------------------------------------------------------
Sys(i).S=[SysA.Si SysA.Si SysA.Si];
Sys(i).g=[SysA.gxy SysA.gz; SysA.gxy SysA.gz; SysA.gxy SysA.gz];
Sys(i).gFrame = [90 90 -45; 30 -90 235; -30 90 -45]*pi/180;
% Sys(i).weight = 1/SysA.iter;
Sys(i).ee = [-2*J12(i)*eye(3) - 2*GG; -2*J13(i)*eye(3) + 2*GG; -2*J23(i)*eye(3) - 2*GG]*cm;
Sys(i).lwpp=[SysA.wG SysA.wL];
[B,calc] = pepper(Sys(i),Exp); % Calculate subspectrum No (i)
calcspc(i,:) = calc; % Put the result in a row of an array
end
y = sum(calcspc); % Sum each subspectrum in the array
figure(1)
plot(B,calcspc)
hold on
plot(B,y)
legend('show')
hold off
figure(2)
plot(phi,J12,phi,J13,phi,J23)
legend('show')
returnpepper in:
Code: Select all
[B,calc] = pepper(Sys(i),Exp); % Calculate subspectrum No (i)
calcspc(i,:) = calc; % Put the result in a row of an arrayThanks!