PostGIS gebruiken voor BAG-data in RD New

PostGIS is een extensie voor de database PostgreSQL. Het project beschrijft zichzelf zo:

PostGIS extends the capabilities of the PostgreSQL relational database by adding support for storing, indexing, and querying geospatial data.

Bron: postgis.net, gelezen op 11 oktober 2026. Een extensie betekent dat u PostgreSQL blijft gebruiken voor beheer, back-ups en rechten, en dat PostGIS daar kolomtypen en functies aan toevoegt.

PostGIS is de aangewezen plek voor ruimtelijke data die meerdere gebruikers raken of die u aan andere systemen koppelt. Voor een eenmalige analyse is een bestand vaak genoeg, zie Python voor GIS.

U leest de BAG OGC API van PDOK met het OAPIF-stuurprogramma van ogr2ogr en schrijft de resultaten naar PostGIS. De BAG slaat zijn geometrie al op in RD New, dus u hoeft niet om te rekenen. Het gaat om een extract: de API is bedoeld voor een gebied, niet voor een heel land.

De volgende commando's laden een vierkante kilometer rond het centrum van Amersfoort. De openingsoptie PAGE_SIZE houdt elk verzoek klein; de standaard is volgens de GDAL-documentatie 1.000.

# Panden en verblijfsobjecten in een vierkante kilometer rond de Onze Lieve Vrouwetoren in Amersfoort
API="OAPIF:https://api.pdok.nl/kadaster/bag/ogc/v2"

ogr2ogr -f PostgreSQL "PG:dbname=bag" "$API" pand   -spat 154500 462700 155500 463700 -spat_srs EPSG:28992   -t_srs EPSG:28992 -nlt PROMOTE_TO_MULTI -lco GEOMETRY_NAME=geom   -overwrite -oo PAGE_SIZE=200

ogr2ogr -f PostgreSQL "PG:dbname=bag" "$API" verblijfsobject   -spat 154500 462700 155500 463700 -spat_srs EPSG:28992   -t_srs EPSG:28992 -lco GEOMETRY_NAME=geom   -overwrite -oo PAGE_SIZE=200

Resultaat op 11 oktober 2026, herhaald op 12 oktober met GDAL 3.13.3: 3.105 panden en 6.070 verblijfsobjecten, voor beide lagen samen tussen 13 en 34 seconden. De panden kwamen als geometry(MultiPolygon,28992), de verblijfsobjecten als geometry(Point,28992), elk met een GiST-index. Alle 3.105 geometrieën waren geldig volgens ST_IsValid. Zet de opties van -spat niet in een shellvariabele: zsh, de standaardshell op macOS, splitst zo'n variabele niet in losse argumenten, en ogr2ogr toonde in de test dan alleen zijn gebruiksinstructie.

Een grotere selectie ging mis: een gebied van 2 bij 1,5 kilometer eindigde met de standaardpaginagrootte in HTTP-status 500 (time-out), en met pagina's van 200 was het na acht minuten nog niet klaar. Voor een hele gemeente of heel Nederland gebruikt u de Atom-download. Meer over de dataset staat op BAG en over PDOK als aanbieder in het overzicht van open geodata. Hoe de API-standaard werkt, leest u bij WFS en OGC API Features.

U maakt een GiST-index op de geometriekolom en draait daarna ANALYZE. Ogr2ogr deed dat hier al: de tabel pand had na het laden de index pand_geom_geom_idx. Handmatig doet u hetzelfde met CREATE INDEX pand_geom_idx ON pand USING GIST (geom);. De PostGIS-handleiding legt het nut uit:

Without indexing, a search for features requires a sequential scan of every record in the database.

Bron: PostGIS-handleiding, hoofdstuk over databasebeheer, gelezen op 11 oktober 2026. Het verschil is op 3.105 rijen klein maar meetbaar. Dezelfde query, alle panden binnen 100 meter van het RD-punt (155000, 463000), gaf de volgende plannen.

-- zonder index (index tijdelijk verwijderd binnen een transactie)
 Aggregate (actual rows=1 loops=1)
   ->  Seq Scan on pand (actual rows=119 loops=1)
         Filter: st_dwithin(geom, '...'::geometry, '100'::double precision)
         Rows Removed by Filter: 2986
 Execution Time: 3.895 ms

-- met de GiST-index
 Aggregate (actual rows=1 loops=1)
   ->  Bitmap Heap Scan on pand (actual rows=119 loops=1)
         Filter: st_dwithin(geom, '...'::geometry, '100'::double precision)
         Rows Removed by Filter: 37
         ->  Bitmap Index Scan on pand_geom_geom_idx (actual rows=156 loops=1)
               Index Cond: (geom && st_expand('...'::geometry, '100'::double precision))
 Execution Time: 0.356 ms

Over vijf herhalingen op 12 oktober 2026 lag de uitvoeringstijd met index tussen 0,23 en 1,1 milliseconde en zonder index tussen 3,9 en 132 milliseconden, op een laptop die ook andere taken draaide. Dat zijn geen benchmarkcijfers; de orde van grootte is het punt. Het patroon telt: de index levert kandidaten op uit een omhullende rechthoek, en ST_DWithin controleert daarna de werkelijke afstand. De documentatie van ST_DWithin zegt het zo:

This function call includes a bounding box comparison that makes use of any indexes that are available on the geometries.

Drie patronen komen telkens terug: dichtstbijzijnde objecten zoeken, objecten aan elkaar koppelen en coördinaten omrekenen. Alle drie draaiden op de tabellen hierboven. Het eerste voorbeeld zoekt panden dichtbij een punt en sorteert op afstand met de operator <->.

