CN114723794A - 基于lsd直线检测的sar图像配准方法 - Google Patents
基于lsd直线检测的sar图像配准方法 Download PDFInfo
- Publication number
- CN114723794A CN114723794A CN202210398024.1A CN202210398024A CN114723794A CN 114723794 A CN114723794 A CN 114723794A CN 202210398024 A CN202210398024 A CN 202210398024A CN 114723794 A CN114723794 A CN 114723794A
- Authority
- CN
- China
- Prior art keywords
- sar
- image
- sate
- mid
- calculating
- 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
Links
- 238000000034 method Methods 0.000 title claims abstract description 41
- 238000001514 detection method Methods 0.000 title claims abstract description 25
- 230000009466 transformation Effects 0.000 claims abstract description 29
- 239000011159 matrix material Substances 0.000 claims abstract description 22
- 238000012545 processing Methods 0.000 claims abstract description 21
- 238000013519 translation Methods 0.000 claims abstract description 16
- 238000007781 pre-processing Methods 0.000 claims abstract description 10
- 238000001914 filtration Methods 0.000 claims description 37
- 238000003708 edge detection Methods 0.000 claims description 9
- 230000000877 morphologic effect Effects 0.000 claims description 9
- 230000001629 suppression Effects 0.000 claims description 8
- 238000012216 screening Methods 0.000 claims description 6
- 238000000605 extraction Methods 0.000 claims description 4
- 238000010606 normalization Methods 0.000 claims description 3
- 230000008859 change Effects 0.000 abstract description 3
- 238000005286 illumination Methods 0.000 abstract description 2
- 230000000007 visual effect Effects 0.000 abstract description 2
- 238000000354 decomposition reaction Methods 0.000 abstract 1
- 230000014616 translation Effects 0.000 description 9
- 238000003384 imaging method Methods 0.000 description 4
- 238000013459 approach Methods 0.000 description 3
- 238000004364 calculation method Methods 0.000 description 3
- 238000004458 analytical method Methods 0.000 description 2
- 238000010586 diagram Methods 0.000 description 2
- 238000005516 engineering process Methods 0.000 description 2
- 230000008569 process Effects 0.000 description 2
- 238000011160 research Methods 0.000 description 2
- 230000009286 beneficial effect Effects 0.000 description 1
- 230000007547 defect Effects 0.000 description 1
- 238000011161 development Methods 0.000 description 1
- 230000000694 effects Effects 0.000 description 1
- 230000006872 improvement Effects 0.000 description 1
- 238000005259 measurement Methods 0.000 description 1
- 230000003287 optical effect Effects 0.000 description 1
- 230000010287 polarization Effects 0.000 description 1
- 230000005855 radiation Effects 0.000 description 1
- 230000003252 repetitive effect Effects 0.000 description 1
- 230000035945 sensitivity Effects 0.000 description 1
- 230000003313 weakening effect Effects 0.000 description 1
Images
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T7/00—Image analysis
- G06T7/30—Determination of transform parameters for the alignment of images, i.e. image registration
- G06T7/33—Determination of transform parameters for the alignment of images, i.e. image registration using feature-based methods
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T5/00—Image enhancement or restoration
- G06T5/70—Denoising; Smoothing
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T7/00—Image analysis
- G06T7/10—Segmentation; Edge detection
- G06T7/13—Edge detection
-
- 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
- G06T2207/10044—Radar image
Landscapes
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- General Physics & Mathematics (AREA)
- Theoretical Computer Science (AREA)
- Computer Vision & Pattern Recognition (AREA)
- Radar Systems Or Details Thereof (AREA)
- Image Processing (AREA)
Abstract
本发明公开了基于LSD（Line Segment Detector）直线检测的SAR图像配准方法，包括：读取SAR图像与卫星数字地图，构建待配准图像对；对待配准图像对进行预处理；采用LSD直线检测方法检测预处理图像结果中的直线特征，以点数据的形式进行存储；计算出角度参数、尺度参数、变换矩阵和平移参数，然后通过尺度变换、平移变换以及变换矩阵的处理，得到SAR图像中的目标点在卫星数字地图中对应的匹配点；载入卫星数字地图配套的经纬度数据库，根据匹配点坐标获取目标点经纬度信息。本发明通过使用基于线特征对SAR图像与卫星数字地图进行配准工作，提升了配准的准确度，对视角和光照变化具有很好的鲁棒性。
Description
技术领域
本发明涉及SAR图像配准技术领域，具体涉及基于LSD直线检测的SAR图像配准方法。
背景技术
合成孔径雷达(Synthetic Aperture Radar，SAR)成像具有全天时、全天候、可穿透性等优点，以及它在不同频段、不同极化下可得到目标的高分辨率雷达图像，为人们提供非常有用的目标信息。因此，一直以来针对SAR图像配准的研究热度经久不息，大量研究成果已经被广泛应用于经济和科技等众多领域，有着广泛的应用前景和发展潜力。
图像配准技术是将相同地区，在不同时刻、不同视角、不同传感器或不同光照条件下拍摄的图像进行空间对准的过程。由于采集位置、空间分辨率、采集时间等不同，使得获取的同一场影的光学与SAR影像间存在较大的辐射与几何差异。因此，要综合利用不同遥感影像必须使它们实现几何意义上的对准，即需要进行图像配准。
图像配准方法一般分为两种：一种是基于灰度的方法，另一种是基于特征的方法。基于灰度的方法使用图像或者图像块的灰度信息，通过相似性度量算法，测度两幅对应图像或对应图像块之间的相似性，搜索参数空间，寻找相似性最优的变换模型的参数值，从而实现图像配准，但是因为其全局优化的计算复杂度太高是不太适用于对时间要求比较严格的配准的。基于特征的方法用对图像中具有稳定性和重复性图像特征的分析，代替对整个图像的分析，降低计算量以及配准算法对传感器变化的敏感性。
对比基于灰度的方法，基于特征的方法因为具有更好的重复性和稳定性，所以它具有更强大的区分能力，对局部有差异的图像对有更好的鲁棒性。而且基于特征的方法使用少量图像特征，对相同图像来说，算法的运算性能普遍高于基于灰度的方法。同时，因为SAR图像分辨率高、数据量大、受噪声影响严重的问题，现有的配准算法处理速度很慢，所以为了高效性和高精度的合成孔径雷达图像的配准，对原有的配准算法改进使其具有高效性和高精度特性是一种比较普遍存在的方式。因为基于特征方法相较于基于灰度方法具有计算量低的特征，所以基于特征的方法是当前SAR图像使用最多的配准方法。但是，由于合成孔径雷达成像方式与可见光成像方式的不同，造成两种不同图像在点特征上的差异，即二者之间能检测出来的角点并不完全是一一对应的，因此基于点特征配准方法降低了匹配精度。
发明内容
发明目的：本发明目的在于针对现有技术的不足，提供一种基于LSD(LineSegment Detector)直线检测的SAR图像配准方法，与有重叠区域的可见光图像进行比对，减弱或者去除因为传感器或者天气变换等原因导致的图像变化的影响，获取更加精准的目标信息。
技术方案：本发明所述基于LSD直线检测的SAR图像配准方法，包括如下步骤：
S1：读取SAR图像与卫星数字地图，构建待配准图像对Msar、Msate；
S2：对待配准图像对Msar、Msate进行预处理，得到预处理后的图像Sarcanny、Satecanny；
S3：采用LSD直线检测方法检测S2预处理后的图像Sarcanny、Satecanny中的直线特征，保存检测到的直线信息，以点数据的形式进行存储，并表示为Vsar、Vsate；
S4：根据Vsar、Vsate计算出角度参数theta、尺度参数s、变换矩阵transMatric和平移参数，然后对图像Msar进行尺度变换、平移变换和变换矩阵处理，得到处理后图像中的目标点在卫星数字地图Msate中对应的匹配点，表示为MP(x，y)；
S5：载入卫星数字地图Msate配套的经纬度数据库，根据匹配点坐标MP(x，y)获取目标点经纬度信息，表示为Lon、Lat。
进一步完善上述技术方案，所述S1中SAR图像与卫星数字地图以多维矩阵的格式进行存储。
进一步地，所述S2中对Msate进行的预处理操作包括灰度化、滤波、二值化、形态学运算、高斯滤波、Canny边缘检测，对Msar进行的预处理操作包括灰度化、滤波、二值化、形态学运算、相干斑抑制、高斯滤波、Canny边缘检测。
进一步地，所述预处理操作过程如下：
S21：将待配准图像对进行灰度化处理
设r、g、b分别表示待处理图像的三个通道，采用平均值法进行灰度化处理，将三通道转换为单通道，得到图中(i，j)处的像素灰度值为：
grayVal_mean(i，j)＝(r(i，j)+g(i，j)+b(i，j))/3 (1)
r(i，j)、g(i，j)、b(i，j)表示三个通道r、g、b在图中(i，j)处的像素灰度值；
将待配准图像对Msar、Msate按公式(1)处理，得到SAR图像与卫星数字地图的灰度图；
S22：对灰度化的待配准图像对进行滤波处理，然后对滤波结果进行二值化处理，再根据二值化图像中非零像素个数的占比重复滤波及二值化操作，迭代至非零像素个数占比达到预设阈值，包括：
对S21中灰度图进行归一化处理，使像素灰度值在0～1范围内，设max[grayVal(i，j)]为图像中最大灰度值，min[grayVal(i，j)]为图像中最小灰度值，通过下式进行图像归一化：
设置卷积核kernel的大小，以核内中心位置对应的像素为目标像素，使用kernel对待配准图像对进行卷积操作，卷积结果为：
按从左到右、从上到下滑动卷积核，直到处理完图像最后一个像素，
对卷积结果进行二值化处理，得到二值化后的结果Sarthrehold、Satethrehold，凸显待配准图像中关键区域的轮廓边缘信息；
S23：对Sarthrehold、Satethrehold进行形态学处理，将断裂区域连接起来，得到处理结果Sardilate、Satedilate；
S24：根据SAR图像Msar的距离向分辨率Dre与方位向分辨率Dis构建滤波核kernelnoise，大小为(9/Dis)×(9/Dre)，以核中心点为原点，通过下式构建：
相干斑抑制步骤如下：
S241、采用构建的滤波核在二值SAR图像Sardilate区域进行滑动操作，从左到右，从上到下，
S242、判断滤波核中元素与对应二值SAR图像Sardilate中元素是否相等且等于255，记录符合条件的元素个数，记为Pixel255num，
S243、判断滤波核中元素值为255，且对应二值SAR图像Sardilate中元素值为0，记录符合条件的元素个数，记为Pixel0/255num，
S244、当Pixel255num与Pixel0/255num的值大于滤波核元素个数的1/3时，判断图像中每行左右邻域或每列上下邻域中像素灰度值为255的像素个数是否超过5，如果超过则设置目标像素的灰度值为0，
S245、重复S241～S244，直到滤波核中心位置与图像右下角像素重合，得到处理结果Sarcoh；
S25：对卫星数字地图Satedilate和相干斑抑制后的SAR图像Sarcoh进行高斯滤波去噪、Canny边缘检测，得到处理后的图像Sarcanny、Satecanny。
进一步地，所述S3采用LSD直线检测方法检测图像Sarcanny与Satecanny中的直线特征，包括：
S31：首先计算图像中所有像素点的梯度大小，记为grd，根据梯度大小grd计算出水平角度agl与梯度幅值|grd|，如下式所示：
agl＝arctan(grdi/-grdj) (6)
利用伪排序得到梯度幅值大的点作为种子点，以该点水平线角度agl作为区域的初始角度θr，在八邻域中寻找与θr的偏差小于预设阈值threshold的点，将该点加入到区域中并更新θr，得到线支持区域，最后通过矩形逼近得到直线特征，检测到的直线特征以点数据的形式进行存储，分别表示为Psar、Psate；
S32：筛选Psar、Psate中最具代表性的直线特征来进行配准工作，筛选步骤如下：
S321、计算Psar、Psate中所有直线段的长度，并按降序排列，表示为Lsar、Lsate，
S322、根据预设直线特征数量Lnum，在Lsar、Lsate中提取前Lnum条直线特征，
S323、提取距离间隔较大的直线特征，
S324、保存直线特征的端点坐标，表示为Vsar、Vsate。
进一步地，所述S4包括：
S41：将直线特征的端点坐标Vsar、Vsate转化为极坐标，表示为Pcsar(angle，solar)、Pcsate(angle，solar)，
S411、计算每条直线段的斜率k与截距b；
S412、计算直线段的中心点坐标(xmid，ymid)、kmid及bmid
kmid＝-1/k (9)
bmid＝ymid-kmid×xmid (10)
S413、计算极坐标的极角与极径：
当kmid＞0时，获得如下三个点：
O(x3，y3)＝O(xmid，ymid)
A(x4，y4)＝A[xmid-1，kmid×(xmid-1)+bmid]
B(x5，y5)＝B(xmid-1，ymid)
当kmid＜0时，获得如下三个点：
O(x3，y3)＝O(xmid，ymid)
A(x4，y4)＝A[xmid+1，kmid×(xmid+1)+bmid]
B(x5，y5)＝B(xmid-1，ymid)
得到极角angle：
dot＝(x3-x4)×(x3-x5)+(y3-y4)×(y3-y5)
angle＝[(a cos(dot/(m1×m2)))×180]/3.1415926
得到极径solar：
令
则有
S42：以卫星数字地图为基准，计算SAR图像中直线与卫星数字地图中直线的角度差theta.
计算Pcsar与Pcsate中angle的差值：angle_diff＝Pcsar(angle)-Pcsate(angle)，将相同的angle-diff值进行计数，出现次数最多的angle_diff值记为theta；
S43：计算尺度参数s，计算步骤如下：
S431、采用Pcsar中的angle统计出SAR图像中平行直线组；
S432、计算SAR图像与卫星数字地图中每组平行直线极坐标的极径差
solarsar＝abs(Pcsar(solarn)-Pcsar(solarn-1))
solarsate＝abs(Pcsate(solarn)-Pcsate(solarn-1))
scale＝solarsate/solarsar
将相同的scale值进行计数，将出现次数最多的scale值记为尺度参数s；
S44、计算变换矩阵transMatric与平移参数
S441、分别计算SAR图像与卫星数字地图中直线特征的交点，分别表示为nodesar、nodesate，
S442、分别计算交点集合nodesar与nodesate中各点间坐标的差分node_diff(xdiff，ydiff)，
S443、将node_diff中的横坐标差值与纵坐标差值分别与阈值比较，保存满足条件的node_diff，同时将相同的node_diff值进行计数；
S444、获取出现次数最多的node_diff中的xdiff、ydiff；
S445、定义与SAR图像尺寸相关参数imgSarx与imgSary，与卫星数字地图尺寸相关参数imgSatex与imgSatey
S446、计算变换矩阵
S447、计算平移参数
S45、结合尺度参数s、变换矩阵transMatric、平移参数计算出SAR图像在卫星数字地图中的匹配点MP(x，y)，设SAR图像中的目标点为C(x，y)，则得到匹配点MP(x，y)：
有益效果：与现有技术相比，本发明的优点在于：本发明通过使用基于线特征对SAR图像与卫星数字地图进行配准工作，提升了配准的准确度，相比于点特征配准方法，对视角和光照变化具有很好的鲁棒性，克服了因不同传感器成像方式差异带来的误差；在特征提取方面，本发明在效率上也高于现有的基于点特征匹配方法；与现有的基于Hough直线特征检测的配准技术相比，本发明具有结果准确、误检可控，同时不需要调节参数的优点。
附图说明
图1是本发明的方法流程图；
图2是本发明中滤波卷积核的示意图；
图3是本发明中待卷积图像区域；
图4是本发明中待配准图像对；
图5是本发明中配准后的效果图；
图6是本发明中SAR图像与卫星数字地图匹配结果；
图7是本发明SAR图像中已匹配目标的经纬度信息结果。
具体实施方式
下面通过附图对本发明技术方案进行详细说明，但是本发明的保护范围不局限于所述实施例。
本发明采用基于LSD(Line Segment Detector)直线检测的SAR图像配准方法：首先SAR图像与卫星数字地图对进行灰度化、滤波、二值化等预处理，然后进行相干斑抑制、高斯滤波、以及边缘检测，获取可以进行直线检测的图像样本，接着使用由直线特征计算得到的尺度参数、角度参数以及变换矩形获取SAR图像中目标点在卫星数字地图中的匹配点，最后根据卫星数字地图中匹配点坐标与卫星数字地图配套经纬度数据库计算SAR图像中目标点的地理位置信息。
如图1所示的基于LSD直线检测的SAR图像配准方法，包括如下步骤：
S1：读取SAR图像与卫星数字地图，使用多通道矩阵存储，构建待配准图像对，并表示为Msar、Msate；
S2：对待配准图像对进行预处理，主要包括灰度化、滤波、二值化、形态学运算，相干斑抑制，高斯滤波、Canny边缘检测；
S3：对S2中的图像结果使用LSD直线检测方法检测图像结果中的直线特征，保存从SAR图像中与卫星数字地图中检测到的直线信息，以点数据的形式进行存储，分别表示为Vsar、Vsate；
S4：根据上述Vsar和Vsate计算出尺度参数、角度参数及变换矩阵，分别表示为s，theta，transMatric，然后通过尺度变换、平移变换以及变换矩阵的处理，得到SAR图像中的目标点在卫星数字地图中对应的匹配点，表示为Mp(x，y)；
S5：载入经纬度数据库，根据匹配点坐标MP(x，y)计算目标点经纬度信息，表示为Lon、Lat。
具体地，S1中使用的SAR图像与卫星数字地图，以多维矩阵的格式进行存储，读取图像关键像素信息。
S2中对图像进行预处理的操作，包括灰度化、滤波、二值化、形态学运算，相干斑抑制(只有SAR图像需要)，高斯滤波、Canny边缘检测，具体地，包括：
S21：设r、g、b分别是SAR图像与卫星数字地图的三个通道，采用平均值法进行灰度化处理，将三通道转换为单通道，即将SAR图像与卫星数字地图的三个通道r、g、b在(i，j)处的像素灰度值r(i，j)、g(i，j)、b(i，j)相加，求平均值grayVal_mean(i，j)，得到在图中(i，j)处的像素灰度值为：
grayVal_mean(i，j)＝(r(i，j)+g(i，j)+b(i，j))/3 (1)
最终将待配准图像对按公式1处理，得到SAR图像与卫星数字地图的灰度图。
S22：对待配准图像进行滤波处理，提高对比度，突显图像边缘等特征，然后对滤波结果进行二值化处理，再根据二值化图像中非零像素个数的占比重复滤波与二值化操作，迭代至非零像素个数占比达到预设阈值，包括
对S21中灰度图进行归一化处理，使像素灰度值在0～1范围内，设max[grayVal(i，j)]为图像中最大灰度值，min[grayVal(i，j)]为图像中最小灰度值，通过下式进行图像归一化：
设卷积核kernel的大小为3×3(kernel的大小根据实际情况可以设置)，如图2所示。以核内中心位置对应的像素为目标像素，使用kernel对待配准图像对进行卷积操作，对应待卷积的图像区域如图3所示，卷积结果为：
按从左到右，从上到下滑动卷积核，直到处理完图像最后一个像素。
接下来，对卷积结果进行二值化处理，得到二值化后的结果Sarthrehold，Satethrehold凸显待配准图像中关键区域的轮廓等边缘信息。
S23：对Sarthrehold与Satethrehold进行形态学处理，将断裂区域连接起来，得到处理结果Sardilate、Satedilate；
S24：对二值化后的SAR图像Sardilate进行滤波处理、，抑制相干斑噪声
根据SAR图像的距离向分辨率与方位向分辨率(分别表示为Dre，Dis)构建滤波核kernelnoise，大小为(9/Dis)×(9/Dre)，以核中心点为原点，通过下式构建：
相干斑噪声抑制步骤如下：
1)采用构建的滤波核在SAR图像Sardilate区域进行滑动操作，从左到右，从上到下；
2)判断滤波核中元素与对应二值SAR图像元素是否相等且等于255，记录符合条件的元素个数，记为Pixel255num；
3)判断滤波核中元素值为255，且对应二值SAR图像元素值为0，记录符合条件的元素个数，记为Pixel0/255num；
4)当Pixel255hum与Pixel0/255num的值大于滤波核元素个数的1/3时，判断图像中每行左右邻域或每列上下邻域中像素灰度值为255的像素个数是否超过5，如果超过则设置目标像素的灰度值为0；
5)重复上述操作，直到滤波核中心位置与图像右下角像素重合，得到结果Sarcoh。
S25：对卫星数字地图Satedilate与相干斑抑制后的SAR图像Sarcoh进行高斯滤波去噪、Canny边缘检测，得到处理后的图像Sarcanny、Satecanny。
S3采用LSD直线检测方法能在线性的时间内得出亚像素级精度的检测结果，且无需调节任何参数。
S31：首先计算图像Sarcanny、Satecanny中所有像素点的梯度大小和水平角度，分别记为grd、agl，根据梯度大小计算出agl与梯度幅值|grd|，如下式所示：
agl＝arctan(grdi/-grdj) (6)
利用伪排序得到梯度幅值大的点作为种子点，以该点水平线角度agl作为区域的初始角度θr，在八邻域中寻找与θr的偏差小于预设阈值threshold的点，将该点加入到区域中并更新θr，得到线支持区域，最后通过矩形逼近得到直线特征，检测到的直线特征以点数据的形式进行存储，分别表示为Psar、Psate。
S32：筛选Psar、Psate中最具代表性的直线特征用来进行配准工作，Psar、Psate中存储了检测到的所有直线段，筛选步骤如下：
1)计算所有直线段的长度，并按降序排列，表示为Lsar、Lsate；
2)根据预设直线特征数量Lnum，在Lsar、Lsate中提取前Lnum条直线特征；
3)提取距离间隔较大的直线特征；
保存直线特征的端点坐标，表示为Vsar、Vsate。
S4中根据上述Vsar、Vsate计算出尺度参数、角度参数及变换矩阵，分别表示为s、theta、transMatric；然后通过尺度变换、平移以及变换矩阵的处理，得到SAR图像中的目标点在卫星数字地图中对应的匹配点，表示为MP(x，y)。
S41：将Vsar、Vsate中直线特征的端点坐标转化为极坐标，表示为PCsar(angle，solar)、Pcsate(angle，solar)
1)计算每条直线段的斜率k与截距b；
2)计算直线段的中心点坐标(xmid，ymid)、kmid及bmid
kmid＝-1/k (9)
bmid＝ymid-kmid×xmid (10)
计算极坐标的极角与极径：
当kmid＞0时，由2中可以获得如下三个点：
O(x3，y3)＝O(xmid，ymid)
A(x4，y4)＝A[xmid-1，kmid×(xmid-1)+bmid]
B(x5，y5)＝B(xmid-1，ymid)
当kmid＜0时，由2中可以获得如下三个点：
O(x3，y3)＝O(xmid，ymid)
A(x4，y4)＝A[xmid+1，kmid×(xmid+1)+bmid]
B(x5，y5)＝B(xmid-1，ymid)
可以得到极角angle：
dot＝(x3-x4)×(x3-x5)+(y3-y4)×(y3-y5)
angle＝[(a cos(dot/(m1×m2)))×180]/3.1415926
可以得到极径solar：
令
则有
S42：以卫星数字地图为基准，计算SAR图像中直线与卫星数字地图中直线的角度差theta。
根据在SAR图像与卫星数字地图中检测到的直线特征计算出SAR图像与卫星数字地图之间的旋转角度，S41中通过直线端点坐标将直线用极坐标，分别表示为Pcsar(angle，solar)与Pcsate(angle，solar)；
计算Pcsar与Pcsate中angle的差值：
angle_diff＝Pcsar(angle)-Pcsate(angle)
将相同的angle_diff值进行计数，出现次数最多的angle_diff值即为theta；
S43：计算尺度参数s
尺度参数计算步骤如下：
(1)采用Pcsar中的angle统计出SAR图像中平行直线组
(2)计算SAR图像与卫星数字地图中每组平行直线极坐标的极径差
solarsar＝abs(Pcsar(solarn)-Pcsar(solarn-1))
solarsate＝abs(Pcsate(solarn)-Pcsate(solarn-1))
scale＝solarsate/solarsar
将相同的scale值进行计数，出现次数最多的scale值即为尺度参数s；
使用s、theta值可以将SAR图像的尺度与朝向变换为与卫星数字地图一致，如图4、图5所示，其中图4为原图，图5为处理结果。
S44：计算变换矩阵transMatric与平移参数
(1)分别计算SAR图像与卫星数字地图中直线特征的交点，分别表示为nodesar与nodesate；
(2)分别计算交点集合nodesar与nodesate中各点间坐标的差分node_diff(xdiff，ydiff)；
(3)将node_diff中的横坐标差值与纵坐标差值分别与阈值比较，保存满足条件的node_diff，同时将相同的node_diff值进行计数；
(4)获取出现次数最多的node_diff中的xdiff与ydiff；
(5)定义与SAR图像尺寸相关参数imgSarx与imgSary，与卫星数字地图尺寸相关参数imgSatex与imgSatey；
(6)计算变换矩阵
则有：
(7)计算平移参数
S45：结合尺度参数s与变换矩阵transMatric计算出SAR图像在卫星数字地图中的匹配点MP(x，y)，设SAR图像中的目标点为C(x，y)，则可以得到匹配点MP(x，y)：
匹配结果如图6所不，其中A点为SAR图像中的目标点，B点为卫星数字地图中与A对应的匹配点。
S5：采用匹配点坐标计算目标点经纬度信息
采用载入的经纬度数据库，结合匹配点坐标MP(x，y)计算出目标点经纬度信息，如图7所示。
如上所述，尽管参照特定的优选实施例已经表示和表述了本发明，但其不得解释为对本发明自身的限制。在不脱离所附权利要求定义的本发明的精神和范围前提下，可对其在形式上和细节上作出各种变化。
Claims (6)
1.基于LSD直线检测的SAR图像配准方法，其特征在于，包括如下步骤：
S1：读取SAR图像与卫星数字地图，构建待配准图像对Msar、Msate；
S2：对待配准图像对Msar、Msate进行预处理，得到预处理后的图像Sarcanny、Satecanny；
S3：采用LSD直线检测方法检测S2预处理后的图像Sarcanny、Satecanny中的直线特征，保存检测到的直线信息，以点数据的形式进行存储，并表示为Vsar、Vsate；
S4：根据Vsar、Vsate计算出角度参数theta、尺度参数s、变换矩阵transMatric和平移参数，然后对图像Msar进行尺度变换、平移变换和变换矩阵处理，得到处理后图像中的目标点在卫星数字地图Msate中对应的匹配点，表示为MP(x，y)；
S5：载入卫星数字地图Msate配套的经纬度数据库，根据匹配点坐标MP(x，y)获取目标点经纬度信息，表示为Lon、Lat。
2.根据权利要求1所述的基于LSD直线检测的SAR图像配准方法，其特征在于：所述S1中SAR图像与卫星数字地图以多维矩阵的格式进行存储。
3.根据权利要求1所述的基于LSD直线检测的SAR图像配准方法，其特征在于：所述S2中对Msate进行的预处理操作包括灰度化、滤波、二值化、形态学运算、高斯滤波、Canny边缘检测，对Msar进行的预处理操作包括灰度化、滤波、二值化、形态学运算、相干斑抑制、高斯滤波、Canny边缘检测。
4.根据权利要求3所述的基于LSD直线检测的SAR图像配准方法，其特征在于，所述预处理操作过程如下：
S21：将待配准图像对进行灰度化处理
设r、g、b分别表示待处理图像的三个通道，采用平均值法进行灰度化处理，将三通道转换为单通道，得到图中(i，j)处的像素灰度值为：
grayVal_mean(i，j)＝(r(i，j)+g(i，j)+b(i，j))/3 (1)
r(i，j)、g(i，j)、b(i，j)表示三个通道r、g、b在图中(i，j)处的像素灰度值；
将待配准图像对Msar、Msate按公式(1)处理，得到SAR图像与卫星数字地图的灰度图；
S22：对灰度化的待配准图像对进行滤波处理，然后对滤波结果进行二值化处理，再根据二值化图像中非零像素个数的占比重复滤波及二值化操作，迭代至非零像素个数占比达到预设阈值，包括：
对S21中灰度图进行归一化处理，使像素灰度值在0～1范围内，设max[grayVal(i，j)]为图像中最大灰度值，min[grayVal(i，j)]为图像中最小灰度值，通过下式进行图像归一化：
设置卷积核kernel的大小，以核内中心位置对应的像素为目标像素，使用kernel对待配准图像对进行卷积操作，卷积结果为：
按从左到右、从上到下滑动卷积核，直到处理完图像最后一个像素，
对卷积结果进行二值化处理，得到二值化后的结果Sarthrehold、Satethrehold，凸显待配准图像中关键区域的轮廓边缘信息；
S23：对Sarthrehold、Satethrehold进行形态学处理，将断裂区域连接起来，得到处理结果Sardilate、Satedilate；
S24：根据SAR图像Msar的距离向分辨率Dre与方位向分辨率Dis构建滤波核kernelnoise，大小为(9/Dis)×(9/Dre)，以核中心点为原点，通过下式构建：
相干斑抑制步骤如下：
S241、采用构建的滤波核在二值SAR图像Sardilate区域进行滑动操作，从左到右、从上到下，
S242、判断滤波核中元素与对应二值SAR图像Sardilate中元素是否相等且等于255，记录符合条件的元素个数，记为Pixel255num，
S243、判断滤波核中元素值为255，且对应二值SAR图像Sardilate中元素值为0，记录符合条件的元素个数，记为Pixel0/255num，
S244、当Pixel255num与Pixel0/255num的值大于滤波核元素个数的1/3时，判断图像中每行左右邻域或每列上下邻域中像素灰度值为255的像素个数是否超过5，如果超过则设置目标像素的灰度值为0，
S245、重复S241～S244，直到滤波核中心位置与图像右下角像素重合，得到处理结果Sarcoh；
S25：对卫星数字地图Satedilate和相干斑抑制后的SAR图像Sarcoh进行高斯滤波去噪、Canny边缘检测，得到处理后的图像Sarcanny、Satecanny。
5.根据权利要求4所述的基于LSD直线检测的SAR图像配准方法，其特征在于，所述S3采用LSD直线检测方法检测图像Sarcanny与Satecanny中的直线特征，包括：
S31：首先计算图像中所有像素点的梯度大小，记为grd，根据梯度大小grd计算出水平角度agl与梯度幅值|grd|，如下式所示：
agl＝arctan(grdi/-grdj) (6)
其中，I(i，j)为图像I在第i行第j个像素上的像素值，
利用伪排序得到梯度幅值大的点作为种子点，以该点水平线角度agl作为区域的初始角度θr，在八邻域中寻找与θr的偏差小于预设阈值threshold的点，将该点加入到区域中并更新θr，得到线支持区域，最后通过矩形逼近得到直线特征，检测到的直线特征以点数据的形式进行存储，分别表示为Psar、Psate；
S32：筛选Psar、Psate中最具代表性的直线特征来进行配准工作，筛选步骤如下：
S321、计算Psar、Psate中所有直线段的长度，并按降序排列，表示为Lsar、Lsate，
S322、根据预设直线特征数量Lnum，在Lsar、Lsate中提取前Lnum条直线特征，
S323、提取距离间隔较大的直线特征，
S324、保存直线特征的端点坐标，表示为Vsar、Vsate。
6.根据权利要求5所述的基于LSD直线检测的SAR图像配准方法，其特征在于，所述S4包括：
S41：将直线特征的端点坐标Vsar、Vsate转化为极坐标，表示为Pcsar(angle，solar)、Pcsate(angle，solar)；
S411、计算每条直线段的斜率k与截距b；
S412、计算直线段的中心点坐标(xmid，ymid)、kmid及bmid
kmid＝-1/k (9)
bmid＝ymid-kmid×xmid (10)
S413、计算极坐标的极角与极径：
当kmid＞0时，获得如下三个点：
O(x3，y3)＝O(xmid，ymid)
A(x4，y4)＝A[xmid-1，kmid×(xmid-1)+bmid]
B(x5，y5)＝B(xmid-1，ymid)
当kmid＜0时，获得如下三个点：
O(x3，y3)＝O(xmid，ymid)
A(x4，y4)＝A[xmid+1，kmid×(xmid+1)+bmid]
B(x5，y5)＝B(xmid-1，ymid)
得到极角angle：
dot＝(x3-x4)×(x3-x5)+(y3-y4)×(y3-y5)
angle＝[(acos(dot/(m1×m2)))×180]/3.1415926
得到极径solar：
令
则有
S42：以卫星数字地图为基准，计算SAR图像中直线与卫星数字地图中直线的角度差theta，
计算Pcsar与Pcsate中angle的差值：angle-diff＝Pcsar(angle)-Pcsate(angle)，
将相同的angle-diff值进行计数，出现次数最多的angle_diff值记为theta；
S43：计算尺度参数S，计算步骤如下：
S431、采用Pcsar中的angle统计出SAR图像中平行直线组；
S432、计算SAR图像与卫星数字地图中每组平行直线极坐标的极径差
solarsar＝abs(Pcsar(solarn)-Pcsar(solarn-1))
solarsate＝abs(Pcsate(solarn)-Pcsate(solarn-1))
scale＝solarsate/solarsar
将相同的scale值进行计数，将出现次数最多的scale值记为尺度参数s；
S44、计算变换矩阵transMatric与平移参数
S441、分别计算SAR图像与卫星数字地图中直线特征的交点，分别表示为nodesar、nodesate，
S442、分别计算交点集合nodesar与nodesate中各点间坐标的差分node_diff(xdiff，ydiff)，
S443、将node_diff中的横坐标差值与纵坐标差值分别与阈值比较，保存满足条件的node_diff，同时将相同的node_diff值进行计数；
S444、获取出现次数最多的node_diff中的xdiff、ydiff；
S445、定义与SAR图像尺寸相关参数imgSarx、imgSary，与卫星数字地图尺寸相关参数imgSatex、imgSatey
S446、计算变换矩阵
则有：
S447、计算平移参数
S45、结合尺度参数s、变换矩阵transMatric、平移参数计算出处理后SAR图像在卫星数字地图中的匹配点MP(x，y)，设处理后SAR图像中的目标点为C(x，y)，则得到匹配点MP(x，y)：
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202210398024.1A CN114723794B (zh) | 2022-04-12 | 2022-04-12 | 基于lsd直线检测的sar图像配准方法 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202210398024.1A CN114723794B (zh) | 2022-04-12 | 2022-04-12 | 基于lsd直线检测的sar图像配准方法 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN114723794A true CN114723794A (zh) | 2022-07-08 |
| CN114723794B CN114723794B (zh) | 2023-03-24 |
Family
ID=82243152
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202210398024.1A Active CN114723794B (zh) | 2022-04-12 | 2022-04-12 | 基于lsd直线检测的sar图像配准方法 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN114723794B (zh) |
Cited By (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN115457222A (zh) * | 2022-09-14 | 2022-12-09 | 北京建筑大学 | 一种在地理信息系统中地理配准三维模型的方法 |
Citations (9)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN102663725A (zh) * | 2012-03-05 | 2012-09-12 | 西北工业大学 | 基于线特征和控制点的可见光和sar图像配准方法 |
| CN103345757A (zh) * | 2013-07-19 | 2013-10-09 | 武汉大学 | 多层次多特征约束下的光学和sar影像自动配准方法 |
| CN103514606A (zh) * | 2013-10-14 | 2014-01-15 | 武汉大学 | 一种异源遥感影像配准方法 |
| CN106886977A (zh) * | 2017-02-08 | 2017-06-23 | 徐州工程学院 | 一种多图自动配准及融合拼接方法 |
| CN107292922A (zh) * | 2017-06-23 | 2017-10-24 | 电子科技大学 | 一种用于光学与合成孔径雷达图像配准的方法 |
| CN109308715A (zh) * | 2018-09-19 | 2019-02-05 | 电子科技大学 | 一种基于点特征和线特征结合的光学图像配准方法 |
| CN111145228A (zh) * | 2019-12-23 | 2020-05-12 | 西安电子科技大学 | 基于局部轮廓点与形状特征融合的异源图像配准方法 |
| CN113077514A (zh) * | 2021-06-03 | 2021-07-06 | 浙江大学 | 一种sar图像船只弱尾迹增强与检测方法 |
| CN113870324A (zh) * | 2020-06-29 | 2021-12-31 | 上海微创卜算子医疗科技有限公司 | 多模态图像的配准方法及其配准装置和计算机可读存储介质 |
-
2022
- 2022-04-12 CN CN202210398024.1A patent/CN114723794B/zh active Active
Patent Citations (9)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN102663725A (zh) * | 2012-03-05 | 2012-09-12 | 西北工业大学 | 基于线特征和控制点的可见光和sar图像配准方法 |
| CN103345757A (zh) * | 2013-07-19 | 2013-10-09 | 武汉大学 | 多层次多特征约束下的光学和sar影像自动配准方法 |
| CN103514606A (zh) * | 2013-10-14 | 2014-01-15 | 武汉大学 | 一种异源遥感影像配准方法 |
| CN106886977A (zh) * | 2017-02-08 | 2017-06-23 | 徐州工程学院 | 一种多图自动配准及融合拼接方法 |
| CN107292922A (zh) * | 2017-06-23 | 2017-10-24 | 电子科技大学 | 一种用于光学与合成孔径雷达图像配准的方法 |
| CN109308715A (zh) * | 2018-09-19 | 2019-02-05 | 电子科技大学 | 一种基于点特征和线特征结合的光学图像配准方法 |
| CN111145228A (zh) * | 2019-12-23 | 2020-05-12 | 西安电子科技大学 | 基于局部轮廓点与形状特征融合的异源图像配准方法 |
| CN113870324A (zh) * | 2020-06-29 | 2021-12-31 | 上海微创卜算子医疗科技有限公司 | 多模态图像的配准方法及其配准装置和计算机可读存储介质 |
| CN113077514A (zh) * | 2021-06-03 | 2021-07-06 | 浙江大学 | 一种sar图像船只弱尾迹增强与检测方法 |
Non-Patent Citations (2)
| Title |
|---|
| XINYUE ZHANG ET AL.: ""Multimodal Remote Sensing Image Registration Methods and Advancements:A Survey"", 《REMOTE SENSING》 * |
| 刘中杰 等: ""基于控制线方法的机载SAR和可见光图像匹配应用研究"", 《航空学报》 * |
Cited By (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN115457222A (zh) * | 2022-09-14 | 2022-12-09 | 北京建筑大学 | 一种在地理信息系统中地理配准三维模型的方法 |
| CN115457222B (zh) * | 2022-09-14 | 2023-06-16 | 北京建筑大学 | 一种在地理信息系统中地理配准三维模型的方法 |
Also Published As
| Publication number | Publication date |
|---|---|
| CN114723794B (zh) | 2023-03-24 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN109409292B (zh) | 基于精细化特征优化提取的异源图像匹配方法 | |
| CN111626190B (zh) | 基于聚类分区进行刻度识别的水位监测方法 | |
| US11403839B2 (en) | Commodity detection terminal, commodity detection method, system, computer device, and computer readable medium | |
| CN112818988B (zh) | 一种指针式仪表自动识别读数方法及系统 | |
| CN106529559A (zh) | 一种指针式圆形多仪表盘实时读数识别方法 | |
| CN101650784B (zh) | 一种利用结构上下文特征进行图像匹配的方法 | |
| Choi et al. | Vehicle detection from aerial images using local shape information | |
| CN114972453B (zh) | 改进的基于lsd与模板匹配的sar图像区域配准方法 | |
| CN108564092A (zh) | 基于sift特征提取算法的向日葵病害识别方法 | |
| CN112085772A (zh) | 一种遥感图像配准方法及装置 | |
| CN112734729B (zh) | 适用于夜间补光条件的水尺水位线图像检测方法、装置及存储介质 | |
| CN107862708A (zh) | 一种sar与可见光图像配准方法 | |
| CN107862319B (zh) | 一种基于邻域投票的异源高分光学影像匹配误差剔除方法 | |
| CN111144250A (zh) | 融合雷达和光学遥感数据的土地覆盖分类方法 | |
| CN114723794B (zh) | 基于lsd直线检测的sar图像配准方法 | |
| CN111476246A (zh) | 应用于复杂环境下指针仪表鲁棒高效智读方法 | |
| CN110390338A (zh) | 一种基于非线性引导滤波与比率梯度的sar高精度匹配方法 | |
| CN111311596A (zh) | 一种基于改进lbp特征的遥感影像变化检测方法 | |
| Changjie et al. | Algorithm of remote sensing image matching based on corner-point | |
| CN110929598A (zh) | 基于轮廓特征的无人机载sar图像匹配方法 | |
| CN115760860B (zh) | 一种基于dxf文件导入的多类型工件尺寸视觉测量方法 | |
| CN104732190B (zh) | 一种基于正交纹理相关分析的合成孔径声纳目标检测方法 | |
| CN116295363A (zh) | 一种星点快速提取与高精度定位方法 | |
| CN115049708B (zh) | 基于lsd直线检测与模板匹配的sar图像配准方法 | |
| CN114898347A (zh) | 一种指针式仪表机器视觉识别方法 |
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 | ||
| CP01 | Change in the name or title of a patent holder | ||
| CP01 | Change in the name or title of a patent holder |
Address after: 210000 Room 301, floor 3, building 75, zone B, entrepreneurship and innovation city, No. 15, Fengji Avenue, Yuhuatai District, Nanjing, Jiangsu Province Patentee after: Nanjing Thunderbolt Information Technology Co.,Ltd. Address before: 210000 Room 301, floor 3, building 75, zone B, entrepreneurship and innovation city, No. 15, Fengji Avenue, Yuhuatai District, Nanjing, Jiangsu Province Patentee before: NANJING LEADING INFORMATION TECHNOLOGY Co.,Ltd. |
|
| PE01 | Entry into force of the registration of the contract for pledge of patent right | ||
| PE01 | Entry into force of the registration of the contract for pledge of patent right |
Denomination of invention: SAR image registration method based on LSD line detection Granted publication date: 20230324 Pledgee: Bank of Nanjing Co.,Ltd. Jiangning sub branch Pledgor: Nanjing Thunderbolt Information Technology Co.,Ltd. Registration number: Y2024980019591 |