Add path coverage options and derive min/max zoom levels

docs/Manifest.rst

@@ -26,6 +26,12 @@ It MAY optionally have:
     identifiers to features. These are additional to any assigned by the ``properties`` file.
 *   for functional connectivity maps, an ``"annotation"`` JSON file assigning anatomical terms to
     features based on their label and anatomical type (System, Organ, FTU).
+*   ``"path-zoom-range"`` to enable path-based minzoom/maxzoom calculation (default ``true``).
+*   ``"initial-zoom"`` for the initial zoom level (default ``4``).
+*   ``"max-zoom"`` for the maximum zoom level (default ``10``).
+*   ``"max-raster-zoom"`` for the maximum zoom level for raster tiles (defaults to ``"max-zoom"``).
+*   ``"path-min-coverage"`` for path zoom range calculation lower coverage bound (default ``0.7``).
+*   ``"path-max-coverage"`` for path zoom range calculation upper coverage bound (default ``5.0``).
 *   a ``"connectivityTerms"`` JSON file specifying equvalences between historical anatomical terms
     used in SCKAN to standard terms (e.g. between FMA and ILX identifiers). **DEPRECATED**
 *   a ``"connectivity"`` JSON file specifying manually defined neuron paths. **DEPRECATED**

mapmaker/__main__.py

@@ -84,14 +84,6 @@ def arg_parser():
     debug_options.add_argument('--tippecanoe', dest='showTippe', action='store_true',
                         help='Show command used to run Tippecanoe')
 
-    zoom_options = parser.add_argument_group('Zoom level')
-    zoom_options.add_argument('--initial-zoom', dest='initialZoom', metavar='N', type=int, default=4,
-                        help='Initial zoom level (defaults to 4)')
-    zoom_options.add_argument('--max-zoom', dest='maxZoom', metavar='N', type=int, default=10,
-                        help='Maximum zoom level (defaults to 10)')
-    zoom_options.add_argument('--max-raster-zoom', dest='maxRasterZoom', metavar='N', type=int,
-                        help='Maximum zoom level of rasterised tiles (defaults to maximum zoom level)')
-
     misc_options = parser.add_argument_group('Miscellaneous')
     misc_options.add_argument('--commit', metavar='GIT_COMMIT',
                         help='The branch/tag/commit to use when the source is a Git repository')

mapmaker/flatmap/manifest.py

@@ -378,6 +378,30 @@ def connectivity_terms(self):
     def description(self):
         return self.__manifest.get('description')
 
+    @property
+    def path_zoom_range(self) -> bool:
+        return bool(self.__manifest.get('path-zoom-range', True))
+
+    @property
+    def initial_zoom(self):
+        return self.__manifest.get('initial-zoom', 4)
+
+    @property
+    def max_zoom(self):
+        return self.__manifest.get('max-zoom', 10)
+
+    @property
+    def max_raster_zoom(self):
+        return self.__manifest.get('max-raster-zoom', self.max_zoom)
+
+    @property
+    def path_min_coverage(self):
+        return self.__manifest.get('path-min-coverage', 0.7)
+
+    @property
+    def path_max_coverage(self):
+        return self.__manifest.get('path-max-coverage', 5.0)
+
     @property
     def exported_properties(self) -> list[str]:
         return self.__manifest.get('exported-properties', [])

mapmaker/maker.py

