From 28f6bf91aadcd9c6e93d2aec97ce3e040e8eae44 Mon Sep 17 00:00:00 2001 From: ThorodanBrom Date: Mon, 16 Jun 2025 17:05:14 +0530 Subject: [PATCH 1/2] Bug-Fix : Fix get tile method for MVT Postgres for WorldCRS84Quad TMS - `ST_TileEnvelope` PostGIS method did not work for WorldCRS84Quad TMS, the MVTs generated with it were not valid tiles for the TMS. This is probably because `ST_TileEnvelope` works properly with WebMercatorQuad, and WebMercatorQuad and WorldCRS84Quad have different tile matrix sizes at the same zoom levels. - Hence instead we manually find the coordinates of the tile envelope based on the zoom level, x and y and use this in `ST_MakeEnvelope` to create the tile envelope for the `ST_AsMVT` query. --- pygeoapi/provider/mvt_postgresql.py | 24 +++++++++++++++++++----- 1 file changed, 19 insertions(+), 5 deletions(-) diff --git a/pygeoapi/provider/mvt_postgresql.py b/pygeoapi/provider/mvt_postgresql.py index c8fa2fc1a..cb36013a2 100644 --- a/pygeoapi/provider/mvt_postgresql.py +++ b/pygeoapi/provider/mvt_postgresql.py @@ -168,11 +168,24 @@ def get_tiles(self, layer=None, tileset=None, if not self.is_in_limits(TileMatrixSetEnum.WORLDCRS84QUAD.value, z, x, y): # noqa raise ProviderTileNotFoundError + # get tile size in degrees based on zoom level. + # Tile size is a function of the zoom level, + # e.g at zoom level 0, tile size is 180, + # at zoom level 1, tile size is 90 and so on. + tile_size_deg = 180/pow(2, int(z)) + + # getting top-left coordinates of the tile + xmin = (tile_size_deg * int(x)) - 180 + ymax = (-tile_size_deg * int(y)) + 90 + + # getting bottom-right coordinates of the tile + xmax = xmin + tile_size_deg + ymin = ymax - tile_size_deg + query = text(""" WITH bounds AS ( - SELECT ST_TileEnvelope(:z, :x, :y, - 'SRID=4326;POLYGON((-180 -90,-180 90,180 90,180 -90,-180 -90))'::geometry) AS boundgeom + SELECT ST_MakeEnvelope(:xmin,:ymin,:xmax,:ymax, 4326) AS boundgeom ), mvtgeom AS ( SELECT ST_AsMVTGeom(ST_CurveToLine({geom}), bounds.boundgeom) AS geom {fields} @@ -184,9 +197,10 @@ def get_tiles(self, layer=None, tileset=None, with self.postgres._engine.connect() as session: result = session.execute(query, { - 'z': z, - 'y': y, - 'x': x + 'xmin': xmin, + 'ymin': ymin, + 'xmax': xmax, + 'ymax': ymax }).fetchone() if len(bytes(result[0])) == 0: From 93e05f0a7a5ad67e4f0cd2ac2337c11fc62ba313 Mon Sep 17 00:00:00 2001 From: ThorodanBrom Date: Wed, 18 Jun 2025 11:51:26 +0530 Subject: [PATCH 2/2] Use format strings instead of parameterized queries since the queries for the different TMSs have different parameters --- pygeoapi/provider/mvt_postgresql.py | 15 +++++---------- 1 file changed, 5 insertions(+), 10 deletions(-) diff --git a/pygeoapi/provider/mvt_postgresql.py b/pygeoapi/provider/mvt_postgresql.py index cb36013a2..72ca719d6 100644 --- a/pygeoapi/provider/mvt_postgresql.py +++ b/pygeoapi/provider/mvt_postgresql.py @@ -154,7 +154,7 @@ def get_tiles(self, layer=None, tileset=None, query = text(""" WITH bounds AS ( - SELECT ST_TileEnvelope(:z, :x, :y) AS boundgeom + SELECT ST_TileEnvelope({z}, {x}, {y}) AS boundgeom ), mvtgeom AS ( SELECT ST_AsMVTGeom(ST_Transform(ST_CurveToLine({geom}), 3857), bounds.boundgeom) AS geom {fields} @@ -162,7 +162,7 @@ def get_tiles(self, layer=None, tileset=None, WHERE ST_Intersects({geom}, ST_Transform(bounds.boundgeom, 4326)) ) SELECT ST_AsMVT(mvtgeom, 'default') FROM mvtgeom; - """.format(geom=self.geom, table=self.table, fields=fields)) # noqa + """.format(geom=self.geom, table=self.table, fields=fields, z=z, x=x, y=y)) # noqa if tileset == TileMatrixSetEnum.WORLDCRS84QUAD.value.tileMatrixSet: if not self.is_in_limits(TileMatrixSetEnum.WORLDCRS84QUAD.value, z, x, y): # noqa @@ -185,7 +185,7 @@ def get_tiles(self, layer=None, tileset=None, query = text(""" WITH bounds AS ( - SELECT ST_MakeEnvelope(:xmin,:ymin,:xmax,:ymax, 4326) AS boundgeom + SELECT ST_MakeEnvelope({xmin},{ymin},{xmax},{ymax}, 4326) AS boundgeom ), mvtgeom AS ( SELECT ST_AsMVTGeom(ST_CurveToLine({geom}), bounds.boundgeom) AS geom {fields} @@ -193,15 +193,10 @@ def get_tiles(self, layer=None, tileset=None, WHERE ST_Intersects({geom}, bounds.boundgeom) ) SELECT ST_AsMVT(mvtgeom, 'default') FROM mvtgeom; - """.format(geom=self.geom, table=self.table, fields=fields)) # noqa + """.format(geom=self.geom, table=self.table, fields=fields, xmin=xmin, ymin=ymin, xmax=xmax, ymax=ymax)) # noqa with self.postgres._engine.connect() as session: - result = session.execute(query, { - 'xmin': xmin, - 'ymin': ymin, - 'xmax': xmax, - 'ymax': ymax - }).fetchone() + result = session.execute(query).fetchone() if len(bytes(result[0])) == 0: return None