grdinterpolate -S is broken for externals passing a memory cube about gmt HOT 14 CLOSED

We could fix the macro to do this

#define gmt_M_file_is_memory(file) (file && !strncmp (file, "@GMTAPI@-", GMTAPI_PREFIX_LEN) && strlen (file) >= GMTAPI_MEMFILE_LEN)

That way it will be recognised even if there are?[layer]at the end (>= vs ==).

Next, we allow the @GMTAPI@-S-I-U-U-U-N-000035?[0] filename to be passed without scrutiny and then grdtrack will call GMT_Read_Data and then gmtapi_import_grid will need a new trick to point to that layer. I think `@GMTAPI@-S-I-U-U-U-N-000035?[layer] should be allowed. I can imagine you also passing layer n of a stack to grdimage etc.

If course, you could probably do this in Julia and pass a grid instead of a cube layer, but what fun is it in that?

from gmt.

How do a load in a cube? Tried

C = gmtread("cube.nc",cube=true);
but your gmtread docs do not mention cube, but have no error.

but grdinterpolate sees it as a grid.

grdinterpolate(C, S=[5.5 4.5])
grdinterpolate [ERROR]: Input memory reference is not a cube!

@GMTAPI@-S-I-G-G-G-N-000002 seems to come from Julia since it is in the opt array.

from gmt.

I think we have to read cubes with gdalread instead of gmtread. The case I debugged was with a cube that was already in memory after creating it from a list of geotiff grids.
Wont be able to confirm more things before the afternoon.

from gmt.

Ok, when back please give me steps so I can run either my or your command.

from gmt.

C = gmtread("cube.nc",cube=true);
but your gmtread docs do not mention cube, but have no error.

Given the way the keyword/values work no kwarg pair is actually an error if the parsing in the program doesn't catch the specific case. And because the number of possible errors is infinite what I do is that I consume all valid options and at the end if there are unconsumed I report that as a warning that is a strong indicator of user mistake. This works fine in pure modules but become quite complicated in modules that call other modules (gmtread became quite complex).

Anyway, it turns out that you can use grdinterpolate to import the cube into memory. i.e.,

C = grdinterpolate("cube.nc")

and then do the

grdinterpolate(C, S=[5.5 4.5])

from gmt.

I did that, but C seems to be a grid. Dont we want it to be a cube? How (in Julia) can I tell the dimensions of C?

from gmt.

No, you did gmtread("cube") and that reads just one layer (needs fix). What I'm saying is to grdinterpolate("cube") and that gives a cube. Possibly not intended but it does what's is needed here. Other option is to update GMT.jl to have the (incomplete) fix that I just committed and do C = gdalread("cube")

from gmt.

Sure, that works, but dont we want passing a cube layer to work down the road?

from gmt.

But passing a cube layer is just like passing a grid to grdtrack, what would be tested here? What I was expecting is to pass the entire cube at once and have all layers be tracked.

To take slices of a cube see ? slicecube on the Julia cmd or https://www.generic-mapping-tools.org/GMTjl_doc/documentation/utilities/slicecube/index.html#slicecube
basically. G = slicecube(cube, n)

from gmt.

If this is all too complicated, leave it. I think I can program the grdinterpolate -S directly in Julia.

from gmt.

I see, so you cannot do something like C{0} to get grid layer 0. I was imagining people with cubes wanting to pass slices like cube[3] to grdimage and we would then need to handle that slicing in gmtapi_import_grid.

from gmt.

Sure we can. I said it above. One can do slicecube(cube, 1) and pass first layer. It has to be a function and not just C[1] because that's just the first element of the 3D array. And it can't neither (well, it can) be C[:,:,1] that looses all grid's metadata and passes just a raw 2D matrix. But when we pass a single layer it's a grid that is being passed and not a cube with only one layer (a needless complicated entity).

from gmt.

So close this or?

from gmt.

Yes better close, though this may strike again in another circumstance. I hope that then we will still remember about this.

from gmt.