# -*- coding: utf-8 -*-
"""Functions for handling CIF files"""
import io
import logging
import math
import CifFile
logger = logging.getLogger(__name__)
bond_order = {1: 'SING', 2: 'DOUB', 3: 'TRIP'}
[docs]class CIFMixin:
"""A mixin for handling CIF files."""
[docs] def to_cif_text(self, configuration=None):
"""Create the text of a CIF file from the confguration.
Parameters
----------
configuration : int = None
The configuration to use, defaults to the current configuration.
Returns
-------
text : str
The text of the file.
"""
lines = []
atoms = self.atoms
if configuration is None:
configuration = self.current_configuration
# Get the chemical formula
formula, empirical_formula, Z = self.formula(configuration)
formula = ' '.join(formula)
empirical_formula = ''.join(empirical_formula)
# And created the file, line-by-line
lines = []
lines.append('# Generated by MolSSI SEAMM')
lines.append(f"data_{empirical_formula}")
# Cell information
if self.periodicity == 3:
cell = self['cell'].cell(configuration)
a, b, c, alpha, beta, gamma = cell.parameters
volume = cell.volume
lines.append("_symmetry_space_group_name_H-M 'P 1'")
lines.append(f'_cell_length_a {a}')
lines.append(f'_cell_length_b {b}')
lines.append(f'_cell_length_c {c}')
lines.append(f'_cell_angle_alpha {alpha}')
lines.append(f'_cell_angle_beta {beta}')
lines.append(f'_cell_angle_gamma {gamma}')
lines.append('_symmetry_Int_Tables_number 1')
lines.append(f'_cell_volume {volume}')
lines.append(f'_cell_formula_units_Z {Z}')
lines.append('loop_')
lines.append(' _symmetry_equiv_pos_site_id')
lines.append(' _symmetry_equiv_pos_as_xyz')
lines.append(" 1 'x, y, z'")
lines.append(f'_chemical_formula_structural {empirical_formula}')
lines.append(f"_chemical_formula_sum '{formula}'")
# The atoms
lines.append('loop_')
lines.append(' _atom_site_type_symbol')
lines.append(' _atom_site_label')
lines.append(' _atom_site_symmetry_multiplicity')
lines.append(' _atom_site_fract_x')
lines.append(' _atom_site_fract_y')
lines.append(' _atom_site_fract_z')
lines.append(' _atom_site_occupancy')
# Need unique names
if 'names' in atoms:
original_names = atoms.get_column('names', configuration)
else:
original_names = atoms.symbols(configuration)
names = []
tmp = {}
for name in original_names:
if name in tmp:
tmp[name] += 1
names.append(name + str(tmp[name]))
else:
tmp[name] = 1
names.append(name)
UVW = atoms.coordinates(
configuration=configuration, fractionals=True, in_cell='molecule'
)
symbols = atoms.symbols(configuration)
for element, name, uvw in zip(symbols, names, UVW):
u, v, w = uvw
lines.append(f'{element} {name} 1 {u:.3f} {v:.3f} {w:.3f} 1')
# And that is it!
return '\n'.join(lines)
[docs] def from_cif_text(self, text, configuration=None):
"""Create the system from a CIF file..
Parameters
----------
text : str
The text from the CIF file
configuration : int = None
The configuration to use, defaults to the current configuration.
Returns
-------
None
"""
cif = CifFile.ReadCif(io.StringIO(text))
data_blocks = [*cif.keys()]
if len(data_blocks) != 1:
raise RuntimeError(
f'There are {len(data_blocks)} data blocks in the cif file.'
)
data_block = cif[data_blocks[0]]
# Reset the system
self.clear(configuration=configuration)
self.periodicity = 3
self.coordinate_system = 'fractional'
# The cell
a = data_block['_cell_length_a']
b = data_block['_cell_length_b']
c = data_block['_cell_length_c']
alpha = data_block['_cell_angle_alpha']
beta = data_block['_cell_angle_beta']
gamma = data_block['_cell_angle_gamma']
self.cell.set_cell(a, b, c, alpha, beta, gamma)
# Add the atoms
# TEMPORARILY lower the symmetry to P1
delta = 1.0e-04
for x, y, z, symbol in zip(
data_block['_atom_site_fract_x'],
data_block['_atom_site_fract_y'],
data_block['_atom_site_fract_z'],
data_block['_atom_site_type_symbol']
): # yapf: disable
xs = []
ys = []
zs = []
symbols = []
# These variables *are* used in the eval below.
x = float(x)
y = float(y)
z = float(z)
for symop in data_block['_space_group_symop_operation_xyz']:
x_eq, y_eq, z_eq = symop.split(',')
x_new = eval(x_eq)
y_new = eval(y_eq)
z_new = eval(z_eq)
# Translate into cell.
x_new = x_new - math.floor(x_new)
y_new = y_new - math.floor(y_new)
z_new = z_new - math.floor(z_new)
# check for almost 1, should be 0
if abs(1 - x_new) < delta:
x_new = 0.0
if abs(1 - y_new) < delta:
y_new = 0.0
if abs(1 - z_new) < delta:
z_new = 0.0
found = False
for x0, y0, z0 in zip(xs, ys, zs):
if (
abs(x_new - x0) < delta and abs(y_new - y0) < delta and
abs(z_new - z0) < delta
):
found = True
break
if not found:
xs.append(x_new)
ys.append(y_new)
zs.append(z_new)
symbols.append(symbol)
self.atoms.append(x=xs, y=ys, z=zs, symbol=symbols)
# self.atoms.append(
# x=data_block['_atom_site_fract_x'],
# y=data_block['_atom_site_fract_y'],
# z=data_block['_atom_site_fract_z'],
# symbol=data_block['_atom_site_type_symbol']
# )
[docs] def to_mmcif_text(self, configuration=None):
"""Create the text of a mmCIF file from the confguration.
Parameters
----------
configuration : int = None
The configuration to use, defaults to the current configuration.
Returns
-------
text : str
The text of the file.
"""
lines = []
atoms = self.atoms
bonds = self.bonds
if configuration is None:
configuration = self.current_configuration
# Get the chemical formula
formula, empirical_formula, Z = self.formula(configuration)
formula = ' '.join(formula)
empirical_formula = ''.join(empirical_formula)
# And created the file, line-by-line
lines = []
lines.append('# Generated by MolSSI SEAMM')
lines.append(f"data_{empirical_formula}")
lines.append(f"_chem_comp.name '{formula}'")
lines.append(f"_chem_comp.id '{empirical_formula}'")
lines.append(f"_chem_comp.formula '{formula}'")
# Cell information
if self.periodicity == 3:
cell = self['cell'].cell(configuration)
a, b, c, alpha, beta, gamma = cell.parameters
volume = cell.volume
lines.append("_symmetry_space_group_name_H-M 'P 1'")
lines.append(f'_cell_length_a {a}')
lines.append(f'_cell_length_b {b}')
lines.append(f'_cell_length_c {c}')
lines.append(f'_cell_angle_alpha {alpha}')
lines.append(f'_cell_angle_beta {beta}')
lines.append(f'_cell_angle_gamma {gamma}')
lines.append('_symmetry_Int_Tables_number 1')
lines.append(f'_cell_volume {volume}')
lines.append(f'_cell_formula_units_Z {Z}')
lines.append('loop_')
lines.append(' _symmetry_equiv_pos_site_id')
lines.append(' _symmetry_equiv_pos_as_xyz')
lines.append(" 1 'x, y, z'")
lines.append(f'_chemical_formula_structural {empirical_formula}')
lines.append(f"_chemical_formula_sum '{formula}'")
# The atoms
lines.append('loop_')
lines.append(' _chem_comp_atom.comp_id')
lines.append(' _chem_comp_atom.atom_id')
lines.append(' _chem_comp_atom.type_symbol')
lines.append(' _chem_comp_atom.model_Cartn_x')
lines.append(' _chem_comp_atom.model_Cartn_y')
lines.append(' _chem_comp_atom.model_Cartn_z')
lines.append(' _chem_comp_atom.pdbx_model_Cartn_x_ideal')
lines.append(' _chem_comp_atom.pdbx_model_Cartn_y_ideal')
lines.append(' _chem_comp_atom.pdbx_model_Cartn_z_ideal')
lines.append(' _chem_comp_atom.pdbx_component_comp_id')
lines.append(' _chem_comp_atom.pdbx_residue_numbering')
# Need unique names
if 'names' in atoms:
original_names = atoms.get_column('names', configuration)
else:
original_names = atoms.symbols(configuration)
names = []
tmp = {}
for name in original_names:
if name in tmp:
tmp[name] += 1
names.append(name + str(tmp[name]))
else:
tmp[name] = 1
names.append(name)
XYZ = atoms.coordinates(
configuration=configuration, in_cell='molecule'
)
XYZa = atoms.coordinates(
configuration=configuration, fractionals=False
)
symbols = atoms.symbols(configuration)
for element, name, xyza, xyz in zip(symbols, names, XYZa, XYZ):
xa, ya, za = xyza
x, y, z = xyz
lines.append(
f'MOL1 {name} {element} {xa:.3f} {ya:.3f} {za:.3f} '
f'{x:.3f} {y:.3f} {z:.3f} HET 1'
)
# The bonds
lines.append('#')
lines.append('loop_')
lines.append(' _chem_comp_bond.comp_id')
lines.append(' _chem_comp_bond.atom_id_1')
lines.append(' _chem_comp_bond.atom_id_2')
lines.append(' _chem_comp_bond.value_order')
for row in bonds.bonds(configuration=configuration):
i = row['i']
j = row['j']
order = bond_order[row['bondorder']]
nm1 = names[i - 1]
nm2 = names[j - 1]
lines.append(f'MOL1 {nm1} {nm2} {order}')
# And that is it!
return '\n'.join(lines)