4.20. 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().

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.

Changed in version 0.9.0: Now allows resids to be given in non sequential order

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:List of Segment instances

Changed in version 0.9.0: Now allows resids in a given Segment to be given in non sequential order.

Changed in version 0.11.0: Returns a list instead of a dict for consistency with the outputs of build_residues()

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.00

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

Return the appropriate topology parser for filename.

Automatic detection is disabled when an explicit format is provided.

Raises:
ValueError

If no appropriate parser could be found.

MDAnalysis.topology.core.guess_angles(bonds)[source]

Given a list of Bonds, find all angles that exist between atoms.

Works by assuming that if atoms 1 & 2 are bonded, and 2 & 3 are bonded, then (1,2,3) must be an angle.

Returns:List of tuples defining the angles. Suitable for use in u._topology

See also

guess_bonds()

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, **kwargs)[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.

Keywords:
fudge_factor

The factor by which atoms must overlap eachother to be considered a bond. Larger values will increase the number of bonds found. [0.72]

vdwradii

To supply custom vdwradii for atoms in the algorithm. Must be a dict of format {type:radii}. The default table of van der Waals radii is hard-coded as MDAnalysis.topology.tables.vdwradii. Any user defined vdwradii passed as an argument will supercede the table values. [None]

lower_bound

The minimum bond length. All bonds found shorter than this length will be ignored. 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. [0.1]

box

Bonds are found using a distance search, if unit cell information is given, periodic boundary conditions will be considered in the distance search. [None]

Returns:

List of tuples suitable for use in Universe topology building.

Warning

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

Raises:ValueError if inputs are malformed or vdwradii data is missing.

New in version 0.7.7.

Changed in version 0.9.0: Updated method internally to use more numpy, should work faster. Should also use less memory, previously scaled as \(O(n^2)\). vdwradii argument now augments table list rather than replacing entirely.

MDAnalysis.topology.core.guess_dihedrals(angles)[source]

Given a list of Angles, find all dihedrals that exist between atoms.

Works by assuming that if (1,2,3) is an angle, and 3 & 4 are bonded, then (1,2,3,4) must be a dihedral.

Returns:List of tuples defining the dihedrals. Suitable for use in u._topology
MDAnalysis.topology.core.guess_improper_dihedrals(angles)[source]

Given a list of Angles, find all improper dihedrals that exist between atoms.

Works by assuming that if (1,2,3) is an angle, and 2 & 4 are bonded, then (2, 1, 3, 4) must be an improper dihedral. ie the improper dihedral is the angle between the planes formed by (1, 2, 3) and (1, 3, 4)

Returns:List of tuples defining the improper dihedrals. Suitable for use in u._topology