Package pipe_control :: Package structure :: Module main
[hide private]
[frames] | no frames]

Module main

source code

Functions [hide private]
 
add_atom(mol_name=None, atom_name=None, res_name=None, res_num=None, pos=[None, None, None], element=None, atom_num=None, chain_id=None, segment_id=None, pdb_record=None)
Add a new atom to the structural data object.
source code
 
add_helix(start=None, end=None, mol_name=None)
Define alpha helical secondary structure for the structural data object.
source code
 
add_model(model_num=None)
Add a new model to the empty structural data object.
source code
 
add_sheet(strand=1, sheet_id='A', strand_count=2, strand_sense=0, start=None, end=None, mol_name=None, current_atom=None, prev_atom=None)
Define beta sheet secondary structure for the structural data object.
source code
numpy rank-3 float64 array, list of str, list of str, list of str, list of int, list of str, list of str
assemble_structural_coordinates(pipes=None, models=None, molecules=None, atom_id=None, lists=False)
Assemble the common atomic coordinates taking sequence alignments into account.
source code
list of lib.structure.internal.object.Internal instances, list of str, list of str
assemble_structural_objects(pipes=None, models=None, molecules=None)
Assemble the atomic coordinates.
source code
 
atomic_fluctuations(pipes=None, models=None, molecules=None, atom_id=None, measure='distance', file=None, format='text', dir=None, force=False)
Write out a correlation matrix of pairwise interatomic distance fluctuations between different structures.
source code
 
average_structure(pipes=None, models=None, molecules=None, atom_id=None, set_mol_name=None, set_model_num=None)
Calculate a mean structure.
source code
 
connect_atom(index1=None, index2=None)
Connect two atoms.
source code
 
com(model=None, atom_id=None)
Calculate the centre of mass (CoM) of all structures.
source code
 
create_diff_tensor_pdb(scale=1.8e-06, file=None, dir=None, force=False)
Create the PDB representation of the diffusion tensor.
source code
 
delete(atom_id=None, model=None, verbosity=1, spin_info=True)
Delete structural data.
source code
 
delete_ss(verbosity=1)
Delete all secondary structure information.
source code
 
displacement(pipes=None, models=None, molecules=None, atom_id=None, centroid=None)
Calculate the rotational and translational displacement between structures or models.
source code
 
find_pivot(pipes=None, models=None, molecules=None, atom_id=None, init_pos=None, func_tol=1e-05, box_limit=200)
Find the pivoted motion of a set of structural models or structures.
source code
 
get_pos(spin_id=None, str_id=None, ave_pos=False)
Load the spins from the structural object into the relax data store.
source code
 
load_spins(spin_id=None, str_id=None, from_mols=None, mol_name_target=None, ave_pos=False, spin_num=True)
Load the spins from the structural object into the relax data store.
source code
 
load_spins_multi_mol(spin_id=None, str_id=None, from_mols=None, mol_name_target=None, ave_pos=False, spin_num=True)
Load the spins from the structural object into the relax data store.
source code
 
pca(pipes=None, models=None, molecules=None, obs_pipes=None, obs_models=None, obs_molecules=None, atom_id=None, algorithm=None, num_modes=4, format='grace', dir=None)
PC analysis of the motions between all the loaded models.
source code
 
read_gaussian(file=None, dir=None, set_mol_name=None, set_model_num=None, verbosity=1, fail=True)
Read structures from a Gaussian log file.
source code
 
read_pdb(file=None, dir=None, read_mol=None, set_mol_name=None, read_model=None, set_model_num=None, alt_loc=None, verbosity=1, merge=False, fail=True)
The PDB loading function.
source code
 
read_xyz(file=None, dir=None, read_mol=None, set_mol_name=None, read_model=None, set_model_num=None, verbosity=1, fail=True)
The XYZ loading function.
source code
float
rmsd(pipes=None, models=None, molecules=None, atom_id=None, atomic=False)
Calculate the RMSD between the loaded models.
source code
 
