search

wdtinc / mapbox-vector-tile-java

Java Mapbox Vector Tile Library for Encoding/Decoding
Apache License 2.0
147 stars 73 forks source link

Polygons crossing the Antimeridian #24

Closed grillorafael closed 6 years ago

grillorafael commented 6 years ago

Hi all,

I noticed I'm having issues displaying polygons crossing Antimeridian. They usually either display totally broken or don't display at all. Is there any specific way to deal with those polygons?

ShibaBandit commented 6 years ago

Hey @grillorafael , can you fill me in on some details for your issue? Are you writing/building MVTs or reading them? What MVT viewer are you using? What data set? The more information I can get the better. Thanks.

grillorafael commented 6 years ago

We are building MVTs from a postgis db and displaying with mapbox-gl. I'll try to get some code samples here

grillorafael commented 6 years ago

So we have a Postgis database and the elements are of type Geometry.

Here is the code to build the tiles:

    val layerParams = new MvtLayerParams()

    val tileBuilder = VectorTile.Tile.newBuilder()
    val decomposedChoices = DecomposedDBChoices.choiceMap(choice)
    val layerBuilder = decomposedChoices
      .map(name => name -> MvtLayerBuild.newLayerBuilder(name.toString, layerParams))
      .toMap
    val layerProps = new MvtLayerProps()

    val userDataConverter = new UserDataKeyValueMapConverter()

    val inTile = new MVTTile(x, y, z)
    val offsetTile = new MVTTile(x + 1, y + 1, z)
    val mercatorInTile = inTile.toMercator
    val mercatorOffsetTile = offsetTile.toMercator

    val lonlatInTile = inTile.toLonLat
    val lonlatOffsetTile = offsetTile.toLonLat
    val fullEnvelope = new Envelope(mercatorInTile.x, mercatorOffsetTile.x, mercatorInTile.y, mercatorOffsetTile.y)
    val clipEnvelope = new Envelope(fullEnvelope)
    val Xpixel = (mercatorInTile.x - mercatorOffsetTile.x).abs
    val Ypixel = (mercatorInTile.y - mercatorOffsetTile.y).abs
    clipEnvelope.expandBy(Xpixel * .01d, Ypixel * .01d)

    val tileBuildFuture = for {
      daoResult <- tiles.getter(choice, lonlatInTile.toCoordinate, lonlatOffsetTile.toCoordinate)
    } yield {
      daoResult.map { result =>
        val tmpGeom =
          List(GeomAdapter.VividSolutiontoLocationTechGeometry(result.geometry))
        tmpGeom.foreach { o =>
          o.setUserData(
            result.toJtsUserData
          )
        }
        val tiles =
          JtsAdapter.createTileGeom(tmpGeom.asJava,
                                    fullEnvelope,
                                    clipEnvelope,
                                    factory,
                                    layerParams,
                                    acceptAllGeomFilter)
        val mvtFeatures =
          tiles.mvtGeoms.asScala.toList.flatMap(g => JtsAdapter.toFeatures(g, layerProps, userDataConverter).asScala)
        mvtFeatures.foreach(y => layerBuilder.filter(f => f._1.toString == result.entityType).head._2.addFeatures(y))
      }
    }
    tileBuildFuture map { _ =>
      layerBuilder.foreach { lb =>
        MvtLayerBuild.writeProps(lb._2, layerProps)
      }
      layerBuilder.foreach { lb =>
        tileBuilder.addLayers(lb._2)
      }
      val protoBuf = tileBuilder.build()
      protoBuf.toByteArray
    }
ShibaBandit commented 6 years ago

Trying my best to follow along with the non-Java-lang code... 

daoResult <- tiles.getter(choice, lonlatInTile.toCoordinate, lonlatOffsetTile.toCoordinate)

What's daoResult? Are you wanting a EPSG:3857 vector tiles? If so, is daoResult geometry in EPSG:3857?

