An addon to download the AHN data in GRASS GIS

dem
dtm
dsm
osgeo
open data
grass gis
data import
Introducing the new GRASS GIS addon to download AHN data. The AHN is the Dutch high-resolution DEM data. It uses the WCS protocol to download the data, ensuring that the extent and resolution of the imported data align with those of the original AHN data.
Author

Paulo van Breugel

Published

March 23, 2024

The r.in.ahn addon

In the previous post, I introduced the ‘Actueel Hoogtebestand Nederland’. To start with a small correction, the version 4 I mentioned isn’t the latest version. Version 5 is being rolled out and is already available for the northern part of the country. But for now, let’s focus on version 4.

In the previous post, I presented some steps to download the AHN for a specific area and import it in GRASS GIS. Downloading was easy using the r.in.wcs addon. However, a few extra steps were required to ensure the imported data would align with the extent and resolution of the original AHN data.

Easy, but why not make it easier yet? So, as it was a rainy day anyway, I used the code presented earlier and wrapped it up in the addon r.in.ahn. Let’s see how to download the DTM for the Land van Cuijk again, but this time using the new addon.


Download the DTM for a selected area

Note that this addon only works in locations with the coordinate reference system RD New (EPSG 28992). This is that the CRS of the original data, and the addon is meant to ensure you download the data as it is. This is akin to how the r.in.gdal import function works. If you want to import and reproject the data on the fly (similar to the r.import function), you can use the r.in.wcs addon.

This example shows the steps to download the 0.5-meter resolution DTM for the Land van Cuijk. You’ll need to install the r.in.wcs and r.in.ahn addons 1 2.

Install the required addons.

First step is to import the Python libraries. Note, this will not be repeated with the next scripts.

import grass.script as gs

Now, you can install the two addons using the g.extension function.

gs.run_command("g.extension", extension="r.in.wcs")
gs.run_command("g.extension", extension="r.in.ahn")

Install the two addons using the g.extension function.

g.extension extension=r.in.wcs
g.extension extension=r.in.ahn

You need the g.extension function to install addons. In the main menu, go to Settings > Addons extension > Install extension from addon. Alternatively, type in g.extension on the command line. The example below is for r.in.wcs. Repeat this step for r.in.ahn.

Figure 1: The g.extension function, started from the command line.
Figure 1: The g.extension function, started from the command line.

Now, download the administrative boundaries of the Dutch municipalities, and extract the boundaries of the “Land van Cuijk”.

Download layer with administrative boundaries of the neighborhood.

gs.run_command(
    "v.in.wfs",
    url="https://service.pdok.nl/cbs/wijkenbuurten/2022/wfs/v1_0?",
    output="municipalities",
    name="gemeenten",
)

Next, extract the boundaries of the municipality of “Land van Cuijk”

gs.run_command(
    "v.extract",
    input="municipalities",
    where="naam = 'Land van Cuijk'",
    output="LandvanCuijk",
)

Get the vector layer with the boundaries of the municipality.

v.in.wfs url="https://service.pdok.nl/cbs/wijkenbuurten/2022/wfs/v1_0?" \
    output=municipalities name=gemeenten

Next, extract the boundaries of the municipality of “Land van Cuijk”

v.extract input=municipalities where="naam = 'Land van Cuijk'" \
    output="LandvanCuijk"

Type in v.in.wfs on the command line or in the console. You can also find the function under main menu > File > Import vector data. This opens the following screen (you need to fill in parameters in two tabs):

Figure 2: Download the vector layer with the municipality boundaries. Define the base URL and the name of the output layer.
Figure 2: Download the vector layer with the municipality boundaries. Define the base URL and the name of the output layer.
Figure 3: Download the vector layer with the municipality boundaries. Fill in the name of the WFS layer to download.
Figure 3: Download the vector layer with the municipality boundaries. Fill in the name of the WFS layer to download.
Figure 4: Extract the boundaries of the municipality of Land van Cuijk. Select the name of the vector layer with municipalities and give the name of the output layer.
Figure 4: Extract the boundaries of the municipality of Land van Cuijk. Select the name of the vector layer with municipalities and give the name of the output layer.
Figure 5: Extract the boundaries of the municipality of Land van Cuijk. Fill in the query, which defines which features you want to select and save. Tip: Use the convenient query builder.
Figure 5: Extract the boundaries of the municipality of Land van Cuijk. Fill in the query, which defines which features you want to select and save. Tip: Use the convenient query builder.

