Cross multiple 2D np arrays to define zones

Question

Using this small reproducible example, I still could not create a new array of integers from 3 arrays that contains unique groups in all three input arrays.

Arrays are associated with topographic properties:

import numpy as np asp = np.array([8,1,1,2,7,8,2,3,7,6,4,3,6,5,5,4]).reshape((4,4)) #aspect slp = np.array([9,10,10,9,9,12,12,9,10,11,11,9,9,9,9,9]).reshape((4,4)) #slope elv = np.array([13,14,14,13,14,15,16,14,14,15,16,14,13,14,14,13]).reshape((4,4)) #elevation 

The idea is that geographic contours are divided into 3 different properties using GIS routines:

• 1-8 for the aspect (1 = north side, 2 = north-east side, etc.).
• 9-12 for tilt (9 = gentle slope ... 12 = steep tilt)
• 13-16 for height (13 = lowest marks ... 16 = highest marks)

The small picture below is trying to portray the result that I got after (the array shown in the lower left corner). Please note that the “answer” given in the graph is just one possible answer. I'm not interested in the final arrangement of integers in the resulting array, as long as the last array contains an integer in each row / column index that identifies unique groups.

For example, the array indices in [0,1] and [0,2] have the same aspect, slope, and height, and therefore get the same integer identifier in the resulting array.

Is there a numpy firmware for this kind of thing?

Answer 1

This can be done using numpy.unique() , and then the following mapping:

The code:

 combined = 10000 * asp + 100 * slp + elv unique = dict(((v, i + 1) for i, v in enumerate(np.unique(combined)))) combined_unique = np.vectorize(unique.get)(combined) 

Security Code:

 import numpy as np asp = np.array([8, 1, 1, 2, 7, 8, 2, 3, 7, 6, 4, 3, 6, 5, 5, 4]).reshape((4, 4)) # aspect slp = np.array([9, 10, 10, 9, 9, 12, 12, 9, 10, 11, 11, 9, 9, 9, 9, 9]).reshape((4, 4)) # slope elv = np.array([13, 14, 14, 13, 14, 15, 16, 14, 14, 15, 16, 14, 13, 14, 14, 13]).reshape((4, 4)) combined = 10000 * asp + 100 * slp + elv unique = dict(((v, i + 1) for i, v in enumerate(np.unique(combined)))) combined_unique = np.vectorize(unique.get)(combined) print(combined_unique) 

Results:

 [[12 1 1 2] [10 13 3 4] [11 9 6 4] [ 8 7 7 5]] 

Answer 2

Each place in the grid is associated with a tuple consisting of one value from asp , slp and elv . For example, the upper left corner has a tuple (8,9,13) . We would like to match this tuple with a number that uniquely identifies this set.

One way to do this is to think of (8,9,13) as an index in the np.arange(9*13*17).reshape(9,13,17) 3D array. This particular array was chosen to hold the largest values ​​in asp , slp and elv :

 In [107]: asp.max()+1 Out[107]: 9 In [108]: slp.max()+1 Out[108]: 13 In [110]: elv.max()+1 Out[110]: 17 

Now we can match the set (8,9,13) with the number 1934:

 In [113]: x = np.arange(9*13*17).reshape(9,13,17) In [114]: x[8,9,13] Out[114]: 1934 

If we do this for each place in the grid, we will get a unique number for each location. We could end up here by allowing these unique numbers to serve as tags.

Or, we can generate smaller labels (starting at 0 and increasing by 1) using np.unique with return_inverse=True :

 uniqs, labels = np.unique(vals, return_inverse=True) labels = labels.reshape(vals.shape) 

So for example

 import numpy as np asp = np.array([8,1,1,2,7,8,2,3,7,6,4,3,6,5,5,4]).reshape((4,4)) #aspect slp = np.array([9,10,10,9,9,12,12,9,10,11,11,9,9,9,9,9]).reshape((4,4)) #slope elv = np.array([13,14,14,13,14,15,16,14,14,15,16,14,13,14,14,13]).reshape((4,4)) #elevation x = np.arange(9*13*17).reshape(9,13,17) vals = x[asp, slp, elv] uniqs, labels = np.unique(vals, return_inverse=True) labels = labels.reshape(vals.shape) 

gives

 array([[11, 0, 0, 1], [ 9, 12, 2, 3], [10, 8, 5, 3], [ 7, 6, 6, 4]]) 

---

This method works fine as long as the values ​​in asp , slp and elv are small integers. If the integers were too large, the product of their maxima could overwhelm the maximum allowable value that can be transferred to np.arange . Moreover, creating such a large array would be inefficient. If the values ​​were float, then they could not be interpreted as indices in the 3D array x .

So, to solve these problems, use np.unique to convert the values ​​in asp , slp and elv to unique integer labels:

 indices = [ np.unique(arr, return_inverse=True)[1].reshape(arr.shape) for arr in [asp, slp, elv] ] M = np.array([item.max()+1 for item in indices]) x = np.arange(M.prod()).reshape(M) vals = x[indices] uniqs, labels = np.unique(vals, return_inverse=True) labels = labels.reshape(vals.shape) 

which gives the same result as above, but works even if asp , slp , elv were float and / or large integers.

---

Finally, we can avoid generating np.arange :

 x = np.arange(M.prod()).reshape(M) vals = x[indices] 

