Pteros  2.0
Molecular modeling library for human beings!
pteros::System Class Reference

The system of atoms. More...

#include <system.h>

Public Member Functions

Constructors and operators
 System ()
 Default constructor.
 
 System (std::string fname)
 Constructor creating system from file.
 
 System (const System &other)
 Copy constructor.
 
Systemoperator= (System other)
 Assignment operator.
 
virtual ~System ()
 Destructor.
 
Adding, deleting and ordering groups of atoms

These methods update resindexes automatically.

Selection append (const System &sys)
 Append other system to this one Returns selection corresponding to appended atoms.
 
Selection append (const Selection &sel)
 Append atoms from selection to this system. More...
 
Selection append (const Atom &at, Vector3f_const_ref coord)
 Append single atom to this system Returns selection corresponding to appended atom.
 
Selection append (const Atom_proxy &at)
 Append Atom_proxy object to this system Returns selection corresponding to appended atom. More...
 
void rearrange (const std::vector< std::string > &sel_strings)
 Rearranges the atoms in the order of provided selection strings. More...
 
void rearrange (const std::vector< Selection > &sel_vec)
 Rearranges the atoms in the order of provided selections. More...
 
void keep (const std::string &sel_str)
 Keep only atoms given by selection string.
 
void keep (const Selection &sel)
 Keep only atoms from given selection.
 
void remove (const std::string &sel_str)
 Remove atoms given by selection string.
 
void remove (Selection &sel)
 Remove atoms of given selection. More...
 
void distribute (const Selection sel, Vector3i_const_ref ncopies, Matrix3f_const_ref shift)
 Creates multiple copies of selection in the system and distributes them in a grid given by 3 translation vectors. More...
 
General properties
int num_frames () const
 Returns the number of frames in selection.
 
int num_atoms () const
 Returns the number of atoms in selection.
 
Selecting atoms.

These functions are just convenience adaptors for constructors of Selection class:

// Conventional way:
sys = System("structure.pdb");
Selection sel(s,"name CA");
// Shorter way:
auto sel = sys.select("name CA");
// Even shorter way using operator ():
auto sel = sys("name CA");

It also allows writing "one-liners" like this:

System("file.pdb").select("name CA").write("ca.pdb");
Selection select (std::string str, int fr=0)
 
Selection operator() (std::string str, int fr=0)
 
Selection select (int ind1, int ind2)
 
Selection operator() (int ind1, int ind2)
 
Selection select (const std::vector< int > &ind)
 
Selection operator() (const std::vector< int > &ind)
 
Selection select (std::vector< int >::iterator it1, std::vector< int >::iterator it2)
 
Selection operator() (std::vector< int >::iterator it1, std::vector< int >::iterator it2)
 
Selection select (const std::function< void(const System &, int, std::vector< int > &)> &callback, int fr=0)
 
Selection operator() (const std::function< void(const System &, int, std::vector< int > &)> &callback, int fr=0)
 
Selection select_all ()
 Convenience functions to select all.
 
Selection operator() ()
 
File IO
void load (std::string fname, int b=0, int e=-1, int skip=0, std::function< bool(System *, int)> on_frame=0)
 Read structure, trajectory or topology from file. More...
 
bool load (const std::unique_ptr< Mol_file > &handler, Mol_file_content what, std::function< bool(System *, int)> on_frame=0)
 Load data into System from the pre-opened file handler. More...
 
Input filtering
void set_filter (std::string str)
 Filters narrow set of atoms and coordinates which are loaded from data files. More...
 
void set_filter (int ind1, int ind2)
 
void set_filter (const std::vector< int > &ind)
 
void set_filter (std::vector< int >::iterator it1, std::vector< int >::iterator it2)
 
Operations with frames
int frame_dup (int fr)
 Duplicates given frame and adds it to the end of frame vector.
 
void frame_append (const Frame &fr)
 Appends provided frame to trajectory.
 
