Cartoon of scanning-SAXS tensor tomography experimental setup 

Matlab code

% Cartoon of the experimental setup used in scanning SAXS tensor tomography

% (SSTT). Data with thanks from the group of Marianne Liebi, Swiss Light

% Source. 


clear

close all

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

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


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

vid.Quality = 100;

vid.FrameRate = 30;

open(vid);


midGrey = [0.5 0.5 0.5];

boneYellow = [224 221 214]/255;

myGold = [0.8 0.64 0];

lp = [-2 -2 1];

pos1 = [-0.1 -0.1 1.2 1.2];


% Dimensions in matlab units

incRadLength = 14; % Distance upstream to which incident radiation extends

FZPpos = -2; % Position upstream of FZP relative to focus @ (0,0,0)

samplePos = 0.0; % Position of sample downstream of focus

DetPos = 20; % Position of detector downstream

halfDetWidth = 4; % Half width of detector

sampleHt = 4.0; % Height of sample

basePlateRad = 2.8; 

imax = 10; % Determines outer radius of FZP


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

aa = 0;

bb = 0.64;

cc = 5; 

% Add roughness to vertices of sphere vertex coordinates to render the bone

% surface to be rough - adds a touch of verisimilitude 

varx = 0.02*(rand(size(x)) - 0.5); 

vary = 0.02*(rand(size(y)) - 0.5);

varz = 0.02*(rand(size(z)) - 0.5);

x = x + varx; 

y = y + vary; 

z = z + varz; 


t = -3.5:0.1:3.5;

tt = (t.^2 + 20)/32;

[x2,y2,z2] = cylinder(tt,200);

% Add roughness to vertices of cylinder vertex coordinates to render the 

% bone surface to be rough - adds a touch of verisimilitude 

varx = 0.01*(rand(size(x2)) - 0.5); 

vary = 0.01*(rand(size(y2)) - 0.5);

varz = 0.01*(rand(size(z2)) - 0.5);

x2 = x2 + varx; 

y2 = y2 + vary; 

z2 = z2 + varz; 


[x3,y3,z3] = cylinder(1,200);


th = 0:pi/180:2*pi;

x4 = basePlateRad*sin(th);

y4 = basePlateRad*cos(th);

z4 = 0.0*x4;


% Set up range of movements of bone sample

xmovement1 = 0:0.05:1;

xmovement2 = 0.95:-0.05:-1;

xmovement3 = -0.95:0.05:0;

xrange = 0.8; 

xmovement = xrange*[xmovement1 xmovement2 xmovement3];


ymovement1 = 0:0.05:1;

ymovement2 = 0.95:-0.05:-1;

ymovement3 = -0.95:0.05:0;

yrange = 1.0; 

ymovement = yrange*[ymovement1 ymovement2 ymovement3];


thetamovement = 0:4:360;


phimax = 44;

phimovement1 = 0:1:phimax;

phimovement2 = phimax:-1:-phimax;

phimovement3 = -phimax:1:0;

phimovement = [phimovement1 phimovement2 phimovement3];


% Set up cones representing x-rays

% Incident x-rays, slightly divergent

t = 1:-0.001:0.8;

[xC1,yC1,zC1] = cylinder(t,800);

% Cone of focussed x-rays and divergent scattered x-rays

t = 0:0.001:1;

[xC2,yC2,zC2] = cylinder(t,800);


% x- and y-coordinates of FZP rings uses this range of angles 

FZPth = 0:pi/100:2*pi; 


% Import SAXS data 

data = h5read('example_dataset.h5','/data');

mask= h5read('example_dataset.h5','/mask'); 

iimax = 137; 


% Colormap for SAXS data 

myCCM = customcolormap([0 0.2 0.5 0.67 0.88 1],...

    {'#880000','#ff0000','#ffff00','#00ff00','#5555bb','#000000'});


% Determine size of SAXS images on detector 

dataXsize = size(data,1);

dataYsize = size(data,2); 

% Crop to central 50% of images from 1/4 to 3/4 in both directions 

cropXstart = round(dataXsize/4); 

cropYstart = round(dataYsize/4); 

