Line data    Source code

       1             : /* +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
       2             :    Copyright (c) 2012-2023 The plumed team
       3             :    (see the PEOPLE file at the root of the distribution for a list of names)
       4             : 
       5             :    See http://www.plumed.org for more information.
       6             : 
       7             :    This file is part of plumed, version 2.
       8             : 
       9             :    plumed is free software: you can redistribute it and/or modify
      10             :    it under the terms of the GNU Lesser General Public License as published by
      11             :    the Free Software Foundation, either version 3 of the License, or
      12             :    (at your option) any later version.
      13             : 
      14             :    plumed is distributed in the hope that it will be useful,
      15             :    but WITHOUT ANY WARRANTY; without even the implied warranty of
      16             :    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
      17             :    GNU Lesser General Public License for more details.
      18             : 
      19             :    You should have received a copy of the GNU Lesser General Public License
      20             :    along with plumed.  If not, see <http://www.gnu.org/licenses/>.
      21             : +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ */
      22             : #include "core/ActionRegister.h"
      23             : #include "core/GenericMolInfo.h"
      24             : 
      25             : namespace PLMD {
      26             : namespace generic {
      27             : 
      28             : //+PLUMEDOC TOPOLOGY MOLINFO
      29             : /*
      30             : This command is used to provide information on the molecules that are present in your system.
      31             : 
      32             : The information on the molecules in your system can either be provided in the form of a pdb file
      33             : or as a set of lists of atoms that describe the various chains in your system. If a pdb file
      34             : is used plumed the MOLINFO command will endeavor to recognize the various chains and residues that
      35             : make up the molecules in your system using the chainIDs and resnumbers from the pdb file. You can
      36             : then use this information in later commands to specify atom lists in terms residues.  For example
      37             : using this command you can find the backbone atoms in your structure automatically.
      38             : Starting with PLUMED 2.7 you can use multiple MOLINFO actions. Every time you perform an atom
      39             : selection, the last available MOLINFO action will be used. This allows you to provide multiple PDB
      40             : files, for instance using different naming conventions (see \issue{134}).
      41             : 
      42             : \warning
      43             : Please be aware that the PDB parser in plumed is far from perfect. You should thus check the log file
      44             : and examine what plumed is actually doing whenever you use the MOLINFO action.
      45             : Also make sure that the atoms are listed in the pdb with the correct order.
      46             : If you are using gromacs, the safest way is to use reference pdb file
      47             : generated with `gmx editconf -f topol.tpr -o reference.pdb`.
      48             : 
      49             : More information of the PDB parser implemented in PLUMED can be found \ref pdbreader "at this page".
      50             : 
      51             : Providing `MOLTYPE=protein`, `MOLTYPE=rna`, or `MOLTYPE=dna` will instruct plumed to look
      52             : for known residues from these three types of molecule. In other words, this is available for
      53             : historical reasons and to allow future extensions where alternative lists will be provided.
      54             : As of now, you can just ignore this keyword.
      55             : 
      56             : Using \ref MOLINFO extends the possibility of atoms selection using the @ special
      57             : symbol. The following shortcuts are available that do not refer to one specific residue:
      58             : 
      59             : \verbatim
      60             : @nucleic : all atoms that are part of a DNA or RNA molecule
      61             : @protein : all atoms that are part of a protein
      62             : @water : all water molecules
      63             : @ions : all the ions
      64             : @hydrogens : all hydrogen atoms (those for which the first non-number in the name is a H)
      65             : @nonhydrogens : all non hydrogen atoms (those for which the first non-number in the name is not a H)
      66             : \endverbatim
      67             : 
      68             : \warning
      69             : Be careful since these choices are based on common names used in PDB files. Always check if
      70             : the selected atoms are correct.
      71             : 
      72             : In addition, atoms from a specific residue can be selected with a symbol in this form:
      73             : 
      74             : \verbatim
      75             : @"definition"-chain_residuenum
      76             : @"definition"-chainresiduenum
      77             : @"definition"-residuenum
      78             : \endverbatim
      79             : 
      80             : So for example
      81             : 
      82             : \verbatim
      83             : @psi-1 will select the atoms defining the psi torsion of residue 1
      84             : @psi-C1  or @psi-C_1 will define the same torsion for residue 1 of chain C.
      85             : @psi-3_1 will define the same torsion for residue 1 of chain 3.
      86             : \endverbatim
      87             : 
      88             : Using the underscore to separate chain and residue is available as of PLUMED 2.5 and allows selecting chains
      89             : with a numeric id.
      90             : 
      91             : In the following are listed the current available definitions:
      92             : 
      93             : For protein residues, the following groups are available:
      94             : 
      95             : \verbatim
      96             : # quadruplets for dihedral angles
      97             : @phi-#
      98             : @psi-#
      99             : @omega-#
     100             : @chi1-#
     101             : @chi2-#
     102             : @chi3-#
     103             : @chi4-#
     104             : @chi5-#
     105             : 
     106             : # all sidechain atoms (excluding glycine, including all hydrogens)
     107             : @sidechain-#
     108             : # all backbone atoms (including hydrogens)
     109             : @back-#
     110             : \endverbatim
     111             : 
     112             : that select the appropriate atoms that define each dihedral angle for residue #.
     113             : 
     114             : For DNA or RNA residues, the following groups are available:
     115             : 
     116             : \verbatim
     117             : # quadruplets for backbone dihedral angles
     118             : @alpha-#
     119             : @beta-#
     120             : @gamma-#
     121             : @delta-#
     122             : @epsilon-#
     123             : @zeta-#
     124             : 
     125             : # quadruplets for sugar dihedral angles
     126             : @v0-#
     127             : @v1-#
     128             : @v2-#
     129             : @v3-#
     130             : @v4-#
     131             : 
     132             : # quadruplet corresponding to the chi torsional angle
     133             : @chi-#
     134             : 
     135             : # backbone, sugar, and base heavy atoms
     136             : @back-#
     137             : @sugar-#
     138             : @base-#
     139             : 
     140             : # ordered triplets of atoms on the 6-membered ring of nucleobases
     141             : # namely:
     142             : #  C2/C4/C6 for pyrimidines
     143             : #  C2/C6/C4 for purines
     144             : @lcs-#
     145             : \endverbatim
     146             : 
     147             : Notice that `zeta` and `epsilon` groups should not be used on 3' end residue and `alpha` and `beta`
     148             : should not be used on 5' end residue.
     149             : 
     150             : Furthermore it is also possible to pick single atoms using the syntax
     151             : `atom-chain_residuenum`, `@atom-chainresiduenum` or `@atom-residuenum`.
     152             : As of PLUMED 2.5, this also works when the residue is not a protein/rna/dna residue.
     153             : For instance, `@OW-100` will select oxygen of water molecule with residue number 100.
     154             : 
     155             : Finally, notice that other shortcuts are available even when not using the \ref MOLINFO command (see \ref atomSpecs).
     156             : 
     157             : \warning If a residue-chain is repeated twice in the reference pdb only the first entry will be selected.
     158             : 
     159             : \bug At the moment the HA1 atoms in a GLY residues are treated as if they are the CB atoms. This may or
     160             : may not be true - GLY is problematic for secondary structure residues as it is achiral.
     161             : 
     162             : \bug If you use WHOLEMOLECULES RESIDUES=1-10 for a 18 amino acid protein
     163             : ( 18 amino acids + 2 terminal groups = 20 residues ) the code will fail as it will not be able to
     164             : interpret terminal residue 1.
     165             : 
     166             : \par Advanced atom selection with mdtraj or MDAnalysis
     167             : 
     168             : Since PLUMED 2.6 it is possible to use the expressive selection syntax of [mdtraj](http://mdtraj.org/latest/atom_selection.html) and/or [MDAnalysis](https://www.mdanalysis.org/docs/documentation_pages/selections.html):
     169             : 
     170             : \plumedfile
     171             : MOLINFO STRUCTURE=helix.pdb PYTHON_BIN=python
     172             : g1: GROUP ATOMS=@mda:backbone
     173             : g2: GROUP ATOMS={@mda:{resnum 1 or resid 3:5}}
     174             : g3: GROUP ATOMS={@mda:{resid 3:5} @mda:{resnum 1}}
     175             : g4: GROUP ATOMS={@mdt:{protein and (backbone or resname ALA)}}
     176             : g5: GROUP ATOMS={@mdt:{mass 5.5 to 20}} # masses guessed by mdtraj based on atom type!
     177             : g6: GROUP ATOMS={@mda:{resid 3:5} @mda:{resnum 1} 1-10}
     178             : \endplumedfile
     179             : 
     180             : Here `@mda:` indicates that `MDAnalysis` language is used, whereas `@mdt:` indicates that `mdtraj` language is used. Notice that these languages typically select atoms in order. If you want to specify a different order, you can chain definitions as in `g3` above (compare with `g2`). Selections can be also chained with standard PLUMED selections (see `g6`).
     181             : 
     182             : The double braces are required due to the way PLUMED parses atom lists. In particular:
     183             : 
     184             : - The outer braces are needed to show PLUMED where the `ATOMS=...` option ends.
     185             : - The inner braces are needed to show PLUMED where each selector ends.
     186             : 
     187             : MDAnalysis also supports geometric selectors based on atomic coordinates. These selectors **are static** and return lists computed using the coordinates stored in the `MOLINFO` pdb file.
     188             : 
     189             : In order to use this syntax you should check the following points at runtime:
     190             : 
     191             : 1. `plumed --no-mpi config has subprocess` prints `subprocess on` (should be ok on most UNIX systems).
     192             : 2. You have a python interpreter with mdtraj and/or MDAnalysis installed. You can check using:
     193             :    - `python -c "import mdtraj"`
     194             :    - `python -c "import MDAnalysis"`
     195             : 
     196             :    In order to install these packages refer to their documentation. Pip or conda install should be ok, provided you make sure the correct python interpreter is in the execution PATH at runtime. Notice that you will only need the package(s) related to the syntax that you want to use.
     197             : 3. In case you installed these modules on a python with a different name (e.g. `python3.6`), the correct check is:
     198             :    - `python3.6 -c "import mdtraj"`
     199             :    - `python3.6 -c "import MDAnalysis"`
     200             : 
     201             :    If this is the case, you should set the environment variable `export PYTHON_BIN=python3.6` or `export PLUMED_PYTHON_BIN=python3.6` (higher priority). Alternatively, directly provide the interpreter in the PLUMED input file
     202             :    using `MOLINFO PYTHON_BIN=python3.6` (even higher priority).
     203             : 4. The PDB file that you provide to `MOLINFO` should have consecutive atom numbers starting from 1. This is currently enforced since reading atom numbers out of order (as PLUMED does) is not supported by other packages.
     204             : 
     205             : \par Advanced atom selection with VMD (experimental)
     206             : 
     207             : Similarly to the `@mda:` and `@mdt:` selectors above, you can use the two following selectors in order to
     208             : access to [VMD](https://www.ks.uiuc.edu/Research/vmd/) syntax for atoms selection:
     209             : - `@vmdexec:`: This selector launches an instance of VMD, so `vmd` executable should be in your execution path.
     210             :   Might be very slow or even crash your simulation. Notice that even if `vmd` executable is used,
     211             :   the implementation is still python based and so a working python interpreter should be provided.
     212             : - `@vmd:`: This selector tries to import the `vmd` python module. Notice that the best way to obtain this module
     213             :   is not within the standard VMD installer but rather by installing the python
     214             :   module that can be found at [this link](http://github.com/Eigenstate/vmd-python).
     215             :   The module is also available on [conda](https://anaconda.org/conda-forge/vmd-python).
     216             :   You should make sure the module is available in the python interpreter used by MOLINFO
     217             :   (check using the command `python -c "import vmd"`).
     218             : 
     219             : These two selectors are experimental and might be removed at some point.
     220             : 
     221             : \par Examples
     222             : 
     223             : In the following example the MOLINFO command is used to provide the information on which atoms
     224             : are in the backbone of a protein to the ALPHARMSD CV.
     225             : 
     226             : \plumedfile
     227             : #SETTINGS MOLFILE=regtest/basic/rt32/helix.pdb
     228             : MOLINFO STRUCTURE=reference.pdb
     229             : ALPHARMSD RESIDUES=all TYPE=DRMSD LESS_THAN={RATIONAL R_0=0.08 NN=8 MM=12} LABEL=a
     230             : \endplumedfile
     231             : 
     232             : The following example prints the distance corresponding to the hydrogen bonds
     233             : in a GC Watson-Crick pair.
     234             : 
     235             : \plumedfile
     236             : #SETTINGS MOLFILE=regtest/basic/rt-ermsd/ref.pdb
     237             : MOLINFO STRUCTURE=reference.pdb MOLTYPE=dna
     238             : hb1: DISTANCE ATOMS=@N2-2,@O2-15
     239             : hb2: DISTANCE ATOMS=@N1-2,@N3-15
     240             : hb3: DISTANCE ATOMS=@O6-2,@N4-15
     241             : PRINT ARG=hb1,hb2,hb3
     242             : \endplumedfile
     243             : 
     244             : This example use MOLINFO to calculate torsion angles
     245             : 
     246             : \plumedfile
     247             : #SETTINGS MOLFILE=regtest/basic/rt32/helix.pdb
     248             : MOLINFO MOLTYPE=protein STRUCTURE=myprotein.pdb
     249             : t1: TORSION ATOMS=@phi-3
     250             : t2: TORSION ATOMS=@psi-4
     251             : PRINT ARG=t1,t2 FILE=colvar STRIDE=10
     252             : \endplumedfile
     253             : 
     254             : */
     255             : //+ENDPLUMEDOC
     256             : 
     257             : 
     258             : /*
     259             : This action is defined in core/ as it is used by other actions.
     260             : Anyway, it is registered here, so that excluding this module from
     261             : compilation will exclude it from plumed.
     262             : */
     263             : 
     264             : typedef PLMD::GenericMolInfo MolInfo;
     265             : 
     266       10587 : PLUMED_REGISTER_ACTION(MolInfo,"MOLINFO")
     267             : 
     268             : }
     269             : }