bwwjohnson commented 4 days ago

This is all pretty low priority and it's a bit extreme to call them bugs, putting this here as a record and for people like me who got confused. May try to do some of the Doc updating myself if I can figure out how.

Suggested fixes (the tldr)

[ ] Update Documentation to add "MGRS" subheading and put current UTM figure in there, making clear those are MGRS latitude codes
[ ] Update Documentation to explicitly state that a "+" is needed to indicate "UTM north of equator"
[ ] Adding a figure from PROJ to replace current UTM Figure
[ ] If EPSG specifications are allowed, document how to use them, as this is useful when working with this data.

--- I can't do the below ----

[ ] Add routines to catch people not adding + or - when specifying UTM projection without latitude code: i.e. -Ju34/1:xxxxx and throw an error along the lines of ERROR: UTM projection insufficiently specified - did you forget to add a latitude code or +/-? because at the moment the error message implies the -R value is wrong.
[ ] Adding routines to catch when people are using EPSG codes wrong for UTM (if this is not supported then a routine to say "this is not supported").

Description of the problems in detail

I have had a few related issues with how UTM coordinates get handled in GMT. I'm happy to handle some of the documentation changes myself if I can, but recording them all here as an issue so it's formalised.

GMT documentation on UTM refers to latitude zones that don't appear in GIS software such as QGIS or packages such as PROJ, or in the EPSG system. It turns out these latitude zones are a feature in the US Army's MGRS grid system, so not really a UTM feature. It might be good to make this clear because UTM zones are split into north (N) and south (S) of the equator in other GIS software, and confusingly MGRS has those letters as latitude zones, but both are north of the equator! Modifying the documentation to point out this difference will save a lot of confusion for newer users. Good article on the difference between UTM and MGRS here.

ISSUE 1

In the GMT documentation on UTM it says

UTM zone (A, B, 1–60, Y, Z). Use negative values for numerical zones in the southern hemisphere or append the latitude modifiers (C–H, J–N, P–X) to specify an exact UTM grid zone.

and

Furthermore, each zone is divided into latitude bands but these are not needed to specify the projection for most cases

This only works in the north if you use a "+" (-Ju+31/1:xxxxx) if you leave it blank (-Ju31/1:xxxx) and assume that counts as a "positive number" you get thrown a misleading -R parsing error:

coast [ERROR]: UTM projection insufficiently specified to auto-determine geographic region
coast [ERROR]: Option -R parsing failure. Correct syntax:
....
....
coast [ERROR]: Offending option -R353932.665723/429763.511095/5618692.59768/5692923.34264+ue

I suggest adding a line in the documentation to explicitly say "you need a plus when not specifying a MGRS latitude zone north of the equator".

ISSUE 2

