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.7773129347647819, 1.0449590000930211]
 [-0.41264525528195983, -0.6911042753586699]
 [0.7519622670966952, 1.1061280870878314]
 [-0.05979725620625776, -0.7489779054709294]
 [1.5578330490893322, 1.010620125702205]
 [-1.5572632100446264, -0.7740916041415213]
 [-1.1021484940386728, 0.9989290830464916]
 [1.1342919190308263, -1.107841110439419]
 [-0.3096660088523369, 2.6711237841834983]
 [1.5714685585059651, -1.0335588048174345]
 ⋮
 [1.1377862006810024, -0.46578642704654144]
 [-0.3932707672066858, -0.9767255283392303]
 [0.0951394065224418, 0.015297271570050995]
 [0.07639690165267327, 1.9388891720587722]
 [-1.2160407092183867, -0.6186041816954277]
 [-0.3707557260615645, -1.187892692044583]
 [0.4979637182154943, 0.8103455408369723]
 [-1.7528722114632016, 2.6170584088128175]
 [0.696934271295315, 0.07381724447824105]

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}}:
 [32.67524753282632, 20.981552667304765]
 [-17.346020188978898, -13.876564296538538]
 [31.609602920302525, 22.2097562813014]
 [-2.5136467707366306, -15.038598996598916]
 [65.48531256490004, 20.29206829375256]
 [-65.46135872210793, -15.542852647967957]
 [-46.33008567076406, 20.057325851996744]
 [47.68127168761489, -22.24415178358279]
 [-13.017168554919735, 53.633036658612895]
 [66.0584969614256, -20.752650100246164]
 ⋮
 [47.82815785502888, -9.352445837513953]
 [-16.53159119215159, -19.611504482497732]
 [3.999294903265289, 0.3071512940550943]
 [3.211432050847673, 38.9306233794787]
 [-51.1176765580625, -12.420847341668125]
 [-15.58514541757918, -23.85149377058924]
 [20.932480378679852, 16.270789229312317]
 [-73.68400915688429, 52.54746725281174]
 [29.296437522392964, 1.4821638003382467]

finally, we can create the histogram.

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

• edges: (-184.0:1.0:177.0, -94.0:1.0:87.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: 33.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.