一. 投影：从不规则梨形球体到平面
我们日常所用的地图是平面地图。但地球实际上是一个赤道略宽两极略扁的不规则的梨形球体，它的表面是个无法直接展平的不规则曲面。将它转换成平面地图要经过一系列投影换算。用下面这张图来直观感受下真实地球和印象里的地球的区别：
真实的地球与印象里的地球，图片来自网络，来源见水印
下面这张动图也很直观：
真实地球动图
1.1 第一重投影：椭球体的规则化
准确地说，地球是个不规则的三轴椭球体。为了将地图展开到平面上，我们首先将真实地球投影成一个规则的椭球。例如著名的 WGS84 坐标体系（1984年世界大地坐标系统），就将地球视为一个规则椭球，其与地图展开相关的参数为：
长半轴 a = 6378137 米
椭球扁率 f = 1/298.257223565
所以对于所有以 WGS84 为基础的平面地图，关于距离、面积等等的计算实际上是真实地球在 WGS84 所定义的椭球面上的计算，这一点会在下一篇博文中详述。
WGS84 参数，图片来自网络
1.2 第二重投影：椭球面的墨卡托投影
为了将地球平面化，历史上人们采用了多种投影方式，其中当今世界上各WEB地图所普遍采用的是墨卡托投影。
想象一个外切于椭球赤道的圆柱面，在球心处放置一个光源，放出的射线将球面各点在圆柱内壁投下影子。沿某条经线剪开此圆柱便可得到平面地图，如下图所示：
墨卡托投影示意图
不难分析，此投影存在以下特征：
- 纵轴方向，纬度的变化是非线性的，越靠近两极，在圆柱面上的投影就越趋向于无穷远，而两极无法被投影到圆柱面上
- 横轴方向，经度的变化是线性的
- 投影区域和都比真实区域大， 处于拉伸状态
所以我们日常所见的 WEB 地图在南北的最大纬度并非 ±90 度，而基本是在 ±85度 左右。
墨卡托投影具有等角性质，即球体上的两点之间的角度方位与平面上的两点之间的角度方位保持不变，因此特别适合用于导航。
墨卡托投影的球面坐标系经纬度与平面坐标系坐标之间的转换公式如下：
墨卡托投影公式
其中x、y是投影展开成平面后以赤道本初子午线交点为原点的平面坐标系的坐标，a 是椭球体长半轴，b是短半轴。L是经度（弧度制），B是纬度（弧度制）。
很容易就能找到开源地图引擎中关于墨卡托投影转换的代码，比如leaflet：
L.Projection.Mercator = {
R: 6378137,
R_MINOR: 6356752.314245179,
bounds: L.bounds([-20037508.34279, -15496570.73972], [20037508.34279, 18764656.23138]),
project: function (latlng) {
var d = Math.PI / 180,
r = this.R,
y = latlng.lat * d,
tmp = this.R_MINOR / r,
e = Math.sqrt(1 - tmp * tmp),
con = e * Math.sin(y);
var ts = Math.tan(Math.PI / 4 - y / 2) / Math.pow((1 - con) / (1 + con), e / 2);
y = -r * Math.log(Math.max(ts, 1E-10));
return new L.Point(latlng.lng * d * r, y);
},
unproject: function (point) {
var d = 180 / Math.PI,
r = this.R,
tmp = this.R_MINOR / r,
e = Math.sqrt(1 - tmp * tmp),
ts = Math.exp(-point.y / r),
phi = Math.PI / 2 - 2 * Math.atan(ts);
for (var i = 0, dphi = 0.1, con; i < 15 && Math.abs(dphi) > 1e-7; i++) {
con = e * Math.sin(phi);
con = Math.pow((1 - con) / (1 + con), e / 2);
dphi = Math.PI / 2 - 2 * Math.atan(ts * con) - phi;
phi += dphi;
}
return new L.LatLng(phi * d, point.x * d / r);
}
};
而在实际应用中为了便于计算，部分厂商采用的是墨卡托投影的变种——Web墨卡托投影，又叫球面墨卡托投影。它接收的输入是 WGS84 的经纬度，但在投影时不再把地球当做椭球而当做半径为6378137米的标准球体，以简化计算。简化过的 Web墨卡托投影，球面坐标系经纬度与平面坐标系坐标之间的转换公式如下：
Web墨卡托投影
其中 R 为标准球体的半径。
公式看起来很累，但同样很容易就能找到已实现的代码：
/*
* @namespace Projection
* @projection L.Projection.SphericalMercator
*
* Spherical Mercator projection — the most common projection for online maps,
* used by almost all free and commercial tile providers. Assumes that Earth is
* a sphere. Used by the `EPSG:3857` CRS.
*/
L.Projection.SphericalMercator = {
R: 6378137,
MAX_LATITUDE: 85.0511287798,
project: function (latlng) {
var d = Math.PI / 180,
max = this.MAX_LATITUDE,
lat = Math.max(Math.min(max, latlng.lat), -max),
sin = Math.sin(lat * d);
return new L.Point(
this.R * latlng.lng * d,
this.R * Math.log((1 + sin) / (1 - sin)) / 2);
},
unproject: function (point) {
var d = 180 / Math.PI;
return new L.LatLng(
(2 * Math.atan(Math.exp(point.y / this.R)) - (Math.PI / 2)) * d,
point.x * d / this.R);
},
bounds: (function () {
var d = 6378137 * Math.PI;
return L.bounds([-d, -d], [d, d]);
})()
};
Web墨卡托投影的计算比椭球面墨卡托投影要简化些，当然这是以一定的精度丢失为代价的。且可以发现代码中已经将纬度范围定死在-85度到+85度之间。我们院里封装的地图API在计算距离和坐标转换时采用的是椭球面墨卡托投影。
上述关于椭球面墨卡托投影和球面墨卡托投影的公式和代码给出了 WGS84 坐标系经纬度与展开后的平面坐标系的相互转换方法。这是第二重投影，也是第一重坐标转换。下面将介绍第二、第三重坐标转换。
二. 坐标转换：从墨卡托投影到像素坐标
上述已实现了 WGS84 经纬度坐标系与平面坐标系之间的换算。而当我们在网页上加载地图时，可以发现每一个像素都对应着一个坐标。这是如何换算的呢？
这就要牵扯出地图学中另一个重要的概念——比例尺（即地图上的一厘米代表着实际上的多少厘米）。到了web地图中我们把比例尺转换成另一个概念——分辨率（Resolution，图上一像素代表实际多少米）。比例尺跟分辨率的换算举个例子：
假设地图的坐标单位是米，整张地图的dpi为96，当前地图在赤道处的比例尺为1:125000000（即图上1米等于实地125000000米）。1英寸=2.54厘米； 1英寸=96像素。
那么计算可得地图赤道上1像素代表实地距离是 125000000*0.0254/96 = 33072.9166666667米。
也就是说，平面地图上一个像素实际上代表的是一个坐标范围，我们取其范围的中心点坐标作为此像素所代表的坐标值。
为什么要强调某一条纬度线（上述例子为赤道）的比例尺？因为根据墨卡托投影的特性，同一张地图中不同纬度线的比例尺是变化的，越靠近两极，图上1米相当于实地的距离越小。
2.1 第二重坐标转换：墨卡托投影到世界平面点
接下来再引入另一个概念：世界平面点。世界平面点本意是指墨卡托投影展开后的平面地图上的一个点，显然这个点应该有自己的坐标。而在上一章中我们已经求得投影坐标 (x,y)，为了计算（x,y）实际所在的像素点，引入了世界平面点这么个中间概念。
常见的设备或浏览器窗口中，像素坐标系是以左上角为原点，以向下向右为正方向的。而上章中的经纬度以及换算得到的平面坐标系坐标（x,y），是以赤道和本初子午线的交点，即地图的中心点为原点的。为了方便计算和展示，世界平面地图同样以左上角为坐标原点，向下向右为正方向，那么就需要进一步换算。
在 leaflet 中，这个坐标系的原点在左上角（0, 0），范围为0~1。对应代码：
/*
* @class Transformation
* @aka L.Transformation
*
* Represents an affine transformation: a set of coefficients `a`, `b`, `c`, `d`
* for transforming a point of a form `(x, y)` into `(a*x + b, c*y + d)` and doing
* the reverse. Used by Leaflet in its projections code.
*
* @example
*
* ```js
* var transformation = new L.Transformation(2, 5, -1, 10),
* p = L.point(1, 2),
* p2 = transformation.transform(p), // L.point(7, 8)
* p3 = transformation.untransform(p2); // L.point(1, 2)
* ```
*/
// factory new L.Transformation(a: Number, b: Number, c: Number, d: Number)
// Creates a `Transformation` object with the given coefficients.
L.Transformation = function (a, b, c, d) {
this._a = a;
this._b = b;
this._c = c;
this._d = d;
};
L.Transformation.prototype = {
// @method transform(point: Point, scale?: Number): Point
// Returns a transformed point, optionally multiplied by the given scale.
// Only accepts actual `L.Point` instances, not arrays.
transform: function (point, scale) { // (Point, Number) -> Point
return this._transform(point.clone(), scale);
},
// destructive transform (faster)
_transform: function (point, scale) {
scale = scale || 1;
point.x = scale * (this._a * point.x + this._b);
point.y = scale * (this._c * point.y + this._d);
return point;
},
// @method untransform(point: Point, scale?: Number): Point
// Returns the reverse transformation of the given point, optionally divided
// by the given scale. Only accepts actual `L.Point` instances, not arrays.
untransform: function (point, scale) {
scale = scale || 1;
return new L.Point(
(point.x / scale - this._b) / this._a,
(point.y / scale - this._d) / this._c);
}
};
对于地球的web墨卡托而言，上述代码中a、b、c、d四个参数为：
/*
* @namespace CRS
* @crs L.CRS.EPSG3857
*
* The most common CRS for online maps, used by almost all free and commercial
* tile providers. Uses Spherical Mercator projection. Set in by default in
* Map's `crs` option.
*/
L.CRS.EPSG3857 = L.extend({}, L.CRS.Earth, {
code: 'EPSG:3857',
projection: L.Projection.SphericalMercator,
transformation: (function () {
var scale = 0.5 / (Math.PI * L.Projection.SphericalMercator.R);
return new L.Transformation(scale, 0.5, -scale, 0.5);
}())
});
其中，scale代表球的周长分之一，b和d都是0.5这代表赤道和本初子午线的交点在世界平面点的位置为(0.5, 0.5)；this._a * point.x + this._b代表x轴方向墨卡托坐标在世界平面点的位置，c=-scale，结合 this._c * point.y + this._d，能计算出y轴方向墨卡托在世界平面点位置。至于c为什么是负的，结合一下经纬度或投影坐标中纬度的性质，以上为正下为负。而到了世界平面坐标中，负的纬度坐标对应的平面点坐标要大于0.5。
2.2 第三重坐标转换：世界平面点到像素坐标
接下来就应该很方便了：只要知道整张平面地图横纵两个方向各自的总像素值，就能很快知道某个点所在的像素点：
设地图横坐标上总像素数为W，纵坐标上总像素数为H。设地图上某点世界平面点坐标为（x,y）,则求其所在像素点坐标：
像素点横坐标：w = Math.floor(x*W);
像素点纵坐标：h = Math.floor(y*H);
那么可能要问，这么简单的计算为什么要分两步？这涉及到地图分辨率分级，或者说地图缩放的问题。一般都会将地图分辨率划分为几个级别，比如百度地图，最高是20级，最低是1级。而在leaflet中，每一级对应的分辨率（赤道上）是：
resolution = 2*PI*R/(256*Math.pow(2,zoom))
对应的整张地图的图片分辨率为：
256*Math.pow(2,zoom) * 256*Math.pow(2,zoom)
其中R是Web墨卡托地球半径，zoom是地图级别。当zoom=0时，整张地图表现为一张256×256的图片。
为什么是256？因为便于计算且大小合适。
现在我们已经实现了从经纬度到像素点坐标的转换过程：
- 经纬度到墨卡托投影坐标
- 墨卡托投影坐标到世界平面点坐标
- 世界平面点坐标到世界地图像素点坐标
当然也同时实现了逆转换代码。
leaflet中对这个过程做了进一步封装：
// @method latLngToPoint(latlng: LatLng, zoom: Number): Point
// Projects geographical coordinates into pixel coordinates for a given zoom.
latLngToPoint: function (latlng, zoom) {
var projectedPoint = this.projection.project(latlng),
scale = this.scale(zoom);
return this.transformation._transform(projectedPoint, scale);
},
// @method scale(zoom: Number): Number
// Returns the scale used when transforming projected coordinates into
// pixel coordinates for a particular zoom. For example, it returns
// `256 * 2^zoom` for Mercator-based CRS.
scale: function (zoom) {
return 256 * Math.pow(2, zoom);
},
接下来的问题是，如何在页面中加载地图？
三. 瓦片：地图的分割与局部加载
我们不难发现，在网页上地图是分割成一个个小方块来加载的，那么为什么要这么做，以及是如何计算现在要加载哪一个图块呢？
3.1 瓦片的分割与定位计算
道理很简单——地图太大了。到达一定级别后，一次性加载整张地图不管是从设备的内存还是网络传输来看都是不可能的。而且用户通常不需要查看整张地图，而是关注一定分辨率下的某一部分。所以局部加载方案应运而生：只加载用户当前关注的那部分地图。
所以大多数厂商把地图等分为一块块 256×256 的图块，称之为瓦片，如下图所示：
地图瓦片
在加载过程中，获取用户当前关注的区域的所有瓦片并将其无缝拼接成地图。每一级瓦片的数量是：
Math.pow(2,zoom*2)
可以计算下第18级地图所拥有的瓦片数量：
Math.pow(2,18*2) = 68719476736
既然已经分割好每一个级别的瓦片，那么地图的缩放就是水到渠成的事了。下面这张图形象地展示了地图如何缩放：
地图瓦片与缩放
地图上每张瓦片有其特定编码，主流厂商大都用x,y,z表示，其中 x 为其横轴方向编码，y为其纵轴方向编码，z为其地图缩放级别。原点为左上角。
现在我们还是以leaflet为例，来计算下已知某点经纬度为（lng，lat），求zoom级别下其所在瓦片：
var lnglat = new L.LatLng(lng,lat);
var point = L.latLngToPoint(lnglat,zoom);
var tileSize = 256;
var xIndex = point.x / tileSize,
yIndex = point.y / tileSize;
其瓦片编码即为 ( xIndex, yIndex, zoom) 。
3.2 瓦片的加载
一般来说在实例化一个地图时，都会给给Map构造函数传入一个zoom和一个center参数。以leaflet为例，已知中心点坐标，要加载一幅地图，我们只需要知道屏幕四个点的经纬度所在范围内的瓦片，再将这些瓦片按照一定的偏移坐标布置即可。
屏幕范围与瓦片加载
上面传入的center代表当前范围的中心点，同时也是屏幕的中心点，那么就可以求出该经纬度对应的像素坐标，这个像素坐标就是屏幕中心点对应的瓦片像素坐标。这里的像素与我们的css像素一一对应，利用屏幕范围可得出屏幕四个角点的瓦片像素坐标。利用这四个点的瓦片坐标，可以求出当前屏幕的瓦片索引范围，加载这些瓦片。代码如下：
_getTiledPixelBounds: function (center) {
var map = this._map,
mapZoom = map._animatingZoom ? Math.max(map._animateToZoom, map.getZoom()) : map.getZoom(),
scale = map.getZoomScale(mapZoom, this._tileZoom),
pixelCenter = map.project(center, this._tileZoom).floor(),
halfSize = map.getSize().divideBy(scale * 2);
return new L.Bounds(pixelCenter.subtract(halfSize), pixelCenter.add(halfSize));
},
_pxBoundsToTileRange: function (bounds) {
var tileSize = this.getTileSize();
return new L.Bounds(
bounds.min.unscaleBy(tileSize).floor(),
bounds.max.unscaleBy(tileSize).ceil().subtract([1, 1]));
},
接下来要注意一些，这时候这些瓦片的坐标范围肯定是大于屏幕的坐标范围，所以要对所有的瓦片做一些偏移。计算过程比较简单，屏幕坐标以左上点为原点，这个点对应的像素坐标已知，只要求出每个瓦片的左上角点的瓦片像素坐标与屏幕左上点的瓦片像素坐标做差值，即可得出在css中的position的偏移值（高级点的用css3的translate）。下面我们可以看看leaflet的处理过程：
_setView: function (center, zoom, noPrune, noUpdate) {
var tileZoom = Math.round(zoom);
if ((this.options.maxZoom !== undefined && tileZoom > this.options.maxZoom) ||
(this.options.minZoom !== undefined && tileZoom < this.options.minZoom)) {
tileZoom = undefined;
}
var tileZoomChanged = this.options.updateWhenZooming && (tileZoom !== this._tileZoom);
if (!noUpdate || tileZoomChanged) {
this._tileZoom = tileZoom;
if (this._abortLoading) { // 如果zoom要发生变化，停止当前所有tiles的加载；通过更改他们的onload onerror事件实现
this._abortLoading();
}
// 1、创建该级别容器瓦片
// 2、 设置zIndex
// 3、设置本级别图层基准点origin
// 4、设置值本级别容器的偏移量
this._updateLevels();
// 1、得到世界的像素bounds
// 2、得通过像素范围除以tileSize得到能够覆盖世界的瓦片范围
// 3、得到坐标系经度和纬度范围内的瓦片范围
this._resetGrid();
if (tileZoom !== undefined) {
// 加载可视范围内的瓦片
// 1、计算可视区域的像素范围
// 2、 将像素范围转变成瓦片格网范围
// 3、计算一个buffer的格网范围
// 4、将不再当前范围内已加载的瓦片打上标签
// 5、如果zoom发生变化重新进行setView
// 6、将格网范围内的tile放入一个数组中
// 7、对数组进行排序，靠近中心点的先加载
// 8、创建瓦片
// (1) 计算瓦片在地图上的偏移量 coords * tileSize - origin
// (2) 加载瓦片数据（图片或者矢量数据）
// (3) 设置图片位置 setPosition
this._update(center);
}
if (!noPrune) {
this._pruneTiles(); // 移除不在范围内的tile; retainParent部分尚没看懂，可能是按照瓦片金字塔保留
}
// Flag to prevent _updateOpacity from pruning tiles during
// a zoom anim or a pinch gesture
this._noPrune = !!noPrune;
}
//将地图的新中心点移到地图中央
this._setZoomTransforms(center, zoom);
},
以上。
网友评论