grdcontour -C & externals is wrong

[joa-quim]

julia> C = makecpt(T=(-7,9,2), Vd=1);
        makecpt  -T-7/9/2

Reading grid from file it's (nearly) OK. Nearly because the aliasing is clearly visible

Now do the same but using data in memory.

julia> Gp = gmtread("peaks.grd");

julia> grdcontour(Gp, C=C, N=true, savefig="p_bad.png", show=1, Vd=1)
        grdcontour  -JX12c/0 -Baf -BWSen -C -N -P -K > C:\TMP\GMTjl_tmp.ps

Finally, another bug is when trying to do the same as a modern one liner

 gmt grdcontour peaks.grd -JX12c/0 -Baf -BWSen -Cc2.cpt -N -png C2

Here it creates two separate figures. One with the colors and another with the contours.

The grid is here but not zipped. Just remove the .zip (needed to cheat github)
peaks.grd.zip

[PaulWessel]

Think the -N predates one-liners...
Will have a look but busy with GMTSAR presentation now and getting tired...

[PaulWessel]

Works in both classic and modern mode but not for one-liner. Problematic with one-liner syntax here since it ends up being two module calls.

[joa-quim]

grdcontour makes a copy of the grid for contouring!!!!!!!!

if ((G_orig = GMT_Duplicate_Data (API, GMT_IS_GRID, GMT_DUPLICATE_DATA, G)) == NULL) Return (GMT_RUNTIME_ERROR); /* Original copy of grid used f

[PaulWessel]

grdcontour clobbers the grid, subtracting out values etc. etc.

[joa-quim]

So grdcontour -C -N calls first grdview, that does a FIRST copy of the grid, and changes the original grid while also tracing contours (though they do not show up in image). Next it calls itself where ANOTHER copy of the grid but this time, when from externals that passed a memory grid, it reads from memory the array that was previously modified by grdview and logically the contours are all wrong. This doesn't happen when one passes a file name because the second time it reads the correct array because it's reading it from disk.

So, even if grdview alone is screwing because it's changing the external array without asking if it can. Not a big deal for me (for the time being) in Julia because I had to make a copy due to the transposition story.

[PaulWessel]

yes, there needs to be a check on the alloc_mode for those grids and if it detects external then make a copy. Sorry, deep in other shits and meetings.

[joa-quim]

But it's already making a copy.

[PaulWessel]

OK, but I think I need to try this in matlab or something to learn what it is doing and if we can be smarter.

[joa-quim]

One clever thing to do is not call grdcontour a second time. It is silly to do grdcontour -C<val> -N<cpt> because the idea of -N is that the contours fall exactly on the color transitions. Then, since grdview already computes the contours, those contours could be sent back to grdcontour that would use them instead of computing them a second time.

[PaulWessel]

Well, if you think it is good use of your time to figure out how grdview module can somehow pass contours back to grdcontour in a form that grdcontour can understand (it traces and does lots of things differently than grdview, not to mention annotations) then knock yourself out and make a branch. I am willing to fix the current problem that fails in Julia and presumably Matlab, but seriously? Rewriting lots of code for what? Saving someone 2 Mb of memory?

[PaulWessel]

If we did not care about contour annotations then there would be no need for the whole -N business in the first place.

[joa-quim]

First example in Matlab's contourf had no annotations
https://www.mathworks.com/help/matlab/ref/contourf.html

[joa-quim]

Rewriting lots of code for what? Saving someone 2 Mb of memory?

That will be our end of times discussion.

Never contested the ratio amount-of-work/memory-saved ratio (specially when it's not my work) but when 2 GB count to make the difference, "Hoops it doesn't work in GMT", I'll still think memory savings are important.

[PaulWessel]

I understand, but it would truly require someone to sit in Julia and waning to make a colored contour-map of earth_relief_30s or something. That only makes since if you are making a poster map at large scale. Is that the typical use?

[joa-quim]

I'll manage to work around with separate grdview + grdcontour calls in Julia so it will work with < 6.1 version. It's the pure GMT solution that I'm arguing about.

[joa-quim]

See, efficiency counts for some people

GenericMappingTools/GMT.jl#396

[stale]

This issue has been automatically marked as stale because it has not had activity in the last 90 days. It will be closed if no further activity occurs within 7 days. Thank you for your contributions.

[stale]

This issue has been automatically marked as stale because it has not had activity in the last 90 days. It will be closed if no further activity occurs within 7 days. Thank you for your contributions.