Cesium中的地形和坐标转换说明

1 Cesium中的地形

Cesium中的地形系统是一种由流式瓦片数据生成地形mesh的技术，厉害指出在于其可以自动模拟出地面、海洋的三维效果。创建地形图层的方式如下：

var terrainProvider = new Cesium.CesiumTerrainProvider({
    url : 'https://assets.agi.com/stk-terrain/v1/tilesets/world/tiles', // 默认立体地表
    // 请求照明
    requestVertexNormals: true,
    // 请求水波纹效果
    requestWaterMask: true
});
viewer.terrainProvider = terrainProvider;

Cesium支持两种类型的地形，STK World Terrain和Small Terrain。

1.1 STK世界地形

STK世界地形（STK World Terrain），其是高分辨率, 基于quantized mesh的地形。这是一种基于网格的地形，可充分利用GL中的Shader来渲染，效果相当逼真。该地形使用了多种数据源，分别适应不同地区和不同精度时的情形，如，美国本土使用美国国家高程数据集（National Elevation Dataset，NED）的高程，精度3-30米；对于欧洲使用EU-DEM高程，精度30米；对于澳洲使用Australia SRTM-derived 1 Second DEM高程，精度30米；对于-60至60纬度段使用CGIAR SRTM高程，精度90米；对于整个地球使用GTOPO30，精度1000米。该地形的生成方式尚未公开，对于封闭的局域网应用时，则需购买AGI的STK terrain server。但是AGI提供了一个webapi可供因特网上调用，并提供了这种地形的格式细节。

1.2 Small Terrain

Small Terrain是中等高分辨率基于heightmap的地形，渲染出的地形效果不如quantized mesh的地形，但也基本能接受。
可以由DEM数据生成这种规范的.terrain文件。生成工具见https://groups.google.com/forum/#!topic/cesium-dev/rBieaEBJHi，需要gdal库和numpy。

2 坐标系

2.1 cesium中的坐标系

Cesium中常用的坐标系主要有两种：

WGS84坐标系和笛卡尔空间直角坐标系。平时我们常见的某个点的经纬度就是在WGS84坐标系下某个点的坐标，它的坐标原点在椭球的质心；

而笛卡尔坐标系主要是用来做空间位置的变化如平移、旋转和缩放等等，它的坐标原点在椭球的中心。

二者的联系如下图, 笛卡尔空间坐标的原点就是椭球的中心。

Pick----屏幕坐标

Cartesian----世界坐标

cartographic-----地理坐标（弧度）

Cesium中的坐标系：

1、平面坐标系（Cartesian2）；

2、笛卡尔空间直角坐标系（Cartesian3）；

3、Cartesian4（unknown，在应用中几乎用不到）

4、Cartographic（地理坐标系下经纬度的弧度表示），通常情况下通过它和WGS84坐标系之间互转。

2.1.1世界坐标

以椭球中心为原点的空间直角坐标系中的一个点的坐标。Cesium中用Cartesian3变量表示，笛卡尔空间直角坐标系，用new Cesium.Cartesian3(x, y, z)创建。

2.1.2 地理坐标

就是测绘中的地理经纬度坐标，地理坐标系，坐标原点在椭球的质心。

经度：参考椭球面上某点的大地子午面与本初子午面间的两面角。东正西负。

纬度：参考椭球面上某点的法线与赤道平面的夹角。北正南负。

Cesuim中没有具体的经纬度对象，要得到经纬度首先需要计算为弧度，再进行转换。

2.1.3 弧度

Cartographic变量表示。

new Cesium.Cartographic(longitude, latitude, height)。是用弧度表示的经纬度，这里的参数也叫做longitude，latitude，即经度和纬度。弧度即角度对应弧长是半径的倍数。

角度转弧度 π/180×角度；
弧度变角度 180/π×弧度。

2.2 坐标转换

Cesium其实是一个封装好的WebGL库，当然这里面就牵扯到好几套坐标问题：屏幕坐标、三维空间坐标、投影坐标。坐标转换肯定是我们在开发任何地理信息系统中经常会碰到的问题，也比较复杂。

“平面坐标系” 和“笛卡尔空间直角坐标系”和“Cartographic”之间的相互转换思路如下所示。

坐标的定义：

1：Cartesian2-----new Cesium.Cartesian2(x, y)

