Introduction to Matlab and producing Digital Elevation Models

Class goals

1. Have a general familiarity with the Matlab programming and scripting environment
2. Understand how a digital elevation model is created from scattered, sparse and dense data

Lessons

• Introduction to Matlab--go through Using Matlab For Geological Applications by Hilley and Arrowsmith. In particular follow through the final example of processing the topographic data from the Tempe Town Lake area.
  
• Look over the Matlab documentation site: http://www.mathworks.com/access/helpdesk/help/techdoc/matlab.html. It is much better and comprehensive than the online help.
• Read these notes on DEM generation: Notes_on_Lidar_interpolation. Note that it applies specifically to high density elevation surveys, but it applies in general to DEMs. The first of the examples below is about interpolation and the second about local binning.

Exercises and review--all due September 20 at the beginning of class

I. Topography along the San Andreas Fault--Interpolation of sparse elevation data

1. Download and review this script:
   topographic_map_script.m
   
2. Download these data:
   vanmatre_data_3col.txt
   Put them in the same directory and make sure that your Matlab session is in that directory.
3. The data come from a topographic survey that I made along the San Andreas Fault in the 1990.
   It is in this part of California: Overview from http://maps.google.com.
   Zoom to the site:
   Van Matre Ranch in Carrizo Plain, California
   Here is a view of the site on the ground:
   
4. The main parts of the topographic_map_script.m that you can change are at the top:
   

   
   %CHANGE THINGS BETWEEN HERE
   %set parameters:
   %number of x and y nodes for grid
   %make these smaller if you get an error about exceeding the max array size
   nx = 100;
   ny = nx;
   %contour interval (this is in the units of the elevation data):
   ci = 2; 
   
   

   Remember that the % sign is the beginning of a comment that the program will ignore
   
   
5. The following are the parts to turn in:
   
   1. Explore the effect of changing resolution of the resultant DEM by changing the grid resolution. Use nx = 10, nx = 25, nx = 100, and nx = 250. Print out figures 2, 3, and 6 for each and write a paragraph describing the changes as you increased the grid resolution. Which do you think is optimal and why (obviously too coarse does not sufficiently represent the landscape, and too fine takes too long to run and is overkill)?
   2. Explore the effect of changing the contour interval. Change ci = 2 meters to others like 0.1 m to 10 m. What do you think is the optimal contour interval for this map? Print out the version of the linear interpolation (figure 3) that you think is best.
   3. Can you tell the difference between the linear (figure 3) and cubic interpolation (figure 5)? Go to http://www.mathworks.com/access/helpdesk/help/techdoc/matlab.html, and type griddata into the function name box. Explain in your own words what our script is doing with this command and why there might be difference between the linear and cubic methods. On what lines is it used in the topographic_map_script.m?

II. Topography in North coastal California--local binning of dense elevation data

1. Download this file of points from northern California: Point_Arena_3col.txt. Note these data are in feet!. Download this new script: local_binning_script.m. Put them in the same directory and make sure that your Matlab session is in that directory. Review the script and follow the comments to see generally what is going on.
   The topography are from north of the Point Arena Lighthouse area and is data gathered by Airborne Laser Swath Mapping (for more information see http://lidar.asu.edu). Here is a link to an amazing set of photographs along the California coastline that includes the Point Arena area and our site: http://www.californiacoastline.org. Note the black car on the road parked above the bend in the river. Here is the location in maps.google.com: location.
2. Note that in this script specifies a grid resolution in feet (dx = 50; %grid resolution would be 50 foot spacing between grid nodes), rather than the number of grid points in the prior case. The search radius is computed as search_radius = sqrt(2).*dx./2, which for a 10 foot grid spacing is 7.071 feet. This figure shows how the search radius (red circle) compares with the grid nodes (plus symbols): searching.gif.
3. Run the script with dx = 50, dx = 10, and dx = 5. That will increase the grid resolution progressively from 50 feet to 10 feet to 5 feet. It may take a while for the 5 foot version to run (took about 8 minutes on my machine. If it goes for longer than about half an hour, push control c and try a different machine). You should be able to see the river (because the laser shots did not return from the water) and its bend as is manifest in the photograph (http://www.californiacoastline.org).
4. The following are the parts to turn in:
   
   1. Print out the 3 sets of figures that go with the increasing resolution. Label them with the grid resolution accordingly.
   2. Using numbers (1, 2, 3, etc.), label at least 5 features on the photograph (print it out): http://www.californiacoastline.org). and label one of your maps (perhaps the highest resolution version) with the same features. On a separate piece of paper, say what those features are.
   3. Answer these questions:
      1. How does the depiction of the landscape improve as the resolution increases?
      2. How does the point density (figure 3) change with increasing resolution? What are typical densities for each resolution? If the computer could handle it, how fine of a resolution do you think it would be appropriate to go to?
      3. To what do you attribute the differences between the different binning approaches (min, max, mean, IDW)?
      4. Comparing the photo and your maps, what parts of the landscape are reasonably well represented by the DEMs you have made and what parts don't seem to show up as well as you expected?