Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Strange polygon buffer behavior on RGeo 2.1.1 #252

Open
durkie opened this issue Apr 1, 2021 · 6 comments
Open

Strange polygon buffer behavior on RGeo 2.1.1 #252

durkie opened this issue Apr 1, 2021 · 6 comments

Comments

@durkie
Copy link

durkie commented Apr 1, 2021

Steps to reproduce

I noticed that a worker of mine was failing after calculating the intersection of a line with a polygon: line.intersection(polygon). The failure would occur with No method "is_empty?" for NilClass, which is strange because even if the polygon and line don't intersect, I should at least always get a geometry response and not nil.

I discovered that the slight amount of buffer that I apply to this particular geometry seems to cause it to generate a very small, zero-area polygon, and intersecting against that returns nil for some reason.

Input polygon:

POLYGON ((115.24384559999999 -8.610341399999996, 115.2439816 -8.610488500000002, 115.2439844 -8.610691900000006, 115.2439875 -8.610914699999995, 115.24411409999999 -8.610929200000001, 115.2441567 -8.611075400000004, 115.2441873 -8.611204700000002, 115.2442387 -8.611217799999991, 115.244305 -8.611229099999989, 115.2443504 -8.611279600000003, 115.244374 -8.611389800000012, 115.2443738 -8.61142559999999, 115.2443951 -8.611495099999999, 115.2444285 -8.611555100000018, 115.24445530000001 -8.611592000000002, 115.24447 -8.611633299999994, 115.24447809999998 -8.611745999999997, 115.24449550000001 -8.611803699999996, 115.2445075 -8.611815000000021, 115.24461529999999 -8.611874100000009, 115.244641 -8.611914599999992, 115.2446801 -8.611980900000006, 115.2447008 -8.612018899999995, 115.2447142 -8.612067999999994, 115.24472249999998 -8.612133, 115.2447272 -8.61218070000001, 115.24474910000002 -8.612236100000004, 115.2447918 -8.612296200000003, 115.24484310000001 -8.6123956, 115.2448663 -8.612417800000003, 115.2449044 -8.6124473, 115.2449313 -8.6124552, 115.2449645 -8.612471200000002, 115.2450308 -8.612515200000004, 115.2450343 -8.612561799999995, 115.24503779999999 -8.612597300000004, 115.24503520000002 -8.612641400000001, 115.2450349 -8.6127015, 115.2450457 -8.612756700000006, 115.245054 -8.612813000000003, 115.245067 -8.612848799999995, 115.24514119999999 -8.612858700000004, 115.2451868 -8.612856500000007, 115.24523599999999 -8.612857899999995, 115.2453135 -8.612859500000013, 115.24537129999999 -8.612862199999995, 115.24538849999999 -8.61287830000002, 115.2454139 -8.612976300000014, 115.2454481 -8.613038899999978, 115.2455229 -8.613085799999993, 115.24554249999998 -8.613105500000003, 115.245578 -8.613148499999994, 115.2455865 -8.613174299999997, 115.2455961 -8.613220900000002, 115.2455994 -8.613306600000001, 115.24561600000001 -8.613439, 115.245723 -8.6135728, 115.2457284 -8.613631800000007, 115.2457279 -8.613750600000003, 115.24573370000002 -8.613815500000015, 115.24574069999998 -8.614159799999996, 115.2457399 -8.614331199999995, 115.2457778 -8.614652300000003, 115.24582659999999 -8.6147651, 115.2458681 -8.614830200000014, 115.2458995 -8.614954100000006, 115.24591619999998 -8.615054600000008, 115.2459413 -8.615218800000008, 115.24599199999999 -8.615449400000017, 115.2460542 -8.615570899999994, 115.2460761 -8.615618799999993, 115.2461383 -8.615742699999998, 115.24617610000001 -8.615820100000008, 115.24619900000002 -8.615913199999994, 115.2462171 -8.616008800000003, 115.24623490000002 -8.616137600000002, 115.24625799999998 -8.616208700000001, 115.2463035 -8.616239000000022, 115.2463205 -8.616262899999995, 115.24645009999999 -8.616433700000002, 115.2465015 -8.616512400000005, 115.246587 -8.61663519999999, 115.2466743 -8.616665000000012, 115.24677410000001 -8.616645899999995, 115.2469036 -8.616597500000012, 115.24697370000001 -8.616599000000008, 115.2470744 -8.616646000000003, 115.24711359999999 -8.616675600000008, 115.24717999999999 -8.616685700000005, 115.2472514 -8.616686000000001, 115.24730059999999 -8.616686200000004, 115.24740890000001 -8.616693999999995, 115.24743589999999 -8.616696599999997, 115.2475159 -8.616709299999997, 115.2476044 -8.616726799999995, 115.2476856 -8.616733299999993, 115.2477914 -8.6167509, 115.2479363 -8.616804299999998, 115.24802469999999 -8.616843799999998, 115.2480862 -8.616873499999983, 115.2481536 -8.61691420000001, 115.2482386 -8.616915899999995, 115.2483273 -8.616880699999982, 115.2483964 -8.616839400000003, 115.24847650000001 -8.6168251, 115.24852360000001 -8.616841100000002, 115.24854510000002 -8.616874699999983, 115.24859049999999 -8.6169395, 115.2486335 -8.616985200000002, 115.2486648 -8.616990099999995, 115.24870319999998 -8.616997400000002, 115.2487633 -8.616985799999995, 115.2488114 -8.616976399999984, 115.24886199999999 -8.616955099999998, 115.2489126 -8.616931399999999, 115.24900360000001 -8.617003600000004, 115.2490347 -8.617025299999995, 115.2491051 -8.617141399999994, 115.24914689999999 -8.617207500000006, 115.2492185 -8.617313199999998, 115.2492926 -8.617401999999998, 115.24938130000001 -8.617450199999993, 115.2495254 -8.617484399999995, 115.24965010000001 -8.617523199999994, 115.24978450000002 -8.6175717, 115.24990699999998 -8.617615299999997, 115.25004360000001 -8.617692399999996, 115.2501995 -8.617762600000006, 115.2503002 -8.617810800000015, 115.250348 -8.617875699999999, 115.2504291 -8.618012399999984, 115.25045750000001 -8.618108199999995, 115.2504907 -8.618208800000005, 115.25053599999998 -8.618287899999999, 115.25061239999998 -8.618391099999997, 115.2506386 -8.618441500000003, 115.25068389999998 -8.618532700000003, 115.25074580000002 -8.618654899999981, 115.2507765 -8.618767500000004, 115.25081940000001 -8.618865799999995, 115.2508575 -8.618942500000003, 115.2508735 -8.619207099999997, 115.25089669999998 -8.619303900000006, 115.25092979999998 -8.619416500000014, 115.2509681 -8.619440600000019, 115.2510692 -8.619509400000013, 115.2511286 -8.61956810000001, 115.25124139999998 -8.61959730000001, 115.2513518 -8.619647999999998, 115.25145270000002 -8.619665299999994, 115.2515535 -8.619689600000001, 115.2516181 -8.619756899999999, 115.25168259999998 -8.619845699999999, 115.2517278 -8.619946400000003, 115.25176580000002 -8.620039899999995, 115.25178775925927 -8.620065100000005, 115.24384559999999 -8.620065100000005, 115.24384559999999 -8.610341399999996))

