pybel - A Python module that simplifies access to the Open Babel API
Note
The openbabel.py module can be accessed through the ob global variable.
Represent an atom.
The underlying Open Babel OBAtom can be accessed using the attribute:
OBAtom
Atomic mass
Atomic number
Coordinates of the atom
Exact mass
Formal charge
Number of non-hydrogen atoms attached
Number of heteroatoms attached
The hybridization of this atom: 1 for sp, 2 for sp2, 3 for sp3, ...
For further details see OBAtom::GetHyb()
The index of the atom in the molecule (starts at 1)
The maximum number of connections expected for this molecule
The isotope for this atom if specified; 0 otherwise.
Partial charge
Spin multiplicity
Atom type
Number of explicit connections
A molecular fingerprint.
The | operator can be used to calculate the Tanimoto coefficient. For example, given two Fingerprints a and b, the Tanimoto coefficient is given by:
tanimoto = a | b
The underlying fingerprint object can be accessed using the attribute fp.
A list of bits set in the fingerprint
Represent a Pybel Molecule.
An iterator (__iter__()) is provided that iterates over the atoms of the molecule. This allows constructions such as the following:
for atom in mymol:
print atom
Add hydrogens.
A list of atoms of the molecule
Calculate descriptor values.
If descnames is not specified, all available descriptors are calculated.
Calculate a molecular fingerprint.
The charge on the molecule
Conformers of the molecule
Access the molecule’s data through a dictionary-like object, MoleculeData.
Are the coordinates 2D, 3D or 0D?
Create a 2D depiction of the molecule.
Optional parameters:
show – display on screen (default is True)
filename – write to file (default is None)
- update – update the coordinates of the atoms
- This sets the atom coordinates to those determined by the structure diagram generator (default is False)
- usecoords – use the current coordinates
- This causes the current coordinates to be used instead of calculating new 2D coordinates (default is False)
OASA is used for 2D coordinate generation and depiction. Tkinter and Python Imaging Library are required for image display.
The molecule’s energy
The exact mass
The molecular formula
Locally optimize the coordinates.
steps – default is 500
If the molecule does not have any coordinates, make3D() is called before the optimization. Note that the molecule needs to have explicit hydrogens. If not, call addh().
Generate 3D coordinates.
steps – default is 50
Once coordinates are generated, hydrogens are added and a quick local optimization is carried out with 50 steps and the MMFF94 forcefield. Call localopt() if you want to improve the coordinates further.
The molecular weight
Remove hydrogens.
The spin multiplicity
The Smallest Set of Smallest Rings (SSSR)
The molecule title
Access any unit cell data
Write the molecule to a file or return a string.
filename – default is None
If a filename is specified, the result is written to a file. Otherwise, a string is returned containing the result.
To write multiple molecules to the same file you should use the Outputfile class.
Store molecule data in a dictionary-type object
Methods and accessor methods are like those of a dictionary except that the data is retrieved on-the-fly from the underlying OBMol.
Example:
>>> mol = readfile("sdf", 'head.sdf').next()
>>> data = mol.data
>>> print data
{'Comment': 'CORINA 2.61 0041 25.10.2001', 'NSC': '1'}
>>> print len(data), data.keys(), data.has_key("NSC")
2 ['Comment', 'NSC'] True
>>> print data['Comment']
CORINA 2.61 0041 25.10.2001
>>> data['Comment'] = 'This is a new comment'
>>> for k,v in data.iteritems():
... print k, "-->", v
Comment --> This is a new comment
NSC --> 1
>>> del data['NSC']
>>> print len(data), data.keys(), data.has_key("NSC")
1 ['Comment'] False
Represent a file to which output is to be sent.
Although it’s possible to write a single molecule to a file by calling the write() method of a Molecule, if multiple molecules are to be written to the same file you should use the Outputfile class.
filename
Close the output file to further writing.
A Smarts Pattern Matcher
Example:
>>> mol = readstring("smi","CCN(CC)CC") # triethylamine
>>> smarts = Smarts("[#6][#6]") # Matches an ethyl group
>>> print smarts.findall(mol)
[(1, 2), (4, 5), (6, 7)]
The numbers returned are the indices (starting from 1) of the atoms that match the SMARTS pattern. In this case, there are three matches for each of the three ethyl groups in the molecule.
A list of supported descriptors
A list of supported forcefields
A list of supported fingerprint types
A dictionary of supported input formats
A list of supported operations
A dictionary of supported output formats
Iterate over the molecules in a file.
filename
You can access the first molecule in a file using the next() method of the iterator:
mol = readfile("smi", "myfile.smi").next()
You can make a list of the molecules in a file using:
mols = list(readfile("smi", "myfile.smi"))
You can iterate over the molecules in a file as shown in the following code snippet:
>>> atomtotal = 0
>>> for mol in readfile("sdf", "head.sdf"):
... atomtotal += len(mol.atoms)
...
>>> print atomtotal
43