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

Class IndexBuilder

source code

object --+
         |
        IndexBuilder

Build an index file with specified groups and the combined group.

This is *not* a full blown selection parser a la Charmm, VMD or
MDAnalysis but a very quick hack.

**Example**

   How to use the :class:`IndexBuilder`::

      G = gromacs.cbook.IndexBuilder('md_posres.pdb', 
                    ['S312:OG','T313:OG1','A38:O','A309:O','@a62549 & r NA'], 
                    offset=-9, out_ndx='selection.ndx')
      groupname, ndx = G.combine()
      del G

   The residue numbers are given with their canonical resids from the
   sequence or pdb. *offset=-9* says that one calculates Gromacs topology
   resids by subtracting 9 from the canonical resid. 

   The combined selection is ``OR`` ed by default and written to
   *selection.ndx*. One can also add all the groups in the initial *ndx*
   file (or the :program:`make_ndx` default groups) to the output (see the
   *defaultgroups* keyword for :meth:`IndexBuilder.combine`).

   Generating an index file always requires calling
   :meth:`~IndexBuilder.combine` even if there is only a single group.

   Deleting the class removes all temporary files associated with it (see
   :attr:`IndexBuilder.indexfiles`).

:Raises:
   If an empty group is detected (which does not always work) then a
   :exc:`gromacs.BadParameterWarning` is issued.

:Bugs:
   If ``make_ndx`` crashes with an unexpected error then this is fairly hard to
   diagnose. For instance, in certain cases it segmentation faults when a tpr
   is provided as a *struct* file and the resulting error messages becomes ::

      GromacsError: [Errno -11] Gromacs tool failed
      Command invocation: make_ndx -o /tmp/tmp_Na1__NK7cT3.ndx -f md_posres.tpr

   In this case run the command invocation manually to see what the problem
   could be.


.. SeeAlso:: In some cases it might be more straightforward to use
             :class:`gromacs.formats.NDX`.

Instance Methods [hide private]
 
__init__(self, struct=None, selections=None, names=None, name_all=None, ndx=None, out_ndx='selection.ndx', offset=0)
Build a index group from the selection arguments.
source code
 
gmx_resid(self, resid)
Returns resid in the Gromacs index by transforming with offset.
source code
 
combine(self, name_all=None, out_ndx=None, operation='|', defaultgroups=False)
Combine individual groups into a single one and write output.
source code
 
cat(self, out_ndx=None)
Concatenate input index files.
source code
 
parse_selection(self, selection, name=None)
Retuns (groupname, filename) with index group.
source code
 
_process_command(self, command, name=None)
Process make_ndx command and return name and temp index file.
source code
 
_process_residue(self, selection, name=None)
Process residue/atom selection and return name and temp index file.
source code
 
_process_range(self, selection, name=None)
Process a range selection.
source code
 
_translate_residue(self, selection, default_atomname='CA')
Translate selection for a single res to make_ndx syntax.
source code
 
check_output(self, make_ndx_output, message=None, err=None)
Simple tests to flag problems with a make_ndx run.
source code
 
_is_empty_group(self, make_ndx_output) source code
 
_has_syntax_error(self, make_ndx_output) source code
 
__del__(self) source code

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

