Cartoon of forward-scattering CXDI for noncrystalline samples

Matlab code

% Cartoon of CXDI in forward scattering direction. Scattered object is an

% ellipsoid with semi-axes of 11, 5, and 5 scattering centers and a random

% electron-density distribution - kind of like a paramecium?

clear

close all

figure('units','pixels','position',[0 0 1920 1080],'ToolBar','none');

set(0,'defaultfigurecolor',[1 1 1]);

vid = VideoWriter('CXDI.mp4','MPEG-4');

vid.Quality = 100;

vid.FrameRate = 30;

open(vid);

myGold = [0.8 0.64 0];

lp = [2 3 1];

lp2 = [1 -2 -1];

% Set up coordinates of front cap of Ewald sphere - I don't want to plot

% the entire sphere, as it is too large and its semiopaque nature would

% compromise the interesting part of the cartoon

[Xe,Ye,Ze] = sphere(2000);

Xe = ESrad*Xe - ESrad; % Shift sphere so front of sphere at origin (0,0,0)

xnot = -80; % Value of x, below which the Ewald sphere is removed

% Only surviving part of sphere for Xe > xnot - a small slice of the sphere

Xe(Xe<xnot) = NaN ; Ye(Xe<xnot) = NaN ; Ze(Xe<xnot) = NaN;

% Collection of random numbers to be sculpted into an ellipsoidal shape

X = rand(23,11,11);

XX = X; % Use XX as replica for FT purposes (FT cannnot handle NaNs)

% Shave off elements outside ellipsoidal surface. For the FT (XX), make

% these zero; for plotting (X), make these NaN

for ii = 1:23

for jj = 1:11

for kk = 1:11

R = (((ii-12)/12)^2 + ((jj-6)/6)^2 + ((kk-6)/6)^2)^0.5;

if (R > 1)

X(ii,jj,kk) = NaN;

XX(ii,jj,kk) = 0;

end

hold on

end

end

end

% Reshape matrices of X for scatter3 plotting purposes

[xR, yR, zR] = meshgrid(1:size(X,2), 1:size(X,1), 1:size(X,3));

xxR = reshape(xR,[],1)-6;

yyR = reshape(yR,[],1)-12;

zzR = reshape(zR,[],1)-6;

X = reshape(X,[11*11*23,1]);

% ________________________________________________________________________

% Shift FT so maximum is in the middle

Y = fftshift(Y);

maxY = max(Y(:));

% Manipulate FT data for plotting purposes

% Only plot out FT to a radius equal to pad/2 and for signal

% stronger than 2% of the FT maximum

if (Y(ii,jj,kk) < 0.02*maxY) || (R > pad/2)

Y(ii,jj,kk) = NaN;

end

end

end

end

% Reshape manipulated FT data for plotting with scatter3

[x, y, z] = meshgrid(1:size(Y,1), 1:size(Y,2), 1:size(Y,3));

cmin = min(log(Y));

cmax = max(log(Y));

% ________________________________________________________________________

% Loop for rotation of object and its FT by 360 degrees

for theta = 0:359

subplot('position',[-0.2 -0.2 1.4 1.55]);

hold off

newplot

hold off

% Coordinates at theta of rotating paramecium ellipsoid

xxRnow = 3*(xxR*cosd(theta) - yyR*sind(theta))-500;

yyRnow = 3*(xxR*sind(theta) + yyR*cosd(theta));

hold on

% Plot ellipsoid paramecium

scatter3(xxRnow,yyRnow,3*zzR,50,X*cmax,'filled','MarkerFaceAlpha',...

0.25,'MarkerEdgeAlpha',0)

colormap jet

% Incident parallel radiation on object

mArrow3([-640 0 0],[-500 0 0],...

'color',myGold,'stemWidth',2.5,...

'tipWidth',7,'facealpha',0.25);

% Plot Ewald sphere

ewaldSurf = surf(Xe,Ye,Ze,'FaceColor',myGold,'LineStyle','none',...

'FaceLighting','flat','DiffuseStrength',1,'SpecularStrength',0.55);

% Vary alpha of Ewald sphere in x-direction with offset of 50 so it

% remains semitransparent for all elements of surf = ewaldSphere

alpha(ewaldSurf,ewaldSurf.XData-50);

% Scattering pattern - transform x- and y-coordinates according to

% rotation angle theta then plot out with scatter3

xxnow = xx*cosd(theta) - yy*sind(theta);

yynow = xx*sind(theta) + yy*cosd(theta);

scatter3(xxnow,yynow,zz,Y.^0.5/5,log(Y(:)),'filled',...

'MarkerFaceAlpha',0.2,'MarkerEdgeAlpha',0)

colormap turbo

caxis([cmin cmax])

% Plot part of scattering pattern that sits on the Ewald sphere

% highlighted in red

clear xxEw yyEw zzEw YEw

mm = 1;

for ii = 1:size(xx)

if (abs(xxnow(ii)) < 0.64)

xxEw(mm) = xxnow(ii);

yyEw(mm) = yynow(ii);

zzEw(mm) = zz(ii);

YEw(mm) = Y(ii);

mm = mm+1;

end

end

scatter3(xxEw,yyEw,zzEw,YEw/10,'r','filled','MarkerFaceAlpha',1,...

'MarkerEdgeAlpha',0)

% Cone of forward-scattered x-rays

t = 0:1;

% Transparency of cone representing scattered radiation proportional to

% the sum of the scattered radiation that sits on the Ewald sphere, YEw

coneAlp = sum(YEw,"omitnan")/2e6;

'FaceColor',myGold,'LineStyle','none','FaceLighting','flat',...

'DiffuseStrength',1);

% Finish off defining plotting appearance

axis off

axis equal

xlim([-640 70]);

ylim([-280 280]);

zlim([-280 280]);

set(gca, 'Projection','perspective');

view(35,22)

light('Position',lp,'Style','infinite');

light('Position',lp2,'Style','infinite');

% Capture frame and write to video

frame = getframe(gcf);

writeVideo(vid,frame);

end

close(vid);