Project 1: First Steps with VTK

 

Handed out: Wednesday, January 16, 2013

Due: Wednesday, January 30, 2013

Preamble

This first project will walk you through some basic aspects of VTK while allowing you to create an interesting visualization along the way. The main objective here is not to create impressive images but to make sure that you have a working installation of VTK and that you have a basic understanding of the general structure of the library. You have two weeks to complete this project (see submission instructions below).

Overview

The topic of this project is geovisualization, in other words the visualization of geographic data. Specifically, you will be visualizing a dataset that combines bathymetry/topography data with satellite imagery of the earth surface. Your visualization will include a graphical user interface that interactively controls the data being displayed.

Description

Bathymetry concerns itself with the measure of the underwater depth of the ocean floor while topography measures land elevation with respect to sea level. Combined, these two pieces of information provide a complete description of the earth surface’s geometry (I use the term bathymetry hereafter to refer to this combined information). NASA’s Blue Marble project provides access to high resolution bathymetry datasets along with satellite imagery of the earth surface acquired at different times of the year. In this project, you will learn how to use VTK to create interactive visualizations of this amazing data.

Practically, the data consists of grayscale images (bathymetry / topography) and color images (satellite views) with resolutions up to 86400 x 43200 (3.7 gigapixels). In each image, the earth surface is sampled at regular angular intervals spanning the domain [-, ] x [-∏ / 2, ∏ / 2] in spherical coordinates (see illustration below).


Specifically, each row of the image corresponds to a parallel (line of constant latitude), while each column corresponds to a meridian (line of constant longitude). In particular, the North and South poles (where the parallels degenerate to points) are each represented by a an entire row of uniform value.

Task 1. Height Map Visualization

A natural way to represent this type of data is through a height map, whereby each (x,y) two-dimensional geolocation is associated with a third z coordinate corresponding to the elevation (or sea depth) at that location, thus creating a surface in 3D that conforms with the elevation profile. VTK provides a convenient algorithm to create such maps called vtkWarpScalar as demonstrated by the program ImageWarp.py (available under Examples/VisualizationAlgorithms/Python in your source directory). Your first task is to write a short program in the language of your choice that takes in input the name of a vtkImageData file containing bathymetry information (in meters) and visualize this data as a height field using vtkWarpScalar. You shall assume that the earth radius is 6,371,000 meters and therefore scale the width of the image (which physically corresponds to earth circumference) to about 40,030,000 meters (the height of the image is half that length).

Note that vtkWarpScalar offers a scale factor (selectable through SetScaleFactor) that controls the amount by which the image is moved up or down for a given scalar value at each point. Include in your program a slider user interface (UI) to control the value of this scaling factor and interactively update the visualization. Graphic UI can be implemented in two main ways in VTK:

  1. The first solution consists in using a general purpose GUI library for which a VTK binding exists (the binding is necessary to allow the GUI to communicate properly with VTK). Tk (available for Python, Java, and C++) and Qt (C++) are both supported by VTK. Corresponding examples are available in the source distribution of VTK under Examples/GUI/Python and Examples/GUI/Qt.

  2. The second solution is to use the vtkWidget mechanism, which offers a limited set of GUI functionalities within VTK. Corresponding examples (in C++) are available under Examples/Widgets/Cxx. For the task at hand, you will find Slider.cxx and Slider2D.cxx most directly relevant.

To improve the appearance of the resulting geometry, you will make use of color coding, in other words, the representation of scalar values through colors. VTK enables this mapping through vtkColorTransferFunction, whereby the user assigns specific colors to select values. A (loosely realistic) color scale for this dataset looks something like the following (it does a better job under water than above sea level).


Refer to the various examples provided in the documentation about the use of vtkColorTransferFunction.

To summarize, your program should be named HeightMap and take the name of the input file from the command line  as follows.

HeightMap[.py] bathymetry.vtk

Furthermore the slider included in your program should support scaling values in the range [0, 50].

Sample Results (the images below link to higher resolution versions)


Task 2. The World is *NOT* Flat

