Ticket #9115: neighbour_grid.py

File neighbour_grid.py, 3.0 KB (added by Owen Arnold, 7 years ago)
Line 
1import numpy as np
2import math
3
4n_dims = 3
5n_bins = 4
6
7def make_nd_array(n_dims, n_bins):
8        base_array = np.arange(0, int(math.pow(n_bins, n_dims)))
9       
10        shape = []
11        for i in range(0, n_dims):
12                shape.append(n_bins)
13        as_grid = base_array.reshape(shape)
14        return as_grid
15       
16
17def get_indexes_from_linear(index, grid):
18        out_indexes = list()
19        shape = grid.shape
20        nBins  = shape[0]
21        # first index
22        ind    = index%nBins
23       
24        out_indexes.append(ind)
25       
26        #what left in linear index after first was removed
27       
28        rest = index/nBins;
29        for d in range(1, len(shape)):
30                nBins = shape[d]
31                ind  = rest%nBins;
32                out_indexes.append( ind )
33                rest = rest/nBins
34               
35        return out_indexes
36
37#as_grid = np.array([[0, 1, 2, 3], [4,5,6,7], [8,9,10,11], [12, 13, 14, 15]]) # 2D
38
39#as_grid = as_grid.flatten() # 1D
40
41as_grid = make_nd_array(n_dims, n_bins)
42
43
44def get_neighbours(grid, i):
45       
46        n_dims = len(grid.shape)
47        offset = 1
48        back =  i - offset
49        forward =i + offset
50       
51       
52        if forward/grid.shape[0] == i/grid.shape[0]:
53                        print forward
54        if back >= 0 and back/grid.shape[0] == i/grid.shape[0]:
55                        print back
56       
57       
58        for j in range(1, n_dims):
59                offset = offset * grid.shape[j-1]
60               
61                offset_next = offset * grid.shape[j]
62               
63                back = i - offset
64                forward =  i + offset
65               
66               
67                if forward/offset_next == i/offset_next:
68                        print forward
69                if back >= 0 and back/offset_next == i/offset_next:
70                        print back
71
72       
73def is_within_boundaries(neigh_index, current_indexes, grid):
74        neigh_indexes  = get_indexes_from_linear(neigh_index, grid)
75        for index in range(len(current_indexes)):
76                i = current_indexes[index]
77                j = neigh_indexes[index]
78                abs_diff = math.fabs(i-j) 
79                if abs_diff < 0 or abs_diff > 1 :
80                        return False
81        return True
82                               
83               
84def get_neighbours_all_touching(grid, i):
85       
86        n_dims = len(grid.shape)
87        offset = 1
88       
89       
90        back =  i - offset
91        forward =i + offset
92
93        indexes = get_indexes_from_linear(i, grid)
94       
95        permitations = list()
96        permitations.append(0)
97       
98        dim_min = [None]
99        dim_max = [None]
100       
101        permitations.append(1)
102        permitations.append(-1)
103
104        if forward/grid.shape[0] != i/grid.shape[0]:
105                        dim_max[0] = indexes[0]
106        if back >= 0 and back/grid.shape[0] != i/grid.shape[0]:
107                        dim_min[0] = indexes[0]
108       
109       
110        for j in range(1, n_dims):
111               
112                dim_min.append(None)
113                dim_max.append(None)
114               
115                offset = offset * grid.shape[j-1]
116                offset_next = offset * grid.shape[j]
117               
118                back = i - offset
119                forward =  i + offset
120               
121                if forward/offset_next != i/offset_next:
122                        dim_max[j] = indexes[j]
123                if back >= 0 and back/offset_next != i/offset_next:
124                        dim_min[j] = indexes[j]
125               
126               
127               
128                extended_perms = list()
129                for perm in permitations:
130                       
131                        extended_perms.append(perm + offset)
132                        extended_perms.append(-1*(perm+offset))
133                       
134                permitations.extend(extended_perms)
135               
136                print "n perms", len(permitations)
137                print permitations
138               
139                           
140        for perm in permitations:
141                if perm == 0:
142                        continue
143               
144                neigh_index = i + perm
145               
146                if is_within_boundaries(neigh_index, indexes, grid):
147                        print neigh_index
148               
149               
150               
151index = 1
152
153print as_grid
154
155get_neighbours_all_touching(as_grid, index)
156
157
158
159
160
161               
162               
163
164