Geospatial data

Overview

Given a table or array containing data, we can georeference these objects onto a geospatial domain with the georef function. The result of the georef function is a GeoTable. If you would like to learn this concept in more depth, check out the recording of our JuliaEO 2024 workshop:

GeoTables.georef — Function

georef(table, domain)

Georeference table on domain from Meshes.jl

Examples

julia> georef((a=rand(100), b=rand(100)), CartesianGrid(10, 10))

source

georef(table, geoms)

Georeference table on vector of geometries geoms from Meshes.jl

Examples

julia> georef((a=rand(10), b=rand(10)), rand(Point, 10))

source

georef(table, coords; [crs])

Georeference table using coordinates coords of points.

Optionally, specify the coordinate reference system crs, which is set by default based on heuristics. Any CRS or EPSG/ESRI code from CoordRefSystems.jl is supported.

Examples

julia> georef((a=[1, 2, 3], b=[4, 5, 6], [(0, 0), (1, 1), (2, 2)])
julia> georef((a=[1, 2, 3], b=[4, 5, 6], [(0, 0), (1, 1), (2, 2)], crs=LatLon)
julia> georef((a=[1, 2, 3], b=[4, 5, 6], [(0, 0), (1, 1), (2, 2)], crs=EPSG{4326})

source

georef(table, names; [crs])

Georeference table using coordinates of points stored in column names.

Optionally, specify the coordinate reference system crs, which is set by default based on heuristics. Any CRS or EPSG/ESRI code from CoordRefSystems.jl is supported.

Examples

georef((a=rand(10), x=rand(10), y=rand(10)), ("x", "y"))
georef((a=rand(10), x=rand(10), y=rand(10)), ("x", "y"), crs=LatLon)
georef((a=rand(10), x=rand(10), y=rand(10)), ("x", "y"), crs=EPSG{4326})

source

georef(tuple)

Georeference a named tuple on CartesianGrid(dims), with dims obtained from the arrays stored in the tuple.

Examples

julia> georef((a=rand(10, 10), b=rand(10, 10))) # 2D grid
julia> georef((a=rand(10, 10, 10), b=rand(10, 10, 10))) # 3D grid

source

The functions values and domain can be used to retrieve the table of attributes and the underlying geospatial domain:

Base.values — Function

values(geotable, [rank])

Return the values of geotable for a given rank as a table.

The rank is a non-negative integer that specifies the parametric dimension of the geometries of interest:

0 - points
1 - segments
2 - triangles, quadrangles, ...
3 - tetrahedrons, hexahedrons, ...

If the rank is not specified, it is assumed to be the rank of the elements of the domain.

source

GeoTables.domain — Function

domain(geotable)

Return underlying domain of the geotable.

source

The GeoIO.jl package can be used to load/save geospatial data from/to various file formats:

GeoIO.load — Function

load(fname, repair=true, layer=0, lenunit=m, kwargs...)

Load geospatial table from file fname stored in any format.

Various repairs are performed on the stored geometries by default, including fixes of orientation in rings of polygons, removal of zero-area triangles, etc.

Some of the repairs can be expensive on large data sets. In that case, we recommend setting repair=false. Custom repairs can be performed with the Repair transform from Meshes.jl.

Optionally, specify the layer to read within the file, and the length unit lenunit of the coordinates when the format does not include units in its specification. Other kwargs are forwarded to the backend packages.

Please use the formats function to list all supported file formats.

Examples

# load coordinates of geojson file as Float64
GeoIO.load("file.geojson", numbertype = Float64)

source

GeoIO.save — Function

save(fname, geotable; kwargs...)

Save geotable to file fname of given format based on the file extension.

Other kwargs are forwarded to the backend packages.

Please use the formats function to list all supported file formats.

Examples

# overwrite an existing shapefile
GeoIO.save("file.shp", force = true)

source

GeoIO.formats — Function

formats([io]; sortby=:format)

Displays in io (defaults to stdout if io is not given) a table with all formats supported by GeoIO.jl and the packages used to load and save each of them.

Optionally, sort the table by the :extension, :load or :save columns using the sortby argument.

source

The GeoArtifacts.jl package provides utility functions to automatically download geospatial data from repositories on the internet.

Examples

using GeoStats
import CairoMakie as Mke

# helper function for plotting two
# variables named T and P side by side
function plot(data)
  fig = Mke.Figure(size = (800, 400))
  viz(fig[1,1], data.geometry, color = data.T)
  viz(fig[1,2], data.geometry, color = data.P)
  fig
end

plot (generic function with 1 method)

Tables

Consider a table (e.g. DataFrame) with 25 samples of temperature T and pressure P:

using DataFrames

table = DataFrame(T=rand(25), P=rand(25))

25×2 DataFrame

Row	T	P
	Float64	Float64
1	0.348329	0.217936
2	0.458625	0.171324
3	0.331716	0.299964
4	0.833794	0.805798
5	0.491734	0.831344
6	0.732932	0.882432
7	0.582012	0.847606
8	0.0634106	0.276897
9	0.154083	0.60053
10	0.128804	0.513346
11	0.68775	0.623933
12	0.729439	0.300338
13	0.969394	0.428641
14	0.737632	0.228709
15	0.915321	0.989978
16	0.684078	0.0203244
17	0.617727	0.726475
18	0.341877	0.579067
19	0.187245	0.940427
20	0.0288956	0.240846
21	0.119165	0.267264
22	0.330471	0.549888
23	0.257002	0.243371
24	0.603894	0.918379
25	0.378081	0.187031

We can georeference this table based on a given set of points:

georef(table, rand(Point, 25)) |> plot

or alternatively, georeference it on a 5x5 regular grid (5x5 = 25 samples):

georef(table, CartesianGrid(5, 5)) |> plot

Another common pattern in geospatial data is when the coordinates of the samples are already part of the table as columns. In this case, we can specify the column names:

table = DataFrame(T=rand(25), P=rand(25), X=rand(25), Y=rand(25), Z=rand(25))

georef(table, ("X", "Y", "Z")) |> plot

Arrays

Consider arrays (e.g. images) with data for various geospatial variables. We can georeference these arrays using a named tuple, and the framework will understand that the shape of the arrays should be preserved in a CartesianGrid:

T, P = rand(5, 5), rand(5, 5)

georef((T=T, P=P)) |> plot

Alternatively, we can interpret the entries of the named tuple as columns in a table:

georef((T=vec(T), P=vec(P)), rand(Point, 25)) |> plot

Files

We can easily load geospatial data from disk without any specific knowledge of file formats:

using GeoIO

zone = GeoIO.load("data/zone.shp")

4×6 GeoTable over 4 GeometrySet
PERIMETER	ACRES	MACROZONA	Hectares	area_m2	geometry
Continuous	Continuous	Categorical	Continuous	Continuous	MultiPolygon
[NoUnits]	[NoUnits]	[NoUnits]	[NoUnits]	[NoUnits]	🖈 GeodeticLatLon{SAD69}
5.8508e6	3.23145e7	Estuario	1.30772e7	1.30772e11	Multi(21×PolyArea)
9.53947e6	2.50594e8	Fronteiras Antigas	1.01412e8	1.01412e12	Multi(1×PolyArea)
1.01743e7	2.75528e8	Fronteiras Intermediarias	1.11502e8	1.11502e12	Multi(1×PolyArea)
7.09612e6	1.61293e8	Fronteiras Novas	6.5273e7	6.5273e11	Multi(2×PolyArea)