The results achieved in the first step allowed you a first look at the data but they do not convey the actual shape of the various continents due to the significant distortion that is caused by the latitude / longitude parameterization of the image. To remedy this problem we need to visualize the data directly on a sphere.

One way to turn the bathymetry dataset (which in an image) into a scalar field defined on a sphere would be to transform the coordinates of its positions from (longitude / latitude) to (x, y, z). If we denote the longitude 𝜃 and the latitude 𝜑, the formula for this transformation is simply:


Doing this transformation in VTK however is made nontrivial by the fact that the input data set (an image) and the resulting sphere correspond to different dataset types: while an image has a regular topology (the connectivity between its pixels) and a regular geometry (the spatial distribution of its pixels), the sphere would have a regular topology but an irregular geometry. In VTK, one dataset is a vtkStructuredPoints while the second is a vtkStructuredGrid. In addition, we need the sphere to be closed, which could not be achieved without changing the topology. Hence, simply modifying the positions’ coordinates would result in an error.

To make your task easier, I am providing you with a sphere dataset that has the correct geometry and topology and supports a one-to-one match with the pixels of the bathymetry dataset. Combining the geometry of the sphere dataset with the scalar (bathymetry) values of the image data will create the desired effect. This combination is made possible in VTK by a filter called vtkMergeFilter.  The usage of this filter is illustrated in the same ImageWarp example discussed above. Once you have the bathymetry mapped onto a sphere, you should apply the same vtkWarpScalar filter used in Task 1 to visualize the relief of the earth. Again, incorporate a slider UI in your program to control the amount of relief emphasis. Note that the sphere that I am giving you is already scaled to the radius of the earth.

Finally, you should incorporate in your visualization the satellite imagery that is available. The solution for that consists in using texture mapping, a mechanism through which an image is “glued” onto a polygonal model. Here the polygonal model is the elevation of the earth surface produced by vtkWarpScalar. The important point here is that texture mapping allows you to use high resolution images at relatively low cost for the rendering pipeline. In particular, the resulting of the satellite images that you will use in this project have much higher resolution than the bathymetry datasets that we use to construct geometry. This is key to create realistic results without overwhelming the rendering engine. Practically, the sphere dataset comes pre-equipped with texture coordinates that tell VTK how to match vertices to locations in the image. To use them however you will need to pass on that information through vtkMergeFilter. You will use vtkTexture to store the satellite image and pass it on to the renderer through SetTexture as demonstrated in examples such as Examples/Rendering/Python/TPlane.py. Note finally that since the images are available in jpeg format, you will need to use a vtkJPEGReader to import them.

Your executable EarthMap for this task should get its input from the command line as follows.

EarthMap[.py] bathymetry.vtk world.jpg sphere.vtk

Sample Results


Task 3. The Inexorable Rise of the Oceans

The steady rise of the sea level is a ongoing trend that represents a major challenge for the world population. The data available in this project gives us all the information necessary to interactively visualize this phenomenon. For a given value of the sea level increase, the goal is to highlight those sections of the continents whose altitude fall below the sea level and are now flooded. Note that this approach is a (coarse) simplification of the actual physical process (in particular, this criterion considers submerged a region of low altitude, even if it is not contiguous to the sea) but we shall accept this limitation for the sake of simplicity.

The idea of selecting regions based on value is closely related to the concept of level set and its topographic equivalent, the contour line. Here we are not interested in selecting a particular elevation value (say the value corresponding to an elevated sea level h) but instead we wish to identify the entire region whose altitude is comprised in the interval [0, h], since those are precisely the regions that are newly submerged (*). VTK provides a filter named vtkClipPolyData that we will use to determine and visualize those regions. A limitation of this filter is that it is only designed to determine which half of a dataset is above or below a certain value. To achieve the intended effect here we need to apply it twice: once to determine the region that is lower than the raised sea level ([-∞, h]) and a second time to remove from this selection the regions that were already below the base sea level ([0, h]) since those regions are already properly visualized as submerged in the satellite image.

