Powder pattern examples

Note

The powder packing features are available only to CSD-Materials and CSD-Enterprise users.

Comparing powder patterns

Note that this script makes use of functionality from the cookbook utility module.

#!/usr/bin/env python
#
# This script can be used for any purpose without limitation subject to the
# conditions at http://www.ccdc.cam.ac.uk/Community/Pages/Licences/v2.aspx
#
# This permission notice and the following statement of attribution must be
# included in all copies or substantial portions of this script.
#
# 2015-06-17: created by the Cambridge Crystallographic Data Centre
#

'''
Generate a HTML report with powder pattern comparisons.

Compare a powder pattern from a .xye file with simulations derived from
crystals in another file.

The resulting HTML report and the powder pattern images are written out to the
specified output directory.

This code makes use of matplotlib ( see https://matplotlib.org/) and so may
require this package to be installed in your python distribution
'''
##############################################################################

import sys
import os
import argparse
import math

import matplotlib
matplotlib.use('Agg')
import matplotlib.pyplot as plot

from ccdc.io import CrystalReader
from ccdc.descriptors import CrystalDescriptors
from utilities import Lineplot


##############################################################################

class Runner(argparse.ArgumentParser):
    '''Class for producing the powder similarity HTML report.'''

    def __init__(self):
        '''Define the command line interface.'''

        super(Runner, self).__init__(description=__doc__)
        self.add_argument(
            'xye_file', type=str,
            help='xye format experimental powder'
        )
        self.add_argument(
            'crystal_file', type=str,
            help='file of crystal structures to compare'
        )
        self.add_argument(
            'output_directory', type=str,
            help='output directory for images and html index'
        )
        self.add_argument(
            '-s', '--start', default=0, type=int,
            help='Crystal at which to start'
        )
        self.add_argument(
            '-e', '--end', default=sys.maxsize, type=int,
            help='Crystal at which to end'
        )
        self.args = self.parse_args()
        if not os.path.exists(self.args.xye_file):
            raise argparse.ArgumentError(
                'xye file does not exist: %s' % self.args.xye_file)
        if not os.path.exists(self.args.crystal_file):
            raise argparse.ArgumentError(
                'crystal file does not exist: %s' % self.args.crystal_file)
        if not os.path.exists(self.args.output_directory):
            os.makedirs(self.args.output_directory)

    def make_title(self):
        '''Return the report title.'''
        return 'Comparison of %s with %s' % (os.path.basename(self.args.xye_file),
                                             os.path.basename(self.args.crystal_file))

    def make_plot(self, patt, name, ref=None, file_name=None, thumb=False):
        '''Generate a plot from the powder pattern and return the file name.'''
        if file_name is None:
            file_name = str(os.path.join(self.args.output_directory,
                                         '%s.png' % name))
        else:
            file_name = str(os.path.join(self.args.output_directory,
                                         '%s.png' % file_name))

        # Create the matplotlib plot.
        pl = Lineplot(
            file_name=file_name,
            xlabel=r'2$\theta$',
            ylabel='Intensity',
            title='PXRD %s' % name
        )

        # Add the first powder pattern.
        pl.add_plot(patt.two_theta, patt.intensity, label='Intensity')

        # Add the second powder pattern if provided.
        if ref is not None:
            pl.add_plot(ref.two_theta, ref.intensity, label='Reference', color='red')

        # If this is a thumbnail, configure, write and return file name.
        if thumb:
            pl.title = pl.xlabel = pl.ylabel = ''
            plot.tick_params(
                axis='x', which='both', bottom=False, top=False, labelbottom=False
            )
            plot.tick_params(
                axis='y', which='both', left=False, right=False, labelleft=False
            )
            pl.fig.set_size_inches(3., 1.)
            pl.write()
            return os.path.basename(file_name)

        # At this point it is not a thumbnail, so add axis, title, etc.
        artists = [
            plot.Line2D((0,0),(0,1), color='b'),
        ]
        labels = ['Intensity']
        if patt.tick_marks is not None:
            absent = [tm.two_theta for tm in patt.tick_marks
                                   if tm.is_systematically_absent]
            present = [tm.two_theta for tm in patt.tick_marks
                                    if not tm.is_systematically_absent]
            if absent:
                pl.add_plot(absent, [-250.0] * len(absent), linestyle=' ',
                            marker='|', markersize=10,
                            label='Systematic absences')
                artists.append(
                    plot.Line2D((0,0),(0,1), color='g')
                )
                labels.append('Systematic absences')
            if present:
                pl.add_plot(present, [-250.0]* len(present), linestyle=' ',
                            marker='|', markersize=10, color='pink',
                            label='Allowed reflections')
                artists.append(
                    plot.Line2D((0,0),(1,0), color='pink')
                )
                labels.append('Allowed reflections')
            if absent or present:
                pl.axes.get_yaxis().set_view_interval(-500, 10000, True)
                pl.axes.get_yaxis().set_ticks([1000*i for i in range(11)])
                pl.axes.legend(artists, labels)
        pl.write()
        return os.path.basename(file_name)

    def make_table(self):
        '''Returns a string containing the HTML table for the report.'''

        def do_one(c):
            '''Return a tuple containing the identifier and the powder pattern.'''
            try:
                p = CrystalDescriptors.PowderPattern.from_crystal(c)
            except RuntimeError:
                p = None
            return c.identifier, p

        l = [
            do_one(self.reader[i]) for i in
            range(self.args.start, min(len(self.reader), self.args.end))
        ]

        # Store the: similarity, identifier and powder pattern
        self.powder_sims = [
            (self.powder.similarity(p), i, p) for i, p in l if p is not None
        ]
        self.powder_sims = [
            (s, i, p) for s, i, p in self.powder_sims if not math.isnan(s)
        ]
        self.powder_sims.sort(reverse=True)

        def row(s, i, p):
            '''Return a HTML table row.

            :param s: similarity
            :param i: identifier
            :param p: pattern
            '''
            refcode = '<TD>%s</TD>' % i
            similarity = '<TD align="right">%.4f</TD>' % s
            thumbnail_img = '<IMG src="%s_thumbnail.png"/>' % i
            thumbnail = '<TD><A HREF="#%s_image">%s</A></TD>' % (i, thumbnail_img)
            overlay_img = '<IMG src="%s_overlay_thumbnail.png">' % i
            overlaid = '<TD><A HREF="#%s_overlay_image">%s</A></TD>' % (i, overlay_img)
            return '<TR>%s%s%s%s</TR>' % (refcode, similarity, thumbnail, overlaid)

        ret = [
            '<TABLE border=1>',
            '<TR>',
            '<TH>Refcode</TH><TH>Similarity</TH><TH>Thumbnail</TH><TH>Overlaid</TH>',
            '</TR>\n'.join(row(*x) for x in self.powder_sims),
            '</TABLE>'
        ]
        for s, i, p in self.powder_sims:
            self.make_plot(p, i, file_name='%s_thumbnail' % i, thumb=True)
        for s, i, p in self.powder_sims:
            self.make_plot(p, i, ref=self.powder,
                           file_name='%s_overlay_thumbnail' % i, thumb=True)
        return '\n'.join(ret)

    def make_plots(self):
        '''Create the plots for the HTML report.

        Returns a HTML string containing references to the plots.
        '''
        def make_one(s, i, p):
            '''Make one pair of plots.

            :param s: similarity
            :param i: identifier
            :param p: pattern
            :returns: HTML string containing containing references to the plots
            '''
            l = [
                '<DIV>',
                '<IMG id="%s_image" src="%s.png"/>' % (i, i),
                '<IMG style="float:right" id="%s_overlay_image" src="%s_overlay.png"/>' % (i, i),
                '</DIV>',
            ]
            self.make_plot(p, i)
            self.make_plot(p, i, self.powder, file_name='%s_overlay' % i)
            return '\n'.join(l)
        l = [
            make_one(*x) for x in self.powder_sims
        ]
        return '\n'.join(l)

    def run(self):
        '''Generate the HTML report.'''

        # Define the reference powder pattern and the crystals to be compared.
        self.powder = CrystalDescriptors.PowderPattern.from_xye_file(self.args.xye_file)
        self.reader = CrystalReader(self.args.crystal_file)

        # Create a plot of the reference powder pattern.
        ref_pattern_plot_file = self.make_plot(self.powder,
            'Reference %s' % os.path.basename(self.args.xye_file),
                                              file_name='reference')

        # Create a list containing the content of the HTML report.
        l = [
            '<HTML><HEAD><TITLE>%s</TITLE></HEAD>' % self.make_title(),
            '<BODY>',
            '<H1>%s</H1>' % self.make_title(),
            '<IMG src="%s"/>' % ref_pattern_plot_file,
            '<HR>',
            self.make_table(),
            '<HR>',
            self.make_plots(),
            '</BODY></HTML>\n'
        ]

        # Write out the content of the list to the output index.html file.
        with open(os.path.join(self.args.output_directory, 'index.html'), 'w') as w:
            w.write('\n'.join(l))


