ADS100空三工程任意范围裁切算法研究

梁 菲

(1.中煤航测遥感集团有限公司，陕西 西安 710199； 2.陕西省地理空间信息工程技术研究中心，陕西 西安 710199)

随着国民经济发展对空间数据获取的现势性要求越来越高，对空间数据产品生产的效率也有了较高的要求，内业数据处理不断在向自动化发展，同时尽量减少外业的工作量，提高作业效率，缩短生产周期。相比传统的框幅式航空摄影测量，ADS(Airborne Digital Sensor)系列三线阵传感器集成了高精度的POS数据，仅需要少量的控制点，就可以得到高精度的空三成果，作业过程可以减少大量的外业工作量；同时获取到前视、下视以及后视三个视角的影像条带，任意两个视角都可以构建立体像对；航线飞行方向可以直接获取到整条航线的影像，可飞行50～60 km的数据量[1]。独特的优势让ADS航摄仪为航空摄影测量自动化开辟了崭新的途径，成为对地观测空间数据获取的主要数据来源方式。

ADS100在数据处理中虽然自动化程度高，外业工作量少，但其巨大的数据量给数据处理也带来了很大的挑战。在生产中，影像数据以航线为单位存储，每条航线的数据量是非常大的，在DEM、DOM、DLG等数据测绘产品实际生产时，模型构建速度比较慢，在交互过程中出现卡顿现象，影响数据生产效率。另外，ADS航摄仪数据处理时所需文件多且组织复杂，格式及各个参数的物理意义不明确，在实际生产应用中，受到数据格式的限制，对数据重新组织比较困难，无法实现自由处理数据。目前现有文献的研究主要分为两个方向，一是对数据处理理论进行研究[2-4]，二是对其生产体系及质量进行评价[5-7]，还未有针对生产实际中的具体问题进行研究实现的。通过研究ADS100航摄仪空三加密成果中工程参数的物理意义(见实例)，提出并实现了一种基于任意范围线的ADS100加密工程裁切方法，并已经应用到数据生产中，极大地提高了构建空三模型的速度，模型影像调度非常流畅，提高了数据生产效率。

1 技术路线

ADS100加密工程裁切以推扫式影像成像原理为基础，充分解析.sup、.ads文件中参数的物理含义，根据范围线对影像进行裁切，重新计算工程参数并更新相应文件。其技术流程如图1所示。

图1 工程裁切完整流程图

(1)基础数据准备

获取测区的空三加密成果、作业范围线数据、WGS84坐标系与范围线坐标系两套坐标系相应的控制点若干。

(2)WGS84与范围线坐标系统转换

ADS100航摄仪集成了高精度GPS和惯性测量装置(IMU)，其初始坐标系统为WGS84，在空三加密完成后，加密成果还是WGS84坐标系统[8]。而在实际生产中作业范围线一般为地方坐标系统，因此就需要在两套空间坐标系统之间进行互相转换。

算法实现时要将ADS数据转换为地方坐标系统，首先将WGS84大地坐标系转换到WGS84空间直角坐标系；然后根据七参数将WGS84空间直角坐标系转换到地方坐标系。不同空间直角坐标系之间的转换七参数计算需要使用原始坐标系(WGS84)、地方坐标系中相应的控制点坐标。其转换过程为：首先将控制点坐标，分别按照相应的椭球及投影参数，将坐标都转换到地心(参心)直角坐标系中，其坐标记为(X,Y,Z)；在地心(参心)直角坐标系中进行两套控制点之间的转换参数计算，即为计算转换七参数，再利用七参数将原始坐标系的数据转换到目标坐标系中。

(3)计算影像与范围线交集

根据工程文件中的影像信息，计算每一块影像在地方坐标系中的地面坐标范围，判断过程为：基于多边形之间的空间关系(包含、相交、外部)，判断多边形之间是否有交集；若有交集，利用点与多边形的空间关系(内部、外部)，采用基于面元的快速检测，计算出交集内的影像范围。

(4)重新计算工程参数

计算要裁切的影像范围在原始影像中的位置，重新计算工程参数中的影像大小、坐标偏移量；输出裁切后的影像；更新ads文件中的影像块信息。

