I'm trying to calculate elevation data for hiking routes, using open data (avoiding licensing constraints like Google).
I was able to read a public DEM of my country (with a 10-metres resolution) using readGDAL (from package RGDAL), and proj4string(mygrid) gives me:
"+proj=utm +zone=32 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0"
The beginning of the .asc file is:
ncols 9000
nrows 8884
xllcorner 323256,181155
yllcorner 4879269,74709
cellsize 10
NODATA_value -9999
978 998 1005 1008 1012 1016 1020 1025 .....
.....
..... 400 Megabytes of elevation values ....
.....
All that I need is to pick up from this grid the elevation data for the specific nodes of a route, so to be able to calculate elevation gain, negative slope, min/max altitude...
I bring the routes data from OpenStreetMap, using the nice package OSMAR, so the data table of my route is something like this:
RouteId NodeId lat lon
1 -13828 -8754 45.36743 7.753664
2 -13828 -8756 45.36762 7.753878
3 -13828 -8758 45.36782 7.754344
4 -13828 -8760 45.36794 7.754541
....
But I've no idea how to transform the latitude/longitude coordinates in the DEM coordinate reference system, and then how to bring the corresponding grid values (doing a sort of average of the nearest points?)
All the documentation I've found googling it is about to render grid maps, not to extract values from them.
Any help will be greatly appreciated!
Cheers, MB
P.S. The second question should be: "Having several grid-tiles, what could I do if a route is across two or more tiles? Merge them, reference both..."
As Simon said , if you are in the same projection you can use :