I have a DEM raster with some negative values. After overlying other shapefiles, I realize that the negative values are for bodies of water (mainly lakes and bays). I want to replace the negative values with the positive values elevation values (lakes). Any suggestion to do this in python? (Con tool will not work: lakes are at different altitude levels).
I would think you could still use the CON tool if you do a little bit of vector work first. If you have a layer of the lakes, you would want them to have an attribute for the base elevation value (as shown in my embedded graphic). I drew this inside a study area polygon because you might want to avoid problems of NoData values by having the study area (green) be 0 or some value to use in the con statement.
So lets say you convert that vector work into a raster and call it studyArea. It has values of 0, 200, or 300 (in my example).
import arcpy from arcpy import env from arcpy.sa import * env.workspace = "C:/myFolder/data" outCon = Con(Raster("studyArea") > 0, "studyArea", "dem") outCon.save("C:/myFolder/output/output.img")
Assuming your lakes are polygons:
- Assign unique IDs to your lakes and do zonal statistics of DEM, maximum. In theory it should pick elevation of individual lake shores
- If output has negatives, create small buffer for lakes polygon and repeat step 1.
- Use Con(IsNull("maxdem"),"dem", "maxdem"), or Con("dem"<0),"maxdem", "dem")
Don't forget to set environment extent=dem extent
I prefer the min instead of the max, as a lake would otherwise overflow BUT beware of dams where the structure could be invisible on coarse resolution DEM's (if you expect large dams in your study area , the mean is safer)
first, extract your lakes with their bank
Con("DEM"<0, 1, 0)to get each individual lake
region group to assign an unique ID to each lake (and background)
shrink the background by one pixel
Then, get the minimum value of the DEM around the lake
Con("DEM"<0, 9999, "DEM")replace negative DEM values with very large values
Zonal statistics "MIN" gets the minimum pixel values on the lake bank
Finally, a combine min values with original DEM
Con("DEM"<0, "zonalMin", "DEM")
By using Numpy, you can get something like this:
x y value
0 0 0
0 1 0
0 2 1
1 0 0
1 1 1
where x and y are row and columns in the raster, and value is the value of that precise cell in the raster (see code below). This way we can easily spot the unwanted values, make a change to the wanted value and put that value back to the raster. Now, the question is: How can we put that value back to the raster???
Code to obtain rows, columns and values of all the raster's cells.
fOut = open('outputFile.txt', w) # Open output file
fOut.write('x' + '\t' + 'y' + '\t' + 'value\n') # Write the header
rstArray = arcpy.RasterToNumPyArray(rasterFile) # Change rasterFile to numpy array
rows, cols = rstArray.shape # Return the rows, columns
for rowNum in xrange(rows): # Loop through the rows
for colNum in xrange(cols): # Loop through the row's columns value = rstArray.item(rowNum, colNum) # Get the value at the cell fOut.write(str(rowNum) + '\t') # Write the row number fOut.write(str(colNum) + '\t') # Write the column number fOut.write(str(value) + '\n') # Write the value and a new line
This code is in the question: "How to get values of each cell in raster attribute table?" in this same website/forum.