1 /*****
2 This file is part of the Babel Program
3 Copyright (C) 1992-96 W. Patrick Walters and Matthew T. Stahl
4 All Rights Reserved
5 All Rights Reserved
6 All Rights Reserved
7 All Rights Reserved
8 
9 For more information please contact :
10 
11 babel@mercury.aichem.arizona.edu
12 ----------------------------------------------------------------------------
13 
14 FILE : orient.c
15 AUTHOR(S) : Pat Walters
16 DATE : 9-1-94
17 PURPOSE : routine to orient a structure along the x and y axes
18 ******/
19 
20 
21 #include "bbltyp.h"
22 #define DELIMS " ,\t\n"
23 #define Primary    1
24 #define Secondary  2
25 #define Tertiary   3
26 
set_std_orientation(ums_type * mol)27 void set_std_orientation(ums_type *mol)
28 {
29   center_at_atom(mol,2);
30   orient_ums(mol,1,2,1,3);
31 }
32 
33 
get_orientation_matrix(ums_type * mol,int x1,int x2,int y1,int y2,matrix_3x3 * m)34 void get_orientation_matrix(ums_type *mol, int x1, int x2, int y1, int y2, matrix_3x3 *m)
35 {
36   vect_type r1,r2,v1,v2,v3;
37 
38   pts_2_vect(mol,&r1,x1,x2);
39   pts_2_vect(mol,&r2,y1,y2);
40   cross_prod(&r1,&r2,&v3);
41   normalize_vect(&v3);
42   cross_prod(&v3,&r1,&v2);
43   normalize_vect(&v2);
44   cross_prod(&v2,&v3,&v1);
45   normalize_vect(&v1);
46 
47   m->a1 = v1.x; m->b1 = v1.y; m->c1 = v1.z;
48   m->a2 = v2.x; m->b2 = v2.y; m->c2 = v2.z;
49   m->a3 = v3.x; m->b3 = v3.y; m->c3 = v3.z;
50 }
51 
ums_plus_vector(ums_type * mol,vect_type * v)52 void ums_plus_vector(ums_type *mol, vect_type *v)
53 {
54   int i;
55 
56   for (i = 1; i <= Atoms; i++)
57   {
58     X(i) += v->x;
59     Y(i) += v->y;
60     Z(i) += v->z;
61   }
62 }
63 
64 
65 
ums_dot_matrix(ums_type * mol,matrix_3x3 * m)66 void ums_dot_matrix(ums_type *mol, matrix_3x3 *m)
67 {
68   int i;
69   vect_type temp;
70 
71   for (i = 1; i <= Atoms; i++)
72   {
73     coord_to_vector(&temp,&Point(i));
74     mat_3x3_dot_vect(m,&temp);
75     vector_to_coord(&Point(i),&temp);
76   }
77 }
78 
79 
orient_ums(ums_type * mol,int x1,int x2,int y1,int y2)80 void orient_ums(ums_type *mol, int x1, int x2, int y1, int y2)
81 {
82   int i;
83   vect_type r1,r2,v1,v2,v3, temp;
84   matrix_3x3 m;
85 
86   pts_2_vect(mol,&r1,x1,x2);
87   pts_2_vect(mol,&r2,y1,y2);
88   cross_prod(&r1,&r2,&v3);
89   normalize_vect(&v3);
90   cross_prod(&v3,&r1,&v2);
91   normalize_vect(&v2);
92   cross_prod(&v2,&v3,&v1);
93   normalize_vect(&v1);
94 
95   m.a1 = v1.x; m.b1 = v1.y; m.c1 = v1.z;
96   m.a2 = v2.x; m.b2 = v2.y; m.c2 = v2.z;
97   m.a3 = v3.x; m.b3 = v3.y; m.c3 = v3.z;
98 
99   for (i = 1; i <= Atoms; i++)
100   {
101     coord_to_vector(&temp,&Point(i));
102     mat_3x3_dot_vect(&m,&temp);
103     vector_to_coord(&Point(i),&temp);
104   }
105 }
106 
coord_to_vector(vect_type * v,coord_type * c)107 void coord_to_vector(vect_type *v, coord_type *c)
108 {
109   v->x = c->x;
110   v->y = c->y;
111   v->z = c->z;
112 }
113 
vector_to_coord(coord_type * c,vect_type * v)114 void vector_to_coord(coord_type *c, vect_type *v)
115 {
116   c->x = v->x;
117   c->y = v->y;
118   c->z = v->z;
119 }
120 
center_at_atom(ums_type * mol,int atom)121 void center_at_atom(ums_type *mol, int atom)
122 {
123   int i;
124   double tempX, tempY, tempZ;
125 
126   tempX = X(atom);
127   tempY = Y(atom);
128   tempZ = Z(atom);
129 
130   for (i = 1; i <= Atoms; i++)
131   {
132     X(i) -= tempX;
133     Y(i) -= tempY;
134     Z(i) -= tempZ;
135   }
136 }
137 
center_at_atoms(ums_type * mol,int atom[],int count)138 void center_at_atoms(ums_type *mol,int atom[],int count)
139 {
140   int i;
141   double tempX, tempY, tempZ;
142 
143   tempX = tempY = tempZ = 0.0;
144   for (i = 0;i < count;i++)
145   {
146     tempX += X(atom[i]);
147     tempY += Y(atom[i]);
148     tempZ += Z(atom[i]);
149   }
150 
151   tempX /= (double)count;
152   tempY /= (double)count;
153   tempZ /= (double)count;
154 
155   for (i = 1; i <= Atoms; i++)
156   {
157     X(i) -= tempX;
158     Y(i) -= tempY;
159     Z(i) -= tempZ;
160   }
161 }
162 
center_at_origin(ums_type * mol,vect_type * v)163 void center_at_origin(ums_type *mol, vect_type *v)
164 {
165   int i;
166   double sumX = 0.0, sumY = 0.0, sumZ = 0.0;
167 
168   for (i = 1; i <= Atoms; i++)
169   {
170     sumX += X(i);
171     sumY += Y(i);
172     sumZ += Z(i);
173   }
174 
175   sumX /= (double)Atoms;
176   sumY /= (double)Atoms;
177   sumZ /= (double)Atoms;
178 
179   for (i = 1; i <= Atoms; i++)
180   {
181     X(i) -= sumX;
182     Y(i) -= sumY;
183     Z(i) -= sumZ;
184   }
185 
186   v->x = sumX;
187   v->y = sumY;
188   v->z = sumZ;
189 }
190 
AlignMol(ums_type * mol)191 void AlignMol(ums_type *mol)
192 {
193   int i;
194   int *cList,*pList,*sList;
195   int cCount,primCount,secCount;
196   vect_type r1,r2,v1,v2,v3;
197   int x,y,z;
198   double primary,secondary,tertiary;
199   char buffer[100];
200   char axis;
201 
202   fprintf(stdout,"Input the atom(s) to be placed at the origin\n");
203   fgets(buffer,sizeof(buffer),stdin);
204   cCount = count_tokens(buffer,DELIMS);
205   if (cCount)
206   {
207     cList = (int *)malloc(sizeof(int) * cCount);
208     if (!cList)
209       fatal_error("Unable to allocate memory for array");
210 
211     for (i = 0;i < cCount;i++)
212       cList[i] = atoi(gettoken(buffer,DELIMS,i+1));
213   }
214   else
215     fatal_error("Atoms to be placed at the origin must be defined");
216 
217   x = y = z = 0;
218 
219   fprintf(stdout,"What is the primary alignment axis (x, y, or z)?\n");
220   fgets(buffer,sizeof(buffer),stdin);
221   lowercase(buffer);
222   axis = '\0';
223   sscanf(buffer,"%c",&axis);
224 
225   if (!axis)
226     fatal_error("You must input x, y, or z!");
227 
228   switch (axis)
229   {
230   case 'x':
231     x = Primary;
232     break;
233   case 'y':
234     y = Primary;
235     break;
236   case 'z':
237     z = Primary;
238     break;
239  }
240 
241   fprintf(stdout,"Enter the atoms to be aligned along the primary axis\n");
242   fgets(buffer,sizeof(buffer),stdin);
243   primCount = count_tokens(buffer,DELIMS);
244   pList = (int *)malloc(sizeof(int) * primCount);
245   if (!pList)
246     fatal_error("Unable to allocate memory for array");
247 
248   for (i = 0;i < primCount;i++)
249     pList[i] = atoi(gettoken(buffer,DELIMS,i+1));
250 
251 
252   fprintf(stdout,"What is the secondary alignment axis (x, y, or z)?\n");
253   fgets(buffer,sizeof(buffer),stdin);
254   lowercase(buffer);
255   axis = '\0';
256   sscanf(buffer,"%c",&axis);
257 
258   if (!axis)
259     fatal_error("You must input x, y, or z!");
260 
261   switch (axis)
262   {
263   case 'x':
264     if (x)
265       fatal_error("X has already been defined");
266     x = Secondary;
267     break;
268   case 'y':
269     if (y)
270       fatal_error("Y has already been defined");
271     y = Secondary;
272     break;
273   case 'z':
274     if (z)
275       fatal_error("Z has already been defined");
276     z = Secondary;
277     break;
278  }
279 
280   fprintf(stdout,"Enter the atoms to be aligned along the secondary axis\n");
281   fgets(buffer,sizeof(buffer),stdin);
282   secCount = count_tokens(buffer,DELIMS);
283   sList = (int *)malloc(sizeof(int) * primCount);
284   if (!sList)
285     fatal_error("Unable to allocate memory for array");
286 
287   for (i = 0;i < secCount;i++)
288     sList[i] = atoi(gettoken(buffer,DELIMS,i+1));
289 
290 
291   if (!(x == Primary || y == Primary || z == Primary))
292     fatal_error("A primary axis must be defined");
293 
294   if (!(x == Secondary || y == Secondary || z == Secondary))
295     fatal_error("A secondary axis must be defined");
296 
297   if (!x)
298     x = Tertiary;
299   if (!y)
300     y = Tertiary;
301   if (!z)
302     z = Tertiary;
303 
304 /*  pts_2_vect(mol,&r1,x1,x2);
305   pts_2_vect(mol,&r2,y1,y2);*/
306 
307   center_at_atoms(mol,cList,cCount);
308 
309   r1.x = r1.y = r1.z = 0;
310   for (i = 0;i < primCount;i++)
311   {
312     r1.x += X(pList[i]);
313     r1.y += Y(pList[i]);
314     r1.z += Z(pList[i]);
315   }
316   r1.x /= (double) primCount;
317   r1.y /= (double) primCount;
318   r1.z /= (double) primCount;
319 
320   r2.x = r2.y = r2.z = 0;
321   for (i = 0;i < secCount;i++)
322   {
323     r2.x += X(sList[i]);
324     r2.y += Y(sList[i]);
325     r2.z += Z(sList[i]);
326   }
327   r2.x /= (double) secCount;
328   r2.y /= (double) secCount;
329   r2.z /= (double) secCount;
330 
331 
332   cross_prod(&r1,&r2,&v3);
333   normalize_vect(&v3);
334   cross_prod(&v3,&r1,&v2);
335   normalize_vect(&v2);
336   cross_prod(&v2,&v3,&v1);
337   normalize_vect(&v1);
338 
339   for (i = 1; i <= Atoms; i++)
340   {
341     primary = v1.x * X(i) + v1.y * Y(i) + v1.z * Z(i);
342     secondary = v2.x * X(i) + v2.y * Y(i) + v2.z * Z(i);
343     tertiary = v3.x * X(i) + v3.y * Y(i) + v3.z * Z(i);
344 
345     switch (x)
346     {
347     case Primary:
348       X(i) = primary;
349       break;
350     case Secondary:
351       X(i) = secondary;
352       break;
353     case Tertiary:
354       X(i) = tertiary;
355       break;
356     }
357 
358     switch (y)
359     {
360     case Primary:
361       Y(i) = primary;
362       break;
363     case Secondary:
364       Y(i) = secondary;
365       break;
366     case Tertiary:
367       Y(i) = tertiary;
368       break;
369     }
370 
371     switch (z)
372     {
373     case Primary:
374       Z(i) = primary;
375       break;
376     case Secondary:
377       Z(i) = secondary;
378       break;
379     case Tertiary:
380       Z(i) = tertiary;
381       break;
382     }
383   }
384 }
385 
find_mol_center(ums_type * mol,vect_type * center,int use_wts)386 void find_mol_center(ums_type *mol, vect_type *center, int use_wts)
387 {
388   int i;
389   double wnorm = 0.0;
390   double maxX = -99e99, maxY = -99e99, maxZ = -99e99;
391   double minX =  99e99, minY =  99e99, minZ =  99e99;
392 
393   center->x = 0.0;
394   center->y = 0.0;
395   center->z = 0.0;
396 
397   if (use_wts)
398     {
399       get_atomic_weights(mol);
400     }
401   else
402     {
403       for (i = 1; i <= Atoms; i++)
404 	{
405 	  Double(i) = 1.0;
406 	}
407     }
408 
409   for (i = 1; i <= Atoms; i++)
410     {
411       center->x += X(i) * Double(i);
412       center->y += Y(i) * Double(i);
413       center->z += Z(i) * Double(i);
414       wnorm += Double(i);
415     }
416 
417   center->x /= wnorm;
418   center->y /= wnorm;
419   center->z /= wnorm;
420 }
421 
422 
423 
424 
425 
426 
427 
428 
429 
430