Package gromacs :: Module cbook :: Class Transformer
[hide private]
[frames] | no frames]

Class Transformer

source code

         object --+    
                  |    
utilities.FileUtils --+
                      |
                     Transformer

Class to handle transformations of trajectories.

1. Center, compact, and fit to reference structure in tpr
   (optionally, only center in the xy plane): :meth:`~Transformer.center_fit`
2. Write compact xtc and tpr with water removed: :meth:`~Transformer.strip_water`
3. Write compact xtc and tpr only with protein: :meth:`~Transformer.keep_protein_only`

Instance Methods [hide private]
 
__init__(self, s='topol.tpr', f='traj.xtc', n=None, force=None, dirname='.')
Set up Transformer with structure and trajectory.
source code
 
rp(self, *args)
Return canonical path to file under dirname with components args
source code
 
center_fit(self, **kwargs)
Write compact xtc that is fitted to the tpr reference structure.
source code
 
fit(self, xy=False, **kwargs)
Write xtc that is fitted to the tpr reference structure.
source code
 
strip_water(self, os=None, o=None, on=None, compact=False, resn='SOL', groupname='notwater', **kwargs)
Write xtc and tpr with water (by resname) removed.
source code
 
keep_protein_only(self, os=None, o=None, on=None, compact=False, groupname='proteinonly', **kwargs)
Write xtc and tpr only containing the protein.
source code

Inherited from utilities.FileUtils: __repr__, check_file_exists, filename, infix_filename

Inherited from utilities.FileUtils (private): _init_filename

Inherited from object: __delattr__, __getattribute__, __hash__, __new__, __reduce__, __reduce_ex__, __setattr__, __str__

Class Variables [hide private]

Inherited from utilities.FileUtils: default_extension

Instance Variables [hide private]

Inherited from utilities.FileUtils: real_filename

Properties [hide private]

Inherited from object: __class__

Method Details [hide private]

__init__(self, s='topol.tpr', f='traj.xtc', n=None, force=None, dirname='.')
(Constructor)

source code 
Set up Transformer with structure and trajectory.

Supply *n* = tpr, *f* = xtc (and *n* = ndx) relative to dirname.

:Keywords:
   *s* 
      tpr file (or similar); note that this should not contain
      position restraints if it is to be used with a reduced
      system (see :meth:`~Transformer.strip_water`)
   *f*
      trajectory (xtc, trr, ...)
   *n*
      index file (it is typically safe to leave this as ``None``; in
      cases where a trajectory needs to be centered on non-standard
      groups this should contain those groups)
   *force*
      Set the default behaviour for handling existing files:
        - ``True``: overwrite existing trajectories
        - ``False``: throw a IOError exception
        - ``None``: skip existing and log a warning [default]

Overrides: object.__init__

rp(self, *args)

source code 

Return canonical path to file under dirname with components args

If args form an absolute path then just return it as the absolute path.

center_fit(self, **kwargs)

source code 

Write compact xtc that is fitted to the tpr reference structure.

See :func:gromacs.cbook.trj_fitandcenter` for details and description of kwargs. The most important ones are listed here but in most cases the defaults should work.

Returns:
dictionary with keys tpr, xtc, which are the names of the the new files

Keywords:

s
Input structure (typically the default tpr file but can be set to some other file with a different conformation for fitting)
n
Alternative index file.
o
Name of the output trajectory.
xy : Boolean
If True then only fit in xy-plane (useful for a membrane normal to z). The default is False.
force
  • True: overwrite existing trajectories
  • False: throw a IOError exception
  • None: skip existing and log a warning [default]

fit(self, xy=False, **kwargs)

source code 
Write xtc that is fitted to the tpr reference structure.

See :func:`gromacs.cbook.trj_xyfitted` for details and
description of *kwargs*. The most important ones are listed
here but in most cases the defaults should work.

:Keywords:
   *s*
     Input structure (typically the default tpr file but can be set to
     some other file with a different conformation for fitting)
   *n*
     Alternative index file.           
   *o*
     Name of the output trajectory. A default name is created. 
     If e.g. *dt* = 100  is one of the *kwargs* then the default name includes
     "_dt100ps".
  *xy* : boolean
     If ``True`` then only do a rot+trans fit in the xy plane
     (good for membrane simulations); default is ``False``.
  *force*
    ``True``: overwrite existing trajectories
    ``False``: throw a IOError exception
    ``None``: skip existing and log a warning [default]
  *kwargs*
     kwargs are passed to :func:`~gromacs.cbook.trj_xyfitted`           

:Returns: 
      dictionary with keys *tpr*, *xtc*, which are the names of the
      the new files              

strip_water(self, os=None, o=None, on=None, compact=False, resn='SOL', groupname='notwater', **kwargs)

source code 
Write xtc and tpr with water (by resname) removed.

:Keywords:
   *os*
      Name of the output tpr file; by default use the original but
      insert "nowater" before suffix.
   *o*
      Name of the output trajectory; by default use the original name but
      insert "nowater" before suffix.
   *on*
      Name of a new index file (without water).
   *compact*
      ``True``: write a compact and centered trajectory
      ``False``: use trajectory as it is [``False``]
   *resn*
      Residue name of the water molecules; all these residues are excluded.
   *groupname*
      Name of the group that is generated by subtracting all waters
      from the system.
   *force* : Boolean
     - ``True``: overwrite existing trajectories
     - ``False``: throw a IOError exception
     - ``None``: skip existing and log a warning [default]           
   *kwargs* 
      are passed on to :func:`gromacs.cbook.trj_compact` (unless the
      values have to be set to certain values such as s, f, n, o
      keywords). The *input* keyword is always mangled: Only the first
      entry (the group to centre the trajectory on) is kept, and as a
      second group (the output group) *groupname* is used.

:Returns: 
      dictionary with keys *tpr*, *xtc*, *ndx* which are the names of the
      the new files
      
.. warning:: The input tpr file should *not* have *any position restraints*;
             otherwise Gromacs will throw a hissy-fit and say 

             *Software inconsistency error: Position restraint coordinates are
             missing*

             (This appears to be a bug in Gromacs 4.x.)

keep_protein_only(self, os=None, o=None, on=None, compact=False, groupname='proteinonly', **kwargs)

source code 
Write xtc and tpr only containing the protein.

:Keywords:
   *os*
      Name of the output tpr file; by default use the original but
      insert "proteinonly" before suffix.
   *o*
      Name of the output trajectory; by default use the original name but
      insert "proteinonly" before suffix.
   *on*
      Name of a new index file.
   *compact*
      ``True``: write a compact and centered trajectory
      ``False``: use trajectory as it is [``False``]
   *groupname*
      Name of the protein-only group.
   *keepalso*
      List of literal make_ndx selections of additional groups that should 
      be kept, e.g. ['resname DRUG', 'atom 6789'].
   *force* : Boolean
     - ``True``: overwrite existing trajectories
     - ``False``: throw a IOError exception
     - ``None``: skip existing and log a warning [default]           
   *kwargs* 
      are passed on to :func:`gromacs.cbook.trj_compact` (unless the
      values have to be set to certain values such as s, f, n, o
      keywords). The *input* keyword is always mangled: Only the first
      entry (the group to centre the trajectory on) is kept, and as a
      second group (the output group) *groupname* is used.

:Returns: 
      dictionary with keys *tpr*, *xtc*, *ndx* which are the names of the
      the new files
      
.. warning:: The input tpr file should *not* have *any position restraints*;
             otherwise Gromacs will throw a hissy-fit and say 

             *Software inconsistency error: Position restraint coordinates are
             missing*

             (This appears to be a bug in Gromacs 4.x.)