From c7b76d411684b329bb29266470e2408aa98e192e Mon Sep 17 00:00:00 2001 From: Joaquim Date: Thu, 25 Apr 2024 23:47:47 +0100 Subject: [PATCH 1/4] Need to cast the UInts floats in UInts GMTimage algebras --- src/grd_operations.jl | 40 +++++++++++++++++++++++++++++++--------- 1 file changed, 31 insertions(+), 9 deletions(-) diff --git a/src/grd_operations.jl b/src/grd_operations.jl index 4c8f9a235..96ffaa1e1 100644 --- a/src/grd_operations.jl +++ b/src/grd_operations.jl @@ -11,7 +11,7 @@ function Base.:+(G1::GMTgrid, G2::GMTgrid) epsg, geog, range, inc, registration, nodata, x, y, v, pad = dup_G_meta(G1) hasnans = (G1.hasnans == 0 && G2.hasnans == 0) ? 0 : ((G1.hasnans + G2.hasnans) == 2) ? 1 : 2 (eltype(G1.z) <: AbstractFloat && eltype(G2.z) <: AbstractFloat) ? (z = G1.z .+ G2.z) : - (z = Matrix{Float32}(undef, size(G1.z)); @inbounds for k = 1:numel(G1.z) z[k] = G1.z[k] + G2.z[k] end) + (z = Matrix{Float32}(undef, size(G1.z)); @inbounds for k = 1:numel(G1.z) z[k] = Float32(G1.z[k]) + Float32(G2.z[k]) end) G3 = GMTgrid(G1.proj4, G1.wkt, epsg, geog, range, inc, registration, nodata, "", "", "", "", G1.names, x, y, v, z, G1.x_unit, G1.y_unit, G1.v_unit, G1.z_unit, G1.layout, 1f0, 0f0, pad, hasnans) setgrdminmax!(G3) # Also take care of NaNs @@ -26,7 +26,7 @@ function Base.:+(G1::GMTgrid, shift::Real) epsg, geog, range, inc, registration, nodata, x, y, v, pad = dup_G_meta(G1) _shift = convert(eltype(G1.z), shift) (eltype(G1.z) <: AbstractFloat) ? (z = G1.z .+ _shift) : - (z = Matrix{Float32}(undef, size(G1.z)); @inbounds for k = 1:numel(G1.z) z[k] = G1.z[k] + _shift end) + (z = Matrix{Float32}(undef, size(G1.z)); @inbounds for k = 1:numel(G1.z) z[k] = Float32(G1.z[k]) + _shift end) G2 = GMTgrid(G1.proj4, G1.wkt, epsg, geog, range, inc, registration, nodata, "", "", "", "", G1.names, x, y, v, z, G1.x_unit, G1.y_unit, G1.v_unit, G1.z_unit, G1.layout, 1f0, 0f0, pad, G1.hasnans) G2.range[5:6] .+= shift @@ -43,7 +43,7 @@ function Base.:-(G1::GMTgrid, G2::GMTgrid) epsg, geog, range, inc, registration, nodata, x, y, v, pad = dup_G_meta(G1) hasnans = (G1.hasnans == 0 && G2.hasnans == 0) ? 0 : ((G1.hasnans + G2.hasnans) == 2) ? 1 : 2 (eltype(G1.z) <: AbstractFloat && eltype(G2.z) <: AbstractFloat) ? (z = G1.z .- G2.z) : - (z = Matrix{Float32}(undef, size(G1.z)); @inbounds for k = 1:numel(G1.z) z[k] = G1.z[k] - G2.z[k] end) + (z = Matrix{Float32}(undef, size(G1.z)); @inbounds for k = 1:numel(G1.z) z[k] = Float32(G1.z[k]) - Float32(G2.z[k]) end) G3 = GMTgrid(G1.proj4, G1.wkt, epsg, geog, range, inc, registration, nodata, "", "", "", "", G1.names, x, y, v, z, G1.x_unit, G1.y_unit, G1.v_unit, G1.z_unit, G1.layout, 1f0, 0f0, pad, hasnans) setgrdminmax!(G3) # Also take care of NaNs @@ -57,7 +57,7 @@ function Base.:-(G1::GMTgrid, shift::Real) _shift = convert(eltype(G1.z), shift) epsg, geog, range, inc, registration, nodata, x, y, v, pad = dup_G_meta(G1) (eltype(G1.z) <: AbstractFloat) ? (z = G1.z .- _shift) : - (z = Matrix{Float32}(undef, size(G1.z)); @inbounds for k = 1:numel(G1.z) z[k] = G1.z[k] - _shift end) + (z = Matrix{Float32}(undef, size(G1.z)); @inbounds for k = 1:numel(G1.z) z[k] = Float32(G1.z[k]) - _shift end) G2 = GMTgrid(G1.proj4, G1.wkt, epsg, geog, range, inc, registration, nodata, "", "", "", "", G1.names, x, y, v, z, G1.x_unit, G1.y_unit, G1.v_unit, G1.z_unit, G1.layout, 1f0, 0f0, pad, G1.hasnans) G2.range[5:6] .-= shift @@ -74,7 +74,7 @@ function Base.:*(G1::GMTgrid, G2::GMTgrid) epsg, geog, range, inc, registration, nodata, x, y, v, pad = dup_G_meta(G1) hasnans = (G1.hasnans == 0 && G2.hasnans == 0) ? 0 : ((G1.hasnans + G2.hasnans) == 2) ? 1 : 2 (eltype(G1.z) <: AbstractFloat && eltype(G2.z) <: AbstractFloat) ? (z = G1.z .* G2.z) : - (z = Matrix{Float32}(undef, size(G1.z)); @inbounds for k = 1:numel(G1.z) z[k] = G1.z[k] * G2.z[k] end) + (z = Matrix{Float32}(undef, size(G1.z)); @inbounds for k = 1:numel(G1.z) z[k] = Float32(G1.z[k]) * Float32(G2.z[k]) end) G3 = GMTgrid(G1.proj4, G1.wkt, epsg, geog, range, inc, registration, nodata, "", "", "", "", G1.names, x, y, v, z, G1.x_unit, G1.y_unit, G1.v_unit, G1.z_unit, G1.layout, 1f0, 0f0, pad, hasnans) setgrdminmax!(G3) # Also take care of NaNs @@ -90,7 +90,7 @@ function Base.:*(G1::GMTgrid, scale::Real) _scale = convert(eltype(G1.z), scale) epsg, geog, range, inc, registration, nodata, x, y, v, pad = dup_G_meta(G1) (eltype(G1.z) <: AbstractFloat) ? (z = G1.z .* _scale) : - (z = Matrix{Float32}(undef, size(G1.z)); @inbounds for k = 1:numel(G1.z) z[k] = G1.z[k] * _scale end) + (z = Matrix{Float32}(undef, size(G1.z)); @inbounds for k = 1:numel(G1.z) z[k] = Float32(G1.z[k]) * _scale end) G2 = GMTgrid(G1.proj4, G1.wkt, epsg, geog, range, inc, registration, nodata, "", "", "", "", G1.names, x, y, v, z, G1.x_unit, G1.y_unit, G1.v_unit, G1.z_unit, G1.layout, 1f0, 0f0, pad, G1.hasnans) G2.range[5:6] .*= scale @@ -106,7 +106,7 @@ function Base.:^(G1::GMTgrid, scale::AbstractFloat) _scale = convert(eltype(G1.z), scale) epsg, geog, range, inc, registration, nodata, x, y, v, pad = dup_G_meta(G1) (eltype(G1.z) <: AbstractFloat) ? (z = G1.z .^ _scale) : - (z = Matrix{Float32}(undef, size(G1.z)); @inbounds for k = 1:numel(G1.z) z[k] = G1.z[k] ^ _scale end) + (z = Matrix{Float32}(undef, size(G1.z)); @inbounds for k = 1:numel(G1.z) z[k] = Float32(G1.z[k]) ^ _scale end) G2 = GMTgrid(G1.proj4, G1.wkt, epsg, geog, range, inc, registration, nodata, "", "", "", "", G1.names, x, y, v, z, G1.x_unit, G1.y_unit, G1.v_unit, G1.z_unit, G1.layout, 1f0, 0f0, pad, G1.hasnans) G2.range[5:6] .^= scale @@ -119,7 +119,7 @@ function Base.:/(G1::GMTgrid, G2::GMTgrid) if (size(G1.z) != size(G2.z)) error("Grids have different sizes, so they cannot be divided.") end epsg, geog, range, inc, registration, nodata, x, y, v, pad = dup_G_meta(G1) (eltype(G1.z) <: AbstractFloat && eltype(G2.z) <: AbstractFloat) ? (z = G1.z ./ G2.z) : - (z = Matrix{Float32}(undef, size(G1.z)); @inbounds for k = 1:numel(G1.z) z[k] = G1.z[k] / G2.z[k] end) + (z = Matrix{Float32}(undef, size(G1.z)); @inbounds for k = 1:numel(G1.z) z[k] = Float32(G1.z[k]) / Float32(G2.z[k]) end) hasnans = any(isnan, z) ? 2 : 1 G3 = GMTgrid(G1.proj4, G1.wkt, epsg, geog, range, inc, registration, nodata, "", "", "", "", G1.names, x, y, v, z, G1.x_unit, G1.y_unit, G1.v_unit, G1.z_unit, G1.layout, 1f0, 0f0, pad, hasnans) @@ -134,7 +134,7 @@ function Base.:/(G1::GMTgrid, scale::Real) _scale = convert(eltype(G1.z), scale) epsg, geog, range, inc, registration, nodata, x, y, v, pad = dup_G_meta(G1) (eltype(G1.z) <: AbstractFloat) ? (z = G1.z ./ _scale) : - (z = Matrix{Float32}(undef, size(G1.z)); @inbounds for k = 1:numel(G1.z) z[k] = G1.z[k] / _scale end) + (z = Matrix{Float32}(undef, size(G1.z)); @inbounds for k = 1:numel(G1.z) z[k] = Float32(G1.z[k]) / _scale end) G2 = GMTgrid(G1.proj4, G1.wkt, epsg, geog, range, inc, registration, nodata, "", "", "", "", G1.names, x, y, v, z, G1.x_unit, G1.y_unit, G1.v_unit, G1.z_unit, G1.layout, 1f0, 0f0, pad, G1.hasnans) G2.range[5:6] ./= scale @@ -218,9 +218,22 @@ Base.:>=(I::GMTimage{T}, val::Real) where T <: Unsigned = mat2img(collect(I.imag """ Subtract two boolean/uint8 mask images. It applies the logical `I1 && !I2` operation. Inherit header parameters from I1 image """ +#= Base.:-(I1::GMTimage{<:Bool}, I2::GMTimage{<:Bool}) = helper_bool_img(I1, collect(I1.image .& .!I2.image)) Base.:-(I1::GMTimage{<:UInt8}, I2::GMTimage{<:UInt8}) = helper_bool_img(I1, collect(reinterpret(Bool, I1.image) .& .!reinterpret(Bool, I2.image))) +=# + +Base.:-(I1::GMTimage{<:Bool}, I2::GMTimage{<:Bool}) = helper_bool_img(I1, sub8_layout(I1, I2)) +Base.:-(I1::GMTimage{<:UInt8}, I2::GMTimage{<:UInt8}) = helper_bool_img(I1, sub8_layout(I1, I2)) +function sub8_layout(I1, I2) + I1.layout[2] == I2.layout[2] && + return eltype(I1) <: Bool ? collect(I1.image .& .!I2.image) : collect(reinterpret(Bool, I1.image) .& .!reinterpret(Bool, I2.image)) + + z1 = (eltype(I1) <: Bool) ? I1.image : reinterpret(Bool, I1.image) + z2 = (eltype(I2) <: Bool) ? I2.image' : reinterpret(Bool, I2.image)' + return collect(z1 .& .!z2) +end # --------------------------------------------------------------------------------------------------- """ @@ -289,6 +302,15 @@ function Base.:!(I::GMTimage{<:UInt8}) return Io end +""" + I = togglemask(I::Union{GMTimage{<:Bool, 2}, GMTimage{<:UInt8, 2}}) -> GMTimage + +Convert between UInt8 and Boolean representations of the mask images. A new object is returned but +the underlying mask matrix is kept unchanged, that is, not copied, but reinterpred in the other type. +""" +togglemask(I::GMTimage{<:Bool, 2}) = return mat2img(reinterpret(UInt8, I.image), I) +togglemask(I::GMTimage{<:UInt8, 2}) = return mat2img(reinterpret(Bool, I.image), I) + # --------------------------------------------------------------------------------------------------- Base.:+(add::T, D1::GMTdataset) where T<:AbstractArray = Base.:+(D1::GMTdataset, add) Base.:+(add::Real, D1::GMTdataset) = Base.:+(D1::GMTdataset, [add;;]) From 0a181fcef7e4ee27e8b956c5ebeb3970ef0a932c Mon Sep 17 00:00:00 2001 From: Joaquim Date: Thu, 25 Apr 2024 23:48:48 +0100 Subject: [PATCH 2/4] Check if GMT types are empty instead of === nothing --- src/utils_types.jl | 41 +++++++++++++++++++++++++++-------------- 1 file changed, 27 insertions(+), 14 deletions(-) diff --git a/src/utils_types.jl b/src/utils_types.jl index 08a8b0e51..dbfb7fa61 100644 --- a/src/utils_types.jl +++ b/src/utils_types.jl @@ -921,7 +921,7 @@ end # --------------------------------------------------------------------------------------------------- """ - I = mat2img(mat::Array{<:Unsigned}; x=[], y=[], hdr=[], proj4="", wkt="", cmap=nothing, kw...) + I = mat2img(mat::Array{<:Unsigned}; x=[], y=[], hdr=[], proj4="", wkt="", cmap=GMTcpt(), kw...) Take a 2D 'mat' array and a `hdr` 1x9 [xmin xmax ymin ymax zmin zmax reg xinc yinc] header descriptor and return a GMTimage type. @@ -943,23 +943,23 @@ If `stretch` is a scalar, scale the values > `stretch` to [0 255] The `kw...` kwargs search for [:layout :mem_layout], [:names] and [:metadata] """ function mat2img(mat::Union{AbstractArray{<:Unsigned}, AbstractArray{<:Bool}, BitMatrix}; x=Float64[], y=Float64[], v=Float64[], hdr=Float64[], - proj4::String="", wkt::String="", cmap=nothing, is_transposed::Bool=false, kw...) + proj4::String="", wkt::String="", cmap=GMTcpt(), is_transposed::Bool=false, kw...) # Take a 2D array of uint8 and turn it into a GMTimage. # Note: if HDR is empty we guess the registration from the sizes of MAT & X,Y - (cmap === nothing && eltype(mat) == Bool) && (cmap = makecpt(T=(0,1), cmap=:gray)) + (cmap === nothing && (eltype(mat) == Bool || eltype(mat) == BitMatrix)) && (cmap = makecpt(T=(0,1), cmap=:gray)) helper_mat2img(mat; x=x, y=y, v=v, hdr=hdr, proj4=proj4, wkt=wkt, cmap=cmap, is_transposed=is_transposed, kw...) end # Special version to desambiguate between UInt8 and UInt16 function mat2img16(mat::AbstractArray{<:Unsigned}; x=Float64[], y=Float64[], v=Float64[], hdr=Float64[], - proj4::String="", wkt::String="", cmap=nothing, is_transposed::Bool=false, kw...) + proj4::String="", wkt::String="", cmap=GMTcpt(), is_transposed::Bool=false, kw...) helper_mat2img(mat; x=x, y=y, v=v, hdr=hdr, proj4=proj4, wkt=wkt, cmap=cmap, is_transposed=is_transposed, kw...) end function helper_mat2img(mat; x=Float64[], y=Float64[], v=Float64[], hdr=Float64[], - proj4::String="", wkt::String="", cmap=nothing, is_transposed::Bool=false, kw...) + proj4::String="", wkt::String="", cmap=GMTcpt(), is_transposed::Bool=false, kw...) color_interp = ""; n_colors = 0; - isa(mat, BitMatrix) && (mat = UInt8.(mat)*UInt8(255)) # Tried but can't find a way to make GMTimage accept also BitMatrix - if (cmap !== nothing) + #isa(mat, BitMatrix) && (mat = collect(mat)) + if (!isempty(cmap)) colormap, labels, n_colors = cpt2cmap(cmap) else (size(mat,3) == 1) && (color_interp = "Gray") @@ -999,9 +999,10 @@ function cpt2cmap(cpt::GMTcpt, start::Float32=NaN32, force_alpha::Bool=true) cmap = zeros(Int32, 256 * nc) (s == 1) && (cmap[1] = 255; cmap[257] = 255; cmap[513] = 255) # nodata pixel color = white n_colors = 256; # Because for GDAL we always send 256 even if they are not all filled + gray_max = (cpt.minmax[2] == 1) ? 1 : 255 # First, carefull, change. Probably should be 'gray_max = cpt.minmax[2]' for n = 1:3 # Write 'cmap' col-wise for m = 1:size(cpt.colormap, 1) - @inbounds cmap[m+s + (n-1)*n_colors] = round(Int32, cpt.colormap[m,n] * 255); + @inbounds cmap[m+s + (n-1)*n_colors] = round(Int32, cpt.colormap[m,n] * gray_max); end end if (have_alpha) # Have alpha color(s) @@ -1125,7 +1126,7 @@ end # --------------------------------------------------------------------------------------------------- function mat2img(mat::Union{GMTgrid,Matrix{<:AbstractFloat}}; x=Float64[], y=Float64[], hdr=Float64[], - proj4::String="", wkt::String="", GI::Union{GItype,Nothing}=nothing, clim=[0,255], cmap=nothing, kw...) + proj4::String="", wkt::String="", GI::Union{GItype,Nothing}=nothing, clim=[0,255], cmap=GMTcpt(), kw...) # This is the same as Matlab's imagesc() ... plus some extras. mi, ma = (isa(mat,GMTgrid)) ? mat.range[5:6] : extrema(mat) (isa(mat,GMTgrid) && mat.hasnans == 2) && (mi = NaN) # Don't know yet so force checking @@ -1147,8 +1148,8 @@ function mat2img(mat::Union{GMTgrid,Matrix{<:AbstractFloat}}; x=Float64[], y=Flo is_transp = isa(mat, GItype) ? ((mat.layout[2] == 'R') ? true : false) : false if (!isa(mat, GMTgrid) && GI !== nothing) I = mat2img(img, GI) - if (cmap !== nothing) I.colormap, I.labels, I.n_colors = cpt2cmap(cmap) - else I.colormap, I.labels, I.n_colors = zeros(Int32,3), String[], 0 # Do not inherit this from GI + if (!isempty(cmap)) I.colormap, I.labels, I.n_colors = cpt2cmap(cmap) + else I.colormap, I.labels, I.n_colors = zeros(Int32,3), String[], 0 # Do not inherit this from GI end elseif (isa(mat, GMTgrid)) I = mat2img(img; x=mat.x, y=mat.y, hdr=hdr, proj4=mat.proj4, wkt=mat.wkt, cmap=cmap, is_transposed=is_transp, kw...) @@ -1184,7 +1185,7 @@ issuing an error. In this case `clim` can be a two elements vector to specify th The default is to let `histogram` guess these values. """ function imagesc(mat::Union{GMTgrid,Matrix{<:AbstractFloat}}; x=Float64[], y=Float64[], hdr=Float64[], - proj4::String="", wkt::String="", GI::Union{GItype,Nothing}=nothing, clim=[0,255], cmap=nothing, kw...) + proj4::String="", wkt::String="", GI::Union{GItype,Nothing}=nothing, clim=[0,255], cmap=GMTcpt(), kw...) # Call 'rescale' and return if the kw 'stretch' is used ((stretch = find_in_kwargs(kw, [:stretch])[1]) !== nothing) && return rescale(mat, stretch=stretch, type=UInt8) @@ -1390,8 +1391,8 @@ function slicecube(GI::GItype; slice::Int=0, α=0.0, angle=0.0, axis="x", cmap=G end end I = mat2img(sc, GI) - if (cmap !== nothing) I.colormap, I.labels, I.n_colors = cpt2cmap(cmap) - else I.colormap, I.labels, I.n_colors = zeros(Int32,3), String[], 0 # Do not inherit this from GI + if (!isempty(cmap)) I.colormap, I.labels, I.n_colors = cpt2cmap(cmap) + else I.colormap, I.labels, I.n_colors = zeros(Int32,3), String[], 0 # Do not inherit this from GI end return mat2grid(z, GI), I else @@ -2103,6 +2104,18 @@ function grid2pix(hdr::Vector{Float64}; pix=true) return hdr end +# --------------------------------------------------------------------------------------------------- +""" + is_stored_transposed(GI::GItype) -> Bool + +Return true if the data in the GMTgrid or GMTimage `GI` is stored transposed or false otherwise. +See more details in the comments of the `gmt2gd` function. +""" +function is_stored_transposed(GI::GItype) + width, height = (length(GI.x), length(GI.y)) .- GI.registration + return width == size(GI, 1) && height == size(GI, 2) +end + #= --------------------------------------------------------------------------------------------------- function mksymbol(f::Function, cmd0::String="", arg1=nothing; kwargs...) # Make a fig and convert it to EPS so it can be used as a custom symbol is plot(3) From 34ac25a5fd30ef0ca65c489824ffdabc600b95f8 Mon Sep 17 00:00:00 2001 From: Joaquim Date: Thu, 25 Apr 2024 23:49:36 +0100 Subject: [PATCH 3/4] Make possible to declare GItypes as empties. --- src/gmt_types.jl | 118 +++++++++++++++++++++++------------------------ 1 file changed, 59 insertions(+), 59 deletions(-) diff --git a/src/gmt_types.jl b/src/gmt_types.jl index c602d9f52..51c208dd7 100644 --- a/src/gmt_types.jl +++ b/src/gmt_types.jl @@ -33,28 +33,28 @@ The fields of this struct are: - `hasnans::Int=0`: 0 -> "don't know"; 1 -> confirmed, "have no NaNs"; 2 -> confirmed, "have NaNs" """ Base.@kwdef mutable struct GMTgrid{T<:Number,N} <: AbstractArray{T,N} - proj4::String - wkt::String - epsg::Int - geog::Int - range::Union{Vector{Float64}, Vector{Any}} - inc::Union{Vector{Float64}, Vector{Any}} - registration::Int - nodata::Union{Float64, Float32} - title::String - remark::String - command::String - cpt::String - names::Vector{String} - x::Array{Float64,1} - y::Array{Float64,1} - v::Union{Vector{<:Real}, Vector{String}, Vector{<:TimeType}} - z::Array{T,N} - x_unit::String - y_unit::String - v_unit::String - z_unit::String - layout::String + proj4::String="" + wkt::String="" + epsg::Int=0 + geog::Int=0 + range::Union{Vector{Float64}, Vector{Any}}=Float64[] + inc::Union{Vector{Float64}, Vector{Any}}=Float64[] + registration::Int=0 + nodata::Union{Float64, Float32}=0.0 + title::String="" + remark::String="" + command::String="" + cpt::String="" + names::Vector{String}=String[] + x::Array{Float64,1}=Float64[] + y::Array{Float64,1}=Float64[] + v::Union{Vector{<:Real}, Vector{String}, Vector{<:TimeType}}=String[] + z::Array{T,N}=Array{Float64,2}(undef,0,0) + x_unit::String="" + y_unit::String="" + v_unit::String="" + z_unit::String="" + layout::String="" scale::Union{Float64, Float32}=1f0 offset::Union{Float64, Float32}=0f0 pad::Int=0 @@ -106,28 +106,28 @@ The fields of this struct are: - `layout::String`: A four character string describing the image memory layout - `pad::Int`: When != 0 means that the array is placed in a padded array of PAD rows/cols """ -mutable struct GMTimage{T<:Union{Unsigned, Bool}, N} <: AbstractArray{T,N} - proj4::String - wkt::String - epsg::Int - geog::Int - range::Vector{Float64} - inc::Vector{Float64} - registration::Int - nodata::Float32 - color_interp::String - metadata::Vector{String} - names::Vector{String} - x::Vector{Float64} - y::Vector{Float64} - v::Vector{Float64} - image::Array{T,N} - colormap::Vector{Int32} - labels::Vector{String} # Labels of a Categorical CPT - n_colors::Int - alpha::Matrix{UInt8} - layout::String - pad::Int +Base.@kwdef mutable struct GMTimage{T<:Union{Unsigned, Bool, BitMatrix}, N} <: AbstractArray{T,N} + proj4::String="" + wkt::String="" + epsg::Int=0 + geog::Int=0 + range::Vector{Float64}=Float64[] + inc::Vector{Float64}=Float64[] + registration::Int=0 + nodata::Float32=0f0 + color_interp::String="" + metadata::Vector{String}=String[] + names::Vector{String}=String[] + x::Vector{Float64}=Float64[] + y::Vector{Float64}=Float64[] + v::Vector{Float64}=Float64[] + image::AbstractArray{T,N}=Array{UInt8,2}(undef,0,0) + colormap::Vector{Int32}=Int32[] + labels::Vector{String}=String[] # Labels of a Categorical CPT + n_colors::Int=0 + alpha::Matrix{UInt8}=Array{UInt8,2}(undef,0,0) + layout::String="" + pad::Int=0 end Base.size(I::GMTimage) = size(I.image) Base.getindex(I::GMTimage{T,N}, inds::Vararg{Int,N}) where {T,N} = I.image[inds...] @@ -216,19 +216,19 @@ The fields of this struct are: - `epsg::Int`: EPSG projection code (Optional) - `geom::Integer`: Geometry type. One of the GDAL's enum (wkbPoint, wkbPolygon, etc...) """ -mutable struct GMTdataset{T<:Real, N} <: AbstractArray{T,N} - data::Array{T,N} - ds_bbox::Vector{Float64} - bbox::Vector{Float64} - attrib::DictSvS - colnames::Vector{String} - text::Vector{String} - header::String - comment::Vector{String} - proj4::String - wkt::String - epsg::Int - geom::Union{UInt32, Int} # 0->Unknown, 1->Point, 2->Line, 3->Polygon, 4->MultiPoint, 5->MultiLine, 6->MultiPolyg +Base.@kwdef mutable struct GMTdataset{T<:Real, N} <: AbstractArray{T,N} + data::Array{T,N}=Array{Float64,2}(undef,0,0) + ds_bbox::Vector{Float64}=Float64[] + bbox::Vector{Float64}=Float64[] + attrib::DictSvS=DictSvS() + colnames::Vector{String}=String[] + text::Vector{String}=String[] + header::String="" + comment::Vector{String}=String[] + proj4::String="" + wkt::String="" + epsg::Int=0 + geom::Union{UInt32, Int}=0 # 0->Unknown, 1->Point, 2->Line, 3->Polygon, 4->MultiPoint, 5->MultiLine, 6->MultiPolyg end Base.size(D::GMTdataset) = size(D.data) Base.getindex(D::GMTdataset{T,N}, inds::Vararg{Int,N}) where {T,N} = D.data[inds...] @@ -286,8 +286,8 @@ GMTdataset(data::Array{Float32,2}, text::String) = GMTdataset(data, Float64[], Float64[], DictSvS(), String[], [text], "", String[], "", "", 0, 0) GMTdataset(data::Array{Float32,2}) = GMTdataset(data, Float64[], Float64[], DictSvS(), String[], String[], "", String[], "", "", 0, 0) -GMTdataset() = - GMTdataset(Array{Float64,2}(undef,0,0), Float64[], Float64[], DictSvS(), String[], String[], "", String[], "", "", 0, 0) +#GMTdataset() = + #GMTdataset(Array{Float64,2}(undef,0,0), Float64[], Float64[], DictSvS(), String[], String[], "", String[], "", "", 0, 0) struct WrapperPluto fname::String end From 7c07ef4898a8b1f318dabb84c1d5fd1cedf262f9 Mon Sep 17 00:00:00 2001 From: Joaquim Date: Thu, 25 Apr 2024 23:50:24 +0100 Subject: [PATCH 4/4] Comments tweaks. --- src/gdal_utils.jl | 4 ++-- src/gmt_main.jl | 12 +++++++----- 2 files changed, 9 insertions(+), 7 deletions(-) diff --git a/src/gdal_utils.jl b/src/gdal_utils.jl index c7e622b30..0132eb2a4 100644 --- a/src/gdal_utils.jl +++ b/src/gdal_utils.jl @@ -366,8 +366,8 @@ function gmt2gd(GI) # Julia show methods. A further problem is that if we use 'GI = gmtread("file", layout="TRB")', GI's # array still has a visible mem layout equal to a column major array. AS an example of this, see: # gmtwrite("z.grd", mat2grid(reshape(1:20, 4,5))); Gr = gmtread("z.grd", layout="TRB"); Gc = gmtread("z.grd"); - # Now both Gr.z and Gc.r are 4x5 matrices, but Gr can't be passed as is to 'gmt2gd(...)' because its - # memory layout is wrong when seen from the GDAL wrapper functions. To disambiguate this, we are not + # Now both Gr.z and Gc.z are 4x5 matrices, but Gr can't be passed as is to 'gmt2gd(...)' because its + # memory layout is wrong when seen from the GDAL wrapper functions. To disambiguate this, we are not # relying in 'size(GI)' but on the length of the x,y coordinates vectors, that are assumed to always be correct. width, height = (GI.layout != "" && GI.layout[2] == 'C') ? (size(GI,2), size(GI,1)) : (length(GI.x), length(GI.y)) .- GI.registration diff --git a/src/gmt_main.jl b/src/gmt_main.jl index b8dc659f0..08dfbca92 100644 --- a/src/gmt_main.jl +++ b/src/gmt_main.jl @@ -483,10 +483,10 @@ function get_image(API::Ptr{Nothing}, object)::GMTimage o = (nz == 1) ? (ny, nx) : (ny, nx, nz) isBRP = startswith(layout, "BRP") (nz == 1 && isBRP) && (layout = "BRPa") # For 1 layer "BRBa" and "BRPa" is actualy the same. - (!isBRP) && @warn("Only 'I' for Images.jl and 'BRP' MEM layouts are allowed.") + (!isBRP && nz > 1) && @warn("Only 'I' for Images.jl and 'BRP' MEM layouts are allowed.") # OK, the above is not always true. Image cubes are not pixel interleaved. But how to detect them? - (I.type > 1 || (nz == 2 || nz > 4)) && (layout = layout[1:2] * "B") # Stil leaves out cases of uint8 cubes. + (I.type > 1 || (nz == 2 || nz > 4)) && (layout = layout[1:2] * "B" * layout[4]) # Stil leaves out cases of uint8 cubes. end t = reshape(unsafe_wrap(Array, data, ny * nx * nz), o) # Apparently the reshape() creates a copy end @@ -923,7 +923,8 @@ function image_init(API::Ptr{Nothing}, Img::GMTimage)::Ptr{GMT_IMAGE} toRP_pad(Img, img_padded, n_rows, n_cols, pad) already_converted = true else - Ib.data = pointer(Img.image) + isa(Img.image, BitMatrix) && (copy_data = collect(Img.image)) # BitMatrx cannot be sent to C and we may need a copy bellow + Ib.data = (!isa(Img.image, BitMatrix)) ? pointer(Img.image) : pointer(copy_data) mem_owned_by_gmt = (pad == 0) ? false : true end @@ -946,9 +947,9 @@ function image_init(API::Ptr{Nothing}, Img::GMTimage)::Ptr{GMT_IMAGE} unsafe_store!(I, Ib) if (!already_converted && !startswith(Img.layout, "BRP")) - img = (mem_owned_by_gmt) ? img_padded : copy(Img.image) + img = (mem_owned_by_gmt) ? img_padded : copy(Img.image) # This copy is a waste when not Change_layout. Needs revisit. (size(img,3) > 2) && GMT_Change_Layout(API, GMT_IS_IMAGE, "BRP", 0, I, img); # Convert to BRP. Not 100% on the > 2 though. - Ib.data = pointer(img) + Ib.data = !isa(Img.image, BitMatrix) ? pointer(img) : pointer(copy_data) unsafe_store!(I, Ib) end @@ -1414,6 +1415,7 @@ Shows information about the `D` GMTdataset (or vector of them). Runs ``show(stdout, "text/plain", any)`` which prints all elements of `any`. Good for printing the entire vector or matrix. """ function info(GI::GItype, showdata::Bool=true; data=true, full=false, crs::Bool=false) + isempty(GI) && return println("Empty object") crs && return print_crs(GI) (data != 1) && (showdata = false) isa(GI, GMTimage) ? println("A GMTimage object with $(size(GI,3)) bands of type $(eltype(GI))") :