@@ -109,20 +109,6 @@ def __init__(self, options: dict[str, Any], logger_port: Optional[int]=None,
         if 'output' not in options:
             options['output'] = './flatmaps'
 
-        # Check zoom settings are valid
-        min_zoom = 0
-        max_zoom = options.get('maxZoom', 10)
-        max_raster_zoom = options.get('maxRasterZoom', max_zoom)
-
-        initial_zoom = options.get('initialZoom', 4)
-        if max_zoom < min_zoom or max_zoom > 15:
-            raise ValueError('Max zoom must be between {} and 15'.format(min_zoom))
-        if max_raster_zoom > max_zoom:
-            raise ValueError(f'Max raster zoom cannot be greater than max zoom ({max_zoom})')
-        if initial_zoom < min_zoom or initial_zoom > max_zoom:
-            raise ValueError(f'Initial zoom cannot be greater than max zoom ({max_zoom})')
-        self.__zoom = (min_zoom, max_zoom, initial_zoom)
-
         if options.get('publish'):
             # Check the given options are compatible with SDS publishing
             errors = False
@@ -153,6 +139,28 @@ def __init__(self, options: dict[str, Any], logger_port: Optional[int]=None,
         if self.__id is None:
             raise ValueError('No id given for map')
 
+        # Check zoom settings are valid
+        min_zoom = 0
+        max_zoom = self.__manifest.max_zoom
+        max_raster_zoom = self.__manifest.max_raster_zoom
+
+        initial_zoom = self.__manifest.initial_zoom
+        if max_zoom < min_zoom or max_zoom > 15:
+            raise ValueError('Max zoom must be between {} and 15'.format(min_zoom))
+        if max_raster_zoom > max_zoom:
+            raise ValueError(f'Max raster zoom cannot be greater than max zoom ({max_zoom})')
+        if initial_zoom < min_zoom or initial_zoom > max_zoom:
+            raise ValueError(f'Initial zoom cannot be greater than max zoom ({max_zoom})')
+
+        path_min_coverage = self.__manifest.path_min_coverage
+        path_max_coverage = self.__manifest.path_max_coverage
+        if path_min_coverage <= 0 or path_max_coverage <= 0:
+            raise ValueError('Path coverage values must be greater than 0')
+        if path_max_coverage < path_min_coverage:
+            raise ValueError('Path max coverage cannot be less than path min coverage')
+
+        self.__zoom = (min_zoom, max_zoom, initial_zoom)
+
         # Publishing requires a ``description.json``
         if options.get('publish') and self.__manifest.description is None:
             raise ValueError('The manifest must specify a JSON `description` file if publishing')

mapmaker/output/geojson.py

@@ -39,6 +39,11 @@
 
 #===============================================================================
 
+# Earth circumference (WGS84 equatorial radius for Web Mercator EPSG:3857)
+EARTH_CIRCUMFERENCE = 2 * math.pi * 6378137.0  # meters
+
+#===============================================================================
+
 class GeoJSONOutput(object):
     def __init__(self, flatmap: FlatMap, layer: MapLayer, output_dir: str):
     #======================================================================
@@ -120,6 +125,11 @@ def __save_features(self, features):
                 if scale > 6 and 'group' not in properties and 'minzoom' not in properties:
                     geojson['tippecanoe']['minzoom'] = 4
             else:
+                if (self.__flatmap.manifest.path_zoom_range
+                and (nodes:=self.__flatmap.connectivity()['paths'].get(feature.models))):
+                    zoom_range = self.__get_path_zoom_range(nodes)
+                    geojson['properties']['minzoom'] = zoom_range[0]
+                    geojson['properties']['maxzoom'] = zoom_range[1]
                 geojson['properties']['scale'] = 10
             geojson['properties'].update(properties)
 
@@ -158,3 +168,43 @@ def __save_features(self, features):
             progress_bar.update(1)
 
         progress_bar.close()
+
+    def __get_path_zoom_range(self, nodes):
+        path_min_coverage = self.__flatmap.manifest.path_min_coverage
+        path_max_coverage = self.__flatmap.manifest.path_max_coverage
+
+        points = [geom.centroid
+                    for node in nodes.get('nodes', [])
+                        if (f := self.__flatmap.get_feature_by_geojson_id(node)) is not None
+                            and (geom := f.geometry) is not None]
+
+        if len(points) > 1:
+            # extent
+            xs = [p.x for p in points]
+            ys = [p.y for p in points]
+            # World-space extents (meters, EPSG:3857)
+            x_extent = max(xs) - min(xs)
+            y_extent = max(ys) - min(ys)
+            max_dist = 0.0
+            for i, p1 in enumerate(points):
+                for p2 in points[i+1:]:
+                    dist = math.hypot(p2.x - p1.x, p2.y - p1.y)
+                    if dist > max_dist:
+                        max_dist = dist
+            extent = max(x_extent, y_extent, max_dist)
+            if not math.isfinite(extent) or extent <= 0:
+                return self.__flatmap.min_zoom, self.__flatmap.max_zoom
+
+            # minzoom and maxzoom
+            z_min = math.ceil(math.log2((path_min_coverage * EARTH_CIRCUMFERENCE) / extent))
+            z_max = math.floor(math.log2((path_max_coverage * EARTH_CIRCUMFERENCE) / extent))
+            if z_max < z_min:
+                z_max = z_min
+            z_min = max(self.__flatmap.min_zoom, z_min)
+            z_max = min(self.__flatmap.max_zoom, z_max)
+            if z_min >= self.__flatmap.max_zoom:
+                z_min = self.__flatmap.max_zoom - 1
+                z_max = self.__flatmap.max_zoom
+            return z_min, z_max
+
+        return self.__flatmap.min_zoom, self.__flatmap.max_zoom