File: phi_from_cncac.c

package info (click to toggle)
garlic 1.6-3
  • links: PTS, VCS
  • area: main
  • in suites: bullseye, buster, sid
  • size: 4,516 kB
  • sloc: ansic: 52,465; makefile: 2,254
file content (134 lines) | stat: -rw-r--r-- 4,493 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
/* Copyright (C) 2000, 2001 Damir Zucic */

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

				phi_from_cncac.c

Purpose:
	Calculate dihedral angle phi, using C[i-1], N[i], CA[i] and C[i].
	Check the distance between C[i-1] and N[i].

Input:
	(1) Pointer to VectorS structure, with C atom coordinates for the
	    previous residue.
	(2) Pointer to VectorS structure, with N atom coordinates.
	(3) Pointer to VectorS structure, with CA atom coordinates.
	(4) Pointer to VectorS structure, with C atom coordinates.
	(5) Pointer to ConfigS structure.

Output:
	Return value.

Return value:
	(1) Dihedral angle phi, on success.
	(2) BADDIHEDANGLE on failure.

Notes:
	(1) This should help you to understand  phi and  psi definitions:

	    ........./................
	    | H * - * N              |
	    |        \    phi1 = 180 |
	    |         * CA           |  residue 1
	    |        /    psi1 = 180 |
	    | O * - * C              |
	    |........\...............|
	    |       N * - * H        |
	    |        /    phi2 = 180 |
	    |    CA *                |  residue 2
	    |        \    psi2 = 180 |
	    |       C * - * O        |
	    |......../...............|

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

#include <stdio.h>

#include <math.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:=================================================*/

void		VectorProduct_ (VectorS *, VectorS *, VectorS *);
double		AbsoluteValue_ (VectorS *);
double		ScalarProduct_ (VectorS *, VectorS *);

/*======calculate phi from C[i-1], N[i], CA[i] and C[i]:=====================*/

double PhiFromCNCAC_ (VectorS *previousC_vectorSP, VectorS *N_vectorSP,
		      VectorS *CA_vectorSP, VectorS *C_vectorSP,
		      ConfigS *configSP)
{
VectorS		N_previousC_vectorS, N_CA_vectorS;
double		bond_length_squared;
VectorS		CA_N_vectorS, CA_C_vectorS;
VectorS		u1S, u2S;
VectorS		v1S, v2S;
double		denom, ratio;
double		alpha, phi;

/* The first pair of auxiliary vectors: */
N_previousC_vectorS.x = previousC_vectorSP->x - N_vectorSP->x;
N_previousC_vectorS.y = previousC_vectorSP->y - N_vectorSP->y;
N_previousC_vectorS.z = previousC_vectorSP->z - N_vectorSP->z;
N_CA_vectorS.x        = CA_vectorSP->x        - N_vectorSP->x;
N_CA_vectorS.y        = CA_vectorSP->y        - N_vectorSP->y;
N_CA_vectorS.z        = CA_vectorSP->z        - N_vectorSP->z;

/* Check the distance between N and previousC: */
bond_length_squared = N_previousC_vectorS.x * N_previousC_vectorS.x +
		      N_previousC_vectorS.y * N_previousC_vectorS.y +
		      N_previousC_vectorS.z * N_previousC_vectorS.z;
if (bond_length_squared < configSP->C_N_min_squared) return BADDIHEDANGLE;
if (bond_length_squared > configSP->C_N_max_squared) return BADDIHEDANGLE;

/* The second pair of auxiliary vectors: */
CA_N_vectorS.x = N_vectorSP->x - CA_vectorSP->x;
CA_N_vectorS.y = N_vectorSP->y - CA_vectorSP->y;
CA_N_vectorS.z = N_vectorSP->z - CA_vectorSP->z;
CA_C_vectorS.x = C_vectorSP->x - CA_vectorSP->x;
CA_C_vectorS.y = C_vectorSP->y - CA_vectorSP->y;
CA_C_vectorS.z = C_vectorSP->z - CA_vectorSP->z;

/* Two vectors perpendicular to N_CA_vectorS,  mutually orthogonal, the */
/* second in the plane defined by N_previousC_vectorS and N_CA_vectorS: */
VectorProduct_ (&u1S, &N_previousC_vectorS, &N_CA_vectorS);
VectorProduct_ (&u2S, &u1S, &N_CA_vectorS);

/* Two vectors perpendicular to  CA_N_vectorS,  mutually orthogonal, */
/* the second in the plane defined by CA_N_vectorS and CA_C_vectorS: */
VectorProduct_ (&v1S, &CA_N_vectorS, &CA_C_vectorS);
VectorProduct_ (&v2S, &CA_N_vectorS, &v1S);

/* Calculate the angle alpha, which will be used to calculate phi: */

/* Avoid division by zero: */
denom = AbsoluteValue_ (&u1S) * AbsoluteValue_ (&v1S);
if (denom == 0.0) return BADDIHEDANGLE;

/* Use the scalar product to calculate the cosine of the angle: */
ratio = ScalarProduct_ (&u1S, &v1S) / denom;

/* Arc cosine is very sensitive to floating point errors: */
if (ratio <= -1.0) alpha = 3.1415927;
else if (ratio >= 1.0) alpha = 0.0;
else alpha = acos (ratio);

/* There are two possible solutions; the right one is resolved here: */
if (ScalarProduct_ (&v2S, &u1S) >= 0) phi = alpha;
else phi = -alpha;

/* Return the angle (in radians): */
return phi;
}

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