GIS/MEA582: Geospatial Modeling and Analysis

This is an unmaintained course material, please see current material at:

Geospatial data models

Resources:

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

Start GRASS GIS

Start GRASS - click on GRASS icon or type

grass70

In startup pannel set GIS Data Directory to path to datasets, for example on MS Windows, C:\Users\myname\grassdata.
For Project location select nc_spm_08_grass7 (North Carolina, State Plane, meters) and for Accessible mapset create a new mapset (called e.g. HW_data_models).
Click Start GRASS.

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 command line.
Hint for running most of the commands in GUI - type or paste the name of the module into the command console 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 higher 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. Here for convenience we use named region which defines both the extent and the resolution. This named region is included in the data (location) we are using but it is possible to create new named regions and use them to bookmark different study areas.

g.region region=swwake_30m -p

The -p flag for g.region is used to print the region we just set.

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 Command console, then Enter to open the GUI panel, 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

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 elev_ned10m_nn aspect=aspect_ned10m_nn

Display the resampled map by adding "aspect_ned10m_nn" to Layer Manager or in command line using:

d.rast aspect_ned10m_nn

Save the displayed map and explain what is going on in the report and how it differs from the aspect computed from the original elevation map?
To save the map, click in Map Display to on the button Save display to graphic file" or alternatively, use the following command:

d.out.file aspect_nn

Now, reinterpolate DEMs using bilinear and bicubic interpolation. Check the 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 elev_ned10m_bil aspect=aspect_ned10m_bil
r.slope.aspect elev_ned10m_bic aspect=aspect_ned10m_bic
d.rast aspect_ned10m_bil
d.rast aspect_ned10m_bic

Save the displayed map and in your report, compare the result 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_bic

Resampling to lower resolution

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

First, display the original elevation and land use maps:

d.rast elev_ned_30m
d.rast landuse96_28m

Then change the region resolution and resample elevation (which is a continuous field) and land use (which has discrete categories).
Explain selection of aggregation method. Can we use average also for landuse? What does mode mean?

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

Before the next computation, remove all map layers from the Layer Manager because we don't need to see them anymore.


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

Remove or switch off the land use, elevation and aspect maps.

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

Display census boundaries using GUI:
Add vector "census_wake2000" Selection > Feature type > boundary (switch off the other types).
Save the displayed map in Map Display click button Save display to graphic file.
Alternatively, use the following commands to control display.

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.vect census_centr icon=basic/circle fcol=green size=10 
d.vect census_wake2000 col=red  type=boundary
d.out.file vect_to_points

Convert contour lines to points (useful for computing DEM from contours):

v.to.points input=elev_ned10m_cont10m output=elev_ned_contpts type=line use=vertex

Display the "elev_ned_contpts" points vector and zoom-in to very small area to see the actual points.

d.vect elev_ned_contpts co=brown icon=basic/point size=3

Convert from vector to raster

Convert vector data to raster for use in raster-based analysis. First, adjust the computational region to resolution 200m:

g.region swwake_30m res=200 -p

Then remove all layers from the Layer Manager.

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

If you haven't done this already, add remove all other map layers from Layer Manager and add "streets_speed_30m" raster layer.
Add legend for "streets_speed_30m" raster using GUI in Map Display:
Add legend > Set Options > Advanced > 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

Convert raster lines to vector lines.

First, set the region and remove map layers from Layer Manager. Then do the conversion.

Explain why we are using r.thin module. You may want to remove all previously used layers from the Layer Manager before you start these new computations.


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

Visually compare the result with streams digitized from airphotos.

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

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.

Use raster value as category number (flag -v) and display vector polygons filled with random colors.
In GUI: Add vector > Colors > Switch on Random colors. You may want to remove all previously used layers from the Layer Manager before you start these new computations.

g.region raster=basin_50K -p

d.rast basin_50K
r.to.vect -sv basin_50K output=basin_50Kval type=area
d.vect -c basin_50Kval
d.vect streams color=blue

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

d.out.file basins