📄 tfw坐标定义.txt
字号:
18 <authority> = AUTHORITY["", ""]
19 <vert cs> = VERT_CS["", <vert datum>, <linear unit>, {<axis>,} {,<authority>}]
20 <vert datum> = VERT_DATUM["", <datum type> {,<authority>}]
21 <datum type> = <number>
22 <compd cs> = COMPD_CS["", <head cs>, <tail cs> {,<authority>}]
23 <head cs> = <coordinate system>
24 <tail cs> = <coordinate system>
25 <twin axes> = <axis>, <axis>
26 <axis> = AXIS["", NORTH | SOUTH | EAST | WEST | UP | DOWN | OTHER]
27 <to wgs84s> = TOWGS84[<seven param>]
28 <seven param> = <dx>, <dy>, <dz>, <ex>, <ey>, <ez>, <ppm>
29 <dx> = <number>
30 <dy> = <number>
31 <dz> = <number>
32 <ex> = <number>
33 <ey> = <number>
34 <ez> = <number>
35 <ppm> = <number>
36 <fitted cs> = FITTED_CS["", <to base>, <base cs>]
37 <to base> = <math transform>
38 <base cs> = <coordinate system>
39 <local cs> = LOCAL_CS["", <local datum>, <unit>, <axis>, {,<axis>}* {,<authority>}]
40 <local datum> = LOCAL_DATUM["", <datum type> {,<authority>}]
这个结构的阅读方法是这样的，所有的尖括号<>中的值都是一个需要被替代的东西（我觉得应该叫做一个类型）。你可以在表的其余中找到对应的类型描述，但是这个类型是不能直接写出来的，也就是说相当于一个链接，也可以看成是一个指针或一个引用，它指向一个对象，但是却不是一个对象，你需要用一个确实的东西来代替它。而且文档的下面有各种类型的解释。在花括号{}中的是可选的项目，表现在实际的坐标表示中，就是可写可不写。等号=是用来表示类型的表达形式。 用 | 分割开的是并列的选项，也就是在实际坐标表示中你只能选用其中的一个。*号表示有数个，多个相同的类型组成。除了上面的其他的符号就代表确实的东西（[]是需要实际写入坐标系表示字符串的一部分，不是这个语法表达表中的表达符号，冒号" 逗号, 同样也是实体的一部分），是需要到坐标系中的。最后注意书写顺序从一条条对应下来，就可以了。比如从开始，它不是一个东西，而它是等类型的其中一个类型。如果我们判断我们需要的坐标系是一个类型，我们就去中找，是等类型其中一个类型。如果我们判断我们需要的坐标系是一个投影坐标类型，我们就到中找，的内容不再是虚的了，我们就把PROJ[.., .., ..]结构写下来。PROJ["", , , {,}* {,}{,}]，其中包括了类型，类型，类型，数个可选的类型，类型，还有可选的和类型组成的。然后到每个类型中去找对应的结构，直到全部内容都写完为止。需要解释的是 是对应的实际的类型的名称和编码。不是在表中有对应。就如上面的例子，基准面名是"North_American_Datum_1927"，对应的类型的类型就用该名称代替。也就是说各种坐标、投影、单位、编码，数值也都算做一个单独的类型。最后直到把那些类型都用实际的东西填满为止。
解释起来很麻烦(汗都说出来了)，如果你还是不理解，算了，就当我没说，其实你可以完全不管这些的（但是如果你要写，可以看看这里，这里来获取更多有关于坐标系统和其参数的信息，你还可以下载postgresql，如果安装了postgis插件，你就拥有了几乎所有常用的坐标系统的wkt表示）。如果我们不亲自写不同坐标间的转换函数。我们只要关心两个坐标系统是否完全一样，其他都可以不去管它。要是说起不同坐标系统之间转换，又有一大堆的东西要说了，还是以后开个专门的页面写写吧。对于我们来说重要的是另外一个部分－－坐标转换。
3.2. 坐标转换
Toggle line numbers Toggle line numbers
1 >>> dataset.GetGCPs()
2 []
3 >>> dataset.GetGeoTransform()
4 (590000.0, 20.0, 0.0, 4928000.0, 0.0, -20.0)
5 >>>
这是我们主要要对付的。只有通过这几个数据，我们才可以知道自己手中这张图该如何铺展到绘制平面上。上面说过了，栅格坐标到地理坐标之间的转换可以通过两种方式来进行：一种是仿射坐标转换，一种是GCPs。一般来说，仿射坐标转换比较常见。仿射坐标转换就是通过几个值来进行栅格框架到地理框架的映射。上面的例子中出现了6个值。这六个值可以用下面的
Toggle line numbers Toggle line numbers
1 Xgeo = GT(0) + Xpixel*GT(1) + Yline*GT(2)
2 Ygeo = GT(3) + Xpixel*GT(4) + Yline*GT(5)
公式来解释，而且数组顺序对应起来了，理解起来也不困难。至于几个参数究竟是什么意思，这里讲得更清楚：
Toggle line numbers Toggle line numbers
1 adfGeoTransform[0] /* top left x */
2 adfGeoTransform[1] /* w-e pixel resolution */
3 adfGeoTransform[2] /* rotation, 0 if image is "north up" */
4 adfGeoTransform[3] /* top left y */
5 adfGeoTransform[4] /* rotation, 0 if image is "north up" */
6 adfGeoTransform[5] /* n-s pixel resolution */
看起来和ESRI为GeoTiff定义的tfw文件格式很像。至于是否真的是完全对应，我还不敢说。
有了上面的公式，我们就可以把栅格的每个点代入公式来求出其在地球上的实际位置（当然，求出的数值的单位是在坐标系统中定义的那个，而非经纬度坐标，要求经纬度坐标，还要通过一步转换）。当然我们可以不用这么麻烦，只要求出整个图像外框上下左右四个边界就可以定下整个图像的位置，因为图像都是矩形的（那些非矩形的图像其实是有另外一层二值的掩膜图像（mask）叠加而成的，这种图像只有两种值，一种是0,一种是1，1代表按照图像原色显示，0代表不显示）。图像中某点的地理位置可以通过外框四点线性计算。
现在看看我们这张图的边界地理范围：
Toggle line numbers Toggle line numbers
1 >>> dataset.GetGeoTransform()
2 (590000.0, 20.0, 0.0, 4928000.0, 0.0, -20.0)
3 >>> dataset.RasterXSize
4 950
5 >>> dataset.RasterYSize
6 700
7 >>>
上边框和左边框都不用算了，是4928000和590000。 求右下角的坐标，代入公式：
右边框 = Xgeo = GT(0) + Xpixel*GT(1) + Yline*GT(2) = 590000.0 + 950*20.0 + 700*0.0 = 609000.0
下边框 = Ygeo = GT(3) + Xpixel*GT(4) + Yline*GT(5) = 4928000.0 + 950*0.0 + 700*(-20.0) = 4914000.0
开ArcGIS验证一下，没错。
Toggle line numbers Toggle line numbers
1 In projected or local coordinates
2 Left: 590000.000000
3 Right: 609000.000000
4 Top: 4928000.000000
5 Bottom: 4914000.000000
好了，有了这些，就可以把图像和地理坐标联系起来了。有了地理坐标，就可以把图像在照屏幕上预先定好的坐标系画出来了。
需要注意的一点是，对于GeoTiff数据源来说，坐标系统和仿射变换参数有两种存储方式，一种是直接存储在.tif文件内部，这种存储方式是按照tif内部的的键值对方式来存储的。具体的一些细节可以参看这里。官方的说明pdf可以在这里下载。也可以把图像仅仅存储到tif文件中，而把坐标系统和仿射参数什么的单独提取出来，创建两个文件，一个是prj文件，存放WKT坐标系字符串，一个是tfw文件，存放仿射转换参数。
3.3. 小插曲：
关于tfw文件格式，网络流传几种说法。一种是这样的:
特别是这几行：
行 说明
1 地图单元中的一个象素在X方向上的X分辨率尺度。
2 平移量。
3 旋转量。
4 地图单元中的一个象素在Y方向上的Y分辨率尺度的负值。
5 象素1，1（左上方）的X地坐标。
6 象素1，1（左上方）的Y地坐标。
这种说法的源头可能是这个，但是他的参数个数是８个 第二种则是在gdal文档中的如下的解释
Toggle line numbers Toggle line numbers
1 WLD -- ESRI World File
2 A world file file is a plain ASCII text file consisting of six values separated by newlines. The format is:
3 pixel X size
4 rotation about the Y axis (usually 0.0)
5 rotation about the X axis (usually 0.0)
6 negative pixel Y size
7 X coordinate of upper left pixel center
8 Y coordinate of upper left pixel center
9 For example:
10 60.0000000000
11 0.0000000000
12 0.0000000000
13 -60.0000000000
14 440750.0000000000
15 3751290.0000000000
16 You can construct that file simply by using your favorite text editor.
17 World file usually has suffix .wld, but sometimes it may has .tfw, tifw, .jgw or other suffixes depending on the image file it comes with.
可以看出参数2和3和第一种的解释不一样。我现在宁可相信GDAL做的第二种解释（如果不是，那著名的grass不是白搞了？），但是这种说法也不一定，因为在gdal中直接解释的是wld文件，别名才是tfw，有没可能tfw有其他官方定义而非ESRI为惟一标准呢？如果去查google，中文和英文的解释会有很多种，越看只能越混乱。幸好大多数情况下是0 -_-! 不管了，先凑合这么用吧。如果谁知道ESRI官方解释，不妨贴上来大家共享）
3.4. GCPs
至于GCPs，那就麻烦了。不过玩过arcview或者mapinfo的配准的应该理解起来都没什么问题。要配准一张栅格图（把一个没有空间信息的图片变成有空间信息的栅格数据集），只要在一个已知的坐标系统中定位几个控制点，然后输入这几个点的地理位置。软件通过读取这几个点在图像中的位置，以及输入的这几个点对应的地理位置，建立一个拟合方程，用方程来描述其他点在这个坐标系统中的地理位置。当然，在不同投影不同坐标系不同软件中拟合方程处理可能不同。所以拟合出来的坐标可能有一些差异。这个过程不好试验也不好验证，手上也没有相关数据，所以在这里就不多说了。
4. 反馈
如果您发现我写的东西中有问题，或者您对我写的东西有意见，请一定要发邮件跟我讲，Email(linux_23@163.com)
⌨️ 快捷键说明
复制代码
Ctrl + C
搜索代码
Ctrl + F
全屏模式
F11
切换主题
Ctrl + Shift + D
显示快捷键
?
增大字号
Ctrl + =
减小字号
Ctrl + -