i trying speed piece of code below vectorization:
[rows,cols] = flow_direction_np.shape elevation_gain = np.zeros((rows,cols), np.float) [i, j], flow in np.ndenumerate(flow_direction_np): try: if flow == 32: elevation_gain[i - 1, j - 1] = elevation_gain[i - 1, j - 1] + sediment_transport_np[i, j] elif flow == 64: elevation_gain[i - 1, j] = elevation_gain[i - 1, j] + sediment_transport_np[i, j] elif flow == 128: elevation_gain[i - 1, j + 1] = elevation_gain[i - 1, j + 1] + sediment_transport_np[i, j] elif flow == 16: elevation_gain[i, j - 1] = elevation_gain[i, j - 1] + sediment_transport_np[i, j] elif flow == 1: elevation_gain[i, j + 1] = elevation_gain[i, j + 1] + sediment_transport_np[i, j] elif flow == 2: elevation_gain[i + 1, j + 1] = elevation_gain[i + 1, j + 1] + sediment_transport_np[i, j] elif flow == 4: elevation_gain[i + 1, j] = elevation_gain[i + 1, j] + sediment_transport_np[i, j] elif flow == 8: elevation_gain[i + 1, j - 1] = elevation_gain[i + 1, j - 1] + sediment_transport_np[i, j] except indexerror: elevation_gain[i, j] = 0 this how code looks @ moment:
elevation_gain = np.zeros_like(sediment_transport_np) nrows, ncols = flow_direction_np.shape lookup = {32: (-1, -1), 16: (0, -1), 8: (+1, -1), 4: (+1, 0), 64: (-1, 0), 128:(-1, +1), 1: (0, +1), 2: (+1, +1)} # initialize array "shifted" mask shifted = np.zeros((nrows+2, ncols+2), dtype=bool) # pad elevation gain zeros tmp = np.zeros((nrows+2, ncols+2), elevation_gain.dtype) tmp[1:-1, 1:-1] = elevation_gain elevation_gain = tmp value, (row, col) in lookup.iteritems(): mask = flow_direction_np == value # reset "shifted" mask shifted.fill(false) shifted[1:-1, 1:-1] = mask # shift mask right amount given value shifted = np.roll(shifted, row, 0) shifted = np.roll(shifted, col, 1) # set values in elevation change offset value in sed_trans elevation_gain[shifted] = elevation_gain[shifted] + sediment_transport_np[mask] the trouble having aren't giving me same result @ end suggestions going wrong?
you can improve performance using np.where indices conditions happen:
ind = np.where( flow_direction_np==32 ) you see ind tuple 2 elements, first indices of first axis , second of second axis of flow_direction_np array.
you can work out indices apply shifts: i-1, j-1 , on...
ind_32 = (ind[0]-1, ind[1]-1) then use fancy indexing update arrays:
elevation_gain[ ind_32 ] += sediment_transport_np[ ind ] edit: applying concept case give this:
lookup = {32: (-1, -1), 16: ( 0, -1), 8: (+1, -1), 4: (+1, 0), 64: (-1, 0), 128: (-1, +1), 1: ( 0, +1), 2: (+1, +1)} num, shift in lookup.iteritems(): ind = np.where( flow_direction_np==num ) ind_num = ind[0] + shift[0], ind[1] + shift[1] elevation_gain[ ind_num] += sediment_transport_np[ ind ] 
Comments
Post a Comment