The Molecular Entity¶
This document describes the EntityHandle
and the related classes.
The Handle Classes¶
Entity
, residue
,
atom
and bond
handles can store
arbitrary generic properties.
- CreateEntity()¶
Creates a new entity. The created entity is empty, that is, it does not contain any atoms, residues, chains, bonds or torsions. To populate the entity, use an entity editor.
- Returns:
The newly created
EntityHandle
Entity Handle¶
- class EntityHandle¶
The entity class represents a molecular structure. Such a structure is in general made up of one or more chains of residues, which in turn are formed by one or more atoms.
The interface of entities is tailored to biological macromolecules, but this does not prevent it to be used for molecules in general: An entity also represent a ligand or a collection of water molecules - hence the rather generic name.
- properties¶
All the
generic properties
are available.
- chains¶
List of all chains of this entity. The chains are in the same order they have been added, for example the order they have been listed in a PDB file.
Also available as
GetChainList()
. To access a single chain, useFindChain()
.This property is read-only.
- Type:
ChainHandleList
(list ofChainHandle
)
- chain_count¶
Number of chains. Read-only. See
GetChainCount()
.- Type:
int
- residues¶
List of all residues of this entity. The returned residues are first sorted by chain and then from N- to C-terminus.
This property is read-only.
Example usage:
for res in ent.residues: print(res.name, res.atom_count)
Also available as
GetResidueList()
. To access a single residue, useFindResidue()
.- Type:
ResidueHandleList
(list ofResidueHandle
)
- residue_count¶
Number of residues. Read-only. See
GetResidueCount()
.- Type:
int
- atoms¶
Get list of all atoms of this entity. To access a single atom, use
FindAtom()
.This property is read-only. Also available as
GetAtomList()
- Type:
AtomHandleList
(list ofAtomHandle
)
- atom_count¶
Number of atoms. Read-only. See
GetAtomCount()
.- Type:
int
- bounds¶
Axis-aligned bounding box of the entity. Read-only
- Type:
- center_of_atoms¶
Center of atoms, that is the average atom position of the entity. Use
center_of_mass
for the the mass-weighted center of the entity. Also available asGetCenterOfAtoms()
.- Type:
- center_of_mass¶
Center of mass of the entity. Also available as
GetCenterOfMass()
- Type:
- geometric_center¶
Mid-point of the axis aligned bounding box of the entity. Also available as
GetGeometricCenter()
- Type:
- valid¶
Validity of handle.
- Type:
bool
- GetName()¶
- Returns:
Name associated to this entity.
- Return type:
str
- SetName(name)¶
- Parameters:
name (
str
) – Sets this as new name to be associated to this entity.
- FindChain(chain_name)¶
Get chain by name. See also
chains
- Parameters:
chain_name (str) – Chain identifier, e.g. “A”
- Returns:
A valid
ChainHandle
, if the entity contains a chain with the given name, an invalidChainHandle
otherwise.
- GetChainCount()¶
See
chain_count
- FindResidue(chain_name, res_num)¶
Get residue by chain name and residue number. See also
GetResidueList()
- Parameters:
chain_name (str) – Chain identifier, e.g. “A”
res_num (
ResNum
) – residue number
- Returns:
A valid
ResidueHandle
if the chain exists and the chain contains a residue of the given residue number, an invalidResidueHandle
otherwise.
- GetResidueCount()¶
See
residue_count
- FindAtom(chain_name, res_num, atom_name)¶
Get atom by chain name, residue number and atom name. See also
atoms
- Parameters:
chain_name (str) – Chain identifier, e.g. “A”
res_num (
ResNum
) – residue numberatom_name (str) – atom name, e.g. CA
- Returns:
A valid
AtomHandle
if an atom matching the parameters could be found, an invalidAtomHandle
otherwise
- GetAtomCount()¶
See
atom_count
- EditXCS([edit_mode=mol.EditMode.UNBUFFERED_EDIT])¶
Request
XCSEditor
for editing the external coordinate system. This call will fail when there are pending changes of the internal coordinate system.- Parameters:
edit_mode (mol.EditMode) – Must be EditMode.BUFFERED_EDIT or EditMode.UNBUFFERED_EDIT. For more details, see the editor documentation.
- Returns:
- EditICS([edit_mode=mol.EditMode.UNBUFFERED_EDIT])¶
Request
ICSEditor
for editing the internal coordinate system, such as torsions, bond lengths and angle between two bonds. This call will fail when there are pending changes of the external coordinate system.- Parameters:
edit_mode (mol.EditMode) – Must be EditMode.BUFFERED_EDIT or EditMode.UNBUFFERED_EDIT. For more details, see the editor documentation.
- Returns:
- Select(query, flags=0)¶
Perform a selection on the entity. The result of the selection is an
EntityView
which (usually) contains only a subset of chains, residues, atoms and bonds of the original entity.Example Usage:
# select calpha atoms of peptides calphas = ent.Select('aname=CA and peptide=true') # select atoms in a box of size 10, centred at the origin in_box = ent.Select('x=-5:5 and y=-5:5 and z=-5:5')
- Parameters:
- Returns:
An
entity view
.- Raises:
QueryError
when the query could not be executed due to syntactic errors.
- CreateFullView()¶
Creates an entity view containing all chains, residues, atoms and bonds of this entity.
# these two lines are identical full=ent.Select('') full=ent.CreateFullView()
- Returns:
- CreateEmptyView()¶
Creates an entity view pointing to this entity, but otherwise empty. This method is usually the starting point for manual entity view creation, e.g.:
view=ent.CreateEmtpyView() for atom in ent.atoms: if ComplicatedPredicate(atom): view.AddAtom(atom)
- Returns:
- Copy()¶
Creates a deep copy of the entity.
- Returns:
A new
EntityHandle
that is an exact copy of this entity.
Note
alternative atom positions are not handled yet.
- GetCenterOfAtoms()¶
See
center_of_atoms
- GetCenterOfMass()¶
See
center_of_mass
- GetGeometricCenter()¶
See
geometric_center
- FindWithin(pos, radius)¶
Find all atoms in sphere of given radius centred at pos. To turn the returned list of atoms into an
EntityView
, useCreateViewFromAtomList()
.- Parameters:
pos (
Vec3
) – Center of sphereradius (float) – The radius of the sphere
- Returns:
AtomHandleList
(list ofAtomHandle
)
Chain Handle¶
- class ChainHandle¶
A chain of one or more
residues
. Chains are always part of an entity (see note on memory management).- properties¶
All the
generic properties
are available.
- atoms¶
Get list of all atoms of this chain. To access a single atom, use
FindAtom()
.This property is read-only. Also available as
GetAtomList()
- Type:
AtomHandleList
(list ofAtomHandle
)
- bounds¶
Axis-aligned bounding box of the chain. Read-only
- Type:
- center_of_atoms¶
Center of atoms (not mass weighted). Also available as
GetCenterOfAtoms()
.- Type:
- geometric_center¶
Mid-point of the axis aligned bounding box of the chain. Also available as
GetGeometricCenter()
- Type:
- description¶
Details about the chain. Not categorised, just text.
- in_sequence¶
Whether the residue numbers are in ascending order. For example:
chain=ent.FindChain("A") print(chain.residues) # [A.GLY1, A.GLY2, A.GLY4A, A.GLY4B] print(chain.in_sequence) # prints true chain=ent.FindChain("B") print(chain.residues) # [B.GLY1, B.GLY4, B.GLY3] print(chain.in_sequence) # prints false
- is_oligosaccharide¶
Indicates if the chain is an oligosaccharide, a branched, non-linear entity of multiple sugars. Also available as
IsOligosaccharide()
.- Type:
bool
- is_polymer¶
Indicates if a chain is a polymer. True for polypeptides, polynucleotides, polysaccharides, oligosaccharides and branched chains. Also available as
IsPolymer()
.- Type:
bool
- is_polynucleotide¶
Indicates if a chain is a nucleic acid. Also available as
IsPolynucleotide()
.- Type:
bool
- is_polypeptide¶
Indicates if a chain is a protein. Also available as
IsPolypeptide()
.- Type:
bool
- is_polysaccharide¶
Indicates if a chain is a polysaccharide. Also available as
IsPolysaccharide()
.- Type:
bool
- mass¶
The total mass of this chain in Dalton. Also available as
GetMass()
- Type:
float
- name¶
The chain name. The name uniquely identifies the chain in the entity. In most cases, the chain name is one character. This is restriction of the PDB file format. However, you are free to use longer names as long as you don’t want to save them as PDB files
This property is read-only. To change the name, use an
XCSEditor
.Also available as
GetName()
- Type:
str
- residue_count¶
Number of residues. Read-only. See
GetResidueCount()
.- Type:
int
- residues¶
List of all residues of this chain. The residues are sorted from N- to C-terminus. Usually the residue numbers are in ascending order (see
in_sequence
).This property is read-only.
Example usage:
chain=ent.FindChain("A") for res in chain.residues: print(res.name, res.atom_count)
Also available as
GetResidueList()
. To access a single residue, useFindResidue()
.- Type:
ResidueHandleList
(list ofResidueHandle
)
- hash_code¶
A unique identifier for this chain. Note that a deep copy of an entity (see
EntityHandle.Copy()
) will have chains with differing identifiers. Shallow copies of the entity preserve the identifier. Also available asGetHashCode()
.- Type:
int
- valid¶
Validity of handle.
- Type:
bool
- FindResidue(res_num)¶
Get residue by residue number. See also
residues
.- Parameters:
res_num (
ResNum
) – residue number- Returns:
A valid
ResidueHandle
if the chain contains a residue with matching residue number, an invalidResidueHandle
otherwise.
- GetResidueCount()¶
See
residue_count
- FindAtom(res_num, atom_name)¶
Get atom by residue number and atom name. See also
atoms
- Parameters:
res_num (
ResNum
) – residue numberatom_name (str) – atom name, e.g. CA
- Returns:
A valid
AtomHandle
if an atom matching the parameters could be found, an invalidAtomHandle
otherwise
- IsOligosaccharide()¶
- IsPolymer()¶
See
is_polymer
- IsPolynucleotide()¶
- IsPolypeptide()¶
See
is_polypeptide
- IsPolysaccharide()¶
- GetDescription()¶
See
description
Residue Handle¶
- class ResidueHandle¶
The residue is either used to represent complete molecules or building blocks in a polymer, e.g. in a protein, DNA or RNA. A residue consists of one or more
atoms
. Residues are always part of aChainHandle
, even if they are ligands or water molecules where the concept of a chain does not apply, and of an entity (see note on memory management).- properties¶
All the
generic properties
are available.
- name¶
The residue name is usually a str of 3 characters, e.g. GLY for glycine or ALA for alanine, but may be shorter, e.g. G for guanosine, or longer for structures loaded from formats other than PDB. Also available as
GetName()
.This property is read-only. To change the name of the residue, use
RenameResidue()
.- Type:
str
- qualified_name¶
The qualified name consists of a residue identifier and chain name. For a glycine with residue number 2 of chain A, the qualified name is “A.GLY2”. Also available as
GetQualifiedName()
.- Type:
str
- number¶
The number of this residue. The residue number has a numeric part and an insertion-code. This property is read-only. Also available as
GetNumber()
.- Type:
- one_letter_code¶
For amino acids, and nucleotides the one_letter_code is an alpha-numeric character. For unknown or more exotic residues, the one letter code is set to ‘?’.
Example
This code-snippet shows how to get the sequence string from a list of residues.
print(''.join([r.one_letter_code for r in chain.residues]))
- Type:
str
- atoms¶
Get list of all atoms of this residue. To access a single atom, use
FindAtom()
.This property is read-only. Also available as
GetAtomList()
- Type:
AtomHandleList
(list ofAtomHandle
)
- atom_count¶
Number of atoms. Read-only. See
GetAtomCount()
.- Type:
int
- bond_count¶
Number of bonds. Read-only. See
GetBondCount()
.- Type:
int
- bounds¶
Axis-aligned bounding box of the residue. Read-only.
- Type:
- center_of_atoms¶
Center of atoms (not mass weighted). Also available as
GetCenterOfAtoms()
.- Type:
- center_of_mass¶
Center of mass. Also available as
GetCenterOfMass()
- Type:
- chain¶
The chain this residue belongs to. Read-only. Also available as
GetChain()
- Type:
- phi_torsion¶
The PHI dihedral angle between this residue and the next. For residues that are not amino acids, residues that do not have all atoms required or residues that do not have bonds between the four atoms involved in the torsion, the PHI torsion is an invalid handle.
Read-only. Also available as
GetPhiTorsion()
- Type:
TorsionHandle
- psi_torsion¶
The PSI dihedral angle between this residue and the previous. For residues that are not amino acids, residues that do not have all atoms required or residues that do not have bonds between the four atoms involved in the torsion, the PSI torsion is an invalid handle.
Read-only. Also available as
GetPsiTorsion()
- omega_torsion¶
The OMEGA dihedral angle between this residue and the previous. For residues that are not amino acids, residues that do not have all atoms required or residues that do not have bonds between the four atoms involved in the torsion, the OMEGA torsion is an invalid handle.
Read-only. Also available as
GetOmegaTorsion()
- Type:
TorsionHandle
- chem_type¶
The chemical type of the residue. The type is only properly set if a compound library is used.
- Type:
- sec_structure¶
The secondary structure of the residue.
- Type:
- is_ligand¶
Whether the residue is a ligand. When loading PDB structures, this property is set based on the HET records. This also means, that this property will most likely not be set properly for all except PDB files coming from pdb.org. When loading MMCIF structures, this property is set based on seqres information and differs from PDB structures. Also available as
IsLigand()
,SetIsLigand()
.- Type:
bool
- is_protein¶
Whether the residue is considered to be part of a protein. This is set when loading a structure if the residue forms a feasible peptide bond to the previous or next residue (see
IsBondFeasible()
). Also available asIsProtein()
,SetIsProtein()
. In contrast toIsPeptideLinking()
this excludes residues which are not connected to neighbouring residues such as CA-only residues or badly positioned ones.- Type:
bool
- peptide_linking¶
Whether residue can form peptide bonds. This is determined based on
chem_class
which is set when loading the structure.- Type:
bool
- nucleotide_linking¶
Whether residue can form nucleotide bonds. This is determined based on
chem_class
which is set when loading the structure.- Type:
bool
- index¶
Residue index (starting at 0) within chain.
- central_atom¶
Central atom used for rendering traces. For peptides, this is usually the CA atom. For nucleotides, this is usually the P atom.
- Type:
- central_normal¶
Normal computed for
central_atom
. Only defined for peptides and nucleotides if all required atoms available. Otherwise, the (1,0,0) vector is returned.- Type:
- hash_code¶
A unique identifier for this residue. Note that a deep copy of an entity (see
EntityHandle.Copy()
) will have residues with differing identifiers. Shallow copies of the entity preserve the identifier. Also available asGetHashCode()
.- Type:
int
- valid¶
Validity of handle.
- Type:
bool
- next¶
Residue after this one in the same chain. Invalid handle returned if there is no next residue. Residues are ordered as in
ChainHandle.residues
independently on whether they are connected or not (seeInSequence()
to check for connected residues).- Type:
- FindAtom(atom_name)¶
Get atom by atom name. See also
atoms
- Parameters:
atom_name (str) – atom name, e.g. CA
- Returns:
A valid
AtomHandle
if an atom with the given name could be found, an invalidAtomHandle
otherwise
- HasAltAtoms()¶
Returns whether the residue has any alternative atom groups
- Return type:
bool
- HasAltAtomGroup(group_name)¶
Returns whether the residue has a certain alternative atom group
- Parameters:
group_name (
str
) – Group of interest- Return type:
bool
- GetAltAtomGroupNames()¶
Returns names of all alternative atom groups in residue
- Return type:
list
ofstr
- GetCurrentAltGroupName()¶
Returns the currently active alternative atom group, empty
str
if residue has no alternative atom groups- Return type:
str
- SwitchAtomPos(group_name)¶
Switches the atoms in residue to the specified alternative atom group
- Parameters:
group_name (
str
) – Group of interest- Return type:
bool
- Returns:
Whether the switch was successful (e.g. False if no such group exists)
- GetQualifiedName()¶
See
qualified_name
- GetOneLetterCode()¶
See
one_letter_code
- GetAtomCount()¶
See
atom_count
- GetBondCount()¶
See
bond_count
- GetCenterOfAtoms()¶
See
center_of_atoms
- GetCenterOfMass()¶
See
center_of_mass
- GetGeometricCenter()¶
See
geometric_center
- GetChemClass()¶
- SetChemClass()¶
See
chem_class
- GetPhiTorsion()¶
See
phi_torsion
- GetPsiTorsion()¶
See
psi_torsion
- GetOmegaTorsion()¶
See
omega_torsion
- GetSecStructure()¶
See
sec_structure
- IsProtein()¶
See
is_protein
- IsPeptideLinking()¶
See
peptide_linking
- IsNucleotideLinking()¶
- GetCentralAtom()¶
- SetCentralAtom()¶
See
central_atom
- GetCentralNormal()¶
See
central_normal
- SetIsLigand()¶
Set the
IsLigand()
flag explicitly.- Parameters:
ligand (bool) – Whether this residue is a ligand or not
- SetIsProtein()¶
Set the
IsProtein()
flag explicitly.- Parameters:
protein (bool) – Whether this residue is a protein or not
Atom Handle¶
- class AtomHandle¶
Represents an atom in a molecular structure. Atoms are always part of a residue and of an entity (see note on memory management)
- properties¶
All the
generic properties
are available.
- name¶
The atom name, usually a capitalized string consisting of the element and an position indicator, e.g. G1. In general, the atom name uniquely identifies an atom within a residue. However, this is not enforced and there may be more than one atoms with the same name.
- Type:
str
- qualified_name¶
The qualified name consists of the atom name as well as a residue identifier and chain name. For CA of a glycine with residue number 2 of chain A, the qualified name is “A.GLY2.CA”.
- Type:
str
- element¶
The atom’s element. By convention in Openstructure, this is the chemical symbol in uppercase, but this is not strictly enforced and may be a non- existing element or an empty string. Also available as
GetElement()
. Read-only.- Type:
str
- pos¶
The atom’s position in global coordinates, transformed by the entity transformation. Also available as
GetPos()
.Read-only.
- Type:
- original_pos¶
The atom’s untransformed position in global coordinates. Also available as
GetOriginalPos()
.Read-only.
- Type:
- radius¶
The van-der-Waals radius of the atom. Also available as
GetRadius()
. Read/write.- Type:
float
- occupancy¶
The atom’s occupancy in the range 0 to 1. Read/write. Also available as
GetOccupancy()
,SetOccupancy()
.- Type:
float
- b_factor¶
The atom’s temperature factor. Read/write. Also available as
GetBFactor()
,SetBFactor()
.- Type:
float
- charge¶
The atom’s charge. Also available as
GetCharge()
,SetCharge()
.- Type:
float
- residue¶
The residue this atom belongs to. Read-only.
- Type:
- is_hetatom¶
Indicates whether this atom is a HETATM. When loading a structure from a PDB file, all non-standard aminoacid and nucleotide atoms as well as ligands are marked as HETATM.
- bonds¶
List of bonds this atom is involved in. The bonds are in no particular order. Also available as
GetBondList()
.- Type:
list of
bond handles
- index¶
Atom index (starting at 0) within entity.
- Type:
int
- hash_code¶
A unique identifier for this atom. Note that a deep copy of an entity (see
EntityHandle.Copy()
) will have atoms with differing identifiers. Shallow copies of the entity preserve the identifier. Atom views on a handle have different identifiers, but the atom view handles (seeAtomView.handle
) have the same identifier. Also available asGetHashCode()
.- Type:
int
- valid¶
Validity of handle.
- Type:
bool
- FindBondToAtom(other_atom)¶
Finds and returns the bond formed between this atom and other_atom. If no bond exists between the two atoms, an invalid
BondHandle
is returned.- Parameters:
other_atom (
AtomHandle
) – The other atom- Return type:
- GetBondCount()¶
- Return type:
int
- GetBondList()¶
See
bonds
- Return type:
BondHandleList
(list ofBondHandle
)
- GetBondPartners()¶
Get list of atoms this atom is bonded with. This method exists as a shortcut to determine all the bonding partners of an atom and avoids code like this:
bond_parters=[] for bond in atom.bonds: if bond.first==atom: bond_partners.append(bond.second) else: bond_partners.append(bond.first)
- Return type:
AtomHandleList
- GetEntity()¶
The entity this atom belongs to
- Return type:
- GetQualifiedName()¶
See
qualified_name
- Return type:
str
- IsHetAtom()¶
See
is_hetatom
- Return type:
bool
Bond Handle¶
- class BondHandle¶
Represents a chemical bond between two atoms (first and second).
- properties¶
All the
generic properties
are available.
- first¶
- second¶
Atoms involved in the bond. No assumptions about the order should be made. With the internal coordinate system enabled, first and second may even be swapped when rebuilding the internal connectivity tree. Also available as
GetFirst()
andGetSecond()
.- Type:
- length¶
Length of the bond. Also available as
GetLength()
.- Type:
float
- bond_order¶
The bond order. Possible values:
1
- single bond2
- double bond3
- triple bond4
- aromatic bond
Also available as
GetBondOrder()
.- Type:
int
- hash_code¶
A unique identifier for this bond handle. Also available as
GetHashCode()
.- Type:
int
- GetBondOrder()¶
See
bond_order
- GetOther(other_atom)¶
Get the other atom. Returns the one of the two atoms that does not match the given one.
- Parameters:
other_atom (
AtomHandle
) – The other atom- Return type:
- SetBondOrder(order)¶
Set the bond order. See
GetBondOrder()
.- Parameters:
order (
int
) – The bond order
See
bond_order
The View Classes¶
Entity
, residue
and
atom
views can store arbitrary
generic properties.
Entity View¶
- class EntityView¶
An entity view represents a structural subset of an
EntityHandle
and containsChainView
s,ResidueView
s,AtomView
s andBondHandle
s. For an introduction, see Introduction to the mol Module.- properties¶
All the
generic properties
are available.
- handle¶
The underlying
handle
of the entity view. Also available asGetHandle()
.- Type:
- chains¶
List of all chains of this entity. The chains are in the same order they have been added.
Also available as
GetChainList()
. To access a single chain, useFindChain()
.This property is read-only.
- Type:
ChainViewList
(list ofChainView
)
- chain_count¶
Number of chains. Read-only. See
GetChainCount()
.- Type:
int
- residues¶
List of all residues of this entity. The returned residues are first sorted by chain and then from N- to C-terminus.
This property is read-only.
Example usage:
for res in ent.residues: print(res.name, res.atom_count)
Also available as
GetResidueList()
. To access a single residue, useFindResidue()
.- Type:
ResidueViewList
(list ofResidueView
)
- residue_count¶
Number of residues. Read-only. See
GetResidueCount()
.- Type:
int
- atoms¶
Get list of all atoms of this entity. To access a single atom, use
FindAtom()
.This property is read-only. Also available as
GetAtomList()
- Type:
AtomViewList
(list ofAtomView
)
- atom_count¶
Number of atoms. Read-only. See
GetAtomCount()
.- Type:
int
- bounds¶
- mass¶
- center_of_atoms¶
- center_of_mass¶
- geometric_center¶
- valid¶
Validity of view.
- Type:
bool
- CreateEmptyView()¶
See
EntityHandle.CreateEmptyView()
- Return type:
- CreateFullView()¶
Returns a copy of this entity. Provided for duck-typing <http://en.wikipedia.org/wiki/Duck_typing> purposes.
- Return type:
- AddChain(chain_handle[, view_add_flags])¶
Add chain to view. By default, only the chain is added to the view, but not its residues and atoms. This behaviour can be changed by passing in an appropriate set of view_add_flags.
The following example just adds a chain without residues and atoms:
pdb = ost.io.LoadPDB(<PDB file name>) v = pdb.CreateEmptyView() v.AddChain(pdb.chains[0])
To get the chain with residues and atoms, go like:
pdb = ost.io.LoadPDB(<PDB file name>) v = pdb.CreateEmptyView() v.AddChain(pdb.chains[0], ost.mol.INCLUDE_ALL)
Note that the view above still lacks bonds which can be added with the
AddAllInclusiveBonds()
method.- Parameters:
chain_handle (
ChainHandle
) – The chain handle to be added.view_add_flags (
int
/ViewAddFlag
) – An ORed together combination ofViewAddFlag
- Return type:
- AddResidue(residue_handle[, view_add_flags])¶
Add residue to view. If the residue’s chain is not already part of the view, it will be added. By default, only the residue is added, but not its atoms. This behaviour can be modified by passing in an appropriate combination of
ViewAddFlag
.- Parameters:
residue_handle (
ResidueHandle
) – The residue handle to be addedview_add_flags (
int
/ViewAddFlag
) – An ORed together combination ofViewAddFlag
- Return type:
- AddAtom(atom_handle[, view_add_flags])¶
Add the atom to view. The chain and residue of the atom are added to the view if they haven’t already been added.
- Parameters:
atom_handle (
AtomHandle
) – The atom handleview_add_flags (
int
/ViewAddFlag
) – An ORed together combination ofViewAddFlag
- Return type:
- AddBond(bond_handle)¶
Add bond to the view.
- Parameters:
bond_handle (BondHandle) – Bond to be added
- Returns:
True if the bond has been added, false if the bond was already present in the view.
- Return type:
bool
- AddAllInclusiveBonds()¶
Convenience method to add all bonds to the view for which both of the bond partners are in the view. This method is useful when manually assembling views.
- RemoveChain(chain)¶
Remove chain and all its residues from the entity view.
- Parameters:
chain (ChainView) – The chain to be removed. May be invalid
- RemoveResidue(residue)¶
Remove residue from view.
- Parameters:
residue (ResidueView) – The residue to be removed. May be invalid
- RemoveAtom(atom)¶
Remove atom from view
- Parameters:
atom (AtomView) – The atom to be removed. May be invalid
- GetAngle(atom1, atom2, atom3)¶
- Parameters:
atom1 (AtomHandle) –
atom2 (AtomHandle) –
atom3 (AtomHandle) –
- Return type:
float
- FindWithin(pos, radius)¶
Find all atoms that are within radius of the given position. See
EntityHandle.FindWithin()
.- Parameters:
pos (
Vec3
) – Center of sphereradius (float) – The radius of the sphere
- Returns:
AtomHandleList
(list ofAtomHandle
)
- FindChain(chain_name)¶
Find chain by name.
- FindResidue(residue)¶
Deprecated. Use
ViewForHandle()
instead.- Parameters:
residue (ResidueHandle) – Residue handle
- Returns:
The residue view pointing the handle, or an invalid handle if the residue is not part of the view
- Return type:
- FindResidue(chain_name, res_num)
Find residue by chain name and residue number.
- Parameters:
chain_name (str) – Chain identifier, e.g. “A”
res_num (
ResNum
) – residue number
- Returns:
The residue if present in the view, an invalid
ResidueView
otherwise- Return type:
- FindAtom(chain_name, res_num, atom_name)¶
- ViewForHandle(handle)¶
Find chain view, residue view or atom view of pointing to the given handle.
- Parameters:
handle – handle to search for
- Returns:
The view pointing the handle, or an invalid handle if the handle is not part of the view
- Return type:
(Chain|Residue|Atom)View
- ExtendViewToResidues()¶
Extend current view to include all atoms of each residue where at least one atom is selected currently.
- Returns:
The extended view
- Return type:
- ExtendViewToSurrounding(handle)¶
Extend current view to include all atoms that are within the sum of their van-der-Waals radius + gap.
This includes all atoms within: at1.GetRadius() + at2.GetRadius() + gap.
- Parameters:
gap (float) – the gap between atoms
- Returns:
The extended view
- Return type:
- Select(query, flags=0)¶
Perform selection on entity view. See
EntityHandle.Select()
.- Parameters:
- Return type:
- Copy()¶
Returns a deep copy of the entity view. The effect is identical to calling
the_copy = view.Select('')
- Return type:
- Dump()¶
Returns a string containing a human-readable summary of the entity view.
- Return type:
str
- GetMass()¶
Get total mass of view.
- Return type:
float
- GetCenterOfMass()¶
Get the center of mass. For a non-mass-weighted center, see
GetCenterOfAtoms()
.- Return type:
- GetGeometricCenter()¶
Get the geometric center (the center of the axis-aligned bounding box).
- Return type:
- GetCenterOfAtoms()¶
Calculates the center of atoms. For a mass-weighted center, use
GetCenterOfMass()
.- Return type:
- GetBondCount()¶
Get number of bonds
- Return type:
int
- GetBondList()¶
See
bonds
- Return type:
BondHandleList
(list ofBondHandle
)
- GetChainCount()¶
See
chain_count
- GetResidueCount()¶
See
residue_count
- GetAtomCount()¶
See
atom_count
Chain View¶
- class ChainView¶
A view representation of a
ChainHandle
. Mostly, the same functionality is provided as for the handle.- properties¶
All the
generic properties
are available.
- handle¶
The chain handle this view points to. Also available as
GetHandle()
.- Type:
- hash_code¶
A unique identifier for this chain view. Note, that this is not the same as for the chain handle (see
ChainHandle.hash_code
). Also available asGetHashCode()
.- Type:
int
- name¶
The chain name. The name uniquely identifies the chain in the entity. In most cases, the chain name is one character. This is a restriction of the PDB file format. However, you are free to use longer names as long as you don’t want to save them as PDB files.
This property is read-only. To change the name, use an
XCSEditor
.Also available as
GetName()
- Type:
str
- residues¶
List of all residues of this chain. The residues are sorted from N- to C-terminus. Usually the residue numbers are in ascending order (see
in_sequence
).This property is read-only.
Example usage:
chain=ent.FindChain("A") for res in chain.residues: print(res.name, res.atom_count)
Also available as
GetResidueList()
. To access a single residue, useFindResidue()
.- Type:
ResidueViewList
(list ofResidueView
)
- in_sequence¶
Whether the residue numbers are in ascending order. For example:
chain=ent.FindChain("A") print(chain.residues) # [A.GLY1, A.GLY2, A.GLY4A, A.GLY4B] print(chain.in_sequence) # prints true chain=ent.FindChain("B") print(chain.residues) # [B.GLY1, B.GLY4, B.GLY3] print(chain.in_sequence) # prints false
Note that the value of in_sequence is independent from the value of
ChainHandle.in_sequence
.- Type:
bool
- atoms¶
Get list of all atoms of this chain. To access a single atom, use
FindAtom()
.This property is read-only. Also available as
GetAtomList()
- Type:
AtomViewList
(list ofAtomView
)
- bounds¶
Axis-aligned bounding box of the chain. Read-only
- Type:
- center_of_atoms¶
- center_of_mass¶
- geometric_center¶
- valid¶
Validity of view.
- Type:
bool
- AddAtom(atom_handle[, view_add_flags])¶
Add atom to the view. If the residue of the atom is not already part of the view, it will be added. If the atom does not belong to chain, the result is undefined.
- Parameters:
atom_handle (
AtomHandle
) – The atom to be addedview_add_flags (
int
/ViewAddFlag
) – An ORed together combination ofViewAddFlag
- Return type:
- AddResidue(residue_handle[, view_add_flags])¶
Add residue to the view. If the atom does not belong to chain, the result is undefined. By default, only the residue, but no atoms are added to the view. To change the behavior, pass in a suitable combination of
ViewAddFlag
.- Parameters:
residue_handle (
ResidueHandle
) – The residue handle to be added.view_add_flags (
int
/ViewAddFlag
) – An ORed together combination ofViewAddFlag
- Return type:
- FindAtom(res_num, atom_name)¶
Find atom with the given residue number and atom name.
- FindResidue(res_num)¶
Find residue by residue number.
- Parameters:
res_num (
ResNum
) –- Return type:
- Returns:
The residue view, or an invalid residue view if no residue with the given residue number is in the view.
- ViewForHandle(handle)¶
Find residue view or atom view of pointing to the given handle.
- Parameters:
handle – handle to search for
- Returns:
The view pointing the handle, or an invalid handle if the handle is not part of the view
- Return type:
(Residue|Atom)View
- GetCenterOfAtoms()¶
See
center_of_atoms
- GetCenterOfMass()¶
See
center_of_mass
- GetEntity()¶
The parent entity.
- GetGeometricCenter()¶
See
geometric_center
- GetResidueByIndex(index)¶
Returns the residue view at index. If index is negative or bigger than the number of residues minus one, an invalid ResidueView is returned.
- Parameters:
index (int) – The index
- Return type:
- InSequence()¶
See
in_sequence
- RemoveResidue(residue)¶
Remove residue from chain. If the residue is not part of the view, the chain is left unchanged.
- Parameters:
residue (
ResidueView
) – The residue view. May be invalid- Return type:
None
- RemoveResidues()¶
Remove all residues from this chain view
- Select(query, flags=0)¶
Perform query on chain view. See
EntityHandle.Select()
.- Parameters:
- Return type:
Residue View¶
- class ResidueView¶
A view representation of a
ResidueHandle
. Mostly, the same functionality is provided as for the handle.- properties¶
All the
generic properties
are available.
- handle¶
The residue handle this view points to. Also available as
GetHandle()
.- Type:
- hash_code¶
A unique identifier for this residue view. Note, that this is not the same as for the residue handle (see
ResidueHandle.hash_code
). Also available asGetHashCode()
.- Type:
int
- name¶
- qualified_name¶
- number¶
- one_letter_code¶
- bounds¶
- mass¶
- center_of_atoms¶
- center_of_mass¶
- geometric_center¶
- phi_torsion¶
- psi_torsion¶
- omega_torsion¶
- chem_class¶
- chem_type¶
- sec_structure¶
- is_ligand¶
- is_protein¶
- peptide_linking¶
- nucleotide_linking¶
- central_atom¶
- central_normal¶
- valid¶
- next¶
- prev¶
See the respective attributes in
ResidueHandle
.
- atoms¶
Get list of all atoms of this residue included in the view. To access a single atom, use
FindAtom()
.This property is read-only. Also available as
GetAtomList()
- Type:
AtomViewList
(list ofAtomView
)
- atom_count¶
Number of atoms included in the view. Read-only. See
GetAtomCount()
.- Type:
int
- chain¶
The chain this residue belongs to. Read-only. Also available as
GetChain()
- Type:
- index¶
Residue index (starting at 0) within chain view.
- RemoveAtom(atom_view)¶
Remove atom from residue and all associated bonds. If the atom is not part of the view, the residue view is left unchanged.
- Parameters:
atom_view (
AtomView
) – The atom to be removed. May be an invalid handle- Return type:
None
- ViewForHandle(handle)¶
Find atom view of pointing to the given handle.
- Parameters:
handle – handle to search for
- Returns:
The view pointing the handle, or an invalid handle if the handle is not part of the view
- Return type:
- GetName()¶
- GetQualifiedName()¶
- GetNumber()¶
- GetOneLetterCode()¶
- GetBounds()¶
- GetMass()¶
- GetCenterOfAtoms()¶
- GetCenterOfMass()¶
- GetGeometricCenter()¶
- GetPhiTorsion()¶
- GetPsiTorsion()¶
- GetOmegaTorsion()¶
- GetChemClass()¶
- SetChemClass()¶
- GetChemType()¶
- SetChemType()¶
- GetSecStructure()¶
- IsLigand()¶
- SetIsLigand()¶
- IsProtein()¶
- IsPeptideLinking()¶
- IsNucleotideLinking()¶
- GetCentralAtom()¶
- SetCentralAtom()¶
- GetCentralNormal()¶
- IsValid()¶
- GetNext()¶
- GetPrev()¶
See the respective methods in
ResidueHandle
.
- FindAtom(atom_name)¶
Find atom by name. Returns the atom, or an invalid atom if no atoms with the given name is part of the view.
- Parameters:
atom_name (str) –
- Return type:
- IsAtomIncluded(atom_handle)¶
Returns true if the given atom is part of the view, false if not.
- Parameters:
atom_handle (
AtomHandle
) –- Return type:
bool
- AddAtom(atom_handle[, flags])¶
Add atom to the view.
- Parameters:
atom_handle (
AtomHandle
) – Atom handle to be addedflags (
int
/ViewAddFlag
) – An ORed together combination ofViewAddFlag
- Return type:
- GetAtomCount()¶
See
atom_count
- Select(query, flags=0)¶
Perform selection on residue view. See
EntityHandle.Select()
.- Parameters:
- Return type:
Atom View¶
- class AtomView¶
A view representation of an
AtomHandle
. Mostly, the same functionality is provided as for the handle.- properties¶
All the
generic properties
are available.
- handle¶
The underlying
AtomHandle
of the atom view. Also available asGetHandle()
.- Type:
- hash_code¶
A unique identifier for this atom view. Note, that this is not the same as for the atom handle (see
AtomHandle.hash_code
). Also available asGetHashCode()
.- Type:
int
Functions¶
Boolean Operators¶
- Intersection(view_a, view_b)¶
Calculates and returns the intersection of view_a and view_b. view_a and view_b must be views of the same entity.
- Parameters:
view_a (EntityView) – first view
view_b (EntityView) – second view
- Difference(view_a, view_b)¶
Calculates and returns the difference between view_a and view_b. view_a and view_b must be views of the same entity.The returned view will contain atoms, residues, chains and bonds that are in view_a and not in view_b.
- Parameters:
view_a (EntityView) – first view
view_b (EntityView) – second view
- Union(view_a, view_b)¶
Calculates and returns the union of view_a and view_b. view_a and view_b must be views of the same entity.The returned view will contain all atoms, residues, chains and bonds that are either in view_a, view_b or part of both views.
- Parameters:
view_a (EntityView) – first view
view_b (EntityView) – second view
Residue Numbering¶
- class ResNum(num, ins_code='\\0')¶
Number for a residue. The residue number has a numeric part and an (optional) insertion-code. You can work with this object as if it was an integer and comparison will look first at the numeric part and then the insertion-code. All access to existing objects is read-only. Openstructure supports a range of (-8388608 to 8388607) for the numeric part. However, the PDB format only supports a range of (-999, 9999). This becomes relevant when a structure is saved in PDB format where an IOException is raised if the PDB range is not respected.
- Parameters:
num (
int
) – Numeric part of residue number.ins_code (
str
) – Alpha-numeric part of residue number (optional insertion code). Only first character kept.
- num¶
Numeric part of residue number.
- Type:
int
- ins_code¶
Alpha-numeric part of residue number (insertion code). Single character.
- Type:
str
ChainType¶
A ChainType fills the ChainHandle.type
attribute. Different types are
described in the ChainType
enum. The type is set with an editor in
EditorBase.SetChainType()
. Further convenience functions are described
here.
- class ChainType¶
The ChainType enum enumerates all types defined by the PDB for the MMCif file format. Following values are supported:
CHAINTYPE_POLY
,CHAINTYPE_NON_POLY
,CHAINTYPE_WATER
,CHAINTYPE_POLY_PEPTIDE_D
,CHAINTYPE_POLY_PEPTIDE_L
,CHAINTYPE_POLY_DN
,CHAINTYPE_POLY_RN
,CHAINTYPE_POLY_SAC_D
,CHAINTYPE_POLY_SAC_L
,CHAINTYPE_POLY_DN_RN
,CHAINTYPE_UNKNOWN
,CHAINTYPE_MACROLIDE
,CHAINTYPE_CYCLIC_PSEUDO_PEPTIDE
,CHAINTYPE_POLY_PEPTIDE_DN_RN
,CHAINTYPE_BRANCHED
,CHAINTYPE_OLIGOSACCHARIDE
,CHAINTYPE_N_CHAINTYPES
Where
CHAINTYPE_N_CHAINTYPES
holds the number of different types available.
- ChainTypeFromString(identifier)¶
Create a ChainType item for a given string.
- Parameters:
identifier (
str
) – String to request a type for.- Raises:
runtime_error
if identifier is unrecognised.- Returns:
ViewAddFlag¶
- class ViewAddFlag¶
Defines flags controlling behaviour of routines adding handles to views:
INCLUDE_ATOMS
- Include all atoms when adding a residue handle to a viewINCLUDE_RESIDUES
- Include all residues when adding a chain to a viewINCLUDE_CHAINS
- Include all chains when creating a new entity viewINCLUDE_ALL
=INCLUDE_ATOMS
|INCLUDE_RESIDUES
|INCLUDE_CHAINS
- Convenience flags to include all substructuresCHECK_DUPLICATES
- If set, it will be checked that no duplicates are created when adding a new handle
Flags can be ORed to combine them.
SecStructure¶
- class SecStructure(type)¶
Defines a secondary structure type following the types defined by DSSP.
- Parameters:
type (
SecStructure.Type
/str
) – Type to be set for this object.
- IsHelical()¶
- Returns:
True, if the set type is any type of helix (i.e. ALPHA_HELIX, PI_HELIX or THREE_TEN_HELIX)
- IsExtended()¶
- Returns:
True, if the set type is any type of beta sheet (i.e. EXTENDED, BETA_BRIDGE)
- IsCoil()¶
- Returns:
True, if the set type is any type of coil (i.e. TURN, BEND or COIL)
- __str__()¶
- Returns:
The character corresponding to the set type (see
SecStructure.Type
)
- class SecStructure.Type¶
Enumerates all popssible secondary structure types distinguished by DSSP. Their values with the corresponding character code are listed here:
ALPHA_HELIX
= ‘H’PI_HELIX
= ‘I’THREE_TEN_HELIX
= ‘G’EXTENDED
= ‘E’
BETA_BRIDGE
= ‘B’TURN
= ‘T’BEND
= ‘S’COIL
= ‘C’
ChemClass¶
- class ChemClass(chem_class)¶
The chemical class is used to broadly categorize residues based on their chemical properties. For example, peptides belong to some PEPTIDE_LINKING class. Possible values as constant variable names and as characters:
PEPTIDE_LINKING
= ‘P’D_PEPTIDE_LINKING
= ‘D’L_PEPTIDE_LINKING
= ‘L’RNA_LINKING
= ‘R’DNA_LINKING
= ‘S’NON_POLYMER
= ‘N’
L_SACCHARIDE
= ‘X’D_SACCHARIDE
= ‘Y’SACCHARIDE
= ‘Z’WATER
= ‘W’UNKNOWN
= ‘U’
The constants are defined directly within the
mol
module. Python can implicitly convert characters to objects of this type. Note however that only the first character of astr
object is considered!- Parameters:
chem_class (
str
) – Chemical class to set.
- IsPeptideLinking()¶
- Returns:
True, if set class is PEPTIDE_LINKING, D_PEPTIDE_LINKING or L_PEPTIDE_LINKING
- IsNucleotideLinking()¶
- Returns:
True, if set class is RNA_LINKING or DNA_LINKING
ChemType¶
- class ChemType¶
The chemical type of a residue is a classification of all compounds obtained from the PDB component dictionary. For example, ions belong to the class IONS, amino acids to AMINOACIDS. Possible values as constant variable names and as characters:
IONS
= ‘I’NONCANONICALMOLS
= ‘M’SACCHARIDES
= ‘S’NUCLEOTIDES
= ‘N’AMINOACIDS
= ‘A’
COENZYMES
= ‘E’WATERCOORDIONS
= ‘C’DRUGS
= ‘D’WATERS
= ‘W’UNKNOWN
= ‘U’
Python can implicitly convert characters to objects of this type. Note however that only the first character of a
str
object is considered!- Parameters:
chem_type (
str
) – Chemical type to set.
- IsIon()¶
- Returns:
True, if set type is IONS or WATERCOORDIONS
- IsNucleotide()¶
- Returns:
True, if set type is NUCLEOTIDES
- IsSaccharide()¶
- Returns:
True, if set type is SACCHARIDES
- IsAminoAcid()¶
- Returns:
True, if set type is AMINOACIDS
- IsCoenzyme()¶
- Returns:
True, if set type is COENZYMES
- IsDrug()¶
- Returns:
True, if set type is DRUGS
- IsNonCanonical()¶
- Returns:
True, if set type is NONCANONICALMOLS
- IsWater()¶
- Returns:
True, if set type is WATERS
- IsKnown()¶
- Returns:
True, if set type is not UNKNOWN