rotate(R=None, origin=None, model=None, atom_id=None, pipe_name=None)
Rotate the structural data about the origin by the specified forwards rotation.
source code
 
sequence_alignment(pipes=None, models=None, molecules=None, msa_algorithm='Central Star', pairwise_algorithm='NW70', matrix='BLOSUM62', gap_open_penalty=1.0, gap_extend_penalty=1.0, end_gap_open_penalty=0.0, end_gap_extend_penalty=0.0)
Superimpose a set of related, but not identical structures.
source code
 
set_vector(spin=None, xh_vect=None)
Place the XH unit vector into the spin container object.
source code
int, int or None, str
structure_loop(pipes=None, molecules=None, models=None, atom_id=None)
Generator function for looping over all internal structural objects, models and molecules.
source code
 
superimpose(pipes=None, models=None, molecules=None, atom_id=None, displace_id=None, method='fit to mean', centre_type='centroid', centroid=None)
Superimpose a set of structures.
source code
 
translate(T=None, model=None, atom_id=None, pipe_name=None)
Shift the structural data by the specified translation vector.
source code
 
vectors(spin_id1=None, spin_id2=None, model=None, verbosity=1, ave=True, unit=True)
Extract the bond vectors from the loaded structures and store them in the spin container.
source code
 
web_of_motion(pipes=None, models=None, molecules=None, atom_id=None, file=None, dir=None, force=False)
Create a PDB representation of the motion between a set of models.
source code
 
write_pdb(file=None, dir=None, model_num=None, compress_type=0, force=False)
The PDB writing function.
source code
Variables [hide private]
  ds = The relax data storage obje...
  status = Status()
  __package__ = 'pipe_control.structure'

Imports: generic_minimise, array, average, concatenate, dot, float64, mean, ones, std, zeros, norm, F_OK, access, getcwd, search, sys, warn, Relax_data_store, Sequence_alignments, is_float, RelaxError, RelaxFileError, vector_angle_atan2, get_file_path, open_write_file, write_data, correlation_matrix, write_xy_data, write_xy_header, tokenise, write_spin_data, msa_general, msa_residue_numbers, msa_residue_skipping, assemble_atomic_coordinates, assemble_coord_array, loop_coord_structures, Displacements, Internal, pca_analysis, diffusion_tensor, atomic_rmsd, per_atom_rmsd, fit_to_first, fit_to_mean, RelaxWarning, RelaxNoPDBFileWarning, RelaxZeroVectorWarning, molmol, pipes, interatomic_loop, check_mol_res_spin_data, create_spin, generate_spin_id_unique, linear_ave, return_spin, spin_loop, cdp_name, check_pipe, get_pipe, check_structure, pipe_centre_of_mass, Status, Pivot_finder


Function Details [hide private]

add_atom(mol_name=None, atom_name=None, res_name=None, res_num=None, pos=[None, None, None], element=None, atom_num=None, chain_id=None, segment_id=None, pdb_record=None)

source code 

Add a new atom to the structural data object.

Parameters:
  • mol_name (str) - The name of the molecule.
  • atom_name (str or None) - The atom name, e.g. 'H1'.
  • res_name (str or None) - The residue name.
  • res_num (int or None) - The residue number.
  • pos (rank-1 or rank-2 array or list of float) - The position vector of coordinates. If a rank-2 array is supplied, the length of the first dimension must match the number of models.
  • element (str or None) - The element symbol.
  • atom_num (int or None) - The atom number.
  • chain_id (str or None) - The chain identifier.
  • segment_id (str or None) - The segment identifier.
  • pdb_record (str or None) - The optional PDB record name, e.g. 'ATOM' or 'HETATM'.

add_helix(start=None, end=None, mol_name=None)

source code 

Define alpha helical secondary structure for the structural data object.

Parameters:
  • start (int) - The residue number for the start of the helix.
  • end (int) - The residue number for the end of the helix.
  • mol_name (str or None) - Define the secondary structure for a specific molecule.

