How to create an elevation profile: Difference between revisions

From GRASS-Wiki
Jump to navigation Jump to search
(wikification)
 
(11 intermediate revisions by 5 users not shown)
Line 1: Line 1:
=== Q: ===
== '''Q:''' How to create an elevation profile from intersecting shapefiles ==
How to create an elevation profile from intersecting shapefiles


=== A: ===
'''A:''' You want to create an elevation profile for a walking track described in one file, when the only available elevation data is contained in another file describing the terrain contours for the same region as the track.
You want to create an elevation profile for a walking track described in one file, when the only available elevation data is contained in another file describing the terrain contours for the same region as the track.
It is assumed both files are available for GRASS GIS in shapefile format.
It is assumed both files are available to GRASS in shapefile format.
In order to ascertain the elevation along the length of the track, you must calculate the elevation at the points where the track intersects with the contour lines, whose elevation is known.
In order to ascertain the elevation along the length of the track, you must calculate the elevation at the points where the track intersects with the contour lines, whose elevation is known.
GRASS can handle data in both vector and raster form. In this instance we will start with the two shapefiles in vector form, one representing the coordinates of the track, the other the coordinates and the elevation of the contour lines. We will then convert them to raster form in order to establish the elevation at the points where the track intersects the contours, and finally produce the track elevation profile data.
GRASS can handle data in both vector and raster form. In this instance we will start with the two shapefiles in vector form, one representing the coordinates of the track, the other the coordinates and the elevation of the contour lines. We will then convert them to raster form in order to establish the elevation at the points where the track intersects the contours, and finally produce the track elevation profile data.


* ''' Importing the shapefiles into GRASS'''
=== Importing the shapefiles into GRASS ===


  {{cmd|v.in.ogr}} -o dsn=/home/george/GRASSDATA/Contour.shp \
  {{cmd|v.in.ogr}} -o dsn=/home/george/GRASSDATA/Contour.shp \
     output=contours_vect min_area=0.0001 type=line snap=-1
     output=contours_vect min_area=0.0001 type=line snap=-1
 
  v.in.ogr -o dsn=/home/george/GRASSDATA/Track_cl.shp \
  v.in.ogr -o dsn=/home/george/GRASSDATA/Track_cl.shp \
     output=tracks_vect min_area=0.0001 type=line snap=-1
     output=tracks_vect min_area=0.0001 type=line snap=-1


At this point you can optionally visualize the map display using the GRASS Display Manager which opens up automatically when you start GRASS.
At this point you can optionally visualize the map display using the GRASS Display Manager which opens up automatically when you start GRASS.
You can also check that the elevation values have been correctly imported using:
{{cmd|v.db.select}}
using elevation as the property to list.


* '''Set the computational region resolution and convert vector to raster'''
=== Set the computational region resolution and convert vector to raster ===


In order to convert the vector files to raster we will use the {{cmd|v.to.rast}} program.
In order to convert the vector files to raster we will use the {{cmd|v.to.rast}} program.
However before you proceed you MUST ensure that the computational region is set to the correct resolution. The v.to.rast documentation has an innocent looking warning which reads as follows:
However before you proceed you MUST ensure that the [[computational region]] is set to the correct resolution. The {{cmd|v.to.rast}} documentation has an innocent looking warning which reads as follows:


''"v.to.rast will only affect data in areas lying inside the boundaries of the current geographic region. Before running v.to.rast, the user should therefore ensure that the current geographic region is correctly set and that the region resolution is at the desired level."''
''"v.to.rast will only affect data in areas lying inside the boundaries of the current geographic region. Before running v.to.rast, the user should therefore ensure that the current geographic region is correctly set and that the region resolution is at the desired level."''


Do not ignore this warning as the most likely outcome will be a useless raster map consisting of large, chunky squares for each pixel.
Do not ignore this warning as the most likely outcome will be a probably useless raster map consisting of large, chunky squares for each pixel.
To change the resolution of the computational region use the g.region command with an appropriate resolution. Some experimentation might be required. For example use:
To change the resolution of the computational region use the {{cmd|g.region}} command with an appropriate resolution. Some experimentation might be required. For example use:


  g.region -s -p res=10
  g.region -s -p res=10
Line 33: Line 35:


  v.to.rast input=contours_vect output=contours_rast use=attr \
  v.to.rast input=contours_vect output=contours_rast use=attr \
     type=line layer=1 column=ELEVATION value=1 rows=4096
     type=line layer=1 column=ELEVATION value=1
 
  v.to.rast input=tracks_vect output=tracks_rast use=val \
  v.to.rast input=tracks_vect output=tracks_rast use=val \
     type=line layer=1 value=1 rows=4096
     type=line layer=1 value=1


* '''Creating the elevation profile'''
=== Creating the elevation profile ===


  r.surf.contour input=contours_rast output=elevmap_raster
  r.surf.contour input=contours_rast output=elevmap_raster


Few seconds later it should have completed its job without errors. Then go ahead with:
Few seconds later it should have completed its job without errors. (See the [[Contour lines to DEM]] wiki page for details.)
 
Then go ahead with:


  r.profile -i input=elevmap_rast output=- null=*
  r.profile -i input=elevmap_rast output=- null=*
Line 50: Line 54:


The *profile* parameter can be set to comma separated geographic
The *profile* parameter can be set to comma separated geographic
coordinates for profile line endpoints. (The interactive flag (*-i*) overrides this option). Alternatively the coordinate pairs can be piped
coordinates for profile line endpoints. (The interactive flag (*-i*) overrides this option). Alternatively the coordinate pairs can be piped from stdin, one comma separated pair per line.
from stdin, one comma separated pair per line.