Side question... is this a server side app sitting between postgis and and a mapbox client? Have you tried the the PostGIS vector tile functionality?

grillorafael commented 6 years ago

Hi @ShibaBandit, This returns coordinates in EPSG:3857 using st_transform function from Postgis.

About using Postgis vector tile, the DB is hosted in RDS and they are not supporting yet (There is a runtime error ERROR: Missing libprotobuf-c)

ShibaBandit commented 6 years ago

Gotcha. Any chance you can post the WKT of one of the problematic geometries and the tile coordinates (x, y, z) so I can build a test case around this? Otherwise I can look into manufacturing a test case and make up some geometry.

grillorafael commented 6 years ago

Yes. I'll send some tomorrow.

grillorafael commented 6 years ago

One of the polygons:

01030000000100000005000000011AF4FEFF7F6640034841644FE329C0070000000080664081126A3D24C32BC039313A79D16F66403F35F9DAF1A72AC0F8947904F7706640E0EE2B8328782AC0011AF4FEFF7F6640034841644FE329C0

As text

POLYGON((179.999999501 -12.943965085,180 -13.8811358634932,179.494320501 -13.3280170849999,179.530153501 -13.234684085,179.999999501 -12.943965085))

EDIT: IGNORE THIS POLYGON

grillorafael commented 6 years ago

Another option

POLYGON((180 -13.8811354749999,179.494320501 -13.3280170849999,179.530153501 -13.234684085,179.044444501 -12.6094450839999,178.7344445 -12.2069450839999,178.6688885 -12.120834084,178.5633335 -11.982778084,178.3080555 -11.6458340829999,178.2463885 -11.5638890829999,177.316944499 -10.3050000819999,177.291666499 -10.2680560819999,177.186110499 -10.1147230819999,177.100833499 -10.088056082,176.229999498 -9.83111108199995,176.086666498 -9.82666708199998,175.944444498 -9.82194508199995,175.864444498 -9.81916708199999,174.949999497 -9.78305608199997,174.718055497 -9.96305608199999,174.56194321 -10.1076860899999,174.502777778 -10.1624999999999,174.446674176 -10.224213961,174.305555496 -10.379445082,174.188020938 -10.5335249819999,174.127777778 -10.6125,173.971111111 -10.8599999999999,173.847860227 -11.0976073739999,173.836388684 -11.11972272,173.725 -11.3902777779999,173.6375 -11.6691666669999,173.574722222 -11.9544444439999,173.525833333 -12.5358333329999,173.540277778 -12.8272222219999,173.580555556 -13.1166666669999,173.646388889 -13.401388889,173.716346921 -13.6153278239999,173.737222496 -13.6791670849999,173.759878285 -13.73181349,173.853055556 -13.948333333,173.9925 -14.206666667,174.154722222 -14.4519444439999,174.338888889 -14.6822222219999,174.453888889 -14.807777778,174.366728386 -14.8960514419999,174.236388496 -15.0280560869999,174.050277496 -15.2572230869999,173.989343302 -15.3471389769999,173.885 -15.5011111109999,173.799372701 -15.654907715,173.741944496 -15.7580560869999,173.622222496 -16.025834087,173.536943827 -16.2729919369999,173.526666667 -16.3027777779999,173.456111111 -16.586388889,173.411388889 -16.8749999999999,173.392777778 -17.166388889,173.400555556 -17.4580555559999,173.412256552 -17.557276197,173.412256555 -17.557276225,173.434722495 -17.7477780889999,173.458309235 -17.8592210039999,173.495277778 -18.033888889,173.581944444 -18.3136111109999,173.56 -18.543333333,172.764876111 -20.0224629729999,174.275709444 -25.0730201499999,174.320277778 -25.0697222219999,174.64 -25.082777778,174.959722222 -25.0697222219999,175.276944444 -25.030833333,175.554269304 -24.973293981,175.588888497 -24.966111096,175.713013011 -24.929497329,175.893055556 -24.8763888889999,176.186944444 -24.7622222219999,176.468333333 -24.6244444439999,176.734722222 -24.464444444,176.984444444 -24.283611111,177.032310278 -24.241915933,177.214999499 -24.082778095,177.309949621 -23.9747088089999,177.309949646 -23.9747087809999,177.613055556 -23.6297222219999,177.741602975 -23.434946825,177.777499499 -23.3805560939999,177.837403156 -23.2682669709999,177.837403183 -23.2682669199999,177.916944444 -23.1191666669999,177.972711347 -22.986330176,178.0311105 -22.8472230939999,178.0627775 -22.723056094,178.062777601 -22.723056248,178.062777778 -22.723055556,178.191944444 -22.918888889,178.378888889 -23.153333333,178.588055556 -23.371944444,178.817777778 -23.5722222219999,179.066111111 -23.7533333329999,179.331111111 -23.913055556,179.610833333 -24.0505555559999,180.000000735694 -24.1924997821196,180 -13.8811358634932,180 -13.8811354749999))
ShibaBandit commented 6 years ago

