Variants to generate an INSPIRE-compliant grid in EPSG:3035 for Germany's federal state Northrhine-Westfalia in PostGIS

nw_3035.md

Variants to generate an INSPIRE-compliant grid in EPSG:3035 for Germany's federal state Northrhine-Westfalia in PostGIS

Get the polygons, which have been released under a zero license: https://www.opengeodata.nrw.de/produkte/geobasis/vkg/dvg/dvg1/dvg1_EPSG25832_Shape.zip

variant 1

• execution time: 2203851,183 ms (36:43,851)
• explain analyze: https://explain.depesz.com/s/5X2O

expand for SQL query

CREATE TEMPORARY TABLE dvg1bld_nw_3035 AS (
    SELECT
        kn::text COLLATE "C",
        ST_Transform(geom, 3035) AS geom
    FROM
        dvg1bld_nw
);
ALTER TABLE dvg1bld_nw_3035 SET (parallel_workers = 32);
CREATE INDEX idx_dvg1bld_nw_3035_spgist ON dvg1bld_nw_3035 USING spgist(geom);
VACUUM ANALYZE dvg1bld_nw_3035;

CREATE TEMPORARY TABLE slices AS ( 
    SELECT
        kn,
        ST_Subdivide(geom) as geom
    FROM
        dvg1bld_nw_3035
);
ALTER TABLE slices SET (parallel_workers = 32);
CREATE INDEX idx_slices_spgist ON slices USING spgist(geom);
CREATE INDEX idx_slices_kn1 ON slices USING btree(kn);
CREATE INDEX idx_slices_kn2 ON slices USING gist(kn, geom);
CREATE INDEX idx_slices_kn3 ON slices USING gist(geom, kn);
VACUUM ANALYZE slices;

CREATE TEMPORARY TABLE grid AS ( 
    SELECT
        (ST_SquareGrid(100, geom)).geom
    FROM
        dvg1bld_nw_3035
);
ALTER TABLE grid SET (parallel_workers = 32);
CREATE INDEX idx_grid_spgist ON grid USING spgist(geom);
VACUUM ANALYZE grid;

SELECT
    CASE WHEN ST_Covers(t1.geom, t2.geom)
     THEN 10000
     ELSE ST_Area(ST_Intersection(t2.geom, t3.geom))
    END
FROM
    slices t1
    JOIN grid t2 ON ST_Intersects(t1.geom, t2.geom)
    JOIN dvg1bld_nw_3035 t3 USING (kn);

variant 2

• execution time: 2921094,539 ms (48:41,095)
• explain analyze: https://explain.depesz.com/s/1gkS

expand for SQL query

WITH
dvg1bld_nw_3035 AS MATERIALIZED (
    SELECT
        kn,
        ST_Transform(geom, 3035) AS geom
    FROM
        dvg1bld_nw
),
slices AS MATERIALIZED (
    SELECT
        kn,
        ST_Subdivide(geom) as geom
    FROM
        dvg1bld_nw_3035
),
grid AS MATERIALIZED (
    SELECT
        (ST_SquareGrid(100, geom)).geom
    FROM
        dvg1bld_nw_3035
)
SELECT
    CASE WHEN ST_Covers(t1.geom, t2.geom)
     THEN 10000
     ELSE ST_Area(ST_Intersection(t2.geom, t3.geom))
    END
FROM
    slices t1
    JOIN grid t2 ON ST_Intersects(t1.geom, t2.geom)
    JOIN dvg1bld_nw_3035 t3 USING (kn);

variant 3

• execution time: 2840714,346 ms (47:20,714)
• explain analyze: https://explain.depesz.com/s/54xp

expand for SQL query

WITH
dvg1bld_nw_3035 AS MATERIALIZED (
    SELECT
        kn,
        ST_Transform(geom, 3035) AS geom
    FROM
        dvg1bld_nw
),
slices AS (
    SELECT
        kn,
        ST_Subdivide(geom) as geom
    FROM
        dvg1bld_nw_3035
),
grid AS (
    SELECT
        (ST_SquareGrid(100, geom)).geom
    FROM
        dvg1bld_nw_3035
)
SELECT
    CASE WHEN ST_Covers(t1.geom, t2.geom)
     THEN 10000
     ELSE ST_Area(ST_Intersection(t2.geom, t3.geom))
    END
