mapbox / mbtiles-extracts

Extract parts of MBTiles into separate files using a GeoJSON with polygons.
ISC License
52 stars 13 forks source link

mbtiles-extracts return nothing #10

Open spatialhast opened 6 years ago

spatialhast commented 6 years ago

For extract parts of an mbtiles I use command: mbtiles-extracts adria.mbtiles polygon3857.geojson admin but it return nothing.

gdalinfo input mbtiles file:

$ gdalinfo adria.mbtiles
Driver: MBTiles/MBTiles
Files: adria.mbtiles
Size is 1218560, 915456
Coordinate System is:
PROJCS["WGS 84 / Pseudo-Mercator",
    GEOGCS["WGS 84",
        DATUM["WGS_1984",
            SPHEROID["WGS 84",6378137,298.257223563,
                AUTHORITY["EPSG","7030"]],
            AUTHORITY["EPSG","6326"]],
        PRIMEM["Greenwich",0,
            AUTHORITY["EPSG","8901"]],
        UNIT["degree",0.0174532925199433,
            AUTHORITY["EPSG","9122"]],
        AUTHORITY["EPSG","4326"]],
    PROJECTION["Mercator_1SP"],
    PARAMETER["central_meridian",0],
    PARAMETER["scale_factor",1],
    PARAMETER["false_easting",0],
    PARAMETER["false_northing",0],
    UNIT["metre",1,
        AUTHORITY["EPSG","9001"]],
    AXIS["X",EAST],
    AXIS["Y",NORTH],
    EXTENSION["PROJ4","+proj=merc +a=6378137 +b=6378137 +lat_ts=0.0 +lon_0=0.0 +x_0=0.0 +y_0=0 +k=1.0 +units=m +nadgrids=@null +wktext  +no_defs"],
    AUTHORITY["EPSG","3857"]]
Origin = (1120261.086547542363405,5743172.557235002517700)
Pixel Size = (2.388657133911758,-2.388657133911758)
Metadata:
  ZOOM_LEVEL=16
Image Structure Metadata:
  INTERLEAVE=PIXEL
Corner Coordinates:
Upper Left  ( 1120261.087, 5743172.557) ( 10d 3'48.52"E, 45d46' 3.08"N)
Lower Left  ( 1120261.087, 3556462.052) ( 10d 3'48.52"E, 30d24'38.81"N)
Upper Right ( 4030983.124, 5743172.557) ( 36d12'39.37"E, 45d46' 3.08"N)
Lower Right ( 4030983.124, 3556462.052) ( 36d12'39.37"E, 30d24'38.81"N)
Center      ( 2575622.105, 4649817.305) ( 23d 8'13.95"E, 38d29'47.74"N)
Band 1 Block=256x256 Type=Byte, ColorInterp=Red
  Overviews: 609280x457728, 304640x228864, 152320x114432, 76160x57216, 38080x28608, 19040x14304, 9520x7152, 4760x3576, 2380x1788, 1190x894, 595x447, 298x223, 149x112
  Mask Flags: PER_DATASET ALPHA
  Overviews of mask band: 609280x457728, 304640x228864, 152320x114432, 76160x57216, 38080x28608, 19040x14304, 9520x7152, 4760x3576, 2380x1788, 1190x894, 595x447, 298x223, 149x112
Band 2 Block=256x256 Type=Byte, ColorInterp=Green
  Overviews: 609280x457728, 304640x228864, 152320x114432, 76160x57216, 38080x28608, 19040x14304, 9520x7152, 4760x3576, 2380x1788, 1190x894, 595x447, 298x223, 149x112
  Mask Flags: PER_DATASET ALPHA
  Overviews of mask band: 609280x457728, 304640x228864, 152320x114432, 76160x57216, 38080x28608, 19040x14304, 9520x7152, 4760x3576, 2380x1788, 1190x894, 595x447, 298x223, 149x112
Band 3 Block=256x256 Type=Byte, ColorInterp=Blue
  Overviews: 609280x457728, 304640x228864, 152320x114432, 76160x57216, 38080x28608, 19040x14304, 9520x7152, 4760x3576, 2380x1788, 1190x894, 595x447, 298x223, 149x112
  Mask Flags: PER_DATASET ALPHA
  Overviews of mask band: 609280x457728, 304640x228864, 152320x114432, 76160x57216, 38080x28608, 19040x14304, 9520x7152, 4760x3576, 2380x1788, 1190x894, 595x447, 298x223, 149x112
Band 4 Block=256x256 Type=Byte, ColorInterp=Alpha
  Overviews: 609280x457728, 304640x228864, 152320x114432, 76160x57216, 38080x28608, 19040x14304, 9520x7152, 4760x3576, 2380x1788, 1190x894, 595x447, 298x223, 149x112

GeoJSON file in EPSG:3857 coordinate system. Source files: mbtiles - https://www.dropbox.com/s/amngyovs64292e2/adria.mbtiles?dl=0, geojson - https://www.dropbox.com/s/fvj6rwaixzy4ah6/polygon3857.geojson?dl=0

There is any spesified requirements for input mbtiles and geojson file?

arunasank commented 5 years ago

Seeing this issue for:

cc/ @mourner @geohacker @batpad

arunasank commented 5 years ago

I did some more digging today. So, the mbtiles-extracts library generates

In several cases, the unprojected [lng, lat] may not be within the polygons in the geojson, even though the polygon intersects part of the tile. Here are a few polygons (all from the geojson referenced above) where the unprojected [lng, lat] (from the tiles in the mbtiles tileset referenced above) is outside the polygon but within the bounding box constricting the polygon:

screen shot 2018-11-09 at 6 43 30 pm

Full set of polygons and unprojected [lng, lat] points that are excluded as a result of a bug somewhere in the data/code here

Why is this the case? Is there a better way to handle this? Is this a bug in the mbtiles-extracts library or the mbtiles file or the geojson file?

arunasank commented 5 years ago

cc also @answerquest, who helped me validate my geojson to make sure that it is following the right hand verification rule, a fix that taught me a lot about geojsons and GIS trivia, but wasn't the solution for the above problem. : )

mourner commented 5 years ago

One way you could try to fix it is replacing single-point whichPolygon query with a bbox one, like described in the readme at the bottom (query.bbox).

arunasank commented 5 years ago

One way you could try to fix it is replacing single-point whichPolygon query with a bbox one, like described in the readme at the bottom (query.bbox).

👋 @mourner Yes I was considering this. So, if I do this, I will be checking the intersection of:

Correct?

mourner commented 5 years ago

@arunasank this line currently queries the center of the tile. Instead, you would unproject 2 corners of the tile ([x, y], [x + 1, y + 1]), make a bounding box out of this ([minLng, minLat, maxLng, maxLat]), and use that as a bbox query (to compare against polygons in the index).

arunasank commented 5 years ago

@mourner got it, that clarifies things. I'll submit a PR. Thanks for the map math! ♥️

arunasank commented 5 years ago

I opened a new issue for the problem I was seeing at https://github.com/mapbox/mbtiles-extracts/issues/12, which seems to be different from the OP's issue.

The solution discussed above does not resolve the OP's issue.