CODEWITHSUNDEEP
pythonnumpy
KDecker
KDecker

Reputation: 7148

Locating indices of points within 3D bounding box using numpy logical_and or where?

I have a large list of points (x, y, z).

points = np.random.rand(999).reshape(333, 3)

I also have two points representing the minimum and maximum corners of a 3D bounding box representing points of interest

min_point = np.random.rand(3)
min_x, min_y, min_z = min_point[0], min_point[1], min_point[2]

max_point = np.random.rand(3)
max_x, max_y, max_z = max_point[0], max_point[1], max_point[2]

I am trying to select indices in points which lie within this bounding box but am having issues. I originally attempted using np.where

poi_inds = np.where(
    points[:, 0] > min_x and points[:, 0] < max_x and
    points[:, 1] > min_y and points[:, 1] < max_y and
    points[:, 2] > min_z and points[:, 2] < max_z
)

Though this results in

ValueError: The truth value of an array with more than one element is ambiguous.

because and'ing the results of each comparison is ambiguous (from my understand numpy can't decide to do element-wise or entire array and).

I found an SO answer which provides a solution in the 2D case and I've attempted to expand it to the 3D case

poi_inds = np.all(np.logical_and.reduce((
    points[:, 0] > min_x, points[:, 0] < max_x,
    points[:, 1] > min_y, points[:, 1] < max_y,
    points[:, 2] > min_z, points[:, 2] < max_z)))

Though this seems to always result in no points being selected (poid_inds = (False), if I remove the np.all I find all indices are False) even after running the SSCCE below multiple times.

import numpy as np

points = np.random.rand(999).reshape(333, 3)

rand_pt0 = np.random.rand(3)
rand_pt1 = np.random.rand(3)
# Ensure the minimum point is always less than the maximum
min_point = np.array([
    rand_pt0[0] if rand_pt0[0] < rand_pt1[0] else rand_pt1[0],
    rand_pt0[1] if rand_pt0[1] < rand_pt1[1] else rand_pt1[1],
    rand_pt0[2] if rand_pt0[2] < rand_pt1[2] else rand_pt1[2],
])
max_point = np.array([
    rand_pt1[0] if rand_pt0[0] < rand_pt1[0] else rand_pt0[0],
    rand_pt1[1] if rand_pt0[1] < rand_pt1[1] else rand_pt0[1],
    rand_pt1[2] if rand_pt0[2] < rand_pt1[2] else rand_pt0[2],
])

min_x, min_y, min_z = min_point[0], min_point[1], min_point[2]
max_x, max_y, max_z = max_point[0], max_point[1], max_point[2]

# poi_inds = np.where(
#     points[:, 0] > min_x and points[:, 0] < max_x and
#     points[:, 1] > min_y and points[:, 1] < max_y and
#     points[:, 2] > min_z and points[:, 2] < max_z
# )

poi_inds = np.logical_and.reduce((
    points[:, 0] > min_x, points[:, 0] < max_x,
    points[:, 1] > min_y, points[:, 1] < max_y,
    points[:, 2] > min_z, points[:, 2] < max_z))

What is the proper use of np.where or np.logical_and in this case to locate the points within the described bounding box?

Upvotes: 0

Views: 1412

Answers (3)

hummat
hummat

Reputation: 166

If points is a Numpy array you can use & instead of np.logical_and to write this as:

  mask_x = (points[:, 0] >= min_x) & (points[:, 0] <= max_x)
  mask_y = (points[:, 1] >= min_y) & (points[:, 1] <= max_y)
  mask_z = (points[:, 2] >= min_z) & (points[:, 2] <= max_z)
  mask = mask_x & mask_y & mask_z
  indices = np.where(mask)

Which, from a quick test, seems to be faster than the solution by KDecker.

Upvotes: 2

Sailesh Panda
Sailesh Panda

Reputation: 21

The problem is that you are comparing an array with an element.

points[:, 0] > min_x

An array cannot be compared with a single element.
Instead, you need to use a loop.

for i in range(333):
    points[i:0] > min_x

Upvotes: 0

KDecker
KDecker

Reputation: 7148

I have found a compact solution which utilizes np.where, functools.reduce, and np.intersect1d

from functools import reduce

poi_inds = reduce(np.intersect1d, (
    np.where(points[:, 0] > min_x),
    np.where(points[:, 0] < max_x),
    np.where(points[:, 1] > min_y),
    np.where(points[:, 1] < max_y),
    np.where(points[:, 2] > min_z),
    np.where(points[:, 2] < max_z)
))

On large lists >100M points this is many times faster than iterating through the points "manually".

poi_inds = []
for ind, pt in enumerate(points):
    if min_x < pt[0] < max_x and min_y < pt[1] < max_y and min_z < pt[2] < max_z:
        poi_inds.append(ind)

Upvotes: 1

Related Questions