11.3.1. Atom selection Hierarchy — MDAnalysis.core.selection

This module contains objects that represent selections. They are constructed and then applied to the group.

In general, Parser.parse() creates a Selection object from a selection string. This Selection object is then passed an AtomGroup through its apply() method to apply the Selection to the AtomGroup.

This is all invisible to the user through the select_atoms() method of an AtomGroup.

class MDAnalysis.core.selection.AltlocSelection(parser, tokens)[source]

Select atoms based on ‘altLoc’ attribute

class MDAnalysis.core.selection.AroundSelection(parser, tokens)[source]
class MDAnalysis.core.selection.AtomICodeSelection(parser, tokens)[source]

Select atoms based on icode attribute

class MDAnalysis.core.selection.AtomNameSelection(parser, tokens)[source]

Select atoms based on ‘names’ attribute

class MDAnalysis.core.selection.AtomTypeSelection(parser, tokens)[source]

Select atoms based on ‘types’ attribute

class MDAnalysis.core.selection.BackboneSelection(parser, tokens)[source]

A BackboneSelection contains all atoms with name ‘N’, ‘CA’, ‘C’, ‘O’.

This excludes OT* on C-termini (which are included by, eg VMD’s backbone selection).

Changed in version 1.0.1: bb_atoms changed to set (from numpy array) performance improved by ~100x on larger systems

class MDAnalysis.core.selection.BaseSelection(parser, tokens)[source]

Selection of atoms in nucleobases.

Recognized atom names (from CHARMM):

‘N9’, ‘N7’, ‘C8’, ‘C5’, ‘C4’, ‘N3’, ‘C2’, ‘N1’, ‘C6’, ‘O6’,’N2’,’N6’, ‘O2’,’N4’,’O4’,’C5M’

Changed in version 1.0.1: base_atoms changed to set (from numpy array) performance improved by ~100x on larger systems

class MDAnalysis.core.selection.ByResSelection(parser, tokens)[source]

Selects all atoms that are in the same segment and residue as selection