I've started a test case with your input geometries. What is the output you are seeing that's invalid on your side for those? Also please give me the zooms you are using to cut tiles. This will help me reproduce and look at the problem you are experiencing.

ShibaBandit commented 6 years ago

Also making my own test case that exaggerates the antimeridian crossing a bit more. What I'm looking for is the inputs to the library and whether it's handling them as advertised. How that st_transform method is handling values like '180.000000735694' in the second example probably matters. I'm not seeing where the first example WKT crosses.

ShibaBandit commented 6 years ago

Here's a WKT that crosses the antimerdian...

POLYGON((
  19037508.00 -1450000.00,
  26000000.00 -1450000.00,
  26000000.00 -1600000.00,
  19037508.00 -1600000.14,
  19037508.00 -1450000.00
  ))
        final WKTReader wktReader = new WKTReader(GEOMETRY_FACTORY);
        final List<Geometry> geom = Collections.singletonList(
                readWkt(new File("src/test/resources/wkt/antimeridian_03.wkt"), wktReader));

        final Envelope tileEnvelope = SphericalMercator.tileXyzToBB(1, 1, 1, DEFAULT_MVT_PARAMS.tileSize);
        final double tileWidth = DEFAULT_MVT_PARAMS.tileSize;
        final double tileHeight = DEFAULT_MVT_PARAMS.tileSize;

        final Envelope clipEnvelope = new Envelope(tileEnvelope);
        final double bufferWidth = tileWidth * .1f;
        final double bufferHeight = tileHeight * .1f;
        clipEnvelope.expandBy(bufferWidth, bufferHeight);

        // Use WKT from these for debug...
        final Geometry tileGeom = GEOMETRY_FACTORY.toGeometry(tileEnvelope);
        final Geometry clipGeom = GEOMETRY_FACTORY.toGeometry(clipEnvelope);

        // Build MVT tile geometry
        final TileGeomResult tileGeomResult = JtsAdapter.createTileGeom(
                geom, tileEnvelope, clipEnvelope, GEOMETRY_FACTORY,
                DEFAULT_MVT_PARAMS, ACCEPT_ALL_FILTER);

JTS Adapter Tile Geometry Output: Brown is the tile boundary, light yellow is the clip boundary, red is the intersection geometry, dark yellow is the leftover clipped geometry of the original.

tile_geom

Don't forget this is all using JTS under the hood... so if you're passing in wrapped coordinates it might not make sense... normalize input geometry and make sure it's valid.

grillorafael commented 6 years ago

Ignore the first polygon I sent you. Sent the wrong one.