The application of vtkClipPolyData will produce a vtkPolyData object that describes the selected region. Visualizing this output directly would lead to visual artifacts since the selected polygons would overlap with the geometry of the textured earth surface (from Task 2). To prevent this issue, we must make sure that the selected region is slightly “lifted” above the earth surface. For that, you will use vtkWarpVector. This filter is similar to vtkWarpScalar in that it deforms the geometry of its input based on its values. However it requires a vector value associated with each vertex and uses the vector to translate the location of the vertex. This translation can be scaled by a user defined factor. Practically you will apply this filter to the radial vector information included in the sphere dataset. Make sure to pass on that information to the vtkMergeFilter used in Task 2 (through SetVectorsConnection()) to make it available in subsequent steps.

Once lifted, the selected region can be color mapped using a color scale that can be seamlessly merged with the underlying texture. The scale that you should use is the following:

RGB(-30) = (0,     0.22,  0.43)

RGB(0)   = (0.125, 0.416, 1.000)

RGB(500) = (0.353, 0.553, 1.000)

RGB(2000)= (0.635, 0.753, 1.000)

RGB(8500)= (0.761, 0.847, 1.000)

Note that the scale starts at altitude -30 meters: this is needed to prevent visible artifacts close to the shore line in the satellite image (where the altitude is zero). Specifically you should first clip below h and then clip above -30 before applying this scale.

For this last part of the project, create a program SeaRise that applies those ideas to create a visualization in which the sea level can be selected interactively (from 0 meter to 9000 meters increase) and the visualization updated accordingly. Use a slider UI to that end. In this case, keep the relief at scale (scale factor = 1) to ensure the consistency of your results. The function call takes the same parameters as in the previous task.

SeaRise[.py] bathymetry.vtk world.jpg sphere.vtk

(*) Practically we need to visualize the region corresponding to altitude [-𝜺, h] (where 𝜺 is a small positive number ) instead of [0, h] to avoid visualization artifacts.

Sample Results


Datasets

As indicated above, the data used in this project is courtesy of NASA Earth Observatory and available on the web site of the Blue Marble project. Since the server can be slow at times, I am providing you with local copies of the relevant files. In addition, I converted the bathymetry files to VTK format to simplify their manipulation in your programs. Finally, I created different versions of the data, corresponding to different resolutions. Make sure to start working with the low resolution datasets as you experiment with your program and move to higher resolution images once your code is working. Also, determine the size of the data that is compatible with your RAM and graphics hardware limitations.

Bathymetry (original URL)


  1. small (1081 x 540)

  2. medium (2161 x 1080)

  3. large (4321 x 2160)

Blue Marble Next Generation (Satellite imaging) (original URL)


  1. small (2160 x 1080)

  2. medium (5400 x 2700)

  3. large (21600 x 10800)

Globe geometry

Note: you must select a geometry file corresponding to the resolution of your bathymetry file.

  1. small (1081 x 540)

  2. medium (2161 x 1080)

  3. large (4321 x 2160)

Submission

Submit your solution before January 30, 2013 11:59pm. turnin will be used to submit all projects this semester. Please follow the procedure described hereafter.

  1. Include the 3 program files along with any other source code you may have

  2. Include Makefile (if applicable). If you are using Visual Studio include the solution file (.sln) and the project files (.vcproj) that are needed to build your project.

  3. Include high resolution sample images showing results for each task

  4. Include a html report summarizing what you have done. The report should include mid-res images with links to high-res ones.

  5. Include README.txt file with compilation / execution instructions (optional).

  6. Include all submitted files in a single directory

  7. Do not include binary files

  8. Do not include data files

  9. Do not use absolute paths in your code


After logging into data.cs.purdue.edu, submit your assignment as follows:

% turnin -c cs530 -p project1 <dir_name>

where <dir_name> is the name of the directory containing all your submission material.

Keep in mind that old submissions are overwritten by new ones whenever you execute this command. You can verify the contents of your submission by executing the following command:

% turnin -v -c cs530 -p project1

Do not forget the -v flag here, as otherwise your submission would be replaced by an empty one.