void frame_copy (int fr1, int fr2)
 Copy all frame data from fr1 to fr2. Fr2 is overwritten!
 
void frame_delete (int b=0, int e=-1)
 Delete specified range of frames. More...
 
void frame_swap (int fr1, int fr2)
 Swaps two specified frames.
 
Inline accessors
Periodic_boxBox (int fr)
 Read/write access for periodic box for given frame.
 
const Periodic_boxBox (int fr) const
 Read only access for periodic box for given frame.
 
float & Time (int fr)
 Read/Write access to the time stamp of given frame.
 
const float & Time (int fr) const
 Read only access to the time stamp of given frame.
 
Eigen::Vector3f & XYZ (int ind, int fr)
 Read/Write access for given coordinate of given frame.
 
const Eigen::Vector3f & XYZ (int ind, int fr) const
 Read only access for given coordinate of given frame.
 
AtomAtom_data (int ind)
 Read/Write access for given atom.
 
const AtomAtom_data (int ind) const
 Read only access for given atom.
 
FrameFrame_data (int fr)
 Get read/write reference for given frame.
 
const FrameFrame_data (int fr) const
 Get read only reference for given frame.
 
Secondary structure functions
void dssp (std::string fname) const
 Determines secondary structure with DSSP algorithm and writes detailed report to file.
 
void dssp (std::ostream &os) const
 Determines secondary structure with DSSP algorithm and writes detailed report to stream.
 
std::string dssp () const
 Determines secondary structure with DSSP algorithm and return it as a code string. More...
 
Manipulating sets of atoms by indexes.

These methods do not update resindexes automatically.

Selection atoms_dup (const std::vector< int > &ind)
 Adds new atoms, which are duplicates of existing ones by index.
 
Selection atoms_add (const std::vector< Atom > &atm, const std::vector< Eigen::Vector3f > &crd)
 Adds new atoms from supplied vectors of atoms and coordinates.
 
void atoms_delete (const std::vector< int > &ind)
 Delete the set of atoms by indexes.
 
Periodicity-related functions
void wrap (int fr, Vector3i_const_ref dims_to_wrap=Eigen::Vector3i::Ones())
 Wrap all system to the periodic box for given frame.
 
Measuring functions
float distance (int i, int j, int fr, bool is_periodic=true, Vector3i_const_ref dims=Eigen::Vector3i::Ones()) const
 Get distance between two atoms for given frame (periodic in given dimensions if needed).
 
float angle (int i, int j, int k, int fr, bool is_periodic=true, Vector3i_const_ref dims=Eigen::Vector3i::Ones()) const
 Get angle in degrees between three atoms for given frame (periodic in given dimensions if needed).
 
float dihedral (int i, int j, int k, int l, int fr, bool is_periodic=true, Vector3i_const_ref dims=Eigen::Vector3i::Ones()) const
 Get dihedral angle in degrees between three atoms for given frame (periodic in given dimensions if needed).
 
Energy functions
Energy_components non_bond_energy (int at1, int at2, int frame, bool is_periodic=true) const
 Compute non-bond energy between two atoms.
 
Energy_components non_bond_energy (const std::vector< Eigen::Vector2i > &nlist, int fr, bool is_periodic=true) const
 Non-bond energy for given list of atom pairs.
 
Utility functions
void clear ()
 Clears the system and prepares for loading completely new structure.
 
bool force_field_ready ()
 Returns true if the force field is set up properly and is able to compute energies.
 
Force_fieldget_force_field ()
 Returns pointer to internal Force_field object for direct manipulation or NULL if force field is not ready.
 
void assign_resindex (int start=0)
 Assign unique resindexes This is usually done automatically upon loading a structure from file.
 
void sort_by_resindex ()
 Sorts atoms by resindex arranging atoms with the same resindexes into contigous pieces. More...
 

Detailed Description

The system of atoms.

