双机载SAR对目标协同成像立体定位的方法与流程

本发明属于sar定位技术领域，具体涉及一种双机载sar对目标协同成像立体定位的方法。

背景技术：

合成孔径雷达(syntheticapertureradar,sar)具有全天候、全天时的对地观测能力，sar在海洋监测方面得到了广泛的应用，并发挥了巨大的社会、经济和军事效益。因成像波段不同，sar系统可以获取不同于光学系统的地物目标信息。

在一般情况下，军事目标在雷达图像上都是由具有较亮的散射中心组成，因而可以采用阈值分割的技术实现目标的检测。然而在具有随机性和复杂性的杂波背景下，sar所成的图像中，背景与目标的对比度随之变化，对固定阈值的图像分割技术往往会产生较多的虚假目标。基于滑动窗口的双参数恒虚警目标检测方法，是一种自适应门限检测方法，能够适应背景杂波的变化，是常用的目标检测算法之一。

sar图像定位与三维信息提取技术可以对较大范围的区域实现全天时、全天候和高精度的目标定位。现有的双基机载sar技术常用的是两部雷达分别作为发射机和接收机的工作模式，这样的方法无法由目标在sar图像的像素位置信息，通过直接解析方程组获得目标的实际位置信息。

技术实现要素：

本发明的目的在于提供一种双机载sar对目标协同成像立体定位的方法。

实现本发明目的的技术方案为：一种双机载sar对目标协同成像立体定位的方法，该方法包括以下步骤：

采用距离多普勒成像算法，对两个机载平台sar1、sar2回波数据进行成像处理；

采用目标检测算法，对sar1、sar2图像同时进行目标检测，输出目标分别在两幅图像中的像点坐标；

将像点坐标代入双sar协同立体定位模型，对目标的实际三维坐标进行牛顿迭代解算。

本发明与现有技术相比，具有以下优点：本发明能够更加有效地对sar图像中目标进行检测及定位；采用双参数cfar检测算法可以自适应背景杂波的变化，得出最适合的门限；将背景窗口分为四个部分进行杂波均值分析可以大大减少计算量；采用两机载sar数据融合定位的三维解算方法，可以不受相干性的限制，得到高精度定位结果。

附图说明

图1是本发明中双机载sar对目标协同成像立体定位方法原理图。

图2是本发明中双机载sar对目标协同成像立体定位方法总流程图。

图3是本发明中目标检测算法流程图。

图4是本发明中地面点三维坐标解算算法流程图。

图5是本发明中双参数cfar检测窗口结构示意图。

图6是本发明中双参数cfar背景窗口分块示意图。

具体实施方式

一种双机载sar对目标协同成像立体定位的方法，该方法包括以下步骤：

步骤1，采用距离多普勒成像算法，对两个机载平台sar1、sar2回波数据进行成像处理；

步骤2，采用目标检测算法，在复杂背景下对sar1、sar2图像同时进行目标检测，输出目标分别在两幅图像中的像点坐标；具体为：

采用双参数cfar算法分别对sar1、sar2所成图像检测指定目标，然后采用形态学滤波算法处理以消除检测结果中的虚假目标，再采用函数regionprops统计被标记的区域的面积分布，显示区域总数，输出检测结果，即目标在图中的像点坐标范围，最后分别取两组像点坐标矩阵的几何中心点作为目标具体像点坐标。

双参数cfar算法的检测过程具体为：

(1a)设定一个检测单元由三个层层嵌套的窗口构成，三个窗口分别是目标窗口t、保护窗口p和背景窗口b，其中目标窗口内是待检测的点，上述3个窗口均正方形；

设定目标窗口边长为a，保护窗口的边长为b，背景窗口环形宽度为c，则得出用于杂波区域计算的像素数量记为numpix，表达式为：numpix＝2c·(2c+2b)；检测器边长记为d，表达式为：d＝b+2c；

(1b)对原sar图像i1、i2的边界扩充，扩充大小为cfar检测器边长的一半，即以对原图像i1、i2进行填充后的图像记为i1'、i2'；

(1c)确定cfar阈值

设置虚警概率为pfa，记cfar阈值为

(1d)利用cfar检测器求解局部阈值，执行单个像素判断，具体为：由杂波区域得到均值和标准偏差，计算双参数cfar检测判别式，遍历所有像素点，判别是否为目标点；

将cfar背景窗口划分为四个部分；

对每个像素(i,j)所对应的背景窗口区域估计，然后分别累加最后求背景区域的均值ub以及标准偏差σb；将该像素(i,j)对应的灰度值x代入双参数cfar检测判别公式进行计算，满足则判别该像素点为目标，赋值为255，否则判别该像素点为背景，赋值为0；由此得到双参数cfar检测后的二值矩阵ia1、ia2。

所述的形态学滤波具体步骤为：

(2a)创建一个以半径为r1的圆形结构元素矩阵b1，以结构元素b1对二值矩阵ia1、ia2进行闭运算，填补轮廓线断裂，得到结果ib1、ib2；

(2b)创建一个以半径为r2的圆形结构元素矩阵b2，以结构元素b2对二值矩阵ib1、ib2进行腐蚀运算，得到结果ic1、ic2；

(2c)以结构元素b2对二值矩阵ic1、ic2进行膨胀运算，得到结果id1、id2。

步骤3，将像点坐标代入双sar协同立体定位模型，对目标的实际三维坐标进行牛顿迭代解算；具体为：

将两组同名点像点坐标代入距离多普勒成像模型，结合sar1、sar2的定向参数及控制点坐标，采用牛顿迭代算法最后迭代出最符合实际条件的目标三维坐标；其中同名点即两幅图像中指代同一目标的点，定向参数包括机载sar的飞行速度和实时坐标。

所述的距离多普勒模型导入牛顿迭代解算三维坐标，具体步骤为：

(4a)由两幅sar图像中同名点的像点t1(i^l,j^l)、t2(i^r,j^r)代入距离公式及多普勒公式，得到同名像点坐标与对应地面点坐标(x,y,z)的关系由如下四个方程构成的方程组表示：

即：

其中，分别代表像点t1、像点t2成像瞬间sar1、sar2天线相位中心位置，分别代表像点t1、像点t2成像瞬间sar1、sar2天线相位中心速度，分别为sar1、sar2成像时的近距延迟，分别为sar1、sar2的斜距向采样间隔；

(4b)计算同名点成像瞬间两天线相位中心位置、速度

若以二次多项式表示天线相位中心位置与成像时刻的关系，则由解算得到的定向参数，利用下式可以分别求得像点t1、像点t2成像瞬间天线相位中心位置和成像瞬间天线相位中心速度

式中，t'为每行间的时间间隔；t^l和t^r分别为像点t1、像点t2在左右图像上的成像时刻；分别为sar1、sar2的天线相位中心加速度初值；分别为sar1、sar2的天线相位中心速度初值；分别为sar1、sar2的天线相位中心位置初值；

(4c)给定地面点三维坐标初值

假设获取测区n个控制点，每个控制点对应的像点坐标(ik,jk)及实际坐标(xk,yk,zk)，分别为:pgc1(i1,j1,x1,y1,z1)，pgc2(i2,j2,x2,y2,z2)，...pgcn(in,jn,xn,yn,zn)，k＝1,2,3...n；

则地面点三维坐标的初值为：

(4d)构建误差方程组

由r-d模型对地面点坐标的线性化形式，可知同名像点t1(i^lj^l)、t2(i^r,j^r)与所对应的地面点p(x,y,z)的线性化关系为：

c·△g-l＝0

其中，c是关于地面点坐标改正量的系数矩阵即：

式中，分别是关于像点t1成像时刻对应的sar1位置的函数；分别是关于像点t1成像时刻对应的sar1速度的函数；分别是关于像点t2成像时刻对应的sar2位置的函数；分别是关于像点t2成像时刻对应的sar2速度的函数。

系数阵c中的各元素分别为：

△g是地面点坐标的改正量向量，△g＝[△x△y△z]^t；

l是r-d模型方程组的初值向量，

(4e)计算三维坐标改正量

c矩阵的法方程式表示为：

c^tc△g-c^tl＝0

解法方程，可求得地面点三维坐标的改正量向量△g：

△g＝(c^tc)^-1c^tl

再在上一次迭代的基础上修正三维坐标初值：

(4f)限差判断

判断改正量是否小于给定限差，若改正量大于限差则返回(4d)，利用改正后的三维坐标重新组误差方程式计算其改正量；若改正量小于等于限差则停止迭代，输出计算出的地面点三维坐标

下面结合附图和实施例对本发明进行详细说明。

实施例

本实施例以杂波海洋环境下定位舰船目标为例，对本发明的示例性实施方式进行详细描述。

