I'm looking for a set of outlines of the river basins in the UK – I don't need a particularly detailed or highly accurate map – just a map saying point X will drain into the Thames or point Y will just go into the sea! Just the main river names would be fine. I have the full Environment Agency rivers set and the OS DTMs if I want… but how do I convert this into a catchment boundary?
[GIS] River drainage basin boundary
catchment-areadatademhydrographyunited-kingdom
Related Solutions
The analysis has already been done in a reply to the antecedent question, but perhaps an illustration will help.
There are two major components of error: the "d8" algorithm, which represents flows in only eight cardinal directions, and the effect of projection (or lack of it). Let's focus on the latter, because this seems to be the principal concern.
The error depends on the distortions in the projection and on the terrain itself. Locally, over a small region, all projection distortions on the earth's surface amount to a stretch in one direction compared a perpendicular direction: this is why a (properly computed) Tissot Indicatrix is a perfect ellipse, because an ellipse is just a stretched circle. The terrain can have any aspect (flow direction). To handle this, let's look at a terrain that indeed has points in all possible directions with simple flow lines: a cone.
Overlaid on this color-shaded contour map of cone elevation is a collection of streamlines showing the directions where water would flow. You can confirm these streamlines are correct by checking that they cross the contours at right angles.
By choosing appropriate units of measurement and an appropriate origin for the coordinate system (at the cone's apex), the equation for the elevation in terms of coordinates (x,y) is simply
z = -Sqrt(x^2 + y^2).
The streamlines are always parallel to the gradient of z (in the reverse direction), computed by differentiating this formula with respect to x and y:
-Grad(z) = (x, y) / Sqrt(x^2 + y^2).
The coefficient 1/Sqrt(x^2 + y^2) doesn't change the direction, so we can ignore it for the purposes of understanding streamlines. Thus, at any location (x,y), the streamline points in the direction (x,y).
The effect of a horizontal stretch in the coordinates (by a factor of 2 in this image) is to stretch all contours (without changing the contour levels: heights are not affected by projections). Although (of course) the contours represent true circles, they no longer look like true circles on the map. Nevertheless, when the streamlines are computed in these coordinates, they must cross the contours at right angles just as before.
The effect of the stretch is to put the elevation at any point of coordinates (x,y) at new coordinates (stretchx, y). Consider this in reverse: the elevation at coordinates (X,Y) = (stretchx,y) must be the value of z computed at (x,y) = (X/stretch,Y). Therefore the equation of the apparent surface in this projection is
z = -Sqrt((x/stretch)^2 + y^2).
Differentiating, we compute
-Grad(z) = (x/stretch^2, y) / Sqrt((x/stretch)^2 + y^2).
Again the common factor matters not; thus, at any location (x,y), the computed streamline points in the direction (x/stretch^2,y). This was the formula used to draw the streamlines in the preceding picture. You can see they correctly cross the contours at right angles.
This third image reprojects the previous picture. The surface is shown once more without distortion. However, the streamlines no longer appear to cross the contours at right angles. This was the case even in the previous picture: due to the distortion therein, the angles only appeared to be right angles. The crossings were incorrect all along. That's why not projecting (or using a nonconformal projection) is a mistake. The question is how big a mistake it might be. Some have claimed it's of little consequence (at least at low to moderate latitudes).
This reprojection (to remove distortion in the map) moves the point at (x*stretch, y) back to (x,y). The stream direction previously computed at this point was stored in a grid (as an angle or a direction code): it does not change. Therefore the computed stream direction at (x,y) is (x/stretch^2, y).
This quantifies the effect of a reprojection on all possible flow directions, as shown by the difference between the first and last graphic. Here's their overlay, without the contour plot for distraction:
Reprojection affects directions differently depending on how the flow is oriented with respect to the major axis of the Tissot Indicatrix. It is a quadratic function of the relative linear distortion in the projection. As such, it exaggerates even slight amounts of distortion. (The factor of two illustrated here is somewhat extreme but realistic: it is the distortion introduced by failing to project--that is, using geographic coordinates as map coordinates--at latitudes of 60 degrees.)
With a little bit of trigonometry one can use these results to compute the angular error in flow direction as a function of the correct direction. Here is a graph of the errors associated with using a geographic (unprojected) coordinate system at latitudes 20, 30, 40, 50, and 60 degrees. (Of course the larger errors are associated with higher latitudes.)
The "true direction" is in degrees east of north. Positive angular differences occur when the apparent direction (computed without projecting lat, lon) is counterclockwise of the true direction.
Remember, you have to superimpose the d8 errors on top of these!
I've been in LiDAR processing for a couple of years now. The best approach we've found is to classify the suspect water points to something other than ground. Should be easy just classifying based on intensity (near nadir points will have high intensity, whereas turbid water will be close to 0) and laser shots are usually absorbed near shore anyway. However, there is still some interpretation and manual modification required
After water points are classified to non-ground, export a ground-only DEM representing intensity. The data voids in the DEM will have 0 values, representing probable water bodies. With some tinkering (i'm not giving away all of my secrets) in model builder, you can convert the raster to polygon, apply some smoothing, and voila...semi-automatic extracted breaklines.
Be warned that this still takes time and a careful eye to do it properly. Nothing is as accurate as the human eye. Good luck. Results may vary.
Best Answer
The trick is to Google "UK river catchment map". (Other combinations that include "watershed" don't work.) You will find
A search involving "watershed" did turn up a historical (1883) map in the Rumsey Map Collection.