polygon.buffer(1e-9) returns:

POLYGON ((115.25178776001319 -8.620065099343037, 115.25178776001319 -8.620065099343039, 115.25178776001319 -8.620065099343035, 115.25178776001319 -8.620065099343037))

and polygon.buffer(1e-9).area is 0.0. line.intersection(polygon.buffer(1e-9)) returns nil

whereas polygon.buffer(1e-8) and polygon.buffer(1e-10) both return the expected polygons, have non-zero areas, and line.intersection(polygon.buffer(1e-8)) or line.intersection(polygon.buffer(1e-10)) both return LINESTRING EMPTY

Ruby version: 2.6.5

OS: Linux

@keithdoggett
Copy link
Member

Are you using a GEOS factory? It's hard to tell exactly what's going on without the intersecting line, but I'm guessing that GEOSIntersection(line, polygon) is returning NULL to indicate an exception which we then return as nil. It's probably raising an exception because the polygon collapses after the buffer operation, but you could try to increase the buffer_resolution of the factory which may help. We are likely going to change this behavior in the next version to explicitly raise an error here.

@durkie
Copy link
Author

durkie commented Apr 6, 2021

Yes, that's using the GEOS factory. The intersecting line in this particular case didn't actually even intersect with the geometry.

But it's strange that I do get different results with different non-intersecting lines:
if line is LINESTRING(0 0,1 1) then you get the expected LINESTRING EMPTY result
if line is LINESTRING(115.25158351287246 -8.625786857591237,115.26185603812337 -8.61876339617126) (some place near the polygon geometry but not intersecting its bounding box) then you get nil as a result

