pythonnumpyopen3d

Optimizing Point Cloud to Voxel Grid with Max Sampling in NumPy


I have two arrays that represent the point coordinates and values respectively. To max sample from this point cloud, I am initializing a grid with the desired size, and looping over each point to assign the max values:

N = 1000000
coords = np.random.randint(0, 256, size=(N, 3))
vals = np.random.rand(N, 3)

grid = np.zeros((3, 256, 256, 256), dtype=np.float16)
for i, pt in enumerate(coords):
    x, y, z = pt

    grid[0, x, y, z] = max(grid[0, x, y, z], vals[i, 0])
    grid[1, x, y, z] = max(grid[1, x, y, z], vals[i, 1])
    grid[2, x, y, z] = max(grid[2, x, y, z], vals[i, 2])

Is there a way I can do this through NumPy without the for loop (which is very slow)?


Solution

  • np.maximum.at can be used to perform the max operation in place at specific indices.

    For example, if we have two arrays arr1 and arr2 to compare and an indices array, we can do this:

    arr1 = np.array([1, 2, 3, 4, 5])
    arr2 = np.array([10, 1, 8, 7, 3])
    indices = np.array([0, 1, 2, 1, 4])
    
    # Perform in-place maximum update
    np.maximum.at(arr1, indices, arr2)
    
    print(arr1)
    # Output: [10  7  8  4  5]
    

    So I believe in your case, you can do something like this:

    indices = coords[:, 0], coords[:, 1], coords[:, 2]
    
    np.maximum.at(grid[0], indices, vals[:, 0])
    np.maximum.at(grid[1], indices, vals[:, 1])
    np.maximum.at(grid[2], indices, vals[:, 2])