Plot of atomic form factor v scattering vector Q
Matlab codes
% Program to plot out atomic form factor of an element v scattering vector
% and track the maximum accessible scattering vector for a given photon
% energy
clear; close all;
myBlue = [0.4 0.44 0.73];
% Prompt user for input file containing 9 parameters to generate form
% factor. If no file is read, the parameters for sulphur are taken
prompt = 'Do you want to input an elemental parameter file (1 x 9)? Y/N [N]: ';
str = input(prompt,'s');
if isempty(str) || (str == 'N') || (str == 'n')
% Parameters needed to generate form factor for sulphur
a1 = 6.9053;
b1 = 1.4679;
a2 = 5.2034;
b2 = 22.2151;
a3 = 1.4379;
b3 = 0.2536;
a4 = 1.5863;
b4 = 56.172;
c = 0.8669;
vid = VideoWriter('fVq_sulphur.mp4','MPEG-4');
else
prompt = 'Name of input file? : ';
filestr = input(prompt,'s');
fileID = fopen(filestr,'r');
formatSpec = '%f';
ffPars = fscanf(fileID,formatSpec);
a1 = ffPars(1);
b1 = ffPars(2);
a2 = ffPars(3);
b2 = ffPars(4);
a3 = ffPars(5);
b3 = ffPars(6);
a4 = ffPars(7);
b4 = ffPars(8);
c = ffPars(9);
prompt = 'Name of output mp4 movie file (including .mp4)? : ';
mp4str = input(prompt,'s');
vid = VideoWriter(mp4str,'MPEG-4');
end
fmax = a1+a2+a3+a4+c;
% Prompt user for input parameters
prompt = 'Minimum photon energy in eV? [default = 500 eV]: ';
hvmin = input(prompt);
if isempty(hvmin)
hvmin = 500;
end
prompt = 'Maximum photon energy in eV? [default = 30000]: ';
hvmax = input(prompt);
if isempty(hvmax)
hvmax = 30000;
end
prompt = 'Step size in photon energy in eV? [default = 50]: ';
hvstep = input(prompt);
if isempty(hvstep)
hvstep = 50;
end
qmax = 4*pi*hvmax/12398.4; % Maximum accessible Q at 2theta = 180 degrees
% at maximum photon energy
q = (0:0.05:qmax);
f1 = a1*exp(-b1*(q/(4*pi)).^2);
f2 = a2*exp(-b2*(q/(4*pi)).^2);
f3 = a3*exp(-b3*(q/(4*pi)).^2);
f4 = a4*exp(-b4*(q/(4*pi)).^2);
f5 = c;
f = f1 + f2 + f3 + f4 + f5; % Scattering form factor
N=(hvmax-hvmin)/hvstep; % Number of frames
open(vid);
figure('ToolBar','none');
for i = hvmin:hvstep:hvmax
lambda = 12398.4/i; % in Angstrom
qnow = (4*pi/lambda); % max accessible Q in reciprocal Angstrom
% for present photon energy
amin = 2*pi/qnow;
f1 = a1*exp(-b1*(qnow/(4*pi))^2);
f2 = a2*exp(-b2*(qnow/(4*pi))^2);
f3 = a3*exp(-b3*(qnow/(4*pi))^2);
f4 = a4*exp(-b4*(qnow/(4*pi))^2);
f5 = c;
fnow = f1 + f2 + f3 + f4 + f5; % Max scattering vector for ith energy
plot(0.0,0.0,qmax,fmax,q,f,'color',myBlue, 'LineWidth', 2.5);
hold on
plot(qnow,fnow,'ro','MarkerFaceColor','r','MarkerSize',10) % marking the ith data point of x and y
hold off
xaxisLabel = get(gca,'XTickLabel');
set(gca,'FontName','Helvetica','fontsize',18);
set(gca,'TickLength',[0.016, 2]);
xlim([0.0 qmax])
ylim([0 fmax+1])
str1 = 'Q [';
str2 = char(197);
str3 = '^{-1}]';
strTot = [str1,str2,str3];
set(gca,'linewidth',2);
xlabel(strTot);
ylabel('f(Q)');
str1 = 'E = ';
str2 = num2str(i,'% 4.0f');
str3 = ' eV; ';
str4 = '{\lambda} = ';
str5 = num2str(lambda,'%4.3f');
str6 = ' ';
str7 = char(197);
strTot = [str1,str2,str3,str4,str5,str6,str7];
str8 = 'a_{min} = ';
str9 = num2str(amin,'%4.3f');
str10 = ' ';
strTot2 = [str8,str9,str10,str7];
hText = text(qmax/16, 0.1*fmax, strTot, 'FontSize',14);
hText2 = text(qmax/1.4, 0.91*fmax, strTot2, 'FontSize',14);
frame = getframe(gcf);
writeVideo(vid,frame);
end
% Output the movie as an mp4 file
close(vid);
% Program to plot out atomic form factor of an element v scattering angle
% in a polar plot as a function of photon energy
clear; close all;
% Coefficients to determine dependence of f(Q) on Q/4pi (sin theta/lambda)
% for sulphur. Change according to element: look up in the International
% Tables of Crystallography
% Prompt user for input file containing 9 parameters to generate form
% factor. If no file is read, the parameters for sulphur are taken
prompt = 'Do you want to input an elemental parameter file (1 x 9)? Y/N [N]: ';
str = input(prompt,'s');
if isempty(str) || (str == 'N') || (str == 'n')
% Parameters needed to generate form factor for sulphur
a1 = 6.9053;
b1 = 1.4679;
a2 = 5.2034;
b2 = 22.2151;
a3 = 1.4379;
b3 = 0.2536;
a4 = 1.5863;
b4 = 56.172;
c = 0.8669;
vid = VideoWriter('PolarDependence_qV2theta_sulphur.mp4','MPEG-4');
else
prompt = 'Name of input file? : ';
filestr = input(prompt,'s');
fileID = fopen(filestr,'r');
formatSpec = '%f';
ffPars = fscanf(fileID,formatSpec);
a1 = ffPars(1);
b1 = ffPars(2);
a2 = ffPars(3);
b2 = ffPars(4);
a3 = ffPars(5);
b3 = ffPars(6);
a4 = ffPars(7);
b4 = ffPars(8);
c = ffPars(9);
prompt = 'Name of output mp4 movie file (including .mp4)? : ';
mp4str = input(prompt,'s');
vid = VideoWriter(mp4str,'MPEG-4');
end
% Prompt user for input parameters
prompt = 'Minimum photon energy in eV? [default = 500 eV]: ';
hvmin = input(prompt);
if isempty(hvmin)
hvmin = 500;
end
prompt = 'Maximum photon energy in eV? [default = 30000]: ';
hvmax = input(prompt);
if isempty(hvmax)
hvmax = 30000;
end
prompt = 'Step size in photon energy in eV? [default = 50]: ';
hvstep = input(prompt);
if isempty(hvstep)
hvstep = 50;
end
figure('ToolBar','none');
open(vid);
set(0,'defaultfigurecolor',[1 1 1]);
theta = (-pi/2:pi/1000:pi/2)';
twoTheta = 2*theta;
fmax = a1+a2+a3+a4+c;
% Create rgb progression through a rainbow
bot = 8:8:128; % 0 to 7f in 16 steps
top = 128:8:248; % 80 to ff
all = 0:16:240; % 0 to ff
botc = 0*all; % vector of zeros
midc = botc + 127; % vector of 127s
topc = botc + 255; % vector of 255s
rcomp = [top topc topc flip(all) botc botc bot flip(bot)]/256; % R in rgb
gcomp = [botc bot top flip(top) midc flip(bot) botc botc]/256; % G in rgb
bcomp = [botc botc botc botc bot midc midc bot]/256; % B in rgb
%N = (hvmax - hvmin)/hvstep;
for i = hvmin:hvstep:hvmax
lambda = 12398.4/i;
q = (4*pi/lambda).*sin(theta);
f1 = a1*exp(-b1*(q/(4*pi)).^2);
f2 = a2*exp(-b2*(q/(4*pi)).^2);
f3 = a3*exp(-b3*(q/(4*pi)).^2);
f4 = a4*exp(-b4*(q/(4*pi)).^2);
f5 = c;
f = f1 + f2 + f3 + f4 + f5; % Scattering vector q
clii = 1+floor(127*(i-hvmin)/(hvmax - hvmin)); % Index to look up rgb components
polarplot(twoTheta,f,'color',[rcomp(clii), gcomp(clii), bcomp(clii)], 'LineWidth', 2.5);
rlimits = [0 fmax*1.025];
rlim(rlimits);
%axis off;
str1 = 'E = ';
str2 = num2str(i,'% 4.0f');
str3 = ' eV; ';
str4 = '{\lambda} = ';
str5 = num2str(lambda,'%4.3f');
str6 = ' ';
str7 = char(197);
strTot = [str1,str2,str3,str4,str5,str6,str7];
hText = text(-pi/2, 1.25*fmax, strTot, 'FontSize',14);
frame = getframe(gcf);
writeVideo(vid,frame);
end
% Output the movie as an mp4 file
close(vid);