2 关键技术

2.1 影像块范围计算

由于ADS100航摄仪影像数据是按照ads文件中的Block存储的，读取每个Block的影像数据，并计算Block影像的四个角的地面坐标，公式中影像列Sample表示为S，影像行Lines表示为L，影像宽度为W。

(1)影像行列号像素坐标计算

计算每个影像块的四个角在整个条带中的像素坐标。

左上角：

Si=Si-1+Wi-1

Li=L

(1)

式中，Si、Li为第i个块的坐标；Si-1为第i-1块的左下角的像素坐标；Wi-1为上一个块的影像宽度；L是整个条带的高度(每个块的影像高度都是相同的)。

左下角：

Si=Si-1+Wi-1

Li=0

(2)

右上角：

Si=Wi+Si-1+Wi-1

Li=L

(3)

式中，Wi为当前块的影像宽度。

右下角：

Si=Wi+Si-1+Wi-1

Li=0

(4)

(2)影像物方坐标计算

步骤1：获取到工程中的基本参数有锚点坐标B0、L0，旋转参数r，缩放比例s，坐标偏移量xt0，yt0。锚点坐标是WGS84坐标系下的经纬度格式，要将其转换到WGS84空间直角坐标系，用X0，Y0，Z0表示。

步骤2：计算物方空间坐标

根据每影像块的四角像素坐标，计算四角的物方坐标Gx，Gy，得到影像的覆盖范围[9]：

Gx=((Si+xt0)cos(r)+(Li+yt0)sin(r))/s

Gy=(-(Si+xt0)sin(r)+(Li+yt0)cos(r))/s

(5)

步骤3：计算地面坐标

根据基本参数信息、物方空间坐标计算每个影像块的地面覆盖范围

X=X0-Gy×sin(B0)×cos(L0)-sin(L0)×Gx

Y=Y0-Gy×sin(B0)×sin(L0)+cos(L0)×Gx

Z=Z0+Gy×cos(B0)

(6)

此时，该地面坐标为WGS84地心直角坐标系，可以直接利用七参数将该地面坐标转换到目标坐标系椭球的地心直角坐标系。将该地心直角坐标转换为经纬度坐标，然后按照目标坐标系的投影方式，投影为空间直角坐标。可以得到转换到目标坐标系下影像的地面覆盖范围。

2.2 基于面元的空间关系快速检测

由于每条航线的ads中会有比较多的影像块，在计算过程中需要循环影像块进行逐块计算影像范围。影像块与范围线关系有以下三种：

(1)当前影像块不包含在范围线内：不计算裁切范围(外部)；

(2)当前影像块完全包含在范围线内：整块影像计算裁切范围(包含)；

(3)当前影像块部分包含在范围线内：计算影像块与范围线交集(相交)。

基于点与多边形的空间关系，采用基于面元的快速检测法计算与范围线相交的影像范围。在实际生产中范围线一般是不规则的多边形，而影像的存储是规则矩形，因此在这里求出的交集是按照影像范围为矩形计算。由于影像数据量比较大，为了提高检测速度，提出了将原始影像抽稀为面元，以面元为单位计算与范围线的交集，面元大小可以设置为50～100个像素之间，如图2所示。其伪代码如下：

FOR(Inti=0,i<影像行,i+=面元大小){

FOR(Intj=0,j<影像列,j+=面元大小){

计算当前像素的地面坐标，判断该地面点是否在范围线内，利用点是否在任意多边形范围内进行判断；

IF(该地面坐标在范围线内){

记录最小的行列号坐标，记录最大的行列号坐标}}}

图2 面元检测算法

按照计算的每个影像块的最大最小坐标，计算影像块的裁切范围，并重新存储裁切后的影像。采用将像素融合为面元的方式，极大地提高了像素与多边形空间关系检测效率。

2.3 重新计算工程参数

