import CairoMakie as Mke
2 Scientific visualization
The visualization ecosystem in Julia is evolving very quickly. Among the various visualization projects, Makie.jl by Danisch and Krumbiegel (2021) is the most advanced for scientific visualization.
Makie.jl is currently organized in backend modules:
- GLMakie.jl is the preferred backend for interactive high-performance visualization.
- WGLMakie.jl is the preferred backend for interactive visualization on the web browser.
- CairoMakie.jl is the preferred backend for publication-quality static visualization.
In this book, we use CairoMakie.jl:
We import the backend as Mke
to avoid polluting the Julia session with names from the visualization stack.
Makie.jl provides a plot recipe system developed after Plots.jl by Breloff (2023), which enables automatic visualization of custom Julia types. The GeoStats.jl framework is integrated with this system, and provides powerful visualization functions for geospatial data.
Julia will automatically trigger the compilation of these visualization functions whenever GeoStats.jl and Makie.jl are loaded in the same session:
using GeoStats
2.1 viz/viz!
The main visualization function that the framework provides is the viz
/viz!
function. The viz
function creates a scene and displays geometries within a geospatial domain. On the other hand, the viz!
function adds more geometries to an existing scene.
Let’s create a small geotable over a Cartesian grid for illustration purposes:
= georef((A=rand(10, 10), B=rand(10, 10))) img
A | B | geometry |
---|---|---|
Continuous | Continuous | Quadrangle |
[NoUnits] | [NoUnits] | 🖈 Cartesian{NoDatum} |
0.720268 | 0.747247 | Quadrangle((x: 0.0 m, y: 0.0 m), ..., (x: 0.0 m, y: 1.0 m)) |
0.651123 | 0.990592 | Quadrangle((x: 1.0 m, y: 0.0 m), ..., (x: 1.0 m, y: 1.0 m)) |
0.965049 | 0.249579 | Quadrangle((x: 2.0 m, y: 0.0 m), ..., (x: 2.0 m, y: 1.0 m)) |
0.182212 | 0.807963 | Quadrangle((x: 3.0 m, y: 0.0 m), ..., (x: 3.0 m, y: 1.0 m)) |
0.549511 | 0.0163637 | Quadrangle((x: 4.0 m, y: 0.0 m), ..., (x: 4.0 m, y: 1.0 m)) |
0.67308 | 0.730661 | Quadrangle((x: 5.0 m, y: 0.0 m), ..., (x: 5.0 m, y: 1.0 m)) |
0.224735 | 0.312531 | Quadrangle((x: 6.0 m, y: 0.0 m), ..., (x: 6.0 m, y: 1.0 m)) |
0.808508 | 0.937077 | Quadrangle((x: 7.0 m, y: 0.0 m), ..., (x: 7.0 m, y: 1.0 m)) |
0.225596 | 0.142237 | Quadrangle((x: 8.0 m, y: 0.0 m), ..., (x: 8.0 m, y: 1.0 m)) |
0.181082 | 0.34248 | Quadrangle((x: 9.0 m, y: 0.0 m), ..., (x: 9.0 m, y: 1.0 m)) |
⋮ | ⋮ | ⋮ |
The georef
function creates a CartesianGrid
starting at the origin whenever the domain is omitted. The size of the grid is taken as the size of the first array in the named tuple:
img.geometry
10×10 CartesianGrid
├─ minimum: Point(x: 0.0 m, y: 0.0 m)
├─ maximum: Point(x: 10.0 m, y: 10.0 m)
└─ spacing: (1.0 m, 1.0 m)
To create a named tuple with a single key in Julia, we need an extra comma after the key/value pair:
=rand(10, 10),) (A
(A = [0.5658614601119324 0.42759266196706325 … 0.013737562681665816 0.21942584209307947; 0.5455953163641568 0.03271726014029375 … 0.9397325748477771 0.050962514548867466; … ; 0.8321079320240431 0.6880013584721631 … 0.1876006736261745 0.35487855179806505; 0.945970490433865 0.34320454824961244 … 0.14454699975064045 0.9704866067352488],)
or a semicolon before the key/value pair:
=rand(10, 10)) (; A
(A = [0.5297225446359806 0.18823426438307245 … 0.9127948392432641 0.9155871004704819; 0.025250888511517466 0.20962085607203285 … 0.8054814744566865 0.7031809230111199; … ; 0.9754621304908608 0.1627062429867605 … 0.12864785650330057 0.17581037377686004; 0.0066567798376847565 0.3188160525889404 … 0.3651655044291209 0.39285592348686615],)
By default, all geometries are displayed with a single color:
viz(img.geometry)
┌ Warning: Found `resolution` in the theme when creating a `Scene`. The `resolution` keyword for `Scene`s and `Figure`s has been deprecated. Use `Figure(; size = ...` or `Scene(; size = ...)` instead, which better reflects that this is a unitless size and not a pixel resolution. The key could also come from `set_theme!` calls or related theming functions.
└ @ Makie ~/.julia/packages/Makie/We6MY/src/scenes.jl:227
We can pass any vector of Colors.jl or numbers (automatically converted to colors) to the function via the color
option. It is common to pass colors from another column of the geotable:
viz(img.geometry, color = img.A)
┌ Warning: Found `resolution` in the theme when creating a `Scene`. The `resolution` keyword for `Scene`s and `Figure`s has been deprecated. Use `Figure(; size = ...` or `Scene(; size = ...)` instead, which better reflects that this is a unitless size and not a pixel resolution. The key could also come from `set_theme!` calls or related theming functions.
└ @ Makie ~/.julia/packages/Makie/We6MY/src/scenes.jl:227
but any vector with the same length can be passed:
viz(img.geometry, color = 1:length(img.A))
┌ Warning: Found `resolution` in the theme when creating a `Scene`. The `resolution` keyword for `Scene`s and `Figure`s has been deprecated. Use `Figure(; size = ...` or `Scene(; size = ...)` instead, which better reflects that this is a unitless size and not a pixel resolution. The key could also come from `set_theme!` calls or related theming functions.
└ @ Makie ~/.julia/packages/Makie/We6MY/src/scenes.jl:227
The alpha
option can be used to control the transparency of each geometry in the domain:
viz(img.geometry, color = img.B, alpha = rand(length(img.B)))
┌ Warning: Found `resolution` in the theme when creating a `Scene`. The `resolution` keyword for `Scene`s and `Figure`s has been deprecated. Use `Figure(; size = ...` or `Scene(; size = ...)` instead, which better reflects that this is a unitless size and not a pixel resolution. The key could also come from `set_theme!` calls or related theming functions.
└ @ Makie ~/.julia/packages/Makie/We6MY/src/scenes.jl:227
Other aesthetic options are available in the official documentation.
As another example, consider the visualization of data over a GeometrySet
:
= GeometrySet([
gset Triangle((12, 12), (15, 15), (12, 15)),
Quadrangle((5, 12), (10, 15), (10, 18), (5, 15))
])
= georef((A=[0.1, 0.2], B=[0.3, 0.4]), gset)
gis
viz(gis.geometry, color = gis.A)
┌ Warning: Found `resolution` in the theme when creating a `Scene`. The `resolution` keyword for `Scene`s and `Figure`s has been deprecated. Use `Figure(; size = ...` or `Scene(; size = ...)` instead, which better reflects that this is a unitless size and not a pixel resolution. The key could also come from `set_theme!` calls or related theming functions.
└ @ Makie ~/.julia/packages/Makie/We6MY/src/scenes.jl:227
We can create a scene with the geometries from the first geotable (“raster data”), and then add the geometries from the second geotable (“vector data”):
viz(img.geometry)
viz!(gis.geometry)
# display current figure
current_figure() Mke.
┌ Warning: Found `resolution` in the theme when creating a `Scene`. The `resolution` keyword for `Scene`s and `Figure`s has been deprecated. Use `Figure(; size = ...` or `Scene(; size = ...)` instead, which better reflects that this is a unitless size and not a pixel resolution. The key could also come from `set_theme!` calls or related theming functions.
└ @ Makie ~/.julia/packages/Makie/We6MY/src/scenes.jl:227
Let’s add an additional set of points and line segments to conclude the example:
viz!([Point(-20, -10), Point(-20, 0), Point(-40, 10)])
viz!([Segment((-40, -10), (0, 0)), Segment((-40, 0), (-20, 10))])
current_figure() Mke.
Makie.jl can set the aspect ratio of the axis after the visualization is created. The following code can be used to adjust the aspect ratio for the data in the scene:
= Mke.current_axis()
ax = Mke.DataAspect()
ax.aspect current_figure() Mke.
Alternatively, it is possible to set the aspect ratio in the viz/viz!
call directly:
viz(CartesianGrid(20, 10), color = 1:200, axis = (; aspect = Mke.DataAspect()))
┌ Warning: Found `resolution` in the theme when creating a `Scene`. The `resolution` keyword for `Scene`s and `Figure`s has been deprecated. Use `Figure(; size = ...` or `Scene(; size = ...)` instead, which better reflects that this is a unitless size and not a pixel resolution. The key could also come from `set_theme!` calls or related theming functions.
└ @ Makie ~/.julia/packages/Makie/We6MY/src/scenes.jl:227
Makie.jl dispatches the viz
and viz!
functions whenever it encounters a geospatial domain, a vector of geometries or a single geometry from Meshes.jl. This means that you can replace viz
with Mke.plot
and viz!
with Mke.plot!
in scripts and the result will be the same.
In the case of Mesh
domains, it is also possible to specify a color for each vertex of the mesh. In this case, the viz
and viz!
functions fill in the domain with using the interpolation routine from the graphics library:
= CartesianGrid(10, 10)
grid
= Mke.Figure()
fig viz(fig[1,1], grid, color = 1:nelements(grid))
viz(fig[1,2], grid, color = 1:nvertices(grid))
fig
┌ Warning: Found `resolution` in the theme when creating a `Scene`. The `resolution` keyword for `Scene`s and `Figure`s has been deprecated. Use `Figure(; size = ...` or `Scene(; size = ...)` instead, which better reflects that this is a unitless size and not a pixel resolution. The key could also come from `set_theme!` calls or related theming functions.
└ @ Makie ~/.julia/packages/Makie/We6MY/src/scenes.jl:227
The viz
and viz!
functions are aware of coordinate reference systems, which is a quite unique feature of the framework, explored in the chapter Map projections.
Below is a preview of this feature in action:
using GeoIO
= GeoIO.load("data/countries.geojson")
world
viz(world.geometry, color = 1:nrow(world))
┌ Warning: Found `resolution` in the theme when creating a `Scene`. The `resolution` keyword for `Scene`s and `Figure`s has been deprecated. Use `Figure(; size = ...` or `Scene(; size = ...)` instead, which better reflects that this is a unitless size and not a pixel resolution. The key could also come from `set_theme!` calls or related theming functions.
└ @ Makie ~/.julia/packages/Makie/We6MY/src/scenes.jl:227
2.2 viewer
As geospatial data scientists we are often interested in quick inspection of intermediate results from multivariate geostatistical analysis. Visualizing all the variables manually with viz
/viz!
can be very time consuming. To address this issue, the framework provides a basic viewer
that displays all variables stored in a geotable:
= georef((A=rand(1000), B=rand(1000)), rand(Point, 1000))
geotable
viewer(geotable)
┌ Warning: Found `resolution` in the theme when creating a `Scene`. The `resolution` keyword for `Scene`s and `Figure`s has been deprecated. Use `Figure(; size = ...` or `Scene(; size = ...)` instead, which better reflects that this is a unitless size and not a pixel resolution. The key could also come from `set_theme!` calls or related theming functions.
└ @ Makie ~/.julia/packages/Makie/We6MY/src/scenes.jl:227
The rand(Point, 1000)
function call creates 1000
random Point
s in 3D space.
It adds interactive elements to the scene, including a menu to select the variable used as color, and a color bar that automatically updates upon menu selection. The viewer
will be particularly useful when we start to work with geospatial transforms in Part II of the book. The pipe operator (|>
) in Julia will be preferred for reasons that will become clear later:
= georef((A=rand(1000), B=rand(1000)), CartesianGrid(10, 10, 10))
geotable
|> viewer geotable
┌ Warning: Found `resolution` in the theme when creating a `Scene`. The `resolution` keyword for `Scene`s and `Figure`s has been deprecated. Use `Figure(; size = ...` or `Scene(; size = ...)` instead, which better reflects that this is a unitless size and not a pixel resolution. The key could also come from `set_theme!` calls or related theming functions.
└ @ Makie ~/.julia/packages/Makie/We6MY/src/scenes.jl:227
The viz/viz!
and the viewer
automatically select color schemes for variables based on data science traits from the DataScienceTraits.jl module. Additionally, the viewer
recognizes units from the Unitful.jl module:
= georef((; A=[1,2,3,4]u"m"), CartesianGrid(2, 2))
geotable
|> viewer geotable
┌ Warning: Found `resolution` in the theme when creating a `Scene`. The `resolution` keyword for `Scene`s and `Figure`s has been deprecated. Use `Figure(; size = ...` or `Scene(; size = ...)` instead, which better reflects that this is a unitless size and not a pixel resolution. The key could also come from `set_theme!` calls or related theming functions.
└ @ Makie ~/.julia/packages/Makie/We6MY/src/scenes.jl:227
2.3 cbar
Any vector of Julia objects implementing the getcolors
function from the Colorfy.jl module can be passed to the color
option of viz
/viz!
. The final vector of colors will be a function of the colormap
and colorrange
options.
To reproduce this behavior in the colorbar, the framework provides the cbar
function. It is similar to Makie’s colorbar, but addresses various practical challenges with missing
values, units, etc.
= CartesianGrid(2, 2)
grid = [1, missing, 3, 4]
vals = "cividis"
cmap
= Mke.Figure()
fig viz(fig[1,1], grid, color = vals, colormap = cmap)
cbar(fig[1,2], vals, colormap = cmap)
fig
┌ Warning: Found `resolution` in the theme when creating a `Scene`. The `resolution` keyword for `Scene`s and `Figure`s has been deprecated. Use `Figure(; size = ...` or `Scene(; size = ...)` instead, which better reflects that this is a unitless size and not a pixel resolution. The key could also come from `set_theme!` calls or related theming functions.
└ @ Makie ~/.julia/packages/Makie/We6MY/src/scenes.jl:227
We are now equipped with a set of visualization functions that can really improve the speed at which we explore and analyze geospatial data. These functions provide a consistent set of aesthetic options that we will cover in more detail with future examples.
Before we start learning the advanced features of the framework, we would like to say a few words about integration with existing GIS technology.