The emergence of phase contrast with sample-detector distance

Matlab code

% Movie cartoon of a phase-contrast tomography experiment

clear all; close all

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

vid.Quality = 100;

vid.FrameRate = 30;

open(vid);

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

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

camSize = 2;

% Set up head and tail parts of zebrafish body as half ellipsoids

[x,y,z] = sphere(270);

x(z<0) = NaN ; y(z<0) = NaN ; z(z<0) = NaN;

% Set up middle part of zebrafish body as modulated cylinder

t = 0:pi/80:pi/2;

r = 0.2+0.8*(sin(t)).^2;

[a,b,c] = cylinder(r,100);

% Set up eyes and otoliths

[x2,y2,z2] = sphere(200);

% Set up backbone of zebrafish as modulated cylinder

t2 = 0:pi/80:4*pi;

r2 = 0.2+0.8*(sin(2*t2)).^6;

[a2,b2,c2] = cylinder(r2,100);

% Set up mirror faces

xm1 = [0 0 0 0 0];

ym1 = [-1 1 1 -1 -1];

zm1 = [-1 -1 1 1 -1];

myBlue = [0.4 0.44 0.73];

skin = [0.945 0.761 0.49];

lp = [-2 -1 1];

view(-34,5);

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

% Set up spherical bowl for lens surfaces

[X,Y,Z] = sphere(250);

X(Z>-0.8) = NaN ; Y(Z>-0.8) = NaN ; Z(Z>-0.8) = NaN;

rL = 1.2;

pos1 = [0.02 0.15 0.75 0.7];

pos2 = [0.8 0.2 0.19 0.6];

% Apply the transform

for ZZ = 0:610 % 150 frames, double the number of projections

hold off

subplot('position',pos1);

newplot

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

% Draw zebrafish embryo

surf(a/8,b/4,c-1,'EdgeColor','none','FaceColor',[0.5 0.5 0.5],...

'FaceLighting','gouraud','DiffuseStrength',1,'FaceAlpha',0.4); % body part 1

hold on

surf(x/8,y/4,z*0.5,'EdgeColor','none','FaceColor',[0.5 0.5 0.5],...

'FaceLighting','gouraud','DiffuseStrength',1,'FaceAlpha',0.4); % body part 2

surf(x/40,y/20,-z*0.5-1,'EdgeColor','none','FaceColor',[0.5 0.5 0.5],...

'FaceLighting','gouraud','DiffuseStrength',1,'FaceAlpha',0.4); % body part 3

% Eyes

surf(x2/10,y2/14+0.14,z2/7+0.25,'EdgeColor','none','FaceColor',[0.5 0.5 0.5],...

'FaceLighting','gouraud','DiffuseStrength',1,'FaceAlpha',0.25);

surf(x2/10,y2/14-0.14,z2/7+0.25,'EdgeColor','none','FaceColor',[0.5 0.5 0.5],...

'FaceLighting','gouraud','DiffuseStrength',1,'FaceAlpha',0.25);

surf(x2/25,y2/25+0.19,z2/25+0.25,'EdgeColor','none','FaceColor',[0.5 0.5 0.5],...

'FaceLighting','gouraud','DiffuseStrength',1,'FaceAlpha',0.4);

surf(x2/25,y2/25-0.19,z2/25+0.25,'EdgeColor','none','FaceColor',[0.5 0.5 0.5],...

'FaceLighting','gouraud','DiffuseStrength',1,'FaceAlpha',0.4);

% Otoliths

s = surf(x2/35,y2/35-0.2,z2/25-0.05,'EdgeColor','none','FaceColor',[0.2 0.2 0.2],...

'FaceLighting','gouraud','DiffuseStrength',1,'FaceAlpha',0.88);

rotate(s,[1 0 0],29,[0,0,0]);

s = surf(x2/35,y2/35+0.2,z2/25-0.05,'EdgeColor','none','FaceColor',[0.2 0.2 0.2],...

'FaceLighting','gouraud','DiffuseStrength',1,'FaceAlpha',0.88);

rotate(s,[1 0 0],-29,[0,0,0]);

% Backbone

s = surf(a2/35,b2/35,c2-0.7,'EdgeColor','none','FaceColor',[0.5 0.5 0.5],...

'FaceLighting','gouraud','DiffuseStrength',1,'FaceAlpha',0.4);

axis off

axis equal