The System is a container for atoms and their coordinates, which are typically loaded from file. All properties of atoms, except the coordinates, are stored in atoms vector. Coordinates are stored as a resizable vector of trajectory frames. The system knows about selections associated with it and sends them signals if selections should adapt to the changes of coordinates of atom properties. Copying and assignment of systems is allowed, but associated selections are not copied.

Member Function Documentation

Selection System::append ( const Selection sel)

Append atoms from selection to this system.

Selection may belong to any system, not necessary current. Returns selection corresponding to appended atoms.

Selection System::append ( const Atom_proxy at)

Append Atom_proxy object to this system Returns selection corresponding to appended atom.

Usage:

Selection sel(s,"name CA");
System new_s;
for(auto& at: sel){
new_s.append(at);
}
void System::rearrange ( const std::vector< std::string > &  sel_strings)

Rearranges the atoms in the order of provided selection strings.

Atom, which are not selected are appended at the end in their previous order.

Note
Selections should not overlap (exception is thrown if they are).
void System::rearrange ( const std::vector< Selection > &  sel_vec)

Rearranges the atoms in the order of provided selections.

Atom, which are not selected are appended at the end in their previous order.

Note
Selections should not overlap (exception is thrown if they are).
void System::remove ( Selection sel)

Remove atoms of given selection.

Warning
Selection becomes invalid after that and is cleared!
void System::distribute ( const Selection  sel,
Vector3i_const_ref  ncopies,
Matrix3f_const_ref  shift 
)

Creates multiple copies of selection in the system and distributes them in a grid given by 3 translation vectors.

Vectors are stored column-wise as shift.col(i)

void System::load ( std::string  fname,
int  b = 0,
int  e = -1,
int  skip = 0,
std::function< bool(System *, int)>  on_frame = 0 
)

Read structure, trajectory or topology from file.

Parameters
bFirst frame to read
eLast frame to read (-1 means up to the end of trajectory)
skipRead only each skip frames
on_frameCallback function, which takes pointer to the System as a first argument and the index of current stored frame as the second. If callback returns false loading stops.
bool System::load ( const std::unique_ptr< Mol_file > &  handler,
Mol_file_content  what,
std::function< bool(System *, int)>  on_frame = 0 
)

Load data into System from the pre-opened file handler.

This is rather low-level method which provides fine control over what should be read. It can be called several times to read trajectory frames one by one from the same pre-opened file.

void System::set_filter ( std::string  str)

Filters narrow set of atoms and coordinates which are loaded from data files.

Only atoms specified by filter are kept in the system. They follow the same rules as ordinary selections but can't be coordinate-dependent and can't be set by callbacs.

Note
Filters do not make loading faster. In fact it becomes slower (up to 2 times for very large frames) and consumes twise as much memory during loading. However when loading finishes the system will contain only selected atoms. If many frames are stored in the system overall memory consumption will be much smaller.

Filter could only be set to empty system (it throws an error otherwise), thus the way of using them is the following:

// Create empty system
System sys;
// Set filter
sys.set_filter("name CA");
// All calls to load will use filter now
sys.load("some_protein.pdb");
sys.load("trajectory.xtc");
Warning
Filters could not be used when loading topologies! load() throws an error if attempting to load topology with filter.
void System::frame_delete ( int  b = 0,
int  e = -1 
)

Delete specified range of frames.

If only

Parameters
bis supplied deletes all frames from b to the end. If only
eis supplied deletes all frames from 0 to e
string System::dssp ( ) const

Determines secondary structure with DSSP algorithm and return it as a code string.

Returns
Code string The code is the same as in DSSP: alphahelix: 'H' betabridge: 'B' strand: 'E' helix_3: 'G' helix_5: 'I' turn: 'T' bend: 'S' loop: ' '
void System::sort_by_resindex ( )

Sorts atoms by resindex arranging atoms with the same resindexes into contigous pieces.

Could be called after atom additions or duplications.


The documentation for this class was generated from the following files: