Given the x-, y-, and z-coordinates of four points (a-b-c-d) in 3-dimensional (3D) space, how to calculate the torsion angle? Overall, this is a well-solved problem in structural biology and chemistry; one can find a description of torsion angle in many text books and on-line documents. The algorithm for its calculation is implementated in virtually every software package in computational structural biology and chemistry.
As basic as the concept is, however, it is important (based on my experience) to have a clear understanding of how torsion angle is defined in order to really get into the 3D world. Here is a worked example using Octave/Matlab of my simplified, geometry-based implementation of calculating torsion angle, including how to determine its sign. No theory or (complicated) mathematical formula, just a step-by-step illustration of how I solve this problem.
- Coordinates of four points are given in variable
abcd = [ 21.350 31.325 22.681 22.409 31.286 21.483 22.840 29.751 21.498 23.543 29.175 22.594 ];
- Two auxiliary functions: norm_vec() to normalize a vector; get_orth_norm_vec() to get the orthogonal component (normalized) of a vector with reference to another vector, which should have already been normalized.
function ovec = norm_vec(vec) ovec = vec / norm(vec); endfunction
function ovec = get_orth_norm_vec(vec, vref) temp = vec - vref * dot(vec, vref); ovec = norm_vec(temp); endfunction
- Get three vectors: b_c is the normalized vector b→c; b_a_orth is the orthogonal component (normalized) of vector b→a with reference to b→c; c_d_orth is similarly defined, as the orthogonal component (normalized) of vector c→d with reference to b→c.
b_c = norm_vec(abcd(3, :) - abcd(2, :)) % [0.2703158 -0.9627257 0.0094077] b_a_orth = get_orth_norm_vec(abcd(1, :) - abcd(2, :), b_c) % [-0.62126 -0.16696 0.76561] c_d_orth = get_orth_norm_vec(abcd(4, :) - abcd(3, :), b_c) % [0.41330 0.12486 0.90199]
- Now the torsion angle is defined as the angle between the two vectors, b_a_orth and c_d_orth, and can be easily calculated by their dot product. The sign of the torsion angle is determined by the relative orientation of the cross product of the same two vectors with reference to the middle vector b→c. Here they are in opposite direction, thus the torsion angle is negative.
angle_deg = acos(dot(b_a_orth, c_d_orth)) * 180 / pi % 65.609 sign = dot(cross(b_a_orth, c_d_orth), b_c) % -0.91075 if (sign < 0) ang_deg = -angle_deg % -65.609 endif
A related concept is the so-called dihedral angle, or more generally the angle between two planes. As long as the normal vectors to the two corresponding planes are defined, the angle between them is easy to work out.
It’s worth noting that the helical twist angle in SCHNAaP and 3DNA is calculated similarly.