by computing vals as the product of indices and steps:

 M = np.r_[1, M[:-1]] strides = M.cumprod() indices = np.stack(indices, axis=-1) vals = (indices * strides).sum(axis=-1) 

So, all together:

 import numpy as np asp = np.array([8,1,1,2,7,8,2,3,7,6,4,3,6,5,5,4]).reshape((4,4)) #aspect slp = np.array([9,10,10,9,9,12,12,9,10,11,11,9,9,9,9,9]).reshape((4,4)) #slope elv = np.array([13,14,14,13,14,15,16,14,14,15,16,14,13,14,14,13]).reshape((4,4)) #elevation def find_labels(*arrs): indices = [np.unique(arr, return_inverse=True)[1] for arr in arrs] M = np.array([item.max()+1 for item in indices]) M = np.r_[1, M[:-1]] strides = M.cumprod() indices = np.stack(indices, axis=-1) vals = (indices * strides).sum(axis=-1) uniqs, labels = np.unique(vals, return_inverse=True) labels = labels.reshape(arrs[0].shape) return labels print(find_labels(asp, slp, elv)) # [[ 3 7 7 0] # [ 6 10 12 4] # [ 8 9 11 4] # [ 2 5 5 1]] 

Answer 3

This is similar to a similar problem with marking unique areas of an image. This is the function I wrote for this, although first you need to combine 3 arrays into one 3D array.

 def labelPix(pix): height, width, _ = pix.shape pixRows = numpy.reshape(pix, (height * width, 3)) unique, counts = numpy.unique(pixRows, return_counts = True, axis = 0) unique = [list(elem) for elem in unique] labeledPix = numpy.zeros((height, width), dtype = int) offset = 0 for index, zoneArray in enumerate(unique): index += offset zone = list(zoneArray) zoneArea = (pix == zone).all(-1) elementsArray, numElements = scipy.ndimage.label(zoneArea) elementsArray[elementsArray!=0] += offset labeledPix[elementsArray!=0] = elementsArray[elementsArray!=0] offset += numElements return labeledPix 

This will indicate unique combinations of three values, as well as assign separate labels for zones that have the same combination of 3 values, but are not in contact with each other.

 asp = numpy.array([8,1,1,2,7,8,2,3,7,6,4,3,6,5,5,4]).reshape((4,4)) #aspect slp = numpy.array([9,10,10,9,9,12,12,9,10,11,11,9,9,9,9,9]).reshape((4,4)) #slope elv = numpy.array([13,14,14,13,14,15,16,14,14,15,16,14,13,14,14,13]).reshape((4,4)) #elevation pix = numpy.zeros((4,4,3)) pix[:,:,0] = asp pix[:,:,1] = slp pix[:,:,2] = elv print(labelPix(pix)) 

returns:

 [[ 0 1 1 2] [10 12 3 4] [11 9 6 4] [ 8 7 7 5]] 

Answer 4

Here's a simple Python method using itertools.groupby . It requires the input to be 1D lists, but this should not be a serious problem. The strategy is to pin lists together with the index number, and then sort the resulting columns. Then we group the same columns together, ignoring the index number when comparing the columns. Then we collect index numbers from each group and use them to create the final output list.

 from itertools import groupby def show(label, seq): print(label, ' '.join(['{:2}'.format(u) for u in seq])) asp = [8, 1, 1, 2, 7, 8, 2, 3, 7, 6, 4, 3, 6, 5, 5, 4] slp = [9, 10, 10, 9, 9, 12, 12, 9, 10, 11, 11, 9, 9, 9, 9, 9] elv = [13, 14, 14, 13, 14, 15, 16, 14, 14, 15, 16, 14, 13, 14, 14, 13] size = len(asp) a = sorted(zip(asp, slp, elv, range(size))) groups = sorted([u[-1] for u in g] for _, g in groupby(a, key=lambda t:t[:-1])) final = [0] * size for i, g in enumerate(groups, 1): for j in g: final[j] = i show('asp', asp) show('slp', slp) show('elv', elv) show('out', final) 

Exit

 asp 8 1 1 2 7 8 2 3 7 6 4 3 6 5 5 4 slp 9 10 10 9 9 12 12 9 10 11 11 9 9 9 9 9 elv 13 14 14 13 14 15 16 14 14 15 16 14 13 14 14 13 out 1 2 2 3 4 5 6 7 8 9 10 7 11 12 12 13 

---

There is no need to do this second sort, we could just use a simple comp list

 groups = [[u[-1] for u in g] for _, g in groupby(a, key=lambda t:t[:-1])] 

or generator expression

 groups = ([u[-1] for u in g] for _, g in groupby(a, key=lambda t:t[:-1])) 

I did this only so that my result matches the output in the question.

Answer 5

Here is one way to solve this problem using dictionary search.

 from collections import defaultdict import itertools group_dict = defaultdict(list) idx_count = 0 for a, s, e in np.nditer((asp, slp, elv)): asp_tuple = (a.tolist(), s.tolist(), e.tolist()) if asp_tuple not in group_dict: group_dict[asp_tuple] = [idx_count+1] idx_count += 1 else: group_dict[asp_tuple].append(group_dict[asp_tuple][-1]) list1d = list(itertools.chain(*list(group_dict.values()))) np.array(list1d).reshape(4, 4) # result array([[ 1, 2, 2, 3], [ 4, 5, 6, 7], [ 7, 8, 9, 10], [11, 12, 12, 13]])