Translating Numpy code into C#

Question

I suspect this might not be considered a great question, but I've completely hit a wall with it and could use some help.

I'm trying to implement this code:

http://www.nathanieltroutman.net/content/calculating-minimum-volume-bounding-box

in C#, with the original being in Python.

This was going fine until I hit this section:

def calcProjections(points, *vectors):
    """Calculates the projection of points (NxD) onto the vectors 
    (MxD) and return the projections p which is a matrix sized (N, M) 
    where N is the number of points and M is the number of vectors.
    p[i][j], is the projection of points[i] onto vectors[j] (which is
    between 0 and 1)."""

    u = np.array(vectors)

    # project the points onto the vectors into on fell swoop
    d = np.dot(points, u.T)

    # this is the dot product of each vector with itself
    v2 = np.diag(np.inner(u, u))

    p = d / v2

    return p

And I'm just struggling to decipher what is actually happening. I'm not sure what the author means by projecting onto specific vectors, or the format of the output (damn you duck typing). The descprition is a bit too vague for me too.

Does anyone have any suggestions or explanations as to what this is doing? Any help greatly appreciated.

Thanks.

Haha, possibly. Regarding the documentation. It's losing me a bit to be honest, in particular what it's doing when dealing with the two arrays that are dotted (dot producted?) together, but even despite that, the general goal of what this part of the code is doing is a little unclear to me. — djcmm476
– djcmm476, Commented Sep 22, 2014 at 15:50
Okay, well let's back up a step. Numpy implementation aside, do you follow the concept of projecting a point onto a vector? If not, give a quick read through this post, hopefully the math doesn't lose you. — Cory Kramer
– Cory Kramer, Commented Sep 22, 2014 at 15:52
Huh, so is this guy just doing the projection shown in that link, but for each point on to each of the three vectors and returning that in a matrix? — djcmm476
– djcmm476, Commented Sep 22, 2014 at 15:58

hpaulj · Accepted Answer · 2014-09-22 17:18:22Z

Here's a sample calculation, done in the interactive Ipython shell:

In [63]: points=np.arange(12,dtype=float).reshape(4,3)
In [64]: vectors=[np.array([1,0,0],dtype=float),np.array([0,1,1],dtype=float)]
In [65]: u=np.array(vectors)
In [66]: points
Out[66]: 
array([[  0.,   1.,   2.],
       [  3.,   4.,   5.],
       [  6.,   7.,   8.],
       [  9.,  10.,  11.]])
In [67]: u
Out[67]: 
array([[ 1.,  0.,  0.],
       [ 0.,  1.,  1.]])
In [68]: d=np.dot(points, u.T)
In [69]: d
Out[69]: 
array([[  0.,   3.],
       [  3.,   9.],
       [  6.,  15.],
       [  9.,  21.]])
In [70]: v2=np.diag(np.inner(u,u))
In [71]: d/v2
Out[71]: 
array([[  0. ,   1.5],
       [  3. ,   4.5],
       [  6. ,   7.5],
       [  9. ,  10.5]])

As specified in the function doc, input is (4,3) points, vectors a list of 2 (3,) vectors, and the output is a (4,2) array, p.

d is the matrix (dot) product of a 4x3 matrix with a 2x3 (or after transpose, 3x2) array, resulting in a 4x2. v2 could also be calculated as np.sum(u,u, axis=1), the magnitude of the 2 'vectors'. p is just that dot product normalized by v2.

If you are familiar with the Einstein summation notation (used in physics) the calculation is also expressed as:

np.einsum('ij,kj->ik',points,u)/np.einsum('ij,ij->i',u,u)

Collectives™ on Stack Overflow

Translating Numpy code into C#

1 Answer 1

Comments

Your Answer

Hot Network Questions

Collectives™ on Stack Overflow

1 Answer 1

Comments

Your Answer

Sign up or log in

Post as a guest

Related