NCSU GIS/MEA582:
Geospatial Modeling and Analysis

Spatial interpolation and approximation II: splines

Resources:

Text files with color rules:

Start GRASS GIS

Start GRASS - click on GRASS icon or type
grass72

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_interpolation_2).
Click Start GRASS.

Change working directory:
Settings > GRASS working environment > Change working directory > select/create any directory
or type cd (stands for change directory) into the GUI Console and hit Enter:

cd
Download all text files with color rules (see above) to the selected directory. Now you can use the commands from the assignment requiring the text file without the need to specify the full path to the file.

Interpolate elevation raster from points using splines with different tension

Compute aspect simultaneously with interpolation and evaluate impact of tension by using tension=40 (default), tension=10 and tension=160.
g.region rural_1m res=1 -p
v.surf.rst input=elev_lid792_randpts elevation=elev_rstdef_1m zcolumn=value aspect=asp_rstdef_1m segmax=30 npmin=140
v.surf.rst input=elev_lid792_randpts elevation=elev_rstt10_1m aspect=asp_rstt10_1m zcolumn=value tension=10 segmax=30 npmin=140
v.surf.rst input=elev_lid792_randpts elevation=elev_rstt160_1m aspect=asp_rstt160_1m zcolumn=value tension=160 segmax=30 npmin=140

Compare the interpolated elevation surfaces using aspect maps.
Change the aspect color table to grey aspect.
Save images for your report.

r.colors asp_rstdef_1m color=aspect
r.colors asp_rstt10_1m color=aspect
r.colors asp_rstt160_1m color=aspect

d.rast elev_rstdef_1m
d.rast asp_rstdef_1m
d.out.file asp_rst_t40
d.rast asp_rstt10_1m
d.out.file asp_rst_t10
d.rast asp_rstt160_1m
d.out.file asp_rst_t160
d.vect elev_lid792_randpts size=1 color=red

Or use 3D views of elev_rstdef_1m, elev_rstt10_1m, elev_rstt160_1m, make sure you switch off the aspect rasters and save the 3 images for your report.

Compute elevation raster and deviations vector point map

For different values of smoothing compare deviation stats for smoothing 0.1 and 10.
Find root mean square deviation rmse.
v.surf.rst input=elev_lid792_randpts elevation=elev_rstdef_1mb zcolumn=value smooth=0.1 deviations=elev_rstdef_devi segmax=30 npmin=140
v.build elev_rstdef_devi
v.surf.rst input=elev_lid792_randpts elevation=elev_rstsm10_1mb zcolumn=value smooth=10 deviations=elev_rstsm10_devi segmax=30 npmin=140
v.build elev_rstsm10_devi
v.info -c elev_rstdef_devi
v.univar elev_rstdef_devi column=flt1 type=point
r.info elev_rstdef_1mb
v.info -c elev_rstsm10_devi
v.univar elev_rstsm10_devi column=flt1 type=point
r.info elev_rstsm10_1mb

Compute and display deviations maps using same color table. You need to use custom color table to see the results well.
Note that we are interpolating here the deviations, not the given elevations.

v.surf.rst input=elev_rstdef_devi elevation=elev_rstdef_devi zcolumn=flt1 segmax=30 npmin=140
v.surf.rst input=elev_rstsm10_devi elevation=elev_rstsm10_devi zcolumn=flt1 segmax=30 npmin=140

Apply the downloaded color table deviations_color.txt to the deviation raster.
Optionally, to view the results in 3D use "elev_rstdef_1mb" for elevation (switch off everything else) and drape the deviations maps as color.

r.colors elev_rstsm10_devi rules=deviations_color.txt
r.colors elev_rstdef_devi raster=elev_rstsm10_devi

d.rast elev_rstdef_devi
d.rast elev_rstsm10_devi
d.legend elev_rstsm10_devi at=2,50,2,6
d.out.file elev_rstsm10_devi

Compute predictive error of interpolation

Compute predictive error of interpolation for each point using cross-validation (no raster output, only points with pred. errors).
v.surf.rst -c input=elev_lid792_randpts zcolumn=value cvdev=elev_rstdef_cv npmin=120 segmax=35
v.build elev_rstdef_cv
v.univar elev_rstdef_cv column=flt1 type=point

Compute raster map of predictive errors and identify locations where the sampling is inadequate.
Optionally, to view the result in 3D use "elev_rstdef_1mb" for elevation (switch off everything else) and drape the crossvalidation map "elev_rstdef_cv" as color.

v.surf.rst input=elev_rstdef_cv elevation=elev_rstdef_cv zcolumn=flt1
r.colors elev_rstdef_cv raster=elev_rstsm10_devi

d.rast elev_rstdef_cv
d.vect elev_rstdef_cv size=2
d.legend elev_rstdef_cv at=2,50,2,6
d.out.file elev_rstdef_cv

Interpolate precipitation with influence of topography

Set the 3D region (read the man page for g.region).
We set tbres to high value - we have just a single level because we are not computing the 3D raster (see lecture for more details).
g.region raster=elev_state_500m -p
g.region t=2000 b=0 tbres=2000 res3=500 -p3

Compute precipitation raster map without influence of elevation (with segmax=700 segmentation is not performed so interpolation function is computed using all points at once).
We will use mask during the interpolation.

r.mask raster=ncmask_500m
v.info -c precip_30ynormals
v.surf.rst input=precip_30ynormals elevation=precip_annual_500m zcolumn=annual segmax=700

Use the downloaded the color table precip_color.txt.
Zoom to computational region when displaying the result.

r.colors precip_annual_500m rules=precip_color.txt

d.rast precip_annual_500m
d.legend precip_annual_500m at=2,30,2,5 range=970,2400

Compute precipitation raster map with elevation.

There is both 3D voxel output and 2D raster output - we want the 2D raster output (cross_output).
Optionally to view the results in 3D, switch off everything except for elev_state_500m and precip_30ynormals_3d,
switch to 3D, set(type in) viewer height at 300000, z-exag at 6, fine res=1,
use precip_anntopo_500m for color, set icon size for points - sphere, 5000.
Display the result and save the image for the report.

v.info -c precip_30ynormals_3d
v.vol.rst input=precip_30ynormals_3d cross_input=elev_state_500m cross_output=precip_anntopo_500m maskmap=elev_state_500m wcolumn=annual zscale=90 segmax=700
r.colors precip_anntopo_500m raster=precip_annual_500m
d.rast precip_anntopo_500m

Try to explain how was elevation used for the precipitation raster interpolation.

After you are finished, remove mask.

r.mask -r