problem with python, error: Never get here.

[automataIA]

Hi.

The problem I have is when I run this line of code

raster_file = 'tif/33TWF.tif'  # GeoTiff
vector_file = 'shp/33TWF.shp' #  Shapefile extracted with pgsql2shp
stats = exact_extract(raster_file, vector_file, ["count", "mean"])`

gives me this error:

............../lib/python3.10/site-packages/exactextract/exact_extract.py", line 320, in exact_extract
    processor.process()
RuntimeError: Never get here.

this is the link to the raster (if it helps): https://we.tl/t-Vn66czZl8w

Context: I'm extracting a shapefile from a postigis DB with ewkb and it had never given me this problem the previous times. I checked both SRIDs (UTM coordinates) 32633 of the two files and they match.

[dbaston]

Are you able to share the shapefile?

[automataIA]

Are you able to share the shapefile?

For now, I can't. But I found that some particular polygons (which I don't know what they are), give me error (perhaps in the presence of cusps in the vertices), even if I checked that there were no empty and invalid polygons via postgis query. Do you have any tips for cleaning up 'crooked' polygons?

[dbaston]

If you run the same process from the command line, the progress output may help you identify the problem polygons. Or in Python, instead of providing the filename to exact_extract, you could provide a GeoPandas data frame with individual polygons until you find the problem. But if you're not able to share the polygon that's causing the problem, I don't think there is much I can do to fix it.

[automataIA]

Now, after applying a postgis buffer of 0 (to smooth the vertices), the error changed, and it gives me this error:

lib/python3.10/site-packages/exactextract/exact_extract.py", line 320, in exact_extract
    processor.process()
RuntimeError: Error calling GEOSCoordSeq_isCCW_r.

[dbaston]

To me that suggests coordinates with a value of NaN or Inf.

[automataIA]

If you run the same process from the command line, the progress output may help you identify the problem polygons.

Or in Python, instead of providing the filename to exact_extract, you could provide a GeoPandas data frame with individual polygons until you find the problem.

1. I tried both methods, and they don't give me extra output.
2. The Windows version doesn't give me any problems. p.s. Could you put a function that skips the problematic lines and puts nan or 0 as the value in that line?

[dbaston]

I tried both methods, and they don't give me extra output.

There is no extra output from the GeoPandas version. What I was suggesting was a way to subset your data to identify the problem polygon.

Alternatively, to get extra output from the CLI, use --progress with --include-col, e.g:

./exactextract --polygons concelhos.gpkg --raster clc2018_v2020_20u1.tif --include-col name --output /tmp/out.csv --progress --stat mode

This will display progress information:

Processing Lagoa.
Processing Nordeste.
Processing Ponta Delgada.
Processing Povoação.
Processing Ribeira Grande.
Processing Vila Franca do Campo.

so the last "Processing" entry before the error can be used to identify the input that is causing the problem.

[automataIA]

./exactextract --polygons concelhos.gpkg --raster clc2018_v2020_20u1.tif --include-col name --output /tmp/out.csv --progress --stat mode

i tried this: exactextract --polygons "shp/zs.shp" --raster "tif/zs.tif" --output /tmp/out.csv --progress --stat mode and it gives me this error: Processing Error: Unexpected field name: .

[dbaston]

You need to use --include-col, see above

[automataIA]

You need to use --include-col, see above

ty. now from this: exactextract --polygons shp/zsd.shp --raster tif/zzs.tif --include-col ID --output /tmp/out.csv --progress --stat mode

it gives me this error:

ERROR 1: SQL Expression Parsing Error: syntax error, unexpected '-', expecting end of string. Occurred around :
SELECT ID FROM zs
                                 ^
ERROR 10: Pointer 'hLayer' is NULL in 'OGR_L_GetSpatialRef'.

ERROR 10: Pointer 'hLayer' is NULL in 'OGR_L_GetLayerDefn'.

[4]    1227422 segmentation fault  exactextract --polygons shp/zs.shp --raster   ID

What should I do?

[dbaston]

What does ogrinfo -al -so shp/zsd.shp return?

[automataIA]

ogrinfo -al -so shp/zsd.shp

this:

INFO: Open of `shp/zs.shp'
      using driver `ESRI Shapefile' successful.

Layer name: zs
Geometry: Polygon
Feature Count: 133542
Extent: (509359.424487, 4499525.694097) - (609779.591119, 4599917.793736)
Layer SRS WKT:
PROJCRS["WGS 84 / UTM zone 33N",
    BASEGEOGCRS["WGS 84",
        DATUM["World Geodetic System 1984",
            ELLIPSOID["WGS 84",6378137,298.257223563,
                LENGTHUNIT["metre",1]]],
        PRIMEM["Greenwich",0,
            ANGLEUNIT["degree",0.0174532925199433]],
        ID["EPSG",4326]],
    CONVERSION["UTM zone 33N",
        METHOD["Transverse Mercator",
            ID["EPSG",9807]],
        PARAMETER["Latitude of natural origin",0,
            ANGLEUNIT["degree",0.0174532925199433],
            ID["EPSG",8801]],
        PARAMETER["Longitude of natural origin",15,
            ANGLEUNIT["degree",0.0174532925199433],
            ID["EPSG",8802]],
        PARAMETER["Scale factor at natural origin",0.9996,
            SCALEUNIT["unity",1],
            ID["EPSG",8805]],
        PARAMETER["False easting",500000,
            LENGTHUNIT["metre",1],
            ID["EPSG",8806]],
        PARAMETER["False northing",0,
            LENGTHUNIT["metre",1],
            ID["EPSG",8807]]],
    CS[Cartesian,2],
        AXIS["easting",east,
            ORDER[1],
            LENGTHUNIT["metre",1]],
        AXIS["northing",north,
            ORDER[2],
            LENGTHUNIT["metre",1]],
    ID["EPSG",32633]]
Data axis to CRS axis mapping: 1,2
ID: Real (32.10)

[automataIA]

Solved. Low quality data, with segments so close together that they created zeros or infinities. Rounding them to the nearest centimeter solves the problem.