add_sheet(strand=1, sheet_id='A', strand_count=2, strand_sense=0, start=None, end=None, mol_name=None, current_atom=None, prev_atom=None)

source code 

Define beta sheet secondary structure for the structural data object.

Parameters:
  • strand (int) - Strand number which starts at 1 for each strand within a sheet and increases by one.
  • sheet_id (str) - The sheet identifier. To match the PDB standard, sheet IDs should range from 'A' to 'Z'.
  • strand_count (int) - The number of strands in the sheet.
  • strand_sense (int) - Sense of strand with respect to previous strand in the sheet. 0 if first strand, 1 if parallel, -1 if anti-parallel.
  • start (int) - The residue number for the start of the sheet.
  • end (int) - The residue number for the end of the sheet.
  • mol_name (str or None) - Define the secondary structure for a specific molecule.
  • current_atom (str or None) - The name of the first atom in the current strand, to link the current back to the previous strand.
  • prev_atom (str or None) - The name of the last atom in the previous strand, to link the current back to the previous strand.

assemble_structural_coordinates(pipes=None, models=None, molecules=None, atom_id=None, lists=False)

source code 

Assemble the common atomic coordinates taking sequence alignments into account.

Parameters:
  • pipes (None or list of str) - The data pipes to assemble the coordinates from.
  • models (None or list of lists of int) - The list of models for each data pipe. The number of elements must match the pipes argument. If set to None, then all models will be used.
  • molecules (None or list of lists of str) - The list of molecules for each data pipe. The number of elements must match the pipes argument.
  • atom_id (str or None) - The molecule, residue, and atom identifier string. This matches the spin ID string format.
  • lists (bool) - A flag which if true will cause the object ID list per molecule, the model number list per molecule, and the molecule name list per molecule to also be returned.
Returns: numpy rank-3 float64 array, list of str, list of str, list of str, list of int, list of str, list of str
The array of atomic coordinates (first dimension is the model and/or molecule, the second are the atoms, and the third are the coordinates); a list of unique IDs for each structural object, model, and molecule; the common list of molecule names; the common list of residue names; the common list of residue numbers; the common list of atom names; the common list of element names.

assemble_structural_objects(pipes=None, models=None, molecules=None)

source code 

Assemble the atomic coordinates.

Parameters:
  • pipes (None or list of str) - The data pipes to assemble the coordinates from.
  • models (None or list of lists of int) - The list of models for each data pipe. The number of elements must match the pipes argument. If set to None, then all models will be used.
  • molecules (None or list of lists of str) - The list of molecules for each data pipe. The number of elements must match the pipes argument.
Returns: list of lib.structure.internal.object.Internal instances, list of str, list of str
The structural objects, structural object names, and data pipes list.

atomic_fluctuations(pipes=None, models=None, molecules=None, atom_id=None, measure='distance', file=None, format='text', dir=None, force=False)

source code 

Write out a correlation matrix of pairwise interatomic distance fluctuations between different structures.

Parameters:
  • pipes (None or list of str) - The data pipes to generate the interatomic distance fluctuation correlation matrix for.
  • models (None or list of lists of int) - The list of models to generate the interatomic distance fluctuation correlation matrix for. The number of elements must match the pipes argument. If set to None, then all models will be used.
  • molecules (None or list of lists of str) - The list of molecules to generate the interatomic distance fluctuation correlation matrix for. The number of elements must match the pipes argument.
  • atom_id (str or None) - The atom identification string of the coordinates of interest. This matches the spin ID string format.
  • measure (str) - The type of fluctuation to measure. This can be either 'distance' or 'angle'.
  • file (str) - The name of the file to write.
  • format (str) - The output format. This can be set to "text" for text file output, or "gnuplot" for creating a gnuplot script.
  • dir (str or None) - The directory where the file will be placed. If set to None, then the file will be placed in the current directory.
  • force (bool) - The force flag which if True will cause the file to be overwritten.

