Line data Source code
1 : /* +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
2 : Copyright (c) 2013-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 "MolDataClass.h"
23 : #include "Exception.h"
24 : #include "Tools.h"
25 : #include "PDB.h"
26 :
27 : namespace PLMD {
28 :
29 38 : unsigned MolDataClass::numberOfAtomsPerResidueInBackbone( const std::string& type ) {
30 38 : if( type=="protein" ) {
31 : return 5;
32 0 : } else if( type=="dna" ) {
33 : return 6;
34 0 : } else if( type=="rna" ) {
35 : return 6;
36 : } else {
37 0 : return 0;
38 : }
39 : }
40 :
41 11809 : bool MolDataClass::allowedResidue( const std::string& type, const std::string& residuename ) {
42 11809 : if( type=="protein" ) {
43 9696 : if(residuename=="ALA") {
44 : return true;
45 9308 : } else if(residuename=="ARG") {
46 : return true;
47 9308 : } else if(residuename=="ASN") {
48 : return true;
49 9308 : } else if(residuename=="ASP") {
50 : return true;
51 9296 : } else if(residuename=="CYS") {
52 : return true;
53 9296 : } else if(residuename=="GLN") {
54 : return true;
55 9296 : } else if(residuename=="GLU") {
56 : return true;
57 9284 : } else if(residuename=="GLY") {
58 : return true;
59 9256 : } else if(residuename=="HIS") {
60 : return true;
61 9256 : } else if(residuename=="ILE") {
62 : return true;
63 9252 : } else if(residuename=="LEU") {
64 : return true;
65 9252 : } else if(residuename=="LYS") {
66 : return true;
67 9246 : } else if(residuename=="MET") {
68 : return true;
69 9246 : } else if(residuename=="PHE") {
70 : return true;
71 9240 : } else if(residuename=="PRO") {
72 : return true;
73 9233 : } else if(residuename=="SER") {
74 : return true;
75 9211 : } else if(residuename=="THR") {
76 : return true;
77 9181 : } else if(residuename=="TRP") {
78 : return true;
79 9149 : } else if(residuename=="TYR") {
80 : return true;
81 9143 : } else if(residuename=="VAL") {
82 : return true;
83 : }
84 : // Terminal groups
85 497 : else if(residuename=="ACE") {
86 : return true;
87 491 : } else if(residuename=="NME") {
88 : return true;
89 485 : } else if(residuename=="NH2") {
90 : return true;
91 : }
92 : // Alternative residue names in common force fields
93 485 : else if(residuename=="GLH") {
94 : return true; // neutral GLU
95 485 : } else if(residuename=="ASH") {
96 : return true; // neutral ASP
97 485 : } else if(residuename=="HID") {
98 : return true; // HIS-D amber
99 485 : } else if(residuename=="HSD") {
100 : return true; // HIS-D charmm
101 485 : } else if(residuename=="HIE") {
102 : return true; // HIS-E amber
103 485 : } else if(residuename=="HSE") {
104 : return true; // HIS-E charmm
105 485 : } else if(residuename=="HIP") {
106 : return true; // HIS-P amber
107 485 : } else if(residuename=="HSP") {
108 : return true; // HIS-P charmm
109 : }
110 : // Weird amino acids
111 485 : else if(residuename=="NLE") {
112 : return true;
113 481 : } else if(residuename=="SFO") {
114 : return true;
115 : } else {
116 481 : return false;
117 : }
118 2113 : } else if( type=="dna" ) {
119 462 : if(residuename=="A") {
120 : return true;
121 397 : } else if(residuename=="A5") {
122 : return true;
123 397 : } else if(residuename=="A3") {
124 : return true;
125 397 : } else if(residuename=="AN") {
126 : return true;
127 397 : } else if(residuename=="G") {
128 : return true;
129 397 : } else if(residuename=="G5") {
130 : return true;
131 397 : } else if(residuename=="G3") {
132 : return true;
133 397 : } else if(residuename=="GN") {
134 : return true;
135 397 : } else if(residuename=="T") {
136 : return true;
137 393 : } else if(residuename=="T5") {
138 : return true;
139 393 : } else if(residuename=="T3") {
140 : return true;
141 393 : } else if(residuename=="TN") {
142 : return true;
143 393 : } else if(residuename=="C") {
144 : return true;
145 393 : } else if(residuename=="C5") {
146 : return true;
147 393 : } else if(residuename=="C3") {
148 : return true;
149 393 : } else if(residuename=="CN") {
150 : return true;
151 393 : } else if(residuename=="DA") {
152 : return true;
153 393 : } else if(residuename=="DA5") {
154 : return true;
155 393 : } else if(residuename=="DA3") {
156 : return true;
157 393 : } else if(residuename=="DAN") {
158 : return true;
159 393 : } else if(residuename=="DG") {
160 : return true;
161 392 : } else if(residuename=="DG5") {
162 : return true;
163 392 : } else if(residuename=="DG3") {
164 : return true;
165 392 : } else if(residuename=="DGN") {
166 : return true;
167 392 : } else if(residuename=="DT") {
168 : return true;
169 392 : } else if(residuename=="DT5") {
170 : return true;
171 392 : } else if(residuename=="DT3") {
172 : return true;
173 392 : } else if(residuename=="DTN") {
174 : return true;
175 392 : } else if(residuename=="DC") {
176 : return true;
177 391 : } else if(residuename=="DC5") {
178 : return true;
179 391 : } else if(residuename=="DC3") {
180 : return true;
181 391 : } else if(residuename=="DCN") {
182 : return true;
183 : } else {
184 391 : return false;
185 : }
186 1651 : } else if( type=="rna" ) {
187 871 : if(residuename=="A") {
188 : return true;
189 808 : } else if(residuename=="A5") {
190 : return true;
191 808 : } else if(residuename=="A3") {
192 : return true;
193 808 : } else if(residuename=="AN") {
194 : return true;
195 808 : } else if(residuename=="G") {
196 : return true;
197 786 : } else if(residuename=="G5") {
198 : return true;
199 782 : } else if(residuename=="G3") {
200 : return true;
201 782 : } else if(residuename=="GN") {
202 : return true;
203 782 : } else if(residuename=="U") {
204 : return true;
205 759 : } else if(residuename=="U5") {
206 : return true;
207 759 : } else if(residuename=="U3") {
208 : return true;
209 759 : } else if(residuename=="UN") {
210 : return true;
211 759 : } else if(residuename=="C") {
212 : return true;
213 744 : } else if(residuename=="C5") {
214 : return true;
215 744 : } else if(residuename=="C3") {
216 : return true;
217 740 : } else if(residuename=="CN") {
218 : return true;
219 740 : } else if(residuename=="RA") {
220 : return true;
221 632 : } else if(residuename=="RA5") {
222 : return true;
223 601 : } else if(residuename=="RA3") {
224 : return true;
225 601 : } else if(residuename=="RAN") {
226 : return true;
227 601 : } else if(residuename=="RG") {
228 : return true;
229 481 : } else if(residuename=="RG5") {
230 : return true;
231 481 : } else if(residuename=="RG3") {
232 : return true;
233 442 : } else if(residuename=="RGN") {
234 : return true;
235 442 : } else if(residuename=="RU") {
236 : return true;
237 342 : } else if(residuename=="RU5") {
238 : return true;
239 342 : } else if(residuename=="RU3") {
240 : return true;
241 311 : } else if(residuename=="RUN") {
242 : return true;
243 311 : } else if(residuename=="RC") {
244 : return true;
245 177 : } else if(residuename=="RC5") {
246 : return true;
247 144 : } else if(residuename=="RC3") {
248 : return true;
249 144 : } else if(residuename=="RCN") {
250 : return true;
251 : } else {
252 141 : return false;
253 : }
254 780 : } else if( type=="water" ) {
255 390 : if(residuename=="SOL") {
256 : return true;
257 : }
258 274 : if(residuename=="WAT") {
259 : return true;
260 : }
261 274 : return false;
262 390 : } else if( type=="ion" ) {
263 390 : if(residuename=="IB+") {
264 : return true;
265 : }
266 390 : if(residuename=="CA") {
267 : return true;
268 : }
269 390 : if(residuename=="CL") {
270 : return true;
271 : }
272 384 : if(residuename=="NA") {
273 : return true;
274 : }
275 372 : if(residuename=="MG") {
276 : return true;
277 : }
278 372 : if(residuename=="K") {
279 : return true;
280 : }
281 372 : if(residuename=="RB") {
282 : return true;
283 : }
284 372 : if(residuename=="CS") {
285 : return true;
286 : }
287 372 : if(residuename=="LI") {
288 : return true;
289 : }
290 372 : if(residuename=="ZN") {
291 : return true;
292 : }
293 372 : return false;
294 : }
295 : return false;
296 : }
297 :
298 4404 : void MolDataClass::getBackboneForResidue( const std::string& type, const unsigned& residuenum, const PDB& mypdb, std::vector<AtomNumber>& atoms ) {
299 4404 : std::string residuename=mypdb.getResidueName( residuenum );
300 4404 : plumed_massert( MolDataClass::allowedResidue( type, residuename ), "residue " + residuename + " unrecognized for molecule type " + type );
301 4404 : if( type=="protein" ) {
302 4404 : if( residuename=="GLY") {
303 0 : atoms.resize(5);
304 0 : atoms[0]=mypdb.getNamedAtomFromResidue("N",residuenum);
305 0 : atoms[1]=mypdb.getNamedAtomFromResidue("CA",residuenum);
306 0 : atoms[2]=mypdb.getNamedAtomFromResidue("HA2",residuenum);
307 0 : atoms[3]=mypdb.getNamedAtomFromResidue("C",residuenum);
308 0 : atoms[4]=mypdb.getNamedAtomFromResidue("O",residuenum);
309 4404 : } else if( residuename=="ACE") {
310 0 : atoms.resize(1);
311 0 : atoms[0]=mypdb.getNamedAtomFromResidue("C",residuenum);
312 8808 : } else if( residuename=="NME"||residuename=="NH2") {
313 0 : atoms.resize(1);
314 0 : atoms[0]=mypdb.getNamedAtomFromResidue("N",residuenum);
315 : } else {
316 4404 : atoms.resize(5);
317 4404 : atoms[0]=mypdb.getNamedAtomFromResidue("N",residuenum);
318 4404 : atoms[1]=mypdb.getNamedAtomFromResidue("CA",residuenum);
319 4404 : atoms[2]=mypdb.getNamedAtomFromResidue("CB",residuenum);
320 4404 : atoms[3]=mypdb.getNamedAtomFromResidue("C",residuenum);
321 4404 : atoms[4]=mypdb.getNamedAtomFromResidue("O",residuenum);
322 : }
323 0 : } else if( type=="dna" || type=="rna" ) {
324 0 : atoms.resize(6);
325 0 : atoms[0]=mypdb.getNamedAtomFromResidue("P",residuenum);
326 0 : atoms[1]=mypdb.getNamedAtomFromResidue("O5\'",residuenum);
327 0 : atoms[2]=mypdb.getNamedAtomFromResidue("C5\'",residuenum);
328 0 : atoms[3]=mypdb.getNamedAtomFromResidue("C4\'",residuenum);
329 0 : atoms[4]=mypdb.getNamedAtomFromResidue("C3\'",residuenum);
330 0 : atoms[5]=mypdb.getNamedAtomFromResidue("O3\'",residuenum);
331 : } else {
332 0 : plumed_merror(type + " is not a valid molecule type");
333 : }
334 4404 : }
335 :
336 5771 : bool MolDataClass::isTerminalGroup( const std::string& type, const std::string& residuename ) {
337 5771 : if( type=="protein" ) {
338 5771 : if( residuename=="ACE" ) {
339 : return true;
340 5225 : } else if( residuename=="NME" ) {
341 : return true;
342 4679 : } else if( residuename=="NH2" ) {
343 : return true;
344 : } else {
345 4679 : return false;
346 : }
347 : } else {
348 0 : plumed_merror(type + " is not a valid molecule type");
349 : }
350 : return false;
351 : }
352 :
353 763 : void MolDataClass::specialSymbol( const std::string& type, const std::string& symbol, const PDB& mypdb, std::vector<AtomNumber>& numbers ) {
354 999 : if(type=="protein" || type=="rna" || type=="dna") {
355 : // symbol should be something like
356 : // phi-123 i.e. phi torsion of residue 123 of first chain
357 : // psi-A321 i.e. psi torsion of residue 321 of chain A
358 : // psi-4_321 is psi torsion of residue 321 of chain 4
359 : // psi-A_321 is equivalent to psi-A321
360 763 : numbers.resize(0);
361 :
362 : // special cases:
363 763 : if(symbol=="protein") {
364 1 : const auto & all = mypdb.getAtomNumbers();
365 133 : for(auto a : all) {
366 132 : auto resname=mypdb.getResidueName(a);
367 132 : Tools::stripLeadingAndTrailingBlanks(resname);
368 264 : if(allowedResidue("protein",resname)) {
369 132 : numbers.push_back(a);
370 : }
371 : }
372 7 : return;
373 : }
374 :
375 762 : if(symbol=="nucleic") {
376 2 : const auto & all = mypdb.getAtomNumbers();
377 457 : for(auto a : all) {
378 455 : auto resname=mypdb.getResidueName(a);
379 455 : Tools::stripLeadingAndTrailingBlanks(resname);
380 1101 : if(allowedResidue("dna",resname) || allowedResidue("rna",resname)) {
381 321 : numbers.push_back(a);
382 : }
383 : }
384 2 : return;
385 : }
386 :
387 760 : if(symbol=="ions") {
388 1 : const auto & all = mypdb.getAtomNumbers();
389 391 : for(auto a : all) {
390 390 : auto resname=mypdb.getResidueName(a);
391 390 : Tools::stripLeadingAndTrailingBlanks(resname);
392 780 : if(allowedResidue("ion",resname)) {
393 18 : numbers.push_back(a);
394 : }
395 : }
396 1 : return;
397 : }
398 :
399 759 : if(symbol=="water") {
400 1 : const auto & all = mypdb.getAtomNumbers();
401 391 : for(auto a : all) {
402 390 : auto resname=mypdb.getResidueName(a);
403 390 : Tools::stripLeadingAndTrailingBlanks(resname);
404 780 : if(allowedResidue("water",resname)) {
405 116 : numbers.push_back(a);
406 : }
407 : }
408 1 : return;
409 : }
410 :
411 758 : if(symbol=="hydrogens") {
412 1 : const auto & all = mypdb.getAtomNumbers();
413 391 : for(auto a : all) {
414 390 : auto atomname=mypdb.getAtomName(a);
415 390 : Tools::stripLeadingAndTrailingBlanks(atomname);
416 390 : auto notnumber=atomname.find_first_not_of("0123456789");
417 390 : if(notnumber!=std::string::npos && atomname[notnumber]=='H') {
418 149 : numbers.push_back(a);
419 : }
420 : }
421 1 : return;
422 : }
423 :
424 757 : if(symbol=="nonhydrogens") {
425 1 : const auto & all = mypdb.getAtomNumbers();
426 391 : for(auto a : all) {
427 390 : auto atomname=mypdb.getAtomName(a);
428 390 : Tools::stripLeadingAndTrailingBlanks(atomname);
429 390 : auto notnumber=atomname.find_first_not_of("0123456789");
430 390 : if(!(notnumber!=std::string::npos && atomname[notnumber]=='H')) {
431 241 : numbers.push_back(a);
432 : }
433 : }
434 1 : return;
435 : }
436 :
437 :
438 : std::size_t dash=symbol.find_first_of('-');
439 756 : if(dash==std::string::npos) {
440 0 : plumed_error() << "Unrecognized symbol "<<symbol;
441 : }
442 :
443 756 : std::size_t firstunderscore=symbol.find_first_of('_',dash+1);
444 756 : std::size_t firstnum=symbol.find_first_of("0123456789",dash+1);
445 756 : std::string name=symbol.substr(0,dash);
446 : unsigned resnum;
447 : std::string resname;
448 : std::string chainid;
449 756 : if(firstunderscore != std::string::npos) {
450 6 : Tools::convert( symbol.substr(firstunderscore+1), resnum );
451 12 : chainid=symbol.substr(dash+1,firstunderscore-(dash+1));
452 750 : } else if(firstnum==dash+1) {
453 1480 : Tools::convert( symbol.substr(dash+1), resnum );
454 : chainid="*"; // this is going to match the first chain
455 : } else {
456 : // if chain id is provided:
457 10 : Tools::convert( symbol.substr(firstnum), resnum );
458 20 : chainid=symbol.substr(dash+1,firstnum-(dash+1)); // this is going to match the proper chain
459 : }
460 756 : resname= mypdb.getResidueName(resnum,chainid);
461 756 : Tools::stripLeadingAndTrailingBlanks(resname);
462 1512 : if(allowedResidue("protein",resname)) {
463 334 : if( name=="phi" && !isTerminalGroup("protein",resname) ) {
464 118 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C",resnum-1,chainid));
465 118 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("N",resnum,chainid));
466 118 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("CA",resnum,chainid));
467 118 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C",resnum,chainid));
468 427 : } else if( name=="psi" && !isTerminalGroup("protein",resname) ) {
469 422 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("N",resnum,chainid));
470 422 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("CA",resnum,chainid));
471 422 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C",resnum,chainid));
472 422 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("N",resnum+1,chainid));
473 6 : } else if( name=="omega" && !isTerminalGroup("protein",resname) ) {
474 2 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("CA",resnum-1,chainid));
475 2 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C",resnum-1,chainid));
476 2 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("N",resnum,chainid));
477 2 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("CA",resnum,chainid));
478 5 : } else if( name=="chi1" && !isTerminalGroup("protein",resname) ) {
479 3 : if ( resname=="GLY" || resname=="ALA" || resname=="SFO" ) {
480 0 : plumed_merror("chi-1 is not defined for ALA, GLY and SFO");
481 : }
482 2 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("N",resnum,chainid));
483 2 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("CA",resnum,chainid));
484 2 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("CB",resnum,chainid));
485 2 : if(resname=="ILE"||resname=="VAL") {
486 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("CG1",resnum,chainid));
487 1 : } else if(resname=="CYS") {
488 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("SG",resnum,chainid));
489 1 : } else if(resname=="THR") {
490 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("OG1",resnum,chainid));
491 1 : } else if(resname=="SER") {
492 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("OG",resnum,chainid));
493 : } else {
494 2 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("CG",resnum,chainid));
495 : }
496 4 : } else if( name=="chi2" && !isTerminalGroup("protein",resname) ) {
497 5 : if ( resname=="GLY" || resname=="ALA" || resname=="SFO" || resname=="CYS" || resname=="SER" ||
498 2 : resname=="THR" || resname=="VAL" ) {
499 0 : plumed_merror("chi-2 is not defined for ALA, GLY, CYS, SER, THR, VAL and SFO");
500 : }
501 2 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("CA",resnum,chainid));
502 2 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("CB",resnum,chainid));
503 :
504 1 : if(resname=="ILE") {
505 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("CG1",resnum,chainid));
506 : } else {
507 2 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("CG",resnum,chainid));
508 : }
509 2 : if(resname=="ASN" || resname=="ASP") {
510 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("OD1",resnum,chainid));
511 1 : } else if(resname=="HIS") {
512 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("ND1",resnum,chainid));
513 4 : } else if(resname=="LEU" || resname=="PHE" || resname=="TRP" || resname=="TYR") {
514 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("CD1",resnum,chainid));
515 1 : } else if(resname=="MET") {
516 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("SD",resnum,chainid));
517 : } else {
518 2 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("CD",resnum,chainid));
519 : }
520 2 : } else if( name=="chi3" && !isTerminalGroup("protein",resname) ) {
521 0 : if (!( resname=="ARG" || resname=="GLN" || resname=="GLU" || resname=="LYS" ||
522 : resname=="MET" )) {
523 0 : plumed_merror("chi-3 is defined only for ARG, GLN, GLU, LYS and MET");
524 : }
525 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("CB",resnum,chainid));
526 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("CG",resnum,chainid));
527 :
528 0 : if(resname=="MET") {
529 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("SD",resnum,chainid));
530 : } else {
531 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("CD",resnum,chainid));
532 : }
533 0 : if(resname=="GLN" || resname=="GLU") {
534 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("OE1",resnum,chainid));
535 0 : } else if(resname=="LYS" || resname=="MET") {
536 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("CE",resnum,chainid));
537 : } else {
538 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("NE",resnum,chainid));
539 : }
540 2 : } else if( name=="chi4" && !isTerminalGroup("protein",resname) ) {
541 0 : if (!( resname=="ARG" || resname=="LYS" )) {
542 0 : plumed_merror("chi-4 is defined only for ARG and LYS");
543 : }
544 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("CG",resnum,chainid));
545 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("CD",resnum,chainid));
546 :
547 0 : if(resname=="ARG") {
548 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("NE",resnum,chainid));
549 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("CZ",resnum,chainid));
550 : } else {
551 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("CE",resnum,chainid));
552 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("NZ",resnum,chainid));
553 : }
554 2 : } else if( name=="chi5" && !isTerminalGroup("protein",resname) ) {
555 0 : if (!( resname=="ARG" )) {
556 0 : plumed_merror("chi-5 is defined only for ARG");
557 : }
558 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("CD",resnum,chainid));
559 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("NE",resnum,chainid));
560 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("CZ",resnum,chainid));
561 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("NH1",resnum,chainid));
562 3 : } else if( name=="sidechain" && !isTerminalGroup("protein",resname) ) {
563 1 : if(resname=="GLY") {
564 0 : plumed_merror("sidechain is not defined for GLY");
565 : }
566 1 : std::vector<AtomNumber> tmpnum1(mypdb.getAtomsInResidue(resnum,chainid));
567 15 : for(unsigned i=0; i<tmpnum1.size(); i++) {
568 66 : if(mypdb.getAtomName(tmpnum1[i])!="N"&&mypdb.getAtomName(tmpnum1[i])!="CA"&&
569 70 : mypdb.getAtomName(tmpnum1[i])!="C"&&mypdb.getAtomName(tmpnum1[i])!="O"&&
570 61 : mypdb.getAtomName(tmpnum1[i])!="HA"&&mypdb.getAtomName(tmpnum1[i])!="H"&&
571 59 : mypdb.getAtomName(tmpnum1[i])!="HN"&&mypdb.getAtomName(tmpnum1[i])!="H1"&&
572 59 : mypdb.getAtomName(tmpnum1[i])!="H2"&&mypdb.getAtomName(tmpnum1[i])!="H3"&&
573 59 : mypdb.getAtomName(tmpnum1[i])!="OXT"&&mypdb.getAtomName(tmpnum1[i])!="OT1"&&
574 73 : mypdb.getAtomName(tmpnum1[i])!="OT2"&&mypdb.getAtomName(tmpnum1[i])!="OC1"&&
575 32 : mypdb.getAtomName(tmpnum1[i])!="OC2") {
576 9 : numbers.push_back( tmpnum1[i] );
577 : }
578 : }
579 2 : } else if( name=="back" && !isTerminalGroup("protein",resname) ) {
580 1 : std::vector<AtomNumber> tmpnum1(mypdb.getAtomsInResidue(resnum,chainid));
581 15 : for(unsigned i=0; i<tmpnum1.size(); i++) {
582 66 : if(mypdb.getAtomName(tmpnum1[i])=="N"||mypdb.getAtomName(tmpnum1[i])=="CA"||
583 70 : mypdb.getAtomName(tmpnum1[i])=="C"||mypdb.getAtomName(tmpnum1[i])=="O"||
584 61 : mypdb.getAtomName(tmpnum1[i])=="HA"||mypdb.getAtomName(tmpnum1[i])=="H"||
585 59 : mypdb.getAtomName(tmpnum1[i])=="HN"||mypdb.getAtomName(tmpnum1[i])=="H1"||
586 59 : mypdb.getAtomName(tmpnum1[i])=="H2"||mypdb.getAtomName(tmpnum1[i])=="H3"||
587 59 : mypdb.getAtomName(tmpnum1[i])=="OXT"||mypdb.getAtomName(tmpnum1[i])=="OT1"||
588 59 : mypdb.getAtomName(tmpnum1[i])=="OT2"||mypdb.getAtomName(tmpnum1[i])=="OC1"||
589 59 : mypdb.getAtomName(tmpnum1[i])=="OC2"||mypdb.getAtomName(tmpnum1[i])=="HA1"||
590 63 : mypdb.getAtomName(tmpnum1[i])=="HA2"||mypdb.getAtomName(tmpnum1[i])=="HA3") {
591 5 : numbers.push_back( tmpnum1[i] );
592 : }
593 : }
594 : } else {
595 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain(name,resnum,chainid));
596 : }
597 :
598 494 : } else if( allowedResidue("rna",resname) || allowedResidue("dna",resname)) {
599 : std::string basetype;
600 480 : if(resname.find_first_of("A")!=std::string::npos) {
601 : basetype+="A";
602 : }
603 480 : if(resname.find_first_of("U")!=std::string::npos) {
604 : basetype+="U";
605 : }
606 480 : if(resname.find_first_of("T")!=std::string::npos) {
607 : basetype+="T";
608 : }
609 480 : if(resname.find_first_of("C")!=std::string::npos) {
610 : basetype+="C";
611 : }
612 480 : if(resname.find_first_of("G")!=std::string::npos) {
613 : basetype+="G";
614 : }
615 480 : plumed_massert(basetype.length()==1,"cannot find type of rna/dna residue "+resname+" "+basetype);
616 480 : if( name=="chi" ) {
617 72 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("O4\'",resnum,chainid));
618 72 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C1\'",resnum,chainid));
619 105 : if(basetype=="T" || basetype=="U" || basetype=="C") {
620 18 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("N1",resnum,chainid));
621 18 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C2",resnum,chainid));
622 51 : } else if(basetype=="G" || basetype=="A") {
623 54 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("N9",resnum,chainid));
624 54 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C4",resnum,chainid));
625 : } else {
626 0 : plumed_error();
627 : }
628 444 : } else if( name=="alpha" ) {
629 14 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("O3\'",resnum-1,chainid));
630 14 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("P",resnum,chainid));
631 14 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("O5\'",resnum,chainid));
632 14 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C5\'",resnum,chainid));
633 437 : } else if( name=="beta" ) {
634 14 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("P",resnum,chainid));
635 14 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("O5\'",resnum,chainid));
636 14 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C5\'",resnum,chainid));
637 14 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C4\'",resnum,chainid));
638 430 : } else if( name=="gamma" ) {
639 58 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("O5\'",resnum,chainid));
640 58 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C5\'",resnum,chainid));
641 58 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C4\'",resnum,chainid));
642 58 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C3\'",resnum,chainid));
643 401 : } else if( name=="delta" ) {
644 2 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C5\'",resnum,chainid));
645 2 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C4\'",resnum,chainid));
646 2 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C3\'",resnum,chainid));
647 2 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("O3\'",resnum,chainid));
648 400 : } else if( name=="epsilon" ) {
649 16 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C4\'",resnum,chainid));
650 16 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C3\'",resnum,chainid));
651 16 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("O3\'",resnum,chainid));
652 16 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("P",resnum+1,chainid));
653 392 : } else if( name=="zeta" ) {
654 14 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C3\'",resnum,chainid));
655 14 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("O3\'",resnum,chainid));
656 14 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("P",resnum+1,chainid));
657 14 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("O5\'",resnum+1,chainid));
658 385 : } else if( name=="v0" ) {
659 6 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C4\'",resnum,chainid));
660 6 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("O4\'",resnum,chainid));
661 6 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C1\'",resnum,chainid));
662 6 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C2\'",resnum,chainid));
663 382 : } else if( name=="v1" ) {
664 6 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("O4\'",resnum,chainid));
665 6 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C1\'",resnum,chainid));
666 6 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C2\'",resnum,chainid));
667 6 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C3\'",resnum,chainid));
668 379 : } else if( name=="v2" ) {
669 6 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C1\'",resnum,chainid));
670 6 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C2\'",resnum,chainid));
671 6 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C3\'",resnum,chainid));
672 9 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C4\'",resnum,chainid));
673 376 : } else if( name=="v3" ) {
674 6 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C2\'",resnum,chainid));
675 6 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C3\'",resnum,chainid));
676 6 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C4\'",resnum,chainid));
677 6 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("O4\'",resnum,chainid));
678 373 : } else if( name=="v4" ) {
679 6 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C3\'",resnum,chainid));
680 6 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C4\'",resnum,chainid));
681 6 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("O4\'",resnum,chainid));
682 6 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C1\'",resnum,chainid));
683 370 : } else if( name=="back" ) {
684 2 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("P",resnum,chainid));
685 2 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("O5\'",resnum,chainid));
686 2 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C5\'",resnum,chainid));
687 2 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C4\'",resnum,chainid));
688 2 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C3\'",resnum,chainid));
689 2 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("O3\'",resnum,chainid));
690 369 : } else if( name=="sugar" ) {
691 108 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C4\'",resnum,chainid));
692 108 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("O4\'",resnum,chainid));
693 108 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C1\'",resnum,chainid));
694 108 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C2\'",resnum,chainid));
695 108 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C3\'",resnum,chainid));
696 315 : } else if( name=="base" ) {
697 25 : if(basetype=="C") {
698 18 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("N1",resnum,chainid));
699 18 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C2",resnum,chainid));
700 18 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("O2",resnum,chainid));
701 18 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("N3",resnum,chainid));
702 18 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C4",resnum,chainid));
703 18 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("N4",resnum,chainid));
704 18 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C5",resnum,chainid));
705 18 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C6",resnum,chainid));
706 16 : } else if(basetype=="U") {
707 4 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("N1",resnum,chainid));
708 4 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C2",resnum,chainid));
709 4 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("O2",resnum,chainid));
710 4 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("N3",resnum,chainid));
711 4 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C4",resnum,chainid));
712 4 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("O4",resnum,chainid));
713 4 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C5",resnum,chainid));
714 4 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C6",resnum,chainid));
715 14 : } else if(basetype=="T") {
716 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("N1",resnum,chainid));
717 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C2",resnum,chainid));
718 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("O2",resnum,chainid));
719 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("N3",resnum,chainid));
720 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C4",resnum,chainid));
721 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("O4",resnum,chainid));
722 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C5",resnum,chainid));
723 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C7",resnum,chainid));
724 0 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C6",resnum,chainid));
725 14 : } else if(basetype=="G") {
726 6 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("N9",resnum,chainid));
727 6 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C4",resnum,chainid));
728 6 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("N3",resnum,chainid));
729 6 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C2",resnum,chainid));
730 6 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("N2",resnum,chainid));
731 6 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("N1",resnum,chainid));
732 6 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C6",resnum,chainid));
733 6 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("O6",resnum,chainid));
734 6 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C5",resnum,chainid));
735 6 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("N7",resnum,chainid));
736 6 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C8",resnum,chainid));
737 11 : } else if(basetype=="A") {
738 22 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("N9",resnum,chainid));
739 22 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C4",resnum,chainid));
740 22 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("N1",resnum,chainid));
741 22 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C2",resnum,chainid));
742 22 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("N3",resnum,chainid));
743 22 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C6",resnum,chainid));
744 22 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("N6",resnum,chainid));
745 22 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C5",resnum,chainid));
746 22 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("N7",resnum,chainid));
747 22 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C8",resnum,chainid));
748 : } else {
749 0 : plumed_error();
750 : }
751 290 : } else if( name=="lcs" ) {
752 568 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C2",resnum,chainid));
753 752 : if(basetype=="T" || basetype=="U" || basetype=="C") {
754 296 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C4",resnum,chainid));
755 296 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C6",resnum,chainid));
756 216 : } else if(basetype=="G" || basetype=="A") {
757 272 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C6",resnum,chainid));
758 272 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain("C4",resnum,chainid));
759 : } else {
760 0 : plumed_error();
761 : }
762 : } else {
763 12 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain(name,resnum,chainid));
764 : }
765 : } else {
766 2 : numbers.push_back(mypdb.getNamedAtomFromResidueAndChain(name,resnum,chainid));
767 : }
768 : } else {
769 0 : plumed_merror(type + " is not a valid molecule type");
770 : }
771 : }
772 :
773 : }
|
