Previous topic

4.1. Topology readers — MDAnalysis.topology

Next topic

4.3. MDAnalysis topology tables

This Page

4.2. Common functions for topology building — MDAnalysis.topology.core

The various topology parsers make use of functions and classes in this module. They are mostly of use to developers.

See also

MDAnalysis.topology.tables for some hard-coded atom information that is used by functions such as guess_atom_type() and guess_atom_mass().

class MDAnalysis.topology.core.Angle(a1, a2, a3)[source]

An angle between three Atom instances. Atom 2 is the apex of the angle

New in version 0.8.

class MDAnalysis.topology.core.Bond(a1, a2, order=None)[source]

A bond between two Atom instances.

Two Bond instances can be compared with the == and != operators. A bond is equal to another if the same atom numbers are connected and they have the same bond order. The ordering of the two atom numbers is ignored as is the fact that a bond was guessed.

The presence of a particular atom can also be queried::
>>> Atom in Bond
True / False
is_guessed[source]

True if the bond was guessed.

See also

guess_bonds()

length()[source]

Length of the bond.

partner(Atom)[source]
Returns :the other Atom in this bond
class MDAnalysis.topology.core.TopologyDict(toptype, members)[source]

A customised dictionary designed for sorting the bonds, angles and torsions present in a group of atoms.

Usage:

topologydict = TopologyDict(topologytype, atoms)
Arguments :
topologytype

one of ‘bond’ ‘angle’ or ‘torsion’; a single TopologyDict can only handle one type of topology.

atoms

a list of MDAnalysis.core.AtomGroup.Atom objects.

Returns :
topologydict

A dictionary of the selected topology type

TopologyDicts are also built lazily from a MDAnalysis.core.AtomGroup.AtomGroup using the MDAnalysis.core.AtomGroup.bondDict() MDAnalysis.core.AtomGroup.angleDict() or MDAnalysis.core.AtomGroup.torsionDict() methods.

The TopologyDict collects all the selected topology type from the atoms and categorises them according to the types of the atoms within. A TopologyGroup containing all of a given bond type can be made by querying with the appropriate key. The keys to the topologyDict are a tuple of the atom types that the bond represents and can be viewed using the keys() method.

For example, from a system containing pure ethanol

>>> td = u.atoms.bondDict
>>> td.keys()
[('C', 'C'),
 ('C', 'H'),
 ('O', 'H'),
 ('C', 'O')]
>>> td['C', 'O']
< TopologyGroup containing 912 bonds >

Note

The key for a bond is taken from the type attribute of the atoms.

Getting and setting types of bonds is done smartly, so a C-C-H angle is considered identical to a H-C-C angle.

Duplicate entries are automatically removed upon creation and combination of different Dicts. This means a bond between atoms 1 and 2 will only ever appear once in a dict despite both atoms 1 and 2 having the bond in their .bond attribute.

Two TopologyDicts can be combined using addition, this will not create any duplicate bonds in the process.

New in version 0.8.

keys()[source]

Returns a list of the different types of available bonds

class MDAnalysis.topology.core.TopologyGroup(bondlist)[source]

A container for a group of bonds (either bonds, angles or torsions):

tg = atomGroup.selectBonds(key)
tg = TopologyDict[key]

key describes the desired bond as a tuple of the involved atom types (as defined by the .type atom attribute). A list of available topology keys can be displayed using the .keys() method.

The TopologyGroup contains MDAnalysis.core.AtomGroup.AtomGroup instances which correspond to the components of the bonds the TopologyGroup contains. Ie a bond has 2 AtomGroups whereas an angle has 3.

The bonds(), angles() and torsions() methods offer a “shortcut” to the Cython distance calculation functions in MDAnalysis.core.distances.

TopologyGroups can be combined with TopologyGroups of the same bond type (ie can combine two angle containing TopologyGroups).

TopologyGroups can be indexed to return a single Bond Angle or Torsion

tg[0], tg[-2]

Or sliced to return a TopologyGroup containing a subset of the original:

tg[4:-4]

New in version 0.8.

angles(result=None, pbc=False)[source]

Calculates the angle in radians formed between a bond between atoms 1 and 2 and a bond between atoms 2 & 3

Keywords :
result

allows a predefined results array to be used, note that this will be overwritten

pbc

apply periodic boundary conditions when calculating angles [False] this is important when connecting vectors between atoms might require minimum image convention

Uses cython implementation

Changed in version 0.8.2: Added pbc option (default False)

bonds(pbc=False, result=None)[source]

Calculates the distance between all bonds in this TopologyGroup

Keywords :
pbc