##############################################################################

if __name__ == '__main__':
    Runner().run()

Calculating many powder patterns

Note that this script makes use of functionality from the cookbook utility module.

#!/usr/bin/env python
#
# This script can be used for any purpose without limitation subject to the
# conditions at http://www.ccdc.cam.ac.uk/Community/Pages/Licences/v2.aspx
#
# This permission notice and the following statement of attribution must be
# included in all copies or substantial portions of this script.
#
# 2015-06-17: created by the Cambridge Crystallographic Data Centre
#

'''
Create and saves powder patterns from all crystals in a file.

Produces a HTML report in the output directory alongside the powder pattern
files.

This code makes use of matplotlib ( see https://matplotlib.org/) and so may
require this package to be installed in your python distribution
'''

import sys
import os
import argparse

import matplotlib
matplotlib.use('Agg')
import matplotlib.pyplot as plot

from ccdc.io import CrystalReader
from ccdc.descriptors import CrystalDescriptors
from ccdc.diagram import DiagramGenerator
from utilities import Lineplot

##############################################################################


class Runner(argparse.ArgumentParser):
    '''Runs the job.'''

    def __init__(self):
        '''Define the command line interface.'''
        super(Runner, self).__init__(description=__doc__)

        self.add_argument(
            'input_file', type=self.valid_file,
            help='Input file containing crystals.'
        )
        self.add_argument(
            'output_directory', type=self.valid_directory,
            help='Directory to place output files.'
        )
        self.add_argument(
            '-d', '--diagram', action='store_true', default=False,
            help='Whether to generate diagrams.'
        )
        self.add_argument(
            '-x', '--xye-file', action='store_false', default=True,
            help='Whether to write xye files.'
        )
        self.add_argument(
            '-r', '--raw-file', action='store_false', default=True,
            help='Whether to write raw files.'
        )
        self.add_argument(
            '-t', '--table', action='store_true', default=False,
            help='Whether to write a cross-similarity table.'
        )
        self.add_argument(
            '-c', '--count', default=sys.maxsize, type=int,
            help='How many structures to process.'
        )
        self.args = self.parse_args()

        html = [
            '<html><head>',
            '<title>PXRD for %s</title>' % os.path.basename(self.args.input_file),
            '</head><body>'
        ]
        self.diag = DiagramGenerator()
        self.all_patts = dict()
        for i, c in enumerate(CrystalReader(self.args.input_file)):
            if i >= self.args.count:
                break
            try:
                patt = CrystalDescriptors.PowderPattern.from_crystal(c)
                self.all_patts[c.identifier] = patt
            except RuntimeError as e:
                patt = None

            if patt is not None and self.args.xye_file:
                patt.write_xye_file(os.path.join(
                    self.args.output_directory, '%s.xye' % c.identifier))
            if patt is not None and self.args.raw_file:
                patt.write_raw_file(os.path.join(
                    self.args.output_directory, '%s.raw' % c.identifier))

            html.append('<div>')
            if self.args.diagram:
                html.append(
                    '<img src="%s"></img>' % self.make_diagram(c)
                )
            if patt is not None:
                html.append('<a name="%s">' % c.identifier)
                html.append(
                    '<img src="%s"></img>' % self.make_plot(patt, c.identifier)
                )
            else:
                html.append('<p>No powder simulation for %s: %s</p>' % (c.identifier, e))
            html.append('</div>')
        if self.args.table:
            html.append(self.make_table())
        html.append('</body></html>')
        f = open(os.path.join(self.args.output_directory, 'index.html'), 'w')
        f.write('\n'.join(html))
        f.close()

    def make_table(self):
        '''Generates a comparison table of all powder simulations.'''
        names = sorted(self.all_patts.keys())

        def do_row(pattern):
            row = [
                '<tr><th>',
                '<a href="#%s">%s</A>' % (pattern, pattern),
                '</th><td>',
                '</td><td>\n'.join('%.3f' % self.all_patts[pattern].similarity(self.all_patts[pattern])
                          for p in names),
                '</td></tr>'
            ]
            return ''.join(row)

        html = [
            '<table border="1" id="comparison_table">',
            '<tr><th>&nbsp;</th><th>',
            '</th><th>'.join('<a href="#%s">%s</a>' % (pattern, pattern) for pattern in names),
            '</th></tr>',
            '\n'.join(do_row(pattern) for pattern in names),
            '</table>'
        ]
        return '\n'.join(html)

    def make_diagram(self, crystal):
        '''Generates a diagram from the crystal.'''
        img = self.diag.image(crystal)
        fname = os.path.join(self.args.output_directory, '%s_diagram.png' % crystal.identifier)
        if img:
            img.save(fname)
        return os.path.basename(fname)

    def make_plot(self, patt, name):
        '''Generate a plot from the powder pattern.'''
        file_name = str(os.path.join(self.args.output_directory, '%s.png' % name))
        pl = Lineplot(
            file_name=file_name,
            xlabel=r'2$\theta$',
            ylabel='Intensity',
            title='PXRD %s' % name
        )
        pl.add_plot(patt.two_theta, patt.intensity, label='Intensity')
        artists = [
            plot.Line2D((0, 0), (0, 1), color='b'),
        ]
        labels = ['Intensity']
        absent = [tm.two_theta for tm in patt.tick_marks if tm.is_systematically_absent]
        present = [tm.two_theta for tm in patt.tick_marks if not tm.is_systematically_absent]
        if absent:
            pl.add_plot(absent, [-250.0] * len(absent), linestyle=' ',
                        marker='|', markersize=10,
                        label='Systematic absences')
            artists.append(
                plot.Line2D((0, 0), (0, 1), color='g')
            )
            labels.append('Systematic absences')
        if present:
            pl.add_plot(present, [-250.0] * len(present), linestyle=' ',
                        marker='|', markersize=10, color='pink',
                        label='Allowed reflections')
            artists.append(
                plot.Line2D((0, 0), (1, 0), color='pink')
            )
            labels.append('Allowed reflections')
        if absent or present:
            pl.axes.get_yaxis().set_view_interval(-500, 10000, True)
            pl.axes.get_yaxis().set_ticks([1000 * i for i in range(11)])
            pl.axes.legend(artists, labels)
        pl.write()
        return os.path.basename(file_name)

    def valid_file(self, file_name):
        '''Return existence of a file.'''
        if os.path.exists(file_name):
            return file_name
        self.error('File %s does not exist' % file_name)

    def valid_directory(self, path):
        '''Return existence of a directory, after trying to make it if necessary.'''
        if os.path.exists(path) and os.path.isdir(path):
            return path
        if not os.path.exists(path):
            os.makedirs(path)
            return path
        self.error('Directory %s does not exist and could not be created' % path)

##############################################################################


if __name__ == '__main__':
    Runner()