Fix the computation of min/max of images with a finite nodata.

src/GMT.jl

@@ -123,8 +123,8 @@ export
 	radar!, stairs, stairs!, stem, stem!,vlines, vlines!, vband, vband!, hspan, hspan!, vspan, vspan!,
 	basemap, basemap!, blockmean, blockmedian, blockmode, clip, clip!,
 	coast, coast!, colorbar, colorbar!, colorscale, colorscale!, contour, contour!, contourf, contourf!, events,
-	filter1d, fitcircle, gmt2kml, gmtbinstats, binstats,
-	gmtconnect, gmtconvert, gmtinfo, gmtmath, gmtregress, gmtread, gmtselect, gmtset, gmtsimplify, gmtspatial,
+	filter1d, fitcircle, gmt2kml, gmtbinstats, binstats, gmtconnect, gmtconvert,
+	gmtinfo, gmtlogo, gmtlogo!, gmtmath, gmtregress, gmtread, gmtselect, gmtset, gmtsimplify, gmtspatial,
 	gmtvector, gmtwrite, gmtwhich, grd2cpt, grd2kml, grd2xyz, grdblend, grdclip, grdcontour, grdcontour!, grdconvert,
 	grdcut, grdedit, grdfft, grdfill, grdfilter, grdgradient, grdhisteq, grdimage, grdimage!, grdinfo, grdinterpolate,
 	grdlandmask, grdmath, grdmask, grdpaste, grdproject, grdsample, grdtrack, grdtrend, grdvector, grdvector!,

src/common_options.jl

@@ -5160,6 +5160,9 @@ macro var"?"(name)
 			         :flatfv,:fv2fv,:icosahedron,:loft,:sphere, :octahedron,:tetrahedron,:torus,:replicant,:revolve
 					]
 					dir = "utilities/"
