Skip to content

Commit

Permalink
Add cubeplot! methods.
Browse files Browse the repository at this point in the history
  • Loading branch information
joa-quim committed Feb 11, 2024
1 parent 869ea35 commit ca6c53e
Show file tree
Hide file tree
Showing 2 changed files with 20 additions and 8 deletions.
3 changes: 3 additions & 0 deletions src/imshow.jl
Original file line number Diff line number Diff line change
Expand Up @@ -96,8 +96,11 @@ end
# If arg1 is a cube:
# - `col cols columns`: set the number of columns in the mosaic
# - `slice`: set the single slice to be plotted. `slice` is the same option as in cubeslice.
# - `title`: The subplot title
# - `titles`: A string vector with pannel titles. Default is the layer names.
# - `cmap=:same`: uses the same CPT (computed from cube's min/max) for all layers.
# - `T, no_interp, tiles`: -T option for grdview
# - `facades, cubeplot`: Call cubeplot.
function imshow(arg1::GMTgrid; kw...)
# Here the default is to show, but if a 'show' was used let it rule
d = KW(kw)
Expand Down
25 changes: 17 additions & 8 deletions src/utils_project.jl
Original file line number Diff line number Diff line change
Expand Up @@ -464,7 +464,7 @@ The `kw...` keyword/value options may be used to pass:
are the axes colorbar labels taking the form: `cmap=(C, "xlabel=Blabla1"[, "ylabel=Blabla2"])`.
"""
function cubeplot(fname1::Union{GMTimage, String}, fname2::Union{GMTimage, String}="", fname3::Union{GMTimage, String}="";
back::Bool=false, notop::Bool=false, xlabel="", ylabel="", zlabel="", title="", show=false, kw...)
back::Bool=false, notop::Bool=false, xlabel="", ylabel="", zlabel="", title="", first=true, show=false, kw...)
# ...
d = KW(kw)
opt_R = ((txt::String = parse_R(d, "")[2]) != "") ? txt[4:end] : "0/9/0/9/-9/0"
Expand All @@ -481,7 +481,9 @@ function cubeplot(fname1::Union{GMTimage, String}, fname2::Union{GMTimage, Strin
else f2, f3 = fname2, fname2 # Two images, repeat second on the vert sides
end

basemap(R=opt_R, J=opt_J, JZ=opt_JZ, p=opt_p, xlabel=xlabel, ylabel=ylabel, zlabel=zlabel, title=title)
opt_UVXY = parse_UVXY(d, "")
first ? basemap(R=opt_R, J=opt_J, JZ=opt_JZ, p=opt_p, xlabel=xlabel, ylabel=ylabel, zlabel=zlabel, title=title, compact=opt_UVXY) :
basemap!(R=opt_R, J=opt_J, JZ=opt_JZ, p=opt_p, xlabel=xlabel, ylabel=ylabel, zlabel=zlabel, title=title, compact=opt_UVXY, B=:autoXYZ)
if ((val = find_in_dict(d, CPTaliases)[1]) !== nothing) # See if want to plot a colorbar.
if isa(val, GMTcpt)
colorbar!(val, position=(outside=true, anchor=:BC, offset="0.9c"), B=:af)
Expand Down Expand Up @@ -523,14 +525,17 @@ function cubeplot(fname1::Union{GMTimage, String}, fname2::Union{GMTimage, Strin
hole_width, hole_depth = hole_width * width, hole_depth * height
xoff = (width - hole_width) * cosd(azim -90)
td = (hole_depth > 0.0) ? "$hole_depth" : "0" # Found that "/0.0" screws up
basemap!(R=opt_R, J="X$(hole_width)/"*td, JZ=opt_JZ, p=opt_p, X=xoff, B=0)
basemap!(R=opt_R, J="X$(hole_width)/"*td, JZ=opt_JZ, p=opt_p, X=xoff, B=:none, title=".") # Need that "." so -B0 plots no axes.
image!(f1, compact=sideplot(plane=:B, vsize=vsz), t=opt_t)
image!(val[1][1], compact=sideplot(plane=:N, vsize=vsz), t=opt_t)
R = image!(val[1][2], compact=sideplot(plane=:W, vsize=vsz), t=opt_t)
end
see && showfig(d...)
return R
end
cubeplot!(fname1::Union{GMTimage, String}, fname2::Union{GMTimage, String}="", fname3::Union{GMTimage, String}="";
back::Bool=false, notop::Bool=false, xlabel="", ylabel="", zlabel="", title="", show=false, kw...) =
cubeplot(fname1, fname2, fname3; back=back, notop=notop, xlabel=xlabel, ylabel=ylabel, zlabel=zlabel, title=title, first=false, show=show, kw...)

# -----------------------------------------------------------------------------------------------
"""
Expand Down Expand Up @@ -654,7 +659,7 @@ cubeplot(C, top="@earth_relief", inset=(lon=110,y=40), topshade=true, zdown=true
```
"""
function cubeplot(G::GMTgrid; top=nothing, topshade=false, zdown::Bool=false, xlabel="", ylabel="", zlabel="", title="",
show=false, interp::Float64=0.0, kw...)
show=false, interp::Float64=0.0, first=true, kw...)
(size(G,3) < 2) && error("First input must be a 3D grid.")
d = KW(kw)
if ((opt_R = parse_R(d, "", false)[2]) != "")
Expand Down Expand Up @@ -724,15 +729,19 @@ function cubeplot(G::GMTgrid; top=nothing, topshade=false, zdown::Bool=false, xl
Ge = interp_vslice(squeeze(slicecube(G, inset[1], y=[G.range[3] inset[2]], axis="x")), inc=interp) # inset East wall
pct_x = (G.range[2] - inset[1]) / (G.range[2] - G.range[1])
pct_y = (inset[2] - G.range[3]) / (G.range[4] - G.range[3])
resetGMT() # To remove "rememberings" left by the grdimage calls (they make annotations on West side??)
gmt_restart() # To remove "rememberings" left by the grdimage calls (they make annotations on West side??)
cubeplot(It, Is, Ie, R=@sprintf("%.12g/%.12g/%.12g/%.12g/%.12g/%.12g", G.range[1:4]..., G.range[7], G.range[8]),
hole=((grdimage(Gs, A=true, B=:none, layout="BRP", C=C), grdimage(Ge, A=true, B=:none, layout="BRP", C=C)), pct_x, pct_y), zsize=zsize, xlabel=xlabel, ylabel=ylabel, zlabel=zlabel, title=title, C=opt_cbar, show=show)
hole=((grdimage(Gs, A=true, B=:none, layout="BRP", C=C), grdimage(Ge, A=true, B=:none, layout="BRP", C=C)), pct_x, pct_y), zsize=zsize, xlabel=xlabel, ylabel=ylabel, zlabel=zlabel, title=title, C=opt_cbar, first=first, show=show; d...)
else
resetGMT()
gmt_restart()
cubeplot(It, Is, Ie, R=@sprintf("%.12g/%.12g/%.12g/%.12g/%.12g/%.12g", G.range[1:4]..., G.range[7], G.range[8]),
zsize=zsize, xlabel=xlabel, ylabel=ylabel, zlabel=zlabel, title=title, C=opt_cbar, show=show)
zsize=zsize, xlabel=xlabel, ylabel=ylabel, zlabel=zlabel, title=title, C=opt_cbar, first=first, show=show; d...)
end
end
cubeplot!(G::GMTgrid; top=nothing, topshade=false, zdown::Bool=false, xlabel="", ylabel="", zlabel="", title="",
interp::Float64=0.0, show=false, kw...) =
cubeplot(G; top=top, topshade=topshade, zdown=zdown, xlabel=xlabel, ylabel=ylabel, zlabel=zlabel, title=title,
interp=interp, first=false, show=show, kw...)

# From https://stackoverflow.com/questions/54879412/get-the-mapping-from-each-element-of-input-to-the-bin-of-the-histogram-in-julia
# get the second output as in Matlab's histc
Expand Down

0 comments on commit ca6c53e

Please sign in to comment.