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