Skip to content

Commit

Permalink
Merge pull request #84 from molssi-seamm/dev
Browse files Browse the repository at this point in the history
Added access to OpenEye mols and SMILES
  • Loading branch information
seamm authored Nov 21, 2024
2 parents 11acf30 + 116b857 commit caec4a8
Show file tree
Hide file tree
Showing 14 changed files with 649 additions and 53 deletions.
6 changes: 6 additions & 0 deletions HISTORY.rst
Original file line number Diff line number Diff line change
@@ -1,6 +1,12 @@
=======
History
=======
2024.11.21 -- Added access to OpenEye mols and SMILES
* Added the ability to create a system from an OpenEye molecule and vice versa.
This gives access to the OpenEye toolkit for generating conformers, etc.
* Added the ability to create a system from a SMILES string and create the SMILES
string from a system using OpenEye.

2024.11.12 -- Added units to properties in Openbabel molecules
* Any properties on the system and configuration are optionally added to the
Openbabel Molecule object for e.g. writing to an SDF file. This adds the units of
Expand Down
2 changes: 1 addition & 1 deletion Makefile
Original file line number Diff line number Diff line change
Expand Up @@ -63,7 +63,7 @@ format: ## reformat with with yapf and isort
black --extend-exclude '_version.py' $(MODULE) tests

test: ## run tests quickly with the default Python
pytest --disable-warnings -rP tests/
pytest --disable-warnings -rP --openeye tests/

