diff --git a/fmtm_splitter/fmtm_algorithm.sql b/fmtm_splitter/fmtm_algorithm.sql index 8b6d116..d8760af 100644 --- a/fmtm_splitter/fmtm_algorithm.sql +++ b/fmtm_splitter/fmtm_algorithm.sql @@ -229,10 +229,10 @@ CREATE TABLE clusteredbuildings AS ( ), -- Cluster the buildings within each splitpolygon. The second term in the - -- call to the ST_ClusterKMeans function is the number of clusters to create, - -- so we're dividing the number of features by a constant (10 in this case) - -- to get the number of clusters required to get close to the right number - -- of features per cluster. + -- call to the ST_ClusterKMeans function is the number of clusters to + -- create, so we're dividing the number of features by a constant + -- (10 in this case) to get the number of clusters required to get close + -- to the right number of features per cluster. -- TODO: This should certainly not be a hardcoded, the number of features -- per cluster should come from a project configuration table buildingstocluster AS ( @@ -380,13 +380,80 @@ CREATE TABLE taskpolygons AS ( ); -- ALTER TABLE taskpolygons ADD PRIMARY KEY(taskid); + +-- Step 1: Identify polygons without any buildings +DROP TABLE IF EXISTS no_building_polygons; +CREATE TABLE no_building_polygons AS ( + SELECT + tp.*, + tp.geom AS no_building_geom + FROM taskpolygons AS tp + LEFT JOIN buildings AS b ON ST_INTERSECTS(tp.geom, b.geom) + WHERE b.geom IS NULL +); + +-- Step 2: Identify neighboring polygons +DROP TABLE IF EXISTS neighboring_polygons; +CREATE TABLE neighboring_polygons AS ( + SELECT + nb.*, + nb.geom AS neighbor_geom, + nb_building_count, + nbp.no_building_geom + FROM taskpolygons AS nb + INNER JOIN no_building_polygons AS nbp + -- Finds polygons that touch each other + ON ST_TOUCHES(nbp.no_building_geom, nb.geom) + LEFT JOIN ( + -- Step 3: Count buildings in the neighboring polygons + SELECT + nb.geom, + COUNT(b.geom) AS nb_building_count + FROM taskpolygons AS nb + LEFT JOIN buildings AS b ON ST_INTERSECTS(nb.geom, b.geom) + GROUP BY nb.geom + ) AS building_counts + ON nb.geom = building_counts.geom +); + +-- Step 4: Find the optimal neighboring polygon to avoid, +-- same polygons with the smallest number of buildings merging into +-- multiple neighboring polygons +DROP TABLE IF EXISTS optimal_neighbors; +CREATE TABLE optimal_neighbors AS ( + SELECT + nbp.no_building_geom, + nbp.neighbor_geom + FROM neighboring_polygons AS nbp + WHERE nbp.nb_building_count = ( + SELECT MIN(nb_building_count) + FROM neighboring_polygons AS np + WHERE np.no_building_geom = nbp.no_building_geom + ) +); + +-- Step 5: Merge the small polygons with their optimal neighboring polygons +UPDATE taskpolygons tp +SET geom = ST_UNION(tp.geom, nbp.no_building_geom) +FROM optimal_neighbors AS nbp +WHERE tp.geom = nbp.neighbor_geom; +DELETE FROM taskpolygons +WHERE geom IN ( + SELECT no_building_geom + FROM no_building_polygons +); + + +DROP TABLE IF EXISTS no_building_polygons; +DROP TABLE IF EXISTS neighboring_polygons; + SELECT POPULATE_GEOMETRY_COLUMNS('public.taskpolygons'::regclass); CREATE INDEX taskpolygons_idx ON taskpolygons USING gist (geom); -- VACUUM ANALYZE taskpolygons; - +-- Generate GeoJSON output SELECT JSONB_BUILD_OBJECT( 'type', 'FeatureCollection', diff --git a/fmtm_splitter/splitter.py b/fmtm_splitter/splitter.py index 1906fb0..9c28e16 100755 --- a/fmtm_splitter/splitter.py +++ b/fmtm_splitter/splitter.py @@ -29,6 +29,7 @@ import geojson import numpy as np from geojson import Feature, FeatureCollection, GeoJSON +from osm_rawdata.postgres import PostgresClient from psycopg2.extensions import connection from shapely.geometry import Polygon, shape from shapely.geometry.geo import mapping @@ -42,7 +43,6 @@ drop_tables, insert_geom, ) -from osm_rawdata.postgres import PostgresClient # Instantiate logger log = logging.getLogger(__name__)