xlim([-4.05 7.1]);

ylim([-1.5 1.5]);

zlim([-1.5 4.1]);

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

% Incident x-rays denoted by arrow

mArrow4([-4 0 0],[-1.4 0 0],'color','y','stemWidth',0.4,...

'tipWidth',0.61,'facealpha',0.7);

% Scintillator screen

scTh = 0.1;

scWi= 1.5;

scHe= 1.5;

scPos = 0.14+ZZ/305;

plotcube([scTh scWi scHe],[scPos -scWi/2 -scHe/2],0.5,'g',0);

view(-34,5);

% Lens

% 1st convex surface

lens1 = surf(2*X+2+ZZ/305,2*Y,2*Z+1.6,'FaceColor',myBlue,'LineStyle','none',...

'FaceLighting','gouraud','DiffuseStrength',1);%,'FaceAlpha',0.35) ;

alpha(lens1,0.5);

rotate(lens1,[0 1 0],90,[2+ZZ/305,0,0]);

% 2nd convex surface

lens2 = surf(2*X+2.1+ZZ/305,2*Y,2*Z+1.6,'FaceColor',myBlue,'LineStyle','none',...

'FaceLighting','gouraud','DiffuseStrength',1);%,'FaceAlpha',0.35) ;

alpha(lens2,0.5);

rotate(lens2,[0 1 0],-90,[2.1+ZZ/305,0,0]);

% Outer rim of lens

[lx, ly, lz] = cylinder(rL,100);

camAp = surf(lx+2+ZZ/305,ly,lz/10,'FaceColor',myBlue,...

'FaceLighting','gouraud','DiffuseStrength',1,'LineStyle','none');

alpha(camAp,0.5);

rotate(camAp,[0 1 0],90,[2+ZZ/305,0,0]);

% 45-degree reflecting mirror

mface1 = fill3(xm1+4+ZZ/305,ym1,zm1,[0.5 0.5 0.5],'LineStyle','none',...

'FaceLighting','gouraud','DiffuseStrength',1,'SpecularStrength',1);

rotate(mface1,[0 1 0],45,[4+ZZ/305,0,0]);

mface2 = fill3(xm1+4.1+ZZ/305,ym1,zm1-0.1,[0.5 0.5 0.5],'LineStyle','none',...

'FaceLighting','gouraud','DiffuseStrength',1,'SpecularStrength',1);

rotate(mface2,[0 1 0],45,[4.1+ZZ/305,0,-0.1]);

mface3 = fill3(ym1*0.05*1.4142+0.05*1.4142+4+ZZ/305,xm1-1,zm1,...

[0.5 0.5 0.5],'LineStyle','none',...

'FaceLighting','gouraud','DiffuseStrength',1,'SpecularStrength',1);

rotate(mface3,[0 1 0],45,[4+ZZ/305,-1,0]);

mface4 = fill3(xm1+4+ZZ/305-1/sqrt(2),ym1,...

zm1*0.05*sqrt(2)-1/sqrt(2)-0.1/sqrt(2),[0.5 0.5 0.5],'LineStyle','none',...

'FaceLighting','gouraud','DiffuseStrength',1,'SpecularStrength',1);

rotate(mface4,[0 1 0],-45,[4+ZZ/305-1/sqrt(2),0,-1/sqrt(2)]);

% Dot-dashed lines

plot3([-2 4+ZZ/305],[0 0],[0 0],'k','LineStyle','-.','LineWidth',1.6);

plot3([4+ZZ/305 4+ZZ/305],[0 0],[0 2],'k','LineStyle','-.','LineWidth',1.6);

% Camera body

% Main body

plotcube([camSize camSize camSize],[4+ZZ/305-camSize/2 -camSize/2 2],...

1,[0.2 0.2 0.2],0);

% Screw ring thread outer surface

rCamAp = 0.9*camSize/2;

[cx, cy, cz] = cylinder(rCamAp,100);

camAp = surf(cx+4+ZZ/305,cy,cz*camSize/16+2-camSize/16,'FaceColor',...

[0.5 0.5 0.5],'LineStyle','none',...

'FaceLighting','gouraud','DiffuseStrength',1);

hold off

subplot('position',pos2);

imshow(frames(:,:,1:3,ZZ+1))

frame = getframe(gcf);

writeVideo(vid,frame);

hold off

end

% Output the movie as an mpg file

close(vid);