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.6866185720920313, 0.29854317292192306]
 [-0.5745471895212764, 1.8227043838752723]
 [0.31764510524281353, -0.4809972437645058]
 [1.2543373769468176, 1.1904486639771292]
 [-0.16577465544186137, 1.004573509818853]
 [-1.1856321421634572, 0.14228292215523003]
 [-0.022791240797673713, 1.4758393415713158]
 [-0.5368067625523009, 1.9694013716560488]
 [1.7935673474305638, 0.6698362599531308]
 [0.9585256518308629, 0.7963661343023223]
 ⋮
 [-1.102884650047637, 0.7349959248036957]
 [0.312162902823187, -1.383884192858594]
 [1.0582293548342525, 0.8537658835061248]
 [-1.0866541508840692, 0.2973115166088481]
 [0.6477211054732739, -0.9242447780868199]
 [0.9027228921020055, 0.10421995256502707]
 [0.3818197800490031, -0.2700404671245361]
 [0.6543876928075286, 0.5379545515110752]
 [0.7214311459367243, 0.36538105832091866]

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}}:
 [-28.37757814468733, 5.7412900387841175]
 [-23.74572787737695, 35.0524663497431]
 [13.128102213782205, -9.250068113364597]
 [51.841092537967995, 22.893543299893885]
 [-6.8513777960837, 19.31897430010147]
 [-49.0015419515726, 2.7362458691110514]
 [-0.9419497855698592, 28.38189742434699]
 [-22.185936227315935, 37.87359920775581]
 [74.12717865232514, 12.881635206212126]
 [39.61535224082436, 15.314933881575248]
 ⋮
 [-45.581632384257766, 14.13472209168881]
 [12.901525721547307, -26.613506025032255]
 [43.736052930106254, 16.418789666007083]
 [-44.91083453948334, 5.717604030317556]
 [26.769966665083768, -17.77417077032222]
 [37.30905404374819, 2.0042561002077557]
 [15.780407180823564, -5.19315390399557]
 [27.045493152023653, 10.345415296824932]
 [29.81636319194507, 7.026650818931667]

finally, we can create the histogram.

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

• edges: ([-182.0, -181.0, -180.0, -179.0, -178.0, -177.0, -176.0, -175.0, -174.0, -173.0 … 170.0, 171.0, 172.0, 173.0, 174.0, 175.0, 176.0, 177.0, 178.0, 179.0], [-86.0, -85.0, -84.0, -83.0, -82.0, -81.0, -80.0, -79.0, -78.0, -77.0 … 86.0, 87.0, 88.0, 89.0, 90.0, 91.0, 92.0, 93.0, 94.0, 95.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: 36.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.