Trace vector contours from a scanned map

Import

• Import a scanned map image. To save processing later you might want to mask out obvious text and other noise in a paint program first.

r.in.gdal

• Georectify map using GUI georectifier or i.points + i.rectify
• Switch to target location
• We now have a registered raster map containing a 300 dpi scan of our A0 sized paper map

Convert to vector map

• Check how many categories (map colors)

r.info chart

8 colors

• Zoom to the map's bounds

g.region rast=chart

• Display map (command line version)

d.mon x0
d.rast chart

• Query map values (command line version)

d.what.rast

black is category 1

• Create a reclass map, with just the black lines

r.reclass in=chart out=chart_black << EOF
1 = 1 black
EOF

• Thin the lines to be 1 raster cell wide

r.thin in=chart_black out=chart_black.thinned

• Sometimes it helps to "grow" the contours before thinning, especially when the scan is not of high quality and contours have "breaks" in them. Note that you may have to repeatedly try different radius values to get the desired effect.

r.grow in=chart out=chart_grown radius=n_cells

• Convert raster lines to vector lines

r.to.vect -s in=chart_black.thinned out=chart_lines

Clean vector map

• Remove any dangles smaller than 1km

v.clean in=chart_lines out=chart_lines_cleaned tool=rmdangle thresh=1000

• Remove known grid lines

v.mkgrid
v.buffer
v.overlay ain=scanned_lines bin=buffered_grid operator=not

• Clean by hand in with the digitizing tool

v.digit

or the QGIS GRASS vector editor

• Connect lines which were broken at dangles

v.build.polylines -q in=chart_lines_cleaned out=chart_polyline

Remove short line segments

• Add category numbers to lines

v.category in=chart_polyline out=chart_polyline_cat type=line

• Create a DB table to hold line length values

v.db.addtable chart_polyline_cat columns="length_km DOUBLE PRECISION"

• Upload line lengths to DB table

v.to.db chart_polyline_cat type=line option=length units=k column=length_km

• Extract line features longer than 5km (cleans out the noise)

v.extract in=chart_polyline_cat out=chart_lines_5km type=line where="length_km > 5"

View

• Display result over raster chart (command line version)

d.vect chart_lines_5km color=red width=2

• Export as a shapefile

v.out.ogr in=chart_lines_5km dsn=chart_lines_gt5km

• View shapefile in QuantumGIS

qgis chart_lines_gt5km/chart_lines_5km_crop.shp

Done!