FROM
    slices t1
    JOIN grid t2 ON ST_Intersects(t1.geom, t2.geom)
    JOIN dvg1bld_nw_3035 t3 USING (kn);

variant 4

• execution time: 4028568,971 ms (01:07:08,569)
• explain analyze: https://explain.depesz.com/s/aPv

expand for SQL query

CREATE TEMPORARY TABLE dvg1bld_nw_3035 AS (
    SELECT
        kn::text COLLATE "C",
        ST_Transform(geom, 3035) AS geom
    FROM
        dvg1bld_nw
);
ALTER TABLE dvg1bld_nw_3035 SET (parallel_workers = 32);
CREATE INDEX idx_dvg1bld_nw_3035_spgist ON dvg1bld_nw_3035 USING spgist(geom);
VACUUM ANALYZE dvg1bld_nw_3035;

WITH
slices AS (
    SELECT
        kn,
        ST_Subdivide(geom) as geom
    FROM
        dvg1bld_nw_3035
),
grid AS (
    SELECT
        (ST_SquareGrid(100, geom)).geom
    FROM
        dvg1bld_nw_3035
)
SELECT
    CASE WHEN ST_Covers(t1.geom, t2.geom)
     THEN 10000
     ELSE ST_Area(ST_Intersection(t2.geom, t3.geom))
    END
FROM
    slices t1
    JOIN grid t2 ON ST_Intersects(t1.geom, t2.geom)
    JOIN dvg1bld_nw_3035 t3 USING (kn);

variant 5

• execution time: 2224610,286 ms (37:04,610)
• explain analyze: https://explain.depesz.com/s/cZGz

expand for SQL query

CREATE TEMPORARY TABLE slices AS ( 
    WITH
    dvg1bld_nw_3035 AS MATERIALIZED (
        SELECT
            kn,
            ST_Transform(geom, 3035) AS geom
        FROM
            dvg1bld_nw
    )
    SELECT
        kn,
        ST_Subdivide(geom) as geom
    FROM
        dvg1bld_nw_3035
);
ALTER TABLE slices SET (parallel_workers = 32);
CREATE INDEX idx_slices_spgist ON slices USING spgist(geom);
CREATE INDEX idx_slices_kn1 ON slices USING btree(kn);
CREATE INDEX idx_slices_kn2 ON slices USING gist(kn, geom);
CREATE INDEX idx_slices_kn3 ON slices USING gist(geom, kn);
VACUUM ANALYZE slices;

WITH
dvg1bld_nw_3035 AS MATERIALIZED (
    SELECT
        kn,
        ST_Transform(geom, 3035) AS geom
    FROM
        dvg1bld_nw
),
grid AS (
    SELECT
        (ST_SquareGrid(100, geom)).geom
    FROM
        dvg1bld_nw_3035
)
SELECT
    CASE WHEN ST_Covers(t1.geom, t2.geom)
     THEN 10000
     ELSE ST_Area(ST_Intersection(t2.geom, t3.geom))
    END
FROM
    slices t1
    JOIN grid t2 ON ST_Intersects(t1.geom, t2.geom)
    JOIN dvg1bld_nw_3035 t3 USING (kn);

variant 6

• execution time: 2151721,188 ms (35:51,721)
• explain analyze: https://explain.depesz.com/s/7lqV

expand for SQL query

CREATE TEMPORARY TABLE grid AS ( 
    WITH
    dvg1bld_nw_3035 AS MATERIALIZED (
        SELECT
            kn,
            ST_Transform(geom, 3035) AS geom
        FROM
            dvg1bld_nw
    )
    SELECT
        (ST_SquareGrid(100, geom)).geom
    FROM
        dvg1bld_nw_3035
);
ALTER TABLE grid SET (parallel_workers = 32);
CREATE INDEX idx_grid_spgist ON grid USING spgist(geom);
VACUUM ANALYZE grid;

WITH
dvg1bld_nw_3035 AS MATERIALIZED (
    SELECT
        kn,
        ST_Transform(geom, 3035) AS geom
    FROM
        dvg1bld_nw
),
slices AS MATERIALIZED (
    SELECT
        kn,
        ST_Subdivide(geom) as geom
    FROM
        dvg1bld_nw_3035
)
SELECT
    CASE WHEN ST_Covers(t1.geom, t2.geom)
     THEN 10000
     ELSE ST_Area(ST_Intersection(t2.geom, t3.geom))
    END
