function [d, a, t] = interncoord(r)

% [d, a, t] = interncoord(r). Compute internal coordinates of a chain.
%
% r: [n 3]: protein chain
% d: [n-1, 1]: bond lengths
% a: [n-2, 1]: bond angles
% t: [n-3, 1]: torsions

dr = diff(r);
d = sqrt(sum(dr.^2, 2));
dr = normr(dr);
dp = dot(-dr(1:end-1,:), dr(2:end,:), 2); % note minus sign
a = acos(dp);
e1 = normr(cross(dr(1:end-2,:), dr(2:end-1,:)));
e2 = normr(cross(dr(2:end-1,:), dr(3:end,:)));
t = acos(dot(e1, e2, 2));
signsin = sign(dot(cross(e1, e2, 2), dr(2:end-1,:), 2));
t(signsin < 0) = 2 * pi - t(signsin < 0);