plotly automatic zooming for "Mapbox maps"

Question

In the plotly website Map Configuration and Styling in Python is described how to automatically zoom a "Geo map":

import plotly.express as px

fig = px.line_geo(lat=[0,15,20,35], lon=[5,10,25,30])            # Creates a "Geo map" figure
fig.update_geos(fitbounds="locations")                           # Automatic Zooming !!!!
fig.show()

and this works, moreover if I try to the same on a "Mapbox map" it does not apply auto zoom:

fig = px.scatter_mapbox(filtered_df, lat="latitude", lon="longitude", color="ID")  # Creates a "Mapbox map" figure
fig.update_layout(mapbox_style="open-street-map")
fig.update_geos(fitbounds="locations")                                             # Automatic Zooming not working!!!

There is not information of how to do it in the Mapbox Map Layers in Python.

Answer 1

Given that I had issues installing planar, I modified one of the answers to only use numpy. Works quite well, and changing the fp values if needed can lead to optimal results.

def get_map_zoom_level(lon=None, lat=None):
"""
Calculate appropriate zoom level for scattermap plotting.
Returns the appropriate zoom-level for plotly-map-graphics.

Args:
    longitudes (list, optional): List of longitude values
    latitudes (list, optional): List of latitude values

Returns:
    float: zoom_level
"""
if (lon is None or lat is None) or (len(lon) != len(lat)):
    return 0

lonlat_pairs = list(zip(lon, lat))

lons, lats = zip(*lonlat_pairs)
min_lon, max_lon = min(lons), max(lons)
min_lat, max_lat = min(lats), max(lats)

width = max_lon - min_lon
height = max_lat - min_lat

# Calculate area (using width * height as a simple approximation)
area = width * height

# 1D-linear interpolation with numpy for zoom level
zoom = np.interp(
    x=area,
    xp=[0, 5**-10, 4**-10, 3**-10, 2**-10, 1**-10, 1**-5],
    fp=[19, 18, 17, 15.5, 9, 6, 4],
)

return zoom

Answer 2

Looking at mapbox's documentation for zoom levels, zoom determines how many tiles mapbox will use to represent the world. world = 2^zoom x 2^zoom tiles.

Given the world is about ~20,000km from pole to pole and ~40,000km across the equator (radius of 6371km * 2 * pi = 40,036km), so at...

• zoom=0, the world will be 1 tile that is 40_000km wide and 20_000km tall.
• zoom=1, 2x2 tiles. Each tile will be 20_000km wide, 10_000km tall
• zoom=3, 8x8 tiles. Each tile will be 40_000 / 8 km wide, and 20_000 / 8 km tall
• and so on...

So we're trying to figure out at what zoom level 1 tile is closest to the width or height.

For charts that are wide, we're solving for log2(40_000 / x-distance). Except the circumference is 40_000 at the equator but shrinks as you near the poles. So instead we have log2(40_000 * cosine(latitude) / x-distance).

For charts that are tall, we're solving for log2(20_000 / y-distance).

So I ended up with something like this (I should probably test this more before posting):

import numpy as np
from geopy import distance

def calculate_mapbox_zoom_center(lons, lats) -> (float, dict):
    """calc zoom and center for plotly mapbox functions

    Temporary solution awaiting official implementation, see:
    https://github.com/plotly/plotly.js/issues/3434
    """

    # because we stick None to make separate lines...
    min_lat = min([l for l in lats if l])
    max_lat = max([l for l in lats if l])
    min_lon = min([l for l in lons if l])
    max_lon = max([l for l in lons if l])
    y_range = distance.geodesic((min_lat, min_lon), (max_lat, min_lon)).kilometers
    x_range = distance.geodesic((min_lat, min_lon), (min_lat, max_lon)).kilometers

    # does this work across the international date line?
    center = dict(lat=(min_lat + max_lat) / 2, lon=(min_lon + max_lon) / 2)

    # figure out zoom
    bound_by_y = y_range > x_range
    if bound_by_y:
        # in mercator the world's height is ~20,000 km
        zoom = np.log2(20_000 / y_range)
    else:
        # in mercator the world's "width" is ~40,000km at the equator but shrinks as you near the poles
        circumference_at_latitude = 40_000 * np.cos(np.abs(center['lat']) * np.pi / 180)
        zoom = np.log2(circumference_at_latitude / x_range)

    return zoom, center

