5

I have two matrices consisting of 3d vectors (numpy 1D arrays) and I need to calculate the angle between the vectors, row-wise, and return the results in a 1d array. I know how to calculate the angle between two 1d vectors. What is the the proper way to do this?

*** The resulting angles are in degrees not radians.

By now I have this:

import numpy as np

A = np.array([[1,0,0],
              [0,1,0],
              [0,0,1]])

B = np.array([[1,0,1],
              [1,1,0],
              [0,1,0]])

def angle(V1,V2):
    """
    angle between vectors V1 and V2 in degrees using
    angle = arccos ( V1 dot V2 / norm(V1) * norm(V2) ) *180/np.pi
    """

    cos_of_angle = V1.dot(V2) / (np.linalg.norm(V1) * np.linalg.norm(V2)) 
    return np.arccos(np.clip(cos_of_angle,-1,1))  * 180/np.pi

Note the scaling term 180/np.pi for the conversion from rad to deg.

I would like to have an array:

C = [ angle(A[0],B[0]) , angle(A[1],B[1])...... and so on]

Really appreciated if someone could help.

RRS
  • 75
  • 1
  • 5

2 Answers2

6

We could use einsum to replace the dot-product computations and axis param for the norm ones to have a vectorized solution, like so -

def angle_rowwise(A, B):
    p1 = np.einsum('ij,ij->i',A,B)
    p2 = np.linalg.norm(A,axis=1)
    p3 = np.linalg.norm(B,axis=1)
    p4 = p1 / (p2*p3)
    return np.arccos(np.clip(p4,-1.0,1.0))

We could optimize further and bring in more of einsum, specifically to compute norms with it. Hence, we could use it like so -

def angle_rowwise_v2(A, B):
    p1 = np.einsum('ij,ij->i',A,B)
    p2 = np.einsum('ij,ij->i',A,A)
    p3 = np.einsum('ij,ij->i',B,B)
    p4 = p1 / np.sqrt(p2*p3)
    return np.arccos(np.clip(p4,-1.0,1.0))

Hence, to solve our case to get output in degrees - 

out = angle_rowwise(A, B)*180/np.pi
Divakar
  • 218,885
  • 19
  • 262
  • 358
  • are these angles relative to the origin (0, 0, 0)? I am trying to translate the idea to using 3D points. – NaN Jun 10 '18 at 02:00
  • @NaN I am not really sure about the theory there. I was just focusing on translating it to a vectorized code. – Divakar Jun 10 '18 at 06:42
  • it must be something different since you get different results than the implementation here https://stackoverflow.com/questions/2827393/angles-between-two-n-dimensional-vectors-in-python a = np.array([[1, 0, 0], [1, 0, 0], [1, 0, 0]]) and b = np.array([[1, 0, 0], [0, 1, 0], [0, -1, 0]]) to represent identical, same and opposite directional vectors pairing – NaN Jun 10 '18 at 07:53
  • @NaN That implementation looks different than the one posted by OP here. I am following this one here. – Divakar Jun 10 '18 at 08:20
  • Dear friends, thank you for the kind help.I edited the original post, as I have noticed I forgot to include the np.arccos in the return statement of the angle function. And as pointed by @Nan, I have included the correction suggested at (https://stackoverflow.com/questions/2827393/angles-between-two-n-dimensional-vectors-in-python) where a Nan would be returned in case of parallel vectors. – RRS Jun 10 '18 at 10:51
  • @RaphaelSantos So, are you sure about the `* 180/np.pi` that you are proposing in this question? I don't see that in the linked question. – Divakar Jun 10 '18 at 11:03
  • @Divakar the original post returned radians, the 180/np.pi is to convert to degrees... preference for explicit reading is simply np.degrees which does the same – NaN Jun 10 '18 at 11:12
  • @RaphaelSantos Edited answer. – Divakar Jun 10 '18 at 11:14
  • @NaN and Divakar I can include np.degrees in the original post if you think it makes it easier for others to understand, otherwise I'll call it a wrap for now. – RRS Jun 10 '18 at 11:29
  • @RaphaelSantos Thin, you can just add in the question that the final output required is in degrees and hence that scaling with `180/np.pi`. Looks good otherwise. – Divakar Jun 10 '18 at 11:31
  • To all.. deja vu https://codereview.stackexchange.com/questions/54347/calculating-element-wise-the-angles-between-two-lists-of-vectors – NaN Jun 11 '18 at 00:51
  • @Divakar, still finding my way around here. Thank you for pointing it out. – RRS Jun 12 '18 at 08:32
2

If you're working with 3D vectors, you can do this concisely using the toolbelt vg. It's a light layer on top of numpy and it works equally well with individual vectors and stacks of vectors.

import numpy as np
import vg

A = np.array([[1,0,0],
              [0,1,0],
              [0,0,1]])

B = np.array([[1,0,1],
              [1,1,0],
              [0,1,0]])

vg.angle(A, B)

I created the library at my last startup, where it was motivated by uses like this: simple ideas which are verbose or opaque in NumPy.

paulmelnikow
  • 16,895
  • 8
  • 63
  • 114