coverage: ## check code coverage quickly with the default Python
pytest -v --cov=$(MODULE) --cov-report term --color=yes tests/
Expand Down
12 changes: 12 additions & 0 deletions molsystem/__init__.py
Original file line number Diff line number Diff line change
Expand Up @@ -5,6 +5,18 @@
A general implementation for molecular and periodic systems.
"""

# Temporary fix for issue with openeye and openbabel conflicting
# Seems to work if we import openeye first
# May not have openeye or openbabel installed, so take care
try:
from openeye import oechem # noqa: F401
except ImportError:
pass
try:
from openbabel import openbabel as ob # noqa: F401
except ImportError:
pass

# Bring up the classes so that they appear to be directly in
# the molsystem package.

Expand Down
2 changes: 2 additions & 0 deletions molsystem/configuration.py
Original file line number Diff line number Diff line change
Expand Up @@ -23,6 +23,7 @@
from .inchi import InChIMixin
from .molfile import MolFileMixin
from .openbabel import OpenBabelMixin
from .openeye import OpenEyeMixin
from .pubchem import PubChemMixin
from .rdkit_ import RDKitMixin
from .pdb import PDBMixin
Expand All @@ -42,6 +43,7 @@ class _Configuration(
SMILESMixin,
TopologyMixin,
OpenBabelMixin,
OpenEyeMixin,
PubChemMixin,
RDKitMixin,
QCSchemaMixin,
Expand Down
93 changes: 71 additions & 22 deletions molsystem/openbabel.py
Original file line number Diff line number Diff line change
Expand Up @@ -61,9 +61,20 @@ def to_OBMol(self, properties=None):
# Set local radical character
ob_mol.AssignSpinMultiplicity(True)

# Add the net charge and spin multiplicity as properties for configurations
pair = ob.OBPairData()

if self.__class__.__name__ == "_Configuration":
pair.SetAttribute("net charge")
pair.SetValue(str(self.charge))
ob_mol.CloneData(pair)

pair.SetAttribute("spin multiplicity")
pair.SetValue(str(self.spin_multiplicity))
ob_mol.CloneData(pair)

if properties == "all":
data = self.properties.get("all", include_system_properties=True)
pair = ob.OBPairData()
for key, value in data.items():
pair.SetAttribute(key)
pair.SetValue(str(value))
Expand All @@ -84,7 +95,29 @@ def to_OBMol(self, properties=None):
def from_OBMol(
self, ob_mol, properties="all", atoms=True, coordinates=True, bonds=True
):
"""Transform an Open Babel molecule into the current object."""
"""Transform an Open Babel molecule into the current object.
Parameters
----------
rdk_mol : rdkit.chem.molecule
The RDKit molecule object
properties : str = "all"
Whether to include all properties or none
atoms : bool = True
Recreate the atoms
coordinates : bool = True
Update the coordinates
bonds : bool = True
Recreate the bonds from the RDKit molecule
Returns
-------
molsystem._Configuration
"""
atnos = []
Xs = []
Ys = []
Expand Down Expand Up @@ -116,9 +149,27 @@ def from_OBMol(
if atoms:
self.clear()

# Get the property data, cast to correct type
data = {}
for item in ob_mol.GetData():
value = item.GetValue()
try:
value = int(value)
except Exception:
try:
value = float(value)
except Exception:
pass
data[item.GetAttribute()] = value

# Check for property items for charge and multiplicity
if self.__class__.__name__ == "_Configuration":
self.charge = ob_mol.GetTotalCharge()
self.spin_multiplicity = ob_mol.GetTotalSpinMultiplicity()
if "net charge" in data:
self.charge = int(data["net charge"])
if "spin multiplicity" in data:
self.spin_multiplicity = int(data["spin multiplicity"])

if atoms:
if any([i != 0.0 for i in qs]):
Expand All @@ -141,26 +192,24 @@ def from_OBMol(

# Record any properties in the database if desired
if properties == "all":
data = ob_mol.GetData()
for item in data:
attribute = item.GetAttribute()
value = item.GetValue()
if not self.properties.exists(attribute):
try:
int(value)
_type = "int"
except Exception:
try:
float(value)
_type = "float"
except Exception:
_type = "str"
self.properties.add(
attribute,
_type,
description="Imported from SDF file",
)
self.properties.put(attribute, value)
for key, value in data.items():
if ",units" not in key and key not in [
"net charge",
"spin multiplicity",
]:
if not self.properties.exists(key):
tmp = key.split("#", maxsplit=1)
if len(tmp) > 1:
units_key = tmp[0] + ",units" + "#" + tmp[1]
else:
units_key = key + ",units"
_type = value.__class__.__name__
if units_key in data:
units = data[units_key]
self.properties.add(key, _type, units=units)
else:
self.properties.add(key, _type)
self.properties.put(key, value)
return self

def coordinates_from_OBMol(self, ob_mol):
Expand Down
184 changes: 184 additions & 0 deletions molsystem/openeye.py
Original file line number Diff line number Diff line change
@@ -0,0 +1,184 @@
# -*- coding: utf-8 -*-

"""Interface to OpenEye OEChem."""

import logging

try:
from openeye import oechem
except ModuleNotFoundError:
oechem_available = False
oechem_licensed = False
else:
oechem_available = True
oechem_licensed = oechem.OEChemIsLicensed()

logger = logging.getLogger(__name__)


def check_openeye_license():
if not oechem_available:
raise RuntimeError(
"OpenEye OEChem is not installed! See "
"https://docs.eyesopen.com/toolkits/python/index.html for detail"
)
if not oechem_licensed:
raise RuntimeError(
"Cannot find a license for OpenEye OEChem! See "
"https://docs.eyesopen.com/toolkits/python/index.html for detail"
)


class OpenEyeMixin:
"""A mixin for handling OpenEye's software via its Python interface."""

def to_OEGraphMol(self, properties=None):
"""Return an OEGraphMol object for the configuration, template, or subset."""
check_openeye_license()
oe_mol = oechem.OEGraphMol()

if self.__class__.__name__ == "_Configuration":
oe_mol.SetIntData("net charge", self.charge)
oe_mol.SetIntData("spin multiplicity", self.spin_multiplicity)

# Create the atoms
oe_atoms = []
if "charge" in self.atoms:
charges = self.atoms.get_column_data("charge")
for atno, xyz, q in zip(
self.atoms.atomic_numbers, self.atoms.coordinates, charges
):
oe_atom = oe_mol.NewAtom(atno)
oe_mol.SetCoords(oe_atom, xyz)
oe_atom.SetPartialCharge(q)
oe_atoms.append(oe_atom)
else:
for atno, xyz in zip(self.atoms.atomic_numbers, self.atoms.coordinates):
oe_atom = oe_mol.NewAtom(atno)
oe_mol.SetCoords(oe_atom, xyz)
oe_atoms.append(oe_atom)

