submarcos / django-vectortiles

Mapbox VectorTiles for django, with PostGIS or Python
https://django-vectortiles.readthedocs.io
MIT License
38 stars 11 forks source link

MultiVectorLayer in single SQL query with PostGIS #43

Open nitrag opened 1 year ago

nitrag commented 1 year ago

I wrote this up for our purposes but perhaps you can incorporate it in a cleaner way than I could:


@dataclass
class VectorLayerQueryset:
    layer_name: str
    queryset: QuerySet
    queryset_limit = None
    fields: tuple = None
    geom_name: str = 'geom'

class MultiLayerVectorTile:
    """
    Base Mixin to handle vector tile generation
    """

    vector_tile_content_type = "application/x-protobuf"
    vector_layers: list[VectorLayerQueryset] = None

    vector_tile_extent = 4096  # define tile extent
    vector_tile_buffer = (
        256  # define buffer around tiles (intersected polygon display without borders)
    )
    vector_tile_clip_geom = (
        True  # define if feature geometries should be clipped in tile
    )

    @classmethod
    def get_bounds(cls, x, y, z):
        return mercantile.xy_bounds(x, y, z)

    def get_tile(self, x, y, z):
        # get tile coordinates from x, y and z
        xmin, ymin, xmax, ymax = self.get_bounds(x, y, z)

        queries = []
        params = []
        for layer in self.vector_layers:
            queryset = layer.queryset

            # keep features intersecting tile
            filters = {
                # GeoFuncMixin implicitly transforms to SRID of geom
                f"{layer.geom_name}__intersects": MakeEnvelope(
                    xmin, ymin, xmax, ymax, 3857
                )
            }
            queryset = queryset.filter(**filters)
            # annotate prepared geometry for MVT
            queryset = queryset.annotate(
                geom_prepared=AsMVTGeom(
                    Transform(layer.geom_name, 3857),
                    MakeEnvelope(xmin, ymin, xmax, ymax, 3857),
                    self.vector_tile_extent,
                    self.vector_tile_buffer,
                    self.vector_tile_clip_geom,
                )
            )
            fields = (
                layer.fields + ("geom_prepared",)
                if layer.fields
                else ("geom_prepared",)
            )
            # limit feature number if limit provided
            limit = layer.queryset_limit
            if limit:
                queryset = queryset[:limit]
            # keep values to include in tile (extra included_fields + geometry)
            queryset = queryset.values(*fields)
            # generate SQL/Params for layer
            sql, p = queryset.query.sql_with_params()
            queries.append(
                f"SELECT ST_ASMVT(subquery.*, '{layer.layer_name}', {self.vector_tile_extent}, 'geom_prepared') as mvt "
                f"FROM ({sql}) as subquery"
            )
            params = params + list(p)

        with connection.cursor() as cursor:
            cursor.execute(
                "SELECT string_agg(mvt, '') as mvt FROM ({}) sub".format(' UNION '.join(queries)),
                params=[*params]
            )
            row = cursor.fetchone()[0]
            # psycopg2 returns memoryview, psycopg returns bytes
            return row.tobytes() if type(row) == memoryview else row or b""
submarcos commented 11 months ago

Ok I undestand that you want to handle multiple layers in single SQL query instead of layer per layer. I will think about this as helper for this. Actually I prefer to make layer per layer because I manage my cache policy per layer. Thank you

Edit: I take this in consideration, maybe we could choose in future, in a view, if we want to make a single query or not