CN114743113A - 光学影像辅助的星载立体sar影像控制点自动生成方法及系统 - Google Patents

光学影像辅助的星载立体sar影像控制点自动生成方法及系统 Download PDF

Info

Publication number
CN114743113A
CN114743113A CN202210439287.2A CN202210439287A CN114743113A CN 114743113 A CN114743113 A CN 114743113A CN 202210439287 A CN202210439287 A CN 202210439287A CN 114743113 A CN114743113 A CN 114743113A
Authority
CN
China
Prior art keywords
image
sar
point
lamp post
coordinates
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.)
Withdrawn
Application number
CN202210439287.2A
Other languages
English (en)
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.)
Anhui University
Original Assignee
Anhui University
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 Anhui University filed Critical Anhui University
Priority to CN202210439287.2A priority Critical patent/CN114743113A/zh
Publication of CN114743113A publication Critical patent/CN114743113A/zh
Withdrawn legal-status Critical Current

Links

Images

Classifications

    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S13/00—Systems using the reflection or reradiation of radio waves, e.g. radar systems; Analogous systems using reflection or reradiation of waves whose nature or wavelength is irrelevant or unspecified
    • G01S13/88—Radar or analogous systems specially adapted for specific applications
    • G01S13/89—Radar or analogous systems specially adapted for specific applications for mapping or imaging
    • G01S13/90—Radar or analogous systems specially adapted for specific applications for mapping or imaging using synthetic aperture techniques, e.g. synthetic aperture radar [SAR] techniques
    • G01S13/9021—SAR image post-processing techniques
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F18/00—Pattern recognition
    • G06F18/20—Analysing
    • G06F18/22—Matching criteria, e.g. proximity measures
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F18/00—Pattern recognition
    • G06F18/20—Analysing
    • G06F18/23—Clustering techniques
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T5/00—Image enhancement or restoration
    • G06T5/73—Deblurring; Sharpening
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00—Image analysis
    • G06T7/70—Determining position or orientation of objects or cameras
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00—Indexing scheme for image analysis or image enhancement
    • G06T2207/10—Image acquisition modality
    • G06T2207/10032—Satellite or aerial image; Remote sensing

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Data Mining & Analysis (AREA)
  • Remote Sensing (AREA)
  • General Physics & Mathematics (AREA)
  • Radar, Positioning & Navigation (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • Bioinformatics & Cheminformatics (AREA)
  • Evolutionary Computation (AREA)
  • General Engineering & Computer Science (AREA)
  • Evolutionary Biology (AREA)
  • Bioinformatics & Computational Biology (AREA)
  • Artificial Intelligence (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Electromagnetism (AREA)
  • Computer Networks & Wireless Communication (AREA)
  • Image Processing (AREA)

Abstract

本发明涉及一种光学影像辅助的星载立体SAR影像控制点生成方法及系统，与现有技术相比实现了从不同升降轨立体SAR影像中自动化的完成灯柱这一类控制点的匹配、提取、误差剔除，并可用的控制点。本发明包括以下步骤：基于高分光学影像上提取灯柱PS点；SAR影像匹配灯柱PS点；SAR影像灯柱PS点坐标精确提取；SAR影像灯柱PS点异常值剔除；SAR影像灯柱PS点坐标误差补偿；立体平差处理生成灯柱控制点三维坐标。本发明实现立体SAR影像同名PS点自动化匹配、提取、误差剔除，生成高精度可用的控制点，为基础测绘以及其他卫星几何校正提供控制数据。

Description

光学影像辅助的星载立体SAR影像控制点自动生成方法及 系统
技术领域
本发明属于星载合成孔径雷达影像地面应用领域，具体来说是基于光学影像辅助的星载立体SAR影像自动生成控制点方法及系统。
背景技术
星载合成孔径雷达SAR(Synthetic Aperture Radar)由于不受云雨影响，可全天时、全天候工作，已成为一种不可或缺的空间对地观测手段。星载SAR 相对光学卫星几何定位不受姿态的影响，全球无控几何定位精度较好。国外TerraSAR卫星经几何定标后，绝对定位精度优于0.3m,国产SAR卫星经几何定标后，绝对定位精度可优于1.5m。因此，可应用星载立体SAR生成高精度控制点，服务于基础测绘。
针对利用星载立体SAR生成地面控制点研究，目前研究主要集中在地面控制点选择和SAR影像匹配两方面。在控制点选取方面，有研究将立体SAR 技术应用于城市区域中少量天然永久散射体PS(Permanent Scatter)点，这些 PS点主要来自建筑物的外墙立面上的窗角或路灯的底部，所有参与立体平差的同名点影像坐标均是通过人工在立体SAR影像中提取。另外也有研究尝试自动化提取和匹配有限数量的建筑物外墙立面上的窗角点，但只适用于相同升降轨且入射角差异不大的立体SAR影像。在SAR影像匹配方面，有研究提出SAR-SIFT方法，该方法通过修改SIFT(Scale-invariant Feature Transform)算法以考虑斑点的统计特性，从而在SAR影像中有效提取局部描述符，但是该方法仅适用于从相同升降轨的SAR影像并且要求入射角差异较小。另外，有研究提出基于Harris的局部特征来检测相同的散射体点，通过使用外部数字高程模型DEM(Digital Elevation Model)和轨道信息对影像进行地理编码来限制搜索空间，并最终使用SIFT进行匹配，尽管该方法在不同升降轨的立体SAR 影像上有一定适用性，但仅适用于孤立的PS点。
现有研究中，从不同观测几何的立体SAR影像中识别和提取同名PS点，目前多为手动提取，应发展一种自动化的提取过程。另外采用不同升降轨的 SAR立体影像配置可获得大的基高比，从而提升立体SAR精度，但要求能够稳健的从不同升降轨配置的立体SAR影像中检测和识别同名PS点。国内还暂无对立体SAR影像自动生成高精度控制点的论文与专利。
发明内容
本发明的目的是为了解决现有技术中无法对具有不同升降轨的立体SAR 影像自动化处理生成可用的控制点，提供一种基于光学影像辅助的星载立体 SAR影像自动生成控制点方法及系统来解决上述问题。
为了实现上述目的，本发明技术方案如下：
一种基于光学影像辅助的星载立体SAR影像自动生成控制点方法，包括以下步骤：
11)基于高分光学影像上提取灯柱点：获取立体SAR影像覆盖区域所对应的高分光学遥感影像I0，基于I0检测到灯柱PS点后并得到其三维坐标，利用 SAR影像的距离多普勒几何定位模型将灯柱的三维坐标反算至SAR影像像方，完成立体SAR影像上同名灯柱PS点的近似像点坐标(τ0,t0)提取。
12)SAR影像匹配灯柱PS点：对SAR影像进行阈值化处理，得到二值图，将光学影像中的路灯点坐标利用SAR的几何定位模型反算至其像方，应用迭代最近点法实现灯柱点集与SAR阈值处理后的亮点点集匹配，从而实现SAR 影像灯柱PS点的匹配，并得到灯柱点对应SAR影像上的粗略坐标。
13)SAR影像灯柱PS点坐标精确提取：以近似坐标(x0,pi,SAR,y0,pi,SAR)为中心截取SAR影像块，进行点目标分析，精确提取灯柱pi点的方位和距离坐标 (xpi,SAR,ySAR)。
14)SAR影像灯柱PS点异常值剔除：潜在灯柱PS点候选对象SCR应足够高，但由于存在误差，进一步通过计算灯柱PI点的相位噪声σφ，进行异常值剔除。
15)SAR影像灯柱PS点坐标误差补偿：对提取的灯柱点pi影像坐标 (τpi,SAR,tpi,SAR)进行系统定位误差和大气延迟误差的校正，得到改正后的时间量纲坐标(τpi,cor,SAR,tpi,cor,SAR)。
16)立体平差处理生成灯柱控制点三维坐标：针对待生成控制点区域内可用的N张(N＞2)立体SAR影像，按照上述步骤1-4提取同名灯柱点对集合，然后进行立体平差处理，生成灯柱控制点的三维坐标。
所述基于高分光学影像上提取灯柱点包括以下步骤：
21)获取对应区域的光学影像I0后，为了使灯柱的阴影更加突出，对光学影像I0进行高频提升滤波处理：
I＝I0+aIm
式中I0为原始图像，a是锐化因子标量，Im是使用非锐化蒙板后的图像，为I0与其对应模糊图像之差，其中a值越高，表示锐化程度越高。
22)在锐化后的光学图像I中提取任意包含灯柱阴影的模板图像T，然后将模板与图像I进行相关运算，得到图像I在像方(u,v)坐标处的相似度：
Figure BDA0003610753240000041
式中I(x,y)和T(x,y)分别表示图像I和模板图像T在(x,y)坐标处的像素值， N1,N2表示模板图像T的大小，
Figure BDA0003610753240000042
和
Figure BDA0003610753240000043
分别表示图像I和模板图像T的平均像素值，ρ(u,v)为归一化互相关值，ρ(u,v)不受图像亮度或对比度的变化影响，因此可提升模板匹配的效果。
23)光学图像I中经模板匹配后得到每个灯柱对象存在多个像素，进一步对模板提取到的像素进行聚类处理，得到表示灯柱唯一的像素点pi，处理如下：
Figure BDA0003610753240000044
式中pi表示计算m(pi)漂移向量的空间坐标点，pj代表pi附近的点，g为具有带宽h的高斯核函数，而|| ||表示欧几里得距离算子。
24)根据光学图像I上灯柱聚类点pi的像素坐标(xpi,I,ypi,I)，获取其平面地理坐标(Latpi,I,Lonpi,I)，并利用平面坐标在开源SRTM DEM数据上内插出对应高程Hpi,I，构成完整的三维坐标(Latpi,I,Lonpi,I,Hpi,I)。
所述SAR影像匹配灯柱PS点包括以下步骤：
31)根据立体SAR影像的元数据信息构建距离多普勒RD几何定位方程：
fR(CS(t),CT,τ)＝|CS(t)-CT|-c/2·τ＝0
Figure BDA0003610753240000045
式中fR、fA分别表示距离方程和多普勒方程，τ、t分别表示距离向时间和方位向时间，CS(t)和
Figure BDA0003610753240000046
分别表示成像t时刻SAR天线相位中心在WGS84坐标系下的位置矢量和速度矢量，CT为观测目标在WGS84坐标系下的位置矢量，fdc(τ)表示在距离双程延时为τ处的多普中心频率；
32)逐个将灯柱聚类点pi的三维坐标(Latpi,I,Lonpi,I,Hpi,I)利用RD定位几何方程反投影至SAR影像的像方，得到灯柱点集。
33)根据一定阈值，对SAR影像上亮点进行二值掩膜处理。
34)利用迭代最近点法，实现反算在SAR像方的灯柱点集与SAR影像上亮点集匹配，SAR影像上只保留匹配到路灯点pi的亮点，从而完成SAR影像灯柱点的匹配，并得到粗略匹配坐标(x0,pi,SAR,y0,pi,SAR)。
所述SAR影像灯柱PS点坐标精确提取包括以下步骤：
41)以近似坐标(x0,pi,SAR,y0,pi,SAR)为中心在SAR影像上截取32*32的影像块，并寻找最大峰值位置；
42)以最大峰值处为中心截取32*32样本窗做二维离散傅里叶变换，32倍采样因子频域补零后进行二维逆离散傅里叶变换，得到升采样后的点目标影像块；
43)对过采样后的影像块，以最大峰值处截取3×3影像子块，进行二维抛物面插值，，椭圆抛物面插值模型如下：
f(x,y)＝a5x2+a4y2+a3xy+a2x+a1y+a0
式中x,y3×3影像子块中的像素坐标，为ai(1≤i≤6)为抛物面插值模型f(x,y)的参数，将3×3子块像素值代入上式中并通过最小二乘求解ai。参数ai求解后，进而获得灯柱pi点精确坐标(xpi,SAR,ySAR)。
所述SAR影像灯柱PS点异常值剔除包括以下步骤：
51)计算灯柱pi点的峰值功率Ppeak和杂波功率Pclutter。峰值功率Ppeak是以峰值(xpi,SAR,ySAR)为中心为2×2分辨率像元的主瓣区域的积分信号能量，杂波功率Pclutter是主瓣和旁瓣区域外的积分信号能量。
52)计算灯柱pi点的信杂比SCR并得到相位噪声
Figure BDA0003610753240000061
计算公式如下：
Figure BDA0003610753240000062
Figure BDA0003610753240000063
53)如果灯柱pi点计算出的σφ大于0.5rad，则予以剔除，否则保留。并将灯柱影像pi在SAR影像上的(xpi,SAR,ySAR)坐标转为时间量纲的坐标 (τpi,SAR,tpi,SAR)。转化公式如下：
t＝tstart+y/PRF
τ＝τmin+x/Fs
式中x、y为待转换的像素坐标，τ、t为转换后的时间量纲的坐标，分别表示距离向双程时延和方位向时间，PRF为方位向脉冲重复频率，Fs为距离向采样频率，τmin为SAR影像近距端双程时延，tstart为方位向起始时间。
所述SAR影像灯柱PS点坐标误差补偿包括以下步骤：
61)从SAR影像元数据中读取系统误差定标参数(Δτcal、Δtcal)，补偿至灯柱pi点影像坐标中，公式如下：
τpi,1＝τpi-Δτcal
tpi,cor＝tpi-Δtcal
上式中τpi,1和tpi,1分别为对灯柱pi点影像初始时间量纲坐标(τpi,SAR,tpi,SAR)进行系统误差补偿后的距离向和方位向时间坐标。
62)进一步对距离向时间τpi,1完成大气对流层延迟误差ΔτTro以及电离层的延迟ΔτIno误差改正，公式如下：
τpi,cor＝τpi,1-ΔτTro-ΔτIno
经上述处理完后，灯柱pi点得到误差补偿后的时间量纲坐标为(τpi,cor,SAR,tpi,cor,SAR)。
所述立体平差处理生成灯柱控制点三维坐标包括以下步骤：
71)按照第一步到第五步，逐个对N张SAR影像提取同名灯柱点对，针对某一灯柱点pi得到同名点对集合
Figure BDA0003610753240000071
2≤m≤N，m表示同名灯柱点pi所在SAR影像的标识，
Figure BDA0003610753240000072
和tpi，cor，SARi分别表示在影像SARi上经过误差改正的距离向和方位向时间坐标。
72)对可用的立体SAR影像构建RD几何定位方程；
73)对某一灯柱点pi的同名点对，构建立体观测方程如下：
Figure BDA0003610753240000074
式中CT＝[Latpi,Lonpi,Hpi]T为待求的三维坐标，其初始值
Figure BDA0003610753240000075
74)对立体观测方程线性化，得到误差方程如下：
Figure BDA0003610753240000076
上式可简记为v＝AdCT-l，。
75)求解立体观测误差方程，得到待求柱点pi点的三维坐标改正值 [dLatpi,dLonpi,dHpi]T。
76)将坐标改正值dCT＝[dLatpi,dLonpi,dHpi]T补偿至初始值
Figure BDA0003610753240000081
中，更新
Figure BDA0003610753240000082
77)重复72)-76)，直到改正值小于阈值10-5为止；
78)重复72)-77)，逐个计算提取出的灯柱点，直到所有灯柱点计算完为止。
基于光学影像辅助的星载立体SAR控制点自动生成系统，包括以下模块：
光学影像灯柱点提取模板，用于对光学影像滤波处理，选择影像上灯柱完成模板匹配和聚类处理，完成灯柱点的平面坐标提取以及从SRTM DEM上提取对应高程。
立体SAR影像灯柱点自动匹配模块，用于首先将提取的灯柱点近似三维坐标反算至SAR影像上，然后利用点目标分析获取灯柱点精确的SAR影像坐标，并计算信杂比及相位噪声，完成异常点的提出。
SAR影像坐标误差补偿模块，用于对匹配的灯柱点SAR影像坐标完成各类误差的补偿，包括系统误差、对流程延迟误差、电离层延迟误差。
立体SAR影像灯柱控制点坐标生成模块，用于构建SAR影像几何定位模型，完成灯柱点立体平差处理，获取三维坐标。
有益效果
本发明的基于光学影像辅助的星载立体SAR影像自动生成控制点方法及系统，与现有技术相比可实现不同升降轨SAR影像中灯柱PS点自动匹配，子像素级精确提取灯柱点在SAR影像上的坐标并完成异常值剔除，立体平差生成大量高精度可用的控制点，本方法可高效率、自动化的完成立体SAR影像上灯柱这一类PS点的匹配、提取，提取出的控制点可更好的服务基础测绘。
附图说明
图1为本发明方法顺序图；
图2为本发明所涉及的方法实施流程图；
图3为光学影像计算灯柱模板相关系数示意图；
图4为SAR影像灯柱PS点匹配示意图；
图5为SAR影像灯柱PS点坐标精确提取示意图。
具体实施方式
为使对本发明的结构特征及所达成的功效有更进一步的了解与认识，用以较佳的实施例及附图配合详细的说明，说明如下：
如图1和图2所示，本发明所述的一种基于光学影像辅助的星载立体SAR 影像自动生成控制点方法，包括以下步骤：
第一步，基于高分光学影像上提取灯柱点：获取立体SAR影像覆盖区域所对应的高分光学遥感影像I0，基于I0检测到灯柱PS点后并得到其三维坐标，利用SAR影像的距离多普勒几何定位模型将灯柱的三维坐标反算至SAR影像像方，提取立体SAR影像上同名灯柱PS点的近似像点坐标(τ0,t0)。其具体步骤如下：
(1)获取对应区域的光学影像I0后，为了使灯柱的阴影更加突出，对光学影像I0进行高频提升滤波处理：
I＝I0+aIm
式中I0为原始图像，a是锐化因子标量，Im是使用非锐化蒙板后的图像，为I0与其对应模糊图像之差，其中a值越高，表示锐化程度越高。
(2)如图3所示，在锐化后的光学图像I中提取任意包含灯柱阴影的模板图像T，然后将模板与图像I进行相关运算，得到图像I在像方(u,v)坐标处的相似度：
Figure BDA0003610753240000101
式中I(x,y)和T(x,y)分别表示图像I和模板图像T在(x,y)坐标处的像素值，N1,N2表示模板图像T的大小，
Figure BDA0003610753240000102
和
Figure BDA0003610753240000103
分别表示图像I和模板图像T的平均像素值，ρ(u,v)为归一化互相关值，ρ(u,v)不受图像亮度或对比度的变化影响，因此可提升模板匹配的效果。
(3)光学图像I中经模板匹配后得到每个灯柱对象存在多个像素，进一步对模板提取到的像素进行聚类处理，得到表示灯柱唯一的像素点pi，处理如下：
Figure BDA0003610753240000104
式中pi表示计算m(pi)漂移向量的空间坐标点，pj代表pi附近的点，g为具有带宽h的高斯核函数，而|| ||表示欧几里得距离算子。
(4)根据光学图像I上灯柱聚类点pi的像素坐标(xpi,I,ypi,I)，获取其平面地理坐标(Latpi,I,Lonpi,I)，并利用平面坐标在开源SRTM DEM数据上内插出对应高程Hpi,I，构成完整的三维坐标(Latpi,I,Lonpi,I,Hpi,I)。
第二步，SAR影像匹配灯柱PS点：对SAR影像进行阈值化处理，得到二值图，将光学影像中的路灯点坐标利用SAR的几何定位模型反算至其像方，应用迭代最近点法实现灯柱点集与SAR阈值处理后的亮点点集匹配，从而实现SAR影像灯柱PS点的匹配，并得到灯柱点对应SAR影像上的粗略坐标。
其具体步骤如下：
(1)根据立体SAR影像的元数据信息构建距离多普勒RD几何定位方程：
fR(CS(t),CT,τ)＝|CS(t)-CT|-c/2·τ＝0
Figure BDA0003610753240000111
式中fR、fA分别表示距离方程和多普勒方程，τ、t分别表示距离向时间和方位向时间，CS(t)和
Figure BDA0003610753240000112
分别表示成像t时刻SAR天线相位中心在WGS84坐标系下的位置矢量和速度矢量，CT为观测目标在WGS84坐标系下的位置矢量，fdc(τ)表示在距离双程延时为τ处的多普中心频率；
(2)逐个将灯柱聚类点pi的三维坐标(Latpi,I,Lonpi,I,Hpi,I)利用RD定位几何方程反投影至SAR影像的像方，得到灯柱点集。
(3)根据一定阈值，对SAR影像上亮点进行二值掩膜处理。
(4)如图4所示，利用迭代最近点法，实现反算在SAR像方的灯柱点集与SAR影像上亮点集匹配，SAR影像上只保留匹配到路灯点pi的亮点，从而完成SAR影像灯柱点的匹配，并得到粗略匹配坐标(x0,pi,SAR,y0,pi,SAR)。
第三步，SAR影像灯柱PS点坐标精确提取：如图5所示，以近似坐标 (x0,pi,SAR,y0,pi,SAR)为中心截取SAR影像块，进行点目标分析，精确提取灯柱pi点的方位和距离坐标(xpi,SAR,ySAR)。其具体步骤如下：
(1)以近似坐标(x0,pi,SAR,y0,pi,SAR)为中心在SAR影像上截取32*32的影像块，并寻找最大峰值位置；
(2)以最大峰值处为中心截取32*32样本窗做二维离散傅里叶变换，32 倍采样因子频域补零后进行二维逆离散傅里叶变换，得到升采样后的点目标影像块；
(3)对过采样后的影像块，以最大峰值处截取3×3影像子块，进行二维抛物面插值，，椭圆抛物面插值模型如下：
f(x,y)＝a5x2+a4y2+a3xy+a2x+a1y+a0
式中x,y3×3影像子块中的像素坐标，为ai(1≤i≤6)为抛物面插值模型f(x,y)的参数，将3×3子块像素值代入上式中并通过最小二乘求解ai。参数ai求解后，进而获得灯柱pi点精确坐标(xpi,SAR,ySAR)。
第四步，SAR影像灯柱PS点异常值剔除：潜在灯柱PS点候选对象SCR 应足够高，但由于存在误差，进一步通过计算灯柱PI点的相位噪声σφ，进行异常值剔除。其具体步骤如下：
(1)计算灯柱pi点的峰值功率Ppeak和杂波功率Pclutter。峰值功率Ppeak是以峰值(xpi,SAR,ySAR)为中心为2×2分辨率像元的主瓣区域的积分信号能量，杂波功率 Pclutter是主瓣和旁瓣区域外的积分信号能量。
(2)计算灯柱pi点的信杂比SCR并得到相位噪声
Figure BDA0003610753240000121
计算公式如下：
Figure BDA0003610753240000122
Figure BDA0003610753240000123
(3)如果灯柱pi点计算出的σφ大于0.5rad，则予以剔除，否则保留。并将灯柱影像pi在SAR影像上的(xpi,SAR,ySAR)坐标转为时间量纲的坐标 (τpi,SAR,tpi,SAR)。转化公式如下：
t＝tstart+y/PRF
τ＝τmin+x/Fs
式中x、y为待转换的像素坐标，τ、t为转换后的时间量纲的坐标，分别表示距离向双程时延和方位向时间，PRF为方位向脉冲重复频率，Fs为距离向采样频率，τmin为SAR影像近距端双程时延，tstart为方位向起始时间。
第五步，SAR影像灯柱PS点坐标误差补偿：对提取的灯柱点pi影像坐标 (τpi,SAR,tpi,SAR)进行系统定位误差和大气延迟误差的校正，得到改正后的时间量纲坐标(τpi,cor,SAR,tpi,cor,SAR)。其具体步骤如下：
(1)从SAR影像元数据中读取系统误差定标参数(Δτcal、Δtcal)，补偿至灯柱pi点影像坐标中，公式如下：
τpi,1＝τpi-Δτcal
tpi,cor＝tpi-Δtcal
上式中τpi,1和tpi,1分别为对灯柱pi点影像初始时间量纲坐标(τpi,SAR,tpi,SAR)进行系统误差补偿后的距离向和方位向时间坐标。
(2)进一步对距离向时间τpi,1完成大气对流层延迟误差ΔτTro以及电离层的延迟ΔτIno误差改正，公式如下：
τpi,cor＝τpi,1-ΔτTro-ΔτIno
经上述处理完后，灯柱pi点得到误差补偿后的时间量纲坐标为 (τpi,cor,SAR,tpi,cor,SAR)。
第六步，立体平差处理生成灯柱控制点三维坐标：针对待生成控制点区域内可用的N张(N＞2)立体SAR影像，按照上述步骤1-4提取同名灯柱点对集合，然后进行立体平差处理，生成灯柱控制点的三维坐标。其具体步骤如下：
(1)按照第一步到第五步，逐个对N张SAR影像提取同名灯柱点对，针对某一灯柱点pi得到同名点对集合
Figure BDA0003610753240000131
2≤m≤N，m表示同名灯柱点pi所在SAR影像的标识，
Figure BDA0003610753240000132
和tpi，cor，SARi分别表示在影像SARi上经过误差改正的距离向和方位向时间坐标。
(2)对可用的立体SAR影像构建RD几何定位方程；
(3)对某一灯柱点pi的同名点对，构建立体观测方程如下：
Figure BDA0003610753240000141
式中CT＝[Latpi,Lonpi,Hpi]T为待求的三维坐标，其初始值
Figure BDA0003610753240000142
(4)对立体观测方程线性化，得到误差方程如下：
Figure BDA0003610753240000143
上式可简记为v＝AdCT-l，。
(5)求解立体观测误差方程，得到待求柱点pi点的三维坐标改正值 [dLatpi,dLonpi,dHpi]T。
(6)将坐标改正值dCT＝[dLatpi,dLonpi,dHpi]T补偿至初始值
Figure BDA0003610753240000144
中，更新
Figure BDA0003610753240000145
(7)重复(2)-(6)，直到改正值小于阈值10-5为止；
(8)重复(2)-(7)，逐个计算提取出的灯柱点，直到所有灯柱点计算完为止。
下面以国产GF-3SAR卫星为例对本发明提出的方法进行说明：选取覆盖嵩山区域5景滑聚模式数据，数据中涵盖升降轨、左右视，不同组合立体角范围为70°到110°。光学影像采用谷歌地图嵩山区域最高分辨率图层的影像，其分辨率约为0.5m。为了使光学影像中灯柱更加突出及后续提高模板匹配精度，首先对光学影像进行滤波处理，选取一类灯柱作为图像模板进行模板匹配后，对匹配出的像素做聚类处理，得到45个灯柱点，通过SRTM DEM获取45个灯柱点的近似三维坐标。通过5景滑聚模式影像的几何定位模型，将45个灯柱点反算至SAR影像像方，并进一步提取精确的像方坐标，再经误差补偿后，进行立体平差处理，解算出灯柱点的三维坐标，经外业GPS测量检查点验证，平面精度优于1m，高程精度优于0.5m。
与现有技术相比，本发明具有如下优点与有益效果：
(1)基于光学影像辅助自动提取灯柱这一类PS点，而不是所有PS点，进而保证了后续控制点的可用性；
(2)可自动从不同升降轨的立体SAR影像中准确匹配灯柱PS点，实现异常值的剔除，从而保证控制点的精度。并且生成的控制点可为基础测绘、其他卫星几何校正提供较好的控制数据补充。
以上显示和描述了本发明的基本原理、主要特征和本发明的优点。本行业的技术人员应该了解，本发明不受上述实施例的限制，上述实施例和说明书中描述的只是本发明的原理，在不脱离本发明精神和范围的前提下本发明还会有各种变化和改进，这些变化和改进都落入要求保护的本发明的范围内。本发明要求的保护范围由所附的权利要求书及其等同物界定。

