OSGeo / gdal

GDAL is an open source MIT licensed translator library for raster and vector geospatial data formats.
https://gdal.org
Other
4.68k stars 2.47k forks source link

ogr2ogr to PGDump with MULTILINESTRINGZMs from FileGDBs drops M value #1327

Closed mojodna closed 5 years ago

mojodna commented 5 years ago

Expected behavior and actual behavior.

Using the nhdflowline layer from https://prd-tnm.s3.amazonaws.com/StagedProducts/Hydrography/NHDPlus/HU4/HighResolution/GDB/NHDPLUS_H_0104_HU4_GDB.zip, M coordinates are dropped when converting to PGDump.

The dimension count is correctly detected (4). However, each vertex only has 3 dimensions (e.g. -70.3263942012797 44.1704549314387 0); Z is 0 (I expect it to be non-zero, but don't yet have a way to confirm, as I don't have an ArcGIS license to extract/inspect features; I'm asking around for help). This discrepancy triggers Dimensions mismatch in lwcollection errors when loading into Postgres unless -dim XYZ is included.

$ ogr2ogr \
  -lco WRITE_EWKT_GEOM=YES \
  -f PGDump \
  /vsistdout/ \
  NHDPLUS_H_0104_HU4_GDB.gdb/ \
  nhdflowline | head -25 | tail -22 | grep -vE "CREATE|ALTER"
SELECT AddGeometryColumn('public','nhdflowline','shape',4269,'MULTILINESTRING',4);
INSERT INTO "public"."nhdflowline" ("shape" , "objectid" , "permanent_identifier", "fdate", "resolution", "gnis_id", "gnis_name", "lengthkm", "reachcode", "flowdir", "wbarea_permanent_identifier", "ftype", "fcode", "mainpath", "innetwork", "visibilityfilter", "shape_length", "nhdplusid", "vpuid", "enabled") VALUES (GeomFromEWKT('SRID=4269;MULTILINESTRING ((-70.3263942012797 44.1704549314387 0,-70.3264092679463 44.1703601981055 0,-70.3262794012799 44.1700241314393 0,-70.3259746012803 44.16955953144 0,-70.3260200679469 44.1693023981071 0,-70.3256998012807 44.1690539981075 0,-70.3256462012808 44.1687889981079 0,-70.3254708012811 44.1683723981085 0,-70.3254326012812 44.1680983314423 0,-70.3252270012815 44.16747039811 0,-70.3251200012817 44.1673185981102 0,-70.3248300679488 44.1664847314448 0,-70.3245936679492 44.1662083314453 0,-70.3243878012828 44.1661071981121 0,-70.3241894012831 44.1659249314457 0,-70.3240904012832 44.1657639981126 0,-70.3239000012836 44.1655861981129 0,-70.3236788012839 44.1652555981134 0,-70.3235336679508 44.1651267314469 0,-70.323373667951 44.165070131447 0,-70.3231680012847 44.165063531447 0,-70.3230382679516 44.1650153981138 0,-70.322748467952 44.1650413314471 0,-70.3224434679525 44.1649323314472 0,-70.3220546012864 44.1649503314472 0,-70.3220242679531 44.164887598114 0,-70.3218718679533 44.1647859314475 0,-70.3217878012868 44.1646743314476 0,-70.3217194012869 44.1643959981147 0))'::TEXT) , 1, '106446373', '2012/03/03 00:34:12', 2, NULL, NULL, 0.848, '01040002004470', 1, NULL, 460, 46006, 0, 1, 50000, 0.00850770588223511, 5000300005409, '0104', 1);

For comparison, ogrinfo sees the M value (e.g. -70.3263942012797 44.1704549314387 0 69.94191), but Z values are still 0:

$ ogrinfo -fid 1 -fields=NO NHDPLUS_H_0104_HU4_GDB.gdb/ nhdflowline
INFO: Open of `NHDPLUS_H_0104_HU4_GDB.gdb/'
      using driver `OpenFileGDB' successful.

Layer name: NHDFlowline
Geometry: 3D Measured Multi Line String
Feature Count: 15872
Extent: (-71.349613, 43.903376) - (-69.850725, 45.338963)
Layer SRS WKT:
GEOGCS["NAD83",
    DATUM["North_American_Datum_1983",
        SPHEROID["GRS 1980",6378137,298.257222101,
            AUTHORITY["EPSG","7019"]],
        TOWGS84[0,0,0,0,0,0,0],
        AUTHORITY["EPSG","6269"]],
    PRIMEM["Greenwich",0,
        AUTHORITY["EPSG","8901"]],
    UNIT["degree",0.0174532925199433,
        AUTHORITY["EPSG","9122"]],
    AUTHORITY["EPSG","4269"]]
FID Column = OBJECTID
Geometry Column = Shape
Permanent_Identifier: String (40.0) NOT NULL
FDate: DateTime (0.0) NOT NULL
Resolution: Integer (0.0) NOT NULL
GNIS_ID: String (10.0)
GNIS_Name: String (65.0)
LengthKM: Real (0.0)
ReachCode: String (14.0)
FlowDir: Integer (0.0) NOT NULL DEFAULT 0
WBArea_Permanent_Identifier: String (40.0)
FType: Integer (0.0) NOT NULL DEFAULT 460
FCode: Integer (0.0) DEFAULT 46003
MainPath: Integer (0.0) NOT NULL DEFAULT 0
InNetwork: Integer(Int16) (0.0) NOT NULL DEFAULT 1
VisibilityFilter: Integer (0.0) NOT NULL DEFAULT 0
Shape_Length: Real (0.0)
NHDPlusID: Real (0.0)
VPUID: String (8.0)
Enabled: Integer(Int16) (0.0)
OGRFeature(NHDFlowline):1
  MULTILINESTRING ZM ((-70.3263942012797 44.1704549314387 0 69.94191,-70.3264092679463 44.1703601981055 0 69.06704,-70.3262794012799 44.1700241314393 0 65.86636,-70.3259746012803 44.16955953144 0 61.15161,-70.3260200679469 44.1693023981071 0 58.77332,-70.3256998012807 44.1690539981075 0 55.66296,-70.3256462012808 44.1687889981079 0 53.20588,-70.3254708012811 44.1683723981085 0 49.21161,-70.3254326012812 44.1680983314423 0 46.68437,-70.3252270012815 44.16747039811 0 40.76493,-70.3251200012817 44.1673185981102 0 39.2029,-70.3248300679488 44.1664847314448 0 31.31574,-70.3245936679492 44.1662083314453 0 28.33707,-70.3243878012828 44.1661071981121 0 26.69017,-70.3241894012831 44.1659249314457 0 24.56506,-70.3240904012832 44.1657639981126 0 22.95,-70.3239000012836 44.1655861981129 0 20.88969,-70.3236788012839 44.1652555981134 0 17.52246,-70.3235336679508 44.1651267314469 0 15.99992,-70.323373667951 44.165070131447 0 14.8218,-70.3231680012847 44.165063531447 0 13.46113,-70.3230382679516 44.1650153981138 0 12.49663,-70.322748467952 44.1650413314471 0 10.56653,-70.3224434679525 44.1649323314472 0 8.31624,-70.3220546012864 44.1649503314472 0 5.74079,-70.3220242679531 44.164887598114 0 5.13128,-70.3218718679533 44.1647859314475 0 3.75844,-70.3217878012868 44.1646743314476 0 2.59345,-70.3217194012869 44.1643959981147 0 0))

Steps to reproduce the problem.

ogr2ogr \
  -lco WRITE_EWKT_GEOM=YES \
  -f PGDump \
  /vsistdout/ \
  NHDPLUS_H_0104_HU4_GDB.gdb/ \
  nhdflowline | head -25 | tail -22 | grep -vE "CREATE|ALTER"

ogrinfo -fid 1 -fields=NO NHDPLUS_H_0104_HU4_GDB.gdb/ nhdflowline

Operating system

MacOS 10.14.3

GDAL version and provenance

2.4.0 from Homebrew.

jratike80 commented 5 years ago

The title is misleading. I believe that error affects only PGDUMP driver because conversion directly into PostGIS creates MULTILINESTRING ZM geometries. Same good result with geopackage.

ogr2ogr -f gpkg test.gpkg NHDPLUS_H_0104_HU4_GDB.gdb nhdflowline

mojodna commented 5 years ago

Against a different FileGDB that I had locally:

spatialite> select AsText(GeomFromGPB(Shape)) from NHDFlowline where OBJECTID=123235;
MULTILINESTRING ZM((-121.385706 43.994888 0 100, -121.385574 43.995108 0 93.9164, -121.385411 43.99544 0 84.95733, -121.385165 43.995945 0 71.3668, -121.384833 43.996498 0 56.046, -121.384578 43.997027 0 41.80569, -121.384407 43.997303 0 34.11488, -121.384222 43.997505 0 27.95647, -121.383994 43.997682 0 21.82034, -121.383628 43.997875 0 13.50255, -121.383081 43.998096 0 1.99649, -121.382986 43.998134 0 0))

It still may be losing the Z value, but I'll open a separate issue for that if/when I'm able to confirm.

rouault commented 5 years ago

I've fixed the issue with PGDump, XYZM and WRITE_EWKT_GEOM=YES per 56f9eca Regarding z==0, I see it also with the proprietary FileGDB driver, and using my low-level dump_gdbtable debugging tool ( https://github.com/rouault/dump_gdbtable ), so I'm pretty confident that this is not a software issue

mojodna commented 5 years ago

Correction: when using EWKB, the error was "ERROR: Column has M dimension but geometry does not"

When generating hex output, it (still) errors with "ERROR: Dimensions mismatch in lwcollection", which corresponds to:

select '01050000A0AD1000000100000001020000C01D0000006C0681A4E39451C048FC9977D115464000000000000000009F71E140487C514034BCB2E3E39451C098B3EB5CCE154640000000000000000087E123624A445140C076FFC2E19451C00808CA59C3154640000000000000000001A4367172775040D8F492C4DC9451C048217120B415464000000000000000008FDFDBF467934E40FC6A4683DD9451C0908073B3AB1546400000000000000000DEAB5626FC624D4020B4FA43D89451C09890B88FA3154640000000000000000030478FDFDBD44B4058242A63D79451C02866BDE09A154640000000000000000043739D465A9A4A4088D37B83D49451C070B60B3A8D1546400000000000000000D7C05609169B484084E342E3D39451C04807023F841546400000000000000000882EA86F99574740EC91E984D09451C0187785AB6F154640000000000000000059DDEA39E9614440F0311FC4CE9451C0A87221B26A1546400000000000000000E02D90A0F89943403CA60D04CA9451C0606B265F4F1546400000000000000000FF092E56D4503F4024158524C69451C0A8F08950461546400000000000000000757632384A563C40C86E0DC5C29451C0C0C22B0043154640000000000000000029AE2AFBAEB03A40200EE784BF9451C0381335073D15464000000000000000001B47ACC5A7903840609DAAE5BD9451C0B85F33C13715464000000000000000003333333333F33640FC2B12C7BA9451C028C8B4ED31154640000000000000000008944DB9C2E33440107B4A27B79451C060B06E18271546400000000000000000D97745F0BF85314038C08EC6B49451C048A16BDF2215464000000000000000008E3BA583F5FF2F40300F7827B29451C0400EA004211546400000000000000000A2B437F8C2A42D406828D7C8AE9451C0E0A942CD201546400000000000000000B2683A3B19EC2A40580DB3A8AC9451C0181E7D391F1546400000000000000000CA6C904946FE284008AC30E9A79451C0D8860813201546400000000000000000E0675C38102225408C6AEDE9A29451C038CFAC801C1546400000000000000000C45F9335EAA120404C03E78A9C9451C0D883AB171D1546400000000000000000BBD05CA791F6164020D8AC0B9C9451C030B56C091B1546400000000000000000BC79AA436E8614402C97768C999451C0303495B41715464000000000000000006FD39FFD48110E40D4B1DC2B989451C070076A0C141546400000000000000000567DAEB662BF044008C1F80C979451C080BF95ED0A15464000000000000000000000000000000000'::geometry;

I'll open another issue for that (and hopefully provide a PR, since your patch provided me with some context to go on).

@themapsmith opened the FileGDB in ArcGIS as well and confirms that the Z value is always 0. Thanks to both of you for looking with different tools!

mojodna commented 5 years ago

Not a bug per se; setting -lco POSTGIS_VERSION=2.2 makes the hex-related error go away.