Csci 210 Lab: Visualizing Terrains

(Laura Toma, developed with the support of NSF award no. 0728780.)

Overview

Test grids:

Vast amounts of terrain data are available from remote sensing technology in the form of aerial photographs. Aerial photographs are a great tool for viewing and are becoming now standard with tools like Google maps or Yahoo maps and others. However, for modeling purposes, a more useful type of data is elevation data: essentially, a terrain elevation model consists of a bunch of points in the terrain and their respective elevation. Terrain elevation data can be used to model and understand a wealth of natural processes in geosciences, and beyond.

Elevation data has traditionally been obtained by interpolating cartographic maps; elevation data can also be estimated based on radar images: Synthetic Aperture Radars (SAR) produce high resolution images of the Earth's surface providing useful information about the physical characteristics of the ground and of the vegetation canopy, such as surface roughness, soil moisture, tree height and bio-mass estimates. By combining two or more SAR images of the same area it is possible to generate elevation maps and surface change maps with unprecedented precision and resolution. This technique is called SAR interferometry. More recently elevation data is available (directly) from technologies like LIDAR.

In this lab you will develop an application to read terrain elevation data from a file and render (visualize) it with various color maps. These is the basic function in a Geographic Information System (GIS).

Representing terrains

Terrains are represented by sampling points on the terrain and recording the geographical coordinates {x, y} of the point and the coresponding elevation {z} of the terrain at that point. Thus a terrain is represented as a cloud of points {x, y, z} stored in a file.

The most common terrain data available comes in form of a grid: this means that the points {x, y} are sampled uniformly with a grid. Grid elevation data is refered to as a digital elevation model or DEM, in short. In a grid terrain the individual coordinates {x, y} of the points need not be specified. It is sufficient if we know the coordinates of the upper-left corner of the grid and the spacing between sample points. If these are specified, the geographical coordinates of the grid points are implicit.

A grid terrain consists of a header, and a set of points {z} that represent the elevation of the terrain sampled with a uniform grid. Here is a grid file, set1.asc, that represents a basin in the Sierra Nevada range in the south-west of the US. Take a look at it.

The format of a terrain grid file is standard. The first 6 lines in the file represent the header of the terrain; it stores the number of rows in the grid, the number of columns, the geographical coordinates of the lower left corner, the spacing in the grid (assumed to be the same both horizontally and vertically), and a nodata value (this is the value that is used for points where the elevation could not be measured).

ncols         391
nrows         472
xllcorner     271845
yllcorner     3875415
cellsize      30
NODATA_value  -9999
Following the header there are nrows lines, each line containing ncols values, for a total of nrows * ncols values. These values represent the elevations sampled from the terrain. Thus, a grid terrain is basically a 2D-array of elevation values. You can assume that the elevations are integers (though you can also use floats or doubles if you want). The dimensions of the grid can be found out from the header. The other information stored in the header (xllcorner, yllcorner, cellsize, NODATA_value) you will happily ignore for this lab since we won't care about geographical coordinates---- but you have to read them anyways in order to get to the point in the file where you can start reading the elevations.

The suggested structure of your program

You will write (at least) two classes: Grid and GridGIS. Class Grid will know to handle grid terrains, that is, read a grid from file and visualize it. Class GridGIS will be the controller that will create and render various grids. Ideally it will have a GUI.

Class Grid

Class Grid will store a grid terrain and will provide useful methods to manipulate grids. Naturally, you will store the grid as a 2D-array.

You will write two main methods: a constructor to load a grid from a file and a method to render a grid on a frame.

The constructor: Reading a grid from a file

Grid (String fileName)
The constructor should take a file name as parameter. It will open the file, read the header, create the array to hold the elevation values, and load it with the values from the file. One thing you will need to figure out is how to work with files in Java. Below is a piece of code that opens a file, creates a scanner on it, and reads integers (while any).

Files in Java

File f = new File(fileName);
Scanner sc = null;
try {
     sc = new Scanner(f);
    } catch (FileNotFoundException e) {
     System.err.println("File not found!");
     System.exit(0);
    }
   
