diff options
author | Bryan Brattlof <hello@bryanbrattlof.com> | 2020-11-12 17:37:25 -0500 |
---|---|---|
committer | Bryan Brattlof <hello@bryanbrattlof.com> | 2020-11-12 17:44:44 -0500 |
commit | 41cb4761013bb76c0ccc56d48d4512e35f5251c6 (patch) | |
tree | 9e1635e653c5eb3ae625368cf0360a42813e5ef6 | |
parent | 6241adbbf78c4ca6beb94883aaa6798f075afd92 (diff) | |
download | norta-41cb4761013bb76c0ccc56d48d4512e35f5251c6.tar.gz norta-41cb4761013bb76c0ccc56d48d4512e35f5251c6.tar.bz2 norta-41cb4761013bb76c0ccc56d48d4512e35f5251c6.zip |
add module to generate OSM basemaps
this is a simple script to ease making visuals for other projects I'm
working on with this dataset.
-rw-r--r-- | basemap.py | 70 | ||||
-rw-r--r-- | readme.rst | 73 |
2 files changed, 137 insertions, 6 deletions
diff --git a/basemap.py b/basemap.py new file mode 100644 index 0000000..86e2ea7 --- /dev/null +++ b/basemap.py @@ -0,0 +1,70 @@ + +from itertools import product +from io import BytesIO +from PIL import Image +import requests +import math + +URL = "https://tile.openstreetmap.org/{z}/{x}/{y}.png" +RADIUS = 6371 # earth's average radius in km +TILE_SIZE = 256 # each tile's size in pixels + + +def point_to_pixels(lon, lat, zoom): + """convert gps coordinates to web mercator""" + r = math.pow(2, zoom) * TILE_SIZE + lat = math.radians(lat) + + x = int((lon + 180.0) / 360.0 * r) + y = int((1.0 - math.log(math.tan(lat) + (1.0 / math.cos(lat))) / + math.pi) / 2.0 * r) + + return x, y + + +def tile(session=None, **url_format_args): + """download tile as PIL.Image from Tile Server API""" + url = url_format_args.get('url', URL) + if not session: + session = requests + + with session.get(url.format(**url_format_args)) as resp: + return Image.open(BytesIO(resp.content)) + + +def image(top, right, bottom, left, zoom=14, **tile_args): + """return an osm map from the given bounding box points""" + + # convert gps coordinates to web mercator + x0, y0 = point_to_pixels(left, top, zoom) + x1, y1 = point_to_pixels(right, bottom, zoom) + + # conver pixel coordinates to tiles + x0_tile, y0_tile = int(x0 / TILE_SIZE), int(y0 / TILE_SIZE) + x1_tile, y1_tile = math.ceil(x1 / TILE_SIZE), math.ceil(y1 / TILE_SIZE) + + assert (x1_tile - x0_tile) * ( + y1_tile - y0_tile) < 50, "That's too many tiles!" + + # crate a single large image to add the tiles to + img = Image.new('RGB', ( + int(x1_tile - x0_tile) * TILE_SIZE, + int(y1_tile - y0_tile) * TILE_SIZE)) + + # download all tiles + with requests.Session() as session: + for x_tile, y_tile in product( + range(x0_tile, x1_tile), range(y0_tile, y1_tile)): + + img.paste( + im=tile(session, x=x_tile, y=y_tile, z=zoom, **tile_args), + box=(int(x_tile - x0_tile) * TILE_SIZE, + int(y_tile - y0_tile) * TILE_SIZE)) + + # crop image to the given bounding box + x, y = x0_tile * TILE_SIZE, y0_tile * TILE_SIZE + return img.crop(( + int(x0 - x), # left + int(y0 - y), # top + int(x1 - x), # right + int(y1 - y))) # bottom @@ -25,8 +25,8 @@ the API. This is on the larger end of what I feel comfortable publishing online. So if you wish for a copy please `send a DM or email me <https://bryanbrattlof.com/connect/>`__ and I'll gladly send it to you. -Preparing The Data -################## +prepare-data.py +############### `prepare-data.py <https://git.bryanbrattlof.com/norta/tree/prepare-data.py>`__: is a small script to convert the ``bus.log.tar.gz`` file into a CSV file named @@ -53,12 +53,73 @@ is a small script to convert the ``bus.log.tar.gz`` file into a CSV file named ) df.set_index('epoch') -Base-map -######## +usage +----- -**TODO** +With the ``bus.log.tar.gz`` inside the ``data`` directory, simply run the python +script to generate the ``bus.csv`` file inside the ``data`` directory like so: -The full write-up is available at +.. code-block:: + + $ python prepare-data.py + +There is no data cleaning involved. This script will only decompress and convert +the log file into a csv file. + +basemap.py +########## + +`basemap.py <https://git.bryanbrattlof.com/norta/tree/basemap.py>`__: is a simple +module that downloads and combines OpenStreetMap tiles into one large image you +can add to your matplotlib visuals. + +dependencies +------------ + +It assumes you have `requests <https://requests.readthedocs.io/en/master/>`__, and +`pillow <https://python-pillow.org/>`__ libraries installed. + +usage +----- + +To use this module, pass the bounding box in GPS coordinates of the area into +the ``top``, ``rgt``, ``bot``, ``lef`` arguments along with the appropriate +``zoom`` level: + +.. code-block:: python + + import basemap + + top, bot = df.lat.max(), df.lat.min() + lef, rgt = df.lon.min(), df.lon.max() + + img = basemap.image(top, rgt, bot, lef, zoom=13) + +The module will return a ``Pillow.Image`` object that you can add to your +matplotlib visuals like this: + +.. code-block:: python + + from matplotlib import pyplot as plt + fig, ax = plt.subplots() + ax.imshow(img, extent=(lef, rgt, bot, top), aspect= 'equal') + plt.show() + +You can also use ``url`` to specify which tile servers you want to use: + +.. code-block:: python + + img = basemap.image(top, rgt, bot, lef, zoom=13, + url="http://c.tile.stamen.com/toner/{z}/{x}/{y}.png") + +Any extra arguments to needed for the ``url`` argument can be passed along as key word arguments in the ``basemap.image()`` function. For example: + +.. code-block:: python + + img = basemap.image(top, rgt, bot, lef, zoom=13, api=API_KEY + url="http://tileserver.example.com/{api}/{z}/{x}/{y}.png") + +The full write-up of how this module works is available at https://bryanbrattlof.com/adding-openstreetmaps-to-matplotlib/ Contributing |