FROM
    slices t1
    JOIN grid t2 ON ST_Intersects(t1.geom, t2.geom)
    JOIN dvg1bld_nw_3035 t3 USING (kn);

variant 1 - modified

• execution time: 2307018,751 ms (38:27,019)
• explain analyze: https://explain.depesz.com/s/D7uv

expand for SQL query

CREATE TEMPORARY TABLE dvg1bld_nw_3035 AS (
    SELECT
        kn::text COLLATE "C",
        ST_Transform(geom, 3035) AS geom
    FROM
        dvg1bld_nw
);
ALTER TABLE dvg1bld_nw_3035 SET (parallel_workers = 32);
CREATE INDEX idx_dvg1bld_nw_3035_spgist ON dvg1bld_nw_3035 USING spgist(geom);
VACUUM ANALYZE dvg1bld_nw_3035;

CREATE TEMPORARY TABLE slices AS ( 
    SELECT
        kn,
        ST_Subdivide(geom) as geom
    FROM
        dvg1bld_nw_3035
);
ALTER TABLE slices SET (parallel_workers = 32);
CREATE INDEX idx_slices_spgist ON slices USING spgist(geom);
CREATE INDEX idx_slices_kn1 ON slices USING btree(kn);
CREATE INDEX idx_slices_kn2 ON slices USING gist(kn, geom);
CREATE INDEX idx_slices_kn3 ON slices USING gist(geom, kn);
VACUUM ANALYZE slices;

CREATE TEMPORARY TABLE grid AS ( 
    SELECT
        (ST_SquareGrid(100, geom)).geom
    FROM
        dvg1bld_nw_3035
);
ALTER TABLE grid SET (parallel_workers = 32);
CREATE INDEX idx_grid_spgist ON grid USING spgist(geom);
VACUUM ANALYZE grid;

SELECT
    CASE WHEN ST_Covers(t3.geom, t2.geom)
     THEN 10000
     ELSE ST_Area(ST_Intersection(t2.geom, t3.geom))
    END
FROM
    slices t1
    JOIN grid t2 ON ST_Intersects(t1.geom, t2.geom)
    JOIN dvg1bld_nw_3035 t3 USING (kn);

variant 2 - modified

• execution time: 821463,687 ms (13:41,464)
• explain analyze: https://explain.depesz.com/s/oaRL

expand for SQL query

WITH
dvg1bld_nw_3035 AS MATERIALIZED (
    SELECT
        kn,
        ST_Transform(geom, 3035) AS geom
    FROM
        dvg1bld_nw
),
slices AS MATERIALIZED (
    SELECT
        kn,
        ST_Subdivide(geom) as geom
    FROM
        dvg1bld_nw_3035
),
grid AS MATERIALIZED (
    SELECT
        (ST_SquareGrid(100, geom)).geom
    FROM
        dvg1bld_nw_3035
)
SELECT
    CASE WHEN ST_Covers(t3.geom, t2.geom)
     THEN 10000
     ELSE ST_Area(ST_Intersection(t2.geom, t3.geom))
    END
FROM
    slices t1
    JOIN grid t2 ON ST_Intersects(t1.geom, t2.geom)
    JOIN dvg1bld_nw_3035 t3 USING (kn);

variant 3 - modified

• execution time: 888230,775 ms (14:48,231)
• explain analyze: https://explain.depesz.com/s/x9YE

expand for SQL query

WITH
dvg1bld_nw_3035 AS MATERIALIZED (
    SELECT
        kn,
        ST_Transform(geom, 3035) AS geom
    FROM
        dvg1bld_nw
),
slices AS (
    SELECT
        kn,
        ST_Subdivide(geom) as geom
    FROM
        dvg1bld_nw_3035
),
grid AS (
    SELECT
        (ST_SquareGrid(100, geom)).geom
    FROM
        dvg1bld_nw_3035
)
SELECT
    CASE WHEN ST_Covers(t3.geom, t2.geom)
     THEN 10000
     ELSE ST_Area(ST_Intersection(t2.geom, t3.geom))
    END
FROM
    slices t1
    JOIN grid t2 ON ST_Intersects(t1.geom, t2.geom)
    JOIN dvg1bld_nw_3035 t3 USING (kn);

variant 4 - modified