Class Variables [hide private]
  RESIDUE = re.compile(r'(?x)(?P<aa>([ACDEFGHIKLMNPQRSTVWY])|(\S...
regular expression to match and parse a residue-atom selection
Instance Variables [hide private]
  offset
*offset* as a number is added to the resids in the first selection scheme; this allows names to be the same as in a crystal structure.
  _command_counter
Auto-labelled groups use this counter.
  make_ndx
Specialized make_ndx that always uses same structure and redirection (can be overridden)
  indexfiles
dict, keyed by group name and pointing to index file for group (Groups are built in separate files because that is more robust as I can clear groups easily.)
Properties [hide private]
  names
Names of all generated index groups.

Inherited from object: __class__

Method Details [hide private]

__init__(self, struct=None, selections=None, names=None, name_all=None, ndx=None, out_ndx='selection.ndx', offset=0)
(Constructor)

source code 
Build a index group from the selection arguments.

If selections and a structure file are supplied then the individual
selections are constructed with separate calls to
:func:`gromacs.make_ndx`. Use :meth:`IndexBuilder.combine` to combine
them into a joint selection.

:Arguments:

   struct : filename
      Structure file (tpr, pdb, ...)

   selections : list
      The list must contain strings or tuples, which must be be one of
      the following constructs:

         "<1-letter aa code><resid>[:<atom name]"

             Selects the CA of the residue or the specified atom
             name.

             example: ``"S312:OA"`` or ``"A22"`` (equivalent to ``"A22:CA"``)
         
         ("<1-letter aa code><resid>", "<1-letter aa code><resid>, ["<atom name>"])

            Selects a *range* of residues. If only two residue
            identifiers are provided then all atoms are
            selected. With an optional third atom identifier,
            only this atom anme is selected for each residue
            in the range. [EXPERIMENTAL]

         "@<make_ndx selection>"

             The ``@`` letter introduces a verbatim ``make_ndx``
             command. It will apply the given selection without any
             further processing or checks.

             example: ``"@a 6234 - 6238"`` or ``'@"SOL"'`` (note the quoting)
             or ``"@r SER & r 312 & t OA"``.

   names : list
      Strings to name the selections; if not supplied or if individuals
      are ``None`` then a default name is created.

   offset : int, dict
      This number is added to the resids in the first selection scheme; this
      allows names to be the same as in a crystal structure. If offset is a 
      dict then it is used to directly look up the resids.

   ndx : filename or list of filenames
      Optional input index file(s).

   out_ndx : filename
      Output index file.
      

Overrides: object.__init__

combine(self, name_all=None, out_ndx=None, operation='|', defaultgroups=False)

source code 
Combine individual groups into a single one and write output.

:Keywords:         
   name_all : string
      Name of the combined group, ``None`` generates a name.  [``None``]
   out_ndx : filename
      Name of the output file that will contain the individual groups
      and the combined group. If ``None`` then default from the class
      constructor is used. [``None``]
   operation : character
      Logical operation that is used to generate the combined group from
      the individual groups: "|" (OR) or "&" (AND) ["|"]
   defaultgroups : bool
      ``True``: append everything to the default groups produced by 
      :program:`make_ndx` (or rather, the groups provided in the ndx file on
      initialization --- if this was ``None`` then these are truly default groups);
      ``False``: only use the generated groups

:Returns:
   ``(combinedgroup_name, output_ndx)``, a tuple showing the
   actual group name and the name of the file; useful when all names are autogenerated.

.. Warning:: The order of the atom numbers in the combined group is
             *not* guaranteed to be the same as the selections on input because
             ``make_ndx`` sorts them ascending. Thus you should be careful when
             using these index files for calculations of angles and dihedrals.
             Use :class:`gromacs.formats.NDX` in these cases.

cat(self, out_ndx=None)

source code 

Concatenate input index files.

Generate a new index file that contains the default Gromacs index groups (if a structure file was defined) and all index groups from the input index files.

Parameters:
  • out_ndx (filename) - Name of the output index file; if None then use the default provided to the constructore. [None].

_process_range(self, selection, name=None)

source code 

Process a range selection.

("S234", "A300", "CA") --> selected all CA in this range ("S234", "A300") --> selected all atoms in this range

Note

Ignores residue type, only cares about the resid (but still required)


Class Variable Details [hide private]

RESIDUE

regular expression to match and parse a residue-atom selection
Value:
re.compile(r'(?x)(?P<aa>([ACDEFGHIKLMNPQRSTVWY])|(\S\S\S))(?P<resid>\d\
+)(:(?P<atom>\w+))?')

Instance Variable Details [hide private]

offset

*offset* as a number is added to the resids in the first selection
scheme; this
allows names to be the same as in a crystal structure. If *offset* is a
dict then it is used to directly look up the resids. Use :meth:`gmx_resid`
to transform a crystal resid to a gromacs resid.

The attribute may be changed directly after init.


Property Details [hide private]

names

Names of all generated index groups.
Get Method:
unreachable.names(self) - Names of all generated index groups.