How can I efficiently bin 3D point cloud data based on point coordinates and desired bin grid size

Question

I have a large point cloud in open3D and I want to basically make a 3D grid and bin the points based on which cube they are in. Other have called it "binning in 3D space."

Example image with grids only in one direction (I want to split into 3D volumes)

Better Image of what I'm trying to do

Example:

import numpy as np

A = np.array([[ 0, -1, 10],
 [ 1, -2 ,11],
 [ 2, -3 ,12],
 [ 3, -4 ,13],
 [ 4, -5 ,14],
 [ 5, -6 ,15],
 [ 6, -7 ,16],
 [ 7, -8 ,17],
 [ 8, -9 ,18]])

#point 1: X,Y,Z
#point 2: X,Y,Z

print(A)

X_segments = np.linspace(0,8,3) #plane at beginning, middle and end - this creates 2 sections where data can be
Y_segments = np.linspace(-9,-1,3)
Z_segments = np.linspace(10,18,3)

#all of these combined form 4 cuboids where data can be 
#its also possible for the data to be outside these cuboids but we can ignore that

bin1 =  A where A[0,:] is > X_segments [0] and < X_segments[1]
    and A where A[1,:] is > Y_segments [0] and < Y_segments[1] 
    and A where A[2,:] is > Z_segments [0] and < Z_segments[1] 

bin2 =  A where A[0,:] is > X_segments [1] and < X_segments[2]
    and A where A[1,:] is > Y_segments [0] and < Y_segments[1] 
    and A where A[2,:] is > Z_segments [0] and < Z_segments[1] 

bin3 =  A where A[0,:] is > X_segments [1] and < X_segments[2]
    and A where A[1,:] is > Y_segments [1] and < Y_segments[2] 
    and A where A[2,:] is > Z_segments [0] and < Z_segments[1] 

bin4 =  A where A[0,:] is > X_segments [1] and < X_segments[2]
    and A where A[1,:] is > Y_segments [1] and < Y_segments[2] 
    and A where A[2,:] is > Z_segments [1] and < Z_segments[2]  

Thanks yall!

Answer 1

You can try the following:

import numpy as np

A = np.array([[ 0, -1, 10],
              [ 1, -2 ,11],
              [ 2, -3 ,12],
              [ 3, -4 ,13],
              [ 4, -5 ,14],
              [ 5, -6 ,15],
              [ 6, -7 ,16],
              [ 7, -8 ,17],
              [ 8, -9 ,18]])

X_segments = np.linspace(0,8,3)
Y_segments = np.linspace(-9,-1,3)
Z_segments = np.linspace(10,18,3)

edges = [X_segments, Y_segments, Z_segments]

print(edges) # just to show edges of the bins

We get:

[[ 0.  4.  8.]
 [-9. -5. -1.]
 [10. 14. 18.]]

Next, apply np.digitize along each coordinate axis:

coords = np.vstack([np.digitize(A.T[i], b, right=True) for i, b in enumerate(edges)]).T
print(coords)

This gives:

[[0 2 0]
 [1 2 1]
 [1 2 1]
 [1 2 1]
 [1 1 1]
 [2 1 2]
 [2 1 2]
 [2 1 2]
 [2 0 2]]

Rows of this array describe positions of the corresponding rows of A in the grid of bins. For example, the second row [1, 2, 1] indicates that the second row of A, i.e. [1, -2, 11] is in the first bin along the X-axis (since 0 < 1 <= 4), in the second bin along the Y-axis (since -5 < -2 <= -1), and in the first bin along the Z-axis (since 10 < 11 <= 14). Then you can pick elements belonging to each cuboid:

# select rows of A that are in the cuboid with coordinates [1, 2, 1]
A[np.all(coords == [1, 2, 1], axis=1)] 

This gives:

array([[ 1, -2, 11],
       [ 2, -3, 12],
       [ 3, -4, 13]])