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.6427821304905709, 0.027817156749768902]
 [0.07132541393147905, 0.7349588837118032]
 [1.3481095642538832, -1.849087298967047]
 [-0.9441642587268025, 2.5712152761096654]
 [-1.854087603449389, 2.199986256417337]
 [-1.4119715950147462, 1.0343238398133665]
 [-0.9160353239480383, 1.1353530632783329]
 [0.6287498282688666, 0.8357737542422743]
 [1.4015985416114123, -0.30291117877094176]
 [-1.3029607796803144, -0.7817085957803981]
 ⋮
 [0.004689799771050458, -0.23509804528589323]
 [-0.1684231274274612, -0.905321822096536]
 [0.5639152799763193, -0.9418042555478809]
 [0.08539984242127531, -0.3623615775199577]
 [-1.5223488222469312, 0.43385090260353476]
 [-1.5785425077965682, 0.5712406488089823]
 [-1.004213868083417, 0.4348680446998973]
 [0.5797457018632493, -0.2250866303570877]
 [1.5675651861789273, -0.39738844730592526]

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}}:
 [-26.692881405463616, 0.5921608795252837]
 [2.961937995717977, 15.645520601141873]
 [55.98308822421131, -39.362655613049114]
 [-39.20840886823892, 54.735036835229174]
 [-76.99489168506184, 46.832456971162216]
 [-58.63509351893036, 22.018286060200605]
 [-38.04029562349372, 24.16895711414532]
 [26.110160509405553, 17.79162859268013]
 [58.20433063493977, -6.44825607637026]
 [-54.108189879879106, -16.640710399478216]
 ⋮
 [0.19475381029724606, -5.004676305472021]
 [-6.994124996796962, -19.272140975751267]
 [23.41776937645828, -20.04876491593705]
 [3.5464082737592877, -7.713813183015912]
 [-63.21874028914163, 9.23565029954694]
 [-65.55230140255901, 12.16035010559197]
 [-41.70209533674295, 9.25730282729456]
 [24.07516097772103, -4.791557174611054]
 [65.09644500862707, -8.459451646574863]

finally, we can create the histogram.

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

• edges: ([-193.0, -192.0, -191.0, -190.0, -189.0, -188.0, -187.0, -186.0, -185.0, -184.0 … 159.0, 160.0, 161.0, 162.0, 163.0, 164.0, 165.0, 166.0, 167.0, 168.0], [-90.0, -89.0, -88.0, -87.0, -86.0, -85.0, -84.0, -83.0, -82.0, -81.0 … 82.0, 83.0, 84.0, 85.0, 86.0, 87.0, 88.0, 89.0, 90.0, 91.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.