图1为本发明提出的一种双机载sar对目标协同成像立体定位的结构示意图，主要由两架固定翼无人机作为平台的两部sar雷达组成，两架无人机以相同方向在目标区域两侧飞行，同时对扫描区域进行成像。

图2为本发明提出的一种双机载sar对目标协同成像立体定位的总流程图，具体实施步骤为：

第一步，采用距离多普勒成像算法，对两个机载平台sar1、sar2接收到的海面扫描区域的回波数据进行高分辨率成像，分辨率达到1米。

所述的距离多普勒成像算法，包括：距离向压缩、距离徙动校正及方位向压缩处理。

距离向压缩的流程为，在频域中通过回波信号与匹配滤波函数相乘完成脉冲压缩，匹配滤波函数为接收频域信号所含的二次相位项的共轭函数。距离徙动校正的解决方法为在频域上乘以一个线性相位。再在方位向上通过与匹配滤波函数相乘完成脉冲压缩。

第二步，采用目标检测算法，流程图如图3所示，对sar1、sar2所成像进行杂波背景下舰船目标检测。所述的目标检测算法，具体步骤为：采用双参数cfar算法分别对sar1、sar2所成图像检测舰船目标，然后采用形态学滤波算法处理以消除检测结果中的虚假目标，再采用函数regionprops统计被标记的区域的面积分布，输出检测结果，即舰船目标在图中的像点坐标范围，再分别取两组像点坐标矩阵的几何中心点作为目标具体像点坐标t1(i^l,j^l)、t2(i^r,j^r)。

所述的双参数cfar检测舰船目标的流程为，确定cfar检测器参数，包括窗口尺寸，保护区宽度，杂波区宽度，遍历扩展后的图像检测得到图像上的高亮像素。

选取合适的滑动窗口，对sar图像中背景窗内参数的高斯分布进行估计，确定阈值将目标估计参数与之做比较，得出图像上的高亮像素区域。

所述的形态学滤波的流程为，从图像中提取连通区域，即可能的舰船目标，对各个区域进行几何特征提取，再基于特征排除不是舰船的目标。

所述的采用函数regionprops对细节进行进一步计算的流程为，采用函数regionprops统计被标记的区域的面积分布，显示区域总数。

所述的双参数cfar检测舰船目标，具体步骤为：

1a)设定一个检测单元由三个层层嵌套的窗口构成，如图5所示。三个窗口分别是目标窗口t、保护窗口p和背景窗口b，其中目标窗口内是待检测的点，海面背景信息从背景窗口获取，保护窗口是为了确保舰船目标部分不会被包含在背景窗口中。上述3个窗口均正方形。

设定目标窗口边长为a，保护窗口的边长为b，背景窗口环形宽度为c，则得出用于杂波区域计算的像素数量记为numpix，表达式为：numpix＝2c·(2c+2b)。检测器边长记为d，表达式为：d＝b+2c。

1b)对原sar图像i1、i2的边界扩充以消除边界对检测目标的影响，扩充大小为cfar检测器边长的一半，即以对原图像i1、i2进行填充后的图像记为i1'、i2'。

1c)确定cfar阈值

设置虚警概率为pfa，记cfar阈值为

1d)利用cfar检测器求解局部阈值，执行单个像素判断(由杂波区域得到均值和标准偏差，计算双参数cfar检测判别式，遍历所有像素点，判别是否为目标点)

将cfar背景窗口划分为四个部分，如图6所示。

对每个像素(i,j)所对应的背景窗口区域估计，然后分别累加最后求背景区域的均值ub以及标准偏差σb。将该像素(i,j)对应的灰度值x代入双参数cfar检测判别公式进行计算，满足则判别该像素点为目标，赋值为255，否则判别该像素点为背景，赋值为0。由此得到双参数cfar检测后的二值矩阵ia1、ia2。

所述的形态学滤波，具体步骤为：

2a)创建一个以半径为r1的圆形结构元素矩阵b1，以结构元素b1对二值矩阵ia1、ia2进行闭运算，填补轮廓线断裂，得到结果ib1、ib2。

2b)创建一个以半径为r2的圆形结构元素矩阵b2，以结构元素b2对二值矩阵ib1、ib2进行腐蚀运算，得到结果ic1、ic2。

2c)以结构元素b2对二值矩阵ic1、ic2进行膨胀运算，得到结果id1、id2。

