File: acceptor_dipole.c

package info (click to toggle)
garlic 1.5-2
  • links: PTS
  • area: main
  • in suites: etch, etch-m68k
  • size: 4,324 kB
  • ctags: 1,378
  • sloc: ansic: 50,306; makefile: 1,088
file content (271 lines) | stat: -rw-r--r-- 9,254 bytes parent folder | download | duplicates (5)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
/* Copyright (C) 2002 Damir Zucic */

/*=============================================================================

				acceptor_dipole.c

Purpose:
	Prepare two vectors which define the hydrogen bond acceptor dipole.
	The first vector is the position of  the dipole center.  The second
	vector is the unit vector  which defines the dipole direction.  Two
	types of acceptor dipoles are recognized by garlic:  C-O dipole and
	C-N-C dipole (which involves the NE2 atom of HIS side chain). 

Input:
	(1) Pointer to VectorS structure,  where the position of the dipole
	    center will be stored.
	(2) Pointer to  VectorS,  where  the unit vector  which defines the
	    direction of dipole will be stored.
	(3) Pointer to AtomS structure,  with data about  the hydrogen bond
	    acceptor.
	(4) Pointer to MolComplexS structure, which contains the data about
	    the macromolecular complex to which  the hydrogen bond acceptor
	    belongs.
	(5) The  array  index  of the  macromolecular complex  to which the
	    hydrogen bond acceptor belongs.

Output:
	(1) The vector which defines the dipole center.
	(2) The unit vector which defines the dipole direction.
	(3) Return value.

Return value:
	(1) Positive on success.
	(2) Negative on failure.

========includes:============================================================*/

#include <stdio.h>

#include <string.h>

#include <X11/Xlib.h>
#include <X11/Xutil.h>
#include <X11/Xos.h>
#include <X11/Xatom.h>

#include "defines.h"
#include "typedefs.h"

/*======function prototypes:=================================================*/

double		AbsoluteValue_ (VectorS *);

/*======prepare vectors which define the acceptor dipole:====================*/