• execution time: (still running)
• explain analyze: (still running)

expand for SQL query

CREATE TEMPORARY TABLE dvg1bld_nw_3035 AS (
    SELECT
        kn::text COLLATE "C",
        ST_Transform(geom, 3035) AS geom
    FROM
        dvg1bld_nw
);
ALTER TABLE dvg1bld_nw_3035 SET (parallel_workers = 32);
CREATE INDEX idx_dvg1bld_nw_3035_spgist ON dvg1bld_nw_3035 USING spgist(geom);
VACUUM ANALYZE dvg1bld_nw_3035;

WITH
slices AS (
    SELECT
        kn,
        ST_Subdivide(geom) as geom
    FROM
        dvg1bld_nw_3035
),
grid AS (
    SELECT
        (ST_SquareGrid(100, geom)).geom
    FROM
        dvg1bld_nw_3035
)
SELECT
    CASE WHEN ST_Covers(t3.geom, t2.geom)
     THEN 10000
     ELSE ST_Area(ST_Intersection(t2.geom, t3.geom))
    END
FROM
    slices t1
    JOIN grid t2 ON ST_Intersects(t1.geom, t2.geom)
    JOIN dvg1bld_nw_3035 t3 USING (kn);

variant 5 - modified

• execution time: 1562439,919 ms (26:02,440)
• explain analyze: https://explain.depesz.com/s/lSiZ

expand for SQL query

CREATE TEMPORARY TABLE slices AS ( 
    WITH
    dvg1bld_nw_3035 AS MATERIALIZED (
        SELECT
            kn,
            ST_Transform(geom, 3035) AS geom
        FROM
            dvg1bld_nw
    )
    SELECT
        kn,
        ST_Subdivide(geom) as geom
    FROM
        dvg1bld_nw_3035
);
ALTER TABLE slices SET (parallel_workers = 32);
CREATE INDEX idx_slices_spgist ON slices USING spgist(geom);
CREATE INDEX idx_slices_kn1 ON slices USING btree(kn);
CREATE INDEX idx_slices_kn2 ON slices USING gist(kn, geom);
CREATE INDEX idx_slices_kn3 ON slices USING gist(geom, kn);
VACUUM ANALYZE slices;

WITH
dvg1bld_nw_3035 AS MATERIALIZED (
    SELECT
        kn,
        ST_Transform(geom, 3035) AS geom
    FROM
        dvg1bld_nw
),
grid AS (
    SELECT
        (ST_SquareGrid(100, geom)).geom
    FROM
        dvg1bld_nw_3035
)
SELECT
    CASE WHEN ST_Covers(t3.geom, t2.geom)
     THEN 10000
     ELSE ST_Area(ST_Intersection(t2.geom, t3.geom))
    END
FROM
    slices t1
    JOIN grid t2 ON ST_Intersects(t1.geom, t2.geom)
    JOIN dvg1bld_nw_3035 t3 USING (kn);

variant 6 - modified

• execution time: 1581953,021 ms (26:21,953)
• explain analyze: https://explain.depesz.com/s/xRxm

expand for SQL query

CREATE TEMPORARY TABLE grid AS ( 
    WITH
    dvg1bld_nw_3035 AS MATERIALIZED (
        SELECT
            kn,
            ST_Transform(geom, 3035) AS geom
        FROM
            dvg1bld_nw
    )
    SELECT
        (ST_SquareGrid(100, geom)).geom
    FROM
        dvg1bld_nw_3035
);
ALTER TABLE grid SET (parallel_workers = 32);
CREATE INDEX idx_grid_spgist ON grid USING spgist(geom);
VACUUM ANALYZE grid;

WITH
dvg1bld_nw_3035 AS MATERIALIZED (
    SELECT
        kn,
        ST_Transform(geom, 3035) AS geom
    FROM
        dvg1bld_nw
),
slices AS MATERIALIZED (
    SELECT
        kn,
        ST_Subdivide(geom) as geom
    FROM
        dvg1bld_nw_3035
)
SELECT
    CASE WHEN ST_Covers(t3.geom, t2.geom)
     THEN 10000
     ELSE ST_Area(ST_Intersection(t2.geom, t3.geom))
    END
FROM
    slices t1
    JOIN grid t2 ON ST_Intersects(t1.geom, t2.geom)
    JOIN dvg1bld_nw_3035 t3 USING (kn);