@keithdoggett
Copy link
Member

Thanks for sending those examples. I tried creating both of those geometries in a C program and got the following error message from GEOS while creating the second linestring (LINESTRING(115.25158351287246 -8.625786857591237,115.26185603812337 -8.61876339617126)).

ERROR: TopologyException: Input geom 1 is invalid: Self-intersection at or near point 115.25178776001319 -8.620065099343039 at 115.25178776001319 -8.620065099343039

This arose from the following statement:

const GEOSGeometry *int2 = GEOSIntersection(line2, buf);

Whereas the other intersection returned LINESTRING EMPTY as expected. So the reason that nil is returned is because the GEOS error handler in the C extension ignores errors and falls back to nil.

My opinion is that this should be handled more explicitly and raise an RGeo::TopologyError or something if we encounter this. @BuonOmo if you have anything to add let us know, but I think this is just another case where we need to be more explicit with error checking. This could also be prevented on the frontend if we don't allow ops with invalid geometries since the buffer in this case is invalid.

I can also share my C program if you think that's helpful.

@durkie
Copy link
Author

durkie commented Apr 12, 2021

Ah cool, glad you were able to get some insight about what's going on. And I would be interested to see the C program if it's no trouble (just to learn a bit more about GEOS). Thanks!

@keithdoggett
Copy link
Member

keithdoggett commented Apr 15, 2021

// main.c
#include <geos_c.h>
#include <stdlib.h>
#include <stdio.h>
#include <stdarg.h>

void notice(const char *fmt, ...)
{
    va_list ap;
    fprintf(stdout, "NOTICE: ");
    va_start(ap, fmt);
    vfprintf(stdout, fmt, ap);
    va_end(ap);
    fprintf(stdout, "\n");
}

void log_and_exit(const char *fmt, ...)
{
    va_list ap;
    fprintf(stdout, "ERROR: ");
    va_start(ap, fmt);
    vfprintf(stdout, fmt, ap);
    va_end(ap);
    fprintf(stdout, "\n");
    exit(1);
}

