1 |
#define EXTERN extern |
2 |
#include "pcwin.h" |
3 |
#include "pcmod.h" |
4 |
#include "field.h" |
5 |
#include "atom_k.h" |
6 |
#include "draw.h" |
7 |
|
8 |
void hadd(void); |
9 |
void hdel(int); |
10 |
void type(void); |
11 |
void nhadd(void); |
12 |
void hcoord(void); |
13 |
static void revec(void); |
14 |
static void xyplan(float *, float *, float *); |
15 |
static void xaxis(float *,float *,float *); |
16 |
static void getvec(void); |
17 |
static void reseq(void); |
18 |
|
19 |
EXTERN struct t_minim_control { |
20 |
int type, method, field, added_const; |
21 |
char added_path[256],added_name[256]; |
22 |
} minim_control; |
23 |
|
24 |
void hadd(void) |
25 |
{ |
26 |
type(); |
27 |
if (field.type == MM3) |
28 |
{ |
29 |
nhadd(); |
30 |
hcoord(); |
31 |
hdel(1); |
32 |
} else if (field.type == AMBER) |
33 |
{ |
34 |
nhadd(); |
35 |
hcoord(); |
36 |
hdel(1); |
37 |
} else if (field.type == OPLSAA) |
38 |
{ |
39 |
nhadd(); |
40 |
hcoord(); |
41 |
hdel(1); |
42 |
} else if (field.type == MMFF94) |
43 |
{ |
44 |
nhadd(); |
45 |
hcoord(); |
46 |
hdel(1); |
47 |
} else |
48 |
{ |
49 |
set_atomtypes(MMX); |
50 |
nhadd(); |
51 |
hcoord(); |
52 |
} |
53 |
type(); |
54 |
generate_bonds(); |
55 |
|
56 |
} |
57 |
|
58 |
void nhadd() |
59 |
{ |
60 |
int i, i1, ik, im, it, itads, j, jk, nhadds, newatom; |
61 |
int jji, bo,ncount; |
62 |
static int iadd[MAXVDW]={ 0,0,0,0,0,3,1,6,0,1, |
63 |
0,0,0,0,3,6,3,0,0,0, |
64 |
0,0,0,0,6,0,0,0,0,0, |
65 |
0,0,0,3,0,0,6,1,1,0, |
66 |
0,4,0,0,0,6,0,6,0,0, |
67 |
3,3,3,0,6,0,0,0,0,0, |
68 |
0,0,0,0,0,1,6,0,0,0, |
69 |
0,0,0,0,0}; |
70 |
|
71 |
/* this routine calculates the #'s of h to be added |
72 |
* to carbon and packs it symbolically into iat(i,j) table |
73 |
* as well as into the connection table via the appropriate |
74 |
* manipulation instructions. hydrogen coordinates |
75 |
* are not added by this routine. lp are also added to oxygen. */ |
76 |
/* 0 to disregard : 1 and 4 for lp's only |
77 |
* 2 and 5 for h's only and 3 and 6 for h and lp |
78 |
* add 3 for trigonal atoms */ |
79 |
|
80 |
ncount = natom; |
81 |
for( i = 1 ; i <= ncount; i++ ) |
82 |
{ |
83 |
if( atom[i].mmx_type > 0 && (atom[i].mmx_type < 300) && atom[i].atomnum != 1 ) |
84 |
{ |
85 |
jji = 0; |
86 |
bo = 0; |
87 |
for (j=0; j < MAXIAT; j++) |
88 |
{ |
89 |
if (atom[i].iat[j] != 0 && atom[i].bo[j] != 9) |
90 |
{ |
91 |
jji ++; |
92 |
bo += atom[i].bo[j]; |
93 |
} |
94 |
} |
95 |
nhadds = atom_def.ligands[atom[i].mmx_type] - jji; |
96 |
|
97 |
if (atom[i].mmx_type >= 11 && atom[i].mmx_type <= 14) |
98 |
nhadds = 0; |
99 |
if (atom[i].mmx_type == 41 || atom[i].mmx_type == 46) |
100 |
nhadds = atom_def.ligands[atom[i].mmx_type] - bo; |
101 |
if (atom[i].mmx_type == 49 || atom[i].mmx_type == 50 || atom[i].mmx_type == 51) |
102 |
nhadds = atom_def.ligands[atom[i].mmx_type] - bo; |
103 |
|
104 |
/* add atom and bond it to atom i */ |
105 |
if( nhadds > 0 ) |
106 |
{ |
107 |
for( im = 1; im <= nhadds; im++ ) |
108 |
{ |
109 |
newatom = make_atom(5,0.F,0.F,0.F,"H"); |
110 |
set_atomdata(newatom,5,5,5,34,29); |
111 |
if (newatom == -1) |
112 |
{ |
113 |
message_alert("Error in nhadd newatom = -1","Error"); |
114 |
hdel(0); |
115 |
return; |
116 |
} |
117 |
make_bond(i,newatom,1); |
118 |
} |
119 |
} |
120 |
} |
121 |
} |
122 |
|
123 |
/* *** check for lp's *** */ |
124 |
if (field.type == MMX) |
125 |
{ |
126 |
ncount = natom; |
127 |
for( i = 1 ; i <= ncount; i++ ) |
128 |
{ |
129 |
it = iadd[atom[i].mmx_type-1]; |
130 |
if( it >= 1 ) |
131 |
{ |
132 |
itads = atom[i].mmx_type ; |
133 |
if (it == 1 || it == 3 || it == 4 || it == 6 ) |
134 |
{ |
135 |
nhadds = 4; |
136 |
for( i1 = 0; i1 < MAXIAT; i1++ ) |
137 |
{ |
138 |
if( atom[i].bo[i1] != 0 && atom[i].bo[i1] != 9 ) |
139 |
nhadds = nhadds - atom[i].bo[i1]; |
140 |
} |
141 |
if( itads == 6 || itads == 15 || itads == 34 || itads == 42 || itads == 53 || itads == 48 ) |
142 |
{ |
143 |
for( ik = 0; ik < MAXIAT; ik++ ) |
144 |
{ |
145 |
if( atom[i].iat[ik] != 0 && atom[i].bo[ik] != 9 ) |
146 |
{ |
147 |
jk = atom[atom[i].iat[ik]].mmx_type; |
148 |
if( jk == 2 || jk == 3 || jk == 4 || jk == 9 |
149 |
|| jk == 10 || jk == 30 || jk == 37 || jk == 40 ) |
150 |
{ |
151 |
nhadds = nhadds - 1; |
152 |
break; |
153 |
} |
154 |
} |
155 |
} |
156 |
} |
157 |
|
158 |
/* correction for sulfoxide */ |
159 |
if( itads == 17 ) |
160 |
nhadds = nhadds + 1; |
161 |
|
162 |
if (itads == 66) |
163 |
nhadds = 2; |
164 |
|
165 |
if( nhadds > 0 ) |
166 |
{ |
167 |
for( i1 = 1; i1 <= nhadds; i1++ ) |
168 |
{ |
169 |
newatom = make_atom(20, 0.F,0.F,0.F,""); |
170 |
set_atomdata(newatom,20,0,0,0,0); |
171 |
make_bond(i,newatom,1); |
172 |
} |
173 |
} |
174 |
} |
175 |
} |
176 |
} |
177 |
} |
178 |
return; |
179 |
} |
180 |
// ================================== |
181 |
void hcoord() |
182 |
{ |
183 |
int i, i1, i3, i4, i5, i6, i7, i8, |
184 |
i9, iadj[10], ii, ii1, ii2, ii3, iii, |
185 |
ij, ijk[3], it, j, j1, j2, j3, j4, j5, j8, |
186 |
j9, k1, k2, l1, na, nhadds, ni; |
187 |
float adj, dis, xpnt, ypnt, zpnt; |
188 |
|
189 |
getvec(); |
190 |
for( i = 1; i <= natom; i++ ) |
191 |
{ |
192 |
it = atom[i].mmx_type; |
193 |
if( it < 300) |
194 |
{ |
195 |
for( j = 0; j < 10; j++ ) |
196 |
iadj[j] = 0; |
197 |
|
198 |
for( j = 0; j < 3; j++ ) |
199 |
ijk[j] = 0; |
200 |
|
201 |
nhadds = 0; |
202 |
l1 = 0; |
203 |
for( j = 0; j < MAXIAT; j++ ) |
204 |
{ |
205 |
if( atom[i].iat[j] != 0 && atom[i].bo[j] != 9 ) |
206 |
{ |
207 |
if( atom[atom[i].iat[j]].mmx_type == 5 || atom[atom[i].iat[j]].mmx_type == 20 ) |
208 |
nhadds++ ; |
209 |
else |
210 |
{ |
211 |
iadj[l1] = atom[i].iat[j]; |
212 |
l1++; |
213 |
} |
214 |
} |
215 |
} |
216 |
if( nhadds != 0 ) |
217 |
{ |
218 |
/* **** calculate the coordinates of the center we |
219 |
* are going to move to the x-axis */ |
220 |
xpnt = atom[i].x; |
221 |
ypnt = atom[i].y; |
222 |
zpnt = atom[i].z; |
223 |
for( ii = 1; ii <= (natom + 4); ii++ ) |
224 |
{ |
225 |
atom[ii].x = atom[ii].x - xpnt; |
226 |
atom[ii].y = atom[ii].y - ypnt; |
227 |
atom[ii].z = atom[ii].z - zpnt; |
228 |
} |
229 |
xpnt = 0.F; |
230 |
ypnt = 0.F; |
231 |
zpnt = 0.F; |
232 |
adj = 0.F; |
233 |
for( i1 = 0; i1 < 10; i1++ ) |
234 |
{ |
235 |
if( iadj[i1] != 0 ) |
236 |
{ |
237 |
na = iadj[i1]; |
238 |
adj += 1.0001F ; |
239 |
xpnt += atom[na].x; |
240 |
ypnt += atom[na].y; |
241 |
zpnt += atom[na].z; |
242 |
} |
243 |
} |
244 |
if( adj > 0.01F ) |
245 |
{ |
246 |
xpnt /= adj; |
247 |
ypnt /= adj; |
248 |
zpnt /= adj; |
249 |
xaxis( &xpnt, &ypnt, &zpnt ); |
250 |
if( xpnt > 0.F ) |
251 |
{ |
252 |
for( i3 = 1; i3 <= (natom + 4); i3++ ) |
253 |
atom[i3].x = -atom[i3].x; |
254 |
} |
255 |
/* **** after moving, check the # of h added to the carbon i */ |
256 |
if( nhadds == 1 ) |
257 |
{ |
258 |
/* **** there is only 1 h added to the carbon i */ |
259 |
for( i4 = 0; i4 < MAXIAT; i4++ ) |
260 |
{ |
261 |
if( atom[i].iat[i4] != 0 && atom[i].bo[i4] != 9 ) |
262 |
{ |
263 |
i5 = atom[i].iat[i4]; |
264 |
if( atom[i5].mmx_type == 5 || atom[i5].mmx_type == 20 ) |
265 |
{ |
266 |
/* i5 is the atom to attach */ |
267 |
if (adj < 2.0) // fix for OH in MMFF |
268 |
{ |
269 |
atom[i5].x = 0.70; |
270 |
atom[i5].y = 0.70; |
271 |
atom[i5].z = 0.00; |
272 |
} else |
273 |
{ |
274 |
atom[i5].x = 1.10F; |
275 |
if( atom[i5].mmx_type == 20 ) |
276 |
atom[i5].x = 0.6F; |
277 |
atom[i5].y = 0.00F; |
278 |
atom[i5].z = 0.00F; |
279 |
} |
280 |
if( atom[i].atomnum == 16 ) |
281 |
{ |
282 |
atom[i5].x = 0.F; |
283 |
atom[i5].y = 0.3796F; |
284 |
atom[i5].z = 1.043F; |
285 |
} |
286 |
} |
287 |
} |
288 |
} |
289 |
} |
290 |
/* **** there are 2 h added to the carbon i */ |
291 |
if( nhadds == 2 ) |
292 |
{ |
293 |
if( it == 2 || it == 3 || it == 29 || it == 30 || it == 9 || it == 41 || it == 46 |
294 |
|| it == 48 || it == 7 || it == 37 || it == 39 |
295 |
|| it == 42 || it == 38 || it == 40 || it == 44 || it == 66) |
296 |
{ |
297 |
for( i8 = 0; i8 < 10; i8++ ) |
298 |
{ |
299 |
if( iadj[i8] != 0 ) |
300 |
{ |
301 |
ii1 = iadj[i8]; |
302 |
for( ii2 = 0; ii2 < MAXIAT; ii2++ ) |
303 |
{ |
304 |
ii3 = 0; |
305 |
if( atom[ii1].iat[ii2] != 0 && atom[ii1].iat[ii2] != i && atom[ii1].bo[ii2] != 9 ) |
306 |
{ |
307 |
ii3 = atom[ii1].iat[ii2]; |
308 |
xpnt = atom[ii3].x; |
309 |
ypnt = atom[ii3].y; |
310 |
zpnt = atom[ii3].z; |
311 |
xyplan( &xpnt, &ypnt, &zpnt ); |
312 |
ij = 0; |
313 |
for( i9 = 0; i9 < MAXIAT; i9++ ) |
314 |
{ |
315 |
if( atom[i].iat[i9] != 0 && atom[i].bo[i9] != 9 ) |
316 |
{ |
317 |
iii = atom[i].iat[i9]; |
318 |
if( atom[iii].mmx_type == 5 || atom[iii].mmx_type == 20 ) |
319 |
{ |
320 |
ijk[ij] = iii; |
321 |
ij++; |
322 |
if( ij == 2 ) |
323 |
{ |
324 |
/* *** =ch2 *** */ |
325 |
atom[ijk[0]].x = 0.5582F; |
326 |
atom[ijk[1]].x = 0.5582F; |
327 |
atom[ijk[0]].y = 0.9478F; |
328 |
atom[ijk[1]].y = -0.9478F; |
329 |
atom[ijk[0]].z = 0.00F; |
330 |
atom[ijk[1]].z = 0.00F; |
331 |
if( it == 7 || it == 38 || it == 42 || it == 39 || it == 66) |
332 |
{ |
333 |
atom[ijk[0]].x = 0.30F; |
334 |
atom[ijk[1]].x = 0.30F; |
335 |
atom[ijk[0]].y = 0.52F; |
336 |
atom[ijk[1]].y = -0.52F; |
337 |
} |
338 |
if( atom[i].mmx_type == 37 ) |
339 |
{ |
340 |
atom[ijk[1]].x = 0.35F; |
341 |
atom[ijk[1]].y = -0.5F; |
342 |
} |
343 |
} |
344 |
} |
345 |
} |
346 |
} |
347 |
break; |
348 |
} |
349 |
if (ii3 == 0) // C=O with no lone pairs |
350 |
{ |
351 |
ij = 0; |
352 |
for( i9 = 0; i9 < MAXIAT; i9++ ) |
353 |
{ |
354 |
if( atom[i].iat[i9] != 0 && atom[i].bo[i9] != 9 ) |
355 |
{ |
356 |
iii = atom[i].iat[i9]; |
357 |
if( atom[iii].mmx_type == 5 || atom[iii].mmx_type == 20 ) |
358 |
{ |
359 |
ijk[ij] = iii; |
360 |
ij++; |
361 |
if( ij == 2 ) |
362 |
{ |
363 |
/* *** =ch2 *** */ |
364 |
atom[ijk[0]].x = 0.5582F; |
365 |
atom[ijk[1]].x = 0.5582F; |
366 |
atom[ijk[0]].y = 0.9478F; |
367 |
atom[ijk[1]].y = -0.9478F; |
368 |
atom[ijk[0]].z = 0.00F; |
369 |
atom[ijk[1]].z = 0.00F; |
370 |
if( it == 7 || it == 38 || it == 42 || it == 39 || it == 66) |
371 |
{ |
372 |
atom[ijk[0]].x = 0.30F; |
373 |
atom[ijk[1]].x = 0.30F; |
374 |
atom[ijk[0]].y = 0.52F; |
375 |
atom[ijk[1]].y = -0.52F; |
376 |
} |
377 |
if( atom[i].mmx_type == 37 ) |
378 |
{ |
379 |
atom[ijk[1]].x = 0.35F; |
380 |
atom[ijk[1]].y = -0.5F; |
381 |
} |
382 |
} |
383 |
} |
384 |
} |
385 |
} |
386 |
break; |
387 |
} |
388 |
} |
389 |
} |
390 |
} |
391 |
} else |
392 |
{ |
393 |
/* ** there are 2 atoms attached to carbon i, |
394 |
* one of which is "ni" and move these two atoms to xy-plane */ |
395 |
for( i6 = 0; i6 < 10; i6++ ) |
396 |
{ |
397 |
if( iadj[i6] != 0 ) |
398 |
{ |
399 |
ni = iadj[i6]; |
400 |
xpnt = atom[ni].x; |
401 |
ypnt = atom[ni].y; |
402 |
zpnt = atom[ni].z; |
403 |
xyplan( &xpnt, &ypnt, &zpnt ); |
404 |
} |
405 |
} |
406 |
ij = 0; |
407 |
for( i7 = 0; i7 < MAXIAT; i7++ ) |
408 |
{ |
409 |
if( atom[i].iat[i7] != 0 && atom[i].bo[i7] != 9 ) |
410 |
{ |
411 |
iii = atom[i].iat[i7]; |
412 |
if( atom[iii].mmx_type == 5 || atom[iii].mmx_type == 20 ) |
413 |
{ |
414 |
ijk[ij] = iii; |
415 |
ij++ ; |
416 |
if( ij == 2 ) |
417 |
{ |
418 |
/* *** -ch2- type 1 atom *** */ |
419 |
atom[ijk[0]].x = 0.6758F; |
420 |
atom[ijk[1]].x = 0.6758F; |
421 |
atom[ijk[0]].y = 0.00F; |
422 |
atom[ijk[1]].y = 0.00F; |
423 |
atom[ijk[0]].z = 0.8807F; |
424 |
atom[ijk[1]].z = -0.8807F; |
425 |
if( atom[i].mmx_type == 25 || atom[i].mmx_type == 46 || atom[i].mmx_type == 67) |
426 |
{ |
427 |
atom[ijk[0]].x = 0.3796F; |
428 |
atom[ijk[1]].x = 0.3796F; |
429 |
atom[ijk[0]].y = 1.043F; |
430 |
atom[ijk[1]].y = -0.5215F; |
431 |
atom[ijk[0]].z = 0.0000F; |
432 |
atom[ijk[1]].z = 0.9032F; |
433 |
} |
434 |
if( atom[i].mmx_type == 22 || atom[i].mmx_type == 30 ) |
435 |
{ |
436 |
/* *** -ch2- type 22 and 3o atoms *** */ |
437 |
atom[ijk[0]].x = 0.2F; |
438 |
atom[ijk[1]].x = 0.2F; |
439 |
atom[ijk[0]].y = 0.0F; |
440 |
atom[ijk[1]].y = 0.0F; |
441 |
atom[ijk[0]].z = 1.0F; |
442 |
atom[ijk[1]].z = -1.0F; |
443 |
} |
444 |
/* *** -o(lp)2- *** */ |
445 |
if (atom[i].mmx_type == 6 && (atom[iadj[0]].mmx_type == 3 || |
446 |
atom[iadj[0]].mmx_type == 2 || atom[iadj[0]].mmx_type == 30 || |
447 |
atom[iadj[0]].mmx_type == 37 || atom[iadj[0]].mmx_type == 40) ) |
448 |
{ |
449 |
atom[ijk[0]].x = 0.450F; |
450 |
atom[ijk[1]].x = 0.2488F; |
451 |
atom[ijk[0]].z = 0.0F; |
452 |
atom[ijk[1]].z = 0.0F; |
453 |
atom[ijk[0]].y = 0.9500F; |
454 |
atom[ijk[1]].y = -0.5460F; |
455 |
} |
456 |
else if( atom[i].mmx_type == 6 || atom[i].mmx_type == 15 || atom[i].mmx_type == |
457 |
42 || atom[i].mmx_type == 34 || atom[i].mmx_type == 35 || atom[i].mmx_type == 53 ) |
458 |
{ |
459 |
atom[ijk[0]].x = 0.2488F; |
460 |
atom[ijk[1]].x = 0.2488F; |
461 |
atom[ijk[0]].y = 0.0F; |
462 |
atom[ijk[1]].y = 0.0F; |
463 |
atom[ijk[0]].z = 0.5460F; |
464 |
atom[ijk[1]].z = -0.5460F; |
465 |
} |
466 |
} |
467 |
} |
468 |
} |
469 |
} |
470 |
} |
471 |
} |
472 |
/* end of nhadds = 2 |
473 |
* **** there are 3 h added to the carbon i *** */ |
474 |
if( nhadds == 3 ) |
475 |
{ |
476 |
for( j1 = 0; j1 < 10; j1++ ) |
477 |
{ |
478 |
if( iadj[j1] != 0 ) |
479 |
{ |
480 |
j2 = iadj[j1]; |
481 |
if( atom[j2].atomnum != 1 && atom[j2].mmx_type != 20 ) |
482 |
na = j2; |
483 |
if( atom[j2].mmx_type == 2 || atom[j2].mmx_type == 3 ) |
484 |
{ |
485 |
for( j3 = 0; j3 < MAXIAT; j3++ ) |
486 |
{ |
487 |
if( (atom[j2].iat[j3] != 0 && atom[j2].iat[j3] != i) && |
488 |
atom[j2].bo[j3] != 9 ) |
489 |
{ |
490 |
j4 = atom[j2].iat[j3]; |
491 |
if( atom[j4].mmx_type == 2 || (atom[j4].mmx_type == 7 && atom[j4].mmx_type == 38) ) |
492 |
{ |
493 |
xpnt = atom[j4].x; |
494 |
ypnt = atom[j4].y; |
495 |
zpnt = atom[j4].z; |
496 |
xyplan( &xpnt, &ypnt, &zpnt ); |
497 |
if( atom[j4].y < 0.F ) |
498 |
{ |
499 |
for( j5 = 1; j5 <= (natom + 4); j5++ ) |
500 |
atom[j5].y = -atom[j5].y; |
501 |
} |
502 |
} |
503 |
} |
504 |
} |
505 |
}else |
506 |
{ |
507 |
/* end of type 2 and type 3 adjacent* */ |
508 |
for( j8 = 0; j8 < MAXIAT; j8++ ) |
509 |
{ |
510 |
if( (atom[j2].iat[j8] != 0 && atom[j2].iat[j8] != i) |
511 |
&& atom[j2].bo[j8] != 9 ) |
512 |
{ |
513 |
j9 = atom[j2].iat[j8]; |
514 |
// if( atom[j9].mmx_type != 5 && atom[j9].mmx_type != 20 ) |
515 |
if(atom[j9].mmx_type != 20 ) |
516 |
{ |
517 |
xpnt = atom[j9].x; |
518 |
ypnt = atom[j9].y; |
519 |
zpnt = atom[j9].z; |
520 |
xyplan( &xpnt, &ypnt, &zpnt ); |
521 |
if( atom[j9].y > 0.F ) |
522 |
{ |
523 |
for( k1 = 1; k1 <= (natom + 4); k1++ ) |
524 |
atom[k1].y = -atom[k1].y; |
525 |
} |
526 |
} |
527 |
} |
528 |
} |
529 |
} |
530 |
/* do type 4 here */ |
531 |
ij = 0; |
532 |
for( k2 = 0; k2 < MAXIAT; k2++ ) |
533 |
{ |
534 |
if( !(atom[i].iat[k2] == 0 || atom[i].bo[k2] == 9) ) |
535 |
{ |
536 |
iii = atom[i].iat[k2]; |
537 |
if( !(atom[iii].mmx_type != 5 && atom[iii].mmx_type != 20) ) |
538 |
{ |
539 |
ijk[ij] = iii; |
540 |
ij++; |
541 |
if( ij == 3 ) |
542 |
{ |
543 |
dis = 0.F; |
544 |
/* *** -ch3 *** |
545 |
* *** add the hydrogens displaced along x enough |
546 |
* to give a normal length c-c bond *** */ |
547 |
atom[ijk[0]].x = 0.3796F + dis; |
548 |
atom[ijk[1]].x = 0.3796F + dis; |
549 |
atom[ijk[2]].x = 0.3796F + dis; |
550 |
atom[ijk[0]].y = 1.043F; |
551 |
atom[ijk[1]].y = -0.5215F; |
552 |
atom[ijk[2]].y = -0.5215F; |
553 |
atom[ijk[0]].z = 0.00F; |
554 |
atom[ijk[1]].z = 0.9032F; |
555 |
atom[ijk[2]].z = -0.9032F; |
556 |
} |
557 |
} |
558 |
} |
559 |
} |
560 |
} |
561 |
} |
562 |
} |
563 |
} |
564 |
} |
565 |
} |
566 |
} |
567 |
/* *** reorient the structure in starting orientation and |
568 |
* turn off atom markers *** */ |
569 |
revec(); |
570 |
return; |
571 |
} |
572 |
// ======================================================== |
573 |
void revec() |
574 |
{ |
575 |
int i; |
576 |
float xpnt, ypnt, zpnt; |
577 |
/* *** this routine used the unit vectors and origin saved |
578 |
* by 'subroutine getvec' to return the molecule to |
579 |
* its original orientation *** */ |
580 |
|
581 |
if( natom <= MAXATOM - 4 ) |
582 |
{ |
583 |
for( i = 1; i <= (natom + 4); i++ ) |
584 |
{ |
585 |
atom[i].x = atom[i].x - atom[natom + 4].x; |
586 |
atom[i].y = atom[i].y - atom[natom + 4].y; |
587 |
atom[i].z = atom[i].z - atom[natom + 4].z; |
588 |
} |
589 |
xpnt = atom[natom + 1].x; |
590 |
ypnt = atom[natom + 1].y; |
591 |
zpnt = atom[natom + 1].z; |
592 |
xaxis( &xpnt, &ypnt, &zpnt ); |
593 |
xpnt = atom[natom + 2].x; |
594 |
ypnt = atom[natom + 2].y; |
595 |
zpnt = atom[natom + 2].z; |
596 |
xyplan( &xpnt, &ypnt, &zpnt ); |
597 |
if( atom[natom + 1].x < 0. ) |
598 |
{ |
599 |
for( i = 1; i <= (natom + 4); i++ ) |
600 |
atom[i].x = -atom[i].x; |
601 |
} |
602 |
if( atom[natom + 2].y < 0. ) |
603 |
{ |
604 |
for( i = 1; i <= (natom + 4); i++ ) |
605 |
atom[i].y = -atom[i].y; |
606 |
} |
607 |
if( atom[natom + 3].z < 0. ) |
608 |
{ |
609 |
for( i = 1; i <= (natom + 4); i++ ) |
610 |
atom[i].z = -atom[i].z; |
611 |
} |
612 |
} |
613 |
return; |
614 |
} |
615 |
/* ------------------------ */ |
616 |
void xyplan(float *xr, float *yr, float *zr) |
617 |
{ |
618 |
int i; |
619 |
float cos3, denom, savex, sin3; |
620 |
|
621 |
/* this subroutine rotates the molecule to place atom (xr,yr,*zr) |
622 |
* in the x,y-plane. note0 y-coordinate will be +. */ |
623 |
|
624 |
denom = sqrt( *yr**yr + *zr**zr ); |
625 |
if( fabs( denom ) >= .00001F ) |
626 |
{ |
627 |
sin3 = *zr/denom; |
628 |
cos3 = *yr/denom; |
629 |
for( i = 1; i <= (natom + 4); i++ ) |
630 |
{ |
631 |
savex = atom[i].y; |
632 |
atom[i].y = atom[i].y*cos3 + atom[i].z*sin3; |
633 |
atom[i].z = atom[i].z*cos3 - savex*sin3; |
634 |
} |
635 |
savex = *yr; |
636 |
*yr = *yr*cos3 + *zr*sin3; |
637 |
*zr = *zr*cos3 - savex*sin3; |
638 |
} |
639 |
return; |
640 |
} |
641 |
/* ------------------------ */ |
642 |
void xaxis(float *xpnt,float *ypnt,float *zpnt) |
643 |
{ |
644 |
int i; |
645 |
float cos1, cos2, denom, save, sin1, sin2; |
646 |
|
647 |
/* this subroutine rotates the molecule about z-axis and then y-axis to |
648 |
* place point with coordinates (i*xpnt,iypnt,izpnt) along x-axis. */ |
649 |
|
650 |
denom = sqrt( *xpnt**xpnt + *ypnt**ypnt ); |
651 |
if( fabs( denom ) < .00001F ) |
652 |
{ |
653 |
sin1 = 0.F; |
654 |
cos1 = 1.F; |
655 |
}else |
656 |
{ |
657 |
sin1 = *ypnt/denom; |
658 |
cos1 = *xpnt/denom; |
659 |
for( i = 1; i <= (natom + 4); i++ ) |
660 |
{ |
661 |
save = atom[i].x; |
662 |
atom[i].x = atom[i].x*cos1 + atom[i].y*sin1; |
663 |
atom[i].y = atom[i].y*cos1 - save*sin1; |
664 |
} |
665 |
save = *xpnt; |
666 |
*xpnt = *xpnt*cos1 + *ypnt*sin1; |
667 |
*ypnt = *ypnt*cos1 - save*sin1; |
668 |
} |
669 |
denom = sqrt( *xpnt**xpnt + *zpnt**zpnt ); |
670 |
if( fabs( denom ) >= .00001F ) |
671 |
{ |
672 |
sin2 = *zpnt/denom; |
673 |
cos2 = *xpnt/denom; |
674 |
for( i = 1; i <= (natom + 4); i++ ) |
675 |
{ |
676 |
save = atom[i].x; |
677 |
atom[i].x = atom[i].x*cos2 + atom[i].z*sin2; |
678 |
atom[i].z = atom[i].z*cos2 - save*sin2; |
679 |
} |
680 |
save = *xpnt; |
681 |
*xpnt = *xpnt*cos2 + *zpnt*sin2; |
682 |
*zpnt = *zpnt*cos2 - save*sin2; |
683 |
} |
684 |
return; |
685 |
} |
686 |
/* ============================================== */ |
687 |
void getvec() |
688 |
{ |
689 |
|
690 |
/* *** this routine saves unit vectors and origin in the |
691 |
* n+1 - n+4 slots of x, y, and z *** */ |
692 |
|
693 |
if(natom > MAXATOM - 4 ) |
694 |
return; |
695 |
atom[natom + 1].x = 1.F; |
696 |
atom[natom + 1].y = 0.F; |
697 |
atom[natom + 1].z = 0.F; |
698 |
atom[natom + 2].x = 0.F; |
699 |
atom[natom + 2].y = 1.F; |
700 |
atom[natom + 2].z = 0.F; |
701 |
atom[natom + 3].x = 0.F; |
702 |
atom[natom + 3].y = 0.F; |
703 |
atom[natom + 3].z = 1.F; |
704 |
atom[natom + 4].x = 0.F; |
705 |
atom[natom + 4].y = 0.F; |
706 |
atom[natom + 4].z = 0.F; |
707 |
return; |
708 |
} |
709 |
/* ============================================== */ |
710 |
void hdel(int lptest) |
711 |
{ |
712 |
/* lptest = 0 to remove H and lone pair */ |
713 |
/* lptest = 1 to remove lone pairs only */ |
714 |
/* lptest = 2 to remove lp & reset various groups for amber */ |
715 |
int i, iatt, iatt1, iatyp, it, it1,ntemp; |
716 |
int j, jatom; |
717 |
/* *** delete hydrogens *** */ |
718 |
ntemp = natom; |
719 |
for( i = 1; i <= ntemp; i++ ) |
720 |
{ |
721 |
if( atom[i].mmx_type != 0 ) |
722 |
{ |
723 |
iatyp = atom[i].mmx_type; |
724 |
if( (lptest == 0 && (iatyp == 5 || iatyp == 20 )) || (lptest >= 1 && iatyp == 20) ) |
725 |
{ |
726 |
iatt = atom[i].iat[0]; |
727 |
iatt1 = atom[i].iat[1]; |
728 |
if( iatt != 0 || iatt1 != 0) |
729 |
{ |
730 |
it = atom[iatt].mmx_type; |
731 |
it1 = atom[iatt1].mmx_type; |
732 |
if (it >= 11 && it <= 14) // halogens |
733 |
goto L_10; |
734 |
if( (it < 300) && (it1 < 300)) // attached metals |
735 |
deleteatom(i); |
736 |
} |
737 |
} |
738 |
} |
739 |
L_10: |
740 |
continue; |
741 |
} |
742 |
reseq(); |
743 |
generate_bonds(); |
744 |
return; |
745 |
} |
746 |
/* -------------------------- */ |
747 |
void reseq() |
748 |
{ |
749 |
int i, ia, ib, iplus, ixx, j, jincr, k; |
750 |
|
751 |
/* ***this subroutine will repack the atom |
752 |
* ***connection table, making sure to fill in |
753 |
* **all "holes" caused by deleting atoms or bonds. |
754 |
* ***the necessary common blocks involve arrays |
755 |
* of the connection table (contb) */ |
756 |
|
757 |
for( i = 1; i < MAXATOM; i++ ) |
758 |
{ |
759 |
if( atom[i].type != 0 ) |
760 |
{ |
761 |
for( ia = 0; ia < (MAXIAT - 1); ia++ ) |
762 |
{ |
763 |
if( atom[i].iat[ia] == 0 ) |
764 |
{ |
765 |
/* *** for every zero entry in the array iat |
766 |
* for atom i, want to move down all remaining |
767 |
* entries of iat by one slot; after doing so |
768 |
* want to store that zero entry in the last |
769 |
* slot--iat *** */ |
770 |
for( ib = ia + 1; ib < MAXIAT; ib++ ) |
771 |
{ |
772 |
atom[i].iat[ib - 1] = atom[i].iat[ib]; |
773 |
atom[i].bo[ib - 1] = atom[i].bo[ib]; |
774 |
} |
775 |
atom[i].iat[MAXIAT - 1] = 0; |
776 |
} |
777 |
} |
778 |
} |
779 |
} |
780 |
|
781 |
/* *** now want to renumber atoms to leave |
782 |
* no "holes" caused by the deleted atoms *** */ |
783 |
jincr = 0; |
784 |
/* i points to first zero entry |
785 |
* kincr points to next real entry */ |
786 |
for( i = 1; i < (MAXATOM - 1); i++ ) |
787 |
{ |
788 |
iplus = 0; |
789 |
if( atom[i].type == 0 ) |
790 |
{ |
791 |
for( jincr = i; jincr < MAXATOM; jincr++ ) |
792 |
{ |
793 |
if( atom[jincr].type != 0 ) |
794 |
{ |
795 |
iplus = jincr; |
796 |
break; |
797 |
} |
798 |
} |
799 |
if (iplus == 0) |
800 |
break; |
801 |
for( j = 0; j < MAXIAT; j++ ) |
802 |
{ |
803 |
atom[i].iat[j] = atom[iplus].iat[j]; |
804 |
atom[i].bo[j] = atom[iplus].bo[j]; |
805 |
} |
806 |
|
807 |
atom[i].x = atom[iplus].x; |
808 |
atom[i].y = atom[iplus].y; |
809 |
atom[i].z = atom[iplus].z; |
810 |
atom[i].type = atom[iplus].type; |
811 |
atom[i].tclass = atom[iplus].tclass; |
812 |
atom[i].type = atom[iplus].type; |
813 |
atom[i].mmx_type = atom[iplus].mmx_type; |
814 |
atom[i].mm3_type = atom[iplus].mm3_type; |
815 |
atom[i].amber_type = atom[iplus].amber_type; |
816 |
atom[i].mmff_type = atom[iplus].mmff_type; |
817 |
atom[i].atomnum = atom[iplus].atomnum; |
818 |
atom[i].formal_charge = atom[iplus].formal_charge; |
819 |
atom[i].atomwt = atom[iplus].atomwt; |
820 |
atom[i].use = atom[iplus].use; |
821 |
strcpy(atom[i].name,atom[iplus].name); |
822 |
atom[i].charge = atom[iplus].charge; |
823 |
atom[i].flags = atom[iplus].flags; |
824 |
atom[i].radius = atom[iplus].radius; |
825 |
atom[i].energy = atom[iplus].energy; |
826 |
/* *** delete the atom whose statistics were just |
827 |
* stored in the i-th atom's connection table *** */ |
828 |
for( ixx = 0; ixx < MAXIAT; ixx++ ) |
829 |
{ |
830 |
atom[iplus].iat[ixx] = 0; |
831 |
atom[iplus].bo[ixx] = 0; |
832 |
} |
833 |
atom[iplus].x = 0.0F; |
834 |
atom[iplus].y = 0.0F; |
835 |
atom[iplus].z = 0.0F; |
836 |
atom[iplus].type = 0; |
837 |
atom[iplus].tclass = 0; |
838 |
atom[iplus].mmx_type = 0; |
839 |
atom[iplus].mm3_type = 0; |
840 |
atom[iplus].amber_type = 0; |
841 |
atom[iplus].mmff_type = 0; |
842 |
atom[iplus].atomnum = 0; |
843 |
atom[iplus].formal_charge = 0; |
844 |
atom[iplus].atomwt = 0.0F; |
845 |
atom[iplus].use = 0; |
846 |
atom[iplus].charge = 0.0F; |
847 |
atom[iplus].flags = 0; |
848 |
atom[iplus].radius = 0.0F; |
849 |
atom[iplus].energy = 0.0F; |
850 |
strcpy(atom[iplus].name,""); |
851 |
|
852 |
for( j = 0; j < MAXIAT; j++ ) |
853 |
{ |
854 |
if( atom[i].iat[j] != 0 ) |
855 |
{ |
856 |
for( k = 0; k < MAXIAT; k++ ) |
857 |
{ |
858 |
if( atom[atom[i].iat[j]].iat[k] == iplus ) |
859 |
atom[atom[i].iat[j]].iat[k] = i; |
860 |
} |
861 |
} |
862 |
} |
863 |
} |
864 |
} |
865 |
/* ***now want to do the same resequencing for iat array */ |
866 |
|
867 |
/* recalc number of atoms */ |
868 |
natom = 0; |
869 |
for (i = 1; i < MAXATOM; i++) |
870 |
{ |
871 |
if (atom[i].type != 0) |
872 |
natom++; |
873 |
} |
874 |
return; |
875 |
} |