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
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.
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.
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
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
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
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
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
> 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.
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.
> > 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.
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
> 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.
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