[gdal-dev] OGR SQL CAST geometry example

[Even Rouault]

Ah, I misunderstood and thought you wanted to extract points from 
linestrings.

Given

$ cat test.csv
id,my_wkt
1,"POINT Z (1 2 3)"

you can do:

$ ogrinfo test.csv -sql "select id, cast(my_wkt as geometry(pointz, 
4326)) from test"
INFO: Open of `test.csv'
       using driver `CSV' successful.

Layer name: test
Geometry: 3D Point
Feature Count: 1
Extent: (1.000000, 2.000000) - (1.000000, 2.000000)
Layer SRS WKT:
GEOGCRS["WGS 84",
     ENSEMBLE["World Geodetic System 1984 ensemble",
         MEMBER["World Geodetic System 1984 (Transit)"],
         MEMBER["World Geodetic System 1984 (G730)"],
         MEMBER["World Geodetic System 1984 (G873)"],
         MEMBER["World Geodetic System 1984 (G1150)"],
         MEMBER["World Geodetic System 1984 (G1674)"],
         MEMBER["World Geodetic System 1984 (G1762)"],
         MEMBER["World Geodetic System 1984 (G2139)"],
         MEMBER["World Geodetic System 1984 (G2296)"],
         ELLIPSOID["WGS 84",6378137,298.257223563,
             LENGTHUNIT["metre",1]],
         ENSEMBLEACCURACY[2.0]],
     PRIMEM["Greenwich",0,
         ANGLEUNIT["degree",0.0174532925199433]],
     CS[ellipsoidal,2],
         AXIS["geodetic latitude (Lat)",north,
             ORDER[1],
             ANGLEUNIT["degree",0.0174532925199433]],
         AXIS["geodetic longitude (Lon)",east,
             ORDER[2],
             ANGLEUNIT["degree",0.0174532925199433]],
     USAGE[
         SCOPE["Horizontal component of 3D system."],
         AREA["World."],
         BBOX[-90,-180,90,180]],
     ID["EPSG",4326]]
Data axis to CRS axis mapping: 2,1
Geometry Column = my_wkt
id: String (0.0)
OGRFeature(test):1
   id (String) = 1
   POINT Z (1 2 3)

Note that probably counterintuitively, if instead of pointz, you try to 
specify point (or even linestring), that will only affect the *declared* 
geometry column type, but will have no consequence on the actual feature 
geometry, and no consistency check is done to ensure they match  
(contrary to PostGIS that is very strict, and would not even allow 
casting that POINT Z (1 2 3) as a geometry(point), but only a 
geometry(pointz))


Le 31/05/2024 à 01:29, Dan Jacobson a écrit :
>>>>>> "ER" == Even Rouault <even.rouault at spatialys.com> writes:
> ER> you can't do operations on geometries with OGR SQL...
>
> OK, I'll try it. Thanks.
> The page should still have an example of actual use of
>    Casting ... POINT[Z], LINESTRING[Z]
>
> Without examples people will try things like POINTZ, "POINTS[Z]"

-- 
http://www.spatialys.com
My software is free, but my time generally not.