Claims (8)

1.一种光学影像辅助的星载立体SAR影像控制点自动生成方法，其特征在于，包括以下步骤：
11)基于高分光学影像上提取灯柱点在SAR影像上近似坐标：获取立体SAR影像覆盖区域所对应的高分光学遥感影像I0，基于I0检测到灯柱点后并得到其近似三维坐标，利用SAR影像的距离多普勒几何定位模型将灯柱的近似三维坐标反算至SAR影像像方，完成立体SAR影像上同名灯柱点的近似像点坐标(τ0,t0)提取。
12)SAR影像匹配灯柱PS点：对SAR影像进行阈值化处理，得到二值图，将光学影像中的路灯点坐标利用SAR的几何定位模型反算至其像方，应用迭代最近点法实现灯柱点集与SAR阈值处理后的亮点点集匹配，从而实现SAR影像灯柱PS点的匹配，并得到灯柱点对应SAR影像上的粗略坐标。
13)SAR影像灯柱PS点坐标精确提取：以近似坐标(x0,pi,SAR,y0,pi,SAR)为中心截取SAR影像块，进行点目标分析，精确提取灯柱pi点的方位和距离坐标(xpi,SAR,ySAR)。
14)SAR影像灯柱PS点异常值剔除：潜在灯柱PS点候选对象SCR应足够高，但由于存在误差，进一步通过计算灯柱pi点的相位噪声σφ，进行异常值剔除。
15)SAR影像灯柱PS点坐标误差补偿：对提取的灯柱点pi影像坐标(τpi,SAR,tpi,SAR)进行系统定位误差和大气延迟误差的校正，得到改正后的时间量纲坐标(τpi,cor,SAR,tpi,cor,SAR)。
16)立体平差处理生成灯柱控制点三维坐标：针对待生成控制点区域内可用的N张(N＞2)立体SAR影像，按照上述步骤11)-15)提取同名灯柱点对集合，然后进行立体平差处理，生成灯柱控制点的三维坐标。
2.根据权利要求1所述的一种光学影像辅助的星载立体SAR影像控制点自动生成方法，其特征在于，所述基于高分光学影像上提取灯柱点包括以下步骤：
21)获取对应区域的光学影像I0后，为了使灯柱的阴影更加突出，对光学影像I0进行高频提升滤波处理：
I＝I0+aIm，
式中，I0为原始图像，a是锐化因子标量，Im是使用非锐化蒙板后的图像，为I0与其对应模糊图像之差，其中a值越高，表示锐化程度越高。
22)在锐化后的光学图像I中提取任意包含灯柱阴影的模板图像T，模板大小本发明设置为256×256，然后将模板与图像I进行相关运算，得到图像I在像方(u,v)坐标处的相似度：
Figure RE-FDA0003663144480000021
式中，I(x,y)和T(x,y)分别表示图像I和模板图像T在(x,y)坐标处的像素值，N1,N2表示模板图像T的大小，
Figure RE-FDA0003663144480000022
和
Figure RE-FDA0003663144480000023
分别表示图像I和模板图像T的平均像素值，ρ(u,v)为归一化互相关值，ρ(u,v)不受图像亮度或对比度的变化影响，因此可提升模板匹配的效果。
23)光学图像I中经模板匹配后得到每个灯柱对象存在多个像素，进一步对模板提取到的像素进行聚类处理，得到表示灯柱唯一的像素点pi，处理如下：
Figure RE-FDA0003663144480000031
式中，pi表示计算m(pi)漂移向量的空间坐标点，pj代表pi附近的点，g为具有带宽h的高斯核函数，而|| ||表示欧几里得距离算子。
24)根据光学图像I上灯柱聚类点pi的像素坐标(xpi,I,ypi,I)，获取其平面地理坐标(Latpi,I,Lonpi,I)，并利用平面坐标在开源SRTM DEM数据上内插出对应高程Hpi,I，构成完整的三维坐标(Latpi,I,Lonpi,I,Hpi,I)。
3.根据权利要求1所述的一种光学影像辅助的星载立体SAR影像控制点自动生成方法，其特征在于，所述SAR影像匹配灯柱点包括以下步骤：
31)根据立体SAR影像的元数据信息构建距离多普勒RD几何定位方程：
Figure RE-FDA0003663144480000032
式中，fR、fA分别表示距离方程和多普勒方程，τ、t分别表示距离向时间和方位向时间，CS(t)和
Figure RE-FDA0003663144480000033
分别表示成像t时刻SAR天线相位中心在WGS84坐标系下的位置矢量和速度矢量，CT为观测目标在WGS84坐标系下的位置矢量，fdc(τ)表示在距离双程延时为τ处的多普中心频率；
32)逐个将灯柱聚类点pi的三维坐标(Latpi,I,Lonpi,I,Hpi,I)利用RD定位几何方程反投影至SAR影像的像方，得到灯柱点集。
33)根据一定阈值，本发明该阈值由对灯柱点像素值统计后获取，对SAR影像上亮点进行二值掩膜处理。
34)利用迭代最近点法，实现反算在SAR像方的灯柱点集与SAR影像上亮点集匹配，SAR影像上只保留匹配到路灯点pi的亮点，从而完成SAR影像灯柱点的匹配，并得到粗略匹配坐标(x0,pi,SAR,y0,pi,SAR)。
4.根据权利要求1所述的一种光学影像辅助的星载立体SAR影像控制点自动生成方法，其特征在于，所述SAR影像灯柱点坐标精确提取包括以下步骤：
41)以近似坐标(x0,pi,SAR,y0,pi,SAR)为中心在SAR影像上截取32*32的影像块，并寻找最大峰值位置；
42)以最大峰值处为中心截取32*32样本窗做二维离散傅里叶变换，32倍采样因子频域补零后进行二维逆离散傅里叶变换，得到升采样后的点目标影像块；
43)对过采样后的影像块，以最大峰值处截取3×3影像子块，进行二维抛物面插值，椭圆抛物面插值模型如下：
f(x,y)＝a5x2+a4y2+a3xy+a2x+a1y+a0，
式中，x,y为在3×3影像子块中的像素坐标，为ai(1≤i≤6)为抛物面插值模型f(x,y)的参数，将3×3子块像素值代入上式中并通过最小二乘求解ai。参数ai求解后，进而获得灯柱pi点精确坐标(xpi,SAR,ySAR)。
5.根据权利要求1所述的一种光学影像辅助的星载立体SAR影像控制点自动生成方法，其特征在于，所述SAR影像灯柱点异常值剔除包括以下步骤：
51)计算灯柱pi点的峰值功率Ppeak和杂波功率Pclutter。峰值功率Ppeak是以峰值(xpi,SAR,ySAR)为中心为2×2分辨率像元的主瓣区域的积分信号能量，杂波功率Pclutter是主瓣和旁瓣区域外的积分信号能量。
52)计算灯柱pi点的信杂比SCR并得到相位噪声
Figure RE-FDA0003663144480000041
计算公式如下：
Figure RE-FDA0003663144480000042
Figure RE-FDA0003663144480000051
53)如果灯柱pi点计算出的σφ大于0.5rad，则予以剔除，否则保留。并将灯柱影像pi在SAR影像上的(xpi,SAR,ySAR)坐标转为时间量纲的坐标(τpi,SAR,tpi,SAR)。转化公式如下：
Figure RE-FDA0003663144480000052
式中，x、y为待转换的像素坐标，τ、t为转换后的时间量纲的坐标，分别表示距离向双程时延和方位向时间，PRF为方位向脉冲重复频率，Fs为距离向采样频率，τmin为SAR影像近距端双程时延，tstart为方位向起始时间。
6.根据权利要求1所述的一种光学影像辅助的星载立体SAR影像控制点自动生成方法，其特征在于，所述SAR影像灯柱点坐标误差补偿包括以下步骤：
61)从SAR影像元数据中读取系统误差定标参数(Δτcal、Δtcal)，补偿至灯柱pi点影像坐标中，公式如下：
Figure RE-FDA0003663144480000053
式中，τpi,1和tpi,1分别为对灯柱pi点影像初始时间量纲坐标(τpi,SAR,tpi,SAR)进行系统误差补偿后的距离向和方位向时间坐标。
62)进一步对距离向时间τpi,1完成大气对流层延迟误差ΔτTro以及电离层的延迟ΔτIno误差改正，公式如下：
τpi,cor＝τpi,1-ΔτTro-ΔτIno，
经上述处理完后，灯柱pi点得到误差补偿后的时间量纲坐标为(τpi,cor,SAR,tpi,cor,SAR)。
7.根据权利要求1所述的一种光学影像辅助的星载立体SAR影像控制点自动生成方法，其特征在于，所述立体平差处理生成灯柱控制点三维坐标包括以下步骤：
71)按照步骤11)-15)，逐个对N张SAR影像提取同名灯柱点对，针对某一灯柱点pi得到同名点对集合
Figure RE-FDA0003663144480000061
m表示同名灯柱点pi所在SAR影像的标识，
Figure RE-FDA0003663144480000062
和tpi,cor,SARi分别表示在影像SARi上经过误差改正的距离向和方位向时间坐标。
72)对可用的立体SAR影像构建距离多普勒几何定位方程；
73)对某一灯柱点pi的同名点对，构建立体观测方程如下：
Figure RE-FDA0003663144480000064
式中，CT＝[Latpi,Lonpi,Hpi]T为待求的三维坐标，其初始值
Figure RE-FDA0003663144480000065
74)对立体观测方程线性化，得到误差方程如下：
Figure RE-FDA0003663144480000066
上式可简记为v＝AdCT-l，。
75)求解立体观测误差方程，得到待求柱点pi点的三维坐标改正值[dLatpi,dLonpi,dHpi]T。
76)将坐标改正值dCT＝[dLatpi,dLonpi,dHpi]T补偿至初始值
Figure RE-FDA0003663144480000072
中，更新
Figure RE-FDA0003663144480000073
77)重复72)-76)，直到改正值小于阈值10-5为止；
78)重复72)-77)，逐个计算提取出的灯柱点，直到所有灯柱点坐标计算完为止。
8.根据权利要求1所述的一种光学影像辅助的星载立体SAR影像控制点自动生成系统，包括以下模块：
光学影像灯柱点提取模板，用于对光学影像滤波处理，选择影像上灯柱完成模板匹配和聚类处理，完成灯柱点的平面坐标提取以及从SRTM DEM上提取对应高程。
立体SAR影像灯柱点自动匹配模块，用于首先将提取的灯柱点近似三维坐标反算至SAR影像上，然后利用点目标分析获取灯柱点精确的SAR影像坐标，并计算信杂比及相位噪声，完成异常点的提出。
SAR影像坐标误差补偿模块，用于对匹配的灯柱点SAR影像坐标完成各类误差的补偿，包括系统误差、对流程延迟误差、电离层延迟误差。
立体SAR影像灯柱控制点坐标生成模块，用于构建SAR影像几何定位模型，完成灯柱点立体平差处理，获取控制点三维坐标。
CN202210439287.2A 2022-04-22 2022-04-22 光学影像辅助的星载立体sar影像控制点自动生成方法及系统 Withdrawn CN114743113A (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202210439287.2A CN114743113A (zh) 2022-04-22 2022-04-22 光学影像辅助的星载立体sar影像控制点自动生成方法及系统

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202210439287.2A CN114743113A (zh) 2022-04-22 2022-04-22 光学影像辅助的星载立体sar影像控制点自动生成方法及系统

Publications (1)

Publication Number Publication Date
CN114743113A true CN114743113A (zh) 2022-07-12

Family

ID=82283741

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202210439287.2A Withdrawn CN114743113A (zh) 2022-04-22 2022-04-22 光学影像辅助的星载立体sar影像控制点自动生成方法及系统

Country Status (1)

Country Link
CN (1) CN114743113A (zh)

Cited By (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN117849795A (zh) * 2024-01-15 2024-04-09 安徽大学 一种顾及成像特性误差补偿的扫描成像模式sar的几何处理方法

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN117849795A (zh) * 2024-01-15 2024-04-09 安徽大学 一种顾及成像特性误差补偿的扫描成像模式sar的几何处理方法
CN117849795B (zh) * 2024-01-15 2024-08-02 安徽大学 一种顾及成像特性误差补偿的扫描成像模式sar的几何处理方法

Similar Documents

Publication Publication Date Title
CN104574347B (zh) 基于多源遥感数据的在轨卫星图像几何定位精度评价方法
WO2019042232A1 (zh) 一种快速鲁棒的多模态遥感影像匹配方法和系统
CN116203562A (zh) 一种光学影像辅助的星载立体sar影像控制点自动生成方法
CN112419380B (zh) 一种基于云掩膜的静止轨道卫星序列影像高精度配准方法
CN109212522B (zh) 一种基于双基星载sar的高精度dem反演方法及装置
CN109523585B (zh) 一种基于方向相位一致性的多源遥感影像特征匹配方法
CN109270527B (zh) 圆迹sar子孔径图像序列联合相关dem提取方法
CN103093459A (zh) 利用机载LiDAR点云数据辅助影像匹配的方法
CN103428408A (zh) 一种适用于帧间的图像稳像方法
CN114743113A (zh) 光学影像辅助的星载立体sar影像控制点自动生成方法及系统
CN115015931A (zh) 无需外部误差校正的实时差分立体sar几何定位方法及系统
CN107463944A (zh) 一种利用多时相高分辨率sar图像的道路信息提取方法
CN112505647B (zh) 一种基于序贯子图像序列的动目标方位速度估计方法
CN108090898B (zh) 基于字典表示的卫星遥感图像典型地标检测方法
CN111127334A (zh) 基于rd平面像素映射的sar图像实时几何校正方法及系统
CN105093222A (zh) 一种sar影像区域网平差连接点自动提取方法
CN115601278A (zh) 基于子图像配准的高精度运动误差补偿的方法
Oh et al. Automated RPCs Bias Compensation for KOMPSAT Imagery Using Orthoimage GCP Chips in Korea
CN113298713A (zh) 一种抗云干扰的在轨快速配准方法
CN117761695B (zh) 基于自适应分区sift的多角度sar三维成像方法
CN110967693B (zh) 一种稳健高效的快速分解投影自动聚焦方法及系统
Toth et al. Recovery of sensor platform trajectory from LiDAR data using reference surfaces
Leitloff et al. Inferring traffic activity from optical satellite images
CN117495933B (zh) 基于视差修正的光电望远镜外挂镜头图像实时配准方法
CN118566927A (zh) 基于光学阴影特征点位预测的sar影像灯杆提取方法及系统

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
WW01 Invention patent application withdrawn after publication

Application publication date: 20220712

WW01 Invention patent application withdrawn after publication