Line data Source code
1 : /* +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
2 : Copyright (c) 2012-2019 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/SetupMolInfo.h"
24 :
25 : namespace PLMD {
26 : namespace setup {
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 :
39 : \warning
40 : Please be aware that the PDB parser in plumed is far from perfect. You should thus check the log file
41 : and examine what plumed is actually doing whenenver you use the MOLINFO action.
42 : Also make sure that the atoms are listed in the pdb with the correct order.
43 : If you are using gromacs, the safest way is to use reference pdb file
44 : generated with `gmx editconf -f topol.tpr -o reference.pdb`.
45 :
46 : More information of the PDB parser implemented in PLUMED can be found \ref pdbreader "at this page".
47 :
48 : Providing `MOLTYPE=protein`, `MOLTYPE=rna`, or `MOLTYPE=dna` will instruct plumed to look
49 : for known residues from these three types of molecule. In other words, this is available for
50 : historical reasons and to allow future extensions where alternative lists will be provided.
51 : As of now, you can just ignore this keyoword.
52 :
53 : Using MOLINFO with a protein's or nucleic acid's pdb extends the possibility of atoms selection using the @ special
54 : symbol in the form
55 :
56 : \verbatim
57 : @"definition"-chainresiduenum
58 : @"definition"-residuenum
59 : \endverbatim
60 :
61 : So for example
62 :
63 : \verbatim
64 : @psi-1 will select the atoms defining the psi torsion of residue 1
65 : @psi-C1 will define the same torsion for residue 1 of chain C.
66 : \endverbatim
67 :
68 : In the following are listed the current available definitions:
69 :
70 : For protein residues, the following groups are available:
71 :
72 : \verbatim
73 : @phi-#
74 : @psi-#
75 : @omega-#
76 : @chi1-#
77 : \endverbatim
78 :
79 : that select the appropriate atoms that define each dihedral angle for residue #.
80 :
81 : For DNA or RNA residues, the following groups are available:
82 :
83 : \verbatim
84 : # quadruplets for backbone dihedral angles
85 : @alpha-#
86 : @beta-#
87 : @gamma-#
88 : @delta-#
89 : @epsilon-#
90 : @zeta-#
91 :
92 : # quadruplets for sugar dihedral angles
93 : @v0-#
94 : @v1-#
95 : @v2-#
96 : @v3-#
97 : @v4-#
98 :
99 : # quadruplet corresponding to the chi torsional angle
100 : @chi-#
101 :
102 : # backbone, sugar, and base heavy atoms
103 : @back-#
104 : @sugar-#
105 : @base-#
106 :
107 : # ordered triplets of atoms on the 6-membered ring of nucleobases
108 : # namely:
109 : # C2/C4/C6 for pyrimidines
110 : # C2/C6/C4 for purines
111 : @lcs-#
112 : \endverbatim
113 :
114 : Notice that `zeta` and `epsilon` groups should not be used on 3' end residue and `alpha` and `beta`
115 : should not be used on 5' end residue.
116 :
117 : Furthermore it is also possible to pick single atoms using the syntax
118 : `@atom-chainresiduenum` or `@atom-residuenum`.
119 :
120 : \warning If a residue-chain is repeated twice in the reference pdb only the first entry will be selected.
121 :
122 : \bug At the moment the HA1 atoms in a GLY residues are treated as if they are the CB atoms. This may or
123 : may not be true - GLY is problematic for secondary structure residues as it is achiral.
124 :
125 : \bug If you use WHOLEMOLECULES RESIDUES=1-10 for a 18 amino acid protein
126 : ( 18 amino acids + 2 terminal groups = 20 residues ) the code will fail as it will not be able to
127 : interpret terminal residue 1.
128 :
129 : \par Examples
130 :
131 : In the following example the MOLINFO command is used to provide the information on which atoms
132 : are in the backbone of a protein to the ALPHARMSD CV.
133 :
134 : \plumedfile
135 : MOLINFO STRUCTURE=reference.pdb
136 : ALPHARMSD RESIDUES=all TYPE=DRMSD LESS_THAN={RATIONAL R_0=0.08 NN=8 MM=12} LABEL=a
137 : \endplumedfile
138 :
139 : The following example prints the distance corresponding to the hydrogen bonds
140 : in a GC Watson-Crick pair.
141 :
142 : \plumedfile
143 : MOLINFO STRUCTURE=reference.pdb
144 : hb1: DISTANCE ATOMS=@N2-1,@O2-14
145 : hb2: DISTANCE ATOMS=@N1-1,@N3-14
146 : hb3: DISTANCE ATOMS=@O6-1,@N4-14
147 : PRINT ARG=hb1,hb2,hb3
148 : \endplumedfile
149 :
150 : This example use MOLINFO to calculate torsions angles
151 :
152 : \plumedfile
153 : MOLINFO MOLTYPE=protein STRUCTURE=myprotein.pdb
154 : t1: TORSION ATOMS=@phi-3
155 : t2: TORSION ATOMS=@psi-4
156 : PRINT ARG=t1,t2 FILE=colvar STRIDE=10
157 : \endplumedfile
158 :
159 : */
160 : //+ENDPLUMEDOC
161 :
162 :
163 : /*
164 : This action is defined in core/ as it is used by other actions.
165 : Anyway, it is registered here, so that excluding this module from
166 : compilation will exclude it from plumed.
167 : */
168 :
169 : typedef PLMD::SetupMolInfo MolInfo;
170 :
171 6507 : PLUMED_REGISTER_ACTION(MolInfo,"MOLINFO")
172 :
173 : }
174 4839 : }
|