average_structure(pipes=None, models=None, molecules=None, atom_id=None, set_mol_name=None, set_model_num=None)

source code 

Calculate a mean structure.

Parameters:
  • pipes (None or list of str) - The data pipes containing structures to average.
  • models (None or list of lists of int) - The list of models for each data pipe. The number of elements must match the pipes argument. If set to None, then all models will be used.
  • molecules (None or list of lists of str) - The list of molecules for each data pipe. The number of elements must match the pipes argument.
  • atom_id (str or None) - The molecule, residue, and atom identifier string. This matches the spin ID string format.
  • set_mol_name (None or str) - The molecule name for the averaged molecule.
  • set_model_num (None or int) - The model number for the averaged molecule.

connect_atom(index1=None, index2=None)

source code 

Connect two atoms.

Parameters:
  • index1 (str) - The global index of the first atom.
  • index2 (str) - The global index of the first atom.

com(model=None, atom_id=None)

source code 

Calculate the centre of mass (CoM) of all structures.

The value will be stored in the current data pipe 'com' variable.

Parameters:
  • model (int or None) - Only use a specific model.
  • atom_id (str or None) - The molecule, residue, and atom identifier string. This matches the spin ID string format. If not given, then all structural data will be used for calculating the CoM.

create_diff_tensor_pdb(scale=1.8e-06, file=None, dir=None, force=False)

source code 

Create the PDB representation of the diffusion tensor.

Parameters:
  • scale (float) - The scaling factor for the diffusion tensor.
  • file (str) - The name of the PDB file to create.
  • dir (str) - The name of the directory to place the PDB file into.
  • force (bool) - Flag which if set to True will overwrite any pre-existing file.

delete(atom_id=None, model=None, verbosity=1, spin_info=True)

source code 

Delete structural data.

Parameters:
  • atom_id (str or None) - The molecule, residue, and atom identifier string. This matches the spin ID string format. If not given, then all structural data will be deleted.
  • model (None or int) - Individual structural models from a loaded ensemble can be deleted by specifying the model number.
  • verbosity (int) - The amount of information to print to screen. Zero corresponds to minimal output while higher values increase the amount of output. The default value is 1.
  • spin_info (bool) - A flag which if True will cause all structural information in the spin containers and interatomic data containers to be deleted as well. If False, then only the 3D structural data will be deleted.

delete_ss(verbosity=1)

source code 

Delete all secondary structure information.

Parameters:
  • verbosity (int) - The amount of information to print to screen. Zero corresponds to minimal output while higher values increase the amount of output. The default value is 1.

displacement(pipes=None, models=None, molecules=None, atom_id=None, centroid=None)

source code 

Calculate the rotational and translational displacement between structures or models.

All results will be placed into the current data pipe cdp.structure.displacements data structure.

Parameters:
  • pipes (None or list of str) - The data pipes to determine the displacements for.
  • models (None or list of lists of int) - The list of models to determine the displacements for. The number of elements must match the pipes argument. If set to None, then all models will be used.
  • molecules (None or list of lists of str) - The list of molecules to determine the displacements for. The number of elements must match the pipes argument.
  • atom_id (str or None) - The atom identification string of the coordinates of interest. This matches the spin ID string format.
  • centroid (list of float or numpy rank-1, 3D array) - An alternative position of the centroid, used for studying pivoted systems.

find_pivot(pipes=None, models=None, molecules=None, atom_id=None, init_pos=None, func_tol=1e-05, box_limit=200)

source code 

Find the pivoted motion of a set of structural models or structures.

The pivot will be placed into the current data pipe cdp.structure.pivot data structure.