所述的采用函数regionprops统计被标记的区域的面积分布，显示区域总数，具体步骤为：采用'area'作为度量数据对图像各个区域中像素总个数进行计算，并确定每个连通区域的坐标范围。

第三步，采用地面点三维坐标解算算法，对舰船目标的实际坐标进行解算。流程图如图4所示，具体步骤为：将两组同名点像点坐标代入距离多普勒成像模型，结合sar1、sar2的定向参数及控制点坐标，采用牛顿迭代算法解算三维坐标的流程为计算同名点成像瞬间天线相位中心位置及速度，代入距离多普勒基本方程并构建误差方程，最后迭代出最符合实际条件的目标三维坐标。

距离多普勒基本方程包括距离公式及多普勒频率方程：

距离公式：rs²＝(x-xs)²+(y-ys)²+(z-zs)²＝(r0+mslant·j)²

可记：f1＝(x-xs)²+(y-ys)²+(z-zs)²-(r0+mslant·j)²

其中，(x,y,z)表示地面点目标坐标，(xs,ys,zs)为成像瞬间天线相位中心位置，r0为近距延迟，mslant为斜距向采样间隔，j为像点的距离向坐标。

多普勒频率方程：

可记：

其中，(vx,vy,vz)表示成像瞬间天线相位中心速度，λ为雷达发射波波长，rs为地面点目标到雷达平台的瞬时位置距离，fdc为多普勒频移参数。

牛顿迭代解算三维坐标的流程为计算同名点成像瞬间天线相位中心位置及速度，代入距离多普勒基本方程并构建误差方程。

所述的距离多普勒模型导入牛顿迭代解算三维坐标，具体步骤为：

4a)由两幅sar图像中同名点的像点t1(i^l,j^l)、t2(i^r,j^r)代入距离公式及多普勒公式，得到同名像点坐标与对应地面点坐标(x,y,z)的关系由如下四个方程构成的方程组表示：

即：

其中，分别代表像点t1、像点t2成像瞬间sar1、sar2天线相位中心位置，分别代表像点t1、像点t2成像瞬间sar1、sar2天线相位中心速度，分别为sar1、sar2成像时的近距延迟，分别为sar1、sar2的斜距向采样间隔。

4b)计算同名点成像瞬间两天线相位中心位置、速度

若以二次多项式表示天线相位中心位置与成像时刻的关系，则由解算得到的定向参数，利用下式可以分别求得像点t1、像点t2成像瞬间天线相位中心位置和成像瞬间天线相位中心速度

式中，t'为每行间的时间间隔；t^l和t^r分别为像点t1、像点t2在左右图像上的成像时刻；分别为sar1、sar2的天线相位中心加速度初值；分别为sar1、sar2的天线相位中心速度初值；分别为sar1、sar2的天线相位中心位置初值。

4c)给定地面点三维坐标初值

假设获取测区n个控制点，每个控制点对应的像点坐标(ik,jk)(k＝1,2,3...n)及实际坐标(xk,yk,zk)，分别为:pgc1(i1,j1,x1,y1,z1)，pgc2(i2,j2,x2,y2,z2)，...pgcn(in,jn,xn,yn,zn)

则地面点三维坐标的初值可取为：

4d)构建误差方程组

由r-d模型对地面点坐标的线性化形式，可知同名像点t1(i^lj^l)、t2(i^r,j^r)与所对应的地面点p(x,y,z)的线性化关系为：

c·△g-l＝0

其中，c是关于地面点坐标改正量的系数矩阵即：

式中，分别是关于像点t1成像时刻对应的sar1位置的函数；分别是关于像点t1成像时刻对应的sar1速度的函数；分别是关于像点t2成像时刻对应的sar2位置的函数；分别是关于像点t2成像时刻对应的sar2速度的函数。

系数阵c中的各元素分别为：

△g是地面点坐标的改正量向量，△g＝[△x△y△z]^t。

l是r-d模型方程组的初值向量，

4e)计算三维坐标改正量

c矩阵的法方程式可表示为：

c^tc△g-c^tl＝0

解法方程，可求得地面点三维坐标的改正量向量△g：

△g＝(c^tc)^-1c^tl

再在上一次迭代的基础上修正三维坐标初值：

4f)限差判断

判断改正量是否小于给定限差，若改正量大于限差则返回(4d)，利用改正后的三维坐标重新组误差方程式计算其改正量；若改正量小于或等于限差则停止迭代，输出计算出的地面点三维坐标。