ssm
ssm

Reputation: 640

numpy extract 3d cubes from 3d array

I would like to extract the 3D cubes (3x3x3) from a boolean 3D array of (180x180x197). This is similar to scipy.ndimage.measurements.label but need to fixed size of (3x3x3). Is there a fast way of doing this than using for loops.

Upvotes: 2

Views: 1439

Answers (2)

ssm
ssm

Reputation: 640

Final solution with help from dnalow

get_starting_point = numpy.vectorize(lambda sl: sl.start)
s = ndimage.generate_binary_structure(3,3)
result = []

while True:
    centers = ndimage.minimum_filter(b, 3, mode="constant")
    labels, num = ndimage.measurements.label(centers, s)
    if not num:
        break
    locations = ndimage.measurements.find_objects(labels)
    sp = get_starting_point(locations)
    for p in sp:
        if numpy.all(a[p[0]-1:p[0]+2, p[1]-1:p[1]+2, p[2]-1:p[2]+2]):
            result.append(p) 
        b[p[0]-1:p[0]+2, p[1]-1:p[1]+2, p[2]-1:p[2]+2] = False

numiter = len(result)
results = numpy.vstack(result)

Upvotes: 1

dnalow
dnalow

Reputation: 984

In your special case, I suggest to use ndimage.minimum_filter Let's say your array is called ``a''. The following:

centers = ndimage.minimum_filter(a, 3, mode="constant")

will only contain ones where your array contained such box of True's Then you can use scipy.ndimage.measurements.label with the default structure to classify the boxes and maybe identify connected boxs. To locate them you can use ndimage.measurements.find_objects

edit:

The way above, you will correctly get the centers of all cubes in the array. To be clear, I think that is the answer to your initial question. In the comments, it turned out, that indeed only non-overlapping cubes are needed. Therefore, one needs to analyze the output of minimum_filter, where I can imagine many approaches.

One can use the following to obtain only one cube per cluster:

s = ndimage.generate_binary_structure(3,3)
labels, num = ndimage.measurements.label(centers, s)
locations = ndimage.measurements.find_objects(labels)
locations = map(lambda slices: [slices[i].start for i in xrange(3)], locations)

Now the problem occurs that cubes are lost which are not overlapping but share a face. Actually one can imagine quite complicated structures of non-overlapping, facesharing cubes. Certainly, there are several solutions (sets of non-overlapping cubes) that can be found for this problem. So it is a completely new task to choose a set of cubes from the found centers and I think you will have to find one which is ideal for you.

One way could be to iterate through all solutions and set each found cube to False:

get_starting_point = numpy.vectorize(lambda sl: sl.start) #to be applied on slices
s = ndimage.generate_binary_structure(3,3)
result = []

while True:
    labels, num = ndimage.measurements.label(centers, s)
    if not num:
        break
    locations = ndimage.measurements.find_objects(labels)
    sp = get_starting_point(locations)
    result.append(sp)
    for p in sp:
        centers[p[0]-1:p[0]+2, p[1]-1:p[1]+2, p[2]-1:p[2]+2] = False

numiter = len(results)
results = numpy.vstack(results)

I guess only very few iterations will be necessary.

I hope this is what you were looking for

Upvotes: 2

Related Questions