遥感图像与地理坐标匹配的快速算法

作者- 段辉明　卫东　葛成辉　孙卫东    来源- 清华大学电子工程系遥感研究室

摘　要　针对在遥感图像处理中, 在进行图像像素坐标与地理坐标相互转换时, 转换速度和转换精度不能同时满足
的矛盾, 引入索引矩阵的概念, 提出并实现了一种新的快速转换算法。该算法能快速、精确地进行地理坐标与图像像
素坐标的相互转换。通过MOD IS 图像验证, 可在接收图像的同时, 实时地将经、纬度网格, 省界(或州界) 等地理信息
叠加到卫星图像上。这对于要求实时处理的应用任务, 如火灾监测的快速定位、军事侦察等具有极其重要的价值。
关键词　　图像坐标　地理坐标　地理信息系统　匹配　MOD IS 图像　转换矩阵
中图法分类号: TP751. 1　　文献标识码: A 　　文章编号: 100628961 (2004) 0120087206
Fa st Coordina t ing Tran sforma t ion Between Image
of Remote and Geodet ic Coordina tes
DUAN Hu i2m ing, W E IDong, GE Cheng2hu i, SUN W ei2dong
(R em ote S ensing L ab, E lectronics E ng ineering D ep artm ent, T sing hua U niversity , B eij ing 100084)
Abstract　Coo rdinat ing t ransfo rmat ion betw een image of GIS data and geodet ic coo rdinates is p roblem frequent ly
met in p rocessing of remo te sensing images, w h ich invo lves a cont radict ion betw een p rocessing speed and conver2
sion accuracy. To so lve th is cont radict ion, a new algo rithm w h ich int roduces a new concep t ion of Index A rray is
developed base on analysis exist ing methods of regist rat ion, w h ich can be used to ret rieve p ixels of co rresponding
coo rdinates direct ly and avo id heavy computat ion need by t radit ional search ing. Th is algo rithm is characterized of
quick accurate t ransfo rmat ion betw een the two coo rdinates, and allow s real2t ime superpo sing geograph ical info,
such as geograph ical coo rdinates and boundaries, on real t ime received remo te sensing images w h ich p roves to be
co rrect by experiments inMOD IS real t ime data receiving.
Keyword　Screen coo rdinates, Geograph ical coo rdinates, MOD IS, Conversion array
1　引　言
随着科技的不断进步, 人们发射了各种各样的遥
感卫星(如NOAA、MOD IS 等) [ 1 ] , 使人们能够从太
空中来观察人类居住的地球。当得到遥感图像后, 其
中一件重要的工作就是将遥感图像与地理坐标进行
精确匹配, 使遥感图像坐标与地理信息相联系, 以利
于提取有用信息。由于遥感图像与地图的坐标系不
同, 因此必须将2 个坐标系统进行匹配, 以便将它们
统一到一个坐标系中。匹配可从如下两个方面进行:
(1) 将遥感图像与地理坐标匹配。这种方法需要
对遥感图像进行重采样, 且所需要的处理时间较长,
因此, 该方法主要用于接收遥感图像的后处理过程
中。当前流行的几个大型遥感数据处理软件, 如En2
vi、A ro Info 等都是采用的这种方法。
(2) 将地图(经纬网格、省界、河流等) 匹配到遥
感图像上, 用于实时观察遥感数据, 但是目前尚无此
类软件, 本文实现算法属于第2 种。
在这里,“匹配”也包括以下2 个方面的意义: 一
是求出图像上每个像素点的地理坐标(经纬度) , 即
在知道图像像素坐标的情况下, 求出地理坐标; 二是
将地图(经纬网格、省界、河流等) 叠加到遥感图像
上, 即在知道地理坐标情况下, 求出图像像素坐标。
对于遥感图像与地理坐标的精确匹配, 虽已有一些
方法, 但是这些方法有一个共同的问题, 那就是不能
同时满足对匹配精度和快速处理的要求, 即, 如要求
匹配精度高, 则处理速度慢; 若要求处理速度快, 则
匹配精度差。
随着遥感图像的应用越来越广, 一些对时间要求
快的应用, 需要一种快速、精确的匹配算法来将遥感图
像与地理信息联系起来, 以便及时获得更多的有用信
息, 如森林、草原火灾的探测, 由于当从遥感图像上发
现火灾情况, 需要的是立即得到火点的具体位置, 即使
不是十分精确, 而目前流行的一些大型软件( 如
Envi) , 对于MOD IS (Moderate Reso lu t ion Im aging
Spect ro radiom eter) 中分辨率成像光谱仪数据, 一般都
需要几个小时后, 才能计算出精确定位的火点位置, 因
而不能满足上述应用的要求。
本文根据遥感图像匹配原理, 通过引入索引矩
阵的概念, 提出了一种能同时满足匹配精度和处理
速度的匹配方法, 并以MOD IS 图像数据为例, 使用
V isualc C 610 为开发工具, 在PC 机上实现了这一
方法。经验证, 用该算法可以在MOD IS 数据接收的
同时, 实时显示MOD IS 图像, 并可实时叠加地理信
息(经纬网格、省界、河流等) , 从而不仅满足了对遥
感图像分析、处理的时效要求, 而且为遥感图像的应
用开辟了新领域。
2　原　理
2. 1　最简单的匹配方法
对于一幅无经纬度网格的遥感图像, 要将其与地
理坐标(经纬度) 相匹配, 最简单的做法是: 首先在图像
上找到2 个或更多的控制点(同时知道图像坐标和地
理坐标的点) , 然后按线性映射关系(每度相当于多少
公里) , 计算出图像上所有点的地理位置(经纬度) [ 2 ]。
如图1 所示: 如果知道A 、B 两点的图像坐标
(X A , YA )、(X B , YB ) , 及地理坐标(UA , V A )、(UB ,
V B ) , 则C 点的地理坐标就可由其图像坐标及A 、B
点地理坐标按线性公式计算出, 反之, C 点的图像坐
标也可由其地理坐标及A 、B 点图像坐标按线性公
式计算出。其计算公式如下:

