Fix some recent regressions. Replace some module calls by their gmt(...) equivalent.

src/GMT.jl

@@ -132,7 +132,7 @@ export
 	stereonet, stereonet!,
 
 	add_opt, isgeog, numel, scan_opt, extrema_nan, parse_RIr, close_PS_file, getsize, meshgrid, parse_B, parse_BJR, parse_I,
-	parse_J, parse_R, ressurectGDAL, CPTaliases, isJLL, POSTMAN,
+	parse_J, parse_R, ressurectGDAL, CPTaliases, isJLL, POSTMAN, PSname,
 
 	mgd77magref, magref,
 

src/contourf.jl

@@ -239,6 +239,6 @@ function auto_contour_interval(range)
 end
 
 # ---------------------------------------------------------------------------------------------------
-contourf!(cmd0::String="", arg1=nothing, arg2=nothing; kw...) = contourf(cmd0, arg1, arg2, false, KW(kw))
+contourf!(cmd0::String="", arg1=nothing, arg2=nothing; kw...) = contourf(cmd0, arg1, arg2; first=false, kw...)
 contourf(arg1, arg2=nothing; kw...) = contourf("", arg1, arg2, true, KW(kw))
 contourf!(arg1, arg2=nothing; kw...) = contourf("", arg1, arg2, false, KW(kw))

src/drawing.jl

@@ -83,7 +83,7 @@ function helper_shapes(x, y, cmd; Vd=0, kw...)
 	((isa(first, Bool) && !first) || (isa(first, Int) && first == 0)) && (CTRLshapes.first[1] = false)
 
 	if (CTRLshapes.first[1])
-		CTRLshapes.fname[1] = GMT.PSname[1]
+		CTRLshapes.fname[1] = PSname[1]
 		if ((val = find_in_dict(d, [:units])[1]) !== nothing)
 			_cmd, opt_B = parse_B(d, "", DEF_FIG_AXES[1])
 			cmd *= " --PROJ_LENGTH_UNIT=p"

src/gadm.jl