[error] c.w.m.a.j.JtsAdapter - found non-noded intersection between LINESTRING ( 1.9353233617348645E7 -1568289.6320478297, 1.9368756501848653E7 -1597937.5938810317 ) and LINESTRING ( 1.7507772912075E7 -1637362.0024390216, 2.0037508342789244E7 -1560583.0925386886 ) [ (1.935996534521406E7, -1581146.9080795648, NaN) ]
org.locationtech.jts.geom.TopologyException: found non-noded intersection between LINESTRING ( 1.9353233617348645E7 -1568289.6320478297, 1.9368756501848653E7 -1597937.5938810317 ) and LINESTRING ( 1.7507772912075E7 -1637362.0024390216, 2.0037508342789244E7 -1560583.0925386886 ) [ (1.935996534521406E7, -1581146.9080795648, NaN) ]
        at org.locationtech.jts.noding.FastNodingValidator.checkValid(FastNodingValidator.java:127)
        at org.locationtech.jts.geomgraph.EdgeNodingValidator.checkValid(EdgeNodingValidator.java:78)
        at org.locationtech.jts.geomgraph.EdgeNodingValidator.checkValid(EdgeNodingValidator.java:43)
        at org.locationtech.jts.operation.overlay.OverlayOp.computeOverlay(OverlayOp.java:231)
        at org.locationtech.jts.operation.overlay.OverlayOp.getResultGeometry(OverlayOp.java:183)
        at org.locationtech.jts.operation.overlay.OverlayOp.overlayOp(OverlayOp.java:86)
        at org.locationtech.jts.operation.overlay.snap.SnapIfNeededOverlayOp.getResultGeometry(SnapIfNeededOverlayOp.java:75)
        at org.locationtech.jts.operation.overlay.snap.SnapIfNeededOverlayOp.overlayOp(SnapIfNeededOverlayOp.java:37)
        at org.locationtech.jts.geom.Geometry.intersection(Geometry.java:1353)
        at com.wdtinc.mapbox_vector_tile.adapt.jts.JtsAdapter.flatIntersection(JtsAdapter.java:162)

So I updated to only filter by this polygon and this is the log message I get.

Mercator polygon:

