Folium Bushfire Map in a Jupyter Notebook
Using Folium to Display Bushfires from Space¶
Introduction¶
I am writing this in November 2019, with a numbers of bushfires raging around Australia. Thankfully, none have come near my home.
I thought it would be a good exercise to retrieve satellite images of the bushfires, and display them in a Jupyter notebook. But first a disclaimer: the Microsoft Edge browser seems to have a problem displaying more than one embedded map; if you are using Edge, and can't see the maps, please switch to Chrome.
Notebook imports¶
import folium
from folium.raster_layers import WmsTileLayer
from folium.raster_layers import TileLayer
import pandas as pd
import sys
import os
import subprocess
import datetime
import platform
import datetime
Notebook magic commands¶
watermark
is used to document the Notebook, lab_black
formats the python code in each cell.
%matplotlib inline
%load_ext watermark
%load_ext lab_black
Folium¶
Introduction to Folium¶
Folium is a python package that allows access to a lot of the functions of the leaflet
Javascript package, but without any of that pesky Javascript.
In a comparison with Cartopy:
- Cartopy has support for a much larger range of projections
- Folium gives you good basemaps out-of-the-box
- Folium allows you to build interactive maps very easily
- Folium allows you to integrate maps into a Flask webserver fairly easily
First Steps¶
The code below shows how to create a map of my local patch. By default, we get an Open Street Map map background (base map). The resulting map is fully pan-and-zoom-able (with a zoom contrl for those of us without a mouse wheel).
start_coords = (-26.52, 153.09)
folium_map = folium.Map(
location=start_coords, zoom_start=13, width='80%'
)
folium_map
Folium Layers¶
Unlike Cartopy, which is oriented towrds producing static cartographic artefacts (maps), Folium fully supports having layers over the basemap, or even turning switching between basemaps.
The example below shows us setting up a map, and then adding a number of the out-of-the-box basemaps. Finally, we add a LayerControl
that allows us the turn these on or off.
start_coords = (-26.52, 153.09)
folium_map = folium.Map(
location=start_coords, zoom_start=13, width='80%'
)
folium.raster_layers.TileLayer(
tiles='OpenStreetMap', name='Open Street Map'
).add_to(folium_map)
folium.raster_layers.TileLayer(
tiles='stamentoner', name='Black/White Map'
).add_to(folium_map)
folium.raster_layers.TileLayer(
tiles='stamenterrain', name='stamenterrain'
).add_to(folium_map)
folium.raster_layers.TileLayer(
tiles='CartoDB dark_matter', name='CartoDB dark_matter'
).add_to(folium_map)
folium.LayerControl().add_to(folium_map)
<folium.map.LayerControl at 0x1fd96ec00f0>
folium_map
External Tile Suppliers¶
The strength of Folium is that it allows you to nominate external Tile providers, to add as layers to your interactive map. In the example below, I nominate two sources of tile from ESRI (World Imagery Tiles, and World Topographic Tiles), and add them to the interactive map. As you alter the pan and zoom of the map, requests are made to get tiles of the appropriate zoom level and extent, to display.
Finally, I add a Layer Control so that I can turn then on or off as needed.
start_coords = (-26.52, 153.09)
folium_map = folium.Map(
location=start_coords, zoom_start=13, width='80%'
)
'''
"http://server.arcgisonline.com/ArcGIS/rest/services/World_Topo_Map/MapServer" // World Topographic Map
"http://server.arcgisonline.com/ArcGIS/rest/services/World_Street_Map/MapServer" // World Street Map
"http://services.arcgisonline.com/ArcGIS/rest/services/Canvas/World_Light_Gray_Base/MapServer" // Light Gray Canvas
"http://services.arcgisonline.com/ArcGIS/rest/services/NatGeo_World_Map/MapServer" // National Geographic World Map
"http://services.arcgisonline.com/ArcGIS/rest/services/Ocean_Basemap/MapServer" // Ocean Basemap
"http://services.arcgisonline.com/ArcGIS/rest/services/World_Terrain_Base/MapServer" // Terrain with Labels
"http://server.arcgisonline.com/ArcGIS/rest/services/World_Imagery/MapServer" // World Imagery
'''
url = (
'http://services.arcgisonline.com/arcgis/rest/services/World_Imagery'
+ '/MapServer/tile/{z}/{y}/{x}'
)
WmsTileLayer(
url=url,
layers=None,
name='ESRI Imagery',
attr='ESRI World Imagery',
).add_to(folium_map)
World_Topo_Map = (
'http://services.arcgisonline.com/arcgis/rest/services/World_Topo_Map'
+ '/MapServer/tile/{z}/{y}/{x}'
)
WmsTileLayer(
url=World_Topo_Map,
layers=None,
name='ESRI Topo Map',
attr='ESRI Topo Map',
).add_to(folium_map)
folium.LayerControl().add_to(folium_map)
folium_map
Bushfires¶
NASA have a service called Global Imagery Browse Services (GIBS), that provides almost up-to-date images from space. This is in contrast to (say) LANDSAT, that only revisits a site in a period measutre in about 10 days. LANDSAT has a goal to image a site with exactly the same sun angle each time, to allow effective image comparison - this precludes timeliness.
NASA have a range of other satellites that do revisit a given site more frequently. I have chosen to use the Moderate Resolution Imaging Spectroradiometer (MODIS). MODIS can be used to give you "True Color" images, that approximate what a human might see from space.
The example below illustrates a number of points:
The range of Coordinate Reference Systems supported by Folium is minimal, but both GIBS and Folium support Web Mercator, identified by the code EPSG:3857. Fortunately, I live relatively close to the Equator, so Mercator distortion is minimal.
The example shows a common Folium pattern: the tile source URL contains a number of placeholders (e.g.
{layer}
,{time}
) that are filled in by Folium from the keyword arguements supplied in theTileLayer
call, or by the current pan and zoom values.The tile size is crucial to get right
If you get tiles from a non-Folium source, you must supply an attribution, to be included in the map (bottom right).
So picking a day a few days ago, and zooming out, the smoke from the bush fires is very clear.
start_coords = (-26.52, 153.09)
folium_map = folium.Map(
tiles='stamenterrain',
location=start_coords,
zoom_start=5,
width='80%',
crs='EPSG3857',
)
TileLayer(
tiles='https://gibs-{s}.earthdata.nasa.gov/wmts/epsg3857/best/'
+ '{layer}/default/{time}/{tileMatrixSet}/{z}/{y}/{x}.jpg',
subdomains='abc',
name='GIBS',
attr='NASA GIBS',
overlay=True,
layer='MODIS_Terra_CorrectedReflectance_TrueColor',
tileMatrixSet='GoogleMapsCompatible_Level9',
time='2019-11-09',
tileSize=256,
).add_to(folium_map)
folium.LayerControl().add_to(folium_map)
folium_map
Visible Infrared Imaging Radiometer Suite (VIIRS)¶
MODIS is not the only system that can provide imagery.
The Visible Infrared Imaging Radiometer Suite (VIIRS) also provide "True Color" images, but can also provide images that "peer through" the smoke.
start_coords = (-26.52, 153.09)
folium_map = folium.Map(
location=start_coords,
zoom_start=5,
width='80%',
crs='EPSG3857',
)
src = (
'http://map1.vis.earthdata.nasa.gov/wmts-webmerc/VIIRS_SNPP_CorrectedReflectance_BandsM3-I3-M11/'
+ 'default/{time}/GoogleMapsCompatible_Level9/{z}/{y}/{x}.jpg'
)
TileLayer(
tiles=src,
subdomains='abc',
name='VIIRS',
attr='NASA VIIRS',
overlay=True,
layer='VIIRS_SNPP_CorrectedReflectance_BandsM11-I2-I1',
tileMatrixSet='GoogleMapsCompatible_Level9',
time='2019-11-09',
tileSize=256,
).add_to(folium_map)
folium.LayerControl().add_to(folium_map)
folium_map
The example below show the VIIRS True Image results.
start_coords = (-26.52, 153.09)
folium_map = folium.Map(
location=start_coords,
zoom_start=5,
width='80%',
crs='EPSG3857',
)
src = (
'http://map1.vis.earthdata.nasa.gov/wmts-webmerc/{layer}/'
+ 'default/{time}/GoogleMapsCompatible_Level9/{z}/{y}/{x}.jpg'
)
TileLayer(
tiles=src,
subdomains='abc',
name='VIIRS',
attr='NASA VIIRS',
overlay=True,
layer='VIIRS_SNPP_CorrectedReflectance_TrueColor',
tileMatrixSet='GoogleMapsCompatible_Level9',
time='2019-11-10',
tileSize=256,
).add_to(folium_map)
folium.LayerControl().add_to(folium_map)
folium_map
Mapping Details¶
Especially in the True Color images, it was sometimes hard to get spatial awareness, because the smoke obscured identifing features (like coastlines). In this example, I have added an Open Street Map layer, with an opacity
set to 0.25. Now, when I zoom in with this layer enabled, I get the coastlines and major towns identified. This is at the cost of losing some detail of the satellite image.
start_coords = (-26.52, 153.09)
folium_map = folium.Map(
location=start_coords,
zoom_start=5,
width='80%',
crs='EPSG3857',
)
src = (
'http://map1.vis.earthdata.nasa.gov/wmts-webmerc/{layer}/'
+ 'default/{time}/GoogleMapsCompatible_Level9/{z}/{y}/{x}.jpg'
)
TileLayer(
tiles=src,
subdomains='abc',
name='True Color',
attr='NASA VIIRS',
overlay=True,
layer='VIIRS_SNPP_CorrectedReflectance_TrueColor',
tileMatrixSet='GoogleMapsCompatible_Level9',
time='2019-11-11',
tileSize=256,
).add_to(folium_map)
TileLayer(
tiles='openstreetmap',
opacity=0.25,
name='Streets',
overlay=True,
).add_to(folium_map)
TileLayer(
tiles=src,
subdomains='abc',
name='Bands M11-12',
attr='NASA VIIRS',
overlay=True,
show=False,
layer='VIIRS_SNPP_CorrectedReflectance_BandsM11-I2-I1',
tileMatrixSet='GoogleMapsCompatible_Level9',
time='2019-11-11',
tileSize=256,
).add_to(folium_map)
TileLayer(
tiles=src,
subdomains='abc',
name='Bands M3-13',
attr='NASA VIIRS',
overlay=True,
show=False,
layer='VIIRS_SNPP_CorrectedReflectance_BandsM3-I3-M11',
tileMatrixSet='GoogleMapsCompatible_Level9',
time='2019-11-11',
tileSize=256,
).add_to(folium_map)
folium.LayerControl().add_to(folium_map)
folium_map
Reproducibility¶
This section of the notebook provides information to support reproducibility.
Notebook version status¶
theNotebook = 'FoliumExamples.ipynb'
# show info to support reproducibility
def python_env_name():
envs = subprocess.check_output(
'conda env list'
).splitlines()
# get unicode version of binary subprocess output
envu = [x.decode('ascii') for x in envs]
active_env = list(
filter(lambda s: '*' in str(s), envu)
)[0]
env_name = str(active_env).split()[0]
return env_name
# end python_env_name
print('python version : ' + sys.version)
print('python environment :', python_env_name())
print('current wkg dir: ' + os.getcwd())
print('Notebook name: ' + theNotebook)
print(
'Notebook run at: '
+ str(datetime.datetime.now())
+ ' local time'
)
print(
'Notebook run at: '
+ str(datetime.datetime.utcnow())
+ ' UTC'
)
print('Notebook run on: ' + platform.platform())
python version : 3.7.1 (default, Dec 10 2018, 22:54:23) [MSC v.1915 64 bit (AMD64)] python environment : ac5-py37 current wkg dir: C:\Users\donrc\Documents\JupyterNotebooks\QLDCrashLocationsNotebookProject\develop Notebook name: FoliumExamples.ipynb Notebook run at: 2019-11-14 12:10:39.987424 local time Notebook run at: 2019-11-14 02:10:39.987424 UTC Notebook run on: Windows-10-10.0.18362-SP0
%watermark
2019-11-14T12:10:40+10:00 CPython 3.7.1 IPython 7.2.0 compiler : MSC v.1915 64 bit (AMD64) system : Windows release : 10 machine : AMD64 processor : Intel64 Family 6 Model 94 Stepping 3, GenuineIntel CPU cores : 8 interpreter: 64bit
%watermark -h -iv
folium 0.10.0 platform 1.0.8 pandas 0.23.4 host name: DESKTOP-SODFUN6