int AcceptorDipole_ (VectorS *center_vectorSP, VectorS *unit_vectorSP,
		     AtomS *acceptor_atomSP,
		     MolComplexS *current_mol_complexSP, int mol_complexI)
{
ResidueS	*current_residueSP;
size_t		first_atomI, last_atomI;
RawAtomS	*acceptor_raw_atomSP;
int		bondI;
TrueBondS	*current_bondSP;
char		alt_location;
AtomS		*partner_atomSP;
double		x1, y1, z1, x2, y2, z2;
double		absolute_value, reciprocal_abs_value;
int		carbon_atomsN;

/* Prepare  the pointer to  the residue to which  the hydrogen */
/* bond acceptor belongs. Only this residue should be scanned: */
current_residueSP = current_mol_complexSP->residueSP +
		    acceptor_atomSP->residue_arrayI;

/* The range of atomic array indices which define this residue: */
first_atomI = current_residueSP->residue_startI;
last_atomI  = current_residueSP->residue_endI;

/* Prepare the pointer to raw atomic data about the acceptor atom: */
acceptor_raw_atomSP = &acceptor_atomSP->raw_atomS;

/*---------------------------------------------------------------------------*/

/* If the acceptor atom is the oxygen atom,  find the only */
/* carbon bound to this oxygen.  It is assumed  that there */
/* is only one carbon atom bound to the given oxygen atom. */
/* This should be true  for standard  amino acid residues. */
if (strcmp (acceptor_raw_atomSP->chemical_symbolA, " O") == 0)
	{
	/* Check all bonds of the acceptor atom  to find the */
	/* partner which forms the dipole with the acceptor: */
	for (bondI = 0; bondI < acceptor_atomSP->bondsN; bondI++)
		{
		/* Pointer to the current bond: */
		current_bondSP = acceptor_atomSP->true_bondSA + bondI;

		/* Only covalently bound partners are valid candidates: */
		if (current_bondSP->bond_typeI != 1) continue;

		/* The bond partner should belong to the same complex: */
		if (current_bondSP->neighbor_mol_complexI != mol_complexI)
			{
			continue;
			}

		/* The bond partner should belong to the same residue: */
		if (current_bondSP->neighbor_arrayI < first_atomI) continue;
		if (current_bondSP->neighbor_arrayI > last_atomI)  continue;

		/* The bond partner  should be the carbon */
		/* atom (see check_dist.c for pair ID's): */
		if (current_bondSP->pairID != 3) continue;

		/* Prepare the pointer to the bond partner (carbon atom): */
		partner_atomSP = current_mol_complexSP->atomSP +
				 current_bondSP->neighbor_arrayI;

		/* Ignore atoms at alternate positions: */
		alt_location = partner_atomSP->raw_atomS.alt_location;
		if ((alt_location != ' ') && (alt_location != 'A')) continue;

		/* If  this point  is reached,  the carbon */
		/* atom bound to the oxygen atom is found. */

		/* Copy atomic coordinates: */
		x1 = partner_atomSP->raw_atomS.x[0];
		y1 = partner_atomSP->raw_atomS.y;
		z1 = partner_atomSP->raw_atomS.z[0];
		x2 = acceptor_atomSP->raw_atomS.x[0];
		y2 = acceptor_atomSP->raw_atomS.y;
		z2 = acceptor_atomSP->raw_atomS.z[0];

		/* Prepare and store the position of dipole center: */
		center_vectorSP->x = 0.5 * (x1 + x2);
		center_vectorSP->y = 0.5 * (y1 + y2);
		center_vectorSP->z = 0.5 * (z1 + z2);

		/* Prepare and store  the unit vector, */
		/* pointing from oxygen to carbon atom */
		/* (don't forget  the normalization!): */
		unit_vectorSP->x = x1 - x2;
		unit_vectorSP->y = y1 - y2;
		unit_vectorSP->z = z1 - z2;
		absolute_value = AbsoluteValue_ (unit_vectorSP);
		if (absolute_value <= 0.0) return -1;
		reciprocal_abs_value = 1.0 / absolute_value;
		unit_vectorSP->x *= reciprocal_abs_value;
		unit_vectorSP->y *= reciprocal_abs_value;
		unit_vectorSP->z *= reciprocal_abs_value;

		/* Return the success indicator: */
		return 1;
		}

	/* If this point  is reached,  all bond partners  were checked */
	/* but no carbon atom was found. Return the failure indicator: */
	return -2;
	}

/*---------------------------------------------------------------------------*/

/* If the acceptor atom is NE2 (HIS side chain), two carbon atoms */
/* bound to the nitrogen should be found. The point in the middle */
/* of these two atoms  will be combined with nitrogen coordinates */
/* to prepare the position of dipole center and  the unit vector. */
/* It is  expected that there are  two carbon atoms  bound to the */
/* nitrogen atom.  This should be true for  standard HIS residue. */
else if (strcmp (acceptor_raw_atomSP->pure_atom_nameA, "NE2") == 0)
	{
	/* Initialize the position of the point */
	/* which lays between two carbon atoms: */
	x1 = 0.0;
	y1 = 0.0;
	z1 = 0.0;

	/* Initialize the number of carbon atoms bound to nitrogen: */
	carbon_atomsN = 0;

	/* Check all bonds of the acceptor atom  to find two */
	/* partners which form the dipole with the acceptor: */
	for (bondI = 0; bondI < acceptor_atomSP->bondsN; bondI++)
		{
		/* Pointer to the current bond: */
		current_bondSP = acceptor_atomSP->true_bondSA + bondI;

		/* Only covalently bound partners are valid candidates: */
		if (current_bondSP->bond_typeI != 1) continue;

		/* The bond partner should belong to the same complex: */
		if (current_bondSP->neighbor_mol_complexI != mol_complexI)
			{
			continue;
			}

		/* The bond partner should belong to the same residue: */
		if (current_bondSP->neighbor_arrayI < first_atomI) continue;
		if (current_bondSP->neighbor_arrayI > last_atomI)  continue;

		/* The bond partner  should be the carbon */
		/* atom (see check_dist.c for pair ID's): */
		if (current_bondSP->pairID != 2) continue;

		/* Prepare the pointer to the bond partner (carbon atom): */
		partner_atomSP = current_mol_complexSP->atomSP +
				 current_bondSP->neighbor_arrayI;

		/* Ignore atoms at alternate positions: */
		alt_location = partner_atomSP->raw_atomS.alt_location;
		if ((alt_location != ' ') && (alt_location != 'A')) continue;

		/* If this point is reached, one of the carbon */
		/* atoms bound to  the nitrogen atom is found. */

		/* Add  the coordinates  of this  carbon atom */
		/* to the coordinates of the point which lays */
		/* in  the middle  between  two carbon atoms: */
		x1 += partner_atomSP->raw_atomS.x[0];
		y1 += partner_atomSP->raw_atomS.y;
		z1 += partner_atomSP->raw_atomS.z[0];

		/* Increment  the counter  which counts the */
		/* carbon atoms bound to the nitrogen atom: */
		carbon_atomsN++;
		}

	/* There should be two carbon atoms bound to the nitrogen atom: */
	if (carbon_atomsN != 2) return -3;

	/* Fix the position of the first point: */
	x1 *= 0.5;
	y1 *= 0.5;
	z1 *= 0.5;

	/* Copy the coordinates of the second point: */
	x2 = acceptor_atomSP->raw_atomS.x[0];
	y2 = acceptor_atomSP->raw_atomS.y;
	z2 = acceptor_atomSP->raw_atomS.z[0];

	/* Prepare and store the position of dipole center: */
	center_vectorSP->x = 0.5 * (x1 + x2);
	center_vectorSP->y = 0.5 * (y1 + y2);
	center_vectorSP->z = 0.5 * (z1 + z2);

	/* Prepare  and store  the unit vector */
	/* pointing from carbon to oxygen atom */
	/* (don't forget  the normalization!): */
	unit_vectorSP->x = x2 - x1;
	unit_vectorSP->y = y2 - y1;
	unit_vectorSP->z = z2 - z1;
	absolute_value = AbsoluteValue_ (unit_vectorSP);
	if (absolute_value <= 0.0) return -4;
	reciprocal_abs_value = 1.0 / absolute_value;
	unit_vectorSP->x *= reciprocal_abs_value;
	unit_vectorSP->y *= reciprocal_abs_value;
	unit_vectorSP->z *= reciprocal_abs_value;

	/* Return the success indicator: */
	return 2;
	}

/*---------------------------------------------------------------------------*/

/* If this point is reached, the acceptor atom was */
/* not recognized - return  the failure indicator: */
return -5;
}

/*===========================================================================*/