Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

grdcut -N on larger region results in no grid extension #31

Open
apasto opened this issue Oct 9, 2022 · 2 comments
Open

grdcut -N on larger region results in no grid extension #31

apasto opened this issue Oct 9, 2022 · 2 comments

Comments

@apasto
Copy link

apasto commented Oct 9, 2022

Hello everyone,

Description of the problem

I am trying to extend a grid onto a new region, including but exceeding the original boundaries, using grdcut with the -N option.
However, the returned grid structure retains the original region unchanged: no nan-padding occurs.

Full script that generated the error

This minimal example reproduces the issue:

region_0 = [-5, 5, -10, 10];
n_el = [21, 41];

x = linspace(region_0(1), region_0(2), n_el(1));
y = linspace(region_0(3), region_0(4), n_el(2));
inc = [x(2) - x(1), y(2) - y(1)];
z = ones(length(y), length(x));
header = [region_0, min(z(:)), max(z(:)), 0, inc];

G_0 = gmt('wrapgrid', z, header);

region_1 = [-6, 6, -11, 11];
G_1 = gmt(['grdcut -Vd -N -R', sprintf('%g/%g/%g/%g', region_1)], G_0);

assert(all(region_1 == G_1.range(1:4)))

This also fails with fill values other than nan, e.g. -N0.

G_1 region is left unchanged, still equal to G_0 extents:

>> G_0.range
ans =    -5     5   -10    10     1     1
>> G_1.range
ans =    -5     5   -10    10     1     1

Writing G_0 to file and calling gmt grdcut -N outside Matlab, from the shell, the nan-fill extension works as expected.

Full error message

grdcut [DEBUG]: Got regular w/e/s/n for region (-6/6/-11/11)
grdcut [INFORMATION]: Processing input grid
grdcut [DEBUG]: Set_Object for family: 1
grdcut [INFORMATION]: Cartesian input grid
grdcut [INFORMATION]: Requested subset exceeds data domain on the left side - nodes in the extra area will be initialized to -nan(ind)
grdcut [INFORMATION]: Requested subset exceeds data domain on the right side - nodes in the extra area will be initialized to -nan(ind)
grdcut [INFORMATION]: Requested subset exceeds data domain on the bottom side - nodes in the extra area will be initialized to -nan(ind)
grdcut [INFORMATION]: Requested subset exceeds data domain on the top side - nodes in the extra area will be initialized to -nan(ind)
grdcut [DEBUG]: Set_Object for family: 1
grdcut [INFORMATION]: Cartesian input grid
grdcut [DEBUG]: gmtapi_expand_headerpad: No pad adjustment needed
grdcut [DEBUG]: Chosen boundary condition for all edges: natural
grdcut [INFORMATION]: gmt_grd_BC_set: Set boundary condition for all edges: natural
grdcut [INFORMATION]: gmt_grd_BC_set: Set boundary condition for left   edge: natural
grdcut [INFORMATION]: gmt_grd_BC_set: Set boundary condition for right  edge: natural
grdcut [INFORMATION]: gmt_grd_BC_set: Set boundary condition for bottom edge: natural
grdcut [INFORMATION]: gmt_grd_BC_set: Set boundary condition for top    edge: natural
grdcut [DEBUG]: Object ID 17 : Registered Grid Memory Reference f3c22b56c0 as an Input resource with geometry Surface [n_objects = 3]
grdcut [DEBUG]: Successfully duplicated a Grid
==> 3 API Objects at end of GMT_Duplicate_Data
--------------------------------------------------------
K.. ID RESOURCE.... FAMILY.... ACTUAL.... DIR... S O M L
--------------------------------------------------------
* 0 15   f3c22bfe30 Grid       Grid       Input  0 Y N 0
* 1 16   f3c22bfbc0 Grid       Grid       Output 0 Y Y 0
* 2 17   f3c22b56c0 Grid       Grid       Input  0 Y N 1
--------------------------------------------------------
grdcut [INFORMATION]: File spec:	W E S N dx dy n_columns n_rows:
grdcut [INFORMATION]: Old:grdcut [INFORMATION]: 	-5	5	-10	10	0.5	0.5	21	41
grdcut [INFORMATION]: New:grdcut [INFORMATION]: 	-6	6	-11	11	0.5	0.5	25	45
grdcut [DEBUG]: GMT_Write_Data: Writing Grid to memory object 16 from object 17 which transfers ownership
grdcut [DEBUG]: gmtapi_begin_io: Output resource access is now enabled [container]
grdcut [DEBUG]: gmtapi_export_data: Messenger dummy output container for object 16 [item 1] freed and set resource=data=NULL
grdcut [DEBUG]: gmtapi_export_grid: Passed ID = 16 and mode = 0
grdcut [INFORMATION]: Referencing grid data to GMT_GRID memory location
grdcut [DEBUG]: Chosen boundary condition for all edges: natural
grdcut [INFORMATION]: gmt_grd_BC_set: Set boundary condition for all edges: natural
grdcut [INFORMATION]: gmt_grd_BC_set: Set boundary condition for left   edge: natural
grdcut [INFORMATION]: gmt_grd_BC_set: Set boundary condition for right  edge: natural
grdcut [INFORMATION]: gmt_grd_BC_set: Set boundary condition for bottom edge: natural
grdcut [INFORMATION]: gmt_grd_BC_set: Set boundary condition for top    edge: natural
grdcut [DEBUG]: GMT_End_IO: Output resource access is now disabled
grdcut [DEBUG]: Set_Object for family: 1
==> 3 API Objects at end of GMT_Write_Data
--------------------------------------------------------
K.. ID RESOURCE.... FAMILY.... ACTUAL.... DIR... S O M L
--------------------------------------------------------
* 0 15   f3c22bfe30 Grid       Grid       Input  0 Y N 0
* 1 16   f3c22b56c0 Grid       Grid       Output 2 Y N 0
* 2 17   f3c22b56c0 Grid       Grid       Input  0 N N 1
--------------------------------------------------------
==> 3 API Objects at end of GMTAPI_Garbage_Collection entry
--------------------------------------------------------
K.. ID RESOURCE.... FAMILY.... ACTUAL.... DIR... S O M L
--------------------------------------------------------
* 0 15   f3c22bfe30 Grid       Grid       Input  0 Y N 0
* 1 16   f3c22b56c0 Grid       Grid       Output 2 Y N 0
* 2 17   f3c22b56c0 Grid       Grid       Input  0 N N 1
--------------------------------------------------------
grdcut [DEBUG]: gmtlib_unregister_io: Unregistering object no 17 [n_objects = 2]
==> 2 API Objects at end of GMTAPI_Garbage_Collection exit
--------------------------------------------------------
K.. ID RESOURCE.... FAMILY.... ACTUAL.... DIR... S O M L
--------------------------------------------------------
* 0 15   f3c22bfe30 Grid       Grid       Input  0 Y N 0
* 1 16   f3c22b56c0 Grid       Grid       Output 2 Y N 0
--------------------------------------------------------

System information

  • Operating system: Windows 10
  • Version of GMT: 6.4.0
  • Version of GMTMex: 2.1.0
@welcome
Copy link

welcome bot commented Oct 9, 2022

👋 Thanks for opening your first issue here! Please make sure you filled out the template with as much detail as possible. We appreciate that you took the time to contribute!

@joa-quim
Copy link
Member

joa-quim commented Oct 9, 2022

Yep, can confirm this behavior in Julia (GMT.jl) too, which functions in a similar way to the gmtmex. This points to the issue being in the GMT lib itself.

Thanks for the report.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

2 participants