int main()
{
    initGEOS(notice, log_and_exit);

    GEOSGeometry *gc1 = GEOSGeomFromWKT("POLYGON ((115.24384559999999 -8.610341399999996, 115.2439816 -8.610488500000002, 115.2439844 -8.610691900000006, 115.2439875 -8.610914699999995, 115.24411409999999 -8.610929200000001, 115.2441567 -8.611075400000004, 115.2441873 -8.611204700000002, 115.2442387 -8.611217799999991, 115.244305 -8.611229099999989, 115.2443504 -8.611279600000003, 115.244374 -8.611389800000012, 115.2443738 -8.61142559999999, 115.2443951 -8.611495099999999, 115.2444285 -8.611555100000018, 115.24445530000001 -8.611592000000002, 115.24447 -8.611633299999994, 115.24447809999998 -8.611745999999997, 115.24449550000001 -8.611803699999996, 115.2445075 -8.611815000000021, 115.24461529999999 -8.611874100000009, 115.244641 -8.611914599999992, 115.2446801 -8.611980900000006, 115.2447008 -8.612018899999995, 115.2447142 -8.612067999999994, 115.24472249999998 -8.612133, 115.2447272 -8.61218070000001, 115.24474910000002 -8.612236100000004, 115.2447918 -8.612296200000003, 115.24484310000001 -8.6123956, 115.2448663 -8.612417800000003, 115.2449044 -8.6124473, 115.2449313 -8.6124552, 115.2449645 -8.612471200000002, 115.2450308 -8.612515200000004, 115.2450343 -8.612561799999995, 115.24503779999999 -8.612597300000004, 115.24503520000002 -8.612641400000001, 115.2450349 -8.6127015, 115.2450457 -8.612756700000006, 115.245054 -8.612813000000003, 115.245067 -8.612848799999995, 115.24514119999999 -8.612858700000004, 115.2451868 -8.612856500000007, 115.24523599999999 -8.612857899999995, 115.2453135 -8.612859500000013, 115.24537129999999 -8.612862199999995, 115.24538849999999 -8.61287830000002, 115.2454139 -8.612976300000014, 115.2454481 -8.613038899999978, 115.2455229 -8.613085799999993, 115.24554249999998 -8.613105500000003, 115.245578 -8.613148499999994, 115.2455865 -8.613174299999997, 115.2455961 -8.613220900000002, 115.2455994 -8.613306600000001, 115.24561600000001 -8.613439, 115.245723 -8.6135728, 115.2457284 -8.613631800000007, 115.2457279 -8.613750600000003, 115.24573370000002 -8.613815500000015, 115.24574069999998 -8.614159799999996, 115.2457399 -8.614331199999995, 115.2457778 -8.614652300000003, 115.24582659999999 -8.6147651, 115.2458681 -8.614830200000014, 115.2458995 -8.614954100000006, 115.24591619999998 -8.615054600000008, 115.2459413 -8.615218800000008, 115.24599199999999 -8.615449400000017, 115.2460542 -8.615570899999994, 115.2460761 -8.615618799999993, 115.2461383 -8.615742699999998, 115.24617610000001 -8.615820100000008, 115.24619900000002 -8.615913199999994, 115.2462171 -8.616008800000003, 115.24623490000002 -8.616137600000002, 115.24625799999998 -8.616208700000001, 115.2463035 -8.616239000000022, 115.2463205 -8.616262899999995, 115.24645009999999 -8.616433700000002, 115.2465015 -8.616512400000005, 115.246587 -8.61663519999999, 115.2466743 -8.616665000000012, 115.24677410000001 -8.616645899999995, 115.2469036 -8.616597500000012, 115.24697370000001 -8.616599000000008, 115.2470744 -8.616646000000003, 115.24711359999999 -8.616675600000008, 115.24717999999999 -8.616685700000005, 115.2472514 -8.616686000000001, 115.24730059999999 -8.616686200000004, 115.24740890000001 -8.616693999999995, 115.24743589999999 -8.616696599999997, 115.2475159 -8.616709299999997, 115.2476044 -8.616726799999995, 115.2476856 -8.616733299999993, 115.2477914 -8.6167509, 115.2479363 -8.616804299999998, 115.24802469999999 -8.616843799999998, 115.2480862 -8.616873499999983, 115.2481536 -8.61691420000001, 115.2482386 -8.616915899999995, 115.2483273 -8.616880699999982, 115.2483964 -8.616839400000003, 115.24847650000001 -8.6168251, 115.24852360000001 -8.616841100000002, 115.24854510000002 -8.616874699999983, 115.24859049999999 -8.6169395, 115.2486335 -8.616985200000002, 115.2486648 -8.616990099999995, 115.24870319999998 -8.616997400000002, 115.2487633 -8.616985799999995, 115.2488114 -8.616976399999984, 115.24886199999999 -8.616955099999998, 115.2489126 -8.616931399999999, 115.24900360000001 -8.617003600000004, 115.2490347 -8.617025299999995, 115.2491051 -8.617141399999994, 115.24914689999999 -8.617207500000006, 115.2492185 -8.617313199999998, 115.2492926 -8.617401999999998, 115.24938130000001 -8.617450199999993, 115.2495254 -8.617484399999995, 115.24965010000001 -8.617523199999994, 115.24978450000002 -8.6175717, 115.24990699999998 -8.617615299999997, 115.25004360000001 -8.617692399999996, 115.2501995 -8.617762600000006, 115.2503002 -8.617810800000015, 115.250348 -8.617875699999999, 115.2504291 -8.618012399999984, 115.25045750000001 -8.618108199999995, 115.2504907 -8.618208800000005, 115.25053599999998 -8.618287899999999, 115.25061239999998 -8.618391099999997, 115.2506386 -8.618441500000003, 115.25068389999998 -8.618532700000003, 115.25074580000002 -8.618654899999981, 115.2507765 -8.618767500000004, 115.25081940000001 -8.618865799999995, 115.2508575 -8.618942500000003, 115.2508735 -8.619207099999997, 115.25089669999998 -8.619303900000006, 115.25092979999998 -8.619416500000014, 115.2509681 -8.619440600000019, 115.2510692 -8.619509400000013, 115.2511286 -8.61956810000001, 115.25124139999998 -8.61959730000001, 115.2513518 -8.619647999999998, 115.25145270000002 -8.619665299999994, 115.2515535 -8.619689600000001, 115.2516181 -8.619756899999999, 115.25168259999998 -8.619845699999999, 115.2517278 -8.619946400000003, 115.25176580000002 -8.620039899999995, 115.25178775925927 -8.620065100000005, 115.24384559999999 -8.620065100000005, 115.24384559999999 -8.610341399999996))");

    GEOSGeometry *buf = GEOSBuffer(gc1, 1e-9, 1);
    if (buf != NULL)
    {
        printf("buf: %s\n", GEOSGeomToWKT(buf));
        printf("buf is_valid: %i\n", GEOSisValid(buf));

        GEOSGeometry *line1 = GEOSGeomFromWKT("LINESTRING(0 0, 1 1)");
        GEOSGeometry *int1 = GEOSIntersection(line1, buf);

        printf("int1: %s\n", GEOSGeomToWKT(int1));
        printf("int1 is_valid: %i\n", GEOSisValid(int1));

        GEOSGeometry *line2 = GEOSGeomFromWKT("LINESTRING(115.25158351287246 -8.625786857591237, 115.26185603812337 -8.61876339617126)");
        GEOSGeometry *int2 = GEOSIntersection(line2, buf);

        // will exit before the program reaches here
        printf("int2: %s\n", GEOSGeomToWKT(int2));
        printf("int2 is_valid: %c\n", GEOSisValid(int2));

        GEOSGeom_destroy(buf);
        GEOSGeom_destroy(line1);
        GEOSGeom_destroy(int1);
        GEOSGeom_destroy(line2);
        GEOSGeom_destroy(int2);
    }

    GEOSGeom_destroy(gc1);

    finishGEOS();

    return 0;
}

