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, 5 months ago) by wdelano
File size: 24089 byte(s)
Log Message:
branch created
Line User Rev File contents
1 tjod 3 #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     }