How to convert from and to UTM coordinates in Mathematica?

Question

How to use the UTM coordinate system (northing, easting) in Mathematica using the built-in functions (GeoPositionENU, GeoGridPosition, etc.). There are UTM Zones implemented in Mathematica, e.g.:

GeoProjectionData["UTMZone33"]

gives

{"TransverseMercator", {"Centering" -> {0, 15},   "CentralScaleFactor" -> 1, "GridOrigin" -> {0, 0},   "ReferenceModel" -> "WGS84"}}

From N45, E15 one should get Easting 500000 and Northing 4982950.4 (zone 33), but I can't reproduce it.

GeoGridPosition[GeoPosition[{45, 15, 0}, "WGS84"], "UTMZone33"]

gives

GeoGridPosition[{0., 4.98494*10^6, 0}, "UTMZone33"]

(Easting offset 500000 is understandable - default, but Northing is way off.)

Any help would be appreciated.

I've wanted to know this for a while too. Right now, I just use the function included in this blog post to go from LatLon to UTM: http://blog.wolfram.com/2009/04/17/mapping-gps-data/ — kale, Jan 08 '13 at 14:56
The output of GeoProjectionData is not correct for the UTM system: the value of GridOrigin is wrong (it should be {500000,0} and so is the value of CentralScaleFactor (it should be 0.9996). — whuber, Jan 08 '13 at 17:04
@kale: the short version of the function would therefore be: ll2utm[coord_] := GeoGridPosition[ GeoPosition[coord, "WGS84Original"], {"UTMZone33", "CentralScaleFactor" -> 0.9996, "GridOrigin" -> {500000, 0}}] [[1, {1, 2}]] and then ll2utm[{45, 15}] // Round gives {500000, 4982950}. — Boocko, Jan 09 '13 at 11:08
@Boocko, Except, you would need to include the logic to figure out what zone you were in and build a string "UTMZone"<>ToString@zone. — kale, Jan 09 '13 at 15:19

score 6 · Accepted Answer · edited Jan 08 '13 at 21:15

6

As in comments by WReach, the correct answer is

GeoGridPosition[ 
  GeoPosition[{45, 15, 0}, "WGS84"], 
  {"UTMZone33", "CentralScaleFactor" -> 0.9996, "GridOrigin" -> {500000, 0}}
]

The problem lies with the scale factor: 0.9996 and the grid origin

If you are familiar with projection systems used in GIS systems, you can check out http://spatialreference.org/ref/epsg/32633/

Details are here:

PROJCS["WGS 84 / UTM zone 33N",
    GEOGCS["WGS 84",
        DATUM["WGS_1984",
            SPHEROID["WGS 84",6378137,298.257223563,
                AUTHORITY["EPSG","7030"]],
            AUTHORITY["EPSG","6326"]],
        PRIMEM["Greenwich",0,
            AUTHORITY["EPSG","8901"]],
        UNIT["degree",0.01745329251994328,
            AUTHORITY["EPSG","9122"]],
        AUTHORITY["EPSG","4326"]],
    UNIT["metre",1,
        AUTHORITY["EPSG","9001"]],
    PROJECTION["Transverse_Mercator"],
    PARAMETER["latitude_of_origin",0],
    PARAMETER["central_meridian",15],
    PARAMETER["scale_factor",0.9996],
    PARAMETER["false_easting",500000],
    PARAMETER["false_northing",0],
    AUTHORITY["EPSG","32633"],
    AXIS["Easting",EAST],
    AXIS["Northing",NORTH]

More explanations of the UTM system can be found here: "The central meridian in each UTM zone has a scale factor of 0.9996"

edited Jan 08 '13 at 21:15

Sjoerd C. de Vries

65,815
14
188
323

answered Jan 08 '13 at 16:01

Tuku

486
2
4

2

I feel this falls a bit short as an answer. Could you elaborate somewhat? – Sjoerd C. de Vries Jan 08 '13 at 16:29
1

I agree with your diagnosis, although the easting should also be 500000 as the OP observes. You might want to mention that the desired result can be obtained using the expression GeoGridPosition[ GeoPosition[{45, 15, 0}, "WGS84"], {"UTMZone33", "CentralScaleFactor" -> 0.9996, "GridOrigin" -> {500000, 0}}]. – WReach Jan 08 '13 at 19:18
Thank you, this works. But to convert back, the simple GeoPosition on previous output should work, but I get errors. Looks like the GeoPosition is not happy with {"UTMZone33",...}. – Boocko Jan 09 '13 at 10:51
1

works in this way: GeoProjectionData["UTMZone33"] gives {"TransverseMercator", {"Centering" -> {0, 15}, "CentralScaleFactor" -> 1, "GridOrigin" -> {0, 0}, "ReferenceModel" -> "WGS84"}}, then use this: GeoPosition@ GeoGridPosition[ GeoPosition[{45, 15, 0}, "WGS84"], {"TransverseMercator", "Centering" -> {0, 15}, "CentralScaleFactor" -> 0.9996, "GridOrigin" -> {500000, 0}, "ReferenceModel" -> "WGS84"}] – Tuku Jan 09 '13 at 15:30
@Tuku that works flawlessly. Thanks. – Boocko Jan 10 '13 at 10:27

score 0 · Answer 2 · answered May 04 '19 at 00:16

The settings appear to have been "fixed" (at least for Mathematica version 12).

Using Los Angeles, CA as an example:

(* Convert from LatLong to UTM *)
(* You need to know what UTM Zone you want *)
lat = 34.052234;
long = -118.243685;
elev = 86.654;
utmZone = "UTMZone11";
utm = GeoGridPosition[GeoPosition[{lat, long, elev}], utmZone]
(* GeoGridPosition[{385215.34335191007`,3.7686452431931966`*^6,86.654`},"UTMZone11"] *)

(* Convert from UTM to LatLong *)
utmZone = "UTMZone11";
utmEasting = 385215.34335191007;
utmNorthing = 3768645.2431931966;
elev = 86.654;
GeoPosition[GeoGridPosition[{utmEasting, utmNorthing, elev}, utmZone]]
(* GeoPosition[{34.052234000000006`,-118.243685`,86.654`}] *)

How to convert from and to UTM coordinates in Mathematica?

2 Answers2

Linked

Related