1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105
| import arcpy from arcpy import env from arcpy.sa import * import copy arcpy.CheckOutExtension("Spatial") env.overwriteOutput = True
def getLocationOfPointInGrid(rasterPath,pointPath,worksapceTempPath): env.workspace = worksapceTempPath cellXSzie = arcpy.GetRasterProperties_management(rasterPath, "CELLSIZEX") cellYSzie = arcpy.GetRasterProperties_management(rasterPath, "CELLSIZEY") leftPosition = arcpy.GetRasterProperties_management(rasterPath, "LEFT") bottomPosition = arcpy.GetRasterProperties_management(rasterPath, "BOTTOM") rstArray = arcpy.RasterToNumPyArray(arcpy.Raster(rasterPath)) rows, cols = rstArray.shape rowcell = {} colcell = {} for rowNum in xrange(rows): for colNum in xrange(cols): rowcell[(rowNum, colNum)] = [rowNum] colcell[(rowNum, colNum)] = [colNum] for rowNum in xrange(rows): for colNum in xrange(cols): rstArray[(rowNum, colNum)] = rowcell[(rowNum, colNum)][0] result = arcpy.NumPyArrayToRaster(rstArray, arcpy.Point(leftPosition[0], bottomPosition[0]), float(cellXSzie[0]),float(cellYSzie[0])) result.save("row.tif") dsc = arcpy.Describe(pointPath) coord_sys = dsc.spatialReference arcpy.DefineProjection_management("row.tif", coord_sys) for rowNum in xrange(rows): for colNum in xrange(cols): rstArray[(rowNum, colNum)] = colcell[(rowNum, colNum)][0] result = arcpy.NumPyArrayToRaster(rstArray, arcpy.Point(leftPosition[0], bottomPosition[0]), float(cellXSzie[0]),float(cellYSzie[0])) result.save("col.tif") dsc = arcpy.Describe(pointPath) coord_sys = dsc.spatialReference arcpy.DefineProjection_management("col.tif", coord_sys) inPointFeatures = pointPath inRasterList = [["row.tif", "row"], ["col.tif", "col"]] ExtractMultiValuesToPoints(inPointFeatures, inRasterList, "") cursor = arcpy.SearchCursor(inPointFeatures) rowN = '' colN = '' for rowcol in cursor: rowN = int(rowcol.row) colN = int(rowcol.col) return [rowN, colN]
def getfloodcell(waterLevel,sourcerow,sourcecol,elevationCell,rows,cols): elevationCellCopy = copy.deepcopy(elevationCell) floodbufferArray=[] floodbufferArray .append(elevationCellCopy[sourcerow,sourcecol]) while(len(floodbufferArray)!=0): rowfirst=floodbufferArray[0][0] colfirst=floodbufferArray[0][1] elevationCellCopy[(rowfirst,colfirst)][3] = True floodbufferArray.pop(0) for row in range(rowfirst - 1, rowfirst + 2): for col in range(colfirst - 1, colfirst + 2): if(row<rows and col<cols and row>=0 and col>=0): elevationinfo=elevationCellCopy[(row,col)][2] if elevationinfo<=waterLevel and elevationCellCopy[(row,col)][3]==False: elevationCellCopy[(row,col)][3] = True floodbufferArray.append((row, col)) return elevationCellCopy
worksapce_path = "D:\\sourceFlood\\temp" dem = "D:\\sourceFlood\\DEM.tif" seed = "D:\\sourceFlood\\seed.shp"
elevationCell={} cellXSzie = arcpy.GetRasterProperties_management(dem, "CELLSIZEX") cellYSzie = arcpy.GetRasterProperties_management(dem, "CELLSIZEY") leftPosition = arcpy.GetRasterProperties_management(dem, "LEFT") topPosition = arcpy.GetRasterProperties_management(dem, "TOP") bottomPosition = arcpy.GetRasterProperties_management(dem, "BOTTOM") maxdem=arcpy.GetRasterProperties_management(dem,"MAXIMUM") mindem=arcpy.GetRasterProperties_management(dem,"MINIMUM") rstArray = arcpy.RasterToNumPyArray(arcpy.Raster(dem))
rows, cols = rstArray.shape print rows,cols for rowNum in xrange(rows): for colNum in xrange(cols): elevationCell[(rowNum, colNum)] = [rowNum,colNum,rstArray.item(rowNum, colNum), False]
RowCol=getLocationOfPointInGrid(dem,seed,worksapce_path)
floodelevationCell=getfloodcell(100,RowCol[0],RowCol[1],elevationCell,rows,cols)
for rowNum in xrange(rows): for colNum in xrange(cols): rstArray[(rowNum, colNum)] = floodelevationCell[(rowNum, colNum)][3] result = arcpy.NumPyArrayToRaster(rstArray, arcpy.Point(leftPosition[0], bottomPosition[0]), float(cellXSzie[0]),float(cellYSzie[0])) result.save("D:\\sourceFlood\\result.tif")
|