ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/freemol/branches/sync4pymol12/src/mengine/src/kangle.c
Revision: 93
Committed: Sat Jan 17 23:24:55 2009 UTC (13 years, 4 months ago) by wdelano
File size: 24089 byte(s)
Log Message:
branch created
Line File contents
1 #define EXTERN extern
2 #include "pcwin.h"
3 #include "pcmod.h"
4
5 #include "angles.h"
6 #include "rings.h"
7 #include "field.h"
8 #include "bonds_ff.h"
9 #include "atom_k.h"
10
11 int is_ring43(int, int, int);
12 int is_ring53(int, int, int);
13 int is_ring63(int, int, int);
14 int is_delocalbond(int, int);
15 int is_TS(int, int, int);
16 int isbond(int,int);
17 int find_bond(int,int);
18 void ts_four(int *, float *, float *,char *, char *,char *,char *);
19 long ipack(int ,int ,int ,int);
20 long kang(int,int,int);
21 void numeral(int,char *,int);
22 float get_general_angle(int);
23
24
25 EXTERN struct t_angk1 {
26 int use_ang3, use_ang4, use_ang5;
27 int nang, nang3, nang4, nang5;
28 int ndel, ndel3, ndel4;
29 char ktype[MAXANGCONST][10],ktype3[MAXANG3CONST][10],ktype4[MAXANG4CONST][10],ktype5[MAXANG5CONST][10];
30 char kdel[MAXANGDEL][10],kdel3[MAXANG3DEL][10],kdel4[MAXANG4DEL][10];
31 float con[MAXANGCONST], ang[MAXANGCONST][3];
32 float con3[MAXANG3CONST], ang3[MAXANG3CONST][3];
33 float con4[MAXANG4CONST], ang4[MAXANG4CONST][3];
34 float con5[MAXANG5CONST], ang5[MAXANG5CONST][3];
35 float condel[MAXANGDEL], angdel[MAXANGDEL][3];
36 float condel3[MAXANG3DEL], angdel3[MAXANG3DEL][3];
37 float condel4[MAXANG4DEL], angdel4[MAXANG4DEL][3];
38 } angk1;
39
40 EXTERN struct t_angkp1 {
41 int npiang;
42 char kpa[MAXPIANGCONST][10];
43 float pacon[MAXPIANGCONST], panat[MAXPIANGCONST][3], piobp[MAXPIANGCONST];
44 } angkp1;
45
46 EXTERN struct t_angf {
47 int use_angf, nfang;
48 char kftype[MAXANGCONST][10];
49 float fcon[MAXANGCONST],fc0[MAXANGCONST],fc1[MAXANGCONST],fc2[MAXANGCONST];
50 float fc3[MAXANGCONST],fc4[MAXANGCONST],fc5[MAXANGCONST],fc6[MAXANGCONST];
51 } angf;
52
53 EXTERN int Missing_constants;
54
55 void kangle()
56 {
57 long int pi_mask, aromatic_mask;
58 int i, j, k, izero, ierr, nh, nRc, nligand, jji;
59 int ia, ib, ic, ie;
60 int iaa, ibb, icc;
61 int cl_a, cl_b, cl_c;
62 int nbo1, nb1, nb2;
63 char pa[4],pb[4],pc[4],pt[10], pt1[10], pt2[10],pt3[10], pt4[10], pt5[10], pt6[10], pt7[10];
64 char zero[4];
65 float bkan, at;
66
67 pi_mask = (1L << PI_MASK);
68 aromatic_mask = (1 << AROMATIC_MASK);
69
70 if( angles.nang != 0 )
71 {
72 for( i = 0; i < angles.nang; i++ )
73 {
74 ia = angles.i13[i][0];
75 ib = angles.i13[i][1];
76 ic = angles.i13[i][2];
77 iaa = atom[ia].type;
78 ibb = atom[ib].type;
79 icc = atom[ic].type;
80 cl_a = atom[ia].tclass;
81 cl_b = atom[ib].tclass;
82 cl_c = atom[ic].tclass;
83
84 /* get bond index each bond for MMFF94 */
85 nb1 = find_bond(ia,ib);
86 nb2 = find_bond(ib,ic);
87
88 /* special metal cases */
89 if( iaa >= 300)
90 cl_a = 300;
91
92 if( ibb >= 300)
93 cl_b = 300;
94
95 if( icc >= 300)
96 cl_c = 300;
97
98 /* special MMX cases */
99 if ( field.type == MMX)
100 {
101 if( iaa == 40 )
102 iaa = 2;
103
104 if( ibb == 40 )
105 ibb = 2;
106
107 if( icc == 40 )
108 icc = 2;
109
110 }
111 numeral(iaa,pa,3);
112 numeral(ibb,pb,3);
113 numeral(icc,pc,3);
114 if( iaa <= icc )
115 {
116 strcpy(pt,pa);
117 strcat(pt,pb);
118 strcat(pt,pc);
119 }
120 if( iaa > icc )
121 {
122 strcpy(pt,pc);
123 strcat(pt,pb);
124 strcat(pt,pa);
125 }
126 izero = 0;
127 strcpy(zero," 0");
128 strcpy(pt1,pa); strcat(pt1,pb); strcat(pt1,zero);
129 strcpy(pt2,pc); strcat(pt2,pb); strcat(pt2,zero);
130 strcpy(pt3,zero); strcat(pt3,pb); strcat(pt3,pc);
131 strcpy(pt4,zero); strcat(pt4,pb); strcat(pt4,pa);
132 strcpy(pt5,zero); strcat(pt5,pb); strcat(pt5,zero);
133
134 numeral(cl_a,pa,3);
135 numeral(cl_b,pb,3);
136 numeral(cl_c,pc,3);
137 if (cl_a < cl_c)
138 {
139 strcpy(pt6,pa); strcat(pt6,pb); strcat(pt6,pc);
140 }else
141 {
142 strcpy(pt6,pc); strcat(pt6,pb); strcat(pt6,pa);
143 }
144
145 if (field.type == MMFF94)
146 {
147 cl_a = atom_k.tclass1[atom[ia].type];
148 cl_b = atom_k.tclass1[atom[ib].type];
149 cl_c = atom_k.tclass1[atom[ic].type];
150 numeral(cl_a,pa,3);
151 numeral(cl_b,pb,3);
152 numeral(cl_c,pc,3);
153 if (cl_a < cl_c)
154 {
155 strcpy(pt7,pa); strcat(pt7,pb); strcat(pt7,pc);
156 }else
157 {
158 strcpy(pt7,pc); strcat(pt7,pb); strcat(pt7,pa);
159 }
160 }
161 angles.acon[i] = 0.0;
162 angles.anat[i] = 0.0;
163 angles.angtype[i] = HARMONIC;
164
165 ierr = FALSE;
166 /* check TS constants */
167 /* check delocalized angles in MMFF94 cyclopropanes */
168 if (field.type == MMFF94 && angk1.ndel3 > 0)
169 {
170 if (isbond(ia,ic) )
171 {
172 nbo1 = bonds_ff.index[nb1] + bonds_ff.index[nb2];
173 if (nbo1 == 1 || nbo1 == 2)
174 {
175 for(j=0; j < angk1.ndel3; j++)
176 {
177 if (strcmp(angk1.kdel3[j],pt) == 0)
178 {
179 angles.acon[i] = angk1.condel3[j];
180 angles.anat[i] = angk1.angdel3[j][0];
181 if (nbo1 == 1)
182 angles.index[i] = 5;
183 else if (nbo1 == 2)
184 angles.index[i] = 6;
185 ierr = TRUE;
186 break;
187 }
188 }
189 }
190
191 if (ierr == TRUE)
192 goto L_10;
193 }
194 }
195 /* check three membered rings */
196 if ( (rings.nring3 > 0) && (angk1.nang3 > 0) )
197 {
198 if (isbond(ia,ic))
199 {
200 /* search three membered ring parameters */
201 for (j = 0; j < angk1.nang3; j++)
202 {
203 if ( (strcmp(angk1.ktype3[j],pt) == 0) || (strcmp(angk1.ktype3[j],pt1) == 0) ||
204 (strcmp(angk1.ktype3[j],pt2) == 0) || (strcmp(angk1.ktype3[j],pt3) == 0) ||
205 (strcmp(angk1.ktype3[j],pt4) == 0) || (strcmp(angk1.ktype3[j],pt5) == 0))
206 {
207 if ( angk1.ang3[j][1] == 0.0 && angk1.ang3[j][2] == 0)
208 {
209 angles.acon[i] = angk1.con3[j];
210 angles.anat[i] = angk1.ang3[j][0];
211 angles.index[i] = 3;
212 ierr = TRUE;
213 break;
214 } else
215 {
216 nh = 0;
217 for (k=0; k < MAXIAT; k++)
218 {
219 ie = atom[ib].iat[k];
220 if ( (atom[ie].type == 5 || atom[ie].type == 36) && (ie != ia) && (ie != ic))
221 nh++;
222 }
223 angles.acon[i] = angk1.con3[j];
224 angles.anat[i] = angk1.ang3[j][nh];
225 angles.index[i] = 3;
226 ierr = TRUE;
227 break;
228 }
229 }
230 }
231 if (ierr == TRUE)
232 goto L_10;
233 else
234 {
235 Missing_constants = TRUE;
236 }
237 }
238 }
239 /* check delocalized angles in MMFF94 cyclobutanes */
240 if (field.type == MMFF94 && angk1.ndel4 > 0)
241 {
242 nRc = is_ring43(ia, ib, ic);
243 if (nRc == TRUE)
244 {
245 nbo1 = bonds_ff.index[nb1] + bonds_ff.index[nb2];
246 if (nbo1 == 1 || nbo1 == 2)
247 {
248 for(j=0; j < angk1.ndel4; j++)
249 {
250 if (strcmp(angk1.kdel4[j],pt) == 0)
251 {
252 angles.acon[i] = angk1.condel4[j];
253 angles.anat[i] = angk1.angdel4[j][0];
254 if (nbo1 == 1)
255 angles.index[i] = 7;
256 else if (nbo1 == 2)
257 angles.index[i] = 8;
258 ierr = TRUE;
259 break;
260 }
261 }
262 }
263 if (ierr == TRUE)
264 goto L_10;
265 }
266 }
267 /* check four memebered rings */
268 if ( (rings.nring4 > 0) && (angk1.nang4 > 0) )
269 {
270 nRc = is_ring43(ia, ib, ic);
271 if (nRc == TRUE)
272 {
273 /* search four membered ring parameters */
274 for (j = 0; j < angk1.nang4; j++)
275 {
276 if ( (strcmp(angk1.ktype4[j],pt) == 0) || (strcmp(angk1.ktype4[j],pt1) == 0) ||
277 (strcmp(angk1.ktype4[j],pt2) == 0) || (strcmp(angk1.ktype4[j],pt3) == 0) ||
278 (strcmp(angk1.ktype4[j],pt4) == 0) || (strcmp(angk1.ktype4[j],pt5) == 0))
279 {
280 if ( angk1.ang4[j][1] == 0.0 && angk1.ang4[j][2] == 0)
281 {
282 angles.acon[i] = angk1.con4[j];
283 angles.anat[i] = angk1.ang4[j][0];
284 angles.index[i] = 4;
285 ierr = TRUE;
286 break;
287 } else
288 {
289 nh = 0;
290 for (k=0; k < MAXIAT; k++)
291 {
292 ie = atom[ib].iat[k];
293 if ( (atom[ie].type == 5 || atom[ie].type == 36) && (ie != ia) && (ie != ic))
294 nh++;
295 }
296 angles.acon[i] = angk1.con4[j];
297 angles.anat[i] = angk1.ang4[j][nh];
298 angles.index[i] = 4;
299 ierr = TRUE;
300 break;
301 }
302 }
303 }
304 if (ierr == TRUE)
305 goto L_10;
306 else
307 {
308 Missing_constants = TRUE;
309 }
310 }
311 }
312 /* check five membered rings */
313 if ( (rings.nring5 > 0) && (angk1.nang5 > 0) )
314 {
315 nRc = is_ring53(ia,ib, ic);
316 if (nRc == TRUE)
317 {
318 /* search five memebered ring parameters */
319 for (j = 0; j < angk1.nang5; j++)
320 {
321 if ( (strcmp(angk1.ktype5[j],pt) == 0) || (strcmp(angk1.ktype5[j],pt1) == 0) ||
322 (strcmp(angk1.ktype5[j],pt2) == 0) || (strcmp(angk1.ktype5[j],pt3) == 0) ||
323 (strcmp(angk1.ktype5[j],pt4) == 0) || (strcmp(angk1.ktype5[j],pt5) == 0) )
324 {
325 if ( angk1.ang5[j][1] == 0.0 && angk1.ang5[j][2] == 0)
326 {
327 angles.acon[i] = angk1.con5[j];
328 angles.anat[i] = angk1.ang5[j][0];
329 angles.index[i] = 0;
330 ierr = TRUE;
331 break;
332 } else
333 {
334 nh = 0;
335 for (k=0; k < MAXIAT; k++)
336 {
337 ie = atom[ib].iat[k];
338 if ( (atom[ie].type == 5 || atom[ie].type == 36) && (ie != ia) && (ie != ic))
339 nh++;
340 }
341 angles.acon[i] = angk1.con5[j];
342 angles.anat[i] = angk1.ang5[j][nh];
343 angles.index[i] = 0;
344 ierr = TRUE;
345 break;
346 }
347 }
348 }
349 if (ierr == TRUE)
350 goto L_10;
351 // don't fail on missing 5 ring parameters
352 }
353 }
354 /* check delocalized angles in MMFF94 */
355 if (field.type == MMFF94 && angk1.ndel > 0)
356 {
357 if ( bonds_ff.index[nb1] == 1 || bonds_ff.index[nb2] == 1 )
358 {
359 nbo1 = bonds_ff.index[nb1] + bonds_ff.index[nb2];
360 for(j=0; j < angk1.ndel; j++)
361 {
362 if (strcmp(angk1.kdel[j],pt) == 0)
363 {
364 angles.acon[i] = angk1.condel[j];
365 angles.anat[i] = angk1.angdel[j][0];
366 if (nbo1 == 1)
367 angles.index[i] = 1;
368 else if (nbo1 == 2)
369 angles.index[i] = 2;
370 ierr = TRUE;
371 break;
372 }
373 }
374 if (ierr == TRUE)
375 goto L_10;
376 }
377 }
378 // look for fourier angles
379 if (angf.nfang > 0)
380 {
381 for (j=0; j < angf.nfang; j++)
382 {
383 if (strcmp(angf.kftype[j],pt) == 0)
384 {
385 angles.fcon[i] = angf.fcon[j];
386 angles.c0[i] = angf.fc0[j];
387 angles.c1[i] = angf.fc1[j];
388 angles.c2[i] = angf.fc2[j];
389 angles.c3[i] = angf.fc3[j];
390 angles.c4[i] = angf.fc4[j];
391 angles.c5[i] = angf.fc5[j];
392 angles.c6[i] = angf.fc6[j];
393 angles.angtype[i] = FOURIER;
394 ierr = TRUE;
395 break;
396 }
397 }
398 if (ierr == TRUE)
399 goto L_10;
400 }
401
402 /* check regular parameters for specific angles*/
403 for (j = 0; j < angk1.nang; j++)
404 {
405 if ( (strcmp(angk1.ktype[j],pt) == 0) )
406 {
407 if ( angk1.ang[j][1] == 0.0 && angk1.ang[j][2] == 0)
408 {
409 angles.acon[i] = angk1.con[j];
410 angles.anat[i] = angk1.ang[j][0];
411 angles.index[i] = 0;
412 ierr = TRUE;
413 break;
414 } else
415 {
416 nh = 0;
417 for (k=0; k < MAXIAT; k++)
418 {
419 ie = atom[ib].iat[k];
420 if ( (atom[ie].atomnum == 1) && (ie != ia) && (ie != ic))
421 nh++;
422 }
423 angles.acon[i] = angk1.con[j];
424 angles.anat[i] = angk1.ang[j][nh];
425 angles.index[i] = 0;
426 ierr = TRUE;
427 break;
428 }
429 }
430 }
431 if (ierr == TRUE)
432 goto L_10;
433 // did not find any specific look for generalized
434 for (j = 0; j < angk1.nang; j++)
435 {
436 if ( (strcmp(angk1.ktype[j],pt) == 0) || (strcmp(angk1.ktype[j],pt1) == 0) ||
437 (strcmp(angk1.ktype[j],pt2) == 0) || (strcmp(angk1.ktype[j],pt3) == 0) ||
438 (strcmp(angk1.ktype[j],pt4) == 0) || (strcmp(angk1.ktype[j],pt5) == 0) ||
439 (strcmp(angk1.ktype[j],pt6) == 0))
440 {
441 if ( angk1.ang[j][1] == 0.0 && angk1.ang[j][2] == 0)
442 {
443 angles.acon[i] = angk1.con[j];
444 angles.anat[i] = angk1.ang[j][0];
445 angles.index[i] = 0;
446 ierr = TRUE;
447 break;
448 } else
449 {
450 nh = 0;
451 for (k=0; k < MAXIAT; k++)
452 {
453 ie = atom[ib].iat[k];
454 if ( (atom[ie].atomnum == 1) && (ie != ia) && (ie != ic))
455 nh++;
456 }
457 angles.acon[i] = angk1.con[j];
458 angles.anat[i] = angk1.ang[j][nh];
459 angles.index[i] = 0;
460 ierr = TRUE;
461 break;
462 }
463 }
464 }
465 if (ierr == TRUE)
466 goto L_10;
467 /* check MMFF for class1 parameters */
468 if (field.type == MMFF94)
469 {
470 for (j = 0; j < angk1.nang; j++)
471 {
472 if ( (strcmp(angk1.ktype[j],pt7) == 0) )
473 {
474 if ( angk1.ang[j][1] == 0.0 && angk1.ang[j][2] == 0)
475 {
476 angles.acon[i] = angk1.con[j];
477 angles.anat[i] = angk1.ang[j][0];
478 angles.index[i] = 0;
479 ierr = TRUE;
480 break;
481 } else
482 {
483 nh = 0;
484 for (k=0; k < MAXIAT; k++)
485 {
486 ie = atom[ib].iat[k];
487 if ( (atom[ie].atomnum == 1) && (ie != ia) && (ie != ic))
488 nh++;
489 }
490 angles.acon[i] = angk1.con[j];
491 angles.anat[i] = angk1.ang[j][nh];
492 angles.index[i] = 0;
493 ierr = TRUE;
494 break;
495 }
496 }
497 }
498 if (ierr == TRUE)
499 goto L_10;
500 }
501 if (ierr != TRUE)
502 {
503 bkan = get_general_angle(ibb);
504 angles.acon[i] = 0.80;
505 angles.anat[i] = bkan;
506 angles.index[i] = 0;
507 // ==== done with angle lookup =============
508 L_10:
509 if (field.type == MMX && ierr != TRUE)
510 {
511 if (ibb >= 300)
512 {
513 nligand = 0;
514 jji = 0;
515 for ( k=0; k < MAXIAT; k++)
516 {
517 if (atom[ib].iat[k] != 0 && atom[ib].bo[k] != 9)
518 jji++;
519 if (atom[ib].bo[k] == 9)
520 nligand++;
521 }
522 if (jji == 2 && nligand == 0)
523 {
524 angles.acon[i] = .1;
525 angles.anat[i] = 180.;
526 }
527 if ((jji == 3 && nligand == 0) || (jji == 2 && nligand == 1) )
528 {
529 angles.acon[i] = .1;
530 angles.anat[i] = 120.;
531 }
532 if ((jji == 4 && nligand == 0) || (jji == 3 && nligand == 1) || (jji == 2 && nligand == 2))
533 {
534 if (jji == 4 && (atom[ib].type >= 300 && (atom[ib].flags & (1L << SATMET_MASK) &&
535 atom[ib].flags & (1L << GT18e_MASK) && atom[ib].flags & (1L << SQPLAN_MASK) &&
536 !(atom[ib].flags & (1L << LOWSPIN_MASK)))) )
537 {
538 angles.acon[i] = .0001;
539 angles.anat[i] = 90.;
540 } else
541 {
542 angles.acon[i] = .1;
543 angles.anat[i] = 109.5;
544 }
545 }
546
547 }
548 }
549 }
550 }
551 }
552 }
553 // ===================================
554 // generalized angles for MMFF with no hydrogen minimization
555 float get_general_angle(int ia)
556 {
557 float angs[80] = {
558 109.0, 120.0, 120.0, 180.0, 0.000, 109.0, 0.000, 109.0, 120.0, 118.0,
559 0.000, 0.000, 0.000, 0.000, 109.0, 109.0, 120.0, 109.0, 109.0, 90.00,
560 0.000, 60.00, 0.000, 0.000, 109.0, 109.0, 0.000, 0.000, 0.000, 90.00,
561 0.000, 0.000, 0.000, 109.0, 0.000, 0.000, 120.0, 120.0, 120.0, 120.0,
562 120.0, 180.0, 120.0, 120.0, 120.0, 120.0, 0.000, 0.000, 109.0, 0.000,
563 109.0, 0.000, 180.0, 120.0, 120.0, 120.0, 120.0, 120.0, 120.0, 180.0,
564 180.0, 120.0, 120.0, 120.0, 120.0, 120.0, 120.0, 109.0, 120.0, 109.0,
565 0.000, 109.0, 109.0, 120.0, 120.0, 109.0, 0.000, 109.0, 109.0, 109.0 };
566
567 if (ia < 80)
568 return (angs[ia-1]);
569 else
570 return 109.0;
571 }
572 /* --------------------------------- */
573 long kang(int i,int j,int k)
574 {
575 long int kang_v;
576 static long i10000 = 10000;
577 static long i100 = 100;
578
579 kang_v = i10000*i + i100*j + k;
580 return( kang_v );
581 }
582
583 /* -------------------- */
584 long ipack(int i,int j,int k,int l)
585 {
586 long int ipack_v;
587 static long itbig = 200000000;
588 static long itmid = 2000000;
589 static long ithou = 1000;
590
591 ipack_v = itbig*i + itmid*j + ithou*k + l;
592 return( ipack_v );
593 }