To compile: gcc -Wall -o geos_demo main.c -lgeos_c

@keithdoggett keithdoggett mentioned this issue Apr 21, 2021
12 tasks
@mjy
Copy link

mjy commented Jul 7, 2023

Leaving a note here. I recently did a brew update and started to see a number of failed tests. I haven't eliminated all potential issues with my machine, but have poked a fair bit.

Setup:

SELECT PostGIS_Version(); -> 3.3 USE_GEOS=1 USE_PROJ=1 USE_STATS=1
Proj: /usr/local/Cellar/proj/9.2.1 (481 files, 708MB) *
M1 - 8GB - Ventura 13.1
Ruby 3.2.1 (RVM)
rgeo (3.0.0)
rgeo-activerecord (7.0.1)
rgeo-geojson (2.1.1)
rgeo-proj4 (4.0.0)
(byebug) ap Gis::FACTORY.marshal_dump
{
    "pref" => "Projected",
    "hasz" => true,
    "hasm" => false,
    "srid" => 4326,
    "wktg" => {
             "tag_format" => "wkt11",
         "emit_ewkt_srid" => nil,
        "square_brackets" => false,
           "convert_case" => "upper"
    },
    "wkbg" => {
           "type_format" => "wkb11",
        "emit_ewkb_srid" => false,
            "hex_format" => true,
         "little_endian" => false
    },
    "wktp" => {
               "support_ewkt" => false,
              "support_wkt12" => false,
               "strict_wkt11" => false,
        "ignore_extra_tokens" => false,
               "default_srid" => nil
    },
    "wkbp" => {
              "support_ewkb" => true,
             "support_wkb12" => false,
        "ignore_extra_bytes" => false,
              "default_srid" => nil
    },
    "bufr" => 1,
      "cs" => "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)\"],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]]",
    "prjc" => "Projector",
    "prjf" => #<RGeo::Geos::CAPIFactory:0x1ce1c srid=4326 bufres=4 flags=10>
}
nil
(byebug) ap k.geo_object
#<RGeo::Geographic::ProjectedPolygonImpl:0x1cffc "POLYGON ((-33.0 -11.0 0.0, -33.0 -23.0 0.0, -21.0 -23.0 0.0, -21.0 -11.0 0.0, -27.0 -13.0 0.0, -33.0 -11.0 0.0))">
nil
(byebug) k.geo_object.buffer(4)
*** RGeo::Error::InvalidGeometry Exception: LinearRings must have 0 or >= 4 points
  • buffer() use on Ubuntu ( kg-config --modversion proj 8.2.1) seems to work fine, our tests are not failing there. I.e. could this be related to 9.2.1?
  • All buffer() calls in our settup are failing with the same Linear Ring warning (regardless of shape type) ... I think
  • I did see the buffer_resolution increase size suggestion (we started at 1, you can see this is 4)

If anyting jumps out I'd be greatful to hear about it, none of this has been confirmed on our Ubuntu setups.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

3 participants