此主题相关图片如下：


此主题相关图片如下：

点纵坐标,U k 表示图像上k 点经度, V k 表示图像上
k 点纬度, k= A ,B , C。
2. 2　简单匹配方法的缺点
虽然上述方法关系简单, 运算速度快, 但是, 有
如下2 个方面的原因限制了它的使用范围:
(1) 根据地球坐标投影原理知道, 地理坐标与图
像坐标并不满足线性关系[ 3 ] , 在不同的纬度, 1°经度所
对应的地面距离不同(图2) , 例如在赤道上, 1°经度等
于60 nm ile, 约111 km , 而随着纬度的增加, 每经度所
对应地面距离相应减少(如图2 所示) , 因此只有在纬
度跨越范围较小的情况下, 才近似满足线性条件。

此主题相关图片如下：

(2) 对于大部分遥感图像, 能很快获得的图像
是“正视图”, 即图像中, 星下点1 个像素所对应的地
面面积与图像边缘1 个像素所对应的地面面积不
同, 如在MOD IS 图像中, 星下点1 个像素对应的地
面面积约1 km 2。而在图像边缘, 1 像素对应的地面
面积约为2×4 km 2。可见, 这也不满足线性条件, 只
有在小范围内才近似满足。
2. 3　改进
上述讨论, 只是知道了简单的匹配方法所存在
的问题, 而改进这种方法, 使之能用于遥感图像匹配
则是本文要解决的问题。
既然简单匹配方法仅适用于小范围遥感图像与
地理坐标的匹配, 那么对于一幅纬度跨度范围大的
遥感图像, 要将其与地理坐标进行精确匹配, 其最直
接的方法是: 首先将图像分割成很多小面积的子图
像, 然后在每幅子图像中找到2 个控制点进行简单
匹配, 换一种说法就是利用增加控制点数量的方法
来实现整幅图像的精确匹配。具体做法是: 对于任意
给定的点, 首先求出与其距离最近的2 个控制点, 再
用这2 个控制点进行简单匹配计算。这样, 控制点越
多, 匹配精度就越高。
但是, 要使用这种方法来实现快速、精确匹配,
还必须首先解决如下2 个问题: 一是怎样快速获取
大量的控制点, 且控制点要均匀分布在遥感图像上。
二是在诸多控制点中, 如何快速找出2 个最佳控制
点(即距给定点距离最近的控制点)。
3　索引矩阵实现的原理与方法
人造卫星的运动轨迹可以由星历表计算出来, 而
每颗卫星的星历表可以从网上获得, 并能不断更新。
由于MOD IS 向下传输的是实时数据, 且数据流中包
含有时间信息, 因此可以使用时间信息, 根据星历表
来计算出卫星在某一时刻的星下点, 再根据卫星的轨
道参数、地球模型、扫描镜夹角、卫星姿态等, 计算出
此时刻卫星遥感图像中对应点的地理坐标。通过使用
这种方式, 就可以计算出遥感图像中每个像素点的地
理坐标(经纬度)。这样就得到了所需要的控制点。从
理论上讲, 控制点越多, 匹配越精确, 虽然匹配最精确
的做法是计算出图像上所有像素的地理坐标, 并将其
作为控制点, 但是在具体实施中, 会遇到一些困难, 如
一幅接收时间为15 m in 的MOD IS 图像, 约有6 000
×1 354个像素, 如果要在程序中计算出如此多的控
制点, 就目前的计算机速度而言, 则时间很长。另外,
由于小范围的遥感图像与地理坐标匹配时, 基本满足
线性条件, 不需要把图像上的每个点作为控制点, 也
能满足精度要求, 因此, 在具体匹配时, 需要根据实际
情况, 均匀选择一些点作为控制点。
3. 1　实现原理及方法
在实现方法中, 在横轴和纵轴方向, 每隔10 个
像素取一个控制点。这样, 就得到一个控制点矩阵
( 即转换矩阵) , 可称它为正向转换(已知图像坐标
x、y , 求出其地理坐标经度、纬度) 矩阵。该矩阵的元
素值表示地理坐标, 而元素在该矩阵中的位置则与
图像中相对应像素的坐标有固定关系, 其关系为

