轻松实现坐标转换
不同地理位置系统转换入门
位置服务 —— 包括基于 GPS 的导航系统和地图站点（如 Google Maps 和 Yahoo! Maps）—— 现在深受客户欢迎。很多企业已经利用了某些位置感知服务，而更多的用户将加入到这个行列中来，因为他们已认识到该服务带来的优势和潜能。在 2006 年，Garter 就曾表示，“位置感知服务在未来两到五年内将成为主流”，并且已经有 “越来越多的组织部署了位置感知移动业务应用程序。”（请参阅 参考资料，获得该报告的链接）。
当企业决定实现某种位置感知应用程序时，编写此类应用程序的任务最终都落在开发人员的身上。构建位置感知服务涉及多种任务，或大或小，其中一项任务（相对较小）可能要将一种系统坐标转换为另一种系统坐标。本文将演示执行此类转换的代码，从而帮助您节省大量的工作。
两种不同的坐标系统
在详细研究本文代码之前，首先需要讨论即将处理的代码所属的坐标系统：较为熟悉的经纬度系统和统一横轴墨卡托投影系统（Universal Transverse Mercator，UTM）。我们还要提到以 UTM 为基础的军事格网参考系 (MGRS)。
经纬度系统
经纬度系统可能是最为人熟知的地理坐标设计方法。它使用两个数值表示位置。纬度 表示从地球中心到地球表面东西方向线之间的角度。经度 指从地球中心到地球表面南北方向线之间的角度。经纬度可以表示为十进制角度（DD），或表示为度、分、和秒（DMS）；后者的格式可表示为诸如
49°30'00" S 12°30'00" E。这是 GPS 设备使用的典型格式。
地球以赤道（0° 纬线）为界，分为南半球和北半球，又以 0° 经线（从南极到北极的假想线，通过英国的格林威治市）为界分为东西半球。北半球的纬度从 0 度到 90 度，而南半球的纬度从 0 度到 -90 度。 东半球的经度范围从 0 度到 180 度，西半球的经度范围为 0 度到 -180 度。
举例说明，坐标 61.44，25.40（使用 DD 单位）或 61°26'24''N，25°23'60''E（使用 DMS 单位）位于芬兰南部。坐标 -47.04, -73.48（使用 DD 单位）或 47°02'24''S，73°28'48''W（使用 DMS 单位）位于智利南部。图 1 展示了表面覆盖经纬线的地球：
图 1. 地球表面覆盖了经纬线
请参阅 参考资料，获得更详细的内容。
统一横轴墨卡托投影
UTM 坐标系统使用基于网格的方法表示坐标。UTM 系统将地球分为 60 个区，每个区基于横轴墨卡托投影。绘图法中的地图投影方法可以在平面中表示一个两维的曲面，例如一个标准地图。图
2 展示了一个横轴墨卡托投影：
图 2. 横轴墨卡托投影
UTM 经度区范围为 1 到 60；其中 58 个区的东西跨度为 6°（稍后详细讨论另外两个区）。经度区涵盖了地球中纬度范围从 80°S 到 84°N 之间的所有区域。
一共有 20 个 UTM 纬度区，每个区的南北跨度为 8°；使用字母 C 到 X 标识（其中没有字母 I 和 O）。A、B、Y、Z 区不在系统范围以内；它们覆盖了南极和北极区。图
3 展示了欧洲的 UTM 区。从图 3 中可看到两个非标准的经度区：32V 区被扩展为覆盖整个挪威的南部，而 31V 区被缩小，所以只覆盖了一片汪洋大海。
图 3. 欧洲地区的 UTM 区
UTM 坐标的表示格式为：经度区纬度区以东以北，其中以东 表示从经度区的中心子午线的投影距离，而以北 表示距离赤道的投影距离。这个两个值的单位均为米。举例来说，使用
UTM 表示经/纬度坐标 61.44，25.40 的结果就是 35 V 414668 6812844；而 经/纬度坐标 -47.04，-73.48 的表示结果为 18 G 615471 4789269。
请参阅 参考资料，获取更多有关 UTM 和 Traverse Mercator
投影的信息。
军事格网参考系
MGRS 是北约（NATO）军事组织使用的标准坐标系统。MGRS 以 UTM 为基础并进一步将每个区划分为 100 km × 100 km 的小方块。这些方块使用两个相连的字母标识：第一个字母表示经度区的东西位置，而第二个字母表示南北位置。
例如，UTM 点 35 V 414668 6812844 等价于 MGRS 点 35VMJ1466812844。该 MGRS 点精度为米，使用 15 个字符表示，其中最后 10 个字符表示指定网格中的以东和以北的值。可以使用 15 个字符表示 MGRS 值（如前例），也可表示为 13、11、9 或 7 个字符；因此，所表示的值的精度分别为 1 米、10 米、100 米、1,000 米或 10,000 米。
本文并未对 MGRS 进行详细说明，但是本文的下载代码包含了经纬度坐标和 MGRS 坐标之间的转换。请参阅 参考资料，获得更多信息。
坐标转换
确定地球上某个位置的经度和纬度坐标的最低需求是，你至少能够看到星星和太阳，并具备一个六分仪和能够显示 GMT 时间的时钟 T。根据空中某个物体与地平线之间的角度可以确定纬度，然后根据地球旋转计算出经度。本文并未详细讨论这些细节（想要了解更多请参阅 参考资料），相反，我们假设您已经具有
DD、DMS 或 UTM 格式的坐标。
在十进制角度和度/分/秒格式之间进行转换
DD 和 DMS 坐标格式之间的转换非常简单。下面给出了 DD 到 DMS 的转换公式：
|
1
2
3
4
5
6
|
DD:
dd.ff
DMS:
dd mm ss
dd=dd
mm.gg=60*ff
ss=60*gg
|
这里的 gg 代表计算的小数部分。负纬度表示位于南半球（S）的位置而负经度表示西半球（W）的位置。例如，假设您具有一个
DD 格式的坐标 61.44，25.40。按照下面的公式将其转换：
|
1
2
3
|
lat
dd=61
lat
mm.gg=60*0.44=26.4
lat
ss=60*0.4=24
|
以及：
|
1
2
3
|
lon
dd=25
lon
mm.gg=60*0.40=24.0
lon
ss=60*0.0=0
|
因此，转换为 DMS 格式的坐标变成了 61°26'24''N 25°24'00''E。
将 DMS 转换为 DD 格式的公式如下所示：
|
1
2
3
4
|
DD:
dd.ff
DMS:
dd mm ss
dd.ff=dd
+ mm/60 + ss/3600
|
注意，南半球（S）的位置为负纬度，西半球（W）位置为负经度。
现在将 DMS 格式坐标 47°02'24''S 和 73°28'48''W 转换为 DD 格式的坐标：
|
1
2
|
lat
dd.ff= - (47 + 2/60 + 24/3600 )=-47.04
lon
dd.ff= - (73 + 28/60 + 48/3600)=-73.48
|
转换后的 DD 格式的坐标为 -47.04 和 -73.48。
在经纬度和 UTM 坐标之间进行转换
十进制坐标可通过一个六分仪和一个记时计确定，与此不同的是，必须通过计算才能确定 UTM 坐标。虽然这些计算无非是最基本的三角形和代数计算，但是所使用的公式非常复杂。如果您阅读了 “The Universal Grids: Universal Transverse Mercator (UTM) and Universal Polar Stereographic (UPS)”（参阅 参考资料 获得链接），就知道它有多复杂了。
本文没有给出 UTM 转换公式，但是可从下面一节中给出的源代码中窥探一二，更多信息请参阅 参考资料 提供的链接。
使用 Java 代码转换坐标
本节介绍了执行坐标转换（十进制角度和 UTM）的库类的源代码。该 Java 类名为 com.ibm.util.CoordinateConversion；其思想是构建一个提供转换方法的类。该类包含实际执行转换的内部类；如果需要的话，可以从 CoordinateConversion 类中重构内部类，从而创建一个库包或向现有包添加类。该类执行的转换精度低于
1 米。
CoordinateConversion 的源代码包含大约 750 行代码，因此本文没有全部显示。以下小节描述了有关方法，本文的 下载 小节中附带了完整的源代码。
CoordinateConversion
CoordinateConversion 是主类，它被实例化为在需要是执行坐标转换。清单 1 展示了相关的公共方法，以及 CoordinateConversion 类中包含的私有内部类：
清单 1. CoordinateConversion
|
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
|
public
class CoordinateConversion
{
public
CoordinateConversion()
{
}
public
double[] utm2LatLon(String UTM)
{
UTM2LatLon
c = new UTM2LatLon();
return
c.convertUTMToLatLong(UTM);
}
public
String latLon2UTM(double latitude, double longitude)
{
LatLon2UTM
c = new LatLon2UTM();
return
c.convertLatLonToUTM(latitude, longitude);
}
//..implementation
omitted
private
class LatLon2UTM
{
public
String convertLatLonToUTM(double latitude, double longitude)
{
//..implementation
omitted
}
//..implementation
omitted
}
private
class LatLon2MGRUTM extends LatLon2UTM
{
public
String convertLatLonToMGRUTM(double latitude, double longitude)
{
//..implementation
omitted
}
//..implementation
omitted
}
private
class MGRUTM2LatLon extends UTM2LatLon
{
public
double[] convertMGRUTMToLatLong(String mgrutm)
{
//..implementation
omitted
}
//..implementation
omitted
}
private
class UTM2LatLon
{
public
double[] convertUTMToLatLong(String UTM)
{
//..implementation
omitted
}
//..implementation
omitted
}
private
class Digraphs
{
//used
to get digraphs when doing conversion between
//lat/long
and MGRS
//..implementation
omitted
}
private
class LatZones
{
//include
methods to determine latitude zones
//..implementation
omitted
}
|
下一节将进一步探讨经纬度与 UTM 之间的转换。
将经纬度转换为 UTM
将经纬度坐标转换为 UTM 坐标需要使用 String latLon2UTM(double latitude, double longitude) 方法。该方法的实现
创建了一个新的内部类 LatLon2UTM c = new LatLon2UTM(); 实例，并将 UTM 坐标返回为由
15 个字符组成的字符串（即精度为 1 米）。LatLon2UTM 方法的实现如清单 2 所示：
清单 2. public String convertLatLonToUTM(double latitude, double longitude)
|
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
|
public
String convertLatLonToUTM(double latitude, double longitude)
{
validate(latitude,
longitude);
String
UTM = "";
setVariables(latitude,
longitude);
String
longZone = getLongZone(longitude);
LatZones
latZones = new LatZones();
String
latZone = latZones.getLatZone(latitude);
double
_easting = getEasting();
double
_northing = getNorthing(latitude);
UTM
= longZone + " " + latZone + " " + ((int) _easting) + " "+ ((int) _northing);
return
UTM;
}
|
该方法执行转换的方法为：调用各种方法获得经纬度区，然后计算以东和以北值，等等。使用 validate() 方法对输入进行验证；如果 (latitude
< -90.0 || latitude > 90.0 || longitude < -180.0 || longitude >= 180.0) 子句为真，将抛出一个 IllegalArgumentException。
清单 3 中的 setVariables() 方法设置计算转换所需的各种变量（请查看 “The Universal
Grids” 获取更多信息；可从 参考资料 获取链接）：
清单 3. protected void setVariables(double latitude, double longitude)
|
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
|
protected
void setVariables(double latitude, double longitude)
{
latitude
= degreeToRadian(latitude);
rho
= equatorialRadius * (1 - e * e) / POW(1 - POW(e * SIN(latitude), 2), 3 / 2.0);
nu
= equatorialRadius / POW(1 - POW(e * SIN(latitude), 2), (1 / 2.0));
double
var1;
if
(longitude < 0.0)
{
var1
= ((int) ((180 + longitude) / 6.0)) + 1;
}
else
{
var1
= ((int) (longitude / 6)) + 31;
}
double
var2 = (6 * var1) - 183;
double
var3 = longitude - var2;
p
= var3 * 3600 / 10000;
S
= A0 * latitude - B0 * SIN(2 * latitude) + C0 * SIN(4 * latitude) - D0
*
SIN(6 * latitude) + E0 * SIN(8 * latitude);
K1
= S * k0;
K2
= nu * SIN(latitude) * COS(latitude) * POW(sin1, 2) * k0 * (100000000) / 2;
K3
= ((POW(sin1, 4) * nu * SIN(latitude) * Math.pow(COS(latitude), 3)) / 24)
*
(5 - POW(TAN(latitude), 2) + 9 * e1sq * POW(COS(latitude), 2) + 4
*
POW(e1sq, 2) * POW(COS(latitude), 4))
*
k0
*
(10000000000000000L);
K4
= nu * COS(latitude) * sin1 * k0 * 10000;
K5
= POW(sin1 * COS(latitude), 3) * (nu / 6)
*
(1 - POW(TAN(latitude), 2) + e1sq * POW(COS(latitude), 2)) * k0
*
1000000000000L;
A6
= (POW(p * sin1, 6) * nu * SIN(latitude) * POW(COS(latitude), 5) / 720)
*
(61 - 58 * POW(TAN(latitude), 2) + POW(TAN(latitude), 4) + 270
*
e1sq * POW(COS(latitude), 2) - 330 * e1sq
*
POW(SIN(latitude), 2)) * k0 * (1E+24);
}
|
清单 4 中的 getLongZone() 方法和 LatZones 类（可从 源代码 获得）用来获得经纬度区。经度区通过
longitude 参数计算而来，而纬度区很难使用 LatZones 类中的数组进行编码。
清单 4. protected String getLongZone(double longitude)
|
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
|
protected
String getLongZone(double longitude)
{
double
longZone = 0;
if
(longitude < 0.0)
{
longZone
= ((180.0 + longitude) / 6) + 1;
}
else
{
longZone
= (longitude / 6) + 31;
}
String
val = String.valueOf((int) longZone);
if
(val.length() == 1)
{
val
= "0" + val;
}
return
val;
}
|
getNorthing() 方法（清单 5）和 getEasting() 方法（清单
6）计算以北和以东的值。两种方法都使用 清单 3 中的 setVariables() 方法设置的变量。
清单 5. protected double getNorthing(double latitude)
|
1
2
3
4
5
6
7
8
9
|
protected
double getNorthing(double latitude)
{
double
northing = K1 + K2 * p * p + K3 * POW(p, 4);
if
(latitude < 0.0)
{
northing
= 10000000 + northing;
}
return
northing;
}
|
清单 6. protected double getEasting()
|
1
2
3
4
|
protected
double getEasting()
{
return
500000 + (K4 * p + K5 * POW(p, 3));
}
|
清单 7 包含了一些示例输出，包括一些经纬度坐标和对应的 UTM 坐标：
清单 7. Latitude/longitude-to-UTM 测试值
|
1
2
3
4
5
6
7
8
9
10
11
|
(
0.0000 0.0000 ) "31 N 166021 0"
(
0.1300 -0.2324 ) "30 N 808084 14385"
(-45.6456
23.3545 ) "34 G 683473 4942631"
(-12.7650
-33.8765 ) "25 L 404859 8588690"
(-80.5434
-170.6540) "02 C 506346 1057742"
(
90.0000 177.0000) "60 Z 500000 9997964"
(-90.0000
-177.0000) "01 A 500000 2035"
(
90.0000 3.0000 ) "31 Z 500000 9997964"
(
23.4578 -135.4545) "08 Q 453580 2594272"
(
77.3450 156.9876) "57 X 450793 8586116"
(-89.3454
-48.9306 ) "22 A 502639 75072"
|
将 UTM 坐标转换为经纬度坐标
UTM 坐标到经纬度坐标的转换要比相反的转换过程容易一些。同样，“The Universal Grids”（请参阅 参考资料）提供了转换公式。清单
8 展示了 convertUTMToLatLong() 方法的代码。该方法返回一个双数组，其中第一个元素（数组索引 [0]）表示纬度，而第二个元素（数组索引 [1]）表示经度。由于
UTM 字符串参数的精度为 1 米，因此经纬度坐标具有与之相同的精度。
清单 8. public double[] convertUTMToLatLong(String UTM)
|
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
|
public
double[] convertUTMToLatLong(String UTM)
{
double[]
latlon = { 0.0, 0.0 };
String[]
utm = UTM.split(" ");
zone
= Integer.parseInt(utm[0]);
String
latZone = utm[1];
easting
= Double.parseDouble(utm[2]);
northing
= Double.parseDouble(utm[3]);
String
hemisphere = getHemisphere(latZone);
double
latitude = 0.0;
double
longitude = 0.0;
if
(hemisphere.equals("S"))
{
northing
= 10000000 - northing;
}
setVariables();
latitude
= 180 * (phi1 - fact1 * (fact2 + fact3 + fact4)) / Math.PI;
if
(zone > 0)
{
zoneCM
= 6 * zone - 183.0;
}
else
{
zoneCM
= 3.0;
}
longitude
= zoneCM - _a3;
if
(hemisphere.equals("S"))
{
latitude
= -latitude;
}
latlon[0]
= latitude;
latlon[1]
= longitude;
return
latlon;
}
|
convertUTMToLatLong() 方法将传入的 UTM 字符串（格式为 34
G 683473 4942631）分解，并使用 getHemisphere() 方法确定字符串表示的位置所在的半球。这种确定非常简单：纬度区 A、C、D、E、F、G、H、J、K、L 和 M 位于南半球，而其余区位于北半球。
清单 9 所示的 setVariables() 方法将设置计算所需的变量，然后立即计算纬度值。经度值则通过经度区计算。
清单 9. protected void setVariables()
|
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
|
protected
void setVariables()
{
arc
= northing / k0;
mu
= arc
/
(a * (1 - POW(e, 2) / 4.0 - 3 * POW(e, 4) / 64.0 - 5 * POW(e, 6) / 256.0));
ei
= (1 - POW((1 - e * e), (1 / 2.0)))
/
(1 + POW((1 - e * e), (1 / 2.0)));
ca
= 3 * ei / 2 - 27 * POW(ei, 3) / 32.0;
cb
= 21 * POW(ei, 2) / 16 - 55 * POW(ei, 4) / 32;
cc
= 151 * POW(ei, 3) / 96;
cd
= 1097 * POW(ei, 4) / 512;
phi1
= mu + ca * SIN(2 * mu) + cb * SIN(4 * mu) + cc * SIN(6 * mu) + cd
*
SIN(8 * mu);
n0
= a / POW((1 - POW((e * SIN(phi1)), 2)), (1 / 2.0));
r0
= a * (1 - e * e) / POW((1 - POW((e * SIN(phi1)), 2)), (3 / 2.0));
fact1
= n0 * TAN(phi1) / r0;
_a1
= 500000 - easting;
dd0
= _a1 / (n0 * k0);
fact2
= dd0 * dd0 / 2;
t0
= POW(TAN(phi1), 2);
Q0
= e1sq * POW(COS(phi1), 2);
fact3
= (5 + 3 * t0 + 10 * Q0 - 4 * Q0 * Q0 - 9 * e1sq) * POW(dd0, 4) / 24;
fact4
= (61 + 90 * t0 + 298 * Q0 + 45 * t0 * t0 - 252 * e1sq - 3 * Q0
*
Q0)
*
POW(dd0, 6) / 720;
lof1
= _a1 / (n0 * k0);
lof2
= (1 + 2 * t0 + Q0) * POW(dd0, 3) / 6.0;
lof3
= (5 - 2 * Q0 + 28 * t0 - 3 * POW(Q0, 2) + 8 * e1sq + 24 * POW(t0, 2))
*
POW(dd0, 5) / 120;
_a2
= (lof1 - lof2 + lof3) / COS(phi1);
_a3
= _a2 * 180 / Math.PI;
}
|
setVariables() 使用以东和以北值设置所需的变量。这些都是类变量并且在 convertUTMToLatLong(String
UTM) 方法中进行设置（参见 清单 8）。
