2

I am trying to identify lakes and other water bodies in a DEM file. The source file is located here. I am using gdaldem aspect to identify contiguous pixels with zero slope change.

gdaldem aspect -co compress=lzw infile_dem.tif outfile_aspect.tif

When I plot the output file, it looks like below (please ignore the white border): enter image description here

Obviously, I am doing the exact opposite of what I need to be doing. Is there any way of inverting the operation of gdaldem aspect? The output of gdalinfo is give below:

Size is 2928, 2285
Coordinate System is:
GEOGCS["WGS 84",
    DATUM["WGS_1984",
        SPHEROID["WGS 84",6378137,298.257223563,
            AUTHORITY["EPSG","7030"]],
        AUTHORITY["EPSG","6326"]],
    PRIMEM["Greenwich",0],
    UNIT["degree",0.0174532925199433],
    AUTHORITY["EPSG","4326"]]
Origin = (81.381111111111110,22.564444444444447)
Pixel Size = (0.002083333333333,-0.002083333333333)
Metadata:
  AREA_OR_POINT=Area
Image Structure Metadata:
  INTERLEAVE=BAND
Corner Coordinates:
Upper Left  (  81.3811111,  22.5644444) ( 81d22'52.00"E, 22d33'52.00"N)
Lower Left  (  81.3811111,  17.8040278) ( 81d22'52.00"E, 17d48'14.50"N)
Upper Right (  87.4811111,  22.5644444) ( 87d28'52.00"E, 22d33'52.00"N)
Lower Right (  87.4811111,  17.8040278) ( 87d28'52.00"E, 17d48'14.50"N)
Center      (  84.4311111,  20.1842361) ( 84d25'52.00"E, 20d11' 3.25"N)
Band 1 Block=2928x1 Type=Int16, ColorInterp=Gray
  NoData Value=-32768
1

You can extract water bodies, i.e. (the most part of) pixels where slope = 0, from a DEM in two steps (rough approximation):

gdaldem slope dem.tif slope.tif
python gdal_calc.py -A slope.tif --outfile=water_bodies.tif --calc="A==0"
0

You aren't doing anything wrong. gdaldem aspect sets the aspect to the no data value when the aspect is undefined. The values you are looking for contain the datasets no data value, which is probably -9999.

  • Hi Kyle... Thanks for the pointer. I added the output of gdalfino to my original post. But how do I get the "data" of the nodota areas? Btw... I checked the documentation of gdaldem; you;re right, the minimum is -9999, but the minimum value of band1 when I look into it is -19. Not sure what's going on. – user3671 Sep 9 '13 at 4:23
  • Aspect is undefined if the slope is 0. So there is no data in those cells. -9999 is a 'flag' for 'not a valid pixel' – user10353 Sep 9 '13 at 17:50

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.