2：Cartesian3---- new Cesium.Cartesian3(x, y, z）

3：Cartographic----new Cesium.Cartographic(longitude, latitude, height) 注：经纬度为弧度单位

转换：

Cartesian3→ Cartesian2，　　　　　　　　Cesium.Cartesian2.fromCartesian3(cartesian, result)→ Cartesian2

经纬度坐标（WGS84）→ Cartesian3，　　 Cesium.Cartesian3.fromDegrees(longitude, latitude, height, ellipsoid, result) → Cartesian3

弧度坐标 → Cartesian3，　　　　　　　　　Cesium.Cartesian3.fromRadians(longitude, latitude, height, ellipsoid, result) → Cartesian3

Cartesian3→ Cartographic，　　　　　　　 Cesium.Cartographic.fromCartesian(cartesian, ellipsoid, result) → Cartographic

经纬度坐标（WGS84）→ Cartographic，　 Cesium.Cartographic.fromDegrees(longitude, latitude, height, result) → Cartographic

单位：
经纬度坐标和弧度坐标也可以通过Cesium.Math来转换。

Cesium.CesiumMath.toDegrees(radians) → Number

Cesium.CesiumMath.toRadians(degrees) → Number

2.2.1 坐标系

new Cesium.Cartesian2(1,1) //表示一个二维笛卡尔坐标系，也就是直角坐标系（屏幕坐标系）
new Cesium.Cartesian3(1,1,1) //表示一个三维笛卡尔坐标系，也是直角坐标系(就是真实世界的坐标系)

2.2.2 二维屏幕坐标系到三维坐标系的转换

var pick1= scene.globe.pick(viewer.camera.getPickRay(pt1), scene) //其中pt1为一个二维屏幕坐标。

2.2.3 三维坐标到地理坐标的转换

var geoPt1= scene.globe.ellipsoid.cartesianToCartographic(pick1) //其中pick1是一个Cesium.Cartesian3对象。

2.2.4 地理坐标到经纬度坐标的转换

var point1=[geoPt1.longitude / Math.PI * 180,geoPt1.latitude / Math.PI * 180]; //其中geoPt1是一个地理坐标。

2.2.5 经纬度坐标转地理坐标（弧度）

var cartographic = Cesium.Cartographic.fromDegree(point) //point是经纬度值
var coord_wgs84 = Cesium.Cartographic.fromDegrees(lng, lat, alt);//单位：度，度，米

2.2.6 经纬度坐标转世界坐标

var cartesian = Cesium.Cartesian3.fromDegree(point)

2.2.7 计算两个三维坐标系之间的距离

var d = Cesium.Cartesian3.distance(
    new Cesium.Cartesian3(pick1.x, pick1.y, pick1.z), 
    new Cesium.Cartesian3(pick3.x, pick3.y, pick3.z)
); //pick1、pick3都是三维坐标系

转换到笛卡尔坐标系后就能运用计算机图形学中的仿射变换知识进行空间位置变换如平移旋转缩放。cesium.js为我们提供了很有用的变换工具类，
Cesium.Cartesian3（相当于Point3D）Cesium.Matrix3（3x3矩阵，用于描述旋转变换）Cesium.Matrix4（4x4矩阵，用于描述旋转
加平移变换），Cesium.Quaternion（四元数，用于描述围绕某个向量旋转一定角度的变换）。下面举个例子：

一个局部坐标为p1(x,y,z)的点，将它的局部坐标原点放置到loc(lng,lat,alt)上，局部坐标的z轴垂直于地表，局部坐标的y轴指向正北，
并围绕这个z轴旋转angle度，求此时p1(x,y,z)变换成全局坐标笛卡尔坐标p2(x1,y1,z1)是多少？

var rotate = Cesium.Math.toRadians(angle);//转成弧度
var quat = Cesium.Quaternion.fromAxisAngle(Cesium.Cartesian3.UNIT_Z, rotate); //quat为围绕这个z轴旋转d度的四元数
var rot_mat3 = Cesium.Matrix3.fromQuaternion(quat);//rot_mat3为根据四元数求得的旋转矩阵

var pt = new Cesium.Cartesian3(x, y, z);//p1的局部坐标
// m2为旋转加平移的4x4变换矩阵，这里平移为(0,0,0)，故填个Cesium.Cartesian3.ZERO
var m = Cesium.Matrix4.fromRotationTranslation(rot_mat3, Cesium.Cartesian3.ZERO);
m = Cesium.Matrix4.multiplyByTranslation(m, pt);//m = m X v

//得到局部坐标原点的全局坐标
var cart3 = ellipsoid.cartographicToCartesian(Cesium.Cartographic.fromDegrees(lng, lat, alt));
//m1为局部坐标的z轴垂直于地表，局部坐标的y轴指向正北的4x4变换矩阵
var m1 = Cesium.Transforms.eastNorthUpToFixedFrame(cart3);
m = Cesium.Matrix4.multiplyTransformation(m, m1);//m = m X m1
var p2 = Cesium.Matrix4.getTranslation(m);//根据最终变换矩阵m得到p2