LandScript - Tutorial

This tutorial will show you how to use LandScript to perform tasks ranging from opening and saving files to sophisticated spatial analysis using map algebra. It is highly recommended that you first familiarlise yourself with LandSerf before completing this tutorial, perhaps by completing the introductory LandSerf tutorial.

For the moment, each tutorial step corresponds approximately to a stage in the LandSerf tutorial. The scripts are commented and should be self-explanatory.

  1. Opening, reprojecting and saving a raster map
  2. Editing raster metadata
  3. Processing a landcover raster
  4. Simple morphometric analysis
  5. Advanced morphometric analysis

Step 1 - Opening, reprojecting and saving a raster map.

Explanation to follow (see the introductory LandSerf tutorial for details). 

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20

# Script to convert binary image file into a UTM LandSerf file.
# v1.2 Jo Wood 21st April, 2009.
# -------------------------------------------------------------
version(1.0);

# Store the base directory in a variable to allow easy changing.
baseDir = "c:/Program Files/LandSerf/data/";

# Open the .bil file and store the raster as a variable.
dem = open(baseDir & "lincolnNED.bil");

# Add the projection metadata to the DEM.
edit(dem,"projection","LatLong");
edit(dem,"ellipsoid","wgs 84");

# Create a new Universal Transverse Mercator (UTM) reprojected DEM.
dem = reproject(dem_,"UTM","true",30,30);

# Save the edited raster as a LandSerf file.
save(dem, baseDir & "landscriptTutorial/step1.srf");

Step 2 - Editing raster metadata.

Explanation to follow (see the introductory LandSerf tutorial for details). 

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29

# Script to edit the metadata of a spatial object.
# v1.2 Jo Wood 21st April, 2009.
# -------------------------------------------------------------
version(1.0);

# Store the base directory in a variable to allow easy changing.
baseDir = "c:/Program Files/LandSerf/data/landscriptTutorial/";

# Open the file from the previous step and store as a variable.
dem = open(baseDir & "step1.srf");

# Give the raster a new colour table.
colouredit(dem,"land3");

# Guess at the water level by finding modal elevation
maxFreq = info(dem,"modef");
if (maxFreq > 1000);
{
  # If there are more than 1000 cells of same height, colour them blue
  waterLevel = info(dem,"mode");
  colouredit(dem,"addRule",waterLevel & " 141 141 166 (D)");
}

# Update title and notes.
edit(dem,"title","Lincoln 30m DEM");
edit(dem,"notes","Elevation data from USGS National Elevation Dataset (NED) available from http://seamless.usgs.gov");

# Save the edited raster as a LandSerf file.
save(dem, baseDir & "step2.srf");

Step 3 - Processing a landcover raster.

Explanation to follow (see the introductory LandSerf tutorial for details). 

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35

# Script to process a landcover raster.
# v1.2 Jo Wood 21st April, 2009.
# -------------------------------------------------------------
version(1.0);

# Store the base directory in a variable to allow easy changing.
baseDir = "c:/Program Files/LandSerf/data/";

# Open the landcover file, reproject it to UTM.
landcover = open(baseDir & "lincolnLandcover.bil");
edit(landcover,"projection","LatLong");
edit(landcover,"ellipsoid","wgs 84");
landcover = reproject(landcover_,"UTM","false",30,30);

# Update the landcover metadata.
edit(landcover,"title","Lincoln landcover");
colouredit(landcover,"file",baseDir & "landcover.ctb");
edit(landcover,"attributes",baseDir & "landcover.atr");
edit(landcover,"attCol",3);
edit(landcover,"type","other");

# Find average height of all forested areas.
dem = open(baseDir & "step2.srf");
forestHeight = new(landcover);
forestHeight = ifelse((landcover >= 40) and (landcover <=49),dem,null());
echo("Average height of forest is " & round(info(forestHeight,"mean")) & "m.");

