1

New Post: Processing satellite conjunctions with numpy efficiently

Original Post: I have a numpy array of shape n x m x r, where the n axis represents an object, the m axis represents a timestep and the r axis represents a position vector in 3-d space. I have an array containing three (x, y and z position values) of m objects at n points in time. This is the format my data is delivered in (from the python module sgp4, specifically a Satrec_Array if anyone's interested) so I can't move this further up the data processing chain.

I want to be able to represent this as an m x n array of position vectors, so essentially "collapsing" the position axis into a single array, so an array containing m x n elements, each of which is a position vector relating to an object at a time.

I'm struggling to find a way to do this efficiently - I can do it with bog standard python loops, but as I scale up the number of objects and timesteps this becomes massively inefficient.

I'm not very well versed in numpy, but none of the solutions I've searched or attempts at using methods such as [h/v/d]stack etc. have given me the right final shape. I also looked into vectorization but as far as I can see that just implements a loop under the hood?

Example with random numbers and an input array of shape (3,2,3)

In[1]: m = 3
       n = 2
       r = 3

In[2]: a = np.random.random((m,n,r))

In[3]: a
Out[3]: 
array([[[0.8416, 0.3694, 0.5708],
        [0.3779, 0.579 , 0.207 ]],

       [[0.7871, 0.6547, 0.0047],
        [0.1115, 0.1445, 0.6147]],

       [[0.8538, 0.2821, 0.8094],
        [0.6214, 0.0147, 0.5852]]])

In[4]: a.shape
Out[4]: (3, 2, 3)
       
In[4]: new_a = np.empty(shape=(m,n), dtype=object)
       for i in range(m):
           for j in range(n):
            new_a[i,j] = a[i,j,:]

In[5]: new_a
Out[5]: 
array([[array([0.8416, 0.3694, 0.5708]), array([0.3779, 0.579 , 0.207 ])],
       [array([0.7871, 0.6547, 0.0047]), array([0.1115, 0.1445, 0.6147])],
       [array([0.8538, 0.2821, 0.8094]), array([0.6214, 0.0147, 0.5852])]],
      dtype=object)

In[6]: new_a.shape
Out[6]: (3, 2)
7
  • 3
    Can you give an example of input and desired output? I can't get it from your description: is your r == 3 and the [x, y, z] coordinates go along the r axis? And you want to change it into what?
    – Headcrab
    Commented Mar 18 at 20:59
  • 2
    Even better, include complete code with your for loop solution. I think you can do this with a moveaxis and reshape, but it's not clear whether the notation is consistent anove. m*n*r won't in general equal m*n, so you can't collapse an array with shape (m, n, r) to an array with shape (m, n) as requested (unless you want an object array in which element is itself an array). Commented Mar 18 at 21:27
  • Good shout, I will add an example tomorrow - thank you!
    – SPZHunter
    Commented Mar 18 at 21:32
  • 1
    So instead of an array of size [n, m, 3] dtype=float you want an array of size [n, m] dtype=array? Why would you want to do that? Commented Mar 19 at 10:46
  • 1
    It would help to show what you would if you had such an array. It doesn't really make sense to use object arrays in numpy.
    – mozway
    Commented Mar 19 at 10:55

1 Answer 1

1

"I'm struggling to find a way to do this efficiently - I can do it with bog standard python loops, but as I scale up the number of objects and timesteps this becomes massively inefficient."

What you want to do makes no sense in numpy. Unless using an object dtype, there is no way to have anything else than numeric data as items. You should keep your nd-array.

In fact, what you have ((m, n, r) shape) is already a (m, n) array of vectors of dimension r. Think of it this way if you like rather than (m, n, r).

Numpy operation can very efficiently operate on a subset of the dimensions. For example:

# add 10/100/1000 to the first dimension
a + np.array([10, 100, 1000])[:, None, None]

array([[[  10.8416,   10.3694,   10.5708],
        [  10.3779,   10.579 ,   10.207 ]],

       [[ 100.7871,  100.6547,  100.0047],
        [ 100.1115,  100.1445,  100.6147]],

       [[1000.8538, 1000.2821, 1000.8094],
        [1000.6214, 1000.0147, 1000.5852]]])
a + np.array([10, 100])[:, None]

# add 10/100 to the second dimension
array([[[ 10.8416,  10.3694,  10.5708],
        [100.3779, 100.579 , 100.207 ]],

       [[ 10.7871,  10.6547,  10.0047],
        [100.1115, 100.1445, 100.6147]],

       [[ 10.8538,  10.2821,  10.8094],
        [100.6214, 100.0147, 100.5852]]])

Your issue is most likely a XY problem. You should keep the (m, n, r) shape and try to solve your ultimate goal with a nd-array.

6
  • I'm used to using MATLAB so I might need to go and properly understand indexing and sub-sets in numpy better! In the end, I want to be able to pass the [X, Y, Z] values of each object at each timestep to a function, for clarity.
    – SPZHunter
    Commented Mar 19 at 11:02
  • @SPZHunter You can access the [X, Y, Z] values through indexing: a[0, 0] returns the [X, Y, Z].
    – tnknepp
    Commented Mar 19 at 11:25
  • Yep - but my issue is I don't know how to pass every [X, Y, Z] values to a function efficiently, currently my only knowledge is doing it through loops and indexing a[i, j] one by one
    – SPZHunter
    Commented Mar 19 at 11:31
  • 1
    Better edit your question to provide a minimal reproducible example of how you would call your function and what it does. However be aware that a function call is often a non-vector operation in python, so this might not be different than a loop. The ideal would be to explain your full goal from A to Z.
    – mozway
    Commented Mar 19 at 11:37
  • @SPZHunter are you trying to perform an operation on all of those [x, y, z]-triplets "at once" instead of performing it sequentially for every triplet at every [row, col] of your array? Like, for example, replacing each [x, y, z] with the length of the corresponding vector by calling np.linalg.norm(a, axis = -1), instead of iterating by m, n of your a and calling norm for every [x, y, z]?
    – Headcrab
    Commented Mar 20 at 1:33

Start asking to get answers

Find the answer to your question by asking.

Ask question

Explore related questions

See similar questions with these tags.