Answer 3

These solutions only kinda work, the problem is that the pixel size changes with the latitude. So if you want a general solution (and cannot talk to mapbox), I used opencv to find the max and min, lat and lon points, on a "test image", then calculate the pixel length in the y axis and the pixel length in the x axis, from max zoom (sometimes 22, or I skip to 20) to min zoom 0.

This is the sample "test image" generated with the max and min lat and long points in pink and orange:

This is an image of the annotated map:

    def get_marker_location(self, color_min, color_max, hsv):
        # cv2.imshow("hsv", hsv)
        # cv2.waitKey(0)
        mask = cv2.inRange(hsv, color_min, color_max)
        # cv2.imshow("Image", mask)
        # cv2.waitKey(0)
        if cv2.countNonZero(mask) == 0:
            return False, None, None
        cnts, hierarchies = cv2.findContours(mask, cv2.RETR_EXTERNAL, cv2.CHAIN_APPROX_SIMPLE)

        for i in cnts:
            M = cv2.moments(i)
            if M['m00'] != 0:
                cX = int(M["m10"] / M["m00"])
                cY = int(M["m01"] / M["m00"])

                return True, cX, cY


    def get_zoom(self, lons, lats, height, width):
        fig = go.Figure()

        # Get max min long lat per owner
        maxlon, minlon = np.amax(lons), np.amin(lons)
        maxlat, minlat = np.amax(lats), np.amin(lats)

        center = {
            'lon': round((maxlon + minlon) / 2, 6),
            'lat': round((maxlat + minlat) / 2, 6)
        }

        # Place each on map in different colors
        # Pink ff00ff hsv(300,100,100)
        pink_min = (147,219,208)
        pink_max = (153,255,255)
        fig.add_trace(go.Scattermapbox(lon=[minlon], lat=[minlat], 
                                        mode='markers', 
                                        marker=dict(color='#FF00FF')))
        # Oraange ff6000 hsv(23,100,100)
        orange_min = np.array([10,230,236])
        orange_max = np.array([12,255,255])
        fig.add_trace(go.Scattermapbox(lon=[maxlon], lat=[maxlat], 
                                        mode='markers', 
                                        marker=dict(color='#FF6000')))

        # Change start number for small area, skip if longer
        start_zoom = 20
        # extra_zoom_threshold = 0.00762 # 25ft in km
        extra_zoom_threshold = 0.01524 # 50ft in km
        max_min_dist = gpd.geodesic(gp.Point(minlat, minlon),gp.Point(maxlat, maxlon)).km
        print(f'max_min_dist: {max_min_dist}')
        print(f'extra_zoom_threshold: {extra_zoom_threshold}')
        if extra_zoom_threshold > max_min_dist:
            start_zoom = 22

        for i in range(start_zoom, -1, -1):
            # Make images
            fig.update_layout(template='plotly_white', 
                                showlegend=False, height=height, width=width, 
                                margin={'l': 0,'r': 0,'b': 0,'t': 0,},
                                mapbox={'style':'carto-positron', 
                                        'zoom':i,
                                        'center': {'lon': center['lon'], 'lat': center['lat']}})
            fig.write_image("check_zoom.png")

            # Check for colors
            img = cv2.imread('check_zoom.png')
            hsv = cv2.cvtColor(img, cv2.COLOR_BGR2HSV)
            pink_loc = dict()
            orange_loc = dict()
            # min
            pink_loc['in_image'], pink_loc['cX'], pink_loc['cY'] = self.get_marker_location(pink_min, pink_max, hsv)
            # max
            orange_loc['in_image'], orange_loc['cX'], orange_loc['cY'] = self.get_marker_location(orange_min, orange_max, hsv)
            # x = 1/0

            if (pink_loc['in_image'] and orange_loc['in_image']) or i == 0:
                # Compute Long and Lat pixel lengths
                # Calc lon lat distance / number of pixels
                pixel_len_lonx = gpd.geodesic(gp.Point(center['lat'], minlon),gp.Point(center['lat'], maxlon)).km / abs(orange_loc['cX'] - pink_loc['cX'])
                pixel_len_laty = gpd.geodesic(gp.Point(minlat, center['lon']),gp.Point(maxlat, center['lon'])).km / abs(orange_loc['cY'] - pink_loc['cY'])
                
                dy = (height/2)*pixel_len_laty # km
                dx = (width/2)*pixel_len_lonx # km
                
                max_lat = gpd.geodesic(kilometers=dy).destination(gp.Point(center['lat'], center['lon']), 0)[0]
                min_lat = gpd.geodesic(kilometers=dy).destination(gp.Point(center['lat'], center['lon']), 180)[0]
                max_lon = gpd.geodesic(kilometers=dx).destination(gp.Point(center['lat'], center['lon']), 90)[1]
                min_lon = gpd.geodesic(kilometers=dx).destination(gp.Point(center['lat'], center['lon']), 270)[1]

                return i, center, max_lat, min_lat, max_lon, min_lon

    def get_mapbox_site_map(self):
        # generate site map
        fig = go.Figure()

        height = 800
        width = 800
        annotation_counter = 0

        if len(self.lines_df):
            l_lat = np.concatenate(self.lines_df['lines'].apply(lambda g: [c[0] for c in g] + [None]).values)
            l_lon = np.concatenate(self.lines_df['lines'].apply(lambda g: [c[1] for c in g] + [None]).values)
        if len(self.gates_df):
            g_lat = np.concatenate(self.gates_df['gates'].apply(lambda g: [c[0] for c in g] + [None]).values)
            g_lon = np.concatenate(self.gates_df['gates'].apply(lambda g: [c[1] for c in g] + [None]).values)
            
        if len(self.gates_df) and len(self.lines_df):
            lg_lat = np.concatenate((l_lat, g_lat))
            lg_lon = np.concatenate((l_lon, g_lon))

            # Get zoom level => Get lat pixel distance, get long pixel distance => Annotate
            zoom, center, max_lat, min_lat, max_lon, min_lon = self.get_zoom(lg_lon[lg_lon != None], lg_lat[lg_lat != None], height, width)
        elif len(self.lines_df):
            # Get zoom level => Get lat pixel distance, get long pixel distance => Annotate
            zoom, center, max_lat, min_lat, max_lon, min_lon = self.get_zoom(l_lon[l_lon != None], l_lat[l_lat != None], height, width)
        elif len(self.gates_df):
            # Get zoom level => Get lat pixel distance, get long pixel distance => Annotate
            zoom, center, max_lat, min_lat, max_lon, min_lon = self.get_zoom(g_lon[g_lat != None], g_lat[g_lat != None], height, width)

        for index in self.lines_df.index:
            line = self.lines_df.at[index, 'lines']

            # 1 is lng x, 0 is lat y
            # x.append(line[0][1])
            # y.append(line[0][0])

            pt1 = line[0]
            pt2 = line[1]
            x_annotation = min(pt1[1], pt2[1]) + (abs(pt2[1] - pt1[1]) / 2)
            y_annotation = min(pt1[0], pt2[0]) + (abs(pt2[0] - pt1[0]) / 2)
            x_annotation = (x_annotation - min_lon)/(max_lon - min_lon)
            y_annotation = (y_annotation - min_lat)/(max_lat - min_lat)
            # print(f'zoom {zoom} x {x_annotation} y {y_annotation}')

            if 2 * (abs(pt2[1] - pt1[1]) / 2) <= (abs(pt2[0] - pt1[0]) / 2):
                xanchor = 'left'
            else:
                xanchor = 'center'

            str_i = str(index).replace('Line ', 'L')
            fig.add_annotation(x=x_annotation,
                                xanchor=xanchor,
                                y=y_annotation,
                                xref='paper',yref='paper',
                                yanchor='bottom',
                                text=f'<b>{str_i}: {str(self.lines_df.at[index, "distance"])}ft</b>',
                                font=dict(size=15, color="darkgreen"),
                                showarrow=False)

            annotation_counter += 1

        if len(self.lines_df):
            fig.add_trace(go.Scattermapbox(lon=l_lon, lat=l_lat, 
                                            mode='markers+lines', 
                                            marker=dict(color='#00B900'), 
                                            line=dict(color='#00B900')))

        annotation_counter = 0

        for index in self.gates_df.index:
            gate = self.gates_df.at[index, 'gates']

            # 1 is lng, 0 is lat
            # x.append(gate[0][1])
            # y.append(gate[0][0])

            pt1 = gate[0]
            pt2 = gate[1]
            x_annotation = min(pt1[1], pt2[1]) + (abs(pt2[1] - pt1[1]) / 2)
            y_annotation = min(pt1[0], pt2[0]) + (abs(pt2[0] - pt1[0]) / 2)
            x_annotation = (x_annotation - min_lon)/(max_lon - min_lon)
            y_annotation = (y_annotation - min_lat)/(max_lat - min_lat)

            if 2 * (abs(pt2[1] - pt1[1]) / 2) <= (abs(pt2[0] - pt1[0]) / 2):
                xanchor = 'left'
            else:
                xanchor = 'center'

            str_i = str(index).replace('Gate ', 'G')
            fig.add_annotation(x=x_annotation,
                                xanchor=xanchor,
                                y=y_annotation,
                                xref='paper',yref='paper',
                                yanchor='bottom',
                                text=f'<b>{str_i}: {str(self.gates_df.at[index, "distance"])}ft</b>',
                                font=dict(size=15, color="#361c00"),
                                showarrow=False)

            annotation_counter += 1

        if len(self.gates_df):
            fig.add_trace(go.Scattermapbox(lon=g_lon, lat=g_lat, 
                                            mode='markers+lines', 
                                            marker=dict(color='#8B4513'), 
                                            line=dict(color='#8B4513')))

        fig.update_yaxes(showticklabels=False)
        fig.update_layout(template='plotly_white', 
                            showlegend=False, height=height, width=width, 
                            margin={'l': 0,'r': 0,'b': 0,'t': 0,},
                            mapbox={'style':'carto-positron', 
                                    'zoom':zoom,
                                    'center': {'lon': center['lon'], 'lat': center['lat']}})

