Histogram

You can use any Makie compatible recipe in GeoMakie, and this includes histograms! Until we get the datashader recipe working, you can also consider this a replacement for that if you use it with the correct nbins.

We use the excellent FHist.jl to create the histogram.

using CairoMakie, GeoMakie
using FHist
import CairoMakie: Point2d

First, we generate random points in a normal distribution:

random_data = randn(Point2d, 100_000)

100000-element Vector{Point{2, Float64}}:
 [-0.1822539443471449, -0.035948734140347174]
 [2.5018591130783743, -0.26317663985245976]
 [-0.07878914258354008, 1.186844985318085]
 [0.20685874306924862, 0.8110986648701664]
 [-0.003974252667858806, -0.20749423040056625]
 [-0.4915499470097861, 1.6906096337105778]
 [2.210843259349297, 0.4958926547033992]
 [-0.47784374516582684, -0.2864602673231118]
 [0.8920878446322862, 0.4367448603863451]
 [0.9209013837830493, -2.2945480221232777]
 ⋮
 [-0.7062547140254513, 0.6236737164043309]
 [-0.6782070587963281, 1.0854292671901782]
 [-0.6332343287824113, -2.0033991176298507]
 [-0.10145408539652316, -0.6748612403556961]
 [0.3300509855838988, 0.020800024059671643]
 [-1.2188663291364965, 0.8483133977995917]
 [0.6226806710429323, -0.5667966243723376]
 [0.9048864027665764, 0.2594629450451621]
 [0.28879119943979725, -0.163020678971824]

then, we rescale them to be within the lat/long bounds of the Earth:

xmin, xmax = extrema(first, random_data)
ymin, ymax = extrema(last, random_data)
latlong_data = random_data .* (Point2d(1/(xmax - xmin), 1/(ymax - ymin)) * Point2d(360, 180),)

100000-element Vector{Point{2, Float64}}:
 [-7.515668471426826, -0.7582899452598816]
 [103.17002314254485, -5.551355412078886]
 [-3.249054921291427, 25.03489039239155]
 [8.53030500323923, 17.10902975842091]
 [-0.16388762163861664, -4.376810265360491]
 [-20.27021391556978, 35.661124578052856]
 [91.16930247561103, 10.46018512144991]
 [-19.705006564745876, -6.042492054922624]
 [36.78733270581761, 9.212542366080038]
 [37.97552651156119, -48.400388378048916]
 ⋮
 [-29.12406810185343, 13.155553862506595]
 [-27.967457314317013, 22.895662928985818]
 [-26.112901407445737, -42.2589958608919]
 [-4.183696949020455, -14.23528547651927]
 [13.61042580012133, 0.4387483866340198]
 [-50.26280925545535, 17.894024236550326]
 [25.67769660013714, -11.955808501930042]
 [37.3151112413888, 5.473020040904381]
 [11.908981834275577, -3.4387008246566126]

finally, we can create the histogram.

h = Hist2D((first.(latlong_data), last.(latlong_data)); nbins = (360, 180))

• edges: ([-195.0, -194.0, -193.0, -192.0, -191.0, -190.0, -189.0, -188.0, -187.0, -186.0 … 157.0, 158.0, 159.0, 160.0, 161.0, 162.0, 163.0, 164.0, 165.0, 166.0], [-96.0, -95.0, -94.0, -93.0, -92.0, -91.0, -90.0, -89.0, -88.0, -87.0 … 76.0, 77.0, 78.0, 79.0, 80.0, 81.0, 82.0, 83.0, 84.0, 85.0])
• bin counts: [0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0; … ; 0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0]
• maximum count: 35.0
• total count: 100000.0

This is what the histogram looks like without any projection,

plot(h)

It's simple to plot to GeoAxis:

plot(h; axis = (; type = GeoAxis))

The projection can also be arbitrary!

plot(h; axis = (; type = GeoAxis, dest = "+proj=tissot +lat_1=60 +lat_2=65"))

---

This page was generated using Literate.jl.