Trouble geo mapping with datashader, holoviews and bokeh

Question

I'm trying to map google phone history locations on to a map using holoviews, datashader and bokeh. Mostly very similar to the examples given in the datashader website. But when I do the map overlay doesn't work as the lat/long gets mangled up.

import datashader as ds
import geoviews as gv
import holoviews as hv
from holoviews.operation.datashader import datashade, dynspread
from datashader import transfer_functions as tf
from colorcet import fire
hv.extension('bokeh')

> df2.head()

lat long
0   -37.7997515 144.9636466
1   -37.7997515 144.9636466
2   -37.7997369 144.9636036
3   -37.7997387 144.9636358
4   -37.7997515 144.9636466

This works to produce an image of the data,

ds_viz = ds.Canvas().points(df2,'lat','long')  
tf.set_background(tf.shade(ds_viz, cmap=fire),"black")

However when I try to overlay it with a map it doesn't work,

from bokeh.models import WMTSTileSource

url = 'https://server.arcgisonline.com/ArcGIS/rest/services/World_Imagery/MapServer/tile/{Z}/{Y}/{X}.jpg'
tile_opts  = dict(width=1000,height=600,bgcolor='black',show_grid=False)
map_tiles  = gv.WMTS(url).opts(style=dict(alpha=0.5), plot=tile_opts)
points     = hv.Points(df2, kdims=['long','lat'])
trips = datashade(points, cmap=fire,width=1000, height=600)

map_tiles * trips

What am I doing wrong?

Urana

Answer 1

It looks like your points are in lon,lat but your map is in Web Mercator coordinates, so you need to project your points into Web Mercator before you overlay them. GeoViews offers comprehensive support for projections, but for this specific case Datashader provides the special-purpose function datashader.utils.lnglat_to_meters. Something like this should work:

df2.loc[:, 'lon'], df.loc[:, 'lat'] = lnglat_to_meters(df2.lon,df2.lat)

Projecting can be slow, so you may want to save the resulting df2 to a Parquet file so that you only have to do it once.

Hartz

I find this works very well, but I was wondering how you know that a map is in Web Mercator coordinates? When I just plot the tiles the axes range seems to show coordinates that range from -180 to 180 (x) and -90 to 90 (y), which seem to be "regular" latitude/longitude numbers in degrees. – Eyeful

Bokeh's maps are normally always in Web Mercator; it doesn't currently support any other projection. Bokeh includes special code to map to and from Web Mercator for the axis labels and hovering, so that most users don't see the underlying values like 9000290, but if you inspect the actual columns of overlaid data using Python or JavaScript you'll see the raw Web Mercator values. – Hartz