NCSU GIS/MEA582:
Geospatial Modeling and Analysis

Geospatial data models

Resources:

• GRASS GIS overview and manual
• Recommendations and tutorial how to use GUI from the first assignment

Start GRASS GIS

grass83

To start working on your Assignment 2 right click on nc_spm_08_grass7 and select Create mapset, type in name for the new mapset, for example, assignment2, press OK. You will get a message that the Current mapset is assignment2 and any new map layer that you will create will be saved here. As mentioned in the previous assignment, the provided data are in the mapset PERMANENT - double click to see the list of available raster and vector map layers.

To save your outputs and to store external data change/create a directory where you can write, e.g. on MS Windows, C:\Users\myname\Documents\GIS582_HW\report2. Then set this directory in GRASS GIS as your working directory: Change working directory:
Settings > GRASS working environment > Change working directory > browse to and select your report2 directory
or type cd (stands for change directory) into the GUI Console and hit Enter:

cd

If you prefer to work in GUI, you should be able to find out yourself the GUI equivalents for the tasks below.
Some hints for GUI are included, but from now on, most of the instructions will be provided as commands for command line.
Hint for running most of the commands in GUI - type or paste the name of the module into the command console in the Console tab and then hit Enter to open the GUI dialog.
Read the manual page for each command you are using for the first time to learn what it is doing and what the parameters mean.

Resampling to a different resolution

Resample the given raster map to higher and lower resolution (30m->10m, 30m->100m) and compare resampling by nearest neighbor with bilinear and bicubic method.

First, set the computation region extent to our study area and set resolution to 30 meters. The computational region (region for short) is set using g.region module, the -p flag is used to print the region coordinates and resolution. Here we use a predefined region which is included in our data set.

g.region region=swwake_30m -p

Zoom to computational region by clicking on the relevant icon on top of the Map Displaay to make sure your display shows the region that we are working with.

Then we display the 30m resolution NED elevation raster.

d.rast elev_ned_30m

To resample it to 10m resolution, first set the computational region to resolution 10m, then resample the raster using the nearest neighbor method.
Hint: To open the r.resamp.interp in GUI, type or paste the module name into the Console tab, then Enter to open the GUI dialog, don't forget to set the method to nearest under Optional tab.

g.region res=10 -p
r.resamp.interp elev_ned_30m out=elev_ned10m_nn method=nearest

The resulting elevation map should display automatically, if not, display the resampled map by adding "elev_ned10m_nn" to Layer Manager in case you don't have it in the Layer Manager already.
Alternatively, use in command line the following:

d.rast elev_ned10m_nn

The elevation map "elev_ned10m_nn" looks the same as the original one, so now check the resampled elevation surface using the aspect map:

r.slope.aspect elevation=elev_ned10m_nn aspect=aspect_ned10m_nn

You should see the aspect map derived from the resampled elevation map in your display, or you can add "aspect_ned10m_nn" to Layer Manager or in command line using:

d.rast aspect_ned10m_nn

To save the map, click in Map Display on the button Save display to graphic file" or alternatively, use the following command to save your map in the current working directory:

d.out.file aspect_nn

Now, reinterpolate DEMs using bilinear and bicubic interpolation. Check the structure of interpolated elevation surfaces using aspect maps.

r.resamp.interp elev_ned_30m out=elev_ned10m_bil meth=bilinear
r.resamp.interp elev_ned_30m out=elev_ned10m_bic meth=bicubic
r.slope.aspect elevation=elev_ned10m_bil aspect=aspect_ned10m_bil
r.slope.aspect elevation=elev_ned10m_bic aspect=aspect_ned10m_bic
d.rast aspect_ned10m_bil
d.rast aspect_ned10m_bic

Save the displayed maps and in your report, compare the results with the previously computed nearest neighbor result.
In Map Display click button Save display to graphic file, or use the following in the command line:

d.out.file aspect_bil

Resample to lower resolution (30m -> 100m).

First, remove all map layers from the Layer Manager because we don't need to see them anymore and then display the original elevation and land use maps:.

d.erase
d.rast elev_ned_30m
d.rast landuse96_28m

g.region res=100 -p
r.resamp.stats elev_ned_30m out=elev_new100m_avg method=average
d.rast elev_new100m_avg
d.out.file elev_100m

d.erase
r.resamp.stats landuse96_28m out=landuse96_100m method=mode
d.rast landuse96_100m
d.out.file lenduse_100m

Converting between vector data types

Convert census blocks polygons to points using their centroids (useful for interpolating a population density trend surface):

v.to.points census_wake2000 type=centroid out=census_centr use=vertex

Note that in both command line and GUI you must either enter the full path to the file you are saving the image in, or you must know the current working directory.

d.erase
d.vect census_centr icon=basic/circle fill_color=green size=10
d.vect census_wake2000 color=red fill_color=none
d.legend.vect
d.out.file vect_to_points

Convert from vector to raster

Convert vector data to raster for use in raster-based analysis. First, remove all previous map layers from Layer Manager and adjust the computational region to resolution 200m:

d.erase
g.region swwake_30m res=200 -p

Convert vector points "schools" to raster. As value for raster use attribute column "CORECAPACI" for core capacity.
To add legend in GUI use Add map elements > Show/hide legend and select "schools_cap_200m".

d.vect schools_wake 
v.info -c schools_wake
v.to.rast schools_wake out=schools_cap_200m use=attr attrcol=CORECAPACI type=point
d.rast schools_cap_200m
d.vect streets_wake co=grey
d.legend schools_cap_200m at=70,30,2,6

Now convert lines in "streets" vector to raster. Set the resolution to 30m and use speed limit attribute.

g.region res=30 -p
v.to.rast streets_wake out=streets_speed_30m use=attr attrcol=SPEED type=line

Raster "streets_speed_30m" will be added to the Layer manager and displayed in the Map Display.
Add legend for "streets_speed_30m" raster using GUI in Map Display:
Add legend > Set Options > Subset > List of discrete cat numbers and type in speed limits 25,35,45,55,65; move legend with mouse as needed.

Alternatively, use the following commands:

d.rast streets_speed_30m 
d.legend streets_speed_30m at=5,30,2,5 use=25,35,45,55,65 

Save the displayed map. In Map Display click button Save display to graphic file, or use the following.

d.out.file vect_to_rast

Convert from raster to vector data

First, remove map layers from Layer Manager and set the region. Hint: right click on the raster legends to select Remove legend. Then do the conversion.

d.erase
g.region raster=streams_derived -p
d.rast streams_derived
r.thin streams_derived output=streams_derived_t
r.to.vect streams_derived_t output=streams_derived_t type=line

Explain why we are using r.thin module.

Visually compare the result with streams digitized from airphotos. To make space for the vector legend you may need to stretch your Map Display window a little wider.

d.vect streams_derived_t color=blue
d.vect streams color=red
d.legend.vect

Save the displayed map (in Map Display click button Save display to graphic file).

d.out.file streams_compare

Convert raster areas representing basins to vector polygons.

First, remove map layers from Layer Manager and set the region. To convert, use raster values as category numbers (flag -v) and display boundaries of vector polygons.
Verify the basin boundaries by displaying them together with streams - the stream networks should "fit" within the basin boundaries.

d.erase
g.region raster=basin_50K -p
r.to.vect -sv basin_50K output=basin_50Kval type=area
d.erase
d.rast basin_50K
d.vect basin_50Kval type=boundary width=2
d.vect streams color=blue

Save the displayed map either using GUI or using the following in case you are working in the command liner:.

d.out.file basins