+			elseif sym in [:gd2gmt, :gmt2gd, 
+						]
+						dir = "gdalfuns/"
 			elseif sym in [:binarize, :bwareaopen, :bwconncomp, :bwdist, :bwhitmiss, :bwperim, :bwskell, :cc2bw, :imdilate,
 						 :imerode, :imopen, :imclose, :imtophat, :imbothat, :immorphgrad, :strel, :imfilter, :imrankfilter,
 						 :imsobel, :graydist, :padarray, :imfill, :imreconstruct, :imsegment, :isodata, :imregionalmax,

src/gdal_utils.jl

@@ -95,10 +95,14 @@ function gd2gmt(_dataset; band::Int=0, bands=Vector{Int}(), sds::Int=0, pad::Int
 	x_min, x_max, y_min, y_max, x_inc, y_inc = getregion(dataset, pad=pad, xSize=xSize, ySize=ySize, gridreg=is_grid, sds=sds)
 
 	if !(eltype(mat) <: Complex)
-		z_min::Float64, z_max::Float64 = (is_grid) ? extrema_nan(mat) : extrema(mat)
+		if (!got_fill_val || !isfinite(fill_val))
+			z_min::Float64, z_max::Float64 = (is_grid) ? extrema_nan(mat) : extrema(mat)
+		else
+			z_min, z_max = extrema_nodataval(mat, fill_val)
+		end
 		isnan(z_min) && (z_min, z_max = extrema_nan(mat))	# Convoluted nc UInt16 layers with scalefactors can result in NaNs
 	else
-		z_min, z_max, z_im_min, z_im_max = extrema(mat)
+		z_min, z_max, z_im_min, z_im_max = extrema(mat)		# SHOULD DO THE SAME AS ABOVE
 	end
 	hdr = [x_min, x_max, y_min, y_max, z_min, z_max, Float64(!is_grid), x_inc, y_inc]
 	prj = getproj(dataset)

src/gmtlogo.jl

@@ -5,8 +5,6 @@ Plots the GMT logo on a map. By default, the GMT logo is 5 cm wide and 2.5 cm hi
 positioned relative to the current plot origin. Use various options to change this and to place
 a transparent or opaque rectangular map panel behind the GMT logo.
 
-See full GMT docs at [`gmtlogo`]($(GMTdoc)gmtlogo.html)
-
 Parameters
 ----------
 
@@ -106,4 +104,7 @@ function jlogo(L::Float64=5.0)
 	s2 = s1 * (1 - 0.06)		# Inner circle diameter. The one that will be filled.
 	t = [L/2 L/2 0 s1; L+L/2 L/2 1 s1; L L/2+H 2 s1; L/2 L/2 3 s2; L+L/2 L/2 4 s2; L L/2+H 5 s2]
 	return " -Sc -C171/43/33,130/83/171,81/143/24,191/101/95,158/122/190,128/171/93 ", t, s2
-end
+end
+
+const gmtlogo  = logo			# Alias
+const gmtlogo! = logo!

src/grd_operations.jl

@@ -32,9 +32,24 @@ function Base.:+(G1::GMTgrid, shift::Real)
 	G2.range[5:6] .+= shift
 	return G2
 end
+
 # Now for images
-Base.:+(G1::GMTimage{T}, shift::Real) where T <: Unsigned = img2grid(G1) + shift
-Base.:+(shift::Real, G1::GMTimage{T}) where T <: Unsigned = img2grid(G1) + shift
+Base.:+(G1::GMTimage{T}, shift::AbstractFloat) where T <: Unsigned = img2grid(G1) + shift
+function Base.:+(G1::GMTimage{T}, shift::Integer) where T <: Unsigned
+	_shift = Int64(shift)
+	shift_typed = convert(eltype(G1.image), _shift)
+	if ((G1.range[5] + _shift) >= typemin(eltype(G1)) && (G1.range[6] + _shift) <= typemax(eltype(G1)))
+		# All values will fit in the same type of G1
+		mat = (!isfinite(G1.nodata)) ? G1.image .+ shift_typed : similar(G1.image)
+		if (isfinite(G1.nodata))
+			nodata = convert(eltype(G1.image), G1.nodata)
+			@inbounds for k in eachindex(mat) mat[k] = ifelse((G1.image[k] != G1.nodata), G1.image[k] + shift_typed, nodata)  end
+		end
+		return mat2img(mat, G1, nodata=G1.nodata)
+	else				# Values will overflow, so make it a Float32 grid
+		return img2grid(G1) + shift
+	end
+end
 
 # ---------------------------------------------------------------------------------------------------
 function Base.:-(G1::GMTgrid, G2::GMTgrid)
@@ -63,9 +78,24 @@ function Base.:-(G1::GMTgrid, shift::Real)
 	G2.range[5:6] .-= shift
 	return G2
 end
+
 # Now for images
-Base.:-(G1::GMTimage{T}, shift::Real) where T <: Unsigned = img2grid(G1) - shift
-Base.:-(shift::Real, G1::GMTimage{T}) where T <: Unsigned = img2grid(G1) - shift
+Base.:-(G1::GMTimage{T}, shift::AbstractFloat) where T <: Unsigned = img2grid(G1) - shift
+function Base.:-(G1::GMTimage{T}, shift::Integer) where T <: Unsigned
+	_shift = Int64(shift)
+	shift_typed = convert(eltype(G1.image), _shift)
+	if ((G1.range[5] - _shift) >= typemin(eltype(G1)) && (G1.range[6] - _shift) <= typemax(eltype(G1)))
+		# All values will fit in the same type of G1
+		mat = (!isfinite(G1.nodata)) ? G1.image .- shift_typed : similar(G1.image)
+		if (isfinite(G1.nodata))
+			nodata = convert(eltype(G1.image), G1.nodata)
+			@inbounds for k in eachindex(mat) mat[k] = ifelse((G1.image[k] != G1.nodata), G1.image[k] - shift_typed, nodata)  end
+		end
+		return mat2img(mat, G1, nodata=G1.nodata)
+	else				# Values will overflow, so make it a Float32 grid
+		return img2grid(G1) - shift
+	end
+end
 
 # ---------------------------------------------------------------------------------------------------
 function Base.:*(G1::GMTgrid, G2::GMTgrid)

src/psxy.jl

@@ -16,7 +16,7 @@ function common_plot_xyz(cmd0::String, arg1, caller::String, first::Bool, is3D::
 end
 function common_plot_xyz(cmd0::String, arg1, caller::String, first::Bool, is3D::Bool, d::Dict{Symbol, Any})
 	(cmd0 != "" && arg1 === nothing && is_in_dict(d, [:groupvar :hue]) !== nothing) && (arg1 = gmtread(cmd0); cmd0 = "")
-	invokelatest(_common_plot_xyz, cmd0, arg1, caller, !first, true, is3D, d)
+	_common_plot_xyz(cmd0, arg1, caller, !first, true, is3D, d)
 	#invokelatest(_common_plot_xyz, cmd0, arg1, caller, !first, true, is3D, d)
 end
 function _common_plot_xyz(cmd0::String, arg1, caller::String, O::Bool, K::Bool, is3D::Bool, d::Dict{Symbol, Any})

src/utils_types.jl

@@ -1379,7 +1379,7 @@ function mat2img(mat::Union{AbstractMatrix{UInt16},AbstractArray{UInt16,3}}; x=F
 	# stretch = [v1 v2 v3 v4 v5 v6] scales firts band >= v1 && <= v2 to [0 255], second >= v3 && <= v4, same for third
 	# If the img8 argument is used, it should contain a pre-allocated UInt8 array with the exact same size as MAT.
 	# The 'scale_only' kw option makes it return just the scaled array an no GMTimage creation (useful when img8 is a view).
-	# Use the keyword NOCONV to return GMTimage UInt16 type. I.e., no conversion to UInt8
+	# Use the keyword `noconv=true` to return GMTimage UInt16 type. I.e., no conversion to UInt8
 
 	d = KW(kw)
 	x::Vector{Float64} = vec(x);	y::Vector{Float64} = vec(y);	v::Vector{Float64} = vec(v);
@@ -1536,8 +1536,8 @@ end
 
 # ---------------------------------------------------------------------------------------------------
 # This method creates a new GMTimage but retains all the header data from the IMG object
-function mat2img(mat, I::GMTimage; names::Vector{String}=String[], metadata::Vector{String}=String[])
-	range = copy(I.range);	(size(mat,3) == 1) && (range[5:6] .= extrema(mat))
+function mat2img(mat, I::GMTimage; names::Vector{String}=String[], metadata::Vector{String}=String[], nodata=NaN)
+	range = copy(I.range);	(size(mat,3) == 1) && (range[5:6] .= !isfinite(nodata) ? extrema(mat) : extrema_nodataval(mat, nodata))
 	GMTimage(I.proj4, I.wkt, I.epsg, I.geog, range, copy(I.inc), I.registration, NaN32, I.color_interp, metadata, names, copy(I.x), copy(I.y), zeros(size(mat,3)), mat, copy(I.colormap), String[], I.n_colors, Array{UInt8,2}(undef,1,1), I.layout, 0)
 end
 function mat2img(mat, G::GMTgrid; names::Vector{String}=String[], metadata::Vector{String}=String[])

test/test_misc.jl

@@ -77,6 +77,16 @@
 	I1 ⊻ I2;
 	!I1;
 
+	I = mat2img(reshape(collect(UInt16(11):UInt16(26)), 4,4), noconv=true);
+	I2 = I - 11;
+	@test eltype(I2) == UInt16
+	I2 = I - 12
+	@test eltype(I2) == Float32
+	I2 = I + 11;
+	@test eltype(I2) == UInt16
+	I2 = I + typemax(UInt16)
+	@test eltype(I2) == Float32
+
 	@test GMT.bin2dec("10111") == 23
 	@test GMT.dec2bin(23) == "10111"
 	@test GMT.sub2ind((3,3), [1 2 3 1], [2 2 2 3]) == [4  5  6  7]