//here is how you would read a file of integers
while (sc.hasNextInt()) {
    int nextOne = sc.nextInt();
  
}
In a nutshell, the code creates a new File, then opens the file with the Scanner class. Opening a scanner from a file must be checked for an exception to make sure everything went ok. Once it has been opened, then the Scanner can basically act as an interator.

If you need to revisit scanners, we talked about scanners in one of the first lectures. The doc pages will be a lot of help, as usual.

Testing your grid constructor

To test that you read the right values, write a printGrid method that prints the grid (header and values); it should print the same information as when you open the grid.asc file. You may want to break this function into two: printHeader and printValues.

In addition to the printGrid method, write a printInfo method that prints the important information about the grid: rows, cols, h_min and h_max.

See some test grids at the top of this page. Some of them (europe.asc) are large, so better not open them in a browser. Your Grid class will not be able to handle large files (try it), but it should handle the rest of them just fine. While testing, you should try getting the small ones to work before you try any of the larger ones.

If you want to play with more data, you can download your own DEMs at USGS geodata site. One of the test DEMs I included, europe.asc, comes from this site under GTOPO30 project, which contains Global 30 Arc-Second Elevation Data Set global DEM with a horizontal grid spacing of 30 arc seconds (approximately 1 kilometer). These files come in a different format, so if you want to use one of them, let me know and I will be very happy to help you convert it. Another site where you can download (more recent) data is: here. This is SRTM data collected by NASA SRTM mission in 2000 and released for free for research and educational purposes. It covers the entire world at 3 arc second resolution (approx. 90m). It comes in tiles. Download a tile in a location that you like.

Rendering a grid

Now that your Grid constructor is working, you are ready to render a grid. Class Grid should provide a function to render the grid in a window. To render a grid you will loop through all the points (i,j) in the grid and draw each point (i,j) at an appropriate position and with an appropriate color.

Position: You will use the position (i,j) of the point in the grid to determine the position (x, y) of the point in the window. Note that we ignore true geographical coordinates, and we use (row, column) of the point--- we don't care about the true geographical location anyways, at least, not for this lab.

One thing to be careful with is the transformation from (row,column) to window (x,y) coordinates. You want the upper-left corner of the grid (row=0,col=0) to map to the upper-left corner of the window (x=0,y=0), i to y and j to x.

One problem you have to solve when rendering a grid is to work out the mapping of grid points to pixels on the screen. When you open a window in Java of size 400 by 300, well, that is 400 pixels wide by 300 pixels long. Now you have a grid that has, say, 1000 columns by 500 rows. How do you map one pixel of the map to a pixel on the screen? If the grid is larger than the window, as in this case, then several points from the grid will map to the same pixel on the window. Otherwise, if the grid is smaller than the window, a point in the grid will map to more than one point on the window.

You will need to figure out is how to convert from grid coordinates to window coordinates. First of all, the window should keep the same height-width ratio as the grid, so that the grid does not look distorted. Your program could handle small grids (therefore stretching the grid to fit on the window); and large grids, therefore finding a way to compress the grid to the window size.

In other words: when you render the grid open a window of default width (say DEFAULT_WIDTH=300) and that has the same aspect ratio as the image. If you render an image that's smaller than the window, you have to stretch it out to the window size. If you render an image that's larger than the window, you have to compress the image to fit in the window. Think of finding a transformation (function) that maps every point (i,j) in the grid to a point (x,y) in the window.

Colors: You will use the elevation of point (i,j) to determine the color with which you draw the point.

First, some background on colors. In Java the color is specified as a triplet (r, g, b), where r, g, b are values in the range [0.0, 1.0] or {0, 1, 2, ..255} specifying the red, green and blue intensity in the color, respectively. 0.0 represents the darkest, and 1.0 the lightest. WHITE corresponds to (1.0, 1.0, 1.0) and BLACK corresponds to (0.0, 0.0, 0.0). Remember that on a graphics object you can change the color using

Color c = new Color (r, g, b);
g.setColor(c);

In order to make the color reflect the elevation, you will have to come up with a formula to set r,g,b as function of the elevation; to draw a point, use drawLine with the same start and end point---this will, in effect, draw a point. Every time you call drawLine, you need to set the color as a function of the elevation of that line. In your case, the line will be a point (and you know the height of the point).

