Sunday, 17 June 2012

Digital Terrain Modelling and Mapping

This is the first post in a new blog series focusing on all the fun you can do with digital terrain data. It’s going to be more related to topographic than thematic mapping, although the two map types can be used in combination. I’ll also try to restrict myself to data and software in the public domain.

My focus area is going to be Jotunheimen in Norway, a mountainous area popular among hikers and climbers. The 29 highest mountains in Norway are all in Jotunheimen, including the very highest - Galdhøpiggen (2469 m). Here's a few pictures I took on a hike back in July 2006 which gives you a glimpse of the beautiful landscape:

Skogadalen valley and Hurrungane mountain range

Fannaråkhytta

Bøverbreen glacier

It’s many different data formats for storing and working with digital terrain data. I’ll start looking at Digital Elevation Models (DEM) represented as rasters. A raster DEM is a grid of squares or pixels where each square represents an elevation at a geographic location. The size of each square (in meters or arc-seconds) tells us about the level of detail or data resolution.

In Norway, it’s possible to get digital terrain data with 10 meters resolution, but unfortunately it’s not yet available to the public. Instead I found a suitable dataset on Viewfinder Panoramas, created from topographic maps. You’ll find  DEMs with global coverage at various resolutions here, here and here.

I downloaded two files from Viewfinder Panoramas (N61E007.hgt and N61E008.hgt), each covering an area of 1 degree latitude and longitude marked in yellow:


Bigger map | KML

You can use the great GDAL library to combine and clip the DEMs to the area of interest. I used gdalbuildvrt to create a combined virtual dataset of the DEM files:

gdalbuildvrt jotunheimen.vrt N61E007.hgt N61E008.hgt

Then you can use gdalwarp to project and clip the DEM to the area of interest:

gdalwarp -t_srs EPSG:32632 -te 432000 6790000 492000 6850000 -r bilinear jotunheimen.vrt jotunheimen.tif

The original DEM data is unprojected and I’m projecting the data using UTM zone 32N (EPSG:32632). UTM coordinates (meters) are used to clip the dataset to an area covering 3,600 km² (60 x 60 km). The area is shown with a red square on the map above. The -r bilinear option is important when projecting elevation data because other resampling methods tend to produce odd stripes in the image.

We finally have a DEM for Jotunheimen (download). You won’t see much if you try to open this image in Windows. The image is almost all black. Use gdalinfo to see the properties of this image:

gdalinfo -mm jotunheimen.tif

Driver: GTiff/GeoTIFF
Files: jotunheimen.tif
Size is 3134, 3134
Coordinate System is:
PROJCS["WGS 84 / UTM zone 32N",
    GEOGCS["WGS 84",
        DATUM["WGS_1984",
            SPHEROID["WGS 84",6378137,298.257223563,
                AUTHORITY["EPSG","7030"]],
            AUTHORITY["EPSG","6326"]],
        PRIMEM["Greenwich",0],
        UNIT["degree",0.0174532925199433],
        AUTHORITY["EPSG","4326"]],
    PROJECTION["Transverse_Mercator"],
    PARAMETER["latitude_of_origin",0],
    PARAMETER["central_meridian",9],
    PARAMETER["scale_factor",0.9996],
    PARAMETER["false_easting",500000],
    PARAMETER["false_northing",0],
    UNIT["metre",1,
        AUTHORITY["EPSG","9001"]],
    AUTHORITY["EPSG","32632"]]
Origin = (432000.000000000000000,6850000.000000000000000)
Pixel Size = (19.144862795149969,-19.144862795149969)
Metadata:
  AREA_OR_POINT=Area
Image Structure Metadata:
  INTERLEAVE=BAND
Corner Coordinates:
Upper Left  (432000.000, 6850000.000) (7d42'39.90"E, 61d46'36.81"N)
Lower Left  (432000.000, 6790000.000) (7d43'59.45"E, 61d14'18.21"N)
Upper Right (492000.000, 6850000.000) (8d50'54.02"E, 61d46'58.29"N)
Lower Right (492000.000, 6790000.000) (8d51'3.39"E, 61d14'39.21"N)
Center      (462000.000, 6820000.000) (8d17'9.19"E, 61d30'42.34"N)
Band 1 Block=3134x1 Type=Int16, ColorInterp=Gray
    Computed Min/Max=15.000,2466.000


What does it tell you?

  • The image format is GeoTIFF, which allows georeferencing information to be embedded within a TIFF file.
  • The image or grid size is 3134 x 3134 pixels or cells.
  • The coordinate system is UTM zone 32N using the WGS 84 spheroid.
  • The coordinate unit is meters.
  • Each pixel or grid cell covers an area of approximately 19 x 19 meters.
  • The extent of the area is defined by four coordinate pairs.
  • The image is single band and the elevation data is stored as 16-bit integers. The minimum elevation is 15 meters and the maximum is 2,466 meters.
Let’s use gdal_translate convert the GeoTiff image to PNG showing elevation with 256 shades of gray (-scale).

gdal_translate -of PNG -ot Byte -scale 15 2466 0 256 jotunheimen.tif jotunheimen.png



This image gives you a glimpse of the elevation in Jotunheimen. Dark areas are valleys with low altitude. Light areas are mountain ridges with high altitude. Areas with uniform shades are lakes or areas with less variation in elevation.

Hillshade is a more effective technique to visualise terrain data. That’s the topic for the next blog post

2 comments:

Vegard Paulsen said...

Veldig bra! Jeg prøver for tide å lære meg litt mer om dette selv, for å frese ut height maps fra Norge i MDF på cnc maskina mi. Har du laget et gråskala heightmap fra Oslofjorden fra oslo ned til Fredrikstad? Det hadde vært kjempetøft å prøve å frese ut i mdf! Jeg sliter med å finne noen kart som både viser dybde under vannlinja og høyde over vannlinja i samme gråskala. Vet du hvordan jeg kan få til det?

Bjørn Sandvik said...

Hei Vegard,

Du trenger en detaljert høydemodell for å få litt futt ut av "flate" oslofjorden. 10 x 10 meter høydedata blir frigitt fra kartverket 1. oktober i år. Tilsvarende data for havdybde er det vanskelig å få tak i.