I have an issue with numpy that I can't solve. I have 3D arrays (x,y,z) filled with 0 and 1. For instance, one slice in the z axis :
array([[1, 0, 1, 0, 1, 1, 0, 0],
[0, 0, 1, 1, 0, 1, 1, 0],
[1, 0, 1, 1, 0, 0, 0, 1],
[0, 0, 0, 0, 0, 0, 0, 0],
[1, 1, 1, 0, 1, 0, 0, 1],
[1, 0, 0, 0, 0, 1, 0, 1],
开发者_JAVA技巧 [0, 0, 0, 0, 1, 0, 0, 0],
[0, 0, 1, 0, 1, 1, 0, 1]])
And I want this result :
array([[1, 1, 1, 1, 1, 1, 0, 0],
[0, 0, 1, 1, 1, 1, 1, 0],
[1, 1, 1, 1, 1, 1, 1, 1],
[0, 0, 0, 0, 0, 0, 0, 0],
[1, 1, 1, 1, 1, 1, 1, 1],
[1, 1, 1, 1, 1, 1, 1, 1],
[0, 0, 0, 0, 1, 0, 0, 0],
[0, 0, 1, 1, 1, 1, 1, 1]])
That is to say, what I want to do for each slice z is to scan line by line right to left and left to right (x axis) and the first time I have a 1 I want to fill the rest of the line with ones.
Is there an efficient way to compute that ?
Thanks a lot.
Nico !
Accessing NumPy array elements one by one is not very efficient. You may do better with just plain Python lists. They also have an index
method which can search for the first entry of the value in the list.
from numpy import *
a = array([[1, 0, 1, 0, 1, 1, 0, 0],
[0, 0, 1, 1, 0, 1, 1, 0],
[1, 0, 1, 1, 0, 0, 0, 1],
[0, 1, 0, 0, 1, 0, 1, 0],
[1, 1, 1, 0, 1, 0, 0, 1],
[1, 0, 0, 0, 0, 1, 0, 1],
[0, 0, 0, 0, 1, 0, 0, 0],
[0, 0, 0, 0, 0, 0, 0, 0],
[0, 0, 1, 0, 1, 1, 0, 1]])
def idx_front(ln):
try:
return list(ln).index(1)
except ValueError:
return len(ln) # an index beyond line end
def idx_back(ln):
try:
return len(ln) - list(reversed(ln)).index(1) - 1
except ValueError:
return len(ln) # an index beyond line end
ranges = [ (idx_front(ln), idx_back(ln)) for ln in a ]
for ln, (lo,hi) in zip(a, ranges):
ln[lo:hi] = 1 # attention: destructive update in-place
print "ranges =", ranges
print a
Output:
ranges = [(0, 5), (2, 6), (0, 7), (1, 6), (0, 7), (0, 7), (4, 4), (8, 8), (2, 7)]
[[1 1 1 1 1 1 0 0]
[0 0 1 1 1 1 1 0]
[1 1 1 1 1 1 1 1]
[0 1 1 1 1 1 1 0]
[1 1 1 1 1 1 1 1]
[1 1 1 1 1 1 1 1]
[0 0 0 0 1 0 0 0]
[0 0 0 0 0 0 0 0]
[0 0 1 1 1 1 1 1]]
Actually, this is a basic binary image morphology operation.
You can do it in one step for the entire 3D array using scipy.ndimage.morphology.binary_fill_holes
You just need a slightly different structure element. In a nutshell, you want a structuring element that looks like this for the 2D case:
[[0, 0, 0],
[1, 1, 1],
[0, 0, 0]]
Here's a quick example:
import numpy as np
import scipy.ndimage as ndimage
a = np.array( [[1, 0, 1, 0, 1, 1, 0, 0],
[0, 0, 1, 1, 0, 1, 1, 0],
[1, 0, 1, 1, 0, 0, 0, 1],
[0, 1, 0, 0, 1, 0, 1, 0],
[1, 1, 1, 0, 1, 0, 0, 1],
[1, 0, 0, 0, 0, 1, 0, 1],
[0, 0, 0, 0, 0, 0, 0, 0],
[0, 0, 0, 0, 1, 0, 0, 0],
[0, 0, 1, 0, 1, 1, 0, 1]])
structure = np.zeros((3,3), dtype=np.int)
structure[1,:] = 1
filled = ndimage.morphology.binary_fill_holes(a, structure)
print filled.astype(np.int)
This yields:
[[1 1 1 1 1 1 0 0]
[0 0 1 1 1 1 1 0]
[1 1 1 1 1 1 1 1]
[0 1 1 1 1 1 1 0]
[1 1 1 1 1 1 1 1]
[1 1 1 1 1 1 1 1]
[0 0 0 0 0 0 0 0]
[0 0 0 0 1 0 0 0]
[0 0 1 1 1 1 1 1]]
The real advantage to this (Other than speed... It will be much faster and more memory efficient than using lists!) is that it will work just as well for 3D, 4D, 5D, etc arrays.
We just need to adjust the structuring element to match the number of dimensions.
import numpy as np
import scipy.ndimage as ndimage
# Generate some random 3D data to match what we want...
x = (np.random.random((10,10,20)) + 0.5).astype(np.int)
# Make the structure (I'm assuming that "z" is the _last_ dimension!)
structure = np.zeros((3,3,3))
structure[1,:,1] = 1
filled = ndimage.morphology.binary_fill_holes(x, structure)
print x[:,:,5]
print filled[:,:,5].astype(np.int)
Here's a slice from the random input 3D array:
[[1 0 1 0 1 1 0 1 0 0]
[1 0 1 1 0 1 0 1 0 0]
[1 0 0 1 0 1 1 1 1 0]
[0 0 0 1 1 0 1 0 0 0]
[1 0 1 0 1 0 0 1 1 0]
[1 0 1 1 0 1 0 0 0 1]
[0 1 0 1 0 0 1 0 1 0]
[0 1 1 0 1 0 0 0 0 1]
[0 0 0 1 1 1 1 1 0 1]
[1 0 1 1 1 1 0 0 0 1]]
And here's the filled version:
[[1 1 1 1 1 1 1 1 0 0]
[1 1 1 1 1 1 1 1 0 0]
[1 1 1 1 1 1 1 1 1 0]
[0 0 0 1 1 1 1 0 0 0]
[1 1 1 1 1 1 1 1 1 0]
[1 1 1 1 1 1 1 1 1 1]
[0 1 1 1 1 1 1 1 1 0]
[0 1 1 1 1 1 1 1 1 1]
[0 0 0 1 1 1 1 1 1 1]
[1 1 1 1 1 1 1 1 1 1]]
The key difference here is that we did this for every slice of the entire 3D array in one step.
After a moments thought, following your description and corner case with all zero rows, this will be still quite straightforward with numpy
like:
In []: A
Out[]:
array([[1, 0, 1, 0, 1, 1, 0, 0],
[0, 0, 1, 1, 0, 1, 1, 0],
[1, 0, 1, 1, 0, 0, 0, 1],
[0, 1, 0, 0, 1, 0, 1, 0],
[1, 1, 1, 0, 1, 0, 0, 1],
[1, 0, 0, 0, 0, 1, 0, 1],
[0, 0, 0, 0, 0, 0, 0, 0],
[0, 0, 0, 0, 1, 0, 0, 0],
[0, 0, 1, 0, 1, 1, 0, 1]])
In []: v= 0< A.sum(1) # work only with rows at least one 1
In []: A_v= A[v, :]
In []: (r, s), a= A_v.nonzero(), arange(v.sum())
In []: se= c_[searchsorted(r, a), searchsorted(r, a, side= 'right')- 1]
In []: for k in a: A_v[k, s[se[k, 0]]: s[se[k, 1]]]= 1
..:
In []: A[v, :]= A_v
In []: A
Out[]:
array([[1, 1, 1, 1, 1, 1, 0, 0],
[0, 0, 1, 1, 1, 1, 1, 0],
[1, 1, 1, 1, 1, 1, 1, 1],
[0, 1, 1, 1, 1, 1, 1, 0],
[1, 1, 1, 1, 1, 1, 1, 1],
[1, 1, 1, 1, 1, 1, 1, 1],
[0, 0, 0, 0, 0, 0, 0, 0],
[0, 0, 0, 0, 1, 0, 0, 0],
[0, 0, 1, 1, 1, 1, 1, 1]])
Update:
After some more tinkering, here is a more 'pythonic' implementation and way much simpler, than the above one. So, the following lines:
for k in xrange(A.shape[0]):
m= A[k].nonzero()[0]
try: A[k, m[0]: m[-1]]= 1
except IndexError: continue
are quite straightforward ones. And they'll perform very well, indeed.
I can't think of a more efficient way than what you describe:
For every line
Scan line from the left until you find a
1
.If no
1
is find continue with next line.Otherwise scan from the right to find the last
1
in the line.Fill everything in the current line between the positions from 1. and 3. with
1
s.
精彩评论