0

Our team is creating a module for performing offline downloads of Leaflet tile layers that utilizes GeoJSON to optimize the process.

The typical approach is to only use a bounding box and to generate tile URLs for every tile within that box. This approach isn't adequate for what we're doing because we are working with agriculture fields which may be non-contiguous and spread far apart physically. Therefore the overwhelming majority of the time is often wasted iterating through and calling tiles that don't exist. We also need to iterate through several zoom levels in order to obtain a fully functional offline layer, so this can easily involve iterating through 100,000+ tiles. This is unacceptable for our case since we want to support mobile devices (mbtiles is an alternative option, but we'd like to avoid it if possible).

We are optimizing this process by accepting GeoJSON data that describes the true shape of our agriculture fields. We are still iterating through every tile in the bounding box of the shape, but we ultimately only call the URLs of those tiles which intersect the bounding box and the GeoJSON shape.

The problem stems from needing to project the coordinates of the shape's bounding box into XYZ, allowing us to iterate within the boundary by tile indices (if there's a way to avoid this please enlighten me). We then need to convert the XYZ point back into lat/lng in order to determine if the tile's coordinates intersect with the GeoJSON shape. When converting back to lat/lng there is a pretty severe loss of accuracy / drift occurring, and this causes the intersection detector to always return false.

This snippet outlines the issue that occurs during project/unproject

for (var j = tileBounds.min.y; j <= tileBounds.max.y; j += 1) {
  for (var i = tileBounds.min.x; i <= tileBounds.max.x; i += 1) {
    var tilePoint = new L.Point(i, j)
    var tileLatLng = map.layerPointToLatLng(tilePoint) // FIXME: this value is experiencing drift!
    var tileShape = { type: 'Point', coordinates: [tileLatLng.lng, tileLatLng.lat] }
    var tileIntersects = gju.pointInPolygon(tileShape, origShape)

    if (tileIntersects) {
      // ...
    }
  }
}

Interestingly the use of map.unproject (instead of map.layerPointToLatLng) produces much worse results. layerPointToLatLng produces data that's "close", but still several KM away from the expected coordinate.

For instance, if the lat/lng bounding box of the shape is

[-120.210483180564, 36.1665531949778, -120.192638904505, 36.1764258657549]

Then calling unproject on the tile coordinates will produce something like

{lat: 85.00333833822032, lng: -179.7664546966553}

which is completely inaccurate (or at least, completely unexpected).

However, using layerPointToLatLng gets us closer to an expected coordinate, but still off:

 {lat: 35.732091125624216, lng: -120.022931098938}

layerPointToLatLng is simply calling unproject, but it also subtracts the pixel origin from the result.

I have tried experimenting with different Coordinate Reference Systems but they make no difference in the results. We are using Leaflet's default CRS, which is EPSG3857.

What else could be causing this drift? If this is somehow expected behavior, which I'm assuming it could be, what is the reasoning behind it?

I put together a Codepen which outlines the problem. You will need to look at the console output in order to get the full picture.

If you zoom out 6 levels or so you will see a marker cluster to the south and slightly east of the expected area. This helps to visualize the high amount of drift taking place.

https://codepen.io/slurmulon/pen/zRVzjX

1

This is an instance of a XY problem. You don't want to fix bounding box projections in Leaflet, what you really want to do is literally...

generate tile URLs for every tile within [...] agriculture fields which may be non-contiguous and spread far apart physically.

This is a textbook use case for the supermercado library. That'll allow you to get a complete list of tile coordinates from geometries.

There are other approaches, like creating vector tiles from your polygon data and iterating through non-empty vector tiles. Or leveraging the spherical-mercator library in javascript.

  • Thanks for sharing that library, it's just what I needed. I definitely don't want to fix bounding projections in Leaflet, I more want to understand why the drift even occurs. I viewed the projection drift as a roadblock to the overall problem. Anyways, I'll be more aware of the XY problem from now on. – Erik Vavro Mar 8 '18 at 19:13

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.