python - Edge of Zone in 2D NumPy Array -
i'm using python 2.7.10 arcpy site package (arcgis 10.4.1) attempt solve following problem. have 2d integer numpy array of zones representing watershed boundaries (watershedarray) , 2d floating point numpy array representing elevation values (elevationarray). attempting find lowest value in elevationarray along outside edge of contiguous group of particular watershedarray values. - in first iteration, working watershed id = 1, need find lowest value in elevationarray masked edge of watershedarray = 1. watershed outlet - now, based on outlet, need find lowest adjacent elevationarray value outside of watershedarray = 1. let's ends being watershedarray = 16. - in next iteration, need find lowest value in elevationarray, masked outer edge of combined contiguous group of watershedarray = 1 or 16.
to date i've been using np.where/zip make lists of appropriate watershedarray indices. loop through list nested statements check edge cells (insidecells < 9), loop through edgelist lowest value. create list of edge cells matching value (there possibility of multiple outlets. seems inefficient , i'm sure using combination of masks, scipy erosion, etc. can things done more effectively, i'm not sure how proceed.
here visual of problem. starting brown watershed in upper right. outlets in red.
def findoutlet(demarray, sinkarray, sinkid, contributingset): # array shape. used ensure adjacent neighbor checks remain within bounds arrayrows = int(demarray.shape[0]) arraycols = int(demarray.shape[1]) # establish set of cells who's neighbors have been scanned, not need rechecked norecheckset = set() controws,contcols = np.where(sinksheds == sinkid) newcontributinglist = zip(controws,contcols) contributinglist = list(contributingset) + newcontributinglist contributingset = set(contributinglist) del newcontributinglist #arcpy.addmessage("final contributing cell count: " + str(len(contributingset))) # initiate z value edge cell candidate larger acceptable legitimate value outletz = 99999.0 #initialize empty edge set edgelist = [] # scan contributing cells in contributingset determine if @ edge of contributing zone counting adjacent cells within contributingset thecontributor in contributinglist: contributingrow = thecontributor[0] contributingcol = thecontributor[1] # elevation value of current contributing cell contributingz = float(demarray[contributingrow][contributingcol]) # initiate counter number of adjacent cells in contributingset insidecells = 0 # iterate through adjacent row , column surrounding current contributing cell adjacentcol in range(contributingcol - 1, contributingcol + 2): adjacentrow in range(contributingrow - 1, contributingrow + 2): # if adjacent cell in contributingset, increment insidecells counter if (adjacentrow, adjacentcol) in contributingset: insidecells += 1 # if insidecells < 9 (9 = self , neighbors), , if contributing cell's elevation lower current edgez elevation (minimum), assign new edgez candidate if insidecells < 9: edgelist.append((thecontributor)) if contributingz < outletz: outletz = contributingz outletlist = [] theedge in edgelist: edgerow = theedge[0] edgecol = theedge[1] if float(demarray[edgerow][edgecol]) == outletz: outletlist.append(theedge) # arcpy.addmessage("potential outlets: " + str(len(outletlist))) return outletlist, contributingset
Comments
Post a Comment