ANN: New GMT interface

Hi,

GMT.jl got a new tag (0.2) that introduces a new higher-lever type of interface to the GMT package.

While this new interface makes it more similar to Plots.jl calls, and indeed there is some overlap, I would like to recall that GMT.jl ads a geographically dimension to your plots.

Contributions, specially regarding on how to create a better documentation without repeating the original GMT docs, are welcome.

4 Likes

Great! This package has always interested me, but the old syntax scared me. I’ll be checking it out more closely, now.

Naive question: why not separate the backend (manipulating geographic and Cartesian data sets) from the plotting frontend, add some nice plot recipes and let the actual plotting be handled by Plots.jl (or Makie.jl)? Won’t that make your life easier?
Sorry if I overstepped…

The idea to split plotting into dedicated tasks/libraries e.g. coordinate transformation is not a new one…
Keep in mind, that GMT is a big-block of SW and written in C. Splitting out something there is not straight forward.

1 Like

Sorry, so the plotting is intertwined with the “backend” calculations and therefore it’s hard to separate the two?

GMT is essentially made of two small libraries (total ~ 4 Mb) and a driver executable (the gmt.exe). The Julia wrapper replicates this command line driver functionality and, in a single shot, gives access to whole GMT. The compact syntax that scares some people is needed given the large level of details one may control combined with a large number of options a module can perform.

The new upper level off access hides away some of this complexity and could in principle be added as a Plots backend but this would be a large work (I can’t even say I understand the Plots code) and furthermore I’m not even convinced I liked the Plots syntax very much with its arguments at times symbols, other numbers and other text strings, sorry.

1 Like

I see, cool. Thanks for the clarification :+1:

GMT is mostly known for its graphics/images quality but it actually is a very powerful data processing tool too and this two features are intimately connected so the notion of splitting them sounds a bit weird to me.

I have to apologize, i was a little bit confused as i was in a parallel discussion about plotting…

For me split means: make internal processing steps of GMT available in the GMT.jl adaptation to use it in other plotting.

make internal processing steps of GMT available in the GMT.jl adaptation to use it in other plotting.

But you can. GMT communicates in-out via its data types that can simply be accessed by other codes. In particular, there is a large potential of cooperation with Images.jl that could benefit enormously from the GMT capacity to read all sort of images formats (which it does via GDAL).

Wow! These examples look great! Thank you for your work, I will definitely give it a try.

I always wanted to be able to plot the results of spatial estimation problems in my GeoStats.jl package on spherical coordinates like in your first example. Suppose I have a 2D Julia array A with values A[i,j] and coordinates in degress lon[i], lat[j], how would you plot the array of values A as a heatmap on the globe?

Reading the docs now…

I would also love to see some sort of integration with Plots.jl and PlotRecipes.jl specially. Ideally I could define recipes to plot on the globe without depending on GMT.jl.

If your 2D array A is already on a regular grid the best thing would be to wrap in a GMTgrid type and than use it in imshow ord grdimage. imshow is currently just an alias+some_guessings of grdimage but the idea is to extend it to also deal with the 3D grdview. If A is not a regular grid you create one with surface like in the contours example. Remember that elaborate GMT plots are made with the layer cake model where you keep adding new layers, likely from different modules, to your plot and at the end display it.

You cannot NOT depend on GMT.jl because it’s the working beast that knows how to communicate with the GMT library.

@joa-quim can you maybe share a minimum working example with the 2D Julia array I mentioned?

OK, a fake example that should give you a Mercator map of Gulf of Angola random colorized ocean displayed from PostScript

x = collect(linspace(0,20,256));     # Longitude vector between 0 and 20 degree East
y = collect(linspace(0,20,256));     # Latitude vector between 0 and 20 degree North
A = rand(Float32, 256, 256);         # Your A array (note the dimensions nx x ny)

G = GMT.GMTgrid("", "", [minimum(x), maximum(x), minimum(y), maximum(y), minimum(A), maximum(A)], [x[2]-x[1], y[2]-y[1]], 
	            0, NaN, "", "", "", "", x, y, A, "x", "y", "z", "");

grdimage(G, proj="M10c", frame="a")    # This will use the default color scale
coast!(region=[0 20 0 20], land="brown", show=true)   # Pay attention to use 'coast!' and not 'coast' because we are appending to an existing plot

This is what I get when I run the code:

Does it look correct?

Nevermind, I fixed the example to something more intuitive.

I need to read the docs to get something out of GMT, it is not easy to catch up the syntax. It would be great to have a set of working examples discussing the features.

Well, everybody has to read the docs. Either for learning or for advanced usage. I’m a GMT developer and have to read the docs all the time for non trivial things. GMT is just too big for anyone to know it by hart.

People complain about the syntax but let me explore a bit just the frame="a" option above. It means draw axes with automatic labeling. Now, if we want we can control the label font, the label font size the spacing from the labeling from the axes, the tick marks length, which axes to plot, which axes to annotate and more. So, for example if one want to plot all 4 axes but annotate only the left and bottom, one would use frame="a WSen". But for plotting only the left axes (the West) it would be frame="a W". And this still by letting the automatic labeling decide what to do for us. Imagine what it would take to spell all this options. Just for the WSen one would have to do something like.

annotate_west=true, annotate_south=true, draw_only_east=true, draw_only_north=true

expand this to other options and it would take a page just to fine control a simple plot. This is the challenge of creating an upper level interface. The other challenge is on how to document it. I can’t just duplicate all the GMT documentation in GMT.jl. So far what I came out with was to put a link after each option directing to the GMT corresponding documentation, See for example the coast module. But, as mentioned in the opening post of this thread, I most than welcome suggestions on how to improve both the package and its docs.

@joa-quim thanks for sharing the full reference, it is helpful.

I have a suggestion for you. Try to concentrate your efforts on writing your own tutorials with GMT.jl for common tasks people always ask you. Tutorials are much more effective and helpful than tons of text with formal specifications or a hand of disconnected examples. Specially when it comes to a plotting package (no one has time to learn the internals of complex transformations producing the graphics on the display).