[R-sig-Geo] contour data

[Roger Bivand]

On Tue, 26 May 2009, John Callahan wrote:

> I have a DEM (2 meter spacing, IMG file format) and I'd like to create 
> contour vector data output, shapefiles would be great.  I just downloaded R 
> 2.9.0.  Looking for more info, I keep finding references to the contour 
> function in the grahics package, but that seems like it only produces plots 
> instead of output vectors. 
> I know the sp and spatstat packages are good for various types of geospatial 
> analysis.  Is one of these 'the next generation' of the other?  Would you 
> recommend one or the other specifically for contour generation from a DEM? 
> (My platform is Windows.)  Thanks.

Just roughly:

download.file("http://geomorphometry.org/data/DEM25m.zip",
   destfile="DEM25m.zip")
fname <- zip.file.extract(file="DEM25m.asc", zipname="DEM25m.zip")
file.copy(fname, "DEM25m.asc")
# to give some data to play with
library(rgdal)
dem <- readGDAL("DEM25m.asc")
# GDAL has many formats, IMG isn't very descriptive, try and see which 
# driver suits for reading to a SpatialGridDataFrame - I'm assuming one
# band only, using the first and only band next
im <- as.image.SpatialGridDataFrame(dem)
cl <- contourLines(im)
# contourLines() takes the same interval arguments as contour()
library(maptools)
SLDF <- ContourLines2SLDF(cl)
# convert to a SpatialLinesDataFrame with the contour labels as 
# attributes and export
writeOGR(SLDF, ".", "my_contours", driver="ESRI Shapefile")
mc <- readOGR(".", "my_contours")
summary(dem)
summary(mc)
image(dem, col=gray.colors(20))
plot(mc, col=terrain.colors(8), add=TRUE)

The colours in the last line are rather sleight of hand, but for 
demonstration they work this time.

Hope this helps,

Roger

>
> - John
>
>
> PS - At first pass, I'm going to give gdal_contour 
> (http://www.gdal.org/gdal_contour.html) a try.  I'd really like to see how 
> this compares with some method within R.
>
> I'm also looking at SAGA for the R interface options.  I just installed the 
> most recent QGIS with the ManageR plugin to see what tools that brings in. 
> And GRASS is always an option (either through QGIS or standalone) but I don't 
> quite understand how mapsets work yet. 
> **************************************************
> John Callahan
> Geospatial Application Developer
> Delaware Geological Survey, University of Delaware
> 227 Academy St, Newark DE 19716-7501
> Tel: (302) 831-3584  Email: john.callahan at udel.edu
> http://www.dgs.udel.edu
>
> _______________________________________________
> R-sig-Geo mailing list
> R-sig-Geo at stat.math.ethz.ch
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
>

-- 
Roger Bivand
Economic Geography Section, Department of Economics, Norwegian School of
Economics and Business Administration, Helleveien 30, N-5045 Bergen,
Norway. voice: +47 55 95 93 55; fax +47 55 95 95 43
e-mail: Roger.Bivand at nhh.no