I have a grid in XYZ which is oblique; so points not aligned with X/Y axis.
I want to use xyz2grd to format this into a netCDF, but how should I assign the -I?
If the grid was aligned with the X/Y axis it would be really easy. In this case -I12.5/25.
But it’s not, so how do you do it in this case?
Can we imagine your dots inside a square pixel (Pixel registration)?
I guess you could, but does it matter?
Not sure if my actual data use gridline or pixel node registration.
Do you know how many degrees the grid is rotated by?
GMT does not support non-NS oriented grids, so basically that means you should not use xyz2grd in this case. Likely, you’ll have to interpolate the data again.
But if the angle of rotation is known, then the points can be rotated and then run through xyz2grd
I don’t know the angle of rotation exactly, but I can give a reasonable estimate.
(Now that you mention it, how do you rotate a set of points with GMT..)
The grid in question is an interpretation from a seismic cube that I have exported. So it’s basically a Z-value for each bin for all inlines and crosslines. Maybe I have to use rows and columns instead of coordinates, ala Format output with grd2xyz - Q&A - GMT Community Forum.
Thanks for all input!
But when you do that you are moving a certain z from one position to another, without interpolation. That dosen’t seem correct to me.
What are the projected units for the suspected oblique projection?
Is it not possible to just make a Cartesian grid in these units of the data and then e.g plot it using -Jx and the corresponding -J(a|b|o) projection for example?
It’s in UTM31, so meter.
I did create(/reformat) a Cartesion grid the normal way (ref. -I12.5/25 in my original post), but this creates a grid with gaps, since my ‘region of data’ is rotated. So as far as xyz2grd is concerned, the spacing of -I12.5/25 is not correct (but it would be if you rotated the data points to align with N-S, ref. timhume’s answer above).
But if they were UTM they would be aligned NS. What you seem to have is a rotated local reference system.
Hmm, you’re probably right. That’s what it’s called.
I’m not sure, but maybe you could try reprojecting your data and then using xyz2grd.
This is an example from gmt project:
To make a Cartesian coordinate transformation of mydata.xy so that the new origin is at 5,3 and the new
axis (p) makes an angle of 20 degrees with the old axis, use:
gmt project mydata.xy -C5/3 -A70 -Fpq > mydata.pq
The idea of rotation isn’t to move the points at all, but to create a new X and Y axis which the data are aligned with. Of course one would want to add appropriate metadata to the grd file so that an unsuspecting user didn’t think the axes were geographic north-south and west-east.
I think my original reply was confusing. I should have said rotate the axes rather than rotate the data.
Still think the easiest is to regrid the data and aply it an outside mask to let it clearly show the rotated nature of the data.
Curiously, when many hears ago I used to work with aeromagnetic surveys most of the grids where not NS and we had software that knew how to dean with local rotated reference systems. Nowadays, I think GeoTIFF has metadata to store angle of rotation (and the GDAL affine matix has place for it) but, again, GMT does not know anything about rotated CS.
Thank you all for very good input. I’ll think and see what I do. I can always just interpolate the data, as @Joaquim suggested. I just hoped there was a more elegant way of doing it.
Thanks!
I thought it might make sense to specify (e g to xyz2grd for gridding) that the input data are not geographical, presumably using -fic like it says in the docs here: https://docs.generic-mapping-tools.org/latest/gmt.html#f-full - unless this gets detected automatically (a not very well educated guess).
For a small local region it however might not make a lot of difference, as far as i can understand.
Maybe. I’ll just have to play around of test a bit. Thanks!