Answer 4

Based on this question on plotly.com with the first version of the function below I came up with the following final solution:

def get_plotting_zoom_level_and_center_coordinates_from_lonlat_tuples(
        longitudes=None, latitudes=None, lonlat_pairs=None):
    """Function documentation:\n
    Basic framework adopted from Krichardson under the following thread:
    https://community.plotly.com/t/dynamic-zoom-for-mapbox/32658/6

    # NOTE:
    # THIS IS A TEMPORARY SOLUTION UNTIL THE DASH TEAM IMPLEMENTS DYNAMIC ZOOM
    # in their plotly-functions associated with mapbox, such as go.Densitymapbox() etc.

    Returns the appropriate zoom-level for these plotly-mapbox-graphics along with
    the center coordinate tuple of all provided coordinate tuples.
    """

    # Check whether the list hasn't already be prepared outside this function
    if lonlat_pairs is None:
        # Check whether both latitudes and longitudes have been passed,
        # or if the list lenghts don't match
        if ((latitudes is None or longitudes is None)
                or (len(latitudes) != len(longitudes))):
            # Otherwise, return the default values of 0 zoom and the coordinate origin as center point
            return 0, (0, 0)

        # Instantiate collator list for all coordinate-tuples
        lonlat_pairs = [(longitudes[i], latitudes[i]) for i in range(len(longitudes))]

    # Get the boundary-box via the planar-module
    b_box = planar.BoundingBox(lonlat_pairs)

    # In case the resulting b_box is empty, return the default 0-values as well
    if b_box.is_empty:
        return 0, (0, 0)

    # Otherwise, get the area of the bounding box in order to calculate a zoom-level
    area = b_box.height * b_box.width

    # * 1D-linear interpolation with numpy:
    # - Pass the area as the only x-value and not as a list, in order to return a scalar as well
    # - The x-points "xp" should be in parts in comparable order of magnitude of the given area
    # - The zoom-levels are adapted to the areas, i.e. start with the smallest area possible of 0
    # which leads to the highest possible zoom value 20, and so forth decreasing with increasing areas
    # as these variables are antiproportional
    zoom = np.interp(x=area,
                     xp=[0, 5**-10, 4**-10, 3**-10, 2**-10, 1**-10, 1**-5],
                     fp=[20, 17, 16, 15, 14, 7, 5])

    # Finally, return the zoom level and the associated boundary-box center coordinates
    return zoom, b_box.center

