In the past couple days, I had to delineate a watershed using the PCRaster tool for some hydrological modeling work. Despite the excellent tutorial of https://github.com/jvdkwast/PCRasterTutorials/tree/main/PCRasterCatchmentDelineation, there were several pitfalls that took me a long while to figure out.
The first is, how do I know what x-y coordinates to put into "location.txt"? I only know the latitude and longitude of my stream gauge, but they did not work when I tried putting them in. It turns out the x-y coordinates are the coordinates of this stream gauge in the projected coordinate system of the flow direction (ldd) map. As such, the easiest way to approach the problem seems to be to use the QGIS PCRaster plugin. First, import the digital elevation map into QGIS, reproject it into a projected coordinate system, and set the layer's coordinate system to be the same as the digital elevation map. Second, calculate the ldd map using the PCRaster plugin's `lddcreate` command, so that the ldd map has the same projected coordinate system. Third, load the text file containing the lat-lon of the stream gauge as a map layer. QGIS automatically does the conversion, so that the gauge is shown in its correct position in the projected coordinate system. Fourth, calculate the stream order from the ldd map using the PCRaster plugin's `streamorder` command and check if the gauge is actually located in the correct stream. This is important because the lat-lon coordinates and the digital elevation maps are not entirely accurate, often causing the gauge to be outside the stream implied by the ldd map. If the gauge is outside, find the nearest in-stream pixel that seems to make sense. Placing the curser on top of the nearest in-stream pixel finally makes QGIS show the desired x-y coordinates in its status bar.
The second is, how to install the PCRaster plugin in QGIS? The plugin does not work on its own. Instead, one must install PCRaster in the same Python environment as QGIS accesses. In the end, I had to create a new conda environment and install PCRaster and QGIS both from the `conda-forge` channel. Moreover, the default `conda create -n ` command now pulls Python 3.11, but this Python version does not support compatible PCRaster and QGIS releases. I set the Python version to 3.9 to finally make things work.
The third is, why is my delineated catchment area incorrect? My stream gauge is on a major river so I know it must cover a large catchment area. Yet I got a zero-sized catchment at my first try. It turned out that in running the `lddcreate` command in PCRaster plugin, one must remove pits - pixels that draw flow from all directions and do not let water flow out - that do not make sense. One can set thresholds in the "Outflow depth", "Core area", "Core volume", and "Catchment precipitation" in the options of the `lddcreate` command to remove pits that are too shallow or too small. The default thresholds were too small for me, and as a result, the pits broke up my major river into many small ones in the stream order map created from the ldd map. After I set the Core area to be five times the default 9999999 (map units, for me it was meters), my delineated catchment finally made sense.