此主题相关图片如下：

其中, i 表示矩阵元素的位置是第i 行; j 表示矩阵
元素的位置是第j 列; x , y 表示图像像素坐标位置。
从原理上讲, 有了上述正向转换矩阵, 对于遥感
图像上任一个像素点, 都可以将其转换为地理坐标,
也可以将任意一个地理坐标点转换到显示图像的图
像坐标系中, 方法如下:
(1) 正向转换　即由图像坐标转换成地理坐
标。对于任意一个已知图像坐标的像素点, 若要将其
图像坐标转换成地理坐标, 则首先根据其图像坐标,
使用上述对应关系, 求出正向转换矩阵中与该像素
距离最近的2 个元素的行、列数( i, j ) , 再根据( i, j )
来得到对应的地理坐标。这样, 仅需要做几次乘ö除
法, 就可求出2 个合适的控制点, 再用线性公式, 就
可求得该像素点的地理坐标。
(2) 反向转换　即由地理坐标转换成图像坐
标。对于任意一个已知地理坐标的地面点, 若要将其
地理坐标转换成图像坐标, 则同样要根据其地理坐
标, 从正向转换矩阵中, 求出距离该地面点最近的2
个控制点。在以前的方法中, 一般都是通过遍历正向
转换矩阵的每个元素来求出该地面点与转换矩阵中
每个元素的距离, 再取出与该地面点距离最小的2
个点作为计算该地面点的控制点。
将正向转换与反向转换过程进行比较, 发现它们
的不同点在于求出2 个最佳控制点的方法不同: 正向
转换是利用图像坐标与转换矩阵元素位置的对应关
系来得到2 个最佳控制点, 其与正向转换矩阵的大小
无关; 反向转换则是用遍历正向转换矩阵的每个元素
的方法来求出2 个最佳控制点, 且矩阵越大(元素越
多) , 求出2 个最佳控制点所消耗的时间越多。正是由
于这个原因, 在遥感图像匹配时, 造成了这样一个现
象, 即在知道图像像素坐标(x , y ) 的情况下, 就能很
快计算出精确的地理坐标; 而在知道某点地理坐标的
情况下, 计算图像坐标时, 其计算精度和计算速度则
形成了一对矛盾。例如在MOD IS 数据中, 一个轨道
数据的接收时间约10～ 15m in, 其控制点矩阵的元素
约有135×600 个。一般做一次正向转换, 求出2 个最
佳控制点, 只需要做2 次乘法(或除法) ; 而在反向转
换中, 则首先要求出转换点与转换矩阵中元素的距
离, 由距离公式d = x 2+ y 2知道, 其需要2 次乘法,
一次加法, 一次开平方, 才能求出该点与转换矩阵中
某个元素的距离。也就是说, 反向转换必须先做135
×600 次这种运算, 然后, 再从求得的135×600 个距
离中, 求出最小的2 个, 这又需要做135×600 次比
较, 因此, 做一次反向转换, 至少要比做一次正向转换
多用135×600×4 倍的时间。
图4　反向转换后,MOD IS 图像与地理坐标匹配图
通过分析形成这种现象的原因发现, 由于正向转
换矩阵中, 暗中包含了一种索引关系, 即该矩阵中元
素序号的顺序与图像像素坐标的大小顺序相同或相
反, 且该矩阵中2 个元素之间的距离, 与图像坐标系
中2 个点的距离成正比, 如在实现的算法中, 该矩阵
中2 个元素位置坐标的行差为5, 列差为4, 则对应图
像上这2 个点的横坐标相差5×10 个像素, 纵坐标相
差4×10 个像素。为说明这个概念, 假设正向转换矩
阵是由将图像中每个像素都作为控制点而产生的, 那
么在这种情况下, 图像像素坐标和该矩阵元素位置坐
标的步进一样, 都等于1, 因此只要将两者零点重合,
那么图像像素坐标就等于该矩阵元素位置坐标。这样
在知道图像像素坐标(x , y ) 的情况下, 直接用图像像
素坐标(x , y ) 作为该矩阵元素的位置坐标, 就可求出
对应的地理坐标(如图3 所示)。

