Angle between 3 points in 3d space
Asked Answered
G

5

21

I have 3 points containing X, Y, Z coordinates:

var A = {x: 100, y: 100, z: 80},
    B = {x: 100, y: 175, z: 80},
    C = {x: 100, y: 100, z: 120};

The coordinates are pixels from a 3d CSS transform. How can I get the angle between vectors BA and BC? A math formula will do, JavaScript code will be better. Thank you.

enter image description here

Gaseous answered 1/11, 2013 at 15:26 Comment(6)
en.wikipedia.org/wiki/Law_of_cosinesAesculapius
see #4143756Nigel
I know that I have to normalize the vectors but do not know howGaseous
Divide a vector by its length. And it's length is simply the root of its squares.Aesculapius
does the picture points to the right angle?Hewet
how you find the coordinates after 3-d transform ?Lunulate
M
40

In pseudo-code, the vector BA (call it v1) is:

v1 = {A.x - B.x, A.y - B.y, A.z - B.z}

Similarly the vector BC (call it v2) is:

v2 = {C.x - B.x, C.y - B.y, C.z - B.z}

The dot product of v1 and v2 is a function of the cosine of the angle between them (it's scaled by the product of their magnitudes). So first normalize v1 and v2:

v1mag = sqrt(v1.x * v1.x + v1.y * v1.y + v1.z * v1.z)
v1norm = {v1.x / v1mag, v1.y / v1mag, v1.z / v1mag}

v2mag = sqrt(v2.x * v2.x + v2.y * v2.y + v2.z * v2.z)
v2norm = {v2.x / v2mag, v2.y / v2mag, v2.z / v2mag}

Then calculate the dot product:

res = v1norm.x * v2norm.x + v1norm.y * v2norm.y + v1norm.z * v2norm.z

And finally, recover the angle:

angle = acos(res)
Merna answered 1/11, 2013 at 15:41 Comment(16)
Nice, just let me put it in JavascriptGaseous
would you have a jsfiddle of the result @Gaseous ?Saucepan
@Saucepan Sorry no fiddle but look trough typefolly.com/interaction as it should be in thereGaseous
man, isn't there something more effective and less precise? acos and two sqrts in a render loop can build up slow quite fastBronwynbronx
@Bronwynbronx well the question wasn't about speed and of course there are fast sqrt and trig algorithms but why would you ever need to calculate the actual angle inside a render loop?Merna
well to me it seems like it's almost always. you manipulate the objects somehow, then want to perform some logic based on their location. i mean, the possibilities when you need to do trigonometry in the render loop are limitless. btw i'm not trying to hate your algorithm or anything, i was just trying to find something fasterBronwynbronx
@Bronwynbronx Well, two points: 1) there are very fast algorithms in shaders that make many optimizations unnecessary, and 2) trigonometry in this context is just a different way of viewing linear algebra. The algebra itself is almost always faster and vectorizes well - knowledge of an actual angle is mostly a step on the way somewhere else - jump that step with the linear algebra and avoid any need for the computation in the first place. It's always quicker to do nothing ;-)Merna
@Bronwynbronx Don't worry, I'm not offended - you can't really "hate" my algorithm, because it isn't mine anyway, it's just an expression of the fact that the dot product of two unit vectors is simply the cosine of the angle between them. The question was clear and the answer appropriate - if you need something faster in a different context, that's a different question with a different solution. Why not ask it?Merna
Do we really need to normalize the vectors first?Scatology
@SnehasishKarmakar the dot product is the cosine of the angle, scaled by the product of the magnitudes. So either normalise before or after.Merna
@ Roger Rowland: Yeah, basically it is the same thing. It was stupid of me not to carefully notice before bothering you! :)Scatology
@RogerRowland would it work to obtain torsion angles?Ramekin
@Ramekin Well, I had to google for torsion angles and I guess you’re a chemist, so if I’ve understood correctly you’re talking about the angle between two planes, each of which is defined by three points? Although the technique handles three dimensions ok, it might be more convenient to use the alternative plane equation, where the plane is described as a combination of a point on the plane and a normal vector. Then the angle between the planes could be derived just from the two normals.Merna
@Ramekin Note that this method would not recover the sign of the angle, which I see may be important in your application.Merna
@RogerRowland thank you and yes, what you say is true. The application to chemistry can be exemplified through this image cbio.bmt.tue.nl/pumma/uploads/Theory/dihedral.png and what I have done is in this post #46652042Ramekin
@Saucepan Here's a gist of it coded it up in Javascript: gist.github.com/andfaulkner/c4ad12a72d29bcd653eb4b8cca2ae476Friseur
V
6
double GetAngleABC( double* a, double* b, double* c )
{
    double ab[3] = { b[0] - a[0], b[1] - a[1], b[2] - a[2] };
    double bc[3] = { c[0] - b[0], c[1] - b[1], c[2] - b[2]  };

    double abVec = sqrt(ab[0] * ab[0] + ab[1] * ab[1] + ab[2] * ab[2]);
    double bcVec = sqrt(bc[0] * bc[0] + bc[1] * bc[1] + bc[2] * bc[2]);

    double abNorm[3] = {ab[0] / abVec, ab[1] / abVec, ab[2] / abVec};
    double bcNorm[3] = {bc[0] / bcVec, bc[1] / bcVec, bc[2] / bcVec};

    double res = abNorm[0] * bcNorm[0] + abNorm[1] * bcNorm[1] + abNorm[2] * bcNorm[2];

    return acos(res)*180.0/ 3.141592653589793;
}


double a[] = {1, 0, 0};

double b[] = {0, 0, 0};

double c[] = {0, 1, 0};

