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 GLMakie, GeoMakie
using FHist
import GLMakie: Point2d

First, we generate random points in a normal distribution:

random_data = randn(Point2d, 100_000)

100000-element Vector{Point{2, Float64}}:
 [0.9611710548410338, 1.3539558849788076]
 [1.546894144996102, 1.4777685598940002]
 [0.6192651589082052, -1.3331887275363703]
 [1.7207215700950595, 1.2141049455602626]
 [-1.824089736188693, -0.31675312211973733]
 [-1.3755161447585458, 0.8251206001491073]
 [-0.7823761088620949, 0.1443903087935883]
 [-0.15507315641436806, 0.008193924347544183]
 [-0.05498719599137963, 1.0111570345744907]
 [0.5653623526371226, 1.5934084522558587]
 ⋮
 [0.8442354494339527, 0.12563176664714457]
 [1.271149215897578, 1.5508830793492443]
 [-0.04963137329937564, -0.1361796985592117]
 [-1.4157917229096568, 0.45735187687617673]
 [1.0573562293503542, 0.5218051873753633]
 [-0.0220598242710306, 0.3667565464423203]
 [1.238655290092802, -0.8221921229511365]
 [-1.7811991239262013, 0.5302684100024248]
 [-1.1446363747119364, -0.4867754087669667]

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}}:
 [41.87141742563803, 26.254712114976076]
 [67.38722533536918, 28.655577735597376]
 [26.976998355498708, -25.852013809849165]
 [74.95965548678072, 23.542846688503687]
 [-79.46267460000718, -6.14219571334389]
 [-59.9214992823675, 16.000007132721446]
 [-34.08258756130476, 2.7998888528244006]
 [-6.7554394517333325, 0.15888931627933575]
 [-2.3954028003904204, 19.607460730673434]
 [24.62883473735576, 30.89796400286851]
 ⋮
 [36.77736104383264, 2.436139834695029]
 [55.374971146964796, 30.07334967412128]
 [-2.1620875268341515, -2.6406759786211853]
 [-61.67602105695734, 8.868562111841591]
 [46.06152445367757, 10.11838356524869]
 [-0.9609903520672586, 7.111817785168599]
 [53.959440868231525, -15.94322069926302]
 [-77.59423430455841, 10.282494875006709]
 [-49.86370241253172, -9.439117155598131]

finally, we can create the histogram.

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

• edges: ([-170.0, -169.0, -168.0, -167.0, -166.0, -165.0, -164.0, -163.0, -162.0, -161.0 … 182.0, 183.0, 184.0, 185.0, 186.0, 187.0, 188.0, 189.0, 190.0, 191.0], [-97.0, -96.0, -95.0, -94.0, -93.0, -92.0, -91.0, -90.0, -89.0, -88.0 … 75.0, 76.0, 77.0, 78.0, 79.0, 80.0, 81.0, 82.0, 83.0, 84.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: 32.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.