Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

fix GeoData export #134

Merged
merged 2 commits into from
Jul 29, 2024
Merged
Changes from 1 commit
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
33 changes: 19 additions & 14 deletions src/utils.jl
Original file line number Diff line number Diff line change
Expand Up @@ -70,7 +70,7 @@

Add `new_fields` fields to a `Q1Data` dataset; set `cellfield` to `true` if the field is a cell field; otherwise it is a vertex field
"""
function addfield(V::Q1Data,new_fields::NamedTuple; cellfield=false)
function addfield(V::Q1Data,new_fields::NamedTuple; cellfield=false)
if cellfield
return Q1Data(V.x.val, V.y.val, V.z.val, V.fields, merge(V.cellfields, new_fields))
else
Expand Down Expand Up @@ -122,7 +122,7 @@



"""
"""
cross_section_volume(Volume::AbstractGeneralGrid; dims=(100,100), Interpolate=false, Depth_level=nothing; Lat_level=nothing; Lon_level=nothing; Start=nothing, End=nothing, Depth_extent=nothing )

Creates a cross-section through a volumetric (3D) `GeoData` object.
Expand Down Expand Up @@ -259,7 +259,7 @@
Creates a cross-section through a surface (2D) `GeoData` object.

- Cross-sections can be horizontal (map view at a given depth), if `Depth_level` is specified
- They can also be vertical, either by specifying `Lon_level` or `Lat_level` (for a fixed lon/lat), or by defining both `Start=(lon,lat)` & `End=(lon,lat)` points.
- They can also be vertical, either by specifying `Lon_level` or `Lat_level` (for a fixed lon/lat), or by defining both `Start=(lon,lat)` & `End=(lon,lat)` points. Start and End points will be in km!

- IMPORTANT: The surface to be extracted has to be given as a gridded GeoData object. It may also contain NaNs where it is not defined. Any points lying outside of the defined surface will be considered NaN.

Expand Down Expand Up @@ -351,9 +351,14 @@

end

# create GeoData structure with the interpolated points
Data_profile = GeoData(Lon, Lat, depth_intp, (fields_new));

# create GeoData/CartData structure with the interpolated points
if isa(S,GeoData)
Data_profile = GeoData(Lon, Lat, depth_intp, fields_new);
elseif isa(S,CartData)
Data_profile = CartData(Lon, Lat, depth_intp, fields_new);

Check warning on line 358 in src/utils.jl

View check run for this annotation

Codecov / codecov/patch

src/utils.jl#L357-L358

Added lines #L357 - L358 were not covered by tests
else
error("still to be implemented")

Check warning on line 360 in src/utils.jl

View check run for this annotation

Codecov / codecov/patch

src/utils.jl#L360

Added line #L360 was not covered by tests
end
return Data_profile
end

Expand Down Expand Up @@ -545,9 +550,9 @@
julia> Lon,Lat,Depth = lonlatdepth_grid(10:20,30:40,(-300:25:0)km);
julia> Data = Depth*2; # some data
julia> Vx,Vy,Vz = ustrip(Data*3),ustrip(Data*4),ustrip(Data*5);
julia> Data_set3D = GeoData(Lon,Lat,Depth,(Depthdata=Data,LonData=Lon, Velocity=(Vx,Vy,Vz)));
julia> Data_cross = cross_section(Data_set3D, Depth_level=-100km)
GeoData
julia> Data_set3D = GeoData(Lon,Lat,Depth,(Depthdata=Data,LonData=Lon, Velocity=(Vx,Vy,Vz)));
julia> Data_cross = cross_section(Data_set3D, Depth_level=-100km)
GeoData
size : (11, 11, 1)
lon ϵ [ 10.0 : 20.0]
lat ϵ [ 30.0 : 40.0]
Expand Down Expand Up @@ -601,7 +606,7 @@
```
"""
function flatten_cross_section(V::CartData)

x_new = sqrt.((V.x.val.-V.x.val[1,1,1]).^2 .+ (V.y.val.-V.y.val[1,1,1]).^2) # NOTE: the result is in km, as V.x and V.y are stored in km


Expand Down Expand Up @@ -1244,11 +1249,11 @@
return GeoData(New.lon.val,New.lat.val,Znew, fields_new)
end

"""
"""
Surf_interp = interpolate_datafields_2D(V::GeoData, x::AbstractRange, y::AbstractRange; Lat::Number, Lon::Number)

Interpolates a 3D data set `V` with a projection point `proj=(Lat, Lon)` on a plane defined by `x` and `y`, where `x` and `y` are uniformly spaced.
Returns the 2D array `Surf_interp`.
Returns the 2D array `Surf_interp`.
"""
function interpolate_datafields_2D(V::GeoData, x::AbstractRange, y::AbstractRange; Lat=49.9929, Lon=8.2473)
# Default: Lat=49.9929, Lon=8.2473 => Mainz (center of universe)
Expand Down Expand Up @@ -1340,7 +1345,7 @@

return depth_new, fields_new
end

# Extracts a sub-data set using indices
function ExtractDataSets(V::AbstractGeneralGrid, iLon, iLat, iDepth)

Expand Down Expand Up @@ -1815,4 +1820,4 @@
end
end
return inside
end
end
Loading