Add gdal vector unnest, or convert gdal vector explode-collections into gdal vector explode?

[dbaston]

Feature description

I implemented a pipeline step to convert a single feature with array field(s) into multiple features with scalar field(s). You can convert multiple array fields, provided that they are all the same length in every feature. An example use case is with gdal raster zonal-stats, to get a table of the pixels covered by a polygon:

echo -e "wkt\n\"POINT (440720 3751320)\"" | gdal pipeline read --input-format CSV --oo GEOM_POSSIBLE_NAMES=wkt /vsistdin/ ! buffer 150 ! zonal-stats --input autotest/gcore/data/byte.tif --zones _PIPE_ --pixels fractional --stat center_x --stat center_y --stat values --stat coverage ! unnest --field center_x --field center_y --field values --field coverage ! write --of CSV /vsistdout/

which outputs

center_x,center_y,values,coverage
440750,3751290,107,1
440810,3751290,123,1
440870,3751290,132,0.423658132553101
440750,3751230,115,1
440810,3751230,132,0.88035124540329
440870,3751230,107,0.07479527592659
440750,3751170,115,0.423658132553101
440810,3751170,132,0.07479527592659

Then I started wondering if gdal vector explode-collections should instead be reworked to allow exploding geometry and/or array attribute fields, subject to the same constraint that everything being exploded must have the same number of parts/elements.

So

`gdal vector explode --geom`
`gdal vector explode --field array_1 --field array_2`
`gdal vector explode --geom --field array_1` 

But this turns out to be at odds from the current behavior of gdal vector explode-collections, where we get the cartesian product of the exploded fields rather than "zipping".

Any thoughts?

Additional context

The name "unnest" comes from the PostgreSQL function that does the same thing.

[rouault]

Then I started wondering if gdal vector explode-collections should instead be reworked to allow exploding geometry and/or array attribute

could potentially be interesting for the S-101 driver that natively produces features like

{
        "type": "Feature",
        "geometry": {
            "type": "GeometryCollection",
            "geometries": [
                {"type": "MultiPoint", "coordinates": [[-2.0, -49.0], [-2.01, -49.1]]},
                {"type": "MultiPoint", "coordinates": [[-2.0, -49.0], [-2.01, -49.1]]},
            ],
        },
        "properties": {
            "recordId": 6,
            "recordVersion": 1,
            "producingAgency": 12345,
            "featureIdentificationNumber": 66,
            "featureIdentificationSubdivision": 1,
            "geometryLayerName": ["MultiPoint2D", "MultiPoint2D"],
            "geometryRecordId": [1, 1],
            "scaleMinimum": [1000, 3000],
            "scaleMaximum": [2000, 4000],
        },
        "id": 1,
    }

where we get the cartesian product of the exploded fields rather than "zipping".

have both mode available: --mode=cartesian-product|zip ?

[dbaston]

have both mode available: --mode=cartesian-product|zip ?

Possible of course, it's just a bit uglier to implement. What is the use case for cartesian product?

[rouault]

No actual use case for cartesian product for explode collections but i didn't see much alternative to be able to process arbitrary datasets without restriction on their content.well an alternative could be to explode thé first geometry column and set other geometry columns to null and then proceed simlilarly with other geometry columns

[dbaston]

I've updated gdal vector unnest to handle any combination of attribute and geometry fields, using the zipping method. I would propose renaming it to gdal vector explode, as a more general replacement for gdal vector explode-collections.

The interface is fairly minimal:

  --field <FIELD>                                      Attribute fields(s) to explode [may be repeated]
  --geom-field <GEOM_FIELD>                            Geometry field(s) to explode [may be repeated]
  --index-field <INDEX-FIELD>                          Name of the output index field (default: )

Any field that is not specified is left as-is. If you need the cartesian product, you can get it by calling unnest twice in a pipeline, e.g.

gdal.alg.vector.pipeline(pipeline=f'read {src_fname} ! unnest --geom-field 0 ! unnest --geom-field 1 ! write {dst_fname}')

(The above accepts an index for geometry fields in addition to the geometry field name -- I would like to apply this to other utilities as well)

gdal vector explode-collections also has a --geometry-type and --skip-on-type-mismatch flags. These complicate the implementation (only slightly) but it's not entirely obvious how they should work with multiple geometry columns (I guess they apply to all of them?). I think it would be more clear to remove these options and use a pipeline to achieve the desired output:

gdal.alg.vector.pipeline(pipeline=f'read {src_fname} ! unnest --geom-field 0 ! set-geom-type --geometry-type LINESTRING --skip ! write {dst_fname}')