# and bonds
index = {j: i for i, j in enumerate(self.atoms.ids)}
for row in self.bonds.bonds():
oe_mol.NewBond(
oe_atoms[index[row["i"]]], oe_atoms[index[row["j"]]], row["bondorder"]
)

if properties == "all":
data = self.properties.get("all", include_system_properties=True)
for key, value in data.items():
_type = self.properties.type(key)
if _type == "int":
oe_mol.SetIntData(key, value)
elif _type == "float":
oe_mol.SetDoubleData(key, value)
elif _type == "str":
oe_mol.SetStringData(key, value)
else:
raise ValueError(f"Can't handle property of type '{_type}'")

# Units, if any
units = self.properties.units(key)
if units is not None and units != "":
tmp = key.split("#", maxsplit=1)
if len(tmp) > 1:
oe_mol.SetStringData(tmp[0] + ",units" + "#" + tmp[1], units)
else:
oe_mol.SetStringData(key + ",units", units)
return oe_mol

def from_OEMol(
self, oe_mol, properties="all", atoms=True, coordinates=True, bonds=True
):
"""Transform an OpenEye molecule into the current object."""
check_openeye_license()

atnos = []
Xs = []
Ys = []
Zs = []
qs = []
for oe_atom in oe_mol.GetAtoms():
atno = oe_atom.GetAtomicNum()
atnos.append(atno)
x, y, z = oe_mol.GetCoords(oe_atom)
Xs.append(x)
Ys.append(y)
Zs.append(z)
qs.append(oe_atom.GetPartialCharge())

index = {at: i for i, at in enumerate(oe_mol.GetAtoms())}
Is = []
Js = []
BondOrders = []
for oe_bond in oe_mol.GetBonds():
oe_i = oe_bond.GetBgn()
oe_j = oe_bond.GetEnd()
i = index[oe_i]
j = index[oe_j]
bondorder = oe_bond.GetOrder()
Is.append(i)
Js.append(j)
BondOrders.append(bondorder)
logger.debug(f"bond {i} - {j} {bondorder}")

if atoms:
self.clear()

# Get the property data, cast to correct type
data = {}
for tmp in oe_mol.GetDataIter():
tag = tmp.GetTag()
attribute = oechem.OEGetTag(tag)
value = oe_mol.GetData(tag)
data[attribute] = value

# Check for property items for charge and multiplicity
if self.__class__.__name__ == "_Configuration":
if "net charge" in data:
self.charge = int(data["net charge"])
if "spin multiplicity" in data:
self.spin_multiplicity = int(data["spin multiplicity"])

if atoms:
if any([i != 0.0 for i in qs]):
if "formal_charge" not in self.atoms:
self.atoms.add_attribute("charge", coltype="float", default=0)
ids = self.atoms.append(x=Xs, y=Ys, z=Zs, atno=atnos, charge=qs)
else:
ids = self.atoms.append(x=Xs, y=Ys, z=Zs, atno=atnos)
else:
ids = self.atoms.ids

if coordinates:
xyz = [[x, y, z] for x, y, z in zip(Xs, Ys, Zs)]
self.atoms.coordinates = xyz

if atoms or bonds:
i = [ids[x] for x in Is]
j = [ids[x] for x in Js]
self.bonds.append(i=i, j=j, bondorder=BondOrders)

# Record any properties in the database if desired
if properties == "all":
for key, value in data.items():
if ",units" not in key and key not in [
"",
"net charge",
"spin multiplicity",
]:
if not self.properties.exists(key):
tmp = key.split("#", maxsplit=1)
if len(tmp) > 1:
units_key = tmp[0] + ",units" + "#" + tmp[1]
else:
units_key = key + ",units"
_type = value.__class__.__name__
if units_key in data:
units = data[units_key]
self.properties.add(key, _type, units=units)
else:
self.properties.add(key, _type)
self.properties.put(key, value)
return self
Loading

0 comments on commit caec4a8

Please sign in to comment.