POLYGON ((20037508.342789244 -1560583.0479908173, 19981216.358455975 -1497232.44588009, 19985205.269769568 -1486557.137581056, 19931136.39121586 -1415147.017975345, 19896627.348958626 -1369269.1870654249, 19889329.688420184 -1359463.1980369468, 19877579.3595695 -1343748.5061339263, 19849161.94259876 -1305428.2478044252, 19842297.20356003 -1296115.808114179, 19738831.970647845 -1153382.553427536, 19736018.036559567 -1149202.7848828451, 19724267.596389394 -1131860.2003064523, 19714774.60417301 -1128844.9030731008, 19617833.80661624 -1099804.248379368, 19601878.050042357 -1099302.1751065387, 19586045.969422758 -1098768.7014008262, 19577140.410159294 -1098454.8570211122, 19475344.858289506 -1094375.45834727, 19449524.970316954 -1114714.1906592846, 19432146.63002155 -1131064.4876201905, 19425560.364258748 -1137263.0696633745, 19419314.93985244 -1144243.2103423364, 19403605.68025342 -1161806.5989473464, 19390521.793106247 -1179248.3552335203, 19383815.555211276 -1188191.6551411022, 19366375.501616538 -1216234.1319285487, 19352655.275969848 -1243177.611701249, 19351378.269644488 -1245686.488245749, 19338978.53806145 -1276395.1510278175, 19329238.082617052 -1308080.3236675279, 19322249.692336965 -1340524.3305145602, 19316807.406108018 -1406751.2716207672, 19318415.35437022 -1439999.6548067757, 19322899.05610746 -1473064.3076378794, 19330227.589214247 -1505627.5096007837, 19338015.281713378 -1530120.6809864347, 19340339.140092406 -1537433.7057614166, 19342861.170987397 -1543466.039109835, 19353233.617348645 -1568289.6320478297, 19368756.501848653 -1597937.5938810317, 19386814.996997055 -1626118.5446400617, 19407316.33658859 -1652604.387740755, 19420118.078029808 -1667057.1268723146, 19410415.415218584 -1677223.3450284004, 19395906.04503374 -1692433.7153839474, 19375188.263282683 -1718862.1226816499, 19368405.099834714 -1729239.417305958, 19356789.656588387 -1747019.9140019412, 19347257.66926569 -1764793.375066282, 19340864.79072793 -1776721.1715876255, 19327537.39865118 -1807714.6918119488, 19318044.22064257 -1836358.4977597166, 19316900.17242457 -1839812.8856359886, 19309045.963858012 -1872730.913210239, 19304067.508877818 -1906279.1451314767, 19301995.729478214 -1940202.314996254, 19302861.547764678 -1974211.2362327129, 19304164.09668116 -1985792.9313489683, 19304164.09701513 -1985792.9346182102, 19306664.99401612 -2008047.4120350352, 19309290.657902393 -2021077.168611503, 19313405.977284513 -2041515.560979306, 19323053.666412394 -2074289.0371811301, 19320610.82208057 -2101243.8896629307, 19232098.035639517 -2275692.164451594, 19400283.232942592 -2884716.1730047823, 19405244.557188973 -2884310.861113917, 19440835.872137308 -2885915.4391782135, 19476427.18708564 -2884310.861113917, 19511740.20330699 -2879532.280017855, 19542611.865506507 -2872464.747964094, 19546465.656442944 -2871582.709753931, 19560283.13413637 -2867087.4410747406, 19580325.3785669 -2860569.404813194, 19613040.93992888 -2846567.120952145, 19644365.007767234 -2829686.109797525, 19674019.283243705 -2810105.6823032047, 19701818.233836494 -2788005.637248687, 19707146.634103782 -2782914.449160216, 19727483.505158935 -2763498.2296491503, 19738053.30439072 -2750326.57034395, 19738053.30717373 -2750326.5669326973, 19771794.90273133 -2708352.6253326936, 19786104.73595721 -2684703.5337785725, 19790100.718730137 -2678105.797117097, 19796769.16332402 -2664493.3901419137, 19796769.16632966 -2664493.3839619756, 19805623.659001224 -2646436.205707584, 19811831.60224631 -2630365.587299248, 19818332.56622103 -2613553.281291495, 19821857.720535982 -2598561.137940095, 19821857.73177925 -2598561.1565258726, 19821857.751482803 -2598561.0730105564, 19836236.518969394 -2622212.563499624, 19857047.079393417 -2650572.2727424167, 19880331.406254783 -2677061.430111272, 19905903.96703171 -2701367.4890559665, 19933548.307208277 -2723379.379467087, 19963047.972268496 -2742817.123656242, 19994186.5075851 -2759569.7351219007, -20037508.26089217 -2776882.651658859, 20037508.342789244 -1560583.0925386886, 20037508.342789244 -1560583.0479908173))
grillorafael commented 6 years ago

It seems to break on all zoom levels for this guy

ShibaBandit commented 6 years ago

The mercator polygon you posted is not valid. JTS clipping won't work with invalid polygons. This is because your output is wrapping the coordinates like I said above.

QGIS 2.18 QuickWKT:

screen shot 2018-03-28 at 10 00 08 am screen shot 2018-03-28 at 10 00 05 am

Check out a similar issue here on the mapbox github: https://github.com/mapbox/mapbox-gl-js/issues/3250

Blog dealing with this problem: https://macwright.org/2016/09/26/the-180th-meridian.html

PostGIS: https://gis.stackexchange.com/questions/29975/postgis-incorrect-interpretation-of-a-polygon-that-intersects-the-180th-meridia

grillorafael commented 6 years ago

I'll have a look. I'll close the issue since it doesn't seem related to this library.

Thanks for helping me with this!

ShibaBandit commented 6 years ago

Let me know if there's anything else I can help with. I might look into adding a README section for this, was a good question.