CN113739767A - 针对国产面阵摆扫成像系统获取的影像生产正射影像的方法 - Google Patents

针对国产面阵摆扫成像系统获取的影像生产正射影像的方法 Download PDF

Info

Publication number
CN113739767A
CN113739767A CN202110971773.4A CN202110971773A CN113739767A CN 113739767 A CN113739767 A CN 113739767A CN 202110971773 A CN202110971773 A CN 202110971773A CN 113739767 A CN113739767 A CN 113739767A
Authority
CN
China
Prior art keywords
image
coordinate system
point
tpc
coordinate
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Granted
Application number
CN202110971773.4A
Other languages
English (en)
Other versions
CN113739767B (zh
Inventor
段延松
赵新博
张祖勋
柯涛
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Wuhan University WHU
Original Assignee
Wuhan University WHU
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Wuhan University WHU filed Critical Wuhan University WHU
Priority to CN202110971773.4A priority Critical patent/CN113739767B/zh
Publication of CN113739767A publication Critical patent/CN113739767A/zh
Application granted granted Critical
Publication of CN113739767B publication Critical patent/CN113739767B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • GPHYSICS
    • G01MEASURING; TESTING
    • G01CMEASURING DISTANCES, LEVELS OR BEARINGS; SURVEYING; NAVIGATION; GYROSCOPIC INSTRUMENTS; PHOTOGRAMMETRY OR VIDEOGRAMMETRY
    • G01C11/00Photogrammetry or videogrammetry, e.g. stereogrammetry; Photographic surveying
    • G01C11/04Interpretation of pictures

Landscapes

  • Engineering & Computer Science (AREA)
  • Multimedia (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Radar, Positioning & Navigation (AREA)
  • Remote Sensing (AREA)
  • Image Processing (AREA)
  • Processing Or Creating Images (AREA)
  • Studio Devices (AREA)

Abstract

本发明涉及一种针对国产面阵摆扫成像系统获取的影像生产正射影像的方法。首先选择合理的局部坐标系作为基准坐标系,将导航数据信息转换到基准坐标系中,然后在DEM辅助下对所有影像进行同名点自动提取,形成影像连接点网,之后使用连接点开展基于DEM的区域网平差,解算每张影像的外方位参数,最后基于最优化理论生产目标区域的整体最优正射影像。本发明引入GPS/IMU参数作为初值,使用DEM数据辅助平差,解决了传统摄影测量无法直接进行平差解算的问题。本发明可以更快速、更准确地完成影像拼接,生产合格的正射影像产品,填补国产面阵摆扫成像系统数据处理的空白。

Description

针对国产面阵摆扫成像系统获取的影像生产正射影像的方法
技术领域
本发明属于遥感测绘技术领域,特别是涉及一种针对国产面阵摆扫成像系统获取的影像生产正射影像的方法。
背景技术
面阵摆扫成像系统是一种快速获取地面影像信息的设备,与传统固定框幅成像系统和多视倾斜成像系统相比,数据获取效率高,整体成本低,具有很强的实用性。早在2008年,国际上就已经推出商用系统,典型代表是以色列的A3面阵摆扫成像系统和日本的SamplX面阵摆扫成像系统。这些系统已经形成商业化产品,对外提供了完整解决方案,方案包含了从系统航拍安装到最终测绘产品(包括空中三角测量、数字高程模型、正射影像、数字线划图)的全流程生产过程操作说明,但是在具体技术细节上,如空中三角测量处理过程和技术细节等方面没有提供相关信息。厂商在设备说明中也仅提供了如相机分辨率、相机焦距、影像格式等设备输出成果参数,对于更核心的参数如成像系统POS系统工作原理、参数等都没有提供,只提供了对系统获取的影像信息开展生产的软件。经过国内相关技术领域的科研院所10年的努力,成功研制出国产面阵摆扫成像系统。国产面阵摆扫成像系统由成像硬件和处理软件组成,成像硬件主要由镜头、成像CCD面板、GPS测位设备、IMU测姿设备、摆扫驱动设备、存储设备等核心部件组成,如图1所示。处理软件主要包含影像信号解码、GPS数据处理、IMU数据处理以及测绘产品生产模块。
在成像原理上,面阵摆扫成像系统与传统成像系统存在较大差异,其最大特点是连续摆扫成像,摆扫周期内,任意两张影像可近似理解为同心成像,可无缝拼接,重叠区域为同光线成像,不能进行立体定位。此外,国产面阵摆扫成像系统所使用的IMU受自主技术限制,精度只能做到5~10角秒,无法实现影像精准定位定姿态。因此,面阵摆扫成像系统的处理软件需要充分考虑以上特点。此外,面阵摆扫成像系统的每张影像都具有GPS信息和IMU信息,可为处理软件提供重要的信息。
目前现有处理技术没有充分考虑摆扫影像的特点,特别是同心成像无法进行立体定位的特性,因此现有系统无法生产合格产品。本发明提出的方法针对国产面阵摆扫成像系统获取影像连续摆扫成像、影像重叠小、IMU精度不高等特点,引入DEM数据辅助处理,在处理效率和成果精度之间寻找平衡,可以更快速、更准确地完成影像拼接,生产合格的正射影像产品,有效填补国产面阵摆扫成像系统数据处理的空白。
发明内容
本发明针对现有技术的不足,提供一种针对国产面阵摆扫成像系统获取的影像生产正射影像的方法。首先选择合理的局部坐标系作为基准坐标系,将导航数据信息转换到基准坐标系中,然后在DEM辅助下对所有影像进行同名点自动提取,形成影像连接点网,之后使用连接点开展基于DEM的区域网平差,解算每张影像的外方位参数,最后基于最优化理论生产目标区域的整体最优正射影像。
为了达到上述目的,本发明提供的技术方案是一种针对国产面阵摆扫成像系统获取的影像生产正射影像的方法,包括以下步骤:
步骤1,建立局部坐标系作为统一的基准坐标系,将导航数据信息转换到基准坐标系中;
步骤2,在概略DEM辅助下进行影像同名点自动提取;
步骤2.1,在概略DEM辅助下,根据摆扫成像周期性和先后顺序获取影像间的相邻关系;
步骤2.2,根据步骤2.1获取的影像相邻关系,采用基于两类膨胀进行立体影像密集匹配的方法开展同名点匹配;
步骤3,基于DEM辅助进行区域网平差,将步骤1转换处理后的导航数据信息作为平差初始值,解算每张影像的外方位元素参数;
步骤3.1,利用所有影像的POS信息计算区域网范围,根据区域网范围从全球SRTM数据中获取对应的DEM数据;
步骤3.2,利用外方位角元素计算每张影像的旋转矩阵,使用步骤1得到的角元素作为计算的初始值;
步骤3.3,以单束同名光线为观测单元,分别计算各连接像点光线与DEM的交点,得到对应的物方坐标;
步骤3.4,计算同名像点对应的物方点间的差异大小,如果大于角元素精度推导的理论物方定位精度α,则将同名像点标记为粗差点,不参与平差;否则,取同名像点对应的物方坐标的平均值作为真值,参与平差;
步骤3.5,逐点构建误差方程,解算外方位元素;
步骤3.6,重复步骤3.2至步骤3.5直到影像的每个角元素的变化不超过10-6弧度且所有像点改正数的平方和小于设定的阈值λ;
步骤4,整体最优正射影像生产;
步骤4.1,计算每张原始影像四个角点对应的地面坐标;
步骤4.2,选择正射影像对应的最佳原始影像进行反投影;
步骤4.3,反投影法生产正射影像。
而且,所述步骤1中局部坐标系采用定点切平面坐标系统,即TPC坐标系(TangentPlane Coordinate System),坐标定义如下:在WGS84坐标框架下,选择飞行区域中几何中心点经度L0,纬度B0作为平面原点,以地球椭球高为0作椭球切平面,定义过原点纬度线的东方向为X坐标轴,过原点经度线北方向为Y坐标轴,过原点背向地心的垂线为Z坐标轴,如图4所示。
TPC坐标系与像空间坐标系之间存在转换关系如下:
Figure BDA0003226085700000031
即从像空间坐标系→载荷本体坐标系→IMU本体坐标系→导航坐标系→地心地固坐标系→TPC坐标系,各坐标系定义及转换如图5所示。
各坐标系间旋转矩阵计算公式如下:
Figure BDA0003226085700000032
Figure BDA0003226085700000033
Figure BDA0003226085700000041
Figure BDA0003226085700000042
Figure BDA0003226085700000043
式中,
Figure BDA0003226085700000044
为TPC坐标系下的影像姿态旋转矩阵,
Figure BDA0003226085700000045
是影像的三个外方位角元素,为待求解参数;
Figure BDA0003226085700000046
是地心地固坐标系与TPC坐标系之间的转换矩阵,由TPC坐标原点的经度L0,纬度B0定义;
Figure BDA0003226085700000047
是导航坐标系与地心地固坐标系之间的转换矩阵,由载荷所在位置的经度Lz、纬度Bz定义;
Figure BDA0003226085700000048
是导航坐标系和IMU本体坐标系之间的旋转矩阵,(Φ,Θ,Ψ)是IMU记录的俯仰、横滚、偏航数值;
Figure BDA00032260857000000412
是载荷本体坐标系与IMU本体坐标系之间的转换矩阵,又称为设备安置矩阵,由IMU安置在载荷上的ex、ey、ez三个角定义;
Figure BDA00032260857000000413
是像空间坐标系与载荷本体坐标系之间的转换矩阵。
求得旋转矩阵
Figure BDA0003226085700000049
后,反解可得到影像在TPC坐标系中的姿态角,即:
Figure BDA00032260857000000410
Figure BDA00032260857000000411
ω=arcsin(-b3) (9)
Figure BDA0003226085700000051
GPS点位置信息可以使用地球坐标定义模型无损地转换到任意需要的坐标系中,从原始的GPS坐标(L,B,H)转换到局部切平面坐标系中的TPC坐标经过如下变换:
首先将原始的GPS数据转换到地心坐标系中:
Xc=(N+H)cos B cos L (11)
Yc=(N+H)cos B sin L (12)
Zc=[N(1-e2)+H]sin B (13)
Figure BDA0003226085700000052
Figure BDA0003226085700000053
式中,N为卯酉圈的半径,a为地球椭球的长半轴,b为地球椭球的短半轴,(Xc,Yc,Zc)为地心坐标系的坐标。
然后将局部切平面坐标系原点坐标(L0,B0,0)带入式(11)~(13),求得局部切平面坐标系原点在地心坐标系中的坐标,记为
Figure BDA0003226085700000057
从地心系坐标到TPC坐标系的计算公式为:
Figure BDA0003226085700000054
式中,
Figure BDA0003226085700000055
是地心地固坐标系与TPC坐标系之间的转换矩阵,由式(2)计算可得;(XTPC,YTPC,ZTPC)为TPC坐标系的坐标。
Figure BDA0003226085700000056
即为转换到基准坐标系下的导航信息。
而且,所述步骤2中根据摆扫成像周期性和先后顺序获取影像间的相邻关系是先将摆扫一个周期得到的N张影像分为一组,每张影像根据其在摆扫周期中的次序将它们编号为1至N,顺序连续的影像必为邻接影像。各组之间的相邻关系可根据GPS和IMU信息,将获取的影像投影到全球DEM表示的地面上,然后根据每张影像的地面范围求影像重叠情况,也可根据重叠区域反算出原始影像的重叠区域,以此确认各组间影像的相邻关系。
而且,所述步骤3.3中以单束光线为观测单元,分别计算各连接点与DEM的交点,得到对应的物方坐标,计算公式如下:
Figure BDA0003226085700000061
Figure BDA0003226085700000062
式中,(Xg,Yg)为交点的物方坐标(即TCP坐标);zdem为交点的DEM高程;Xs、Ys、Zs是影像的外方位线元素,步骤1得到的(XTPC,YTPC,ZTPC)可作为外方位线元素的初始值;(x,y)是像点坐标;f是相机焦距。
而且,所述步骤3.4中理论物方定位精度α的计算方式如下:
Figure BDA0003226085700000063
式中,H为摄影位置到地面的高度,ε为成像系统的IMU角元素精度,θ为主光轴与重力方向的夹角。
而且,所述步骤3.5中逐点构建误差方程为:
Figure BDA0003226085700000064
式中,(x,y)是像点坐标,(x)、(y)是用未知数初值获取的像点坐标,vx、vy是误差变量,
Figure BDA0003226085700000065
表示求偏导,d代表误差改正,求解出外方位元素的改正数之后即可更新每张影像的外方位元素。
而且,所述步骤4.3中反投影法生产正射影像的方法如下:设正射影像上任一点P的坐标为(X′,Y′),由正射影像左下角图廓点地面坐标(X0,Y0)与正射影像比例尺分母M,计算P点所对应的地面坐标(X,Y)如下式所示:
X=X0+M×X′ (22)
Y=Y0+M×Y′ (23)
计算正射影像上点P对应原始影像上相应像点坐标p(x,y)为:
Figure BDA0003226085700000071
Figure BDA0003226085700000072
式中,X0、Y0为所有原始影像的角点对应的地面坐标中的最小值,Z为P点的高程,由DEM内插求得;ai,bi,ci(i=1,2,3)是影像的旋转矩阵元素,Xs、Ys、Zs是影像的线元素,f是原始影像的焦距。
由于所得的像点坐标不一定正好落在像元素中心,为此需进行灰度内插,本发明采用双线性内插的方法求得像点p的灰度值g(x,y),将像点的灰度值赋给正射影像上对应的像素点P(X′,Y′),即:
G(X′,Y′)=g(x,y) (26)
依次对每个像元进行上述运算,即能获得正射影像。
与现有技术相比,本发明具有如下优点:
(1)国产面阵摆扫成像系统连续摆扫成像,利用其成像周期内的影像连续即相邻的特点,使用GPS/IMU与DEM数据求交获取影像邻接矩阵和重叠区域,针对性开展影像匹配,快速提取影像连接点。
(2)国产面阵摆扫成像系统的POS精度欠佳,用于直接地理定位精度有限,不能直接生产正射影像,仍需要平差处理;但该系统的影像之间重叠度小、交会条件弱,本发明引入GPS/IMU参数作为初值,使用DEM数据辅助平差解决了传统摄影测量无法直接进行平差解算的问题。
(3)本发明可以更快速、更准确地完成影像拼接,生产合格的正射影像产品,填补国产面阵摆扫成像系统数据处理的空白。
附图说明
图1为本发明实施例国产面阵摆扫成像系统核心部件示意图。
图2为本发明实施例技术流程图。
图3为本发明实施例IMU旋转角的国际标准定义示意图。
图4为本发明实施例中局部切平面坐标系(TPC坐标系)定义示意图。
图5为本发明实施例中IMU坐标系转换为TPC坐标系的示意图。
图6为本发明实施例国产面阵摆扫成像系统摆扫成像的GPS/IMU参数转换结果示意图。
图7为本发明实施例影像邻接矩阵示意图。
图8是本发明实施例的正射影像数据成果。
图9是本发明实施例的正射影像数据拼接部位的局部放大图。
具体实施方式
本发明提供涉及一种针对国产面阵摆扫成像系统获取的影像生产正射影像的方法,首先选择合理的局部坐标系作为基准坐标系,将导航数据信息转换到基准坐标系中,然后在DEM辅助下对所有影像进行同名点自动提取,形成影像连接点网,之后使用连接点开展基于DEM的区域网平差,解算每张影像的外方位参数,最后基于最优化理论生产目标区域的整体最优正射影像。
下面结合附图和实施例对本发明的技术方案作进一步说明。
如图2所示,本发明实施例的流程包括以下步骤:
步骤1,建立局部坐标系作为统一的基准坐标系,将导航数据信息转换到基准坐标系中。
国产面阵摆扫成像系统安装有GPS设备和IMU设备,GPS设备记录飞行器位置信息,包含经度、纬度和高程,IMU记录飞行器姿态信息,包含俯仰、横滚、偏航3个角度值,IMU坐标系按国际标准进行定义,如图3所示。
局部坐标系采用定点切平面坐标系统,即TPC坐标系(Tangent Plane CoordinateSystem),坐标定义如下:在WGS84坐标框架下,选择飞行区域中几何中心点经度L0,纬度B0作为平面原点,以地球椭球高为0作椭球切平面,定义过原点纬度线的东方向为X坐标轴,过原点经度线北方向为Y坐标轴,过原点背向地心的垂线为Z坐标轴,如图4所示。
TPC坐标系与像空间坐标系之间存在转换关系如下:
Figure BDA0003226085700000091
即从像空间坐标系→载荷本体坐标系→IMU本体坐标系→导航坐标系→地心地固坐标系→TPC坐标系,各坐标系定义及转换如图5所示。
各坐标系间旋转矩阵计算公式如下:
Figure BDA0003226085700000092
Figure BDA0003226085700000093
Figure BDA0003226085700000094
Figure BDA0003226085700000095
Figure BDA0003226085700000096
式中,
Figure BDA0003226085700000097
为TPC坐标系下的影像姿态旋转矩阵,
Figure BDA0003226085700000098
是影像的三个外方位角元素,为待求解参数;
Figure BDA0003226085700000101
是地心地固坐标系与TPC坐标系之间的转换矩阵,由TPC坐标原点的经度L0,纬度B0定义;
Figure BDA0003226085700000102
是导航坐标系与地心地固坐标系之间的转换矩阵,由载荷所在位置的经度Lz、纬度Bz定义;
Figure BDA0003226085700000103
是导航坐标系和IMU本体坐标系之间的旋转矩阵,(Φ,Θ,Ψ)是IMU记录的俯仰、横滚、偏航数值;
Figure BDA0003226085700000104
是载荷本体坐标系与IMU本体坐标系之间的转换矩阵,又称为设备安置矩阵,由IMU安置在载荷上的ex、ey、ez三个角定义;
Figure BDA0003226085700000105
是像空间坐标系与载荷本体坐标系之间的转换矩阵。
求得旋转矩阵
Figure BDA0003226085700000106
后,反解可得到影像在TPC坐标系中的姿态角,即:
Figure BDA0003226085700000107
Figure BDA0003226085700000108
ω=arcsin(-b3) (9)
Figure BDA0003226085700000109
GPS点位置信息可以使用地球坐标定义模型无损地转换到任意需要的坐标系中,从原始的GPS坐标(L,B,H)转换到局部切平面坐标系中的TPC坐标经过如下变换:
首先将原始的GPS数据转换到地心坐标系中:
Xc=(N+H)cos B cos L (11)
Yc=(N+H)cos B sin L (12)
Zc=[N(1-e2)+H]sin B (13)
Figure BDA00032260857000001010
Figure BDA0003226085700000111
式中,N为卯酉圈的半径,a为地球椭球的长半轴,b为地球椭球的短半轴,(Xc,Yc,Zc)为地心坐标系的坐标。
然后将局部切平面坐标系原点坐标(L0,B0,0)带入式(11)~(13),求得局部切平面坐标系原点在地心坐标系中的坐标,记为
Figure BDA0003226085700000114
从地心系坐标到TPC坐标系的计算公式为:
Figure BDA0003226085700000112
式中,
Figure BDA0003226085700000115
是地心地固坐标系与TPC坐标系之间的转换矩阵,由式(2)计算可得;(XTPC,YTPC,ZTPC)为TPC坐标系的坐标。
Figure BDA0003226085700000113
即为转换到基准坐标系下的导航信息。
步骤2,在概略DEM辅助下进行影像同名点自动提取。
针对国产面阵摆扫成像系统获取的影像,可以充分利用其周期摆扫规律和GPS/IMU参数,配合全球DEM数据(主要是SRTM数据),实现快速有效的影像匹配。面阵摆扫成像系统的影像和其GPS和IMU参数投影到全球DEM数据处理原理和效果示意图如图6所示。
步骤2.1,在概略DEM辅助下,根据摆扫成像周期性和先后顺序获取影像间的相邻关系。
首先,将摆扫一个周期得到的N(本实施例中,N取36)张影像分为一组,每张影像根据其在摆扫周期中的次序将它们编号为1至N,顺序连续的影像必为邻接影像。各组之间的相邻关系可根据GPS和IMU信息,将获取的影像投影到全球DEM表示的地面上,然后根据每张影像的地面范围求影像重叠情况,也可根据重叠区域反算出原始影像的重叠区域,以此确认各组间影像的相邻关系,最终确认的邻接矩阵如图7所示。
步骤2.2,根据步骤2.1获取的影像相邻关系,采用基于两类膨胀进行立体影像密集匹配的方法开展同名点匹配。
步骤3,基于DEM辅助进行区域网平差,将步骤1转换处理后的导航数据信息作为平差初始值,解算每张影像的外方位元素参数。
步骤3.1,利用所有影像的POS信息计算区域网范围,根据区域网范围从全球SRTM数据中获取对应的DEM数据。
步骤3.2,利用外方位角元素计算每张影像的旋转矩阵,使用步骤1得到的角元素作为计算的初始值。
外方位角元素计算影像的旋转矩阵如下式所示:
Figure BDA0003226085700000121
步骤3.3,以单束同名光线为观测单元,分别计算各连接像点光线与DEM的交点,得到对应的物方坐标,计算公式如下:
Figure BDA0003226085700000122
Figure BDA0003226085700000123
式中,(Xg,Yg)为交点的物方坐标(即TCP坐标);zdem为交点的DEM高程;Xs、Ys、Zs是影像的外方位线元素,步骤1得到的(XTPC,YTPC,ZTPC)可作为外方位线元素的初始值;(x,y)是像点坐标;f是相机焦距。
步骤3.4,计算同名像点对应的物方点间的差异大小,如果大于角元素精度推导的理论物方定位精度α,则将同名像点标记为粗差点,不参与平差;否则,取同名像点对应的物方坐标的平均值作为真值,参与平差。
设成像系统的IMU角元素精度为ε,主光轴与重力方向的夹角为θ,则理论物方定位精度α为:
Figure BDA0003226085700000131
式中,H为摄影位置到地面的高度。
步骤3.5,逐点构建误差方程,解算未知数改正数,误差方程为:
Figure BDA0003226085700000132
式中,(x,y)是像点坐标,(x)、(y)是用未知数初值获取的像点坐标,vx、vy是误差变量,
Figure BDA0003226085700000133
表示求偏导,d代表误差改正,求解出外方位元素的改正数之后即可更新每张影像的外方位元素。
步骤3.6,重复步骤3.2至步骤3.5直到影像的每个角元素的变化不超过10-6弧度且所有像点改正数的平方和小于设定的阈值λ。
步骤4,整体最优正射影像生产。
步骤4.1,计算每张原始影像四个角点对应的地面坐标。
利用公式(18)和(19)计算每张原始影像四个角点对应的地面坐标构成四边形,并求得每个四边形的重心。
步骤4.2,选择正射影像对应的最佳原始影像进行反投影。
设正射影像上任一点P对应的地面坐标为(X,Y),选择与点P(X,Y)距离最近的四边形重心对应的原始影像作为最佳原始影像进行反投影。
步骤4.3,反投影法生产正射影像。
设正射影像上任一点P的坐标为(X′,Y′),由正射影像左下角图廓点地面坐标(X0,Y0)与正射影像比例尺分母M,计算P点所对应的地面坐标(X,Y)如下式所示:
X=X0+M×X′ (22)
Y=Y0+M×Y′ (23)
计算正射影像上点P对应原始影像上相应像点坐标p(x,y)为:
Figure BDA0003226085700000141
Figure BDA0003226085700000142
式中,X0、Y0为所有原始影像的角点对应的地面坐标中的最小值,Z为P点的高程,由DEM内插求得;ai,bi,ci(i=1,2,3)是影像的旋转矩阵元素,Xs、Ys、Zs是影像的线元素,f是原始影像的焦距。
由于所得的像点坐标不一定正好落在像元素中心,为此需进行灰度内插,本发明采用双线性内插的方法求得像点p的灰度值g(x,y),将像点的灰度值赋给正射影像上对应的像素点P(X′,Y′),即:
G(X′,Y′)=g(x,y) (26)
依次对每个像元进行上述运算,即能获得正射影像。
采用步骤4.1-步骤4.3生产整体正射影像,每个像素仅计算一次,避免了冗余计算,大幅提高了处理效率,由于每个像素都选择了最优原始影像,最终得到的整体正射影像也是最优化的。
具体实施时,以上流程可采用计算机软件技术实现自动运行流程。
本文中所描述的具体实施例仅仅是对本发明精神作举例说明。本发明所属技术领域的技术人员可以对所描述的具体实施例做各种各样的修改或补充或采用类似的方式替代,但并不会偏离本发明的精神或者超越所附权利要求书所定义的范围。

Claims (8)

1.一种针对国产面阵摆扫成像系统获取的影像生产正射影像的方法,其特征在于,包括如下步骤:
步骤1,建立局部坐标系作为统一的基准坐标系,将导航数据信息转换到基准坐标系中;
步骤2,在概略DEM辅助下进行影像同名点自动提取;
步骤2.1,在概略DEM辅助下,根据摆扫成像周期性和先后顺序获取影像间的相邻关系;
步骤2.2,根据步骤2.1获取的影像相邻关系,采用基于两类膨胀进行立体影像密集匹配的方法开展同名点匹配;
步骤3,基于DEM辅助进行区域网平差,将步骤1转换处理后的导航数据信息作为平差初始值,解算每张影像的外方位元素参数;
步骤4,整体最优正射影像生产;
步骤4.1,计算每张原始影像四个角点对应的地面坐标;
步骤4.2,选择正射影像对应的最佳原始影像进行反投影;
步骤4.3,反投影法生产正射影像。
2.如权利要求1所述的一种针对国产面阵摆扫成像系统获取的影像生产正射影像的方法,其特征在于:所述步骤1中局部坐标系采用定点切平面坐标系统,即TPC坐标系,坐标定义如下:在WGS84坐标框架下,选择飞行区域中几何中心点经度L0,纬度B0作为平面原点,以地球椭球高为0作椭球切平面,定义过原点纬度线的东方向为X坐标轴,过原点经度线北方向为Y坐标轴,过原点背向地心的垂线为Z坐标轴;TPC坐标系与像空间坐标系之间存在转换关系如下:
Figure FDA0003226085690000011
即从像空间坐标系→载荷本体坐标系→IMU本体坐标系→导航坐标系→地心地固坐标系→TPC坐标系;
各坐标系间旋转矩阵计算公式如下:
Figure FDA0003226085690000012
Figure FDA0003226085690000021
Figure FDA0003226085690000022
Figure FDA0003226085690000023
Figure FDA0003226085690000024
式中,
Figure FDA0003226085690000025
为TPC坐标系下的影像姿态旋转矩阵,
Figure FDA0003226085690000026
是影像的三个外方位角元素,为待求解参数;
Figure FDA0003226085690000027
是地心地固坐标系与TPC坐标系之间的转换矩阵,由TPC坐标原点的经度L0,纬度B0定义;
Figure FDA0003226085690000028
是导航坐标系与地心地固坐标系之间的转换矩阵,由载荷所在位置的经度Lz、纬度Bz定义;
Figure FDA0003226085690000029
是导航坐标系和IMU本体坐标系之间的旋转矩阵,(Φ,Θ,Ψ)是IMU记录的俯仰、横滚、偏航数值;
Figure FDA00032260856900000210
是载荷本体坐标系与IMU本体坐标系之间的转换矩阵,又称设备安置矩阵,由IMU安置在载荷上的ex、ey、ez三个角定义;
Figure FDA00032260856900000211
是像空间坐标系与载荷本体坐标系之间的转换矩阵;
求得旋转矩阵
Figure FDA0003226085690000031
后,反解可得到影像在TPC坐标系中的姿态角,即:
Figure FDA0003226085690000032
Figure FDA0003226085690000033
ω=arcsin(-b3) (9)
Figure FDA0003226085690000034
GPS点位置信息可以使用地球坐标定义模型无损地转换到任意需要的坐标系中,从原始的GPS坐标(L,B,H)转换到局部切平面坐标系中的TPC坐标经过如下变换:
首先将原始的GPS数据转换到地心坐标系中:
Xc=(N+H)cos B cos L (11)
Yc=(N+H)cos B sin L (12)
Zc=[N(1-e2)+H]sin B (13)
Figure FDA0003226085690000035
Figure FDA0003226085690000036
式中,N为卯酉圈的半径,a为地球椭球的长半轴,b为地球椭球的短半轴,(Xc,Yc,Zc)为地心坐标系的坐标;
然后将局部切平面坐标系原点坐标(L0,B0,0)带入式(11)~(13),求得局部切平面坐标系原点在地心坐标系中的坐标,记为
Figure FDA0003226085690000037
从地心系坐标到TPC坐标系的计算公式为:
Figure FDA0003226085690000038
式中,
Figure FDA0003226085690000041
是地心地固坐标系与TPC坐标系之间的转换矩阵,由式(2)计算可得;(XTPC,YTPC,ZTPC)为TPC坐标系的坐标;
Figure FDA0003226085690000042
即为转换到基准坐标系下的导航信息。
3.如权利要求1所述的一种针对国产面阵摆扫成像系统获取的影像生产正射影像的方法,其特征在于:所述步骤2中根据摆扫成像周期性和先后顺序获取影像间的相邻关系是先将摆扫一个周期得到的N张影像分为一组,每张影像根据其在摆扫周期中的次序将它们编号为1至N,顺序连续的影像必为邻接影像;各组之间的相邻关系根据GPS和IMU信息,将获取的影像投影到全球DEM表示的地面上,然后根据每张影像的地面范围求影像重叠情况,或根据重叠区域反算出原始影像的重叠区域,以此确认各组间影像的相邻关系。
4.如权利要求2所述的一种针对国产面阵摆扫成像系统获取的影像生产正射影像的方法,其特征在于:所述步骤3基于DEM辅助进行区域网平差包括以下几个步骤:
步骤3.1,利用所有影像的POS信息计算区域网范围,根据区域网范围从全球SRTM数据中获取对应的DEM数据;
步骤3.2,利用外方位角元素计算每张影像的旋转矩阵,使用步骤1得到的角元素作为计算的初始值;
步骤3.3,以单束同名光线为观测单元,分别计算各连接像点光线与DEM的交点,得到对应的物方坐标;
步骤3.4,计算同名像点对应的物方点间的差异大小,如果大于角元素精度推导的理论物方定位精度α,则将同名像点标记为粗差点,不参与平差;否则,取同名像点对应的物方坐标的平均值作为真值,参与平差;
步骤3.5,逐点构建误差方程,解算外方位元素;
步骤3.6,重复步骤3.2至步骤3.5直到影像的每个角元素的变化不超过10-6弧度且所有像点改正数的平方和小于设定的阈值λ。
5.如权利要求4所述的一种针对国产面阵摆扫成像系统获取的影像生产正射影像的方法,其特征在于:所述步骤3.3中以单束光线为观测单元,分别计算各连接点与DEM的交点,得到对应的物方坐标计算公式如下:
Figure FDA0003226085690000051
Figure FDA0003226085690000052
式中,(Xg,Yg)为交点的物方坐标,即TCP坐标;zdem为交点的DEM高程;Xs、Ys、Zs是影像的外方位线元素,步骤1得到的(XTPC,YTPC,ZTPC)可作为外方位线元素的初始值;(x,y)是像点坐标;f是相机焦距。
6.如权利要求5所述的一种针对国产面阵摆扫成像系统获取的影像生产正射影像的方法,其特征在于:所述步骤3.4中理论物方定位精度α的计算方式如下:
Figure FDA0003226085690000053
式中,H为摄影位置到地面的高度,ε为成像系统的IMU角元素精度,θ为主光轴与重力方向的夹角。
7.如权利要求6所述的一种针对国产面阵摆扫成像系统获取的影像生产正射影像的方法,其特征在于:所述步骤3.5中逐点构建误差方程为:
Figure FDA0003226085690000054
式中,(x,y)是像点坐标,(x)、(y)是用未知数初值获取的像点坐标,vx、vy是误差变量,
Figure FDA0003226085690000055
表示求偏导,d代表误差改正,求解出外方位元素的改正数之后即可更新每张影像的外方位元素。
8.如权利要求1所述的一种针对国产面阵摆扫成像系统获取的影像生产正射影像的方法,其特征在于:所述步骤4.3中反投影法生产正射影像的方法如下:设正射影像上任一点P的坐标为(X′,Y′),由正射影像左下角图廓点地面坐标(X0,Y0)与正射影像比例尺分母M,计算P点所对应的地面坐标(X,Y)如下式所示:
X=X0+M×X′ (22)
Y=Y0+M×Y′ (23)
计算正射影像上点P对应原始影像上相应像点坐标p(x,y)为:
Figure FDA0003226085690000061
Figure FDA0003226085690000062
式中,X0、Y0为所有原始影像的角点对应的地面坐标中的最小值,Z为P点的高程,由DEM内插求得;ai,bi,ci(i=1,2,3)是影像的旋转矩阵元素,Xs、Ys、Zs是影像的线元素,f是原始影像的焦距;
由于所得的像点坐标不一定正好落在像元素中心,为此需进行灰度内插,采用双线性内插的方法求得像点p的灰度值g(x,y),将像点的灰度值赋给正射影像上对应的像素点P(X′,Y′),即:
G(X′,Y′)=g(x,y) (26)
依次对每个像元进行上述运算,即能获得正射影像。
CN202110971773.4A 2021-08-24 2021-08-24 针对国产面阵摆扫成像系统获取的影像生产正射影像的方法 Active CN113739767B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202110971773.4A CN113739767B (zh) 2021-08-24 2021-08-24 针对国产面阵摆扫成像系统获取的影像生产正射影像的方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202110971773.4A CN113739767B (zh) 2021-08-24 2021-08-24 针对国产面阵摆扫成像系统获取的影像生产正射影像的方法

Publications (2)

Publication Number Publication Date
CN113739767A true CN113739767A (zh) 2021-12-03
CN113739767B CN113739767B (zh) 2022-09-13

Family

ID=78732418

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202110971773.4A Active CN113739767B (zh) 2021-08-24 2021-08-24 针对国产面阵摆扫成像系统获取的影像生产正射影像的方法

Country Status (1)

Country Link
CN (1) CN113739767B (zh)

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN114485574A (zh) * 2021-12-21 2022-05-13 武汉大学 基于卡尔曼滤波模型的三线阵影像pos辅助对地定位方法
CN114972536A (zh) * 2022-05-26 2022-08-30 中国人民解放军战略支援部队信息工程大学 一种航空面阵摆扫式相机定位和标定方法

Citations (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
EP0385705A2 (en) * 1989-02-27 1990-09-05 Texas Instruments Incorporated Apparatus and method for digitized 3D video system
WO2008075335A1 (en) * 2006-12-20 2008-06-26 Elbit Systems Electro-Optics Elop Ltd. Airborne photogrammetric imaging system and method
US20100246971A1 (en) * 2007-10-30 2010-09-30 Pasco Corporation House change determining method, house change determining program, house change determining image generating method, and house change determining image
CN105352509A (zh) * 2015-10-27 2016-02-24 武汉大学 地理信息时空约束下的无人机运动目标跟踪与定位方法
CN106780729A (zh) * 2016-11-10 2017-05-31 中国人民解放军理工大学 一种无人机序列影像批处理三维重建方法
CN109903352A (zh) * 2018-12-24 2019-06-18 中国科学院遥感与数字地球研究所 一种卫星遥感影像大区域无缝正射影像制作方法
CN110763205A (zh) * 2019-11-05 2020-02-07 新疆维吾尔自治区测绘科学研究院 数字摄影测量系统生成边境狭长区域正射影像图的方法

Patent Citations (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
EP0385705A2 (en) * 1989-02-27 1990-09-05 Texas Instruments Incorporated Apparatus and method for digitized 3D video system
WO2008075335A1 (en) * 2006-12-20 2008-06-26 Elbit Systems Electro-Optics Elop Ltd. Airborne photogrammetric imaging system and method
US20100246971A1 (en) * 2007-10-30 2010-09-30 Pasco Corporation House change determining method, house change determining program, house change determining image generating method, and house change determining image
CN105352509A (zh) * 2015-10-27 2016-02-24 武汉大学 地理信息时空约束下的无人机运动目标跟踪与定位方法
CN106780729A (zh) * 2016-11-10 2017-05-31 中国人民解放军理工大学 一种无人机序列影像批处理三维重建方法
CN109903352A (zh) * 2018-12-24 2019-06-18 中国科学院遥感与数字地球研究所 一种卫星遥感影像大区域无缝正射影像制作方法
CN110763205A (zh) * 2019-11-05 2020-02-07 新疆维吾尔自治区测绘科学研究院 数字摄影测量系统生成边境狭长区域正射影像图的方法

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
杨国鹏 等: "面阵摆扫航空相机序列图像的大区域无缝拼接", 《测绘科学》 *
莫德林: "航空线阵摆扫式相机成像仿真", 《光学学报》 *

Cited By (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN114485574A (zh) * 2021-12-21 2022-05-13 武汉大学 基于卡尔曼滤波模型的三线阵影像pos辅助对地定位方法
CN114972536A (zh) * 2022-05-26 2022-08-30 中国人民解放军战略支援部队信息工程大学 一种航空面阵摆扫式相机定位和标定方法
CN114972536B (zh) * 2022-05-26 2023-05-09 中国人民解放军战略支援部队信息工程大学 一种航空面阵摆扫式相机定位和标定方法

Also Published As

Publication number Publication date
CN113739767B (zh) 2022-09-13

Similar Documents

Publication Publication Date Title
CN109903352B (zh) 一种卫星遥感影像大区域无缝正射影像制作方法
CN110675450B (zh) 基于slam技术的正射影像实时生成方法及系统
CN110648398B (zh) 基于无人机航摄数据的正射影像实时生成方法及系统
CN111458720B (zh) 复杂山区基于机载激光雷达数据的倾斜摄影建模方法
Li Potential of high-resolution satellite imagery for national mapping products
CN108305237B (zh) 考虑不同光照成像条件的多立体影像融合制图方法
CN102506824B (zh) 一种城市低空无人机系统生成数字正射影像图的方法
US5596494A (en) Method and apparatus for acquiring digital maps
CN113192193B (zh) 基于Cesium三维地球框架的高压输电线路走廊三维重建方法
CN113739767B (zh) 针对国产面阵摆扫成像系统获取的影像生产正射影像的方法
CN112465732A (zh) 一种车载激光点云与序列全景影像的配准方法
CN110706273B (zh) 一种基于无人机的实时塌方区域面积的测量方法
CN112598608A (zh) 一种基于目标区域的光学卫星快速融合产品制作方法
CN108253942B (zh) 一种提高倾斜摄影测量空三质量的方法
CN110986888A (zh) 一种航空摄影一体化方法
CN108801225A (zh) 一种无人机倾斜影像定位方法、系统、介质及设备
CN115471619A (zh) 基于立体成像高分辨率卫星影像的城市三维模型构建方法
CN110780313A (zh) 一种无人机可见光立体测量采集建模方法
CN110853140A (zh) 一种dem辅助的光学视频卫星影像稳像方法
CN107705272A (zh) 一种空间影像的高精度几何校正方法
CN116824079A (zh) 基于全信息摄影测量的三维实体模型构建方法和装置
Jacobsen Orientation of high resolution optical space images
CN113483739B (zh) 海上目标位置测量方法
CN114777745A (zh) 一种基于无迹卡尔曼滤波的倾斜取证建模方法
CN111044076B (zh) 基于参考底图的高分一号b卫星几何检校方法

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant