Back from a journey or a holiday, you can use the crumpled tourist maps and street atlas pages at the bottom of your handbag or in the back seat of your car to ell the story of your travels. Digital holiday snapshots add accents and highlights to your stories; GPS traces provide a narrative structure for them. To get started, though, you need a base map which is georeferencedone for which you know the geospatial position of the top left corner, the space represented by each pixel in the image, and the rotation or skew of the image away from True North.
How can we georeference a random tourist map? If you know the precise location of several points on the ground and can identify them on the map, you can align them to fit on the map. It's possible to do this by hand, adjusting the projection of the points until they fit, but after we'd done this manually a few times, we were yearning for an automated solutionso we wrote one.
3.13.1. The Web Interface
Before starting, you should have at least three waypoints that you've collected "on the ground," with an east and a north coordinate for each point. The waypoints should be well separated from each other on the map you want to georeference. This will provide a more accurate approximation of the map's extents, rotation, and scale.
Next, visit http://mappinghacks.com/georeference/ to load up the interactive georeferencer. Hit the "Browse" button to view an image from your local hard drive. Now start adding ground control points. The program doesn't care about the kind of values you use for the spatial positions, so they can be in WGS84 latitude and longitude, the default for GPS units, or in UTM, which can be a lot easier for doing spatial calculations and plotting points.
Figure 3-39 shows the work in progress. For each ground control point, click the place on the map where your waypoint was taken. The x- and y-coordinates, in pixels relative to the image, will appear on the lefthand side. Type in your longitude or UTM easting, and latitude or UTM northing. Now click "Add GCP" to add this point as a ground control point. You should see a number appear on the map, indicating your new point. After three or more points have been added, the map is ready to "Rectify."
If successful, you should see a series of six long numbers in the box on the lefthand side. What are they? They represent the coefficients of a matrix transformation that can be applied to the digital image, to reproject it into cartographic space. Now what?
3.13.2. The World File
The simplest way to make use of this transformation is to put it in a separate file that accompanies the image. The format is just the list of six numbers, separated with line breaks. This is called a world file, often seen in the company of TIFFs as a TIFF world file (.tfw). A world file can accompany any image format with GIS raster data in it. A common convention is to take the first and last letters of the file extension, followed by a "w," for the world file extensionso a JPEG would travel with a .jgw, or an ESRI BIL raster with a .blw.
Once you have a world file describing an image, you can load it into your favorite GIS program. QGIS is a good choice for quick-and-easy viewing of GIS data. Also, the GDAL toolkit understands world files and can use them to convert raster data between different formats and different coordinate systems, as described in [Hack #68] .
3.13.3. The Matrix Transformation
The six numbers (coefficients) in the world file are values from a 3x3 affine transformation matrix. The matrix describes a linear transformation between one 2-D coordinate space and anothersay, from pixels to latitude and longitude. For each dimension, x and y, there's a scale between one coordinate system and the other, the skew between one and the other, and the offset of the imagee.g., where it begins, at its northwestern corner. Figure 3-40 shows the linear algebra used in the affine transform.
Figure 3-40. Transforming coordinates via an affine matrix
The last three values are constants there to make the matrix algebra work out, and in an affine transformation, they are always [0 0 1]. So the matrix can be represented as a vector of numbers [a b c d e f]. This is the order that the world file lists the coefficients in. If a map is oriented to True North, then the rotational coefficients b and c from its world file will be zero.
This transformation matrix might be used in the following fashion. Suppose that we have a map in a rectangular projection. If a user clicks on the map at (x, y), what is the geographic location corresponding to that point? If we have a world file expressed, for example, in UTM coordinates, then we can calculate the easting and northing as (x', y') with the following formulas:
x' = ax + cy + e y' = bx + dy + f
The actual target coordinate system doesn't matter; one can use the same technique to turn (x, y) into latitude and longitude, if the world file uses the same coordinates.
3.13.4. Projecting Arbitrary Maps in SVG
The affine transformation matrix doesn't have to be stored in a world file. SVG, the Scalable Vector Graphics format, has native support for matrix transformation. To take a raster map and use it directly as a base map in an SVG document, we can apply a transformation to the image like this:
The content of the transform parameter, you'll notice, is a matrix transform in exactly the form as the one in the world file! Now you can plot points in their original geographic coordinates directly onto this map. Unfortunately, the transformation process can produce a lossy effect in the SVG viewer, and the map can look jagged, blurry, or warpedespecially at extreme latitudes, if you're using unprojected coordinates. The other downside is that SVG viewers typically only display points with integer accuracy, so you have to multiply latitude and longitude by 100,000 to scale them up to something viewable.
Another approach to the same problem might be to transform the geometry of the points and leave the map alone. We can do this by taking the inverse of the original affine matrix, and applying it to the points, instead of the base map.