@@ -47,7 +47,7 @@ function gadm(country, subregions...; children::Bool=false, names::Bool=false, c
 	reportlevels && return nlayers
 
 	function _filterlayer(layer, layer_name, level, value, all::Bool=false)
-		filtered = Vector{GMT.Gdal.Feature}(undef,0)
+		filtered = Vector{Gdal.Feature}(undef,0)
 		# They Fck changed this between version 3.6 and 4.5. Before, it was always NAME_0, Name_1, ...
 		key = "NAME_$(level)"
 		if (level == 0)

src/imgtiles.jl

@@ -1006,21 +1006,21 @@ function geocoder(address::String; options=String[])::GDtype
 	# Get the geocoder info for a given address. Adapted from https://www.itopen.it/geocoding-with-gdal/
 
 	_ops = isempty(options) ? C_NULL : options	# The default is ["SERVICE", "OSM_NOMINATIM"]
-	hSession = GMT.Gdal.OGRGeocodeCreateSession(_ops)
-	hLayer = GMT.Gdal.OGRGeocode(hSession, address, C_NULL, [""])
-	hFDefn = GMT.Gdal.OGR_L_GetLayerDefn(hLayer)
-	hFeature = GMT.Gdal.OGR_L_GetNextFeature(hLayer)
-	count = GMT.Gdal.OGR_FD_GetFieldCount(hFDefn)
+	hSession = Gdal.OGRGeocodeCreateSession(_ops)
+	hLayer = Gdal.OGRGeocode(hSession, address, C_NULL, [""])
+	hFDefn = Gdal.OGR_L_GetLayerDefn(hLayer)
+	hFeature = Gdal.OGR_L_GetNextFeature(hLayer)
+	count = Gdal.OGR_FD_GetFieldCount(hFDefn)
 	(count == 0) && (@warn("No result found for the address $address"); return GMTdataset())
 	dic = Dict{String,String}()
 	for k = 0:count-1
-		hFieldDefn = GMT.Gdal.OGR_FD_GetFieldDefn(hFDefn,k)
-		((val = GMT.Gdal.OGR_F_GetFieldAsString(hFeature, k)) != "") && (dic[GMT.Gdal.OGR_Fld_GetNameRef(hFieldDefn)] = val)
+		hFieldDefn = Gdal.OGR_FD_GetFieldDefn(hFDefn,k)
+		((val = Gdal.OGR_F_GetFieldAsString(hFeature, k)) != "") && (dic[Gdal.OGR_Fld_GetNameRef(hFieldDefn)] = val)
 	end
 
 	BB = parse.(Float64, split(dic["boundingbox"], ","))
-	GMT.Gdal.OGRGeocodeFreeResult(hLayer)
-	GMT.Gdal.OGRGeocodeDestroySession(hSession)
+	Gdal.OGRGeocodeFreeResult(hLayer)
+	Gdal.OGRGeocodeDestroySession(hSession)
 	D::GDtype = mat2ds([parse(Float64, dic["lon"]) parse(Float64, dic["lat"])], attrib=dic, proj4=prj4WGS84, geom=wkbPoint)
 	D.ds_bbox = [BB[3], BB[4], BB[1], BB[2]]
 	return D

src/imshow.jl

@@ -79,8 +79,8 @@ function imshow(arg1, x::AbstractVector{Float64}=Float64[], y::AbstractVector{Fl
 			isFV3D = (isa(arg1, GMTfv) && size(arg1.verts, 2) > 2)
 		end
 		CTRL.limits[1:4] = ginfo[1:4];		CTRL.limits[7:10] = ginfo[1:4]
-		call_plot3 = ((isa(arg1, GMTdataset) && arg1.geom == Gdal.wkbPolygonZM) ||
-		              (isa(arg1, Vector{<:GMTdataset}) && arg1[1].geom == Gdal.wkbPolygonZM) ||
+		call_plot3 = ((isa(arg1, GMTdataset) && arg1.geom == wkbPolygonZM) ||
+		              (isa(arg1, Vector{<:GMTdataset}) && arg1[1].geom == wkbPolygonZM) ||
 					  isFV3D) ? true : false			# Should evolve into a fun that detects the several plot3d cases.
 		!call_plot3 && (call_plot3 = isplot3(kw))
 		return (call_plot3) ? plot3d(arg1; show=see, kw...) : plot(arg1; show=see, kw...)

src/maregrams.jl

@@ -112,7 +112,7 @@ function maregrams(x::Real, y::Real; days=2, starttime::String="", printurl::Boo
 	# Find the closest station to input x,y
 	@assert -180 <= x <= 360 && -90 <= y <= 90  "Coordinates must be between -180 and 360 and -90 and 90."
 	d = read_maregrams()
-	dists = mapproject(d[:pos], G="$x/$y", o=2)
+	dists = gmt("mapproject -G$(x)/$(y) -o2", d[:pos])
 	ind = argmin(dists)
 	maregrams(code=d[:code][ind], days=days, starttime=starttime, printurl=printurl)
 end

src/pcolor.jl

@@ -42,7 +42,7 @@ This form takes a grid (or the file name of one) as input an paints it's cell wi
 ---
     pcolor(GorD; kwargs...)
 
-If `GorD` is either a GMTgrid or a GMTdataset containing a Pearson correlation matrix obtained with ``GMT.cor()``,
+If `GorD` is either a GMTgrid or a GMTdataset containing a Pearson correlation matrix obtained with ``cor()``,
 the processing recieves a special treatment. In this case, other than the `labels` keyword, user is also
 interested in seing if the automatic choice of x-annotaions angle is correct. If not, one can force it
 by setting the `rotx` (ot `slanted`) keywords.
@@ -73,7 +73,7 @@ by setting the `rotx` (ot `slanted`) keywords.
 Display a Pearson's correlation matrix
 
 ```julia
-pcolor(GMT.cor(rand(4,4)), labels=:y, colorbar=1, show=true)
+pcolor(cor(rand(4,4)), labels=:y, colorbar=1, show=true)
 ```
 """
 function pcolor(X_::VMr, Y_::VMr, C::Union{Nothing, AbstractMatrix{<:Real}}=nothing; first::Bool=true, kwargs...)
@@ -155,7 +155,7 @@ function pcolor(cmd0::String="", arg1=nothing; first=true, kwargs...)
 	(cmd0 != "") && (arg1 = gmtread(cmd0))
 	(isa(arg1, Matrix)) && (arg1 = mat2grid(Float32.(arg1)))
 	isdataframe(arg1) && (arg1 = df2ds(arg1))
-	if ((arg1 == arg1' && arg1[1] == 1 && arg1[end] == 1))		# A corr matrix (computed with GMT.cor())
+	if ((arg1 == arg1' && arg1[1] == 1 && arg1[end] == 1))		# A corr matrix (computed with cor())
 		return pcolor(mat2ds(arg1.z); first=first, kwargs...)
 	elseif (isa(arg1, GMTdataset))		# Arrive here when arg1 was originally a DataFrame
 		return pcolor(arg1; first=first, kwargs...)
@@ -177,7 +177,7 @@ function pcolor(D::GMTdataset; first=true, kwargs...)
 	end
 	d[:xticks] = (colnames, ang)
 	d[:yticks] = colnames[end:-1:1]
-	if (D == D' && D[1] == 1 && D[end] == 1)		# A correlation matrix (that computed with GMT.cor())
+	if (D == D' && D[1] == 1 && D[end] == 1)		# A correlation matrix (that computed with cor())
 		(is_in_dict(d, CPTaliases) === nothing) && (d[:C] = makecpt(T=(-1,1.0,0.1), C="tomato,azure1,dodgerblue4", Z=true))
 		for n = 1:size(z,2), m = 1:size(z,1)
 			m < n && (z[m,n] = NaN)					# Since the matrix is symmetric, remove the upper triangle

src/plot.jl

@@ -472,7 +472,7 @@ end
 function _bar(arg, first::Bool, d::Dict{Symbol,Any})
 
 	do_cat = ((haskey(d, :stack) || haskey(d, :stacked)) && isvector(arg) && length(arg) > 2) ? false : true
-	is_waterfall = ((val = hlp_desnany_str(d, [:stack :stacked])) !== "" && startswith(val, "water"))
+	is_waterfall = ((val = hlp_desnany_str(d, [:stack :stacked], false)) !== "" && startswith(val, "water"))
 	if (is_waterfall)
 		isa(arg, Vector) && (arg = reshape(arg, 1, length(arg)))		# Waterfall stacks must be matrices
 		(arg[1] != 0) && (arg = hcat(repeat([1.0],size(arg,1)), arg))	# If first el != 0 assume coord is missing

src/proj_utils.jl

@@ -359,7 +359,8 @@ function buffergeo(line::Matrix{<:Real}; width=0, unit=:m, np=120, flatstart=fal
 	# be able to use GMT's mapproject to find only the points that are at the exact 'width' distance from the
 	# line, but at this moment there seems to be an issue in GMT and distances are shorter. However we can
 	# do some cheap statistics to get read of the more inner points a get an almost perfec buffer.
-	Ddist2line = mapproject(D, L=(line=line, unit=:e))		# Find the distance from buffer points to input line
+	#Ddist2line = mapproject(D, L=(line=line, unit=:e))		# Find the distance from buffer points to input line
+	Ddist2line = gmt("mapproject -L+ue", D, line)		# Find the distance from buffer points to input line
 	if (isa(Ddist2line, GMTdataset))
 		d_mean = mean(view(Ddist2line.data, :, 3))
 		ind = view(Ddist2line.data, :, 3) .>= d_mean
@@ -486,15 +487,15 @@ function geodesic_long(lonlat1::VMr, lonlat2::VMr; step=0.0, np=180, proj::Strin
 	dest2, = geod(lonlat1, azA, 40050, unit=:km, proj=proj, epsg=epsg)	# Destination of a point after the perimeter
 	t = geodesic([dest1[1] dest1[2]; dest2[1] dest2[2]], step=1000, proj=proj, epsg=epsg)	# Temp geodesic arround the perimeter point
 	gmtwrite(_name, t);				# Workaround a bug in <= 6.4.0
-	x = mapproject([lonlat1[1] lonlat1[2]], L=_name)	# Assume the closest point to A corresponds to full perimeter.
+	x = gmt("mapproject -L$_name", [lonlat1[1] lonlat1[2]])	# Assume the closest point to A corresponds to full perimeter.
 	d, = invgeod(dest1, [x[4] x[5]], proj=proj, epsg=epsg)	# Distance from dest1 to the point where the geodesic does a full perimeter.
 	dtot = 39950000 + d				# Total length of this geodesic (A bit in excess probably due to geodesics not closing)
 	(step != 0) && (np = round(Int, dtot / step) + 1)	# Takes precedence over 'np' because that one has a default value.
 
 	function get_mindist(az)
 		D = geod(lonlat1, az, linspace(0,dtot,360), proj=proj, epsg=epsg)[1]
 		gmtwrite(_name, D);
-		mapproject([lonlat2[1] lonlat2[2]], L=_name)	# Find closest point to B in the passing by geodesic
+		gmt("mapproject -L$_name", [lonlat2[1] lonlat2[2]])		# Find closest point to B in the passing by geodesic.
 	end
 
 	# Using the azimuth computed from the shortest distance between A & B will always give only an aproximation
@@ -512,7 +513,7 @@ function geodesic_long(lonlat1::VMr, lonlat2::VMr; step=0.0, np=180, proj::Strin
 	az = (dt < d0) ? az + (d0 / dd) * daz : az - (d0 / dd) * daz
 
 	d = geod(lonlat1, az, dtot-distAB)[1]				# First estimate of the arrival point
-	D = mapproject(lonlat2, G=(d[1], d[2]))				# Distance between estimation and true point
+	D = gmt("mapproject -G$(d[1])/$(d[2])", (isa(lonlat2, Vector{<:Real})) ? reshape(lonlat2, 1, :) : lonlat2)	# Distance between estimation and true point
 	total_len = dtot - distAB - D[end]					# Correct 'dtot' that was estimated in excess.
 	p,a = geod(lonlat1, az, linspace(0,total_len,np))	# Possibly because the path is slightly hellicoidal.
 	prj = (proj != "" && (contains(proj,"=lon") || contains(proj,"=lat"))) ? proj : prj4WGS84

src/psxy.jl

@@ -36,7 +36,7 @@ function _common_plot_xyz(cmd0::String, arg1, caller::String, O::Bool, K::Bool,
 	cmd, is_gridtri, arg1 = parse_grid2tri_case(d, cmd, caller, is3D, isFV, O, arg1)
 
 	isa(arg1, GMTdataset) && (arg1 = with_xyvar(d, arg1))		# See if we have a column request based on column names
-	if ((val = hlp_desnany_int(d, [:decimate])) !== -999)		# Asked for a clever data decimation?
+	if ((val = hlp_desnany_int(d, [:decimate], del=true)) !== -999)		# Asked for a clever data decimation?
 		dec_factor::Int = (val == 1) ? 10 : val
 		arg1 = lttb(arg1, dec_factor)
 	end

src/spatial_funs.jl

@@ -646,7 +646,7 @@ function Base.:in(D1::GMTdataset, D2::GMTdataset)::Union{Bool, Int, Vector{Int}}
 	(!is2Pol && indPoint == 0 && indPol == 0) && error("Input arguments must have a Point and Polygon geometries, or both be Polygons.")
 	if     (indPoint == 1 && indPol == 2)  inwhichpolygon([D2], D1.data)
 	elseif (indPoint == 2 && indPol == 1)  inwhichpolygon([D1], D2.data)
-	else   GMT.contains(D1, D2)		# Returns `true` if D1 contains D2.
+	else   contains(D1, D2)		# Returns `true` if D1 contains D2.
 	end
 end
 

src/utils_types.jl

@@ -361,8 +361,8 @@ function _mat2ds(mat::Array{T,N}, txt::Union{String,Vector{String}}, hdr::Vector
 
 	# By default convert to Doubles, except if instructed to NOT to do it.
 	#(find_in_dict(d, [:datatype])[1] === nothing) && (eltype(mat) != Float64) && (mat = Float64.(mat))
-	_geom::Int = Int((geom == 0 && (2 <= length(mat) <= 3)) ? Gdal.wkbPoint : (geom == 0 ? Gdal.wkbUnknown : geom))	# Guess geom
-	(multi && _geom == 0 && size(mat,1) == 1) && (_geom = Int(Gdal.wkbPoint))	# One row with many columns and MULTI => Points
+	_geom::Int = Int((geom == 0 && (2 <= length(mat) <= 3)) ? wkbPoint : (geom == 0 ? wkbUnknown : geom))	# Guess geom
+	(multi && _geom == 0 && size(mat,1) == 1) && (_geom = Int(wkbPoint))	# One row with many columns and MULTI => Points
 	if (isempty(xx))				# No coordinates transmitted
 		if (ndims(mat) == 3)
 			coln = fill_colnames(coln, size(mat,2)-2, is_geog)
@@ -959,7 +959,7 @@ function tabletypes2ds(arg, interp=0)
 	#(arg === nothing || isa(arg, GDtype) || isa(arg, Matrix{<:Real})) && return arg
 	isdataframe(arg) && return df2ds(arg)				# DataFrames are(?) easier to deal with.
 	isODE(arg) && return ODE2ds(arg, interp=interp)		# DifferentialEquations type is a complex beast.
-	(isa(arg, GMT.Gdal.AbstractDataset) || isa(arg, GMT.Gdal.AbstractGeometry)) && return gd2gmt(arg)
+	(isa(arg, Gdal.AbstractDataset) || isa(arg, Gdal.AbstractGeometry)) && return gd2gmt(arg)
 
 	# Guesswork, it may easily screw.
 	colnames = [i for i in fields(arg) if Base.nonmissingtype(eltype(getproperty(arg, i))) <: AbstractFloat]