-- panden binnen 50 meter van een punt, dichtstbijzijnde eerst
SELECT identificatie, bouwjaar, round(ST_Area(geom)::numeric, 1) AS m2
FROM pand
WHERE ST_DWithin(geom, ST_SetSRID(ST_MakePoint(155000, 463000), 28992), 50)
ORDER BY geom <-> ST_SetSRID(ST_MakePoint(155000, 463000), 28992)
LIMIT 5;

  identificatie   | bouwjaar |  m2
------------------+----------+-------
 0307100000333887 |     1445 | 211.0
 0307100000379020 |     1920 |  90.1
 0307100000371329 |     1800 |  29.1
 0307100000536538 |     1922 |   8.6
 0307100000372777 |     2004 |  54.4

Het tweede telt per pand de verblijfsobjecten in de extract, en vergelijkt dat met het aantal dat de BAG zelf bij het pand vermeldt.

-- verblijfsobjecten per pand, via een ruimtelijke join
SELECT p.identificatie, p.aantal_verblijfsobjecten AS bag_aantal, count(v.*) AS in_extract
FROM pand p
JOIN verblijfsobject v ON ST_Within(v.geom, p.geom)
WHERE p.status = 'Pand in gebruik' AND p.aantal_verblijfsobjecten > 20
GROUP BY 1, 2
ORDER BY 2 DESC
LIMIT 5;

  identificatie   | bag_aantal | in_extract
------------------+------------+------------
 0307100000314600 |        272 |        272
 0307100000351286 |        128 |        128
 0307100000374052 |         97 |         97
 0307100000343270 |         86 |         18
 0307100000376537 |         77 |         77

Het vierde pand heeft 86 verblijfsobjecten in de BAG maar 18 in de extract. De rest ligt buiten het gekozen vierkant. Een extract met een rechthoekig gebied knipt dus objecten door, en tellingen per pand zijn alleen geldig voor panden die er helemaal in liggen. Dat is geen fout in de query maar een eigenschap van de selectie.

Het derde voorbeeld rekent een WGS84-punt om naar RD. Het resultaat wijkt minder dan een millimeter af van (155000, 463000) en is gelijk aan wat pyproj en proj4js met de juiste definitie geven.

-- WGS84 naar RD, met dezelfde PROJ-bibliotheek als pyproj
SELECT ST_AsText(ST_Transform(ST_SetSRID(ST_MakePoint(5.3872035, 52.1551723), 4326), 28992));

POINT(154999.9995416017 463000.00031318853)

De adresgegevens zitten als tekstkolommen in verblijfsobject, zoals postcode, huisnummer en openbare_ruimte_naam. Een snelle check op een enkel adres kan ook zonder database met de BAG-adrescheck.

Vier valkuilen zijn gemeten, de rest volgt uit de documentatie.

• ST_Distance in de WHERE-clausule. ST_Distance(geom, punt) < 100 gebruikte de index niet en las de hele tabel. Gebruik ST_DWithin. De PostGIS-documentatie noemt ST_Distance ook zelf niet-indexeerbaar in haar voorbeeld.
• Een functie om de kolom. ST_DWithin(ST_Transform(geom, 4326), ...) leidde ook tot een volledige tabelscan. Reken liever het punt om naar het stelsel van de kolom dan de kolom naar het stelsel van het punt.
• Het verkeerde SRID. Een punt met SRID 4326, of een punt zonder SRID maar met graden als coördinaten, gaf tegen de RD-kolom in de test met de GiST-index stil nul rijen. Een RD-punt zonder SRID gaf wel meteen de fout Operation on mixed SRID geometries, en met de index uitgeschakeld gaven ook de eerste twee die fout (PostgreSQL 17, PostGIS 3.6.4, herhaald op 12 oktober 2026). Of u een fout ziet, hangt dus af van het queryplan. Zet altijd SRID 28992 op uw punten.
• Eenheden. Met geometry in RD zijn afstanden meters, in WGS84 graden. De handleiding noemt voor geography meters als eenheid.
• Kolomnamen. De BAG-API levert kolommen met een punt in de naam: pand.href in verblijfsobject en verblijfsobject.href in pand. Zet die tussen dubbele aanhalingstekens, of laat ze bij het laden weg.

Na een grote import draait u ANALYZE zodat de planner de juiste keuzes maakt. Hoe u zo'n database beheert en hostelt, staat bij GIS-beheer en hosting.

PostGIS loont niet als er één gebruiker is die een eenmalige analyse draait. Dan volstaat een GeoPackage met Python voor GIS, en bespaart u het beheer van een databaseserver. Over het bestandsformaat leest u op GeoPackage, en over de conversieprogramma's op GDAL.

Ook als gegevens alleen worden getoond op een kaart en nooit veranderen, is een database niet nodig. Vooraf gebouwde tegels of bestanden achter een webserver zijn eenvoudiger te beheren en goedkoper te hosten. PostGIS komt in beeld bij bewerkingen door meerdere gebruikers, bij koppelingen met andere systemen en bij query's over grote bestanden.

Wie een open data-platform of dataplatform bouwt, heeft meestal beide nodig: PostGIS voor het beheer en een API voor het delen. Dat laatste doet bijvoorbeeld pygeoapi. TechinFlow levert zelf dataplatforms en heeft dus belang bij deze keuze, daarom zijn de afwegingen hier met bronnen onderbouwd.