apply periodic boundary conditions when calculating distance [False]

result

allows a predefined results array to be used, note that this will be overwritten

Uses cython implementation

torsions(result=None, pbc=False)[source]

Calculate the torsional angle in radians for this topology group.

Defined as the angle between a plane formed by atoms 1, 2 and 3 and a plane formed by atoms 2, 3 and 4.

Keywords :
result

allows a predefined results array to be used, note that this will be overwritten

pbc

apply periodic boundary conditions when calculating angles [False] this is important when connecting vectors between atoms might require minimum image convention

Uses cython implementation.

Changed in version 0.8.2: Added pbc option (default False)

class MDAnalysis.topology.core.Torsion(a1, a2, a3, a4)[source]

Torsion (dihedral angle) between four Atom instances.

The torsion is defined as the angle between the planes formed by (1, 2, 3) and (2, 3, 4).

New in version 0.8.

MDAnalysis.topology.core.build_bondlists(atoms, bonds=None, angles=None, torsions=None)[source]

Construct the topology lists of each Atom.

The lists are stored in the attributes MDAnalysis.core.AtomGroup.Atom.bonds MDAnalysis.core.AtomGroup.Atom.angles MDAnalysis.core.AtomGroup.Atom.torsions and consist of a list of Bond Angle and Torsion instances respectively

Changed in version 0.8.

MDAnalysis.topology.core.build_residues(atoms)[source]

Create a list Residue instances from a list of Atom instances.

Updating residues also changes the underlying Atom instances, which record to which residue an atom belongs.

Returns :List of Residue instances

New in version 0.8.

MDAnalysis.topology.core.build_segments(atoms)[source]

Create all Segment instancess from a list of Atom instances.

The function also builds the Residue instances by tracking residue numbers.

Updating segments also changes the underlying Atom instances, which record to which residue and segment an atom belongs.

Returns :structure dict, which associates a segname with a Segment
MDAnalysis.topology.core.get_atom_mass(element)[source]

Return the atomic mass in u for element.

Masses are looked up in MDAnalysis.topology.tables.masses.

Warning

Unknown masses are set to 0.

MDAnalysis.topology.core.get_parser_for(filename, permissive=False, bonds=False, format=None)[source]

Return the appropriate topology parser for filename.

Automatic detection is disabled when an explicit format is provided.

MDAnalysis.topology.core.guess_atom_charge(atomname)[source]

Guess atom charge from the name.

Warning

Not implemented; simply returns 0.

MDAnalysis.topology.core.guess_atom_element(atomname)[source]

Guess the element of the atom from the name.

Looks in dict to see if element is found, otherwise it uses the first character in the atomname. The table comes from CHARMM and AMBER atom types, where the first character is not sufficient to determine the atom type. Some GROMOS ions have also been added.

See also

guess_atom_type() and MDAnalysis.topology.tables (where the data are stored)

MDAnalysis.topology.core.guess_atom_mass(atomname)[source]

Guess a mass based on the atom name.

guess_atom_element() is used to determine the kind of atom.

Warning

Anything not recognized is simply set to 0; if you rely on the masses you might want to double check.

MDAnalysis.topology.core.guess_atom_type(atomname)[source]

Guess atom type from the name.

At the moment, this function simply returns the element, as guessed by guess_atom_element().

MDAnalysis.topology.core.guess_bonds(atoms, coords, fudge_factor=0.72, vdwradii=None, lower_bound=0.1)[source]

Guess if bonds exist between two atoms based on their distance.

Bond between two atoms is created, if the two atoms are within

\[\begin{split}d < f * (R_1 + R_2)\end{split}\]

of each other, where \(R_1\) and \(R_2\) are the VdW radii of the atoms and \(f\) is an ad-hoc fudge_factor. This is the same algorithm that VMD uses.

The table of van der Waals radii is hard-coded as MDAnalysis.topology.tables.vdwradii or a user-specified table can be provided as a dictionary in the keyword argument vdwradii. Atoms are found by their Atom.type.

lower_bound defines a heuristic cutoff below which a bond is too short to exist. This is useful for parsing PDB with altloc records where atoms with altloc A and B maybe very close together and there should be no chemical bond between them.

Warning

No check is done after the bonds are guessed to see if Lewis structure is correct. This is wrong and will burn somebody.

The code is also in pure python now, so it’s slow.

New in version 0.7.7.

MDAnalysis.topology.core.guess_format(filename, format=None)[source]

Returns the type of topology file filename.

The current heuristic simply looks at the filename extension but more complicated probes could be implemented here or in the individual packages (e.g. as static methods).

If format is supplied then it overrides the auto detection.