GMT mapproject: Corner geographic coordinates of a skewed Local XY grid

Q&A
1

I have an ASCII comma-separated grid file in Local XY coordinates (km units) that is 800km long in the Y dimension and 500 km wide in the X dimension. It is centered upon 33.535°N and 134.675°E and is rotated -20°N. The top few lines looks like this (Y in the first column, X in second column):

% head -5 test_file
-400.000,-250.000,1.6
-400.000,-248.000,1.6129
-400.000,-246.000,1.62602
-400.000,-244.000,1.63934
-400.000,-242.000,1.65289

I would like to obtain the lat/lon geographical coordinates in decimal degrees for each of the four corners of this grid. For that, a colleague suggested using GMT mapproject like this:

% mapproject test_file -JoA134.675/33.535/70/1:1 -Fk -Rk-400/400/-250/250 -C -V -:

That command produces a header that contains the extreme Local XY min/max values (“Input extreme values”) and the corresponding extreme geographic coordinates (“Output extreme values”) , followed by three-column output with numbers that I do not understand. So, I’m obviously making a mistake.

Presumably, once I’ve got mapproject working properly and can actually get decimal degrees to be output, the output records could be pushed through “gmtinfo” to find the corner coordinates of the grid.

Sorry for being a dummy but could someone point me in the right direction with the mapproject command? Thank you very much for any help that is given!
Andrew

2

Are you sure your data is in oblique Mercator projection? Because -JoA134.675/33.535/70/1:1 is that. Your mapproject command project all the points but lacks the -I to tell it to do the inverse projection. And I think you should use -Joa instead of -JoA.
Anyway, if you run gmt gmtinfo in your data it will report only the min/max of each column and that’s what you want to pass to mapproject to just get the corner coordinates.

3

Thank you, Joaquim - your suggestions worked! Yes, I believe that the local XY values are in oblique Mercator projection and, yes, I had not realised that I needed the -I inverse option. I made a new file containing just the corner coordinates. The -JoA and -Joa options give identical output. Anyway, here’s what I got, with a resultant map image!

% seafloor% cat nankai_corners_localxy.txt
-250,400
-250,-400
250,-400
250,400

% mapproject nankai_corners_localxy.txt -JoA134.675/33.535/70/1:1 -Fk -Rk-400/400/-250/250 -C -V -I > nankai_corners_degrees.txt

% cat nankai_corners_degrees.txt
130.542891311 36.0779857363
133.665842859 29.3839321025
138.565161289 30.8627777068
135.786023803 37.6770565271

And, when plotted in GeoMapApp, the grid location is correct.nankai_grid

Thank you for answering my question!
Andrew