If using data projected in EPSG:32631 (UTM 31N) you get some weird effects -Ju31(N-Z)-- correct for all zone specifications -Ju31(A-M) -- incorrect for all zone specifications (thinks you're in the southern hemisphere)

This is probably not a priority because if someone just looks at the chart and uses a MGRS zone that is approximately right, then it works! But it is confusing.

I suggest updating documentation to state the difference between UTM and MGRS and that MGRS codes have nothing to do with UTM as understood by most GIS software.

ISSUE 3

Specifying your projection using the EPSG code you projected it into doesn't work, but it's difficult to tell if this is because EPSG codes are not supported, or if the -R argument is wrong:

Using: -J"EPSG:32631/1:650000"

Results in a similar misleading -R parsing error for Issue 1 (you get a similar error if you try to do it with -R${gridfile} like in Example 28)

coast [ERROR]: UTM projection insufficiently specified to auto-determine geographic region
coast [ERROR]: Option -R parsing failure. Correct syntax:
...
...
coast [ERROR]: Offending option -R353932.665723/429763.511095/5618692.59768/5692923.34264+ue

I suggest routines that either catch people using EPSG codes and telling them it's not allowed, or updating documentation to reflect how they actually work, because I don't understand how this is insufficiently specified.

Code to replicate these issues:

Here is some code that replicates these issues, with a small DEM attached, projected into EPSG:32631 (UTM 31N) test_dem.zip


#!/bin/bash

# gmt version 6.4.0

file='bash-gmt-utm-test'

bounds="-R353932.665723/429763.511095/5618692.59768/5692923.34264+ue"

gmt grdconvert test_dem1.tif -Gtest_dem.grd
grid=test_dem.grd
gmt begin $file

    gmt grdimage $grid -Jx1:650000 
    gmt basemap $bounds -Jx1:650000 -B20000g20000+u"@:8:m@::" -BesNE --FONT_ANNOT_PRIMARY=10p   --MAP_GRID_CROSS_SIZE_PRIMARY=0.25c --FONT_LABEL=10p

    # "correct" MGRS ZONE -- plots correct frame
    gmt coast  $bounds -Ju31U/1:650000 -Df -B -Ia -W

    # "incorrect" MGRS ZONES (above the equator) -- plots correct frame
    gmt coast  $bounds -Ju31S/1:650000 -Df -B -Ia -W
    gmt coast $bounds -J"EPSG:32631/1:650000" -Df -B -Ia -W

    # "incorrect" MGRS ZONE (below the equator) -- plots incorrect frame
    gmt coast  $bounds -Ju31D/1:650000 -Df -B -Ia -W

    # UTM without any latitude zone -- throws a -R parsing error
    gmt coast  $bounds -Ju31/1:650000 -Df -B -Ia -W

    # UTM specified using the EPSG code
    gmt coast $bounds -J"EPSG:32631/1:650000" -Df -B -Ia -W

    # Region specified using grid
    gmt coast -R$grid -J"EPSG:32631/1:650000" -Df -B -Ia -W

gmt end show

welcome[bot] commented 4 days ago

👋 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!

Please make sure you read our Contributing Guide and abide by our Code of Conduct.

joa-quim commented 4 days ago

Hello Ben,

Thank you very much for this complete bug report.

We are a bit busy working on materials for the AGU, so advances on this may be slow. But I'm going to ask you to split this issue in parts otherwise it will quickly get confusing to address its sub-items. Another thing is, could update to GMT6.5? Or even better, can you build from source code?

Two quick answers.

EPSG's are allowed and documented at the end of -J but likely their implementation for -R...+ue is probably faulty or even missing at all.

-Ju|U doesn't need a + in mapproject, so when it complains for its absence the problem is certainly something else. See:

C:\v>echo -9 40 | gmt mapproject -Ju29/1:1 -C -F
500000  4427757.21888

C:\v>echo -9 40 | gmt mapproject -Ju+29/1:1 -C -F
500000  4427757.21888

C:\v>echo -9 40 | gmt mapproject -Ju-29/1:1 -C -F
500000  14427757.2189

C:\v>echo -9 40 | gmt mapproject -J32629
500000  4427757.21874

bwwjohnson commented 3 days ago

Hi Joaquim, thanks for your fast reply and hope AGU goes well. Thanks for pointing out that -JEPSG works in mapproject, that makes more sense that it threw a -R error then.

I will raise each of the issues I listed above as new issues in github and link the issue numbers to this report. I was lazy and just got the latest version from apt-get prior to making this report, but I will do as you suggest and build from source before creating the new issues.

Maybe you could add a documentation tag to this issue thread so the UTM/MGRS stuff distinction can be made? (I don't see a way for me to add these tags)

joa-quim commented 3 days ago

Don't know if the tags are permission dependent, but that is done on the page's right side

GenericMappingTools / gmt

UTM zones in GMT Documentation are actually MGRS zones + associated confusion #8623

Description of the problems in detail

ISSUE 1

ISSUE 2

ISSUE 3

Code to replicate these issues: