Line data Source code
1 : /* +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
2 : Copyright (c) 2013-2020 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 "MolDataClass.h"
23 : #include "Exception.h"
24 : #include "Tools.h"
25 : #include "PDB.h"
26 :
27 : namespace PLMD {
28 :
29 25 : unsigned MolDataClass::numberOfAtomsPerResidueInBackbone( const std::string& type ) {
30 25 : if( type=="protein" ) return 5;
31 0 : else if( type=="dna" ) return 6;
32 0 : else if( type=="rna" ) return 6;
33 0 : else return 0;
34 : }
35 :
36 8121 : bool MolDataClass::allowedResidue( const std::string& type, const std::string& residuename ) {
37 8121 : if( type=="protein" ) {
38 6014 : if(residuename=="ALA") return true;
39 5754 : else if(residuename=="ARG") return true;
40 5754 : else if(residuename=="ASN") return true;
41 5754 : else if(residuename=="ASP") return true;
42 5746 : else if(residuename=="CYS") return true;
43 5746 : else if(residuename=="GLN") return true;
44 5746 : else if(residuename=="GLU") return true;
45 5738 : else if(residuename=="GLY") return true;
46 5730 : else if(residuename=="HIS") return true;
47 5730 : else if(residuename=="ILE") return true;
48 5728 : else if(residuename=="LEU") return true;
49 5728 : else if(residuename=="LYS") return true;
50 5724 : else if(residuename=="MET") return true;
51 5724 : else if(residuename=="PHE") return true;
52 5720 : else if(residuename=="PRO") return true;
53 5715 : else if(residuename=="SER") return true;
54 5711 : else if(residuename=="THR") return true;
55 5691 : else if(residuename=="TRP") return true;
56 5681 : else if(residuename=="TYR") return true;
57 5677 : else if(residuename=="VAL") return true;
58 : // Terminal groups
59 489 : else if(residuename=="ACE") return true;
60 483 : else if(residuename=="NME") return true;
61 477 : else if(residuename=="NH2") return true;
62 : // Alternative residue names in common force fiels
63 477 : else if(residuename=="GLH") return true; // neutral GLU
64 477 : else if(residuename=="ASH") return true; // neutral ASP
65 477 : else if(residuename=="HID") return true; // HIS-D amber
66 477 : else if(residuename=="HSD") return true; // HIS-D charmm
67 477 : else if(residuename=="HIE") return true; // HIS-E amber
68 477 : else if(residuename=="HSE") return true; // HIS-E charmm
69 477 : else if(residuename=="HIP") return true; // HIS-P amber
70 477 : else if(residuename=="HSP") return true; // HIS-P charmm
71 : // Weird amino acids
72 477 : else if(residuename=="NLE") return true;
73 475 : else if(residuename=="SFO") return true;
74 475 : else return false;
75 2107 : } else if( type=="dna" ) {
76 462 : if(residuename=="A") return true;
77 397 : else if(residuename=="A5") return true;
78 397 : else if(residuename=="A3") return true;
79 397 : else if(residuename=="AN") return true;
80 397 : else if(residuename=="G") return true;
81 397 : else if(residuename=="G5") return true;
82 397 : else if(residuename=="G3") return true;
83 397 : else if(residuename=="GN") return true;
84 397 : else if(residuename=="T") return true;
85 393 : else if(residuename=="T5") return true;
86 393 : else if(residuename=="T3") return true;
87 393 : else if(residuename=="TN") return true;
88 393 : else if(residuename=="C") return true;
89 393 : else if(residuename=="C5") return true;
90 393 : else if(residuename=="C3") return true;
91 393 : else if(residuename=="CN") return true;
92 393 : else if(residuename=="DA") return true;
93 393 : else if(residuename=="DA5") return true;
94 393 : else if(residuename=="DA3") return true;
95 393 : else if(residuename=="DAN") return true;
96 393 : else if(residuename=="DG") return true;
97 392 : else if(residuename=="DG5") return true;
98 392 : else if(residuename=="DG3") return true;
99 392 : else if(residuename=="DGN") return true;
100 392 : else if(residuename=="DT") return true;
101 392 : else if(residuename=="DT5") return true;
102 392 : else if(residuename=="DT3") return true;
103 392 : else if(residuename=="DTN") return true;
104 392 : else if(residuename=="DC") return true;
105 391 : else if(residuename=="DC5") return true;
106 391 : else if(residuename=="DC3") return true;
107 391 : else if(residuename=="DCN") return true;
108 391 : else return false;
109 1645 : } else if( type=="rna" ) {
110 865 : if(residuename=="A") return true;
111 808 : else if(residuename=="A5") return true;
112 808 : else if(residuename=="A3") return true;
113 808 : else if(residuename=="AN") return true;
114 808 : else if(residuename=="G") return true;
115 786 : else if(residuename=="G5") return true;
116 782 : else if(residuename=="G3") return true;
117 782 : else if(residuename=="GN") return true;
118 782 : else if(residuename=="U") return true;
119 759 : else if(residuename=="U5") return true;
120 759 : else if(residuename=="U3") return true;
121 759 : else if(residuename=="UN") return true;
122 759 : else if(residuename=="C") return true;
123 744 : else if(residuename=="C5") return true;
124 744 : else if(residuename=="C3") return true;
125 740 : else if(residuename=="CN") return true;
126 740 : else if(residuename=="RA") return true;
127 632 : else if(residuename=="RA5") return true;
128 601 : else if(residuename=="RA3") return true;
129 601 : else if(residuename=="RAN") return true;
130 601 : else if(residuename=="RG") return true;
131 481 : else if(residuename=="RG5") return true;
132 481 : else if(residuename=="RG3") return true;
133 442 : else if(residuename=="RGN") return true;
134 442 : else if(residuename=="RU") return true;
135 342 : else if(residuename=="RU5") return true;
136 342 : else if(residuename=="RU3") return true;
137 311 : else if(residuename=="RUN") return true;
138 311 : else if(residuename=="RC") return true;
139 177 : else if(residuename=="RC5") return true;
140 144 : else if(residuename=="RC3") return true;
141 144 : else if(residuename=="RCN") return true;
142 141 : else return false;
143 780 : } else if( type=="water" ) {
144 390 : if(residuename=="SOL") return true;
145 274 : if(residuename=="WAT") return true;
146 274 : return false;
147 390 : } else if( type=="ion" ) {
148 390 : if(residuename=="IB+") return true;
149 390 : if(residuename=="CA") return true;
150 390 : if(residuename=="CL") return true;
151 384 : if(residuename=="NA") return true;
152 372 : if(residuename=="MG") return true;
153 372 : if(residuename=="K") return true;
154 372 : if(residuename=="RB") return true;
155 372 : if(residuename=="CS") return true;
156 372 : if(residuename=="LI") return true;
157 372 : if(residuename=="ZN") return true;
158 372 : return false;
159 : }
160 : return false;
161 : }
162 :
163 2652 : void MolDataClass::getBackboneForResidue( const std::string& type, const unsigned& residuenum, const PDB& mypdb, std::vector<AtomNumber>& atoms ) {
164 2652 : std::string residuename=mypdb.getResidueName( residuenum );
165 2652 : plumed_massert( MolDataClass::allowedResidue( type, residuename ), "residue " + residuename + " unrecognized for molecule type " + type );
166 2652 : if( type=="protein" ) {
167 2652 : if( residuename=="GLY") {
168 0 : atoms.resize(5);
169 0 : atoms[0]=mypdb.getNamedAtomFromResidue("N",residuenum);
170 0 : atoms[1]=mypdb.getNamedAtomFromResidue("CA",residuenum);
171 0 : atoms[2]=mypdb.getNamedAtomFromResidue("HA1",residuenum);
172 0 : atoms[3]=mypdb.getNamedAtomFromResidue("C",residuenum);
173 0 : atoms[4]=mypdb.getNamedAtomFromResidue("O",residuenum);
174 2652 : } else if( residuename=="ACE") {
175 0 : atoms.resize(1);
176 0 : atoms[0]=mypdb.getNamedAtomFromResidue("C",residuenum);
177 5304 : } else if( residuename=="NME"||residuename=="NH2") {
178 0 : atoms.resize(1);
179 0 : atoms[0]=mypdb.getNamedAtomFromResidue("N",residuenum);
180 : } else {
181 2652 : atoms.resize(5);
182 7956 : atoms[0]=mypdb.getNamedAtomFromResidue("N",residuenum);
183 7956 : atoms[1]=mypdb.getNamedAtomFromResidue("CA",residuenum);
184 7956 : atoms[2]=mypdb.getNamedAtomFromResidue("CB",residuenum);
185 7956 : atoms[3]=mypdb.getNamedAtomFromResidue("C",residuenum);
186 7956 : atoms[4]=mypdb.getNamedAtomFromResidue("O",residuenum);
187 : }
188 0 : } else if( type=="dna" || type=="rna" ) {
189 0 : atoms.resize(6);
190 0 : atoms[0]=mypdb.getNamedAtomFromResidue("P",residuenum);
191 0 : atoms[1]=mypdb.getNamedAtomFromResidue("O5\'",residuenum);
192 0 : atoms[2]=mypdb.getNamedAtomFromResidue("C5\'",residuenum);
193 0 : atoms[3]=mypdb.getNamedAtomFromResidue("C4\'",residuenum);
194 0 : atoms[4]=mypdb.getNamedAtomFromResidue("C3\'",residuenum);
195 0 : atoms[5]=mypdb.getNamedAtomFromResidue("O3\'",residuenum);
196 : }
197 : else {
198 0 : plumed_merror(type + " is not a valid molecule type");
199 : }
200 2652 : }
201 :
202 3409 : bool MolDataClass::isTerminalGroup( const std::string& type, const std::string& residuename ) {
203 3409 : if( type=="protein" ) {
204 3409 : if( residuename=="ACE" ) return true;
205 3082 : else if( residuename=="NME" ) return true;
206 2755 : else if( residuename=="NH2" ) return true;
207 2755 : else return false;
208 : } else {
209 0 : plumed_merror(type + " is not a valid molecule type");
210 : }
211 : return false;
212 : }
213 :
214 585 : void MolDataClass::specialSymbol( const std::string& type, const std::string& symbol, const PDB& mypdb, std::vector<AtomNumber>& numbers ) {
215 813 : if(type=="protein" || type=="rna" || type=="dna") {
216 : // symbol should be something like
217 : // phi-123 i.e. phi torsion of residue 123 of first chain
218 : // psi-A321 i.e. psi torsion of residue 321 of chain A
219 : // psi-4_321 is psi torsion of residue 321 of chain 4
220 : // psi-A_321 is equivalent to psi-A321
221 585 : numbers.resize(0);
222 :
223 : // special cases:
224 585 : if(symbol=="protein") {
225 1 : const auto & all = mypdb.getAtomNumbers();
226 133 : for(auto a : all) {
227 132 : auto resname=mypdb.getResidueName(a);
228 132 : Tools::stripLeadingAndTrailingBlanks(resname);
229 264 : if(allowedResidue("protein",resname)) numbers.push_back(a);
230 : }
231 8 : return;
232 : }
233 :
234 584 : if(symbol=="nucleic") {
235 2 : const auto & all = mypdb.getAtomNumbers();
236 457 : for(auto a : all) {
237 455 : auto resname=mypdb.getResidueName(a);
238 455 : Tools::stripLeadingAndTrailingBlanks(resname);
239 1365 : if(allowedResidue("dna",resname) || allowedResidue("rna",resname)) numbers.push_back(a);
240 : }
241 2 : return;
242 : }
243 :
244 582 : if(symbol=="ions") {
245 1 : const auto & all = mypdb.getAtomNumbers();
246 391 : for(auto a : all) {
247 390 : auto resname=mypdb.getResidueName(a);
248 390 : Tools::stripLeadingAndTrailingBlanks(resname);
249 780 : if(allowedResidue("ion",resname)) numbers.push_back(a);
250 : }
251 1 : return;
252 : }
253 :
254 581 : if(symbol=="water") {
255 1 : const auto & all = mypdb.getAtomNumbers();
256 391 : for(auto a : all) {
257 390 : auto resname=mypdb.getResidueName(a);
258 390 : Tools::stripLeadingAndTrailingBlanks(resname);
259 780 : if(allowedResidue("water",resname)) numbers.push_back(a);
260 : }
261 1 : return;
262 : }
263 :
264 580 : if(symbol=="hydrogens") {
265 1 : const auto & all = mypdb.getAtomNumbers();
266 391 : for(auto a : all) {
267 390 : auto atomname=mypdb.getAtomName(a);
268 390 : Tools::stripLeadingAndTrailingBlanks(atomname);
269 : auto notnumber=atomname.find_first_not_of("0123456789");
270 780 : if(notnumber!=std::string::npos && atomname[notnumber]=='H') numbers.push_back(a);
271 : }
272 1 : return;
273 : }
274 :
275 579 : if(symbol=="nonhydrogens") {
276 1 : const auto & all = mypdb.getAtomNumbers();
277 391 : for(auto a : all) {
278 390 : auto atomname=mypdb.getAtomName(a);
279 390 : Tools::stripLeadingAndTrailingBlanks(atomname);
280 : auto notnumber=atomname.find_first_not_of("0123456789");
281 780 : if(!(notnumber!=std::string::npos && atomname[notnumber]=='H')) numbers.push_back(a);
282 : }
283 1 : return;
284 : }
285 :
286 :
287 : std::size_t dash=symbol.find_first_of('-');
288 578 : if(dash==std::string::npos) plumed_error() << "Unrecognized symbol "<<symbol;
289 :
290 578 : std::size_t firstunderscore=symbol.find_first_of('_',dash+1);
291 : std::size_t firstnum=symbol.find_first_of("0123456789",dash+1);
292 578 : std::string name=symbol.substr(0,dash);
293 : unsigned resnum;
294 : std::string resname;
295 : std::string chainid;
296 578 : if(firstunderscore != std::string::npos) {
297 12 : Tools::convert( symbol.substr(firstunderscore+1), resnum );
298 12 : chainid=symbol.substr(dash+1,firstunderscore-(dash+1));
299 572 : } else if(firstnum==dash+1) {
300 1124 : Tools::convert( symbol.substr(dash+1), resnum );
301 : chainid="*"; // this is going to match the first chain
302 : } else {
303 : // if chain id is provided:
304 20 : Tools::convert( symbol.substr(firstnum), resnum );
305 20 : chainid=symbol.substr(dash+1,firstnum-(dash+1)); // this is going to match the proper chain
306 : }
307 1156 : resname= mypdb.getResidueName(resnum,chainid);
308 578 : Tools::stripLeadingAndTrailingBlanks(resname);
309 1156 : if(allowedResidue("protein",resname)) {
310 206 : if( name=="phi" && !isTerminalGroup("protein",resname) ) {
311 111 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C",resnum-1,chainid));
312 111 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("N",resnum,chainid));
313 111 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("CA",resnum,chainid));
314 111 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C",resnum,chainid));
315 132 : } else if( name=="psi" && !isTerminalGroup("protein",resname) ) {
316 189 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("N",resnum,chainid));
317 189 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("CA",resnum,chainid));
318 189 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C",resnum,chainid));
319 189 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("N",resnum+1,chainid));
320 6 : } else if( name=="omega" && !isTerminalGroup("protein",resname) ) {
321 3 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("CA",resnum,chainid));
322 3 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C",resnum,chainid));
323 3 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("N",resnum+1,chainid));
324 3 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("CA",resnum+1,chainid));
325 4 : } else if( name=="chi1" && !isTerminalGroup("protein",resname) ) {
326 3 : if ( resname=="GLY" || resname=="ALA" || resname=="SFO" ) plumed_merror("chi-1 is not defined for ALA, GLY and SFO");
327 4 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("N",resnum,chainid));
328 3 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("CA",resnum,chainid));
329 3 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("CB",resnum,chainid));
330 2 : if(resname=="ILE"||resname=="VAL") {
331 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("CG1",resnum,chainid));
332 1 : } else if(resname=="CYS") {
333 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("SG",resnum,chainid));
334 1 : } else if(resname=="THR") {
335 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("OG1",resnum,chainid));
336 1 : } else if(resname=="SER") {
337 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("OG",resnum,chainid));
338 : } else {
339 3 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("CG",resnum,chainid));
340 : }
341 2 : } else if( name=="chi2" && !isTerminalGroup("protein",resname) ) {
342 5 : if ( resname=="GLY" || resname=="ALA" || resname=="SFO" || resname=="CYS" || resname=="SER" ||
343 2 : resname=="THR" || resname=="VAL" ) plumed_merror("chi-2 is not defined for ALA, GLY, CYS, SER, THR, VAL and SFO");
344 4 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("CA",resnum,chainid));
345 3 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("CB",resnum,chainid));
346 :
347 1 : if(resname=="ILE") {
348 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("CG1",resnum,chainid));
349 : } else {
350 3 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("CG",resnum,chainid));
351 : }
352 2 : if(resname=="ASN" || resname=="ASP") {
353 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("OD1",resnum,chainid));
354 1 : } else if(resname=="HIS") {
355 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("ND1",resnum,chainid));
356 4 : } else if(resname=="LEU" || resname=="PHE" || resname=="TRP" || resname=="TYR") {
357 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("CD1",resnum,chainid));
358 1 : } else if(resname=="MET") {
359 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("SD",resnum,chainid));
360 : } else {
361 3 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("CD",resnum,chainid));
362 : }
363 0 : } else if( name=="chi3" && !isTerminalGroup("protein",resname) ) {
364 0 : if (!( resname=="ARG" || resname=="GLN" || resname=="GLU" || resname=="LYS" ||
365 0 : resname=="MET" )) plumed_merror("chi-3 is defined only for ARG, GLN, GLU, LYS and MET");
366 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("CB",resnum,chainid));
367 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("CG",resnum,chainid));
368 :
369 0 : if(resname=="MET") {
370 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("SD",resnum,chainid));
371 : } else {
372 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("CD",resnum,chainid));
373 : }
374 0 : if(resname=="GLN" || resname=="GLU") {
375 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("OE1",resnum,chainid));
376 0 : } else if(resname=="LYS" || resname=="MET") {
377 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("CE",resnum,chainid));
378 : } else {
379 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("NE",resnum,chainid));
380 : }
381 0 : } else if( name=="chi4" && !isTerminalGroup("protein",resname) ) {
382 0 : if (!( resname=="ARG" || resname=="LYS" )) plumed_merror("chi-4 is defined only for ARG and LYS");
383 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("CG",resnum,chainid));
384 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("CD",resnum,chainid));
385 :
386 0 : if(resname=="ARG") {
387 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("NE",resnum,chainid));
388 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("CZ",resnum,chainid));
389 : } else {
390 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("CE",resnum,chainid));
391 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("NZ",resnum,chainid));
392 : }
393 0 : } else if( name=="chi5" && !isTerminalGroup("protein",resname) ) {
394 0 : if (!( resname=="ARG" )) plumed_merror("chi-5 is defined only for ARG");
395 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("CD",resnum,chainid));
396 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("NE",resnum,chainid));
397 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("CZ",resnum,chainid));
398 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("NH1",resnum,chainid));
399 0 : } else numbers.push_back(mypdb.getNamedAtomFromResidueAndChain(name,resnum,chainid));
400 :
401 1425 : } else if( allowedResidue("rna",resname) || allowedResidue("dna",resname)) {
402 : std::string basetype;
403 474 : if(resname.find_first_of("A")!=std::string::npos) basetype+="A";
404 474 : if(resname.find_first_of("U")!=std::string::npos) basetype+="U";
405 474 : if(resname.find_first_of("T")!=std::string::npos) basetype+="T";
406 474 : if(resname.find_first_of("C")!=std::string::npos) basetype+="C";
407 474 : if(resname.find_first_of("G")!=std::string::npos) basetype+="G";
408 474 : plumed_massert(basetype.length()==1,"cannot find type of rna/dna residue "+resname+" "+basetype);
409 474 : if( name=="chi" ) {
410 102 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("O4\'",resnum,chainid));
411 102 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C1\'",resnum,chainid));
412 99 : if(basetype=="T" || basetype=="U" || basetype=="C") {
413 27 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("N1",resnum,chainid));
414 27 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C2",resnum,chainid));
415 47 : } else if(basetype=="G" || basetype=="A") {
416 75 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("N9",resnum,chainid));
417 75 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C4",resnum,chainid));
418 0 : } else plumed_error();
419 440 : } else if( name=="alpha" ) {
420 21 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("O3\'",resnum-1,chainid));
421 21 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("P",resnum,chainid));
422 21 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("O5\'",resnum,chainid));
423 21 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C5\'",resnum,chainid));
424 433 : } else if( name=="beta" ) {
425 21 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("P",resnum,chainid));
426 21 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("O5\'",resnum,chainid));
427 21 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C5\'",resnum,chainid));
428 21 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C4\'",resnum,chainid));
429 426 : } else if( name=="gamma" ) {
430 81 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("O5\'",resnum,chainid));
431 81 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C5\'",resnum,chainid));
432 81 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C4\'",resnum,chainid));
433 81 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C3\'",resnum,chainid));
434 399 : } else if( name=="delta" ) {
435 3 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C5\'",resnum,chainid));
436 3 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C4\'",resnum,chainid));
437 3 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C3\'",resnum,chainid));
438 3 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("O3\'",resnum,chainid));
439 398 : } else if( name=="epsilon" ) {
440 24 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C4\'",resnum,chainid));
441 24 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C3\'",resnum,chainid));
442 24 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("O3\'",resnum,chainid));
443 24 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("P",resnum+1,chainid));
444 390 : } else if( name=="zeta" ) {
445 21 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C3\'",resnum,chainid));
446 21 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("O3\'",resnum,chainid));
447 21 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("P",resnum+1,chainid));
448 21 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("O5\'",resnum+1,chainid));
449 383 : } else if( name=="v0" ) {
450 9 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C4\'",resnum,chainid));
451 9 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("O4\'",resnum,chainid));
452 9 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C1\'",resnum,chainid));
453 9 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C2\'",resnum,chainid));
454 380 : } else if( name=="v1" ) {
455 9 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("O4\'",resnum,chainid));
456 9 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C1\'",resnum,chainid));
457 9 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C2\'",resnum,chainid));
458 9 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C3\'",resnum,chainid));
459 377 : } else if( name=="v2" ) {
460 9 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C1\'",resnum,chainid));
461 9 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C2\'",resnum,chainid));
462 9 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C3\'",resnum,chainid));
463 9 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C4\'",resnum,chainid));
464 374 : } else if( name=="v3" ) {
465 9 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C2\'",resnum,chainid));
466 9 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C3\'",resnum,chainid));
467 9 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C4\'",resnum,chainid));
468 9 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("O4\'",resnum,chainid));
469 371 : } else if( name=="v4" ) {
470 9 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C3\'",resnum,chainid));
471 9 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C4\'",resnum,chainid));
472 9 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("O4\'",resnum,chainid));
473 9 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C1\'",resnum,chainid));
474 368 : } else if( name=="back" ) {
475 3 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("P",resnum,chainid));
476 3 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("O5\'",resnum,chainid));
477 3 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C5\'",resnum,chainid));
478 3 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C4\'",resnum,chainid));
479 3 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C3\'",resnum,chainid));
480 3 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("O3\'",resnum,chainid));
481 367 : } else if( name=="sugar" ) {
482 156 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C4\'",resnum,chainid));
483 156 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("O4\'",resnum,chainid));
484 156 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C1\'",resnum,chainid));
485 156 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C2\'",resnum,chainid));
486 156 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C3\'",resnum,chainid));
487 315 : } else if( name=="base" ) {
488 25 : if(basetype=="C") {
489 27 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("N1",resnum,chainid));
490 27 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C2",resnum,chainid));
491 27 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("O2",resnum,chainid));
492 27 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("N3",resnum,chainid));
493 27 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C4",resnum,chainid));
494 27 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("N4",resnum,chainid));
495 27 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C5",resnum,chainid));
496 27 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C6",resnum,chainid));
497 16 : } else if(basetype=="U") {
498 6 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("N1",resnum,chainid));
499 6 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C2",resnum,chainid));
500 6 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("O2",resnum,chainid));
501 6 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("N3",resnum,chainid));
502 6 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C4",resnum,chainid));
503 6 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("O4",resnum,chainid));
504 6 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C5",resnum,chainid));
505 6 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C6",resnum,chainid));
506 14 : } else if(basetype=="T") {
507 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("N1",resnum,chainid));
508 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C2",resnum,chainid));
509 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("O2",resnum,chainid));
510 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("N3",resnum,chainid));
511 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C4",resnum,chainid));
512 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("O4",resnum,chainid));
513 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C5",resnum,chainid));
514 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C7",resnum,chainid));
515 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C6",resnum,chainid));
516 14 : } else if(basetype=="G") {
517 9 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("N9",resnum,chainid));
518 9 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C4",resnum,chainid));
519 9 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("N3",resnum,chainid));
520 9 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C2",resnum,chainid));
521 9 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("N2",resnum,chainid));
522 9 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("N1",resnum,chainid));
523 9 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C6",resnum,chainid));
524 9 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("O6",resnum,chainid));
525 9 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C5",resnum,chainid));
526 9 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("N7",resnum,chainid));
527 9 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C8",resnum,chainid));
528 11 : } else if(basetype=="A") {
529 33 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("N9",resnum,chainid));
530 33 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C4",resnum,chainid));
531 33 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("N1",resnum,chainid));
532 33 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C2",resnum,chainid));
533 33 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("N3",resnum,chainid));
534 33 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C6",resnum,chainid));
535 33 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("N6",resnum,chainid));
536 33 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C5",resnum,chainid));
537 33 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("N7",resnum,chainid));
538 33 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C8",resnum,chainid));
539 0 : } else plumed_error();
540 290 : } else if( name=="lcs" ) {
541 852 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C2",resnum,chainid));
542 752 : if(basetype=="T" || basetype=="U" || basetype=="C") {
543 444 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C4",resnum,chainid));
544 444 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C6",resnum,chainid));
545 216 : } else if(basetype=="G" || basetype=="A") {
546 408 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C6",resnum,chainid));
547 408 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C4",resnum,chainid));
548 0 : } else plumed_error();
549 12 : } else numbers.push_back(mypdb.getNamedAtomFromResidueAndChain(name,resnum,chainid));
550 2 : } else numbers.push_back(mypdb.getNamedAtomFromResidueAndChain(name,resnum,chainid));
551 : }
552 : else {
553 0 : plumed_merror(type + " is not a valid molecule type");
554 : }
555 : }
556 :
557 : }
|
