diff --git a/requirements.txt b/requirements.txt index f31556f6..f3d9dfb3 100644 --- a/requirements.txt +++ b/requirements.txt @@ -25,3 +25,4 @@ wsgiref==0.1.2 zope.dottedname==4.1.0 edtf==0.9.3 mapbox-vector-tile==1.2.0 +git+https://github.com/tilezen/raw_tiles@master#egg=raw_tiles diff --git a/setup.py b/setup.py index 6087f149..16c4a799 100644 --- a/setup.py +++ b/setup.py @@ -44,6 +44,7 @@ 'pyproj', 'python-dateutil', 'PyYAML', + 'raw_tiles', 'redis', 'requests', 'Shapely', diff --git a/tests/test_query_fixture.py b/tests/test_query_fixture.py index 7dbb04b8..767d6d7f 100644 --- a/tests/test_query_fixture.py +++ b/tests/test_query_fixture.py @@ -1,14 +1,19 @@ import unittest -class TestQueryFixture(unittest.TestCase): +class FixtureTestCase(unittest.TestCase): def _make(self, rows, min_zoom_fn, props_fn, relations=[], - layer_name='testlayer'): - from tilequeue.query.fixture import LayerInfo + layer_name='testlayer', label_placement_layers={}): + from tilequeue.query.common import LayerInfo from tilequeue.query.fixture import make_fixture_data_fetcher layers = {layer_name: LayerInfo(min_zoom_fn, props_fn)} - return make_fixture_data_fetcher(layers, rows, relations=relations) + return make_fixture_data_fetcher( + layers, rows, label_placement_layers=label_placement_layers, + relations=relations) + + +class TestQueryFixture(FixtureTestCase): def test_query_simple(self): from shapely.geometry import Point @@ -141,3 +146,173 @@ def _rel(id, nodes=None, ways=None, rels=None): _rel(4, rels=[3]), _rel(5, rels=[2, 4]), ], 5) + + +class TestLabelPlacement(FixtureTestCase): + + def _test(self, layer_name, props): + from ModestMaps.Core import Coordinate + from tilequeue.tile import coord_to_mercator_bounds + from shapely.geometry import box + + def min_zoom_fn(shape, props, fid, meta): + return 0 + + tile = Coordinate(zoom=15, column=0, row=0) + bounds = coord_to_mercator_bounds(tile) + shape = box(*bounds) + + rows = [ + (1, shape, props), + ] + + label_placement_layers = { + 'polygon': set([layer_name]), + } + fetch = self._make( + rows, min_zoom_fn, None, relations=[], layer_name=layer_name, + label_placement_layers=label_placement_layers) + + read_rows = fetch(16, bounds) + return read_rows + + def test_named_item(self): + from shapely import wkb + + layer_name = 'testlayer' + read_rows = self._test(layer_name, {'name': 'Foo'}) + + self.assertEquals(1, len(read_rows)) + + label_prop = '__label__' + self.assertTrue(label_prop in read_rows[0]) + point = wkb.loads(read_rows[0][label_prop]) + self.assertEqual(point.geom_type, 'Point') + + +class TestGeometryClipping(FixtureTestCase): + + def _test(self, layer_name, bounds, factor): + from shapely.geometry import box + + def min_zoom_fn(shape, props, fid, meta): + return 0 + + boxwidth = bounds[2] - bounds[0] + boxheight = bounds[3] - bounds[1] + # make shape overlap the edges of the bounds. that way we can check to + # see if the shape gets clipped. + shape = box(bounds[0] - factor * boxwidth, + bounds[1] - factor * boxheight, + bounds[2] + factor * boxwidth, + bounds[3] + factor * boxheight) + + props = {'name': 'Foo'} + + rows = [ + (1, shape, props), + ] + + fetch = self._make( + rows, min_zoom_fn, None, relations=[], layer_name=layer_name) + + read_rows = fetch(16, bounds) + self.assertEqual(1, len(read_rows)) + return read_rows[0] + + def test_normal_layer(self): + from ModestMaps.Core import Coordinate + from tilequeue.tile import coord_to_mercator_bounds + from shapely import wkb + + tile = Coordinate(zoom=15, column=10, row=10) + bounds = coord_to_mercator_bounds(tile) + + read_row = self._test('testlayer', bounds, 1.0) + clipped_shape = wkb.loads(read_row['__geometry__']) + # for normal layers, clipped shape is inside the bounds of the tile. + x_factor = ((clipped_shape.bounds[2] - clipped_shape.bounds[0]) / + (bounds[2] - bounds[0])) + y_factor = ((clipped_shape.bounds[2] - clipped_shape.bounds[0]) / + (bounds[2] - bounds[0])) + self.assertAlmostEqual(1.0, x_factor) + self.assertAlmostEqual(1.0, y_factor) + + def test_water_layer(self): + # water layer should be expanded by 10% on each side. + from ModestMaps.Core import Coordinate + from tilequeue.tile import coord_to_mercator_bounds + from shapely import wkb + + tile = Coordinate(zoom=15, column=10, row=10) + bounds = coord_to_mercator_bounds(tile) + + read_row = self._test('water', bounds, 1.0) + clipped_shape = wkb.loads(read_row['__geometry__']) + # for water layer, the geometry should be 10% larger than the tile + # bounds. + x_factor = ((clipped_shape.bounds[2] - clipped_shape.bounds[0]) / + (bounds[2] - bounds[0])) + y_factor = ((clipped_shape.bounds[2] - clipped_shape.bounds[0]) / + (bounds[2] - bounds[0])) + self.assertAlmostEqual(1.1, x_factor) + self.assertAlmostEqual(1.1, y_factor) + + +class TestNameHandling(FixtureTestCase): + + def _test(self, input_layer_names, expected_layer_names): + from shapely.geometry import Point + from tilequeue.query.common import LayerInfo + from tilequeue.query.fixture import make_fixture_data_fetcher + from tilequeue.tile import coord_to_mercator_bounds + from tilequeue.tile import mercator_point_to_coord + + def min_zoom_fn(shape, props, fid, meta): + return 0 + + def props_fn(shape, props, fid, meta): + return {} + + shape = Point(0, 0) + props = {'name': 'Foo'} + + rows = [ + (1, shape, props), + ] + + layers = {} + for name in input_layer_names: + layers[name] = LayerInfo(min_zoom_fn, props_fn) + fetch = make_fixture_data_fetcher(layers, rows) + + feature_coord = mercator_point_to_coord(16, shape.x, shape.y) + read_rows = fetch(16, coord_to_mercator_bounds(feature_coord)) + self.assertEqual(1, len(read_rows)) + + all_names = set(expected_layer_names) | set(input_layer_names) + for name in all_names: + properties_name = '__%s_properties__' % name + self.assertTrue(properties_name in read_rows[0]) + actual_name = read_rows[0][properties_name].get('name') + if name in expected_layer_names: + expected_name = props.get('name') + self.assertEquals(expected_name, actual_name) + else: + # check the name doesn't appear anywhere else + self.assertEquals(None, actual_name) + + def test_name_single_layer(self): + # in any oone of the pois, landuse or buildings layers, a name + # by itself will be output in the same layer. + for layer_name in ('pois', 'landuse', 'buildings'): + self._test([layer_name], [layer_name]) + + def test_precedence(self): + # if the feature is in the pois layer, then that should get the name + # and the other layers should not. + self._test(['pois', 'landuse'], ['pois']) + self._test(['pois', 'buildings'], ['pois']) + self._test(['pois', 'landuse', 'buildings'], ['pois']) + # otherwise, landuse should take precedence over buildings. + self._test(['landuse', 'buildings'], ['landuse']) diff --git a/tests/test_query_rawr.py b/tests/test_query_rawr.py new file mode 100644 index 00000000..14fdaa69 --- /dev/null +++ b/tests/test_query_rawr.py @@ -0,0 +1,582 @@ +import unittest + + +class TestGetTable(object): + """ + Mocks the interface expected by raw_tiles.index.index.index_table, + which provides "table lookup". Here, we just return static stuff + previously set up in the test. + """ + + def __init__(self, tables): + self.tables = tables + + def __call__(self, table_name): + return self.tables.get(table_name, []) + + +class RawrTestCase(unittest.TestCase): + """ + Base layer of the tests, providing a utility function to create a data + fetcher with a set of mocked data. + """ + + def _make(self, min_zoom_fn, props_fn, tables, tile_pyramid, + layer_name='testlayer', label_placement_layers={}, + source='test'): + from tilequeue.query.common import LayerInfo + from tilequeue.query.rawr import make_rawr_data_fetcher + + layers = {layer_name: LayerInfo(min_zoom_fn, props_fn)} + return make_rawr_data_fetcher( + layers, tables, tile_pyramid, source, + label_placement_layers=label_placement_layers) + + +class TestQueryRawr(RawrTestCase): + + def test_query_simple(self): + # just check that we can get back the mock data we put into a tile, + # and that the indexing/fetching code respects the tile boundary and + # min_zoom function. + + from shapely.geometry import Point + from tilequeue.query.rawr import TilePyramid + from tilequeue.tile import coord_to_mercator_bounds + from tilequeue.tile import mercator_point_to_coord + + feature_min_zoom = 11 + + def min_zoom_fn(shape, props, fid, meta): + return feature_min_zoom + + shape = Point(0, 0) + # get_table(table_name) should return a generator of rows. + tables = TestGetTable({ + 'planet_osm_point': [(0, shape.wkb, {})], + }) + + zoom = 10 + max_zoom = zoom + 5 + coord = mercator_point_to_coord(zoom, shape.x, shape.y) + tile_pyramid = TilePyramid(zoom, coord.column, coord.row, max_zoom) + + fetch = self._make(min_zoom_fn, None, tables, tile_pyramid) + + # first, check that it can get the original item back when both the + # min zoom filter and geometry filter are okay. + feature_coord = mercator_point_to_coord( + feature_min_zoom, shape.x, shape.y) + read_rows = fetch( + feature_min_zoom, coord_to_mercator_bounds(feature_coord)) + + self.assertEquals(1, len(read_rows)) + read_row = read_rows[0] + self.assertEquals(0, read_row.get('__id__')) + # query processing code expects WKB bytes in the __geometry__ column + self.assertEquals(shape.wkb, read_row.get('__geometry__')) + self.assertEquals({'min_zoom': 11}, + read_row.get('__testlayer_properties__')) + + # now, check that if the min zoom or geometry filters would exclude + # the feature then it isn't returned. + read_rows = fetch(zoom, coord_to_mercator_bounds(coord)) + self.assertEquals(0, len(read_rows)) + + read_rows = fetch( + feature_min_zoom, coord_to_mercator_bounds(feature_coord.left())) + self.assertEquals(0, len(read_rows)) + + def test_query_min_zoom_fraction(self): + # test that fractional min zooms are included in their "floor" zoom + # tile. this is to allow over-zooming of a zoom N tile until N+1, + # where the next zoom tile kicks in. + + from shapely.geometry import Point + from tilequeue.query.rawr import TilePyramid + from tilequeue.tile import coord_to_mercator_bounds + from tilequeue.tile import mercator_point_to_coord + + def min_zoom_fn(shape, props, fid, meta): + return 11.999 + + shape = Point(0, 0) + tables = TestGetTable({ + 'planet_osm_point': [(0, shape.wkb, {})] + }) + + zoom = 10 + max_zoom = zoom + 5 + coord = mercator_point_to_coord(zoom, shape.x, shape.y) + tile_pyramid = TilePyramid(zoom, coord.column, coord.row, max_zoom) + + fetch = self._make(min_zoom_fn, None, tables, tile_pyramid) + + # check that the fractional zoom of 11.999 means that it's included in + # the zoom 11 tile, but not the zoom 10 one. + feature_coord = mercator_point_to_coord(11, shape.x, shape.y) + read_rows = fetch(11, coord_to_mercator_bounds(feature_coord)) + self.assertEquals(1, len(read_rows)) + + feature_coord = feature_coord.zoomBy(-1).container() + read_rows = fetch(10, coord_to_mercator_bounds(feature_coord)) + self.assertEquals(0, len(read_rows)) + + def test_query_past_max_zoom(self): + # check that features with a min_zoom beyond the maximum zoom are still + # included at the maximum zoom. since this is the last zoom level we + # generate, it must include everything. + + from shapely.geometry import Point + from tilequeue.query.rawr import TilePyramid + from tilequeue.tile import coord_to_mercator_bounds + from tilequeue.tile import mercator_point_to_coord + + def min_zoom_fn(shape, props, fid, meta): + return 20 + + shape = Point(0, 0) + tables = TestGetTable({ + 'planet_osm_point': [(0, shape.wkb, {})] + }) + + zoom = 10 + max_zoom = zoom + 6 + coord = mercator_point_to_coord(zoom, shape.x, shape.y) + tile_pyramid = TilePyramid(zoom, coord.column, coord.row, max_zoom) + + fetch = self._make(min_zoom_fn, None, tables, tile_pyramid) + + # the min_zoom of 20 should mean that the feature is included at zoom + # 16, even though 16<20, because 16 is the "max zoom" at which all the + # data is included. + feature_coord = mercator_point_to_coord(16, shape.x, shape.y) + read_rows = fetch(16, coord_to_mercator_bounds(feature_coord)) + self.assertEquals(1, len(read_rows)) + + # but it should not exist at zoom 15 + feature_coord = feature_coord.zoomBy(-1).container() + read_rows = fetch(10, coord_to_mercator_bounds(feature_coord)) + self.assertEquals(0, len(read_rows)) + + def test_root_relation_id(self): + # check the logic for finding a root relation ID for station complexes. + + from shapely.geometry import Point + from tilequeue.query.rawr import TilePyramid + from tilequeue.tile import coord_to_mercator_bounds + from tilequeue.tile import mercator_point_to_coord + + def min_zoom_fn(shape, props, fid, meta): + return 10 + + def _test(rels, expected_root_id): + shape = Point(0, 0) + props = { + 'railway': 'station', + 'name': 'Foo Station', + } + tables = TestGetTable({ + 'planet_osm_point': [(1, shape.wkb, props)], + 'planet_osm_rels': rels, + }) + + zoom = 10 + max_zoom = zoom + 6 + coord = mercator_point_to_coord(zoom, shape.x, shape.y) + tile_pyramid = TilePyramid(zoom, coord.column, coord.row, max_zoom) + + fetch = self._make(min_zoom_fn, None, tables, tile_pyramid, + layer_name='pois') + + feature_coord = mercator_point_to_coord(16, shape.x, shape.y) + read_rows = fetch(16, coord_to_mercator_bounds(feature_coord)) + self.assertEquals(1, len(read_rows)) + + props = read_rows[0]['__pois_properties__'] + self.assertEquals(expected_root_id, + props.get('mz_transit_root_relation_id')) + + # the fixture code expects "raw" relations as if they come straight + # from osm2pgsql. the structure is a little cumbersome, so this + # utility function constructs it from a more readable function call. + def _rel(id, nodes=None, ways=None, rels=None): + way_off = len(nodes) if nodes else 0 + rel_off = way_off + (len(ways) if ways else 0) + parts = (nodes or []) + (ways or []) + (rels or []) + members = [""] * len(parts) + tags = ['type', 'site'] + return (id, way_off, rel_off, parts, members, tags) + + # one level of relations - this one directly contains the station + # node. + _test([_rel(2, nodes=[1])], 2) + + # two levels of relations r3 contains r2 contains n1. + _test([_rel(2, nodes=[1]), _rel(3, rels=[2])], 3) + + # asymmetric diamond pattern. r2 and r3 both contain n1, r4 contains + # r3 and r5 contains both r4 and r2, making it the "top" relation. + _test([ + _rel(2, nodes=[1]), + _rel(3, nodes=[1]), + _rel(4, rels=[3]), + _rel(5, rels=[2, 4]), + ], 5) + + +class TestLabelPlacement(RawrTestCase): + + def _test(self, layer_name, props): + from ModestMaps.Core import Coordinate + from tilequeue.tile import coord_to_mercator_bounds + from shapely.geometry import box + from tilequeue.query.rawr import TilePyramid + + top_zoom = 10 + max_zoom = top_zoom + 6 + + def min_zoom_fn(shape, props, fid, meta): + return top_zoom + + tile = Coordinate(zoom=15, column=0, row=0) + top_tile = tile.zoomTo(top_zoom).container() + tile_pyramid = TilePyramid( + top_zoom, top_tile.column, top_tile.row, max_zoom) + + bounds = coord_to_mercator_bounds(tile) + shape = box(*bounds) + tables = TestGetTable({ + 'planet_osm_polygon': [ + (1, shape.wkb, props), + ] + }) + + label_placement_layers = { + 'polygon': set([layer_name]), + } + fetch = self._make( + min_zoom_fn, None, tables, tile_pyramid, layer_name=layer_name, + label_placement_layers=label_placement_layers) + + read_rows = fetch(tile.zoom, bounds) + return read_rows + + def test_named_item(self): + # check that a label is generated for features in label placement + # layers. + + from shapely import wkb + + layer_name = 'testlayer' + read_rows = self._test(layer_name, {'name': 'Foo'}) + + self.assertEquals(1, len(read_rows)) + + label_prop = '__label__' + self.assertTrue(label_prop in read_rows[0]) + point = wkb.loads(read_rows[0][label_prop]) + self.assertEqual(point.geom_type, 'Point') + + +class TestGeometryClipping(RawrTestCase): + + def _test(self, layer_name, tile, factor): + from shapely.geometry import box + from tilequeue.query.rawr import TilePyramid + from tilequeue.tile import coord_to_mercator_bounds + + top_zoom = 10 + max_zoom = top_zoom + 6 + + def min_zoom_fn(shape, props, fid, meta): + return top_zoom + + top_tile = tile.zoomTo(top_zoom).container() + tile_pyramid = TilePyramid( + top_zoom, top_tile.column, top_tile.row, max_zoom) + + bounds = coord_to_mercator_bounds(tile) + boxwidth = bounds[2] - bounds[0] + boxheight = bounds[3] - bounds[1] + # make shape overlap the edges of the bounds. that way we can check to + # see if the shape gets clipped. + shape = box(bounds[0] - factor * boxwidth, + bounds[1] - factor * boxheight, + bounds[2] + factor * boxwidth, + bounds[3] + factor * boxheight) + + props = {'name': 'Foo'} + + tables = TestGetTable({ + 'planet_osm_polygon': [ + (1, shape.wkb, props), + ], + }) + + fetch = self._make( + min_zoom_fn, None, tables, tile_pyramid, layer_name=layer_name) + + read_rows = fetch(tile.zoom, bounds) + self.assertEqual(1, len(read_rows)) + return read_rows[0] + + def test_normal_layer(self): + # check that normal layer geometries are clipped to the bounding box of + # the tile. + + from ModestMaps.Core import Coordinate + from tilequeue.tile import coord_to_mercator_bounds + from shapely import wkb + + tile = Coordinate(zoom=15, column=10, row=10) + bounds = coord_to_mercator_bounds(tile) + + read_row = self._test('testlayer', tile, 1.0) + clipped_shape = wkb.loads(read_row['__geometry__']) + # for normal layers, clipped shape is inside the bounds of the tile. + x_factor = ((clipped_shape.bounds[2] - clipped_shape.bounds[0]) / + (bounds[2] - bounds[0])) + y_factor = ((clipped_shape.bounds[2] - clipped_shape.bounds[0]) / + (bounds[2] - bounds[0])) + self.assertAlmostEqual(1.0, x_factor) + self.assertAlmostEqual(1.0, y_factor) + + def test_water_layer(self): + # water layer should be clipped to the tile bounds expanded by 10%. + + from ModestMaps.Core import Coordinate + from tilequeue.tile import coord_to_mercator_bounds + from shapely import wkb + + tile = Coordinate(zoom=15, column=10, row=10) + bounds = coord_to_mercator_bounds(tile) + + read_row = self._test('water', tile, 1.0) + clipped_shape = wkb.loads(read_row['__geometry__']) + # for water layer, the geometry should be 10% larger than the tile + # bounds. + x_factor = ((clipped_shape.bounds[2] - clipped_shape.bounds[0]) / + (bounds[2] - bounds[0])) + y_factor = ((clipped_shape.bounds[2] - clipped_shape.bounds[0]) / + (bounds[2] - bounds[0])) + self.assertAlmostEqual(1.1, x_factor) + self.assertAlmostEqual(1.1, y_factor) + + +class TestNameHandling(RawrTestCase): + + def _test(self, input_layer_names, expected_layer_names): + from shapely.geometry import Point + from tilequeue.query.common import LayerInfo + from tilequeue.query.rawr import TilePyramid + from tilequeue.query.rawr import make_rawr_data_fetcher + from tilequeue.tile import coord_to_mercator_bounds + from tilequeue.tile import mercator_point_to_coord + + top_zoom = 10 + max_zoom = top_zoom + 6 + + def min_zoom_fn(shape, props, fid, meta): + return top_zoom + + def props_fn(shape, props, fid, meta): + return {} + + shape = Point(0, 0) + props = {'name': 'Foo'} + + tables = TestGetTable({ + 'planet_osm_point': [ + (1, shape.wkb, props), + ], + }) + + tile = mercator_point_to_coord(16, shape.x, shape.y) + top_tile = tile.zoomTo(top_zoom).container() + tile_pyramid = TilePyramid( + top_zoom, top_tile.column, top_tile.row, max_zoom) + + layers = {} + for name in input_layer_names: + layers[name] = LayerInfo(min_zoom_fn, props_fn) + source = 'test' + fetch = make_rawr_data_fetcher(layers, tables, tile_pyramid, source) + + read_rows = fetch(tile.zoom, coord_to_mercator_bounds(tile)) + # the RAWR query goes over features multiple times because of the + # indexing, so we can't rely on all the properties for one feature to + # be all together in the same place. this loops over all the features, + # checking that there's only really one of them and gathering together + # all the __%s_properties__ from all the rows for further testing. + all_props = {} + for row in read_rows: + self.assertEquals(1, row['__id__']) + self.assertEquals(shape.wkb, row['__geometry__']) + for key, val in row.items(): + if key.endswith('_properties__'): + self.assertFalse(key in all_props) + all_props[key] = val + + all_names = set(expected_layer_names) | set(input_layer_names) + for name in all_names: + properties_name = '__%s_properties__' % name + self.assertTrue(properties_name in all_props) + actual_name = all_props[properties_name].get('name') + if name in expected_layer_names: + expected_name = props.get('name') + self.assertEquals(expected_name, actual_name) + else: + # check the name doesn't appear anywhere else + self.assertEquals(None, actual_name) + + def test_name_single_layer(self): + # in any oone of the pois, landuse or buildings layers, a name + # by itself will be output in the same layer. + for layer_name in ('pois', 'landuse', 'buildings'): + self._test([layer_name], [layer_name]) + + def test_precedence(self): + # if the feature is in the pois layer, then that should get the name + # and the other layers should not. + self._test(['pois', 'landuse'], ['pois']) + self._test(['pois', 'buildings'], ['pois']) + self._test(['pois', 'landuse', 'buildings'], ['pois']) + # otherwise, landuse should take precedence over buildings. + self._test(['landuse', 'buildings'], ['landuse']) + + +class TestMeta(RawrTestCase): + + def test_meta_gate(self): + # test that the meta is passed to the min zoom function, and that we + # can use it to get information about the highways that a gate is + # part of. + + from shapely.geometry import Point, LineString + from tilequeue.query.rawr import TilePyramid + from tilequeue.tile import coord_to_mercator_bounds + from tilequeue.tile import mercator_point_to_coord + + feature_min_zoom = 11 + + def min_zoom_fn(shape, props, fid, meta): + self.assertIsNotNone(meta) + + # expect meta to have a source, which is a string name for the + # source of the data. + self.assertEquals('test', meta.source) + + # expect meta to have a list of relations, which is empty for this + # test. + self.assertEquals(0, len(meta.relations)) + + # only do this for the node + if fid == 0: + # expect meta to have a list of ways, each of which is a (fid, + # shape, props) tuple, of which only props is used. + self.assertEquals(1, len(meta.ways)) + way_fid, way_shape, way_props = meta.ways[0] + self.assertEquals(1, way_fid) + self.assertEquals({'highway': 'secondary'}, way_props) + + # only set a min zoom for the node - this just simplifies the + # checking later, as there'll only be one feature. + return feature_min_zoom if fid == 0 else None + + shape = Point(0, 0) + way_shape = LineString([[0, 0], [1, 1]]) + # get_table(table_name) should return a generator of rows. + tables = TestGetTable({ + 'planet_osm_point': [(0, shape.wkb, {'barrier': 'gate'})], + 'planet_osm_line': [(1, way_shape.wkb, {'highway': 'secondary'})], + 'planet_osm_ways': [(1, [0], ['highway', 'secondary'])], + }) + + zoom = 10 + max_zoom = zoom + 5 + coord = mercator_point_to_coord(zoom, shape.x, shape.y) + tile_pyramid = TilePyramid(zoom, coord.column, coord.row, max_zoom) + + fetch = self._make(min_zoom_fn, None, tables, tile_pyramid) + + # first, check that it can get the original item back when both the + # min zoom filter and geometry filter are okay. + feature_coord = mercator_point_to_coord( + feature_min_zoom, shape.x, shape.y) + read_rows = fetch( + feature_min_zoom, coord_to_mercator_bounds(feature_coord)) + + self.assertEquals(1, len(read_rows)) + read_row = read_rows[0] + self.assertEquals(0, read_row.get('__id__')) + # query processing code expects WKB bytes in the __geometry__ column + self.assertEquals(shape.wkb, read_row.get('__geometry__')) + self.assertEquals({'min_zoom': 11, 'barrier': 'gate'}, + read_row.get('__testlayer_properties__')) + + def test_meta_route(self): + # test that we can use meta in the min zoom function to find out which + # route(s) a road is part of. + + from shapely.geometry import LineString + from tilequeue.query.rawr import TilePyramid + from tilequeue.tile import coord_to_mercator_bounds + from tilequeue.tile import mercator_point_to_coord + from tilequeue.query.common import deassoc + + feature_min_zoom = 11 + + rel_tags = [ + 'type', 'route', + 'route', 'road', + 'ref', '101', + ] + + def min_zoom_fn(shape, props, fid, meta): + self.assertIsNotNone(meta) + + # expect meta to have a source, which is a string name for the + # source of the data. + self.assertEquals('test', meta.source) + + # expect meta to have a list of ways, but empty for this test. + self.assertEquals(0, len(meta.ways)) + + # expect meta to have a list of relations, each of which is a dict + # containing at least the key 'tags' mapped to a list of + # alternating k, v suitable for passing into deassoc(). + self.assertEquals(1, len(meta.relations)) + rel = meta.relations[0] + self.assertIsInstance(rel, dict) + self.assertIn('tags', rel) + self.assertEquals(deassoc(rel_tags), deassoc(rel['tags'])) + + return feature_min_zoom + + shape = LineString([[0, 0], [1, 1]]) + # get_table(table_name) should return a generator of rows. + tables = TestGetTable({ + 'planet_osm_line': [(1, shape.wkb, {'highway': 'secondary'})], + 'planet_osm_rels': [(2, 0, 1, [1], [''], rel_tags)], + }) + + zoom = 10 + max_zoom = zoom + 5 + coord = mercator_point_to_coord(zoom, *shape.coords[0]) + tile_pyramid = TilePyramid(zoom, coord.column, coord.row, max_zoom) + + fetch = self._make(min_zoom_fn, None, tables, tile_pyramid) + + # first, check that it can get the original item back when both the + # min zoom filter and geometry filter are okay. + feature_coord = mercator_point_to_coord( + feature_min_zoom, *shape.coords[0]) + read_rows = fetch( + feature_min_zoom, coord_to_mercator_bounds(feature_coord)) + + self.assertEquals(1, len(read_rows)) + read_row = read_rows[0] + self.assertEquals(1, read_row.get('__id__')) + self.assertEquals({'min_zoom': 11, 'highway': 'secondary'}, + read_row.get('__testlayer_properties__')) diff --git a/tilequeue/query/common.py b/tilequeue/query/common.py new file mode 100644 index 00000000..bbbc29ab --- /dev/null +++ b/tilequeue/query/common.py @@ -0,0 +1,334 @@ +from collections import namedtuple +from collections import defaultdict +from itertools import izip + + +def namedtuple_with_defaults(name, props, defaults): + t = namedtuple(name, props) + t.__new__.__defaults__ = defaults + return t + + +class LayerInfo(namedtuple_with_defaults( + 'LayerInfo', 'min_zoom_fn props_fn shape_types', (None,))): + + def allows_shape_type(self, shape): + if self.shape_types is None: + return True + typ = shape_type_lookup(shape) + return typ in self.shape_types + + +def deassoc(x): + """ + Turns an array consisting of alternating key-value pairs into a + dictionary. + + Osm2pgsql stores the tags for ways and relations in the planet_osm_ways and + planet_osm_rels tables in this format. Hstore would make more sense now, + but this encoding pre-dates the common availability of hstore. + + Example: + >>> from raw_tiles.index.util import deassoc + >>> deassoc(['a', 1, 'b', 'B', 'c', 3.14]) + {'a': 1, 'c': 3.14, 'b': 'B'} + """ + + pairs = [iter(x)] * 2 + return dict(izip(*pairs)) + + +# fixtures extend metadata to include ways and relations for the feature. +# this is unnecessary for SQL, as the ways and relations tables are +# "ambiently available" and do not need to be passed in arguments. +Metadata = namedtuple('Metadata', 'source ways relations') + + +def shape_type_lookup(shape): + typ = shape.geom_type + if typ.startswith('Multi'): + typ = typ[len('Multi'):] + return typ.lower() + + +# list of road types which are likely to have buses on them. used to cut +# down the number of queries the SQL used to do for relations. although this +# isn't necessary for fixtures, we replicate the logic to keep the behaviour +# the same. +BUS_ROADS = set([ + 'motorway', 'motorway_link', 'trunk', 'trunk_link', 'primary', + 'primary_link', 'secondary', 'secondary_link', 'tertiary', + 'tertiary_link', 'residential', 'unclassified', 'road', 'living_street', +]) + + +class Relation(object): + def __init__(self, obj): + self.id = obj['id'] + self.tags = deassoc(obj['tags']) + way_off = obj['way_off'] + rel_off = obj['rel_off'] + self.node_ids = obj['parts'][0:way_off] + self.way_ids = obj['parts'][way_off:rel_off] + self.rel_ids = obj['parts'][rel_off:] + + +def mz_is_interesting_transit_relation(tags): + public_transport = tags.get('public_transport') + typ = tags.get('type') + return public_transport in ('stop_area', 'stop_area_group') or \ + typ in ('stop_area', 'stop_area_group', 'site') + + +# starting with the IDs in seed_relations, recurse up the transit relations +# of which they are members. returns the set of all the relation IDs seen +# and the "root" relation ID, which was the "furthest" relation from any +# leaf relation. +def mz_recurse_up_transit_relations(seed_relations, osm): + root_relation_ids = set() + root_relation_level = 0 + all_relations = set() + + for rel_id in seed_relations: + front = set([rel_id]) + seen = set([rel_id]) + level = 0 + + if root_relation_level == 0: + root_relation_ids.add(rel_id) + + while front: + new_rels = set() + for r in front: + new_rels |= osm.transit_relations(r) + new_rels -= seen + level += 1 + if new_rels and level > root_relation_level: + root_relation_ids = new_rels + root_relation_level = level + elif new_rels and level == root_relation_level: + root_relation_ids |= new_rels + front = new_rels + seen |= front + + all_relations |= seen + + root_relation_id = min(root_relation_ids) if root_relation_ids else None + return all_relations, root_relation_id + + +# extract a name for a transit route relation. this can expand comma +# separated lists and prefers to use the ref rather than the name. +def mz_transit_route_name(tags): + # prefer ref as it's less likely to contain the destination name + name = tags.get('ref') + if not name: + name = tags.get('name') + if name: + name = name.strip() + return name + + +def is_station_or_stop(fid, shape, props): + "Returns true if the given (point) feature is a station or stop." + return ( + props.get('railway') in ('station', 'stop', 'tram_stop') or + props.get('public_transport') in ('stop', 'stop_position', 'tram_stop') + ) + + +def is_station_or_line(fid, shape, props): + """ + Returns true if the given (line or polygon from way) feature is a station + or transit line. + """ + + railway = props.get('railway') + return railway in ('subway', 'light_rail', 'tram', 'rail') + + +Transit = namedtuple( + 'Transit', 'score root_relation_id ' + 'trains subways light_rails trams railways') + + +def mz_calculate_transit_routes_and_score(osm, node_id, way_id): + candidate_relations = set() + if node_id: + candidate_relations.update(osm.relations_using_node(node_id)) + if way_id: + candidate_relations.update(osm.relations_using_way(way_id)) + + seed_relations = set() + for rel_id in candidate_relations: + rel = osm.relation(rel_id) + if mz_is_interesting_transit_relation(rel.tags): + seed_relations.add(rel_id) + del candidate_relations + + # TODO: if the station is also a multipolygon relation? + + # this complex query does two recursive sweeps of the relations + # table starting from a seed set of relations which are or contain + # the original station. + # + # the first sweep goes "upwards" from relations to "parent" relations. if + # a relation R1 is a member of relation R2, then R2 will be included in + # this sweep as long as it has "interesting" tags, as defined by the + # function mz_is_interesting_transit_relation. + # + # the second sweep goes "downwards" from relations to "child" relations. + # if a relation R1 has a member R2 which is also a relation, then R2 will + # be included in this sweep as long as it also has "interesting" tags. + all_relations, root_relation_id = mz_recurse_up_transit_relations( + seed_relations, osm) + del seed_relations + + # collect all the interesting nodes - this includes the station node (if + # any) and any nodes which are members of found relations which have + # public transport tags indicating that they're stations or stops. + stations_and_stops = set() + for rel_id in all_relations: + rel = osm.relation(rel_id) + for node_id in rel.node_ids: + if is_station_or_stop(*osm.node(node_id)): + stations_and_stops.add(node_id) + + if node_id: + stations_and_stops.add(node_id) + + # collect any physical railway which includes any of the above + # nodes. + stations_and_lines = set() + for node_id in stations_and_stops: + for way_id in osm.ways_using_node(node_id): + if is_station_or_line(*osm.way(way_id)): + stations_and_lines.add(way_id) + + if way_id: + stations_and_lines.add(way_id) + + # collect all IDs together in one array to intersect with the parts arrays + # of route relations which may include them. + all_routes = set() + for lookup, ids in ((osm.relations_using_node, stations_and_stops), + (osm.relations_using_way, stations_and_lines), + (osm.relations_using_rel, all_relations)): + for i in ids: + for rel_id in lookup(i): + rel = osm.relation(rel_id) + if rel.tags.get('type') == 'route' and \ + rel.tags.get('route') in ('subway', 'light_rail', 'tram', + 'train', 'railway'): + all_routes.add(rel_id) + + routes_lookup = defaultdict(set) + for rel_id in all_routes: + rel = osm.relation(rel_id) + route = rel.tags.get('route') + if route: + route_name = mz_transit_route_name(rel.tags) + routes_lookup[route].add(route_name) + trains = routes_lookup['train'] + subways = routes_lookup['subway'] + light_rails = routes_lookup['light_rail'] + trams = routes_lookup['tram'] + railways = routes_lookup['railway'] + del routes_lookup + + # if a station is an interchange between mainline rail and subway or + # light rail, then give it a "bonus" boost of importance. + bonus = 2 if trains and (subways or light_rails) else 1 + + score = (100 * min(9, bonus * len(trains)) + + 10 * min(9, bonus * (len(subways) + len(light_rails))) + + min(9, len(trams) + len(railways))) + + return Transit(score=score, root_relation_id=root_relation_id, + trains=trains, subways=subways, light_rails=light_rails, + railways=railways, trams=trams) + + +# properties for a feature (fid, shape, props) in layer `layer_name` at zoom +# level `zoom`. also takes an `osm` parameter, which is an object which can +# be used to look up nodes, ways and relations and the relationships between +# them. +def layer_properties(fid, shape, props, layer_name, zoom, osm): + layer_props = props.copy() + + # need to make sure that the name is only applied to one of + # the pois, landuse or buildings layers - in that order of + # priority. + # + # TODO: do this for all name variants & translations + if layer_name in ('pois', 'landuse', 'buildings'): + layer_props.pop('name', None) + + # urgh, hack! + if layer_name == 'water' and shape.geom_type == 'Point': + layer_props['label_placement'] = True + + if shape.geom_type in ('Polygon', 'MultiPolygon'): + layer_props['area'] = shape.area + + if layer_name == 'roads' and \ + shape.geom_type in ('LineString', 'MultiLineString') and \ + fid >= 0: + mz_networks = [] + mz_cycling_networks = set() + mz_is_bus_route = False + for rel_id in osm.relations_using_way(fid): + rel = osm.relation(rel_id) + typ, route, network, ref = [rel.tags.get(k) for k in ( + 'type', 'route', 'network', 'ref')] + if route and (network or ref): + mz_networks.extend([route, network, ref]) + if typ == 'route' and \ + route in ('hiking', 'foot', 'bicycle') and \ + network in ('icn', 'ncn', 'rcn', 'lcn'): + mz_cycling_networks.add(network) + if typ == 'route' and route in ('bus', 'trolleybus'): + mz_is_bus_route = True + + mz_cycling_network = None + for cn in ('icn', 'ncn', 'rcn', 'lcn'): + if layer_props.get(cn) == 'yes' or \ + ('%s_ref' % cn) in layer_props or \ + cn in mz_cycling_networks: + mz_cycling_network = cn + break + + if mz_is_bus_route and \ + zoom >= 12 and \ + layer_props.get('highway') in BUS_ROADS: + layer_props['is_bus_route'] = True + + layer_props['mz_networks'] = mz_networks + if mz_cycling_network: + layer_props['mz_cycling_network'] = mz_cycling_network + + is_poi = layer_name == 'pois' + is_railway_station = props.get('railway') == 'station' + is_point_or_poly = shape.geom_type in ( + 'Point', 'MultiPoint', 'Polygon', 'MultiPolygon') + + if is_poi and is_railway_station and \ + is_point_or_poly and fid >= 0: + node_id = None + way_id = None + if shape.geom_type in ('Point', 'MultiPoint'): + node_id = fid + else: + way_id = fid + + transit = mz_calculate_transit_routes_and_score( + osm, node_id, way_id) + layer_props['mz_transit_score'] = transit.score + layer_props['mz_transit_root_relation_id'] = ( + transit.root_relation_id) + layer_props['train_routes'] = transit.trains + layer_props['subway_routes'] = transit.subways + layer_props['light_rail_routes'] = transit.light_rails + layer_props['tram_routes'] = transit.trams + + return layer_props diff --git a/tilequeue/query/fixture.py b/tilequeue/query/fixture.py index 14f9fc41..0cc77c8d 100644 --- a/tilequeue/query/fixture.py +++ b/tilequeue/query/fixture.py @@ -1,281 +1,118 @@ -from collections import namedtuple -from collections import defaultdict from shapely.geometry import box from tilequeue.process import lookup_source -from itertools import izip from tilequeue.transform import calculate_padded_bounds +from tilequeue.query.common import Metadata +from tilequeue.query.common import Relation +from tilequeue.query.common import layer_properties +from tilequeue.query.common import shape_type_lookup +from tilequeue.query.common import mz_is_interesting_transit_relation +from collections import defaultdict + + +class OsmFixtureLookup(object): + + def __init__(self, rows, rels): + # extract out all relations and index by ID. this is helpful when + # looking them up later. + relations = {} + nodes = {} + ways = {} + ways_using_node = {} + + for (fid, shape, props) in rows: + if fid >= 0: + if shape.geom_type in ('Point', 'MultiPoint'): + nodes[fid] = (fid, shape, props) + features = props.get('__ways__', []) + ways_using_node[fid] = [f[0] for f in features] + else: + ways[fid] = (fid, shape, props) + + for r in rels: + r = Relation(r) + assert r.id not in relations + relations[r.id] = r + + relations_using_node = defaultdict(list) + relations_using_way = defaultdict(list) + relations_using_rel = defaultdict(list) + + for rel_id, rel in relations.items(): + for (ids, index) in ((rel.node_ids, relations_using_node), + (rel.way_ids, relations_using_way), + (rel.rel_ids, relations_using_rel)): + for osm_id in ids: + index[osm_id].append(rel_id) + + transit_relations = defaultdict(set) + for rel_id, rel in relations.items(): + if mz_is_interesting_transit_relation(rel.tags): + for member in rel.rel_ids: + transit_relations[member].add(rel_id) + + # looks up relation IDs + self._relations_using_node = relations_using_node + self._relations_using_way = relations_using_way + self._relations_using_rel = relations_using_rel + # looks up way IDs + self._ways_using_node = ways_using_node + # looks up Relation objects + self._relations = relations + # looks up (fid, shape, props) feature objects + self._ways = ways + self._nodes = nodes + # returns the set of transit relation IDs that contain the given + # relation IDs + self._transit_relations = transit_relations + + def relations_using_node(self, node_id): + "Returns a list of relation IDs which contain the node with that ID." + + return self._relations_using_node.get(node_id, []) + + def relations_using_way(self, way_id): + "Returns a list of relation IDs which contain the way with that ID." + + return self._relations_using_way.get(way_id, []) + + def relations_using_rel(self, rel_id): + """ + Returns a list of relation IDs which contain the relation with that + ID. + """ + + return self._relations_using_rel.get(rel_id, []) + + def ways_using_node(self, node_id): + "Returns a list of way IDs which contain the node with that ID." + + return self._ways_using_node.get(node_id, []) + + def relation(self, rel_id): + "Returns the Relation object with the given ID." + + return self._relations[rel_id] + def way(self, way_id): + """ + Returns the feature (fid, shape, props) which was generated from the + given way. + """ + + return self._ways[way_id] + + def node(self, node_id): + """ + Returns the feature (fid, shape, props) which was generated from the + given node. + """ -def namedtuple_with_defaults(name, props, defaults): - t = namedtuple(name, props) - t.__new__.__defaults__ = defaults - return t + return self._nodes[node_id] + def transit_relations(self, rel_id): + "Return transit relations containing the relation with the given ID." -class LayerInfo(namedtuple_with_defaults( - 'LayerInfo', 'min_zoom_fn props_fn shape_types', (None,))): - - def allows_shape_type(self, shape): - if self.shape_types is None: - return True - typ = _shape_type_lookup(shape) - return typ in self.shape_types - - -def deassoc(x): - """ - Turns an array consisting of alternating key-value pairs into a - dictionary. - - Osm2pgsql stores the tags for ways and relations in the planet_osm_ways and - planet_osm_rels tables in this format. Hstore would make more sense now, - but this encoding pre-dates the common availability of hstore. - - Example: - >>> from raw_tiles.index.util import deassoc - >>> deassoc(['a', 1, 'b', 'B', 'c', 3.14]) - {'a': 1, 'c': 3.14, 'b': 'B'} - """ - - pairs = [iter(x)] * 2 - return dict(izip(*pairs)) - - -# fixtures extend metadata to include ways and relations for the feature. -# this is unnecessary for SQL, as the ways and relations tables are -# "ambiently available" and do not need to be passed in arguments. -Metadata = namedtuple('Metadata', 'source ways relations') - - -def _shape_type_lookup(shape): - typ = shape.geom_type - if typ.startswith('Multi'): - typ = typ[len('Multi'):] - return typ.lower() - - -# list of road types which are likely to have buses on them. used to cut -# down the number of queries the SQL used to do for relations. although this -# isn't necessary for fixtures, we replicate the logic to keep the behaviour -# the same. -BUS_ROADS = set([ - 'motorway', 'motorway_link', 'trunk', 'trunk_link', 'primary', - 'primary_link', 'secondary', 'secondary_link', 'tertiary', - 'tertiary_link', 'residential', 'unclassified', 'road', 'living_street', -]) - - -class Relation(object): - def __init__(self, obj): - self.id = obj['id'] - self.tags = deassoc(obj['tags']) - way_off = obj['way_off'] - rel_off = obj['rel_off'] - self.node_ids = obj['parts'][0:way_off] - self.way_ids = obj['parts'][way_off:rel_off] - self.rel_ids = obj['parts'][rel_off:] - - -def mz_is_interesting_transit_relation(tags): - public_transport = tags.get('public_transport') - typ = tags.get('type') - return public_transport in ('stop_area', 'stop_area_group') or \ - typ in ('stop_area', 'stop_area_group', 'site') - - -# starting with the IDs in seed_relations, recurse up the transit relations -# of which they are members. returns the set of all the relation IDs seen -# and the "root" relation ID, which was the "furthest" relation from any -# leaf relation. -def mz_recurse_up_transit_relations(seed_relations, relations): - transit_relations = defaultdict(set) - for rel_id, rel in relations.items(): - if mz_is_interesting_transit_relation(rel.tags): - for member in rel.rel_ids: - transit_relations[member].add(rel_id) - - root_relation_ids = set() - root_relation_level = 0 - all_relations = set() - - for rel_id in seed_relations: - front = set([rel_id]) - seen = set([rel_id]) - level = 0 - - if root_relation_level == 0: - root_relation_ids.add(rel_id) - - while front: - new_rels = set() - for r in front: - new_rels |= transit_relations[r] - new_rels -= seen - level += 1 - if new_rels and level > root_relation_level: - root_relation_ids = new_rels - root_relation_level = level - elif new_rels and level == root_relation_level: - root_relation_ids |= new_rels - front = new_rels - seen |= front - - all_relations |= seen - - root_relation_id = min(root_relation_ids) if root_relation_ids else None - return all_relations, root_relation_id - - -# extract a name for a transit route relation. this can expand comma -# separated lists and prefers to use the ref rather than the name. -def mz_transit_route_name(tags): - # prefer ref as it's less likely to contain the destination name - name = tags.get('ref') - if not name: - name = tags.get('name') - if name: - name = name.strip() - return name - - -Transit = namedtuple( - 'Transit', 'score root_relation_id ' - 'trains subways light_rails trams railways') - - -def mz_calculate_transit_routes_and_score(rows, rels, node_id, way_id): - # extract out all relations and index by ID. this is helpful when - # looking them up later. - relations = {} - nodes = {} - ways = {} - ways_using_node = {} - - for (fid, shape, props) in rows: - if fid >= 0: - if shape.geom_type in ('Point', 'MultiPoint'): - nodes[fid] = (fid, shape, props) - features = props.get('__ways__', []) - ways_using_node[fid] = [f[0] for f in features] - else: - ways[fid] = (fid, shape, props) - - for r in rels: - r = Relation(r) - assert r.id not in relations - relations[r.id] = r - - relations_using_node = defaultdict(list) - relations_using_way = defaultdict(list) - relations_using_rel = defaultdict(list) - - for rel_id, rel in relations.items(): - for (ids, index) in ((rel.node_ids, relations_using_node), - (rel.way_ids, relations_using_way), - (rel.rel_ids, relations_using_rel)): - for osm_id in ids: - index[osm_id].append(rel_id) - - candidate_relations = set() - if node_id: - candidate_relations.update(relations_using_node.get(node_id, [])) - if way_id: - candidate_relations.update(relations_using_way.get(way_id, [])) - - seed_relations = set() - for rel_id in candidate_relations: - rel = relations[rel_id] - if mz_is_interesting_transit_relation(rel.tags): - seed_relations.add(rel_id) - del candidate_relations - - # TODO: if the station is also a multipolygon relation? - - # this complex query does two recursive sweeps of the relations - # table starting from a seed set of relations which are or contain - # the original station. - # - # the first sweep goes "upwards" from relations to "parent" relations. if - # a relation R1 is a member of relation R2, then R2 will be included in - # this sweep as long as it has "interesting" tags, as defined by the - # function mz_is_interesting_transit_relation. - # - # the second sweep goes "downwards" from relations to "child" relations. - # if a relation R1 has a member R2 which is also a relation, then R2 will - # be included in this sweep as long as it also has "interesting" tags. - all_relations, root_relation_id = mz_recurse_up_transit_relations( - seed_relations, relations) - del seed_relations - - # collect all the interesting nodes - this includes the station node (if - # any) and any nodes which are members of found relations which have - # public transport tags indicating that they're stations or stops. - stations_and_stops = set() - for rel_id in all_relations: - rel = relations[rel_id] - for node_id in rel.node_ids: - fid, shape, props = nodes[node_id] - railway = props.get('railway') in ('station', 'stop', 'tram_stop') - public_transport = props.get('public_transport') in \ - ('stop', 'stop_position', 'tram_stop') - if railway or public_transport: - stations_and_stops.add(fid) - - if node_id: - stations_and_stops.add(node_id) - - # collect any physical railway which includes any of the above - # nodes. - stations_and_lines = set() - for node_id in stations_and_stops: - for way_id in ways_using_node[node_id]: - fid, shape, props = ways[way_id] - railway = props.get('railway') - if railway in ('subway', 'light_rail', 'tram', 'rail'): - stations_and_lines.add(way_id) - - if way_id: - stations_and_lines.add(way_id) - - # collect all IDs together in one array to intersect with the parts arrays - # of route relations which may include them. - all_routes = set() - for lookup, ids in ((relations_using_node, stations_and_stops), - (relations_using_way, stations_and_lines), - (relations_using_rel, all_relations)): - for i in ids: - for rel_id in lookup.get(i, []): - rel = relations[rel_id] - if rel.tags.get('type') == 'route' and \ - rel.tags.get('route') in ('subway', 'light_rail', 'tram', - 'train', 'railway'): - all_routes.add(rel_id) - - routes_lookup = defaultdict(set) - for rel_id in all_routes: - rel = relations[rel_id] - route = rel.tags.get('route') - if route: - route_name = mz_transit_route_name(rel.tags) - routes_lookup[route].add(route_name) - trains = routes_lookup['train'] - subways = routes_lookup['subway'] - light_rails = routes_lookup['light_rail'] - trams = routes_lookup['tram'] - railways = routes_lookup['railway'] - del routes_lookup - - # if a station is an interchange between mainline rail and subway or - # light rail, then give it a "bonus" boost of importance. - bonus = 2 if trains and (subways or light_rails) else 1 - - score = (100 * min(9, bonus * len(trains)) + - 10 * min(9, bonus * (len(subways) + len(light_rails))) + - min(9, len(trams) + len(railways))) - - return Transit(score=score, root_relation_id=root_relation_id, - trains=trains, subways=subways, light_rails=light_rails, - railways=railways, trams=trams) + return self._transit_relations.get(rel_id, set()) class DataFetcher(object): @@ -284,13 +121,16 @@ def __init__(self, layers, rows, rels, label_placement_layers): """ Expect layers to be a dict of layer name to LayerInfo. Expect rows to be a list of (fid, shape, properties). Label placement layers should - be a set of layer names for which to generate label placement points. + be a dict of geometry type ('point', 'linestring', 'polygon') to set + of layer names, meaning that each feature of the given type in any of + the named layers should additionally get a generated label point. """ self.layers = layers self.rows = rows self.rels = rels self.label_placement_layers = label_placement_layers + self.osm = OsmFixtureLookup(self.rows, self.rels) def __call__(self, zoom, unpadded_bounds): read_rows = [] @@ -354,92 +194,18 @@ def __call__(self, zoom, unpadded_bounds): # if the feature exists in any label placement layer, then we # should consider generating a centroid label_layers = self.label_placement_layers.get( - _shape_type_lookup(shape), {}) + shape_type_lookup(shape), {}) if layer_name in label_layers: generate_label_placement = True - layer_props = props.copy() - layer_props['min_zoom'] = min_zoom + layer_props = layer_properties( + fid, shape, props, layer_name, zoom, self.osm) - # need to make sure that the name is only applied to one of - # the pois, landuse or buildings layers - in that order of - # priority. - # - # TODO: do this for all name variants & translations - if layer_name in ('pois', 'landuse', 'buildings'): - layer_props.pop('name', None) - - # urgh, hack! - if layer_name == 'water' and shape.geom_type == 'Point': - layer_props['label_placement'] = True - - if shape.geom_type in ('Polygon', 'MultiPolygon'): - layer_props['area'] = shape.area - - if layer_name == 'roads' and \ - shape.geom_type in ('LineString', 'MultiLineString'): - mz_networks = [] - mz_cycling_networks = set() - mz_is_bus_route = False - for rel in rels: - rel_tags = deassoc(rel['tags']) - typ, route, network, ref = [rel_tags.get(k) for k in ( - 'type', 'route', 'network', 'ref')] - if route and (network or ref): - mz_networks.extend([route, network, ref]) - if typ == 'route' and \ - route in ('hiking', 'foot', 'bicycle') and \ - network in ('icn', 'ncn', 'rcn', 'lcn'): - mz_cycling_networks.add(network) - if typ == 'route' and route in ('bus', 'trolleybus'): - mz_is_bus_route = True - - mz_cycling_network = None - for cn in ('icn', 'ncn', 'rcn', 'lcn'): - if layer_props.get(cn) == 'yes' or \ - ('%s_ref' % cn) in layer_props or \ - cn in mz_cycling_networks: - mz_cycling_network = cn - break - - if mz_is_bus_route and \ - zoom >= 12 and \ - layer_props.get('highway') in BUS_ROADS: - layer_props['is_bus_route'] = True - - layer_props['mz_networks'] = mz_networks - if mz_cycling_network: - layer_props['mz_cycling_network'] = mz_cycling_network - - is_poi = layer_name == 'pois' - is_railway_station = props.get('railway') == 'station' - is_point_or_poly = shape.geom_type in ( - 'Point', 'MultiPoint', 'Polygon', 'MultiPolygon') - - if is_poi and is_railway_station and \ - is_point_or_poly and fid >= 0: - node_id = None - way_id = None - if shape.geom_type in ('Point', 'MultiPoint'): - node_id = fid - else: - way_id = fid - - transit = mz_calculate_transit_routes_and_score( - self.rows, self.rels, node_id, way_id) - layer_props['mz_transit_score'] = transit.score - layer_props['mz_transit_root_relation_id'] = ( - transit.root_relation_id) - layer_props['train_routes'] = transit.trains - layer_props['subway_routes'] = transit.subways - layer_props['light_rail_routes'] = transit.light_rails - layer_props['tram_routes'] = transit.trams - - if layer_props: - props_name = '__%s_properties__' % layer_name - read_row[props_name] = layer_props - if layer_name == 'water': - has_water_layer = True + layer_props['min_zoom'] = min_zoom + props_name = '__%s_properties__' % layer_name + read_row[props_name] = layer_props + if layer_name == 'water': + has_water_layer = True # if at least one min_zoom / properties match if read_row: diff --git a/tilequeue/query/rawr.py b/tilequeue/query/rawr.py new file mode 100644 index 00000000..927315ca --- /dev/null +++ b/tilequeue/query/rawr.py @@ -0,0 +1,578 @@ +from collections import namedtuple, defaultdict +from shapely.geometry import box +from shapely.wkb import loads as wkb_loads +from tilequeue.query.common import layer_properties +from tilequeue.query.common import is_station_or_stop +from tilequeue.query.common import is_station_or_line +from tilequeue.query.common import deassoc +from tilequeue.query.common import mz_is_interesting_transit_relation +from tilequeue.query.common import shape_type_lookup +from tilequeue.transform import calculate_padded_bounds +from raw_tiles.tile import shape_tile_coverage +from math import floor + + +class Relation(object): + """ + Relation object holds data about a relation and provides a nicer interface + than the raw tuple by turning the tags array into a dict, and separating + out the "parts" array of IDs into separate lists for nodes, ways and other + relations. + """ + + def __init__(self, rel_id, way_off, rel_off, parts, members, tags): + self.id = rel_id + self.tags = deassoc(tags) + self.node_ids = parts[0:way_off] + self.way_ids = parts[way_off:rel_off] + self.rel_ids = parts[rel_off:] + + +class TilePyramid(namedtuple('TilePyramid', 'z x y max_z')): + """ + Represents a "tile pyramid" of all tiles which are geographically + contained within the tile `z/x/y` up to a maximum zoom of `max_z`. This is + the set of tiles corresponding to one RAWR tile. + """ + + def tile(self): + from raw_tiles.tile import Tile + return Tile(self.z, self.x, self.y) + + def bounds(self): + from ModestMaps.Core import Coordinate + from tilequeue.tile import coord_to_mercator_bounds + + coord = Coordinate(zoom=self.z, column=self.x, row=self.y) + bounds = coord_to_mercator_bounds(coord) + + return bounds + + def bbox(self): + return box(*self.bounds()) + + +# weak enum type used to represent shape type, rather than use a string. +# (or, where we have to use a string, at least use a shared single instance +# of a string) +class ShapeType(object): + point = 1 + line = 2 + polygon = 3 + + _LOOKUP = ['point', 'line', 'polygon'] + + @staticmethod + def lookup(typ): + "turn the enum into a string" + return ShapeType._LOOKUP[typ-1] + + +# determine the shape type from the raw WKB bytes. this means we don't have to +# parse the WKB, which can be an expensive operation for large polygons. +def _wkb_shape(wkb): + reverse = ord(wkb[0]) == 1 + type_bytes = map(ord, wkb[1:5]) + if reverse: + type_bytes.reverse() + typ = type_bytes[3] + if typ == 1 or typ == 4: + return ShapeType.point + elif typ == 2 or typ == 5: + return ShapeType.line + elif typ == 3 or typ == 6: + return ShapeType.polygon + else: + assert False, "WKB shape type %d not understood." % (typ,) + + +# return true if the tuple of values corresponds to, and each is an instance +# of, the tuple of types. this is used to make sure that argument lists are +# the right "shape" before destructuring (splatting?) them in a function call. +def _match_type(values, types): + if len(values) != len(types): + return False + for val, typ in zip(values, types): + if not isinstance(val, typ): + return False + return True + + +# return true if the tags indicate that this is a gate +def _is_gate(props): + return props.get('barrier') == 'gate' + + +# return true if the tags indicate that this is a highway, cycleway or footway +# which might be part of a route relation. note that this is pretty loose, and +# might return true for things we don't eventually render as roads, but is just +# aimed at cutting down the number of items we need in our index. +def _is_routeable(props): + return props.get('whitewater') == 'portage_way' or 'highway' in props + + +class OsmRawrLookup(object): + """ + Implements the interface needed by the common code (e.g: layer_properties) + to look up information about node, way and relation IDs. For database + lookups, we previously did this with a JOIN, and the fixture data source + just iterates over the (small) number of items. + + For RAWR tiles, we index the data to provide faster lookup, and are more + selective about what goes into the index. + """ + + def __init__(self): + self.nodes = {} + self.ways = {} + self.relations = {} + + self._ways_using_node = defaultdict(list) + self._relations_using_node = defaultdict(list) + self._relations_using_way = defaultdict(list) + self._relations_using_rel = defaultdict(list) + + def add_row(self, *args): + # there's only a single dispatch from the indexing function, which + # passes row data from the table. we have to figure out here what + # kind of row it was, and send the data on to the right function. + + # IDs can be either ints or longs, and generally we don't care which, + # so we accept either as the type for that position in the function. + num = (int, long) + + if _match_type(args, (num, (str, bytes), dict)): + self.add_feature(*args) + + elif _match_type(args, (num, list, list)): + self.add_way(*args) + + elif _match_type(args, (num, num, num, list, list, list)): + self.add_relation(*args) + + else: + raise Exception("Unknown row shape for OsmRawrLookup.add_row: %s" % + (repr(map(type, args)),)) + + def add_feature(self, fid, shape_wkb, props): + if fid < 0: + return + + shape_type = _wkb_shape(shape_wkb) + if is_station_or_stop(fid, None, props) and \ + shape_type == ShapeType.point: + # must be a station or stop node + self.nodes[fid] = (fid, shape_wkb, props) + + elif _is_gate(props) and shape_type == ShapeType.point: + # index the highways that use gates to influence min zoom + self.nodes[fid] = (fid, shape_wkb, props) + + elif (is_station_or_line(fid, None, props) and + shape_type != ShapeType.point): + # must be a station polygon or stop line + self.ways[fid] = (fid, shape_wkb, props) + + elif _is_routeable(props) and shape_type == ShapeType.line: + # index routable items (highways, cycleways, footpaths) to + # get the relations using them. + self.ways[fid] = (fid, shape_wkb, props) + + def add_way(self, way_id, nodes, tags): + for node_id in nodes: + if node_id in self.nodes: + assert way_id in self.ways + self._ways_using_node[node_id].append(way_id) + + def add_relation(self, rel_id, way_off, rel_off, parts, members, tags): + r = Relation(rel_id, way_off, rel_off, parts, members, tags) + is_transit_relation = mz_is_interesting_transit_relation(r.tags) + is_route = 'route' in r.tags and \ + ('network' in r.tags or 'ref' in r.tags) + if is_route or is_transit_relation: + self.relations[r.id] = r + for node_id in r.node_ids: + if node_id in self.nodes: + self._relations_using_node[node_id].append(rel_id) + for way_id in r.way_ids: + if way_id in self.ways: + self._relations_using_way[way_id].append(rel_id) + for member_rel_id in r.rel_ids: + self._relations_using_rel[member_rel_id].append(rel_id) + + def relations_using_node(self, node_id): + "Returns a list of relation IDs which contain the node with that ID." + + return self._relations_using_node.get(node_id, []) + + def relations_using_way(self, way_id): + "Returns a list of relation IDs which contain the way with that ID." + + return self._relations_using_way.get(way_id, []) + + def relations_using_rel(self, rel_id): + """ + Returns a list of relation IDs which contain the relation with that + ID. + """ + + return self._relations_using_rel.get(rel_id, []) + + def ways_using_node(self, node_id): + "Returns a list of way IDs which contain the node with that ID." + + return self._ways_using_node.get(node_id, []) + + def relation(self, rel_id): + "Returns the Relation object with the given ID." + + return self.relations[rel_id] + + def way(self, way_id): + """ + Returns the feature (fid, shape, props) which was generated from the + given way. + """ + + return self.ways[way_id] + + def node(self, node_id): + """ + Returns the feature (fid, shape, props) which was generated from the + given node. + """ + + return self.nodes[node_id] + + def transit_relations(self, rel_id): + "Return transit relations containing the relation with the given ID." + + return set(self.relations_using_rel(rel_id)) + + +# yield all the tiles at the given zoom level which intersect the given bounds. +def _tiles(zoom, unpadded_bounds): + from tilequeue.tile import mercator_point_to_coord + from raw_tiles.tile import Tile + + minx, miny, maxx, maxy = unpadded_bounds + topleft = mercator_point_to_coord(zoom, minx, miny) + bottomright = mercator_point_to_coord(zoom, maxx, maxy) + + # make sure that the bottom right coordinate is below and to the right + # of the top left coordinate. it can happen that the coordinates are + # mixed up due to small numerical precision artefacts being enlarged + # by the conversion to integer and y-coordinate flip. + assert topleft.zoom == bottomright.zoom + bottomright.column = max(bottomright.column, topleft.column) + bottomright.row = max(bottomright.row, topleft.row) + + for x in range(int(topleft.column), int(bottomright.column) + 1): + for y in range(int(topleft.row), int(bottomright.row) + 1): + tile = Tile(zoom, x, y) + yield tile + + +# the object which gets indexed. this is a normal (fid, shape, props) tuple +# expanded to include a dict of layer name to min zoom in `layer_min_zooms`. +# this means that the properties don't have to be copied and altered to +# include the min zoom for each layer, reducing the memory footprint. +_Feature = namedtuple('_Feature', 'fid shape properties layer_min_zooms') + + +class _LazyShape(object): + """ + This proxy exists so that we can avoid parsing the WKB for a shape unless + it is actually needed. Parsing WKB is pretty fast, but multiplied over + many thousands of objects, it can become the slowest part of the indexing + process. Given that we reject many features on the basis of their + properties alone, lazily parsing the WKB can provide a significant saving. + """ + + def __init__(self, wkb): + self.wkb = wkb + self.obj = None + + def __getattr__(self, name): + if self.obj is None: + self.obj = wkb_loads(self.wkb) + return getattr(self.obj, name) + + +_Metadata = namedtuple('_Metadata', 'source ways relations') + + +def _make_meta(source, fid, shape_type, osm): + ways = [] + rels = [] + + # fetch ways and relations for any node + if fid >= 0 and shape_type == ShapeType.point: + for way_id in osm.ways_using_node(fid): + ways.append(osm.way(way_id)) + for rel_id in osm.relations_using_node(fid): + rels.append(osm.relation(rel_id)) + + # and relations for any way + if fid >= 0 and shape_type == ShapeType.line: + for rel_id in osm.relations_using_way(fid): + rels.append(osm.relation(rel_id)) + + # have to transform the Relation object into a dict, which is + # what the functions called on this data expect. + # TODO: reusing the Relation object would be better. + rel_dicts = [] + for r in rels: + tags = [] + for k, v in r.tags.items(): + tags.append(k) + tags.append(v) + rel_dicts.append(dict(tags=tags)) + + return _Metadata(source, ways, rel_dicts) + + +class _LayersIndex(object): + """ + Index features by the tile(s) that they appear in. + + This is done by calculating a min-min-zoom, the lowest min_zoom for that + feature across all layers, and then adding that feature to a list for each + tile it appears in from the min-min-zoom up to the max zoom for the tile + pyramid. + """ + + def __init__(self, layers, tile_pyramid, source): + self.layers = layers + self.tile_pyramid = tile_pyramid + self.tile_index = defaultdict(list) + self.source = source + self.delayed_features = [] + + def add_row(self, fid, shape_wkb, props): + shape = _LazyShape(shape_wkb) + # single object (hence single id()) will be shared amongst all layers. + # this allows us to easily and quickly de-duplicate at later layers in + # the stack. + feature = _Feature(fid, shape, props, {}) + + # delay min zoom calculation in order to collect more information about + # the ways and relations using a particular feature. + self.delayed_features.append(feature) + + def index(self, osm): + for feature in self.delayed_features: + self._index_feature(feature, osm) + del self.delayed_features + + def _index_feature(self, feature, osm): + fid = feature.fid + shape = feature.shape + props = feature.properties + layer_min_zooms = feature.layer_min_zooms + + # grab the shape type without decoding the WKB to save time. + shape_type = _wkb_shape(shape.wkb) + + meta = _make_meta(self.source, fid, shape_type, osm) + for layer_name, info in self.layers.items(): + shape_type_str = ShapeType.lookup(shape_type) + if info.shape_types and shape_type_str not in info.shape_types: + continue + min_zoom = info.min_zoom_fn(shape, props, fid, meta) + if min_zoom is not None: + layer_min_zooms[layer_name] = min_zoom + + # quick exit if the feature didn't have a min zoom in any layer. + if not layer_min_zooms: + return + + # lowest zoom that this feature appears in any layer. note that this + # is clamped to the max zoom, so that all features that appear at some + # zoom level appear at the max zoom. this is different from the min + # zoom in layer_min_zooms, which is a property that will be injected + # for each layer and is used by the _client_ to determine feature + # visibility. + min_zoom = min(self.tile_pyramid.max_z, min(layer_min_zooms.values())) + + # take the minimum integer zoom - this is the min zoom tile that the + # feature should appear in, and a feature with min_zoom = 1.9 should + # appear in a tile at z=1, not 2, since the tile at z=N is used for + # the zoom range N to N+1. + # + # we cut this off at this index's min zoom, as we aren't interested + # in any tiles outside of that. + floor_zoom = max(self.tile_pyramid.z, int(floor(min_zoom))) + + # seed initial set of tiles at maximum zoom. all features appear at + # least at the max zoom, even if the min_zoom function returns a + # value larger than the max zoom. + zoom = self.tile_pyramid.max_z + tiles = shape_tile_coverage(shape, zoom, self.tile_pyramid.tile()) + + while zoom >= floor_zoom: + parent_tiles = set() + for tile in tiles: + self.tile_index[tile].append(feature) + parent_tiles.add(tile.parent()) + + zoom -= 1 + tiles = parent_tiles + + def __call__(self, tile): + return self.tile_index.get(tile, []) + + +class DataFetcher(object): + + def __init__(self, layers, tables, tile_pyramid, label_placement_layers, + source): + """ + Expect layers to be a dict of layer name to LayerInfo (see fixture.py). + Tables should be a callable which returns a generator over the rows in + the table when called with that table's name. + """ + + from raw_tiles.index.index import index_table + + self.layers = layers + self.tile_pyramid = tile_pyramid + self.label_placement_layers = label_placement_layers + self.source = source + self.layer_indexes = {} + + table_indexes = defaultdict(list) + + self.layers_index = _LayersIndex( + self.layers, self.tile_pyramid, self.source) + for shape_type in ('point', 'line', 'polygon'): + table_name = 'planet_osm_' + shape_type + table_indexes[table_name].append(self.layers_index) + + self.osm = OsmRawrLookup() + # NOTE: order here is different from that in raw_tiles index() + # function. this is because here we want to gather up some + # "interesting" feature IDs before we look at the ways/rels tables. + for typ in ('point', 'line', 'polygon', 'ways', 'rels'): + table_name = 'planet_osm_' + typ + source = tables(table_name) + extra_indexes = table_indexes[table_name] + index_table(source, self.osm, *extra_indexes) + + # there's a chicken and egg problem with the indexes: we want to know + # which features to index, but also calculate the feature's min zoom, + # which might depend on ways and relations not seen yet. one solution + # would be to do this in two passes, but that might mean paying a cost + # to decompress or deserialize the data twice. instead, the index + # buffers the features and indexes them in the following step. this + # might mean we buffer more information in memory than we technically + # need if many of the features are not visible, but means we get one + # single set of _Feature objects. + self.layers_index.index(self.osm) + + def _named_layer(self, layer_min_zooms): + # we want only one layer from ('pois', 'landuse', 'buildings') for + # each feature to be assigned a name. therefore, we use the presence + # or absence of a min zoom to check whether these features as in these + # layers, and therefore which should be assigned the name. handily, + # the min zooms are already pre-calculated as layer_min_zooms from the + # index. + for layer_name in ('pois', 'landuse', 'buildings'): + if layer_name in layer_min_zooms: + return layer_name + return None + + def _lookup(self, zoom, unpadded_bounds): + features = [] + seen_ids = set() + + for tile in _tiles(zoom, unpadded_bounds): + tile_features = self.layers_index(tile) + for feature in tile_features: + feature_id = id(feature) + if feature_id not in seen_ids: + seen_ids.add(feature_id) + features.append(feature) + + return features + + def __call__(self, zoom, unpadded_bounds): + read_rows = [] + bbox = box(*unpadded_bounds) + + # check that the call is fetching data which is actually within the + # bounds of the tile pyramid. we don't have data outside of that, so + # can't fulfil requests. if these assertions are tripping, it probably + # indicates a programming error - has the wrong DataFetcher been + # loaded? + assert zoom <= self.tile_pyramid.max_z + assert zoom >= self.tile_pyramid.z + assert bbox.within(self.tile_pyramid.bbox()) + + for (fid, shape, props, layer_min_zooms) in self._lookup( + zoom, unpadded_bounds): + # reject any feature which doesn't intersect the given bounds + if bbox.disjoint(shape): + continue + + # place for assembing the read row as if from postgres + read_row = {} + generate_label_placement = False + + # add name into whichever of the pois, landuse or buildings + # layers has claimed this feature. + name = props.get('name', None) + named_layer = self._named_layer(layer_min_zooms) + + for layer_name, min_zoom in layer_min_zooms.items(): + # we need to keep fractional zooms, e.g: 4.999 should appear + # in tiles at zoom level 4, but not 3. also, tiles at zooms + # past the max zoom should be clamped to the max zoom. + tile_zoom = min(self.tile_pyramid.max_z, floor(min_zoom)) + if tile_zoom > zoom: + continue + + layer_props = layer_properties( + fid, shape, props, layer_name, zoom, self.osm) + layer_props['min_zoom'] = min_zoom + + if name and named_layer == layer_name: + layer_props['name'] = name + + read_row['__' + layer_name + '_properties__'] = layer_props + + # if the feature exists in any label placement layer, then we + # should consider generating a centroid + label_layers = self.label_placement_layers.get( + shape_type_lookup(shape), {}) + if layer_name in label_layers: + generate_label_placement = True + + if read_row: + read_row['__id__'] = fid + + # if this is a water layer feature, then clip to an expanded + # bounding box to avoid tile-edge artefacts. + clip_box = bbox + if layer_name == 'water': + pad_factor = 1.1 + clip_box = calculate_padded_bounds( + pad_factor, unpadded_bounds) + clip_shape = clip_box.intersection(shape) + read_row['__geometry__'] = bytes(clip_shape.wkb) + + if generate_label_placement: + read_row['__label__'] = bytes( + shape.representative_point().wkb) + + read_rows.append(read_row) + + return read_rows + + +# tables is a callable which should return a generator over the rows of the +# table when called with the table name. +def make_rawr_data_fetcher(layers, tables, tile_pyramid, source, + label_placement_layers={}): + return DataFetcher(layers, tables, tile_pyramid, label_placement_layers, + source)