Reciprocal-space 3D plot of SSX data for first 1000 patterns  
Matlab code 

% Program to plot out successive randomly oriented diffraction patterns

% recorded in an SSX experiment on lysozyme using a data file exptable.dat

% that was parsed from the original data file all.dat using the matlab

% routine parseDataFile.m

 

clear; close all;

 

vid = VideoWriter('serialMXdataCollationMovie.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]);

set(gca,'linewidth',7);

 

% Open data file. This has the format:

%

% Chunk 1   h    k    l          I

% -38    7    1     -27.51

% -37    0    4      24.75  ...

% ...

% Chunk 2   h    k    l          I

% -32  -20   -4      47.25

% -30  -24    0      46.00  ...

% etc etc

%

% The program finds lines including the character 'h' and uses the data in

% between for each diffraction pattern

 

fid = fopen('exptable.dat','r');

S = textscan(fid,'%s','delimiter','\n');

SS = S{1};

fclose(fid);

% Get lines with 'h' in them - these are headers of data chunks

idx = find(contains(SS,'h')); % A list of line numbers for the headers

lengthFile = size(SS); % Number of lines in data file

ff = lengthFile(end,1);

numChunks = str2num(cell2mat(SS(ff)));

 

boxSize = 39.9999;

colormap jet;

% Quasitetragonal u.c. lattice constants of lysozyme in nm

abax = 7.86; % a = b

cax = 3.79; % c

cdab = cax/abax; % a/c

 

for i = 1:1000 % Number of diffraction patterns to plot. Gets REALLY slow for large i!!

    clear Pmat

    clear dumArray

    startInt = idx(i);

    endInt = idx(i+1);

    

    Pmat = zeros(endInt-startInt-1,4); % e.g. 5 to 500, 495 elements

    for j = endInt-startInt-1:-1:1

        Pmat(j,:) = str2num(cell2mat(SS(j+startInt))); % After much todo, at last a line of numbers

        % Determine if data point lies outside Q = 38 in a* and b* units

        h = Pmat(j,1);

        k = Pmat(j,2);

        l = Pmat(j,3);

        hkl = (h^2 + k^2 + (l/cdab)^2)^0.5;

        if (hkl >= 38) % Remove data above maximum plotted Q-value

            Pmat(j,:) = [];

        elseif (Pmat(j,4) <= 0) % Negative intensity data removed

            Pmat(j,:) = [];

        end

    end

    

    dumArray = [boxSize-2 0 0 0.00001; -(boxSize-2) 0 0 0.00001; ...

        0 boxSize-2 0 0.00001; 0 -(boxSize-2) 0 0.00001; ...

        0 0 (boxSize-2)*cdab 0.00001; 0 0 -(boxSize-2)*cdab 0.00001];

    circSize = 0.16*(0.000001+(Pmat(:,4)).^0.5);

    colorGrad = (Pmat(:,4)).^0.1; % Map colormap to intensity range of data

    if (i==1) % Start with an empty volume in reciprocal space

        scatter3(dumArray(:,1), dumArray(:,2), dumArray(:,3), dumArray(:,4), ...

            [0 0 0], 'MarkerFaceColor','Flat');

        xlim([-boxSize boxSize]);

        ylim([-boxSize boxSize]);

        zlim(cdab*[-boxSize boxSize]);

        % Drawn axes through origin. Matlab doesn't support this natively.

        % Comment out next seven lines if you want box axes

        line([0,0], [0,0], cdab*[-1.05*boxSize,1.05*boxSize], 'LineWidth', 2, 'Color', 'k');

        line([0,0], [-1.05*boxSize,1.05*boxSize], [0,0], 'LineWidth', 2, 'Color', 'k');

        line([-1.05*boxSize,1.05*boxSize], [0,0], [0,0], 'LineWidth', 2, 'Color', 'k');

        htext = text(-1.07*boxSize,0,0,'h','FontSize',14);

        ktext = text(0,-1.07*boxSize,0,'k','FontSize',14);

        ltext = text(-0.35,0,1.1*boxSize*cdab,'l','FontSize',14);

        

        axis square

        axis off % Comment out if you want box axes

        set(gca,'XTickLabel',[]);

        set(gca,'YTickLabel',[]);

        set(gca,'ZTickLabel',[]);

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

        hold on

        

        pauseMov(vid,100);

    end

    scatter3(Pmat(:,1), Pmat(:,2), Pmat(:,3), circSize, colorGrad, 'MarkerFaceColor','Flat');

    LineHandle.NodeChildren.LineWidth = 0; % Set linewidth of outline of circles in scatter plot to zero

    

    cb = colorbar;

    set(cb,'position',[.35 .25 .01 .2],'TickLength',0.0); % Size and position of colorbar

    cb.Color = [1 1 1];

    

    % Comment next seven lines out for box axes

    line([0,0], [0,0], cdab*[-1.05*boxSize,1.05*boxSize], 'LineWidth', 2, 'Color', 'k');

    line([0,0], [-1.05*boxSize,1.05*boxSize], [0,0], 'LineWidth', 2, 'Color', 'k');

    line([-1.05*boxSize,1.05*boxSize], [0,0], [0,0], 'LineWidth', 2, 'Color', 'k');

    delete(htext); delete(ktext); delete(ltext)

    htext = text(-1.07*boxSize,0,0,'h','FontSize',14);

    ktext = text(0,-1.07*boxSize,0,'k','FontSize',14);

    ltext = text(-0.35,0,1.1*boxSize*cdab,'l','FontSize',14);

    

    axis square % Keep scaling 1:1:1 for x,y,z

    axis off % Comment out for box axes

    set(gca,'XTickLabel',[]);

    set(gca,'YTickLabel',[]);

    set(gca,'ZTickLabel',[]);

    

    if (i > 1)

        delete(hText);

    end

    str1 = 'Capture ';

    str2 = num2str(i,'% 3i');

    strTot = [str1,str2];

    hText = text(-1*boxSize, 0*boxSize, -0.82*boxSize/2,strTot,'FontSize',14);

    

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

    hold on

    xlim([-boxSize boxSize]);

    ylim([-boxSize boxSize]);

    zlim(cdab*[-boxSize boxSize]);

    

    %Store the frame

    frame = getframe(gcf);

    writeVideo(vid,frame);

    if (i <= 5)

        pauseMov(vid,15);

    end

    if (i >= 6) && (i <= 15)

        pauseMov(vid,7);

    end

    if (i >= 16) && (i <= 30)

        pauseMov(vid,4);

    end

    if (i >= 31) && (i <= 60)

        pauseMov(vid,3);

    end

    clear Pmat

    clear dumArray

end

delete(hText)

delete(htext)

delete(ktext)

delete(ltext)

[caz,cel] = view;

for i = caz:-caz/120:0 %120 default

    azimuth = i;

    elevation = cel*(i/caz);

    set(gca,'View',[azimuth,elevation]);

    frame = getframe(gcf);

    writeVideo(vid,frame);

end

 

pauseMov(vid,45);

% Output the movie as an mpg file

close(vid);


% Parsing program of huge Crystfel data file to extract only h k l

% and intensity up to a certain h k l, ignoring negative intensities.

% Note: this program will produced a parsed file with a size in the range

% of 100s of MB!! It takes a couple of minutes to generate on my 2018

% macbook pro, too!

 

clear; close all;

Str = fileread('all.dat');

CStr = strsplit(Str, newline);

N = size(CStr,2);

 

fileID = fopen('exptable.dat','w');

fileChunk=0;

printToFile = 0;

for i = 1:N

    % Find start of data chunk

    lineId = contains(CStr(i), 'h    k    l'); % Logical identifier of lines that start new data chunk

    if (lineId == logical(1)) % i.e. YES! this line contains 'h k l'!!

        fileChunk = fileChunk+1;

        %CStr(i)

        fprintf(fileID,'Chunk %i',fileChunk);

        printToFile = 1; % Allows next line to be printed to file

    end

    % Find end of data chunk

    lineId2 = contains(CStr(i), 'End of reflections'); % Line below the end of latest data chunk

    if (lineId2 == logical(1)) % i.e. YES! this line contains 'End of reflections'!!

        %CStr(i-1)

        printToFile = 0; % Stops next line to be printed to file

    end

    % Writing to file part of loop

    if (printToFile == 1)

        thisLine = cell2mat(CStr(i));

        fprintf(fileID,'%28.28s \n',thisLine);

    end

end

fprintf(fileID,'%i \n',fileChunk);fileChunk % Number of diffraction patterns

fclose(fileID);