Changed in version 1.0.0: Use "resindices" instead of "resids" (see #2669 and #2672)

class MDAnalysis.core.selection.DistanceSelection[source]

Base class for distance search based selections

validate_dimensions(dimensions)[source]

Check if the system is periodic in all three-dimensions.

Parameters:dimensions (numpy.ndarray) – 6-item array denoting system size and angles
Returns:Returns argument dimensions if system is periodic in all three-dimensions, otherwise returns None
Return type:None or numpy.ndarray
class MDAnalysis.core.selection.MoleculeTypeSelection(parser, tokens)[source]

Select atoms based on ‘moltypes’ attribute

class MDAnalysis.core.selection.NucleicBackboneSelection(parser, tokens)[source]

Contains all atoms with name “P”, “C5’”, C3’”, “O3’”, “O5’”.

These atoms are only recognized if they are in a residue matched by the NucleicSelection.

Changed in version 1.0.1: bb_atoms changed to set (from numpy array) performance improved by ~100x on larger systems

class MDAnalysis.core.selection.NucleicSelection(parser, tokens)[source]

All atoms in nucleic acid residues with recognized residue names.

Recognized residue names:

  • from the CHARMM force field ::
    awk ‘/RESI/ {printf “’”’”%s”’”’,”,$2 }’ top_all27_prot_na.rtf
  • recognized: ‘ADE’, ‘URA’, ‘CYT’, ‘GUA’, ‘THY’
  • recognized (CHARMM in Gromacs): ‘DA’, ‘DU’, ‘DC’, ‘DG’, ‘DT’

Changed in version 0.8: additional Gromacs selections

Changed in version 1.0.1: nucl_res changed to set (from numpy array) performance improved by ~100x on larger systems

class MDAnalysis.core.selection.NucleicSugarSelection(parser, tokens)[source]

Contains all atoms with name C1’, C2’, C3’, C4’, O2’, O4’, O3’.

Changed in version 1.0.1: sug_atoms changed to set (from numpy array) performance improved by ~100x on larger systems

class MDAnalysis.core.selection.PointSelection(parser, tokens)[source]
class MDAnalysis.core.selection.PropertySelection(parser, tokens)[source]

Some of the possible properties: x, y, z, radius, mass,

Possible splitting around operator:

prop x < 5 prop x< 5 prop x <5 prop x<5

class MDAnalysis.core.selection.ProteinSelection(parser, tokens)[source]

Consists of all residues with recognized residue names.

Recognized residue names in ProteinSelection.prot_res.

  • from the CHARMM force field::
    awk ‘/RESI/ {printf “’”’”%s”’”’,”,$2 }’ top_all27_prot_lipid.rtf
  • manually added special CHARMM, OPLS/AA and Amber residue names.

Changed in version 1.0.1: prot_res changed to set (from numpy array) performance improved by ~100x on larger systems

class MDAnalysis.core.selection.RecordTypeSelection(parser, tokens)[source]

Select atoms based on ‘record_type’ attribute

class MDAnalysis.core.selection.ResidSelection(parser, tokens)[source]

Select atoms based on numerical fields

Allows the use of ‘:’ and ‘-‘ to specify a range of values For example

resid 1:10
class MDAnalysis.core.selection.ResidueNameSelection(parser, tokens)[source]

Select atoms based on ‘resnames’ attribute

class MDAnalysis.core.selection.SameSelection(parser, tokens)[source]

Selects all atoms that have the same subkeyword value as any atom in selection

Changed in version 1.0.0: Map "residue" to "resindices" and "segment" to "segindices" (see #2669 and #2672)

class MDAnalysis.core.selection.SegmentNameSelection(parser, tokens)[source]

Select atoms based on ‘segids’ attribute

class MDAnalysis.core.selection.SelectionParser[source]

A small parser for selection expressions. Demonstration of recursive descent parsing using Precedence climbing (see http://www.engr.mun.ca/~theo/Misc/exp_parsing.htm). Transforms expressions into nested Selection tree.

For reference, the grammar that we parse is

E(xpression)--> Exp(0)
Exp(p) -->      P {B Exp(q)}
P -->           U Exp(q) | "(" E ")" | v
B(inary) -->    "and" | "or"
U(nary) -->     "not"
T(erms) -->     segid [value]
                | resname [value]
                | resid [value]
                | name [value]
                | type [value]
expect(token)[source]

Anticipate and remove a given token

parse(selectstr, selgroups, periodic=None)[source]

Create a Selection object from a string.

Parameters:
  • selectstr (str) – The string that describes the selection
  • selgroups (AtomGroups) – AtomGroups to be used in group selections
  • periodic (bool, optional) – for distance based selections, whether to consider periodic boundary conditions
Returns:

  • The appropriate Selection object. Use the .apply method on
  • this to perform the selection.

Raises:

SelectionError – If anything goes wrong in creating the Selection object.

class MDAnalysis.core.selection.SphericalLayerSelection(parser, tokens)[source]
class MDAnalysis.core.selection.SphericalZoneSelection(parser, tokens)[source]
class MDAnalysis.core.selection.StringSelection(parser, tokens)[source]
MDAnalysis.core.selection.grab_not_keywords(tokens)[source]

Pop tokens from the left until you hit a keyword

Parameters:tokens (collections.deque) – deque of strings, some tokens some not
Returns:values – All non keywords found until a keyword was hit
Return type:list of strings

Note

This function pops the values from the deque

Examples

grab_not_keywords([‘H’, ‘and’,’resname’, ‘MET’]) >>> [‘H’]

grab_not_keywords([‘H’, ‘Ca’, ‘N’, ‘and’,’resname’, ‘MET’]) >>> [‘H’, ‘Ca’ ,’N’]

grab_not_keywords([‘and’,’resname’, ‘MET’]) >>> []

MDAnalysis.core.selection.is_keyword(val)[source]

Is val a selection keyword?

Returns False on any of the following strings:
  • keys in SELECTIONDICT (tokens from Selection objects)
  • keys in OPERATIONS (tokens from LogicOperations)
  • (Parentheses)
  • The value None (used as EOF in selection strings)
MDAnalysis.core.selection.return_empty_on_apply(func)[source]

Decorator to return empty AtomGroups from the apply() function without evaluating it