> I have a jpg map of altitudinal contours in my study area which I > would like to map and use in R to calculate the altitude and slope of > GPS coordinates.
Estimating elevation from contours is not as common as it used to be because of the availability of near-global high resolution raster elevation data. What kind of spatial resolution are you after; and are you aware of the SRTM90/30 and ASTER30 data? Here are some examples using techniques that Barry alluded to. library(raster) library(maptools) # setting up some data: # contour lines for volcano cl <- ContourLines2SLDF(contourLines(volcano)) # level should be numeric cl@data$level <- as.numeric(as.character(cl@data$level)) # RasterLayer for volcano r <- flip(t(raster(volcano)), 'y') plot(r) plot(cl, add=T) # I would like to interpolate from points, by sampling from the lines, e.g. p <- spsample(cl, 1000, 'regular') plot(p) # However, as Barry mentioned spplot does not keep the attributes of the sampled lines. # But in this example it seems reasonable to do: p <- as(cl, 'SpatialPointsDataFrame') spplot(p, 'level') # now there are many ways to interpolate. See, e.g., the 'gstat' and 'automap' packages. See raster::interpolate for an example with splines. # using gstat and inverse distrance weighted interpolation: library(gstat) g <- gstat(id="level", formula = level~1, data=p, nmax=7, set=list(idp = .5)) x1 <- interpolate(r, g) par(mfrow=c(1,2)) plot(x1) plot(cl, add=T) # compare values with original: plot(x1, r) # You _can_ also first rasterize the lines; which would be better if your cells are very small relative to the inter-node distance of the lines rr <- rasterize(cl, r, field='level', progress='window') v <- data.frame(rasterToPoints(rr)) colnames(v)[3] = 'level' g2 <- gstat(id="level", formula=level~1, locations=~x+y, data=v, nmax=7, set=list(idp = .5)) x2 <- interpolate(r, g2) par(mfrow=c(1,2)) plot(x2) plot(cl, add=T) plot(x2-r) plot(cl, add=T) # for slope and aspect # RasterLayer must have a CRS. Making something up for this example: projection(x1) <- "+proj=longlat +datum=WGS84" slope <- slopeAspect(x1, out='slope') gpspoints <- cbind(0.5, 0.5) extract(slope, gpspoints) -- View this message in context: http://r-sig-geo.2731867.n2.nabble.com/Mapping-contours-from-jpg-map-tp6558992p6562497.html Sent from the R-sig-geo mailing list archive at Nabble.com. _______________________________________________ R-sig-Geo mailing list R-sig-Geo@r-project.org https://stat.ethz.ch/mailman/listinfo/r-sig-geo