Google Groups no longer supports new Usenet posts or subscriptions. Historical content remains viewable.
Dismiss

How to calculate angle between two 3D vectors with sign

505 views
Skip to first unread message

Dhritarashtra

unread,
May 26, 2009, 3:15:18 AM5/26/09
to
I have two 3D vectors, say, a = (ax ay az) and b = (bx by bz). The easiest way to calculate the angle is:

theta = atan2 ( magnitude ( cross (a, b) ), dot (a, b) )

However, in 3D, i get values only between 0 and pi, because clockwise and counterclockwise directions make no sense in 3D. If i have to get an angle between -pi and pi, the only way is to rotate the coordinate system so that both vectors lie in one plane, say x-y. Is that correct? (if there is any other way, please let me know.)

In other words, i need to rotate the coordinate system so that my new z-axis (say z') points in the direction of cross (a, b). In matlab, how do i rotate the coordinate system in this fashion, and get (ax' ay' 0) and (bx' by' 0)?

Thanks in advance,
Regards,
Arvind

Matt

unread,
May 26, 2009, 3:55:03 AM5/26/09
to
"Dhritarashtra" <marv...@rediffmail.com> wrote in message <gvg4u6$s6r$1...@fred.mathworks.com>...

> I have two 3D vectors, say, a = (ax ay az) and b = (bx by bz). The easiest way to calculate the angle is:
>
> theta = atan2 ( magnitude ( cross (a, b) ), dot (a, b) )
>
> However, in 3D, i get values only between 0 and pi, because clockwise and counterclockwise directions make no sense in 3D. If i have to get an angle between -pi and pi, the only way is to rotate the coordinate system so that both vectors lie in one plane, say x-y. Is that correct? (if there is any other way, please let me know.)
>
> In other words, i need to rotate the coordinate system so that my new z-axis (say z') points in the direction of cross (a, b).
------

You need to come up with a fixed normal vector, n, specifying the plane in in which a and b will always lie. Then, you modify your formula as follows

theta = atan2 ( abs( cross (a, b) ), dot (a, b) ) * sign( dot( cross(a,b), n) )


How you come up with the normal vector depends on details of problem that you haven't yet given us. Clearly a and b are not arbitrary 3D vectors in your case, if there is a distinction to be made between theta and -theta.

Matt

unread,
May 26, 2009, 4:28:01 AM5/26/09
to
"Matt " <x...@whatever.com> wrote in message <gvg78n$242$1...@fred.mathworks.com>...

> "Dhritarashtra" <marv...@rediffmail.com> wrote in message <gvg4u6$s6r$1...@fred.mathworks.com>...
> > I have two 3D vectors, say, a = (ax ay az) and b = (bx by bz). The easiest way to calculate the angle is:
> >
> > theta = atan2 ( magnitude ( cross (a, b) ), dot (a, b) )
> >
> > However, in 3D, i get values only between 0 and pi, because clockwise and counterclockwise directions make no sense in 3D. If i have to get an angle between -pi and pi, the only way is to rotate the coordinate system so that both vectors lie in one plane, say x-y. Is that correct? (if there is any other way, please let me know.)
------

Another possibility is to define the convention that

theta(a,b)=-theta(b,a)

in which case you can use the formula

theta=acos( dot(a,b) /(norm(a)*norm(b) ) *sign(cross(a,b))

Again, though, we would have to know more about the underlying problem to know if this is an acceptable approach.

Dhritarashtra

unread,
May 26, 2009, 4:56:01 AM5/26/09
to
"Matt " <x...@whatever.com> wrote in message <gvg78n$242$1...@fred.mathworks.com>...
> "Dhritarashtra" <marv...@rediffmail.com> wrote in message <gvg4u6$s6r$1...@fred.mathworks.com>...
> > I have two 3D vectors, say, a = (ax ay az) and b = (bx by bz). The easiest way to calculate the angle is:
> >
> > theta = atan2 ( magnitude ( cross (a, b) ), dot (a, b) )
> >
> > However, in 3D, i get values only between 0 and pi, because clockwise and counterclockwise directions make no sense in 3D. If i have to get an angle between -pi and pi, the only way is to rotate the coordinate system so that both vectors lie in one plane, say x-y. Is that correct? (if there is any other way, please let me know.)
> >
> > In other words, i need to rotate the coordinate system so that my new z-axis (say z') points in the direction of cross (a, b).
> ------
>
> You need to come up with a fixed normal vector, n, specifying the plane in in which a and b will always lie. Then, you modify your formula as follows
>
> theta = atan2 ( abs( cross (a, b) ), dot (a, b) ) * sign( dot( cross(a,b), n) )

But the vector in the direction of cross (a, b) obviously defines the plane of a and b, right? However, in that case, cross (a, b) and n will always be parallel, so sign( dot( cross(a,b), n) will be always 1. Or did i misunderstand something?

>
>
> How you come up with the normal vector depends on details of problem that you haven't yet given us. Clearly a and b are not arbitrary 3D vectors in your case, if there is a distinction to be made between theta and -theta.

I am not sure this will help, but i have some DNA molecules with one kink separating two straight sections of each molecule. I am fitting lines to the two straight sections, and want to calculate the angle between them. Roughly in the same plane, the kink can bend the molecule in the clockwise or counterclockwise direction. I want to calculate both the magnitude and direction of this bend.

Thanks for your help.
Regards,
Arvind

Dhritarashtra

unread,
May 26, 2009, 5:14:02 AM5/26/09
to
"Matt " <x...@whatever.com> wrote in message <gvg96h$i5j$1...@fred.mathworks.com>...

> "Matt " <x...@whatever.com> wrote in message <gvg78n$242$1...@fred.mathworks.com>...
> > "Dhritarashtra" <marv...@rediffmail.com> wrote in message <gvg4u6$s6r$1...@fred.mathworks.com>...
> > > I have two 3D vectors, say, a = (ax ay az) and b = (bx by bz). The easiest way to calculate the angle is:
> > >
> > > theta = atan2 ( magnitude ( cross (a, b) ), dot (a, b) )
> > >
> > > However, in 3D, i get values only between 0 and pi, because clockwise and counterclockwise directions make no sense in 3D. If i have to get an angle between -pi and pi, the only way is to rotate the coordinate system so that both vectors lie in one plane, say x-y. Is that correct? (if there is any other way, please let me know.)
> ------
>
> Another possibility is to define the convention that
>
> theta(a,b)=-theta(b,a)

That holds true, for the case i am studying.

>
> in which case you can use the formula
>
> theta=acos( dot(a,b) /(norm(a)*norm(b) ) *sign(cross(a,b))

Sorry, i don't understand. Sign(cross(a,b)) will be a 3x1 matrix. Acos( dot(a,b) /(norm(a)*norm(b) ) is just a number. Multiplying the two will give me a 3X1 matrix. How does that tell me the sign of angle theta?

Thanks and regards,
Arvind

Bruno Luong

unread,
May 26, 2009, 7:22:01 AM5/26/09
to
"Dhritarashtra" <marv...@rediffmail.com> wrote in message <gvg4u6$s6r$1...@fred.mathworks.com>...

In order to define the sign of the angle in a consistent manner in R^3, you MUST restrict the normal vector of a and b in one looking direction (called it forward), and prevent the normal vector to look backward. Let's call the directive forward vector F (you must select it). Consistent angles in R^3 requires to discard all the pairs of vectors that have normal vectors co-linear with F, i.e., rank({a, b, F})<=2. Now you just need to to add to the atan2() angle by PI if the cross product dot(cross(a,b),F) or det([a b F]) < 0. This defines a consistent angles (excepted for (a,b) co linear with F).

It is not possible to define consistent angle for all pairs in general.

Example: if you pick F as (0,0,1).', then you decide to look from above of (x,y) plane and compute the angle from this direction.

(In the mathematical term, a normal vector belongs to the sphere S2, one can define the equivalent relation (E) n1 ~ n2 <=> (n1 = n2) or (n1 = -n2), and take the class of equivalences by (E). Now *you* must select a fixed directive vector called forward F in R^3. For any class N = {n, -n, n in S2} such that <n,F>!=0, we select a representation n such that <n,F> > 0).

Hope it is not confusing.

Bruno

Dhritarashtra

unread,
May 26, 2009, 8:11:01 AM5/26/09
to
"Bruno Luong" <b.l...@fogale.findmycountry> wrote in message <gvgjcp$62d$1...@fred.mathworks.com>...

> "Dhritarashtra" <marv...@rediffmail.com> wrote in message <gvg4u6$s6r$1...@fred.mathworks.com>...
> > I have two 3D vectors, say, a = (ax ay az) and b = (bx by bz). The easiest way to calculate the angle is:
> >
> > theta = atan2 ( magnitude ( cross (a, b) ), dot (a, b) )
> >
> > However, in 3D, i get values only between 0 and pi, because clockwise and counterclockwise directions make no sense in 3D. If i have to get an angle between -pi and pi, the only way is to rotate the coordinate system so that both vectors lie in one plane, say x-y. Is that correct? (if there is any other way, please let me know.)
> >
> > In other words, i need to rotate the coordinate system so that my new z-axis (say z') points in the direction of cross (a, b). In matlab, how do i rotate the coordinate system in this fashion, and get (ax' ay' 0) and (bx' by' 0)?
>
> In order to define the sign of the angle in a consistent manner in R^3, you MUST restrict the normal vector of a and b in one looking direction (called it forward), and prevent the normal vector to look backward. Let's call the directive forward vector F (you must select it). Consistent angles in R^3 requires to discard all the pairs of vectors that have normal vectors co-linear with F, i.e., rank({a, b, F})<=2. Now you just need to to add to the atan2() angle by PI if the cross product dot(cross(a,b),F) or det([a b F]) < 0. This defines a consistent angles (excepted for (a,b) co linear with F).
>
> It is not possible to define consistent angle for all pairs in general.

I think conceptually i understand. You are saying the angle between a and b would be +ve or -ve depending upon whether i look up or look down on the plane containing a and b. Is that correct?

However, it is not clear to me how to go about selecting/quantifying the vector F you suggested. Can you explain what would F be in relation to cross (a, b)? Can i simply select the unit vector in the direction of cross (a, b) as F? But in that case, dot(cross(a,b),F) would always be +ve??? As you see, i am confused about selecting F.

Thanks a lot for your answer.
Regards,
Arvind

Bruno Luong

unread,
May 26, 2009, 8:37:01 AM5/26/09
to
"Dhritarashtra" <marv...@rediffmail.com> wrote in message <gvgm8l$g5d$1...@fred.mathworks.com>...

>
> I think conceptually i understand. You are saying the angle between a and b would be +ve or -ve depending upon whether i look up or look down on the plane containing a and b. Is that correct?

Yes.

>
> However, it is not clear to me how to go about selecting/quantifying the vector F you suggested. Can you explain what would F be in relation to cross (a, b)? Can i simply select the unit vector in the direction of cross (a, b) as F? But in that case, dot(cross(a,b),F) would always be +ve??? As you see, i am confused about selecting F.

Nobody can answer this question at your place. Essentially you have to pick a forward direction that will be likely *never* belong on the same plane the two half DNA a and b. If you can do that, it's all win; if you can't then it screwed, you can never define a consistent angle as wished.

Bruno

Matt

unread,
May 26, 2009, 8:53:01 AM5/26/09
to
"Dhritarashtra" <marv...@rediffmail.com> wrote in message <gvgar1$4d3$1...@fred.mathworks.com>...

> I am not sure this will help, but i have some DNA molecules with one kink separating two straight sections of each molecule. I am fitting lines to the two straight sections, and want to calculate the angle between them. Roughly in the same plane, the kink can bend the molecule in the clockwise or counterclockwise direction. I want to calculate both the magnitude and direction of this bend.

------

In what way does the molecule look different when it bends by +theta as opposed to -theta? Is it equivalent to looking at the molecule from two different sides of the plane containing it? If so, it's not clear why you need this distinction between +/- theta. Conversely, if it is not equivalent, that difference should provide the basis for selecting the normal vector, F that Bruno is talking about.

Matt

unread,
May 26, 2009, 9:01:01 AM5/26/09
to
"Dhritarashtra" <marv...@rediffmail.com> wrote in message <gvgar1$4d3$1...@fred.mathworks.com>...

> "Matt " <x...@whatever.com> wrote in message <gvg78n$242$1...@fred.mathworks.com>...
> > "Dhritarashtra" <marv...@rediffmail.com> wrote in message <gvg4u6$s6r$1...@fred.mathworks.com>...
> > > I have two 3D vectors, say, a = (ax ay az) and b = (bx by bz). The easiest way to calculate the angle is:
> > >
> > > theta = atan2 ( magnitude ( cross (a, b) ), dot (a, b) )
> > >
> > > However, in 3D, i get values only between 0 and pi, because clockwise and counterclockwise directions make no sense in 3D. If i have to get an angle between -pi and pi, the only way is to rotate the coordinate system so that both vectors lie in one plane, say x-y. Is that correct? (if there is any other way, please let me know.)
> > >
> > > In other words, i need to rotate the coordinate system so that my new z-axis (say z') points in the direction of cross (a, b).
> > ------
> >
> > You need to come up with a fixed normal vector, n, specifying the plane in in which a and b will always lie. Then, you modify your formula as follows
> >
> > theta = atan2 ( abs( cross (a, b) ), dot (a, b) ) * sign( dot( cross(a,b), n) )
>
> But the vector in the direction of cross (a, b) obviously defines the plane of a and b, right? However, in that case, cross (a, b) and n will always be parallel, so sign( dot( cross(a,b), n) will be always 1. Or did i misunderstand something?
---------

No, cross(a,b) will be parallel to n but facing possibly in opposite directions, so the sign can also be -1. It is equivalent to det(a,b,n), so it boils down to the approach Bruno was talking about with n=F.

Matt

unread,
May 26, 2009, 10:36:02 AM5/26/09
to
"Dhritarashtra" <marv...@rediffmail.com> wrote in message <gvgbsq$bpp$1...@fred.mathworks.com>...

> > in which case you can use the formula
> >
> > theta=acos( dot(a,b) /(norm(a)*norm(b) ) *sign(cross(a,b))
>
> Sorry, i don't understand. Sign(cross(a,b)) will be a 3x1 matrix. Acos( dot(a,b) /(norm(a)*norm(b) ) is just a number. Multiplying the two will give me a 3X1 matrix. How does that tell me the sign of angle theta?

------

Forget it.

Dhritarashtra

unread,
May 26, 2009, 1:57:01 PM5/26/09
to
"Matt " <x...@whatever.com> wrote in message <gvgond$j90$1...@fred.mathworks.com>...

> "Dhritarashtra" <marv...@rediffmail.com> wrote in message <gvgar1$4d3$1...@fred.mathworks.com>...
>
> > I am not sure this will help, but i have some DNA molecules with one kink separating two straight sections of each molecule. I am fitting lines to the two straight sections, and want to calculate the angle between them. Roughly in the same plane, the kink can bend the molecule in the clockwise or counterclockwise direction. I want to calculate both the magnitude and direction of this bend.
> ------
>
> In what way does the molecule look different when it bends by +theta as opposed to -theta? Is it equivalent to looking at the molecule from two different sides of the plane containing it?

This explanation might be very confusing to visualise for somebody not familiar with DNA structure but i will try. Essentially the straight section of the DNA molecule is like a cylinder. Now if you take any horizontal cross-section of this cylinder, it has two empty volume elements, one a narrower one and the other a wider one, separated by DNA atoms (Imagine a semi-open book held vertically. For a horizontal cross-section, the inner part of the book is the narrower volume element, the rest is the wider element. For those familiar with DNA structure, these are called the minor and the major grooves). I am basically interested in finding if at the kink, the DNA bends towards the narrower volume element and away from the wider element or vice-versa.

The difficulty is: these volume elements are not pointing in the same direction as you move along the radius of the cylinder (the fitted straight line), they rotate about the radius. So if i assume the x-y plane as perpendicular to the radius, the narrow volume element would be pointing towards +ve x at one point, +ve y after some radial distance and so on. So it is very difficult to fix the direction for the entire length and make sure that at the point of kink, say +ve x points in the direction of wider volume element or something like that. This is particularly problematic because i have to do this for many molecules, and in some cases, different sections of the same molecule, and the methodology has to be consistent.

Sorry for the lengthy explanation. I understand you might not be able to suggest a solution to my specific case. But the explanations by you and Bruno have helped a lot in clarifying the problem to myself. Many thanks to both of you. I will see what i can do.

Regards,
Arvind

Matt

unread,
May 26, 2009, 3:16:01 PM5/26/09
to
"Dhritarashtra" <marv...@rediffmail.com> wrote in message <gvhahd$ckc$1...@fred.mathworks.com>...

> The difficulty is: these volume elements are not pointing in the same direction as you move along the radius of the cylinder (the fitted straight line), they rotate about the radius.

------

I don't think you mean the "radius" of the cylinder. I think you mean the cylinder's "long axis". The radius of the cylinder is a line segment lying in the plane of its circular cross sections, connecting the circle's center to its perimeter. The long axis of the cylinder traverses these cross sections.

-----------------


So if i assume the x-y plane as perpendicular to the radius, the narrow volume element would be pointing towards +ve x at one point, +ve y after some radial distance and so on. So it is very difficult to fix the direction for the entire length and make sure that at the point of kink, say +ve x points in the direction of wider volume element or something like that.

--------------------
It doesn't sound like you need to fix the direction for the entire length of the molecule. It sounds like you need to look only at the cross section at the kink and set up an axis which in that particular cross section points from the wide element toward the narrow element.

Dhritarashtra

unread,
May 27, 2009, 1:16:01 AM5/27/09
to
"Matt " <x...@whatever.com> wrote in message <gvhf5h$fqr$1...@fred.mathworks.com>...

> "Dhritarashtra" <marv...@rediffmail.com> wrote in message <gvhahd$ckc$1...@fred.mathworks.com>...
>
> > The difficulty is: these volume elements are not pointing in the same direction as you move along the radius of the cylinder (the fitted straight line), they rotate about the radius.
> ------
>
> I don't think you mean the "radius" of the cylinder. I think you mean the cylinder's "long axis". The radius of the cylinder is a line segment lying in the plane of its circular cross sections, connecting the circle's center to its perimeter. The long axis of the cylinder traverses these cross sections.

Yes, i meant long axis of the cylinder. Sorry for the mistake. When i wrote the mail, it was late night here, and i was half asleep.

>
> -----------------
> So if i assume the x-y plane as perpendicular to the radius, the narrow volume element would be pointing towards +ve x at one point, +ve y after some radial distance and so on. So it is very difficult to fix the direction for the entire length and make sure that at the point of kink, say +ve x points in the direction of wider volume element or something like that.
> --------------------
> It doesn't sound like you need to fix the direction for the entire length of the molecule. It sounds like you need to look only at the cross section at the kink and set up an axis which in that particular cross section points from the wide element toward the narrow element.

Yes. That is what i will have to do.

Thanks and regards,
Arvind

0 new messages