To come up with a logic for setting r, g, b, remember that r, g, b must be values in [0.0, 1.0] (or 0..255), with 0.0 corresponding to highest elevation, and 1.0 to lowest elevation. You'll implement at least 3 different color modes: gray, bucket, and gradient.

Gray scale: The easiest way is to set all color components to the same value. Equal values for r,g and b give a gray color, so it will result in a gray-scale image. Use BLACK for nodata points.

r = g = b = ...
c = new Color(r, g, b);
g.setColor(c)
The value of r,g,b has to be a function of the height of the point, so that:
  1. when the point is NODATA, the color is set to BLACK (r=g=b=0)
  2. when the point has height h_min, r=g=b=1 (white)
  3. when the point has height h_max, r=g=b=0 (gray)
  4. when the point has height h between h_min and h_max, the color should be set as a linear function f(h), such that f(h_min) =1 and f(h_max)=0

Buckets:Another easy solution is to split the height range between h_min and h_max into a number of color buckets, and assign to all heights falling in a bucket the same color. If you chose the number of buckets large (5 or more), this will look pretty good. When you compute the height range, be careful not to include NODATA points. In other words, h_min should not be -9999.

Color gradient: The best color map is to get a smooth coloring, which, instead of assigning to all heights in an interval the same height, it interpolates smoothly between the color at the two ends.

The core of the problem is the following: you have a height range [h1, h2] and you want to map it continuously to a color range [c1=(r1,g1,b1, c2=(r2,g2,b2)]. Let h be the height of a point in an interval [h1,h2]. Let's first consider the red component of the color. You want to find out a linear function r(h) such that:

Similarly, you find g(h) and b(h) for the other components blue and green.

Then, given a height h, you create a color:

  c = new Color(r(h), g(h), b(h))
Note that the three component functions r(h), g(h), b(h) depend on h1, h2, r1,r2, g1,g2, b1,b2, and thus there will be different functions for each height/color interval.
//the height ranges are stored in h[0..4], such that h[0]=h_min and h[4]=h_max
Color getColor_Gradient(height h)
     if h== NODATA ...
     if h in [h[0], h[1]]
             //compute the functions r(h), g(h), b(h)  for this interval
             c = new Color(r(h), g(h), b(h))
     if h in [h[1], h[2]]...
             //compute the functions r(h), g(h), b(h)  for this interval
             c = new Color(r(h), g(h), b(h))
     if h in [h[2], h[3]]...
              //compute the functions r(h), g(h), b(h)  for this interval
             c = new Color(r(h), g(h), b(h))
     if h in [h[3], h[4]]...
              //compute the functions r(h), g(h), b(h)  for this interval
             c = new Color(r(h), g(h), b(h))
     return c; 

For the rendering part, I am most interested in the looks of your rendering. Play with colors, and make it look like a real map. On a real map elevations are coded using colors: high elevations are brown, low elevations are green. Here are some good ones (bucket coloring).
[Brunswick DEM] [set1 DEM]

Extensions

So far you render a terrain by rendering each point separately. In general rendering of a surface is of better quality if the surface is represented as a triangular mesh rather than a discrete set of points (the graphics engine can compute normals to the surface to determine lighting, reflection and all sorts of other stuff). An improved renderer would view the grid as a set of triangles (see below). Instead of drawing points, draw lines between neighboring points so that overall the lines form a triangulation as below. To decide the color of a line, consider its height to be the average height of its endpoints.

Class GridGIS

The interface of GridGIS is open-ended. I imagine that class Grid will use a JFrame to render the window, in a similar way to the image lab from last week. Your GridGIS class should be able to load and display DEMs. If you don't include a GUI, the simplest way is to create a grid from a given file in the main method of GridGIS:
Grid elev = new Grid("set1.asc");
elev.printInfo();
elev.render();


Overview

Read the lab enough times before you start coding, so that what you have to do is clear and you can see how everything fits together. I would even go as far as saying that "the sooner you start coding...the longer it takes". Well, you get the point. Create your code one step at a time, compile often, test and debug before you move to the next step. Encapsulate, use small functions, write nice code. Don't focus on writing code that works, write nice code. In case your are interested, here is a tentative grading key.
Have fun!