[gdal-dev] Spatial view in GeoPackage

[Even Rouault]

Bo Victor,

there was indeed a bug in the OGR GPKG driver preventing you to do what 
you wanted. Fixed per https://github.com/OSGeo/gdal/pull/7848

The issue was specific to the fact that the geometry column of the view 
is computed, and not just the simple selection of one of a table.

Note that such GeoPackage will probably only be readable by 
implementations that install Spatialite (or Spatialite compatible 
functions). It might also be prudent to wrap the result of ST_Multi() 
with AsGPB() to force creating a GeoPackage geometry blob instead of a 
Spatialite one. The OGR GPKG driver is tolerant to receiving both types 
of geometries, but technically the content should be a GeoPackage one.

Inserting into gpkg_ogr_contents is optional. Its purpose is to 
accelerate counting the number of features. But for layers < let's say 
100 k features, it probably doesn't impact performance much.

Even

Le 27/05/2023 à 22:08, Bo Victor Thomsen a écrit :
>
> Hi list -
>
> I'm working with making spatial views "visible" for QGIS, i.e. so they 
> look like a read-only table from the data source manager.
>
> However, I have one problem:  QGIS doesn't recognize the crs for the 
> view. I've tracked the problem down to that *OGR* dosn't recognize the 
> crs for the view. Hence this mail...
>
> What I've done:
>
> *- creating the view: A union of parcels with the same owner (CVR-number)
> CREATE VIEW marker_2019_cvr AS
>   SELECT "CVR", SUM(COALESCE(st_area(geom),0.0)) as area_m2, count(*) 
> as antal, st_multi(st_union(geom)) as geom from marker_2019 group by 1;
>
> -- insert or update gpkg_contents
> INSERT OR REPLACE INTO 
> gpkg_contents(table_name,data_type,identifier,min_x,min_y,max_x,max_y,srs_id)
>   SELECT
>     'marker_2019_cvr' AS table_name,
>     'features' AS data_type,
>     'marker_2019_cvr' AS identifier,
>     MIN(ST_MinX(geom)) AS min_x,
>     MIN(ST_MinY(geom)) AS min_y,
>     MAX(ST_MaxX(geom)) AS max_x,
>     MAX(ST_MaxY(geom)) AS max_y,
>     25832 AS srs_id -- same srs as the base layer for the view
>   FROM marker_2019_cvr;
>
> -- insert or update gpkg_geometry_columns
> INSERT OR REPLACE INTO gpkg_geometry_columns ( table_name, column_name 
> , geometry_type_name, srs_id, z, m)
>   VALUES ('marker_2019_cvr', 'geom', 'MULTIPOLYGON', 25832, 0, 0);
>
> -- insert or update gpkg_ogr_columns (necessary ?)
> INSERT OR REPLACE INTO gpkg_ogr_contents(table_name,feature_count)
>   SELECT
>     'marker_2019_cvr' AS table_name,
>     count(*) AS feature_count
>   FROM marker_2019_cvr;
> *
>
> AFAIK, this is the necessary steps to make a spatial view visible in a 
> GeoPackage. And it works; sort of...
>
> The QGIS datasource manager shows the view as a spatial table. I can 
> add it to the QGIS map, but... it can't find the CRS for the layer and 
> defaults to the projects CRS.
>
> Digging a litle bit deeper: Using ogrinfo like this:
>
> *ogrinfo -ro -so data_marker_vandloeb_soe.gpkg marker_2019_cvr*
>
> shows the layers SRS WKT as (unknown).
>
> Using the same command on the base layer *marker_2019 *gives the 
> correct WKT.
>
> 25832 is the correct EPSG number for the CRS. And the definition 
> exists in the *gpkg_spatial_ref_sys* table.
>
> What am I doing wrong ? Or is it just "out of scope" for OGR ??
>
> -- 
> Med venlig hilsen / Best regards
>
> Bo Victor Thomsen
>
> _______________________________________________
> gdal-dev mailing list
> gdal-dev at lists.osgeo.org
> https://lists.osgeo.org/mailman/listinfo/gdal-dev

-- 
http://www.spatialys.com
My software is free, but my time generally not.
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.osgeo.org/pipermail/gdal-dev/attachments/20230527/4065e43e/attachment.htm>