Answer 5

The Mapbox API documentation shows that zooms are essentially on a log scale. So after some trial and error the following function worked for me:

max_bound = max(abs(x1-x2), abs(y1-y2)) * 111
zoom = 11.5 - np.log(max_bound)

Notes:

• In this example, the xy (lon/lat) coordinates are in decimal degrees
• The 111 is a constant to convert decimal degrees to kilometers
• The value of 11.5 worked for my desired level of zoom/cropping, but I first experimented with values between 10-12

Answer 6

I wrote my own function along with other geojson compatible functions in rv_geojson.py

It takes a list of locations and finds the geometric height and width of the rectangular binding box, good for using with mercator projection. It returns zoom and center.

def zoom_center(lons: tuple=None, lats: tuple=None, lonlats: tuple=None,
        format: str='lonlat', projection: str='mercator',
        width_to_height: float=2.0) -> (float, dict):
    """Finds optimal zoom and centering for a plotly mapbox.
    Must be passed (lons & lats) or lonlats.
    Temporary solution awaiting official implementation, see:
    https://github.com/plotly/plotly.js/issues/3434
    
    Parameters
    --------
    lons: tuple, optional, longitude component of each location
    lats: tuple, optional, latitude component of each location
    lonlats: tuple, optional, gps locations
    format: str, specifying the order of longitud and latitude dimensions,
        expected values: 'lonlat' or 'latlon', only used if passed lonlats
    projection: str, only accepting 'mercator' at the moment,
        raises `NotImplementedError` if other is passed
    width_to_height: float, expected ratio of final graph's with to height,
        used to select the constrained axis.
    
    Returns
    --------
    zoom: float, from 1 to 20
    center: dict, gps position with 'lon' and 'lat' keys

    >>> print(zoom_center((-109.031387, -103.385460),
    ...     (25.587101, 31.784620)))
    (5.75, {'lon': -106.208423, 'lat': 28.685861})
    """
    if lons is None and lats is None:
        if isinstance(lonlats, tuple):
            lons, lats = zip(*lonlats)
        else:
            raise ValueError(
                'Must pass lons & lats or lonlats'
            )
    
    maxlon, minlon = max(lons), min(lons)
    maxlat, minlat = max(lats), min(lats)
    center = {
        'lon': round((maxlon + minlon) / 2, 6),
        'lat': round((maxlat + minlat) / 2, 6)
    }
    
    # longitudinal range by zoom level (20 to 1)
    # in degrees, if centered at equator
    lon_zoom_range = np.array([
        0.0007, 0.0014, 0.003, 0.006, 0.012, 0.024, 0.048, 0.096,
        0.192, 0.3712, 0.768, 1.536, 3.072, 6.144, 11.8784, 23.7568,
        47.5136, 98.304, 190.0544, 360.0
    ])
    
    if projection == 'mercator':
        margin = 1.2
        height = (maxlat - minlat) * margin * width_to_height
        width = (maxlon - minlon) * margin
        lon_zoom = np.interp(width , lon_zoom_range, range(20, 0, -1))
        lat_zoom = np.interp(height, lon_zoom_range, range(20, 0, -1))
        zoom = round(min(lon_zoom, lat_zoom), 2)
    else:
        raise NotImplementedError(
            f'{projection} projection is not implemented'
        )
    
    return zoom, center

Use it as

zoom, center = zoom_center(
    lons=[5, 10, 25, 30],
    lats=[0, 15, 20, 35]
)
fig = px.scatter_mapbox(
    filtered_df, lat="latitude", lon="longitude", color="ID",
    zoom=zoom, center=center
)  # Creates a "Mapbox map" figure