此主题相关图片如下：

其对应关系为, F (X k , Y k ) Z f (U k ,V k ) 。
根据上面的分析, 如果存在这样一个矩阵, 即矩
阵元素的位置坐标与地理坐标有一定的对应关系, 也
就是矩阵元素值即表示图像的像素坐标, 那么, 在知
道某点地理坐标的情况下, 计算图像像素坐标, 就不
必用遍历矩阵的方法来求出2 个最佳控制点, 而只需
根据这种对应关系, 就能按正向转换的方法求出2 个
最佳控制点。这样, 计算速度就与该矩阵的大小无关,
从而解决了转换精度与处理速度的矛盾。这样的矩阵
称为反向转换矩阵, 即产生如下一种新的对应关系
F (U k ,V k ) Z f (X k , Y k )
　　如何得到这个反向转换的索引矩阵? 可利用正
向转换矩阵来产生反向转换矩阵: 首先确定一个位
置坐标为(0, 0) 的元素, 并将其与一个地理坐标关
联, 再计算出其图像坐标, 并将坐标值存入索引矩阵
元素中; 然后, 每间隔1°, 计算一个图像坐标值存入
索引矩阵对应元素中, 即1 个转换索引矩阵元素的
位置坐标差值, 对应地理坐标中1°(为了提高精度,
也可以是其他更小值, 如015°)。
至此, 就得到了2 个转换索引矩阵(或称为查找
表) , F (X k , Y k ) Z f (U k , V k ) 和F (U kV k ) Z f (X k ,
Y k )。其中一个用于正向转换(已知图像坐标, 求地理
坐标) , 另一个用于反向转换(已知地理坐标, 求图像
坐标)。
使用这2 个索引矩阵进行坐标转换, 虽然转换
速度可以极大提高, 但是, 如果仅做以上的处理, 反
向转换出的结果误差较大, 且有明显的几何畸变, 其
原因如图4 所示, 由图4 可见, 当经纬度网格叠加到
卫星图像上时, 其网格呈菱形状。在计算中, 如将其

此主题相关图片如下：

当成矩形计算(图4 中白线矩形) , 就会产生较大的
误差(如图4 所示) , 可是如果完全按菱形修正, 则不
仅计算复杂, 而且耗费时间, 为此在本文的算法中,
则是将菱形近似成旋转了一个角度的矩形, 但仔细
观察, 会发现水平方向和垂直方向的旋转角度不同。
为了进一步提高精度, 可以通过分别计算出水平方
向和垂直方向的夹角来进行修正[ 4 ] , 其最后结果图
5 所示。

此主题相关图片如下：

3. 2　误差分析
经MOD IS 图像匹配验证, 使用本文算法, 其匹
配误差平均在3 个像素左右。该误差产生主要有如
下2 个来源:
(1) 产生正向转换矩阵时, 计算出来的控制点
地理坐标与实际地理坐标有一些误差, 其最大误差
为5 个像素点, 一般在0～ 5 个像素之间随机产生。
这是算法误差产生的主要来源, 今后如要改进此算
法, 应从这里着手。
(2) 算法假定图像中每个像素对应的地面面积
相等, 而实际情况并非如此, 因为图像中每个像素对
应的地面面积随扫描角变化而变化, 如在MOD IS
图像中, 星下点1 个像素对应地面面积约1 km2; 而
在图像边缘, 1 个像素对应地面面积约2×4 km 2。但
是这个匹配误差, 随控制点的增加而减少。
产生匹配误差还有一些其他因素, 如地面高程
的影响、星历表产生时间与数据接收时间的间隔等,
但是由于这些因素对匹配结果的影响不是太大, 因
此算法中没有考虑。
4　结论与说明
本文通过引入索引矩阵的概念, 实现了一种能
同时满足匹配精度和处理速度的匹配算法。此算法
从原理上属于查表法, 与以前算法不同之处在于其
动态产生了2 个查找表, 一个用于正向查找, 一个用
于反向查找, 而以前的方法则不论正向查找或是方
向查找, 都使用同一张查找表。本文算法的特点是:
在保证基本精度的情况下, 坐标匹配速度快, 且运算
速度与图像大小无关, 尤其适用数据量大的遥感图
像。此特点对于火点探测、军事侦察等应用有及其重
要的价值。该算法虽然本文只是对MOD IS 图像处
理, 但是可以推广到其他中、低分辨率卫星图像, 因
为此算法唯一与MOD IS 数据有直接关系的就是星
历表, 只要将星历表换为其他卫星的星历表, 则算法
仍然成立。