(1)ads文件中存储的是当前航线所有影像块的参数信息，根据当前航线每个影像块计算的裁切范围，要重新计算LINES、SAMPLES、LINES_PER_BLOCK、SAMPLES_PER_BLOCK、BLOCK_DATA。要注意的是，ads文件中默认每个影像块的行列数都与参数相同，但是不包含最后一个影像块，最后一个影像块就按照裁切后的实际尺寸生成；另外，第一个影像块尺寸如果与参数不一致时，需要将尺寸调整到与参数一致。最后按照裁切后的所有影像块，生成相应的BLOCK_DATA参数，并生成新的ads文件。

(2)在SUP文件中要重新输出ads、相机、外方位元素的存储路径，重新计算坐标偏移量RECT_XOFFSET、RECT_YOFFSET参数，公式中表示为xt0，yt0。首先利用公式(5)计算裁切后第一个影像块的物方空间坐标，然后利用公式(7)计算出xt0，yt0。

yt0=(sin(r)×s×Gx+cos(r)×s×Gy)/(sin(r)2

+cos(r)2)

xt0=(yt0×cos(r)-s×Gy)/sin(r)

(7)

更新工程文件中的RECT_XOFFSET，RECT_YOFFSET为xt0，yt0。SUP中其他的参数按照裁切后的相应参数更新并写入文件。

3 实例分析

3.1 试验数据

ADS100航摄仪空三加密完成后，为每条航线都生成了加密后的SUP工程文件、ads文件以及odf.adj文件。通过研究ADS100成像原理，解析其加密工程文件，获取到工程中的基本参数信息，如表1所示。

表1 基本参数信息及其物理含义

ADS100获取到的是整条条带的影像，其数据量可达到10 G左右，因此在处理时都会将整条条带影像分块存储，并将其分块信息存储到ads文件中。读取加密工程相对应的ads文件，获取当前航线影像所包含的分块信息blockdata(block0、block1、block2、……、blockn)。

试验选取某测区ADS100航测数据，坐标系统为当地地方独立坐标系，需要进行WGS84坐标与地方独立坐标系的坐标转换。试验测区提供了349个控制点，且在测区范围内均匀分布，根据每个控制点相应的两套坐标，即可求解出七参数。选取测区中一条航线的两个视角(前视和下视)为试验数据，按照两种范围线对当前航线的空三加密工程进行裁切。航线数据基本情况如表2所示。

表2 数据基本情况/km

3.2 裁切实现

按照上述原理在VS2010的平台上实现了ADS100加密工程裁切，并对试验数据进行裁切试验，试验选择一小部分矢量数据作为范围线区域，将裁切前和裁切后的结果对比，如表3所示。

表3 数据裁切前后对比情况

采用上述方法裁切后的加密工程能够无缝兼容入市场上主流的ADS100的数据处理软件，将其导入软件中，裁切前与裁切后的数据套合情况结果如图3、图4所示。

图中，红色线框为地方独立坐标系要裁切的范围线，图3为前视航线裁切前后与范围线的套合情况，图4为下视航线裁切前后与范围线的套合情况，其(b)图中，将范围线外的数据裁切之后，采集的矢量数据与裁切后影像完全套合。

3.3 精度分析

精度验证采用模型检测点、外业检测点对裁切后的模型精度进行验证。模型检测点是在原始空三模型上人工量测的点位，在测区范围内选择了25个均匀分布的检测点；外业检测点在裁切范围内有3个点；人工在裁切后的模型上量测了裁切后的模型点位，并计算了相应检测点的最大差值、中误差，如表4所示。

图3 范围线与前视影像裁切前后套合图

图4 范围线与下视影像裁切前后套合图

表4 模型裁切后精度验证

从计算结果来看，按照1:1 000的模型精度要求，丘陵地区平面小于0.35 m，高程小于0.35 m，外业检测点与模型检测点的检测结果均满足要求。

4 结 语

ADS航摄仪独特的成像方式、不同于常规面阵相机的数据组织与存储，海量的影像数据给后期数据处理带来了很大挑战[10]。基于ADS100的成像原理，通过研究其数据组织及参数的物理意义，实现了基于任意范围线对ADS100加密工程进行裁切，从而提高实际生产效率。目前算法仅实现空三工程的裁切功能，但其研究基础可以扩展到数据处理的其它环节中，对提高实际生产效率具有现实的指导意义。