[PYTHON] Get all XYZ tiles that contain a particular geometry

Abstract abstract

GIS data is roughly divided into rasters (images) and vectors (points, lines, polygons). These are often huge in size and the concept of XYZ tiles exists for efficient data delivery. I will omit a detailed explanation of the concept of XYZ tiles, but it is a mechanism that divides huge data and delivers only the necessary parts. Since it is a standard, the entire globe, except for some polar regions, is always divided into the same tiles, regardless of framework or application. Recently, I had a need for "where does this feature belong to XYZ tiles" and "I want a list of XYZ tiles covering the Japanese archipelago", so the Node.js module and Python Each pattern of modules is explained in this article. The XYZ tiles divide the earth dropped on a plane according to certain rules, so it is possible to calculate based on the parameters (Reference: [Qiita --Map tile calculation summary](https: // qiita) .com / MALORGIS / items / 1a9114dd090e5b891bf7)).

Sample Data Subject

The theme is GeoJSON, which dissolves polygons from Hokkaido. According to Geographical Survey Institute --tile coordinate confirmation page, the tiles including Hokkaido are as shown in the image below. (Zoom level 7 as an example). Visually list the tiles for verification.

All tiles that include Hokkaido at zoom level 7
[x, y, z]
[114, 45, 7]
[116, 45, 7]
[113, 46, 7]
[114, 46, 7]
[115, 46, 7]
[116, 46, 7]
[113, 47, 7]
[114, 47, 7]
[115, 47, 7]

Envs execution environment

• macOS 10.15.4
• Node.js v12.16.1
• Python 3.8.2

JavaScript Node.js module / tile-cover

GitHub - Mapbox / tile-cover License:MIT

Install Install

npm install @mapbox/tile-cover

Usage How to use

This is a sample that reads a local GeoJSON file and gets the index of all tiles that contain the feature. As mentioned above, the polygons in Hokkaido are treated as one feature by dissolve (the CRS is EPSG: 4326 because tile-cover is calculated by latitude and longitude).

tileindex.js

//Library Import
var fs = require('fs')
var cover = require('@mapbox/tile-cover');

//Load local GeoJSON file
var geojson_str = fs.readFileSync('./dissolved_hokkaido.geojson')
var geojson_obj = JSON.parse(geojson_str);

//use tile-cover
var limits = {
    min_zoom: 7,
    max_zoom: 7
};
var tile_index = []
tile_index = cover.tiles(geojson_obj.features[0].geometry, limits)

console.log(tile_index)

Execute this script on the console and the tile list will be output.

node tileindex.js

[
  [ 114, 45, 7 ],
  [ 116, 45, 7 ],
  [ 113, 46, 7 ],
  [ 114, 46, 7 ],
  [ 115, 46, 7 ],
  [ 116, 46, 7 ],
  [ 113, 47, 7 ],
  [ 114, 47, 7 ],
  [ 115, 47, 7 ]
]

All 9 photos have been listed properly.

Tips Reference

--The argument of cover.tiles () is GeoJSON geometry Object and Min_zoom and max_zoom Object. --It is recommended to enter the same value for min_zoom and max_zoom (the behavior is a little peculiar)

Python module / tiletanic GitHub - DigitalGlobe / tiletanic License:MIT

Install Install

pip install tiletanic

Usage How to use

Similar to Node.js, the theme is dissolved Hokkaido, but tiletanic uses the coordinates of Web Mercator instead of latitude and longitude, so CRS is set to EPSG: 3857 (CRS can be converted internally). However, it takes time, so I converted the data itself with QGIS in advance).

tileindex.py

import json

import shapely
import tiletanic

if __name__ == "__main__":
    zoomlevel = 7

    #geojson read
    geojson_dict = {}
    with open('./dissolved_hokkaido_3857.geojson', 'r') as f:
        geojson_dict = json.load(f)

    #The number of features is 1 due to dissolve
    feature = geojson_dict['features'][0]['geometry']
    feature_shape = shapely.geometry.shape(feature)

    #Initialization of tile scheme
    tiler = tiletanic.tileschemes.WebMercator()

    covering_tiles_itr = tiletanic.tilecover.cover_geometry(tiler, feature_shape, zoomlevel)
    covering_tiles = []
    for tile in covering_tiles_itr:
        tile_xyz = [tile[0], tile[1], tile[2]]
        covering_tiles.append(tile_xyz)

    print(covering_tiles)

Run in python

python tileindex.py

[[114, 45, 7], [113, 46, 7], [113, 47, 7], [114, 46, 7], [115, 46, 7], [114, 47, 7], [115, 47, 7], [116, 45, 7], [116, 46, 7]]

The output order is different from the tile-cover of Node.js, but the index of 9 tiles was output properly.

Tips Reference

--The arguments for tiletanic.cover_geometry are an instance of tiletanic.tileschemes, an instance of shapely.geometry, and a zoom level value. --When you install tiletanic with pip, the dependent shapely is also installed at the same time.

Sample codes Sample codes

GitHub - Kanahiro / tileindex-sample All the code, including the GeoJSON files, is in the repository above.