std::cout<< "The angle of ABC is " << GetAngleABC(a,b,c)<< "º " << std::endl;
Volitant answered 3/12, 2015 at 15:14 Comment(1)
This code is viable. I note a dangerous condition: you assume the length of the double* is always three without checking this. Where the user/coder passes an array of different length, this code will produce an exception, but you are not handling exceptions. It's a nitpicking thing, really, since the code works, however I prefer to use try/catch exception handling with ANY array operation (or similar operation, to include any collection with a length/size such as vector) AND I prefer to confirm that the length/size of the array/container is what the function expects (here being 3). Up-vote.Vancevancleave
R
3

The same in python (with output in degrees):

import numpy as np
import math 
import time

def angle_2p_3d(a, b, c):       

    v1 = np.array([ a[0] - b[0], a[1] - b[1], a[2] - b[2] ])
    v2 = np.array([ c[0] - b[0], c[1] - b[1], c[2] - b[2] ])

    v1mag = np.sqrt([ v1[0] * v1[0] + v1[1] * v1[1] + v1[2] * v1[2] ])
    v1norm = np.array([ v1[0] / v1mag, v1[1] / v1mag, v1[2] / v1mag ])

    v2mag = np.sqrt(v2[0] * v2[0] + v2[1] * v2[1] + v2[2] * v2[2])
    v2norm = np.array([ v2[0] / v2mag, v2[1] / v2mag, v2[2] / v2mag ])
    res = v1norm[0] * v2norm[0] + v1norm[1] * v2norm[1] + v1norm[2] * v2norm[2]
    angle_rad = np.arccos(res)

    return math.degrees(angle_rad)


p1 = np.array([1,0,0])
p2 = np.array([0,0,0])
p3 = np.array([0,0,1])

start = time.time()
angle= angle_2p_3d(p1, p2, p3)
end = time.time()

print("angle: ", angle)
print("elapsed in: ", end - start)

Output:

angle: 90.0

elapsed in: 8.392333984375e-05

Repp answered 8/1, 2020 at 16:53 Comment(0)
K
2

@Roger algorithm in swift

func SCNVector3Angle(start: SCNVector3, mid: SCNVector3, end: SCNVector3) -> Double {
    let v1 = (start - mid)
    let v2 = (end - mid)
    let v1norm = v1.normalized()
    let v2norm = v2.normalized()

    let res = v1norm.x * v2norm.x + v1norm.y * v2norm.y + v1norm.z * v2norm.z
    let angle: Double = Double(GLKMathRadiansToDegrees(acos(res)))
    return angle
}

/**
* Subtracts two SCNVector3 vectors and returns the result as a new SCNVector3.
*/
func - (left: SCNVector3, right: SCNVector3) -> SCNVector3 {
    return SCNVector3Make(left.x - right.x, left.y - right.y, left.z - right.z)
}

extension SCNVector3
{
    /**
    * Returns the length (magnitude) of the vector described by the SCNVector3
    */
    func length() -> Float {
        return sqrtf(x*x + y*y + z*z)
    }

    /**
    * Normalizes the vector described by the SCNVector3 to length 1.0 and returns
    * the result as a new SCNVector3.
    */
    func normalized() -> SCNVector3 {
        return self / length()
    }
}
Kirwan answered 13/10, 2016 at 3:25 Comment(0)
A
0

And now for a MATLAB (R2018a) version:

%% Calculate the Anlge Between 3 Points in 3D Space
clear, clc, close all
format compact

%% Test calculating the angle between 3 points
%
% Find angle (∠BC) between vector AB and vector BC
%
%      *  <--- C (100, 100, 120)
%     /
%    / 
%   /
%  *  <------- A (100, 100, 80)
%  |
%  |
%  |
%  *  <------- B (100, 175, 80)
aa = [100, 100, 80];
bb = [100, 175, 80];
cc = [100, 100, 120];
theta = angle_3points(aa, bb, cc);
fprintf("Angle: %.4f degrees\n", theta);  % Result: Angle: 28.0725 degrees


%% Function definition for calculating the angle between 3 points in 3D space
function theta = angle_3points(a, b, c)
    % THETA = ANGLE_3POINTS(A, B, C)
    %  Function to calculate angle between 3 points
    %  
    %  Inputs:
    %    a: The first point. 1x3 vector.
    %    b: The second point. 1x3 vector.
    %    c: The third point. 1x3 vector.
    %  Outputs:
    %    theta: The angle between vector ab and vector bc in degrees.
    if nargin < 3 || isempty(a) || isempty(b) || isempty(c)
        error 'Parameters `a`, `b`, and `c` are all required inputs.';
    end
    
    % Check the sizes of vectors a, b, and c
    [n, m] = size(a);
    if n ~= 1 || m ~= 3
        error 'Parameter `a` must be a 1x3 vector.';
    end
    [n, m] = size(b);
    if n ~= 1 || m ~= 3
        error 'Parameter `b` must be a 1x3 vector.';
    end
    [n, m] = size(c);
    if n ~= 1 || m ~= 3
        error 'Parameter `c` must be a 1x3 vector.';
    end
    clear n m;

    v1 = a - b;
    v1_mag = sqrt(sum(v1.^2));
    v1_norm = v1 ./ v1_mag;
    
    v2 = c - b;
    v2_mag = sqrt(sum(v2.^2));
    v2_norm = v2 ./ v2_mag;
    
    theta = v1_norm(1)*v2_norm(1) + v1_norm(2)*v2_norm(2) + ...
            v1_norm(3)*v2_norm(3);
    theta = acos(theta) * 180 / pi();
end
Achaemenid answered 19/4 at 11:37 Comment(0)

© 2022 - 2024 — McMap. All rights reserved.