Now, set the region to match the extent of the municipality. Note, you do not need to be concerned with the resolution; r.in.ahn will adjust the resolution and extent to make sure the imported data aligns perfectly with the original AHN data. This is done by setting the resolution to 0.5 meters and subsequently extending the region’s extent until it aligns with the original AHN data layer.

Set the region to match the vector layer LandvanCuijk.

gs.run_command("g.region", vector="LandvanCuijk")

Get the vector layer with the boundaries of the municipality.

g.region vector=LandvanCuijk

Type in g.region on the command line or console, or open the function’s window via menu | Settings | Computational region | Set region.

Figure 6: Set the region to match the extent of the vector layer LandvanCuijk.
Figure 6: Set the region to match the extent of the vector layer LandvanCuijk.

Now, you can run the r.in.ahn function to import the layer. Note that by default, the addon will change the region you just defined. It will set the resolution to 0.5 meters (this is the resolution of the AHN data you are about to download). It will furthermore expand the extent so that it aligns perfectly with the AHN data.

Import the DTM using the r.in.ahn function.

gs.run_command("r.in.ahn", product="dtm", output="dtm_05")

Get the vector layer with the boundaries of the municipality.

r.in.ahn product=dtm output=dtm_05

Type in r.in.ahn on the command line or console. A third way is shown in the image below. Open the Tools tab and go to Addons. The function should be available there.

Figure 7: Open the r.in.ahn plugin, select the product to download (dtm or dsm) and provide the output layer name.
Figure 7: Open the r.in.ahn plugin, select the product to download (dtm or dsm) and provide the output layer name.

Download whole tiles

If you set the -t flag, the function will import the DTM or DSM for all 6.5×5 km AHN tiles that overlap with the current region. In addition, a vector polygon layer will be created with the tile boundaries.

Define the region for which you want to download the data, and import the DSM using the r.in.ahn function. Set the -t flag to download the DTM for the area covered by the tiles that overlap with the region.

# Set the region
gs.run_command("g.region", n=412572, s=411280, w=188911, e=190085)

# Import the tile(s) that include the selected region 
gs.run_command("r.in.ahn", product="dtm", flags="t",
                output="dtm_05_subset")
# Set the region
g.region n=412572 s=411280 w=188911 e=190085)

# Import the tile(s) that include the selected region 
r.in.ahn -t product=dtm output=dtm_05_subset)

Open the g.region addon by typing g.region on the command line, in the console, or using the menu. Fill in the northern, southern, western, and eastern limits of the area for which you want to download the data.

Figure 8: Set the region bounds.
Figure 8: Set the region bounds.

Type in r.in.ahn on the command line or console. Or go to the Tools tab and go to Addons. The function should be available there. Fill in the required fields, like in the previous examples. Next, go to the Optional tab and select Download whole tiles.

Figure 9: Open the r.in.ahn plugin, select the product to download (DTM or DSM) and provide the output layer name.
Figure 9: Open the r.in.ahn plugin, select the product to download (DTM or DSM) and provide the output layer name.

The extent of the imported layer covers the 6.5×5 km AHN tile. In addition to the raster layer, you’ll have a vector layer with the boundary of the downloaded tile(s). This vector layer has the same name but with the suffix *_tiles*.

You should be aware that running the function will adjust the computational region so that it aligns with the imported data. You can avoid this by setting the -g flag, as illustrated in Figure 10.

Figure 10: A: The downloaded 6.5 x 5 km AHN tiles that intersect with the user-defined region, here indicated by the orange outline. In addition, a vector layer with the boundaries of the tiles is created. The red outline shows the adjusted region extent after running the function. The blue outlines show the boundaries of the downloaded tiles. B: The same, but with the -g flag set. With this flag set, the user-defined region (red outline) will not be altered .
Figure 10: A: The downloaded 6.5 x 5 km AHN tiles that intersect with the user-defined region, here indicated by the orange outline. In addition, a vector layer with the boundaries of the tiles is created. The red outline shows the adjusted region extent after running the function. The blue outlines show the boundaries of the downloaded tiles. B: The same, but with the -g flag set. With this flag set, the user-defined region (red outline) will not be altered .

Footnotes

  1. You are expected to be familiar with GRASS GIS and the concept of region used in GRASS GIS. If you are new to GRASS GIS, you are warmly recommended to first check out the GRASS GIS Quickstart and the explanation about the GRASS GIS database.↩︎

  2. Downloading the DTM for the whole municipality will take a while. If you want to speed things up, you can work with a smaller area by using your own vector data.↩︎