Data structures

The main data structure used in GeophysicalModelGenerator.jl is GeoData, which contains info about the longitude,latitude, and depth of a data set, as well as several data sets itself. We also provide a UTMData, which is essentially the same but with UTM coordinates, and a CartData structure, which has Cartesian coordinates in kilometers (as used in many geodynamic codes). If one wishes to transfer GeoData to CartData, one needs to provide a ProjectionPoint. For plotting, we transfer this into the ParaviewData structure, which has cartesian coordinates centered around the center of the Earth. We employ the wgs84 reference ellipsoid as provided by the Geodesy.jl package to perform this transformation.

GeophysicalModelGenerator.GeoData — Type

GeoData(lon::Any, lat:Any, depth::GeoUnit, fields::NamedTuple)

Data structure that holds one or several fields with longitude, latitude and depth information.

depth can have units of meter, kilometer or be unitless; it will be converted to km.
fields should ideally be a NamedTuple which allows you to specify the names of each of the fields.
In case you only pass one array we will convert it to a NamedTuple with default name.
A single field should be added as (DataFieldName=Data,) (don't forget the comma at the end).
Multiple fields can be added as well. lon,lat,depth should all have the same size as each of the fields.
In case you want to display a vector field in paraview, add it as a tuple: (Velocity=(Veast,Vnorth,Vup), Veast=Veast, Vnorth=Vnorth, Vup=Vup); we automatically apply a vector transformation when transforming this to a ParaviewData structure from which we generate Paraview output. As this changes the magnitude of the arrows, you will no longer see the [Veast,Vnorth,Vup] components in Paraview which is why it is a good ideas to store them as separate Fields.
Yet, there is one exception: if the name of the 3-component field is colors, we do not apply this vector transformation as this field is regarded to contain RGB colors.
Lat,Lon,Depth should have the same size as the Data array. The ordering of the arrays is important. If they are 3D arrays, as in the example below, we assume that the first dimension corresponds to lon, second dimension to lat and third dimension to depth (which should be in km). See below for an example.

Example

julia> Lat         =   1.0:3:10.0;
julia> Lon         =   11.0:4:20.0;
julia> Depth       =   (-20:5:-10)*km;
julia> Lon3D,Lat3D,Depth3D = lonlatdepth_grid(Lon, Lat, Depth);
julia> Lon3D
3×4×3 Array{Float64, 3}:
[:, :, 1] =
 11.0  11.0  11.0  11.0
 15.0  15.0  15.0  15.0
 19.0  19.0  19.0  19.0

[:, :, 2] =
 11.0  11.0  11.0  11.0
 15.0  15.0  15.0  15.0
 19.0  19.0  19.0  19.0

[:, :, 3] =
 11.0  11.0  11.0  11.0
 15.0  15.0  15.0  15.0
 19.0  19.0  19.0  19.0
julia> Lat3D
 3×4×3 Array{Float64, 3}:
 [:, :, 1] =
  1.0  4.0  7.0  10.0
  1.0  4.0  7.0  10.0
  1.0  4.0  7.0  10.0

 [:, :, 2] =
  1.0  4.0  7.0  10.0
  1.0  4.0  7.0  10.0
  1.0  4.0  7.0  10.0

 [:, :, 3] =
  1.0  4.0  7.0  10.0
  1.0  4.0  7.0  10.0
  1.0  4.0  7.0  10.0
julia> Depth3D
  3×4×3 Array{Unitful.Quantity{Float64, 𝐋, Unitful.FreeUnits{(km,), 𝐋, nothing}}, 3}:
  [:, :, 1] =
   -20.0 km  -20.0 km  -20.0 km  -20.0 km
   -20.0 km  -20.0 km  -20.0 km  -20.0 km
   -20.0 km  -20.0 km  -20.0 km  -20.0 km

  [:, :, 2] =
   -15.0 km  -15.0 km  -15.0 km  -15.0 km
   -15.0 km  -15.0 km  -15.0 km  -15.0 km
   -15.0 km  -15.0 km  -15.0 km  -15.0 km

  [:, :, 3] =
   -10.0 km  -10.0 km  -10.0 km  -10.0 km
   -10.0 km  -10.0 km  -10.0 km  -10.0 km
   -10.0 km  -10.0 km  -10.0 km  -10.0 km
julia> Data        =   zeros(size(Lon3D));
julia> Data_set    =   GeophysicalModelGenerator.GeoData(Lon3D,Lat3D,Depth3D,(DataFieldName=Data,))
GeoData
  size      : (3, 4, 3)
  lon       ϵ [ 11.0 : 19.0]
  lat       ϵ [ 1.0 : 10.0]
  depth     ϵ [ -20.0 km : -10.0 km]
  fields    : (:DataFieldName,)
  attributes: ["note"]

GeophysicalModelGenerator.CartData — Type

CartData(x::Any, y::Any, z::GeoUnit, fields::NamedTuple)

Data structure that holds one or several fields with with Cartesian x/y/z coordinates. Distances are in kilometers

x,y,z can have units of meters, kilometer or be unitless; they will be converted to kilometers
fields should ideally be a NamedTuple which allows you to specify the names of each of the fields.
In case you only pass one array we will convert it to a NamedTuple with default name.
A single field should be added as (DataFieldName=Data,) (don't forget the comma at the end).
Multiple fields can be added as well.
In case you want to display a vector field in paraview, add it as a tuple: (Velocity=(Vx,Vnorth,Vup), Veast=Veast, Vnorth=Vnorth, Vup=Vup); we automatically apply a vector transformation when transforming this to a ParaviewData structure from which we generate Paraview output. As this changes the magnitude of the arrows, you will no longer see the [Veast,Vnorth,Vup] components in Paraview which is why it is a good ideas to store them as separate Fields.
Yet, there is one exception: if the name of the 3-component field is colors, we do not apply this vector transformation as this field is regarded to contain RGB colors.
x,y,z should have the same size as the Data array. The ordering of the arrays is important. If they are 3D arrays, as in the example below, we assume that the first dimension corresponds to x, second dimension to y and third dimension to z (which should be in km). See below for an example.

Example

julia> x        =   0:2:10
julia> y        =   -5:5
julia> z        =   -10:2:2
julia> X,Y,Z    =   xyz_grid(x, y, z);
julia> Data     =   Z
julia> Data_set =   CartData(X,Y,Z, (FakeData=Data,Data2=Data.+1.))
CartData
    size    : (6, 11, 7)
    x       ϵ [ 0.0 km : 10.0 km]
    y       ϵ [ -5.0 km : 5.0 km]
    z       ϵ [ -10.0 km : 2.0 km]
    fields  : (:FakeData, :Data2)
  attributes: ["note"]

CartData is particularly useful in combination with cartesian geodynamic codes, such as LaMEM, which require cartesian grids. You can directly save your data to Paraview with

julia> write_paraview(Data_set, "Data_set")
1-element Vector{String}:
 "Data_set.vts"

If you wish, you can convert this to UTMData (which will simply convert the )

julia> Data_set1 =  convert(GeoData, Data_set)
GeoData
  size  : (116, 96, 25)
  lon   ϵ [ 14.075969111533457 : 14.213417764154963]
  lat   ϵ [ 40.77452227533946 : 40.86110443583479]
  depth ϵ [ -5.4 km : 0.6 km]
  fields: (:FakeData, :Data2)

which would allow visualizing this in paraview in the usual manner:

GeophysicalModelGenerator.UTMData — Type

UTMData(EW::Any, NS:Any, depth::GeoUnit, UTMZone::Int, NorthernHemisphere=true, fields::NamedTuple)

Data structure that holds one or several fields with UTM coordinates (east-west), (north-south) and depth information.

depth can have units of meters, kilometer or be unitless; it will be converted to meters (as UTMZ is usually in meters)
fields should ideally be a NamedTuple which allows you to specify the names of each of the fields.
In case you only pass one array we will convert it to a NamedTuple with default name.
A single field should be added as (DataFieldName=Data,) (don't forget the comma at the end).
Multiple fields can be added as well.
In case you want to display a vector field in paraview, add it as a tuple: (Velocity=(Veast,Vnorth,Vup), Veast=Veast, Vnorth=Vnorth, Vup=Vup); we automatically apply a vector transformation when transforming this to a ParaviewData structure from which we generate Paraview output. As this changes the magnitude of the arrows, you will no longer see the [Veast,Vnorth,Vup] components in Paraview which is why it is a good ideas to store them as separate Fields.
Yet, there is one exception: if the name of the 3-component field is colors, we do not apply this vector transformation as this field is regarded to contain RGB colors.
Lat,Lon,Depth should have the same size as the Data array. The ordering of the arrays is important. If they are 3D arrays, as in the example below, we assume that the first dimension corresponds to lon, second dimension to lat and third dimension to depth (which should be in km). See below for an example.

Example

julia> ew          =   422123.0:100:433623.0
julia> ns          =   4.514137e6:100:4.523637e6
julia> depth       =   -5400:250:600
julia> EW,NS,Depth =   xyz_grid(ew, ns, depth);
julia> Data        =   ustrip.(Depth);
julia> Data_set    =   UTMData(EW,NS,Depth,33, true, (FakeData=Data,Data2=Data.+1.))
UTMData
  UTM zone : 33-33 North
    size    : (116, 96, 25)
    EW      ϵ [ 422123.0 : 433623.0]
    NS      ϵ [ 4.514137e6 : 4.523637e6]
    depth   ϵ [ -5400.0 m : 600.0 m]
    fields  : (:FakeData, :Data2)
  attributes: ["note"]

If you wish, you can convert this from UTMData to GeoData with

julia> Data_set1 =  convert(GeoData, Data_set)
GeoData
  size      : (116, 96, 25)
  lon       ϵ [ 14.075969111533457 : 14.213417764154963]
  lat       ϵ [ 40.77452227533946 : 40.86110443583479]
  depth     ϵ [ -5.4 km : 0.6 km]
  fields    : (:FakeData, :Data2)
  attributes: ["note"]

which would allow visualizing this in paraview in the usual manner:

julia> write_paraview(Data_set1, "Data_set1")
1-element Vector{String}:
 "Data_set1.vts"

GeophysicalModelGenerator.Q1Data — Type

Holds a Q1 Finite Element Data set with vertex and cell data. The specified coordinates are the ones of the vertices.

GeophysicalModelGenerator.FEData — Type

FEData{dim, points_per_cell}

Structure that holds Finite Element info with vertex and cell data. Works in 2D/3D for arbitrary elements

Parameters

vertices with the points on the mesh (dim x Npoints)
connectivity with the connectivity of the mesh (points_per_cell x Ncells)
fields with the fields on the vertices
cellfields with the fields of the cells

GeophysicalModelGenerator.ParaviewData — Type

ParaviewData(x::GeoUnit, y::GeoUnit, z::GeoUnit, values::NamedTuple)

Cartesian data in x/y/z coordinates to be used with Paraview. This is usually generated automatically from the GeoData structure, but you can also invoke do this manually:

julia> Data_set    =   GeophysicalModelGenerator.GeoData(1.0:10.0,11.0:20.0,(-20:-11)*km,(DataFieldName=(-20:-11),))
julia> Data_cart = convert(ParaviewData, Data_set)

GeophysicalModelGenerator.lonlatdepth_grid — Function

Lon, Lat, Depth = lonlatdepth_grid(Lon::Any, Lat::Any, Depth:Any)

Creates 3D arrays of Lon, Lat, Depth from 1D vectors or numbers

Example 1: Create 3D grid

julia> Lon,Lat,Depth =  lonlatdepth_grid(10:20,30:40,(-10:-1)km);
julia> size(Lon)
(11, 11, 10)

Example 2: Create 2D lon/lat grid @ a given depth

julia> Lon,Lat,Depth =  lonlatdepth_grid(10:20,30:40,-50km);
julia> size(Lon)
(11, 11)

Example 3: Create 2D lon/depth grid @ a given lat

julia> Lon,Lat,Depth =  lonlatdepth_grid(10:20,30,(-10:-1)km);
julia> size(Lon)
(11, 11)

Example 4: Create 1D vertical line @ a given lon/lat point

julia> Lon,Lat,Depth =  lonlatdepth_grid(10,30,(-10:-1)km);
julia> size(Lon)
(10, )

GeophysicalModelGenerator.xyz_grid — Function

X,Y,Z = xyz_grid(X_vec::Any, Y_vec::Any, Z_vec::Any)

Creates a X,Y,Z grid. It works just as lonlatdepth_grid apart from the better suited name.

Example 1: Create 3D grid

julia> X,Y,Z =  xyz_grid(10:20,30:40,(-10:-1)km);
julia> size(X)
(11, 11, 10)

See lonlatdepth_grid for more examples.

GeophysicalModelGenerator.ProjectionPoint — Type

struct ProjectionPoint
    Lon     :: Float64
    Lat     :: Float64
    EW      :: Float64
    NS      :: Float64
    zone    :: Integer
    isnorth :: Bool
end

Structure that holds the coordinates of a point that is used to project a data set from Lon/Lat to a Cartesian grid and vice-versa.