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}}:
 [1.0691826533903463, 1.1843784508155402]
 [-1.0179912082623859, -0.0991419113934644]
 [-0.1069605292678584, -0.6835346965648788]
 [-1.6702585044142932, -0.25240248019280226]
 [1.1225699704503154, -0.24235285744539978]
 [0.13228539817968407, 0.22969465396128383]
 [0.20188997016111954, 0.8420038241624407]
 [0.05832036989163219, -2.458614896831935]
 [1.1021155104695395, -0.4851042522820332]
 [-0.942029730398816, 0.673807854879266]
 ⋮
 [2.239385826283209, 0.5267653052643719]
 [-1.2866492163211756, -1.5591499667765447]
 [-0.41888186712520525, -0.32379148978801175]
 [-0.050195454656026865, 0.7835987204183452]
 [-0.5599777626857717, -0.8420140703495678]
 [-1.2651413885743987, 0.07902309441412769]
 [-0.1688207628411023, -0.7714799899878415]
 [-1.7290721670107996, -1.1691961685687178]
 [0.4961611288835051, -1.0467568241353409]

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}}:
 [42.637912322076, 25.53850509142735]
 [-40.59644976926982, -2.1377763223844433]
 [-4.265476674527356, -14.73891585613469]
 [-66.60820341650096, -5.442501947797338]
 [44.766944005012185, -5.225804032089115]
 [5.275406583889876, 4.952857834946791]
 [8.051165831338652, 18.15595254672768]
 [2.325756792018027, -53.014599359998456]
 [43.95124102994968, -10.460201642681906]
 [-37.567183607186145, 14.529178001025556]
 ⋮
 [89.30441979542995, 11.358530224201655]
 [-51.31025676558681, -33.619624991828765]
 [-16.704581080844957, -6.981848246916353]
 [-2.0017434699336087, 16.896575496843962]
 [-22.33134130262704, -18.156173483120043]
 [-50.45254655976047, 1.703958475135965]
 [-6.732399615107615, -16.63526184444853]
 [-68.95363221781771, -25.211132711262064]
 [19.78639911887309, -22.570998707599898]

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], [-93.0, -92.0, -91.0, -90.0, -89.0, -88.0, -87.0, -86.0, -85.0, -84.0 … 79.0, 80.0, 81.0, 82.0, 83.0, 84.0, 85.0, 86.0, 87.0, 88.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.