Parameters:
  • pipes (None or list of str) - The data pipes to use in the motional pivot algorithm.
  • models (None or list of lists of int) - The list of models to use. The number of elements must match the pipes argument. If set to None, then all models will be used.
  • molecules (None or list of lists of str) - The list of molecules to find the pivoted motion for. The number of elements must match the pipes argument.
  • atom_id (str or None) - The molecule, residue, and atom identifier string. This matches the spin ID string format.
  • init_pos (list of float or numpy rank-1, 3D array) - The starting pivot position for the pivot point optimisation.
  • func_tol (None or float) - The function tolerance which, when reached, terminates optimisation. Setting this to None turns of the check.
  • box_limit (int) - The simplex optimisation used in this function is constrained withing a box of +/- x Angstrom containing the pivot point using the logarithmic barrier function. This argument is the value of x.

get_pos(spin_id=None, str_id=None, ave_pos=False)

source code 

Load the spins from the structural object into the relax data store.

Parameters:
  • spin_id (str) - The molecule, residue, and spin identifier string.
  • str_id (int or str) - The structure identifier. This can be the file name, model number, or structure number.
  • ave_pos (bool) - A flag specifying if the average atom position or the atom position from all loaded structures is loaded into the SpinContainer.

load_spins(spin_id=None, str_id=None, from_mols=None, mol_name_target=None, ave_pos=False, spin_num=True)

source code 

Load the spins from the structural object into the relax data store.

Parameters:
  • spin_id (str) - The molecule, residue, and spin identifier string.
  • str_id (int or str) - The structure identifier. This can be the file name, model number, or structure number.
  • from_mols (list of str or None) - The list of similar, but not necessarily identical molecules to load spin information from.
  • mol_name_target (str or None) - The name of target molecule container, overriding the name of the loaded structures
  • ave_pos (bool) - A flag specifying if the average atom position or the atom position from all loaded structures is loaded into the SpinContainer.
  • spin_num (bool) - A flag specifying if the spin number should be loaded.

load_spins_multi_mol(spin_id=None, str_id=None, from_mols=None, mol_name_target=None, ave_pos=False, spin_num=True)

source code 

Load the spins from the structural object into the relax data store.

Parameters:
  • spin_id (str) - The molecule, residue, and spin identifier string.
  • str_id (int or str) - The structure identifier. This can be the file name, model number, or structure number.
  • from_mols (list of str or None) - The list of similar, but not necessarily identical molecules to load spin information from.
  • mol_name_target (str or None) - The name of target molecule container, overriding the name of the loaded structures
  • ave_pos (bool) - A flag specifying if the average atom position or the atom position from all loaded structures is loaded into the SpinContainer.
  • spin_num (bool) - A flag specifying if the spin number should be loaded.

pca(pipes=None, models=None, molecules=None, obs_pipes=None, obs_models=None, obs_molecules=None, atom_id=None, algorithm=None, num_modes=4, format='grace', dir=None)

source code 

PC analysis of the motions between all the loaded models.

Parameters:
  • pipes (None or list of str) - The data pipes to perform the PC analysis on.
  • models (None or list of lists of int) - The list of models to perform the PC analysis on. The number of elements must match the pipes argument. If set to None, then all models will be used.
  • molecules (None or list of lists of str) - The list of molecules to perform the PC analysis on. The number of elements must match the pipes argument.
  • obs_pipes (None or list of str) - The data pipes in the PC analysis which will have zero weight. These structures are for comparison.
  • obs_models (None or list of lists of int) - The list of models in the PC analysis which will have zero weight. These structures are for comparison. The number of elements must match the pipes argument. If set to None, then all models will be used.
  • obs_molecules (None or list of lists of str) - The list of molecules in the PC analysis which will have zero weight. These structures are for comparison. The number of elements must match the pipes argument.
  • atom_id (str or None) - The atom identification string of the coordinates of interest. This matches the spin ID string format.
  • algorithm (str) - The PCA algorithm to use (either 'eigen' or 'svd').
  • num_modes (int) - The number of PCA modes to calculate.
  • format (str) - The graph format to use.
  • dir (str) - The optional directory to place the graphs into.