You can get the coordinate pairs from your track vector with
You can get the coordinate pairs from your track vector with
''v.out.ascii'', and use the ''fs='', option to make the output file
''v.out.ascii'', and use the ''fs='', option to make the output file
"comma-separated".  
"comma-separated".


* '''An alternative''' has been suggested (not tested yet):
=== '''An alternative''' has been suggested (not tested yet) ===
Use the command ''v.drape'' to get profiles on raster maps along lines from vector maps. The command to be executed is something like.
Use the command ''v.drape'' to get profiles on raster maps along lines from vector maps. The command to be executed is something like (note that type=point is required to export to ascii).


  v.drape input=perfil out=perfil_n_mg type=point rast=eigen \
  v.drape input=perfil out=perfil_n_mg type=point rast=eigen \
Line 66: Line 69:


  v.out.ascii perfil_n_mg > perfil_n_mg.txt
  v.out.ascii perfil_n_mg > perfil_n_mg.txt
For lines:
v.drape input=perfil out=perfil_n_mg type=line rast=eigen \
    method=cubic
v.to.points -i input=perfil_n_mg output=perfil_n_mg_point
v.out.ascii perfil_n_mg_point > perfil_n_mg.txt
An example using this method is described here: http://casoilresource.lawr.ucdavis.edu/drupal/node/375


[[Category:FAQ]]
[[Category:FAQ]]
[[Category:HowTo]]

Latest revision as of 01:59, 11 February 2016

Q: How to create an elevation profile from intersecting shapefiles

A: You want to create an elevation profile for a walking track described in one file, when the only available elevation data is contained in another file describing the terrain contours for the same region as the track. It is assumed both files are available for GRASS GIS in shapefile format. In order to ascertain the elevation along the length of the track, you must calculate the elevation at the points where the track intersects with the contour lines, whose elevation is known.

GRASS can handle data in both vector and raster form. In this instance we will start with the two shapefiles in vector form, one representing the coordinates of the track, the other the coordinates and the elevation of the contour lines. We will then convert them to raster form in order to establish the elevation at the points where the track intersects the contours, and finally produce the track elevation profile data.

Importing the shapefiles into GRASS

v.in.ogr -o dsn=/home/george/GRASSDATA/Contour.shp \
   output=contours_vect min_area=0.0001 type=line snap=-1

v.in.ogr -o dsn=/home/george/GRASSDATA/Track_cl.shp \
   output=tracks_vect min_area=0.0001 type=line snap=-1

At this point you can optionally visualize the map display using the GRASS Display Manager which opens up automatically when you start GRASS. You can also check that the elevation values have been correctly imported using: v.db.select using elevation as the property to list.

Set the computational region resolution and convert vector to raster

In order to convert the vector files to raster we will use the v.to.rast program. However before you proceed you MUST ensure that the computational region is set to the correct resolution. The v.to.rast documentation has an innocent looking warning which reads as follows:

"v.to.rast will only affect data in areas lying inside the boundaries of the current geographic region. Before running v.to.rast, the user should therefore ensure that the current geographic region is correctly set and that the region resolution is at the desired level."

Do not ignore this warning as the most likely outcome will be a probably useless raster map consisting of large, chunky squares for each pixel. To change the resolution of the computational region use the g.region command with an appropriate resolution. Some experimentation might be required. For example use:

g.region -s -p res=10

which will create a resolution where every pixel on the map represents 10 metres on the ground. Then from the Map Display window icons bar click the "zoom to..." icon (a magnifying glass over a small map) and from the drop menu select "Set computational region extents to match display". That's it. You can now safely convert your vector files to raster, as follows:

v.to.rast input=contours_vect output=contours_rast use=attr \
   type=line layer=1 column=ELEVATION value=1

v.to.rast input=tracks_vect output=tracks_rast use=val \
   type=line layer=1 value=1

Creating the elevation profile

r.surf.contour input=contours_rast output=elevmap_raster

Few seconds later it should have completed its job without errors. (See the Contour lines to DEM wiki page for details.)

Then go ahead with:

r.profile -i input=elevmap_rast output=- null=*

The -i options allows interactive following of the track on the display map by clicking with the left mouse button along the section of the track you are interested in, and a two column list of values will be printed on stdout, containing the distance between points and the elevation. Use any spreadsheet application to plot the values. Alternatively, you can omit the -i option and pipe the list of coordinates to r.profile as follows:

The *profile* parameter can be set to comma separated geographic coordinates for profile line endpoints. (The interactive flag (*-i*) overrides this option). Alternatively the coordinate pairs can be piped from stdin, one comma separated pair per line.

You can get the coordinate pairs from your track vector with v.out.ascii, and use the fs=, option to make the output file "comma-separated".

An alternative has been suggested (not tested yet)

Use the command v.drape to get profiles on raster maps along lines from vector maps. The command to be executed is something like (note that type=point is required to export to ascii).

v.drape input=perfil out=perfil_n_mg type=point rast=eigen \
   method=cubic

And to convert the profiles to ascii files, like tables, you can run:

v.out.ascii perfil_n_mg > perfil_n_mg.txt

For lines:

v.drape input=perfil out=perfil_n_mg type=line rast=eigen \
   method=cubic

v.to.points -i input=perfil_n_mg output=perfil_n_mg_point

v.out.ascii perfil_n_mg_point > perfil_n_mg.txt

An example using this method is described here: http://casoilresource.lawr.ucdavis.edu/drupal/node/375