# Find average height of scrubland areas.
dem = open(baseDir & "step2.srf");
scrubHeight = new(landcover);
scrubHeight = ifelse(landcover == 51,dem,null());
echo("Average height of scrubland is " & round(info(scrubHeight,"mean")) & "m.");

# Save the edited landcover as a LandSerf file.
save(landcover, baseDir & "landscriptTutorial/step3.srf");

Step 4 - Simple morphometric analysis.

Explanation to follow (see the introductory LandSerf tutorial for details). 

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48

# Script to perform some simple morphometric analysis on a DEM.
# v1.2 Jo Wood 21st April, 2009.
# -------------------------------------------------------------
version(1.0);

# Store the base directory in a variable to allow easy changing.
baseDir = "c:/Program Files/LandSerf/data/landscriptTutorial/";

# Open the DEM and landcover raster from the previous steps.
dem = open(baseDir & "step2.srf");
landcover = open(baseDir & "step3.srf");

# Create a shaded relief image saved as a file.
edit(dem,"vExag",3);
draw(baseDir & "relief.jpg",dem,landcover,"null","null","relief");

# Output some relevant DEM metadata
showMetadata(dem,"elevation");

# Calculate and output slope of DEM
slope = surfparam(dem,"slope");
draw(baseDir & "slope.jpg",dem,slope,"null","null","relief");
showMetadata(slope,"slope");

# Create images of the other surface parameters.
draw(baseDir & "aspect.jpg",dem,surfparam(dem,"aspect"),"null","null","relief");
draw(baseDir & "profc.jpg", dem,surfparam(dem,"profc"),"null","null","relief");
draw(baseDir & "planc.jpg", dem,surfparam(dem,"planc"),"null","null","relief");

# Finally create a surface feature map and save as a LandSerf file.
features = surfparam(dem,"feature");
save(features, baseDir & "LincolnFeatures.srf");


# Function to display some spatial metadata.
# Params: raster Raster from which to extract metadata.
#         name Name of raster.
function showMetadata(raster,name)
{
  echo("\n" & name);
  echo("East-west extent: "& (info(raster,"E")-info(raster,"W"))/1000 & "km.");
  echo("North-south extent: "& (info(raster,"N")-info(raster,"S"))/1000 & "km.");
  echo(name&" range from "& info(raster,"min")&" to "&info(raster,"max")&".");
  echo("Mean "&name&": "& info(raster,"mean"));
  echo("Modal "&name&" value: "& info(raster,"mode"));
  echo("Standard deviation of "&name&": "& info(raster,"stdev"));
  echo("Skewness of "&name&": "& info(raster,"skew") & "\n");
}

Step 5 - Advanced morphometric analysis.

Explanation to follow (see the introductory LandSerf tutorial for details). 

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28

# Script to perform more advanced morphometric analysis on a DEM.
# v1.2 Jo Wood 21st April, 2009.
# -------------------------------------------------------------
version(1.0);

# Store the base directory in a variable to allow easy changing.
baseDir = "c:/Program Files/LandSerf/data/landscriptTutorial/";

# Open the DEM from the earlier step.
dem = open(baseDir & "step2.srf");

# Calculate 11x11 cell smothed elevation and draw output.
winSize = 11*info(dem,"xRes");
echo("Smoothing elevation using "&winSize&"m x "&winSize&"m window");
draw(baseDir & "elev11.jpg", surfparam(dem,"elev",11),"null","null","null","relief");

# Calculate 65x65 cell smothed elevation and draw and save output.
winSize = 65*info(dem,"xRes")/1000;
echo("Smoothing elevation using "&winSize&"km x "&winSize&"km window");
elev65 = surfparam(dem,"elev",65);
draw(baseDir & "elev65.jpg",elev65,"null","null","null","relief");
save(elev65,baseDir&"elev65.srf");

# Create a 5x5 feature network and metric surface network
feat5 = surfparam(dem,"netfeature",5,1);
msn = surfnetwork(dem,feat5);
draw(baseDir & "featnetwork.jpg",dem,"null",msn,"null","relief");
save(msn,baseDir&"msn.vec");