read_gaussian(file=None, dir=None, set_mol_name=None, set_model_num=None, verbosity=1, fail=True)

source code 

Read structures from a Gaussian log file.

Parameters:
  • file (str) - The name of the Gaussian log file to read.
  • dir (str or None) - The directory where the Gaussian log file is located. If set to None, then the file will be searched for in the current directory.
  • set_mol_name (None, str, or list of str) - Set the names of the molecules which are loaded. If set to None, then the molecules will be automatically labelled based on the file name or other information.
  • set_model_num (None, int, or list of int) - Set the model number of the loaded molecule.
  • fail (bool) - A flag which, if True, will cause a RelaxError to be raised if the Gaussian log file does not exist. If False, then a RelaxWarning will be trown instead.
  • verbosity (int) - The amount of information to print to screen. Zero corresponds to minimal output while higher values increase the amount of output. The default value is 1.
Raises:
  • RelaxFileError - If the fail flag is set, then a RelaxError is raised if the Gaussian log file does not exist.

read_pdb(file=None, dir=None, read_mol=None, set_mol_name=None, read_model=None, set_model_num=None, alt_loc=None, verbosity=1, merge=False, fail=True)

source code 

The PDB loading function.

Parameters:
  • file (str) - The name of the PDB file to read.
  • dir (str or None) - The directory where the PDB file is located. If set to None, then the file will be searched for in the current directory.
  • read_mol (None, int, or list of int) - The molecule(s) to read from the file, independent of model. The molecules are numbered consecutively from 1. If set to None, then all molecules will be loaded.
  • set_mol_name (None, str, or list of str) - Set the names of the molecules which are loaded. If set to None, then the molecules will be automatically labelled based on the file name or other information.
  • read_model (None, int, or list of int) - The PDB model to extract from the file. If set to None, then all models will be loaded.
  • set_model_num (None, int, or list of int) - Set the model number of the loaded molecule. If set to None, then the PDB model numbers will be preserved, if they exist.
  • alt_loc (str or None) - The PDB ATOM record 'Alternate location indicator' field value to select which coordinates to use.
  • verbosity (int) - The amount of information to print to screen. Zero corresponds to minimal output while higher values increase the amount of output. The default value is 1.
  • merge (bool) - A flag which if set to True will try to merge the PDB structure into the currently loaded structures.
  • fail (bool) - A flag which, if True, will cause a RelaxError to be raised if the PDB file does not exist. If False, then a RelaxWarning will be trown instead.
Raises:
  • RelaxFileError - If the fail flag is set, then a RelaxError is raised if the PDB file does not exist.

read_xyz(file=None, dir=None, read_mol=None, set_mol_name=None, read_model=None, set_model_num=None, verbosity=1, fail=True)

source code 

The XYZ loading function.

Parameters:
  • file (str) - The name of the XYZ file to read.
  • dir (str or None) - The directory where the XYZ file is located. If set to None, then the file will be searched for in the current directory.
  • read_mol (None, int, or list of int) - The molecule(s) to read from the file, independent of model. If set to None, then all molecules will be loaded.
  • set_mol_name (None, str, or list of str) - Set the names of the molecules which are loaded. If set to None, then the molecules will be automatically labelled based on the file name or other information.
  • read_model (None, int, or list of int) - The XYZ model to extract from the file. If set to None, then all models will be loaded.
  • set_model_num (None, int, or list of int) - Set the model number of the loaded molecule. If set to None, then the XYZ model numbers will be preserved, if they exist.
  • fail (bool) - A flag which, if True, will cause a RelaxError to be raised if the XYZ file does not exist. If False, then a RelaxWarning will be trown instead.
  • verbosity (int) - The amount of information to print to screen. Zero corresponds to minimal output while higher values increase the amount of output. The default value is 1.
Raises:
  • RelaxFileError - If the fail flag is set, then a RelaxError is raised if the XYZ file does not exist.

rmsd(pipes=None, models=None, molecules=None, atom_id=None, atomic=False)