cropXend = round(3*dataXsize/4); 

cropYend = round(3*dataYsize/4); 


subplot('position',pos1);


saxsFrame = 0; 

saxsFrameIndex1 = 1:1:136; 

saxsFrameIndex2 = 137:-1:2; 

% Frame index for SAXS images going back and forth between 1 and 137 

saxsFrameIndex = [saxsFrameIndex1 saxsFrameIndex2 saxsFrameIndex1 saxsFrameIndex2 saxsFrameIndex1 saxsFrameIndex2];

% For loop scanning x-direction

for ii = xmovement 

    hold off 

    sph1 = surf(x+aa,1.1*y+bb+ii-0.35,z+cc+0.1,'FaceAlpha',1,'FaceColor',...

        boneYellow,'LineStyle','none','FaceLighting','flat','DiffuseStrength',0.5);

    material dull

    hold on

    sph2 = surf(0.8*x+aa,0.8*y-bb+ii-0.14,0.88*z+cc+1,'FaceAlpha',1,'FaceColor',...

        boneYellow,'LineStyle','none','FaceLighting','flat','DiffuseStrength',0.5);

    material dull

    sph3 = surf(x+aa,y+bb+ii,z-cc,'FaceAlpha',1,'FaceColor',...

        boneYellow,'LineStyle','none','FaceLighting','flat','DiffuseStrength',0.5);

    material dull

    sph4 = surf(x+aa,y-bb+ii,z-cc,'FaceAlpha',1,'FaceColor',...

        boneYellow,'LineStyle','none','FaceLighting','flat','DiffuseStrength',0.5);

    material dull

    

    cyl1 = surf(0.64*x2,1.25*y2+ii,10.4*z2-5,'FaceAlpha',1,'FaceColor',...

        boneYellow,'LineStyle','none','FaceLighting','flat','DiffuseStrength',0.5);

    material dull

    

    cyl2 = surf(basePlateRad*x3,basePlateRad*y3+ii,z3-6,'FaceAlpha',1,'FaceColor',...

        midGrey,'LineStyle','none','FaceLighting','flat','DiffuseStrength',0.5);

    

    ring = fill3(x4, y4+ii, z4-5,midGrey,'LineStyle','none',...

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

    

    % FZP

    for jj = imax:-2:1

        r1 = 0.125*(jj)^0.5;

        r2 = 0.125*(jj-1)^0.5;

        z_circle1 = r1*cos(FZPth);

        y_circle1 = r1*sin(FZPth);

        x_circle1 = 0.0*y_circle1;

        z_circle2 = r2*cos(FZPth);

        y_circle2 = r2*sin(FZPth);

        x_circle2 = 0.0*y_circle2;

        z_circle = [z_circle1 z_circle2];

        y_circle = [y_circle1 y_circle2];

        x_circle = [x_circle1 x_circle2];

        zone = fill3(x_circle+FZPpos,y_circle,z_circle,'k','LineStyle','none');

    end 

    

    % Incident, slightly divergent, radiation cone

    maxFZPrad = imax^0.5*0.125;

    incRad = surf(maxFZPrad*xC1,maxFZPrad*yC1,-FZPpos+incRadLength*zC1,...

        'FaceColor',myGold,'LineStyle',...

        'none','FaceLighting','flat','DiffuseStrength',1);

    rotate(incRad,[0 1 0],-90,[0,0,0]);

    % Vary alpha of divergent radiation in x-direction with offset so it

    % remains semitransparent for all elements of surf = divRad

    alpha(incRad,0.6*incRad.XData/halfDetWidth-2);

    

    % Focussed radiation cone after FZP

    focRad = surf(maxFZPrad*xC2,maxFZPrad*yC2,-FZPpos*zC2,...

        'FaceAlpha',0.25,'FaceColor',myGold,'LineStyle',...

        'none','FaceLighting','flat','DiffuseStrength',1);

    rotate(focRad,[0 1 0],-90,[0,0,0]);

    

    % Divergent radiation cone after focus

    maxFZPrad = imax^0.5*0.125;

    divRad = surf((1+halfDetWidth)*xC2,(1+halfDetWidth)*yC2,(DetPos)*zC2,...

        'FaceColor',myGold,'FaceAlpha',0.28,'LineStyle',...

        'none','FaceLighting','flat','DiffuseStrength',1);

    rotate(divRad,[0 1 0],90,[0,0,0]);

    % Vary alpha of divergent radiation in x-direction with offset so it

    % remains semitransparent for all elements of surf = divRad

    alpha(divRad,-0.5*divRad.XData/halfDetWidth-1.6);

    

    % Detector

    [V1,F1] = platonic_solid(2,0.88); % Cube

    V1(:,1) = 0.1*V1(:,1)+20.06; % Thin detector 

    V1(:,2) = 2.*halfDetWidth*V1(:,2); % Thin detector 

    V1(:,3) = 2.*halfDetWidth*V1(:,3); % Thin detector 

    patch('Faces',F1,'Vertices',V1,'FaceColor',[0.1 0.1 0.1],...

        'FaceAlpha',1,'EdgeColor','none','FaceLighting','flat',...

        'DiffuseStrength',1,'AmbientStrength',1,'SpecularStrength',0);

    [V1,F1] = platonic_solid(2,0.88); % Cube

    V1(:,1) = V1(:,1)+20.61; % Thin detector 

    V1(:,2) = 2.1*halfDetWidth*V1(:,2); % Thin detector 

    V1(:,3) = 2.1*halfDetWidth*V1(:,3); % Thin detector 

    patch('Faces',F1,'Vertices',V1,'FaceColor',[0.25 0.25 0.25],...

        'FaceAlpha',1,'EdgeColor','none','FaceLighting','flat',...

        'DiffuseStrength',1,'AmbientStrength',1,'SpecularStrength',0);

    

    % Dummy vertical and horizontal lines to stop image resizing

    plot3([22 22],[0 0],[-8 8],...

        'color','w', 'LineWidth', 0.1);

    plot3([22 22],[-4 4],[0 0],...

        'color','w', 'LineWidth', 0.1); 

    

    % SAXS data plotted on detector front 

    saxsFrame = saxsFrame + 1; 

    iii = saxsFrameIndex(saxsFrame); 

    imMax = double(max(max(data(:,:,iii)))); 

    imMin = double(min(min(data(:,:,iii)))); 

    imAdd = 0.0001; 

    im = double(data(cropXstart:cropXend,cropYstart:cropYend,iii)); 

    %im = double(data(:,:,iii)); 

    imNowLin = (im - imMin)./(imMax - imMin);

    %imNow = imAdd+imNowLin.*mask; 

    imNow = imAdd+imNowLin.*mask(cropXstart:cropXend,cropYstart:cropYend); 

    imNowLog = log10(imNow)+255; 

    imMin = min(min(imNowLog));


    XX = [20 20; 20 20];

    YY = [halfDetWidth -halfDetWidth; halfDetWidth -halfDetWidth];

    ZZ = [halfDetWidth halfDetWidth;-halfDetWidth -halfDetWidth]; 

    surf(XX,YY,ZZ,imNowLog,'FaceColor','texturemap','EdgeColor','none');

    colormap(myCCM)

    caxis([imMin 252])

    

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

    

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

    axis equal

    axis off

    frame = getframe(gcf);

    writeVideo(vid,frame);

end % For loop for x-movements



% For loop scanning y-direction

for ii = ymovement 

    hold off     

    sph1 = surf(x+aa,1.1*y+bb-0.35,z+cc+ii+0.1,'FaceAlpha',1,'FaceColor',...

        boneYellow,'LineStyle','none','FaceLighting','flat','DiffuseStrength',0.5);

    material dull

    hold on

    sph2 = surf(0.8*x+aa,0.8*y-bb-0.14,0.88*z+cc+ii+1,'FaceAlpha',1,'FaceColor',...

        boneYellow,'LineStyle','none','FaceLighting','flat','DiffuseStrength',0.5);

    material dull

    sph3 = surf(x+aa,y+bb,z-cc+ii,'FaceAlpha',1,'FaceColor',...

        boneYellow,'LineStyle','none','FaceLighting','flat','DiffuseStrength',0.5);

    material dull

    sph4 = surf(x+aa,y-bb,z-cc+ii,'FaceAlpha',1,'FaceColor',...

        boneYellow,'LineStyle','none','FaceLighting','flat','DiffuseStrength',0.5);

    material dull

    

    cyl1 = surf(0.64*x2,1.25*y2,10.4*z2-5+ii,'FaceAlpha',1,'FaceColor',...

        boneYellow,'LineStyle','none','FaceLighting','flat','DiffuseStrength',0.5);

    material dull

    

    cyl2 = surf(basePlateRad*x3,basePlateRad*y3,z3-6+ii,'FaceAlpha',1,'FaceColor',...

        midGrey,'LineStyle','none','FaceLighting','flat','DiffuseStrength',0.5);

    

    ring = fill3(x4, y4, z4-5+ii,midGrey,'LineStyle','none',...

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

    

    % FZP

    for jj = imax:-2:1

        r1 = 0.125*(jj)^0.5;

        r2 = 0.125*(jj-1)^0.5;

        z_circle1 = r1*cos(FZPth);

        y_circle1 = r1*sin(FZPth);

        x_circle1 = 0.0*y_circle1;

        z_circle2 = r2*cos(FZPth);

        y_circle2 = r2*sin(FZPth);

        x_circle2 = 0.0*y_circle2;

        z_circle = [z_circle1 z_circle2];

        y_circle = [y_circle1 y_circle2];

        x_circle = [x_circle1 x_circle2];

        zone = fill3(x_circle+FZPpos,y_circle,z_circle,'k','LineStyle','none');

    end 

    

    % Incident, slightly divergent, radiation cone

    maxFZPrad = imax^0.5*0.125;

    incRad = surf(maxFZPrad*xC1,maxFZPrad*yC1,-FZPpos+incRadLength*zC1,...

        'FaceColor',myGold,'LineStyle',...

        'none','FaceLighting','flat','DiffuseStrength',1);

    rotate(incRad,[0 1 0],-90,[0,0,0]);

    % Vary alpha of divergent radiation in x-direction with offset so it

    % remains semitransparent for all elements of surf = divRad

    alpha(incRad,0.6*incRad.XData/halfDetWidth-2);

    

    % Focussed radiation cone after FZP

    focRad = surf(maxFZPrad*xC2,maxFZPrad*yC2,-FZPpos*zC2,...

        'FaceAlpha',0.25,'FaceColor',myGold,'LineStyle',...

        'none','FaceLighting','flat','DiffuseStrength',1);

    rotate(focRad,[0 1 0],-90,[0,0,0]);

    

    % Divergent radiation cone after focus

    maxFZPrad = imax^0.5*0.125;

    divRad = surf((1+halfDetWidth)*xC2,(1+halfDetWidth)*yC2,(DetPos)*zC2,...

        'FaceColor',myGold,'FaceAlpha',0.28,'LineStyle',...

        'none','FaceLighting','flat','DiffuseStrength',1);

    rotate(divRad,[0 1 0],90,[0,0,0]);

    % Vary alpha of divergent radiation in x-direction with offset so it

    % remains semitransparent for all elements of surf = divRad

    alpha(divRad,-0.5*divRad.XData/halfDetWidth-1.6);

    

    % Detector

    [V1,F1] = platonic_solid(2,0.88); % Cube

    V1(:,1) = 0.1*V1(:,1)+20.06; % Thin detector 

    V1(:,2) = 2.*halfDetWidth*V1(:,2); % Thin detector 

    V1(:,3) = 2.*halfDetWidth*V1(:,3); % Thin detector 

    patch('Faces',F1,'Vertices',V1,'FaceColor',[0.1 0.1 0.1],...

        'FaceAlpha',1,'EdgeColor','none','FaceLighting','flat',...

        'DiffuseStrength',1,'AmbientStrength',1,'SpecularStrength',0);

    [V1,F1] = platonic_solid(2,0.88); % Cube

    V1(:,1) = V1(:,1)+20.61; % Thin detector 

    V1(:,2) = 2.1*halfDetWidth*V1(:,2); % Thin detector 

    V1(:,3) = 2.1*halfDetWidth*V1(:,3); % Thin detector 

    patch('Faces',F1,'Vertices',V1,'FaceColor',[0.25 0.25 0.25],...

        'FaceAlpha',1,'EdgeColor','none','FaceLighting','flat',...

        'DiffuseStrength',1,'AmbientStrength',1,'SpecularStrength',0);

    

    % Dummy vertical and horizontal lines to stop image resizing

    plot3([22 22],[0 0],[-8 8],...

        'color','w', 'LineWidth', 0.1);

    plot3([22 22],[-4 4],[0 0],...

        'color','w', 'LineWidth', 0.1);

    

    % SAXS data plotted on detector front 

    saxsFrame = saxsFrame + 1; 

    iii = saxsFrameIndex(saxsFrame); 

    imMax = double(max(max(data(:,:,iii)))); 

    imMin = double(min(min(data(:,:,iii)))); 

    imAdd = 0.0001; 

    im = double(data(cropXstart:cropXend,cropYstart:cropYend,iii)); 

    imNowLin = (im - imMin)./(imMax - imMin);

    imNow = imAdd+imNowLin.*mask(cropXstart:cropXend,cropYstart:cropYend); 

    imNowLog = log10(imNow)+255; 

    imMin = min(min(imNowLog));


    XX = [20 20; 20 20];

    YY = [halfDetWidth -halfDetWidth; halfDetWidth -halfDetWidth];

    ZZ = [halfDetWidth halfDetWidth;-halfDetWidth -halfDetWidth]; 

    surf(XX,YY,ZZ,imNowLog,'FaceColor','texturemap','EdgeColor','none');

    colormap(myCCM)

    caxis([imMin 252])


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

    

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

    axis equal

    axis off

    frame = getframe(gcf);

    writeVideo(vid,frame);

end % For loop for y-movements


% For loop scanning theta rotation

for ii = thetamovement 

    hold off 

    sph1 = surf(x+aa,1.1*y+bb-0.35,z+cc+0.1,'FaceAlpha',1,'FaceColor',...

        boneYellow,'LineStyle','none','FaceLighting','flat','DiffuseStrength',0.5); 

    material dull

    rotate(sph1,[0 0 1],ii,[0,0,0]);

    hold on

    sph2 = surf(0.8*x+aa,0.8*y-bb-0.14,0.88*z+cc+1,'FaceAlpha',1,'FaceColor',...

        boneYellow,'LineStyle','none','FaceLighting','flat','DiffuseStrength',0.5);

    material dull

    rotate(sph2,[0 0 1],ii,[0,0,0]);

    sph3 = surf(x+aa,y+bb,z-cc,'FaceAlpha',1,'FaceColor',...

        boneYellow,'LineStyle','none','FaceLighting','flat','DiffuseStrength',0.5);

    material dull

    rotate(sph3,[0 0 1],ii,[0,0,0]);

    sph4 = surf(x+aa,y-bb,z-cc,'FaceAlpha',1,'FaceColor',...

        boneYellow,'LineStyle','none','FaceLighting','flat','DiffuseStrength',0.5);

    material dull

    rotate(sph4,[0 0 1],ii,[0,0,0]);

    

    cyl1 = surf(0.64*x2,1.25*y2,10.4*z2-5,'FaceAlpha',1,'FaceColor',...

        boneYellow,'LineStyle','none','FaceLighting','flat','DiffuseStrength',0.5);

    material dull

    rotate(cyl1,[0 0 1],ii,[0,0,0]);

    

    cyl2 = surf(basePlateRad*x3,basePlateRad*y3,z3-6,'FaceAlpha',1,'FaceColor',...

        midGrey,'LineStyle','none','FaceLighting','flat','DiffuseStrength',0.5);

    

    ring = fill3(x4, y4, z4-5,midGrey,'LineStyle','none',...

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

    

    % FZP

    for jj = imax:-2:1

        hold on

        r1 = 0.125*(jj)^0.5;

        r2 = 0.125*(jj-1)^0.5;

        z_circle1 = r1*cos(FZPth);

        y_circle1 = r1*sin(FZPth);

        x_circle1 = 0.0*y_circle1;

        z_circle2 = r2*cos(FZPth);

        y_circle2 = r2*sin(FZPth);

        x_circle2 = 0.0*y_circle2;

        z_circle = [z_circle1 z_circle2];

        y_circle = [y_circle1 y_circle2];

        x_circle = [x_circle1 x_circle2];

        zone = fill3(x_circle+FZPpos,y_circle,z_circle,'k','LineStyle','none');

    end 

    

    % Incident, slightly divergent, radiation cone

    maxFZPrad = imax^0.5*0.125;

    incRad = surf(maxFZPrad*xC1,maxFZPrad*yC1,-FZPpos+incRadLength*zC1,...

        'FaceColor',myGold,'LineStyle',...

        'none','FaceLighting','flat','DiffuseStrength',1);

    rotate(incRad,[0 1 0],-90,[0,0,0]);

    % Vary alpha of divergent radiation in x-direction with offset so it

    % remains semitransparent for all elements of surf = divRad

    alpha(incRad,0.6*incRad.XData/halfDetWidth-2);

    

    % Focussed radiation cone after FZP

    focRad = surf(maxFZPrad*xC2,maxFZPrad*yC2,-FZPpos*zC2,...

        'FaceAlpha',0.25,'FaceColor',myGold,'LineStyle',...

        'none','FaceLighting','flat','DiffuseStrength',1);

    rotate(focRad,[0 1 0],-90,[0,0,0]);

    

    % Divergent radiation cone after focus

    maxFZPrad = imax^0.5*0.125;

    divRad = surf((1+halfDetWidth)*xC2,(1+halfDetWidth)*yC2,(DetPos)*zC2,...

        'FaceColor',myGold,'FaceAlpha',0.28,'LineStyle',...

        'none','FaceLighting','flat','DiffuseStrength',1);

    rotate(divRad,[0 1 0],90,[0,0,0]);

    % Vary alpha of divergent radiation in x-direction with offset so it

    % remains semitransparent for all elements of surf = divRad

    alpha(divRad,-0.5*divRad.XData/halfDetWidth-1.6);

    

    % Detector

    [V1,F1] = platonic_solid(2,0.88); % Cube

    V1(:,1) = 0.1*V1(:,1)+20.06; % Thin detector 

    V1(:,2) = 2.*halfDetWidth*V1(:,2); % Thin detector 

    V1(:,3) = 2.*halfDetWidth*V1(:,3); % Thin detector 

    patch('Faces',F1,'Vertices',V1,'FaceColor',[0.1 0.1 0.1],...

        'FaceAlpha',1,'EdgeColor','none','FaceLighting','flat',...

        'DiffuseStrength',1,'AmbientStrength',1,'SpecularStrength',0);

    [V1,F1] = platonic_solid(2,0.88); % Cube

    V1(:,1) = V1(:,1)+20.61; % Thin detector 

    V1(:,2) = 2.1*halfDetWidth*V1(:,2); % Thin detector 

    V1(:,3) = 2.1*halfDetWidth*V1(:,3); % Thin detector 

    patch('Faces',F1,'Vertices',V1,'FaceColor',[0.25 0.25 0.25],...

        'FaceAlpha',1,'EdgeColor','none','FaceLighting','flat',...

        'DiffuseStrength',1,'AmbientStrength',1,'SpecularStrength',0);


    % Dummy vertical and horizontal lines to stop image resizing

    plot3([22 22],[0 0],[-8 8],...

        'color','w', 'LineWidth', 0.1);

    plot3([22 22],[-4 4],[0 0],...

        'color','w', 'LineWidth', 0.1);

    

    % SAXS data plotted on detector front 

    saxsFrame = saxsFrame + 1; 

    iii = saxsFrameIndex(saxsFrame); 

    imMax = double(max(max(data(:,:,iii)))); 

    imMin = double(min(min(data(:,:,iii)))); 

    imAdd = 0.0001; 

    im = double(data(cropXstart:cropXend,cropYstart:cropYend,iii)); 

    imNowLin = (im - imMin)./(imMax - imMin);

    imNow = imAdd+imNowLin.*mask(cropXstart:cropXend,cropYstart:cropYend); 

    imNowLog = log10(imNow)+255; 

    imMin = min(min(imNowLog));


    XX = [20 20; 20 20];

    YY = [halfDetWidth -halfDetWidth; halfDetWidth -halfDetWidth];

    ZZ = [halfDetWidth halfDetWidth;-halfDetWidth -halfDetWidth]; 

    surf(XX,YY,ZZ,imNowLog,'FaceColor','texturemap','EdgeColor','none');

    colormap(myCCM)

    caxis([imMin 252])


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

    

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

    axis equal

    axis off

    frame = getframe(gcf);

    writeVideo(vid,frame);

end % For loop for theta rotations 


% For loop scanning phi rotations 

for ii = phimovement 

    hold off 

    sph1 = surf(x+aa,1.1*y+bb-0.35,z+cc+0.1,'FaceAlpha',1,'FaceColor',...

        boneYellow,'LineStyle','none','FaceLighting','flat','DiffuseStrength',0.5); 

    material dull

    rotate(sph1,[0 1 0],ii,[0,0,0]);

    hold on

    sph2 = surf(0.8*x+aa,0.8*y-bb-0.14,0.88*z+cc+1,'FaceAlpha',1,'FaceColor',...

        boneYellow,'LineStyle','none','FaceLighting','flat','DiffuseStrength',0.5);

    material dull

    rotate(sph2,[0 1 0],ii,[0,0,0]);

    sph3 = surf(x+aa,y+bb,z-cc,'FaceAlpha',1,'FaceColor',...

        boneYellow,'LineStyle','none','FaceLighting','flat','DiffuseStrength',0.5);

    material dull

    rotate(sph3,[0 1 0],ii,[0,0,0]);

    sph4 = surf(x+aa,y-bb,z-cc,'FaceAlpha',1,'FaceColor',...

        boneYellow,'LineStyle','none','FaceLighting','flat','DiffuseStrength',0.5);

    material dull

    rotate(sph4,[0 1 0],ii,[0,0,0]);

    

    cyl1 = surf(0.64*x2,1.25*y2,10.4*z2-5,'FaceAlpha',1,'FaceColor',...

        boneYellow,'LineStyle','none','FaceLighting','flat','DiffuseStrength',0.5);

    material dull

    rotate(cyl1,[0 1 0],ii,[0,0,0]);

    

    cyl2 = surf(basePlateRad*x3,basePlateRad*y3,z3-6,'FaceAlpha',1,'FaceColor',...

        midGrey,'LineStyle','none','FaceLighting','flat','DiffuseStrength',0.5);

    rotate(cyl2,[0 1 0],ii,[0,0,0]);

    

    ring = fill3(x4, y4, z4-5,midGrey,'LineStyle','none',...

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

    rotate(ring,[0 1 0],ii,[0,0,0]);

    

    % FZP

    for jj = imax:-2:1

        hold on

        r1 = 0.125*(jj)^0.5;

        r2 = 0.125*(jj-1)^0.5;

        z_circle1 = r1*cos(FZPth);

        y_circle1 = r1*sin(FZPth);

        x_circle1 = 0.0*y_circle1;

        z_circle2 = r2*cos(FZPth);

        y_circle2 = r2*sin(FZPth);

        x_circle2 = 0.0*y_circle2;

        z_circle = [z_circle1 z_circle2];

        y_circle = [y_circle1 y_circle2];

        x_circle = [x_circle1 x_circle2];

        zone = fill3(x_circle+FZPpos,y_circle,z_circle,'k','LineStyle','none');

    end 

    

    % Incident, slightly divergent, radiation cone

    maxFZPrad = imax^0.5*0.125;

    incRad = surf(maxFZPrad*xC1,maxFZPrad*yC1,-FZPpos+incRadLength*zC1,...

        'FaceColor',myGold,'LineStyle',...

        'none','FaceLighting','flat','DiffuseStrength',1);

    rotate(incRad,[0 1 0],-90,[0,0,0]);

    % Vary alpha of divergent radiation in x-direction with offset so it

    % remains semitransparent for all elements of surf = divRad

    alpha(incRad,0.6*incRad.XData/halfDetWidth-2);

    

    % Focussed radiation cone after FZP

    focRad = surf(maxFZPrad*xC2,maxFZPrad*yC2,-FZPpos*zC2,...

        'FaceAlpha',0.25,'FaceColor',myGold,'LineStyle',...

        'none','FaceLighting','flat','DiffuseStrength',1);

    rotate(focRad,[0 1 0],-90,[0,0,0]);

    

    % Divergent radiation cone after focus

    maxFZPrad = imax^0.5*0.125;

    divRad = surf((1+halfDetWidth)*xC2,(1+halfDetWidth)*yC2,(DetPos)*zC2,...

        'FaceColor',myGold,'FaceAlpha',0.28,'LineStyle',...

        'none','FaceLighting','flat','DiffuseStrength',1);

    rotate(divRad,[0 1 0],90,[0,0,0]);

    % Vary alpha of divergent radiation in x-direction with offset so it

    % remains semitransparent for all elements of surf = divRad

    alpha(divRad,-0.5*divRad.XData/halfDetWidth-1.6);

    

    % Detector

    [V1,F1] = platonic_solid(2,0.88); % Cube

    V1(:,1) = 0.1*V1(:,1)+20.06; % Thin detector 

    V1(:,2) = 2.*halfDetWidth*V1(:,2); % Thin detector 

    V1(:,3) = 2.*halfDetWidth*V1(:,3); % Thin detector 

    patch('Faces',F1,'Vertices',V1,'FaceColor',[0.1 0.1 0.1],...

        'FaceAlpha',1,'EdgeColor','none','FaceLighting','flat',...

        'DiffuseStrength',1,'AmbientStrength',1,'SpecularStrength',0);

    [V1,F1] = platonic_solid(2,0.88); % Cube

    V1(:,1) = V1(:,1)+20.61; % Thin detector 

    V1(:,2) = 2.1*halfDetWidth*V1(:,2); % Thin detector 

    V1(:,3) = 2.1*halfDetWidth*V1(:,3); % Thin detector 

    patch('Faces',F1,'Vertices',V1,'FaceColor',[0.25 0.25 0.25],...

        'FaceAlpha',1,'EdgeColor','none','FaceLighting','flat',...

        'DiffuseStrength',1,'AmbientStrength',1,'SpecularStrength',0);

    % Dummy vertical and horizontal lines to stop image resizing

    plot3([22 22],[0 0],[-8 8],...

        'color','w', 'LineWidth', 0.1);

    plot3([22 22],[-4 4],[0 0],...

        'color','w', 'LineWidth', 0.1);  

    

    % SAXS data plotted on detector front 

    saxsFrame = saxsFrame + 1; 

    iii = saxsFrameIndex(saxsFrame); 

    imMax = double(max(max(data(:,:,iii)))); 

    imMin = double(min(min(data(:,:,iii)))); 

    imAdd = 0.0001; 

    im = double(data(cropXstart:cropXend,cropYstart:cropYend,iii)); 

    imNowLin = (im - imMin)./(imMax - imMin);

    imNow = imAdd+imNowLin.*mask(cropXstart:cropXend,cropYstart:cropYend); 

    imNowLog = log10(imNow)+255; 

    imMin = min(min(imNowLog));


    XX = [20 20; 20 20];

    YY = [halfDetWidth -halfDetWidth; halfDetWidth -halfDetWidth];

    ZZ = [halfDetWidth halfDetWidth;-halfDetWidth -halfDetWidth]; 

    surf(XX,YY,ZZ,imNowLog,'FaceColor','texturemap','EdgeColor','none');

    colormap(myCCM)

    caxis([imMin 252])


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

    

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

    axis equal

    axis off

    frame = getframe(gcf);

    writeVideo(vid,frame);

end % For loop for phi rotations 


close(vid);