source code 

Calculate the RMSD between the loaded models.

The RMSD value will be placed into the current data pipe cdp.structure.rmsd data structure.

Parameters:
  • pipes (None or list of str) - The data pipes to determine the RMSD for.
  • models (None or list of lists of int) - The list of models to determine the RMSD for. The number of elements must match the pipes argument. If set to None, then all models will be used.
  • molecules (None or list of lists of str) - The list of molecules to determine the RMSD for. The number of elements must match the pipes argument.
  • atom_id (str or None) - The atom identification string of the coordinates of interest. This matches the spin ID string format.
  • atomic (bool) - A flag which if True will allow for per-atom RMSDs to be additionally calculated.
Returns: float
The RMSD value.

rotate(R=None, origin=None, model=None, atom_id=None, pipe_name=None)

source code 

Rotate the structural data about the origin by the specified forwards rotation.

Parameters:
  • R (numpy 3D, rank-2 array or a 3x3 list of floats) - The forwards rotation matrix.
  • origin (numpy 3D, rank-1 array or list of len 3 or None) - The origin of the rotation. If not supplied, the origin will be set to [0, 0, 0].
  • model (int) - The model to rotate. If None, all models will be rotated.
  • atom_id (str or None) - The molecule, residue, and atom identifier string. Only atoms matching this selection will be used.
  • pipe_name (None or str) - The name of the data pipe containing the structures to translate. This defaults to the current data pipe.

sequence_alignment(pipes=None, models=None, molecules=None, msa_algorithm='Central Star', pairwise_algorithm='NW70', matrix='BLOSUM62', gap_open_penalty=1.0, gap_extend_penalty=1.0, end_gap_open_penalty=0.0, end_gap_extend_penalty=0.0)

source code 

Superimpose a set of related, but not identical structures.

Parameters:
  • pipes (None or list of str) - The data pipes to include in the alignment and superimposition.
  • models (list of lists of int or None) - The list of models to for each data pipe superimpose. The number of elements must match the pipes argument. If set to None, then all models will be used.
  • molecules (None or list of str) - The molecule names to include in the alignment and superimposition. The number of elements must match the pipes argument.
  • msa_algorithm (str) - The multiple sequence alignment (MSA) algorithm to use.
  • pairwise_algorithm (str) - The pairwise sequence alignment algorithm to use.
  • matrix (str) - The substitution matrix to use.
  • gap_open_penalty (float) - The penalty for introducing gaps, as a positive number.
  • gap_extend_penalty (float) - The penalty for extending a gap, as a positive number.
  • end_gap_open_penalty (float) - The optional penalty for opening a gap at the end of a sequence.
  • end_gap_extend_penalty (float) - The optional penalty for extending a gap at the end of a sequence.

set_vector(spin=None, xh_vect=None)

source code 

Place the XH unit vector into the spin container object.

Parameters:
  • spin (SpinContainer instance) - The spin container object.
  • xh_vect (array of len 3) - The unit vector parallel to the XH bond.

structure_loop(pipes=None, molecules=None, models=None, atom_id=None)

source code 

Generator function for looping over all internal structural objects, models and molecules.

Parameters:
  • pipes (None or list of str) - The data pipes to loop over.
  • models (None or list of lists of int) - The list of models for each data pipe. The number of elements must match the pipes argument. If set to None, then all models will be used.
  • molecules (None or list of lists of str) - The list of molecules for each data pipe. The number of elements must match the pipes argument.
  • atom_id (None or str) - The molecule, residue, and atom identifier string of the coordinates of interest. This matches the spin ID string format.
Returns: int, int or None, str
The data pipe index, model number, and molecule name.

superimpose(pipes=None, models=None, molecules=None, atom_id=None, displace_id=None, method='fit to mean', centre_type='centroid', centroid=None)

source code 

Superimpose a set of structures.

Parameters:
  • pipes (None or list of str) - The data pipes to include in the alignment and superimposition.
  • models (list of lists of int or None) - The list of models to for each data pipe superimpose. The number of elements must match the pipes argument. If set to None, then all models will be used.
  • molecules (None or list of str) - The molecule names to include in the alignment and superimposition. The number of elements must match the pipes argument.
  • atom_id (str or None) - The molecule, residue, and atom identifier string. This matches the spin ID string format.
  • displace_id (None, str, or list of str) - The atom ID string for restricting the displacement to a subset of all atoms. If not set, then all atoms will be translated and rotated. This can be a list of atom IDs with each element corresponding to one of the structures.
  • method (str) - The superimposition method. It must be one of 'fit to mean' or 'fit to first'.
  • centre_type (str) - The type of centre to superimpose over. This can either be the standard centroid superimposition or the CoM could be used instead.
  • centroid (list of float or numpy rank-1, 3D array) - An alternative position of the centroid to allow for different superpositions, for example of pivot point motions.

translate(T=None, model=None, atom_id=None, pipe_name=None)

source code 

Shift the structural data by the specified translation vector.

Parameters:
  • T (numpy rank-1, 3D array or list of float) - The translation vector.
  • model (int or None) - The model to translate. If None, all models will be rotated.
  • atom_id (str or None) - The molecule, residue, and atom identifier string. Only atoms matching this selection will be used.
  • pipe_name (None or str) - The name of the data pipe containing the structures to translate. This defaults to the current data pipe.

vectors(spin_id1=None, spin_id2=None, model=None, verbosity=1, ave=True, unit=True)

source code 

Extract the bond vectors from the loaded structures and store them in the spin container.

Parameters:
  • spin_id1 (str) - The spin identifier string of the first spin of the pair.
  • spin_id2 (str) - The spin identifier string of the second spin of the pair.
  • model (str) - The model to extract the vector from. If None, all vectors will be extracted.
  • verbosity (int) - The higher the value, the more information is printed to screen.
  • ave (bool) - A flag which if True will cause the average of all vectors to be extracted.
  • unit (bool) - A flag which if True will cause the function to calculate the unit vectors.

web_of_motion(pipes=None, models=None, molecules=None, atom_id=None, file=None, dir=None, force=False)

source code 

Create a PDB representation of the motion between a set of models.

This will create a PDB file containing the atoms of all models, with identical atoms links using CONECT records. This function only supports the internal structural object.

Parameters:
  • pipes (None or list of str) - The data pipes to generate the web between.
  • models (None or list of lists of int) - The list of models to generate the web between. The number of elements must match the pipes argument. If set to None, then all models will be used.
  • molecules (None or list of lists of str) - The list of molecules to generate the web between. The number of elements must match the pipes argument.
  • atom_id (str or None) - The atom identification string of the coordinates of interest. This matches the spin ID string format.
  • file (str) - The name of the PDB file to write.
  • dir (str or None) - The directory where the PDB file will be placed. If set to None, then the file will be placed in the current directory.
  • force (bool) - The force flag which if True will cause the file to be overwritten.

write_pdb(file=None, dir=None, model_num=None, compress_type=0, force=False)

source code 

The PDB writing function.

Parameters:
  • file (str or file instance) - The name of the PDB file to write. This can also be a file instance.
  • dir (str or None) - The directory where the PDB file will be placed. If set to None, then the file will be placed in the current directory.
  • model_num (None or int) - The model to place into the PDB file. If not supplied, then all models will be placed into the file.
  • compress_type (int) - The compression type. The integer values correspond to the compression type: 0, no compression; 1, Bzip2 compression; 2, Gzip compression.
  • force (bool) - The force flag which if True will cause the file to be overwritten.

Variables Details [hide private]

ds

Value:
The relax data storage object.

Data pipes:
  None

Data store objects:
  __dict__ <type 'dict'>: dict() -> new empty dictionary
  __doc__ <type 'str'>: The relax data storage object.
...