基于先验DEM数据的星载SAR图像几何精校正方法与流程

本发明属于信号处理技术领域，特别涉及一种基于先验dem数据的星载sar图像几何精校正方法，适用于减小星载sar因为系统误差造成的sar图像场景中的像素点成像位置偏离真实位置点的问题，从而改善星载sar的成像质量，有助于用户对sar图像的判读和解译。

背景技术：

星载sar是从空间对地观测的一种有效手段，具有全天时、全天候、多波段及高分辨率的优点，能够详细、准确地测绘地形、地貌，获取地球表面的信息，最终生成目标场景的高分辨地图。

但是星载sar在从雷达脉冲的发射一直到最后通过成像算法生成星载sar图像的过程中，会有很多非理想因素误差引入其中，这样会使得实际得到的星载sar图像中每个像素点对应的场景中位置和实际场景中的位置都会有一定的偏移，使得实际得到的星载sar图像与真实场景并不完全对应，进而无法详细、准确地测绘地形、地貌，获取地球表面的信息。

技术实现要素：

针对上述现有技术存在的不足，本发明的目的在于提出一种基于先验dem数据的星载sar图像几何精校正方法，该种基于先验dem数据的星载sar图像几何精校正方法使用星载sar回波数据进行成像处理，得到实际星载sar图像，然后由星载sar的pos数据结合rd模型定位星载sar所在场景中每个像素点的坐标，利用ulaby模型产生后向散射系数矩阵，反演出理想的星载sar图像，理想的星载sar图像和实际sar图像进行比对，求出两者之间的偏移量，从而修正实际sar图像因全链路非理想因素而引起的系统误差，最终完成星载sar图像的几何精校正。

为达到上述技术目的，本发明采用如下技术方案予以实现。

一种基于先验dem数据的星载sar图像几何精校正方法，包括以下步骤：

步骤1，获取待校正的实际星载sar图像，所述待校正的实际星载sar图像为p×n维矩阵，p表示待校正的实际星载sar图像的方位向包含的像素点个数，n表示待校正的实际星载sar图像的距离向包含的像素点个数，p、n分别为大于0的正整数；

步骤2，建立待校正的实际星载sar图像的距离多普勒模型，记为成像模型，然后利用所述成像模型对待校正的实际星载sar图像中的p×n个像素点分别进行定位，分别得到待校正的实际星载sar图像中p×n个像素点各自的定位坐标；

步骤3，根据待校正的实际星载sar图像中p×n个像素点各自的定位坐标，计算得到待校正的实际星载sar图像中p×n个像素点各自的电磁波入射角；

步骤4，根据待校正的实际星载sar图像中p×n个像素点各自的电磁波入射角对待校正的实际星载sar图像进行反演，得到最终反演的理想星载sar图像；

步骤5，分别计算待校正的实际星载sar图像相对于最终反演的理想星载sar图像在方位向的像素偏差量na，以及待校正的实际星载sar图像相对于最终反演的理想星载sar图像在距离向的像素偏差量nr；

步骤6，根据待校正的实际星载sar图像相对于最终反演的理想星载sar图像在方位向的像素偏差量na，以及待校正的实际星载sar图像相对于最终反演的理想星载sar图像在距离向的像素偏差量nr，对待校正的实际星载sar图像进行几何精校正，得到几何精校正后的实际星载sar图像。

本发明与现有的技术相比具有以下优点：

第一，本发明利用先验dem信息，可以实现无控制点的图像几何精校正，不需要在在场景中布置角反射器，减少了工作量，提高了效率。

第二，本发明中利用了srtm库中的dem来提高定位精度，从而使星载sar图像的反演更加准确，提高了几何精校正的准确性。

附图说明

下面结合附图和具体实施方式对本发明作进一步详细说明。

图1为本发明的一种基于先验dem数据的星载sar图像几何精校正方法流程图；

图2为待校正的实际星载sar图像的距离多普勒模型示意图；

图3为求解待校正的实际星载sar图像中j个像素点各自的定位坐标的流程图；

图4为计算待校正的实际星载sar图像中的像素点电磁波入射角过程示意图；

图5为将几何精校正后的实际星载sar图像中p×n个像素点分别转换到经纬度坐标系中的过程示意图；

图6为高分三号星载sar图像示意图；

图7a为图6中地形起伏较大的山地部分示意图；

图7b为根据本发明方法对几何精校正后的图7a中每一个像素点分别进行反演后得到的理想星载sar图像示意图；

图8为对图6进行几何精校正后得到的星载sar图像示意图；

图9为图8经过地理编码后得到的结果示意图；

图10为对图9中的每一个像素点进行评估后的结果示意图。

具体实施方式

参照图1，为本发明的一种基于先验dem数据的星载sar图像几何精校正方法流程图；其中所述基于先验dem数据的星载sar图像几何精校正方法，包括以下步骤：

步骤1，获取待校正的实际星载sar图像，所述待校正的实际星载sar图像为p×n维矩阵，p表示待校正的实际星载sar图像的方位向包含的像素点个数，n表示待校正的实际星载sar图像的距离向包含的像素点个数，p、n分别为大于0的正整数。

具体地，确定星载sar，并获取星载sar回波数据，然后通过sar成像算法对获取的星载sar回波数据进行成像处理，得到待校正的实际星载sar图像；也可以直接下载网上提供的仅仅进行了成像处理的初级sar图像。

步骤2，建立待校正的实际星载sar图像的距离多普勒模型，记为成像模型，然后利用所述成像模型对待校正的实际星载sar图像中的p×n个像素点分别进行定位，分别得到待校正的实际星载sar图像中p×n个像素点各自的定位坐标。

2.1建立待校正的实际星载sar图像的距离多普勒模型，如图2所示，在所述距离多普勒模型中，卫星沿着自身轨道飞行，且卫星自带的gps会对星载雷达的空间位置坐标进行记录，星载雷达位于卫星上，星载雷达的空间位置包括m个点，m个点各自的坐标都是wgs84坐标系中的坐标，wgs84坐标系的原点是地球的中心；赤道和0度经线的交点，与地球中心的连线构成wgs84坐标系的x轴；北极点和地球的中心连线构成wgs84坐标系的z轴，再由右手坐标系原则确定出wgs84坐标系的y轴；将星载雷达的空间位置包括的m个点中第i个点记为si，i＝1,2,3,...,m，第i个点的坐标为(xi,yi,zi)，xi表示第i个点si在wgs84坐标系的x轴坐标值，yi表示第i个点si在wgs84坐标系的y轴坐标值，zi表示第i个点si在wgs84坐标系的z轴坐标值；将m个点各自的坐标记为星载雷达pos数据，星载雷达pos数据对应包括m个pos数据；将m个点各自的坐标进行逐点连接后形成星载雷达航迹。

根据星载雷达pos数据，并由速度和位移关系，计算得到星载雷达的空间位置中第l个点处的瞬时速度prf表示为脉冲重复频率,sl+1表示第l+1个点,sl表示第l个点,m表示星载雷达的空间位置包含的点个数,与pos数据个数取值相等。

星载雷达在卫星轨道上飞行的过程中，会向地面发射电磁波，就像黑夜里墙头上的手电筒倾斜的照射到地面上，在地面上形成一个椭圆形的光照区域一样，星载雷达发射的电磁波照射到地面上时也会形成一个椭圆的照射区域，随着星载雷达从第1个pos数据到第m个pos数据运动，该椭圆的照射区域会慢慢在地面上移动，直到星载雷达运动到第m个pos数据的位置,对应椭圆的照射区域形成一个长方形条带，将所述长方形条带记为合成孔径雷达(sar)图像成像区域,sar图像成像区域对应一个长方形条带。

将与星载雷达航迹平行的sar图像成像区域的对应方向,记作sar图像成像区域的方位向；将与星载雷达航迹垂直的sar图像成像区域的对应方向,记作sar图像成像区域的距离向，且sar图像成像区域对应待校正的实际星载sar图像，待校正的实际星载sar图像的方位向包含p个像素点，待校正的实际星载sar图像的距离向包含n个像素点，因此待校正的实际星载sar图像为p×n维矩阵。

2.2rd模型即距离多普勒模型，距离多普勒模型的提出是为了对待校正的实际星载sar图像中的p×n个像素点分别进行定位，即利用m个pos数据结合成像模型对sar图像成像区域中的p×n个像素点分别进行定位。

将待校正的实际星载sar图像中第m行、第n列处的像素点记作tmn，tmn的坐标记作(xmn,ymn,zmn),xmn表示待校正的实际星载sar图像中第m行、第n列处的像素点在wgs84坐标系的x轴坐标值，ymn表示待校正的实际星载sar图像中第m行、第n列处的像素点在wgs84坐标系的y轴坐标值，zmn表示待校正的实际星载sar图像中第m行、第n列处的像素点在wgs84坐标系的z轴坐标值。

星载雷达空间位置包含的点个数与待校正的实际星载sar图像的方位向包含的像素数个数取值相同，且一一对应；将第i个点si和待校正的实际星载sar图像中第m行、第n列处的像素点tmn之间的斜距记为n＝1,2,3,…,n，m＝1,2,3,…,p，i和m取值相同，且通过星载雷达开始接收回波的时间和星载雷达接收的回波在距离向的采样率求出,且tstarti表示星载雷达的空间位置中第i个点si开始接收回波的时间,fs表示星载雷达的距离向采样率,c为光速；由此可以看出,第i个点si和待处理的实际星载sar图像中第m行、第n列处的像素点tmn之间的斜距与m无关。

对于待校正的实际星载sar图像中p×n个像素点各自的坐标，通过对应得到p×n个斜距均求解出来。

2.3以下通过三个非线性方程组去求解待校正的实际星载sar图像中第m行、第n列处的像素点tmn的坐标。

方程一：距离方程

星载雷达开始接收回波时间和星载雷达接收的回波距离向采样率都是已知的，对于待校正的实际星载sar图像中的每个像素点，通过星载雷达开始接收回波的时间和星载雷达接收的回波距离向采样率准确确定每个像素点反射回波到星载雷达的时间，因而可以知道星载雷达与该像素点的距离，如图1中所示；也就是说，通过距离方程得到第i个点si和待校正的实际星载sar图像中第m行、第n列处的像素点tmn之间的斜距记为其满足方程如下：

显然只用方程(1)无法求出待校正的实际星载sar图像中第m行、第n列处的像素点tmn，要解出tmn的坐标(xmn,ymn,zmn),必须寻找其他方程。

方程二：多普勒方程：

计算得到星载雷达的多普勒频率fd求解方法如下：

其中表示星载雷达的径向速度，→表示矢量，λ表示为星载雷达载频对应的波长；对于本发明要研究的问题，计算得到第i个点si和待校正的实际星载sar图像中第m行、第n列处的像素点tmn之间的斜距对应的多普勒频率其表达式为：

其中，表示星载雷达的空间位置中第l个点处的瞬时速度矢量，通过pos数据可以求解；表示第i个点si和待校正的实际星载sar图像中第m行、第n列处的像素点tmn之间的斜距矢量，→表示矢量，λ为星载雷达载频对应的波长，||表示求绝对值操作。

星载雷达采取正侧视观测模式时，星载雷达的空间位置中第l个点处的瞬时速度矢量和垂直，因此第i个点si和待校正的实际星载sar图像中第m行、第n列处的像素点tmn之间的斜距对应的多普勒频率取值为0，所以

方程三：待校正的实际星载sar图像中第m行、第n列处的像素点tmn满足的地球椭球方程为：

其中，re为地球的赤道半径，rp为地球的极地半径，re和rp都是已知量；h表示待校正的实际星载sar图像中第m行、第n列处的像素点tmn的海拔高度，且h是未知量；xmn表示待校正的实际星载sar图像中第m行、第n列处的像素点在wgs84坐标系的x轴坐标值，ymn表示待校正的实际星载sar图像中第m行、第n列处的像素点在wgs84坐标系的y轴坐标值，zmn表示待校正的实际星载sar图像中第m行、第n列处的像素点在wgs84坐标系的z轴坐标值；n＝1,2,3,…,n，m＝1,2,3,…,p，p表示待校正的实际星载sar图像的方位向包含的像素点个数，n表示待校正的实际星载sar图像的距离向包含的像素点个数。

因为待校正的实际星载sar图像所在场景是地球上的某一区域，所以将该区域内像素点的坐标分别代入地球椭球方程中，在h已知的情况下，将满足等式(4)。

联立方程(1)、(3)、(4)三个方程，三个未知数(假设h已知)，分别计算得到待校正的星载实际sar图像中第m行、第n列处的像素点tmn的坐标(xmn,ymn,zmn)，求解完成也就意味着完成待校正的实际星载sar图像中第m行、第n列处的像素点tmn的定位。

2.4使用迭代法求解一组方程，来对待校正的实际星载sar图像中第m行、第n列处的像素点tmn进行定位，tmn的坐标为(xmn,ymn,zmn)，xmn表示待校正的实际星载sar图像中第m行、第n列处的像素点在wgs84坐标系的x轴坐标值，ymn表示待校正的实际星载sar图像中第m行、第n列处的像素点在wgs84坐标系的y轴坐标值，zmn表示待校正的实际星载sar图像中第m行、第n列处的像素点在wgs84坐标系的z轴坐标值；n＝1,2,3,…,n，m＝1,2,3,…,p，p表示待校正的实际星载sar图像的方位向包含的像素点个数，n表示待校正的实际星载sar图像的距离向包含的像素点个数，m和n的初始值分别为1。

假设h已知，并联立式(1)、式(3)和式(4)三个方程，三个未知数，通过求解方程组得出待校正的实际星载sar图像中第m行、第n列处的像素点坐标；而实际情况中，一般是不知道sar图像场景中每一个像素点的海拔高度的，也就是说式(1)、式(3)和式(4)三个方程中包含四个未知数，三个方程，显然是无法求解的；这种情况下，本发明解决策略是先假设待校正的实际星载sar图像中第m行、第n列处的像素点海拔高度为0。

参照图3，为求解待校正的实际星载sar图像中j个像素点各自的定位坐标的流程图；初始化：令f表示第f次迭代，f的初始值为1，令h0表示待校正的实际星载sar图像中每个像素点的海拔高度预设值，其中每个像素点分别为sar图像场景中每个像素点；设定待校正的实际星载sar图像中包含j个像素点，令j表示第j个像素点，j∈{1,2,…,j}，j的初始值为1，j＝p×n。

2.5将第f-1次迭代后第j个像素点的海拔高度预设值hf-1带入方程(4)，并根据式(1)和式(3)，得到如下方程组：

其中，xi表示第i个点si在wgs84坐标系的x轴坐标值，xjf表示第f次迭代后第i个点si在wgs84坐标系的x轴坐标值，yi表示第i个点si在wgs84坐标系的y轴坐标值，yjf表示第f次迭代后第i个点si在wgs84坐标系的y轴坐标值，zi表示第i个点si在wgs84坐标系的z轴坐标值，zjf表示第f次迭代后第i个点si在wgs84坐标系的z轴坐标值，表示第f次迭代后第i个点si和待校正的实际星载sar图像中第j个像素点之间的斜距，fdjf表示第f次迭代后第i个点si和待校正的实际星载sar图像中第j个像素点之间的斜距对应的多普勒频率，表示星载雷达的空间位置中第i个点处的瞬时速度矢量，表示第f次迭代后第i个点si和待校正的实际星载sar图像中第j个像素点之间的斜距矢量，→表示矢量，λ表示星载雷达载频对应的波长，||表示求绝对值操作，re表示地球的赤道半径，hj(f-1)表示第f-1次迭代后第j个像素点的海拔高度预设值，rp表示地球的极地半径。

求解以上方程组，三个方程，三个未知数，计算得到第f次迭代后第j个像素点的坐标(xjf,yjf,zjf)，所述第f次迭代后第j个像素点的坐标(xjf,yjf,zjf)是在wgs84坐标系中的，将(xjf,yjf,zjf)转换到经纬度坐标系中，其转换过程直接调用matlab中的坐标转换函数ecef2lla，即将第f次迭代后第j个像素点的坐标(xjf,yjf,zjf)转换到了经纬度坐标系中，得到第f次迭代后第j个像素点的经纬度坐标(ljf,bjf,hjf)，其中ljf为第f次迭代后第j个像素点的经度，bjf为第f次迭代后第j个像素点的纬度，hjf为第f次迭代后第j个像素点的海拔高度。

然后利用航天飞机雷达地形测绘使命(srtm)数据库中的数字高程模型(dem)数据去找出地球上经纬度(ljf,bjf)对应的海拔高度hjf'，hjf'为经纬度(ljf,bjf)对应的标准海拔高度,且此时的第f-1次迭代后第j个像素点的海拔高度预设值hf-1和hjf'不相等；其中，所述srtm数据库中的dem数据为先验dem数据。

2.6如果|hf-1-hjf'|≥ε0,ε0表示设定的门限阈值，本实施例中ε0＝10^-2，则令f加1，将经纬度(ljf,bjf)对应的海拔高度hjf'作为第f-1次迭代后第j个像素点的海拔高度预设值hf-1，返回2.5。

如果|hf-1-hjf'|<ε0，则迭代停止，并将迭代停止时得到的第f次迭代后第j个像素点的坐标(xjf,yjf,zjf)，作为待校正的实际星载sar图像中第j个像素点的定位坐标。

2.7令j加1，返回2.5，直到得到待校正的实际星载sar图像中第1个像素点的定位坐标至待校正的实际星载sar图像中第j个像素点的定位坐标，此时求得了待校正的实际星载sar图像中j个像素点各自的定位坐标，完成了待校正的实际星载sar图像中j个像素点各自的准确定位，j＝p×n。

步骤3，根据待校正的实际星载sar图像中p×n个像素点各自的定位坐标，计算得到待校正的实际星载sar图像中p×n个像素点各自的电磁波入射角。

具体地，步骤2结束后可得到待校正的实际星载sar图像中p×n个像素点各自的准确定位，在星载雷达的pos数据已知的情况下，根据空间几何关系得到每个像素点处电磁波的局部入射角，为下一步的反演星载sar图像作准备。

参照图4，为计算待校正的实际星载sar图像中的像素点电磁波入射角过程示意图；步骤3的子步骤为：

3.1根据向量夹角计算公式，计算得到待校正的实际星载sar图像中第m行、第n列处像素点的电磁波入射角θmn：

其中，tmn表示待校正的实际星载sar图像中第m行、第n列处的像素点，tm(n+2)表示待校正的实际星载sar图像中第m行、第n+2列处的像素点，当n＝n-1时，待校正的实际星载sar图像中第m行、第n+1列处的像素点tm(n+1)不存在，直接舍去，即记为0；当n＝n时，待校正的实际星载sar图像中第m行、第n+2列处的像素点tm(n+2)不存在，直接舍去，即记为0；si表示第i个点，n＝1,2,3,…,n，n+2＝1,2,3,…,n，m＝1,2,3,…,p，i＝1,2,3,…,m，p表示待处理的实际星载sar图像的方位向包含的像素数个数，m表示星载雷达的空间位置包含的点个数，i和m取值相等且一一对应；·表示点乘,*表示乘法，n表示待校正的实际星载sar图像的距离向包含的像素点个数，→表示矢量。

3.2令m不变，且n分别取1至n，重复执行3.1，进而分别得到待校正的实际星载sar图像中第m行、第1列处像素点的电磁波入射角θm1至待校正的实际星载sar图像中第m行、第n列处像素点的电磁波入射角θmn，记为待校正的实际星载sar图像中第m行n个像素点的电磁波入射角θm。

3.3令m分别取1至p，重复执行3.2，进而分别得到待校正的实际星载sar图像中第1行n个像素点的电磁波入射角θ1至待校正的实际星载sar图像中第p行n个像素点的电磁波入射角θp，记为待校正的实际星载sar图像中p×n个像素点各自的电磁波入射角。

步骤4，根据待校正的实际星载sar图像中p×n个像素点各自的电磁波入射角，以及ulaby地物散射系数经验模型，对待校正的实际星载sar图像进行反演，得到最终反演的理想星载sar图像sima。

步骤4的子步骤为：

4.1根据大量的实验和研究，待校正的实际星载sar图像的幅度与其电磁波入射角存在近似关系；即根据待校正的实际星载sar图像中第m行、第n列处像素点的电磁波入射角θmn，计算得到理想星载sar图像中第m行、第n列处像素点的幅度amn：

amn＝1-sin(θmn)(7)

其中，sin表示求正弦操作。

4.2令m不变，且n分别取1至n，重复执行4.1，进而分别得到理想星载sar图像中第m行、第1列处像素点的幅度am1至理想星载sar图像中第m行、第n列处像素点的幅度amn，记为理想星载星载sar图像中第m行n个像素点的幅度am。

4.3令m分别取1至p，重复执行4.2，进而分别得到理想sar图像中第1行n个像素点的幅度a1至理想星载sar图像中第p行n个像素点的幅度ap，记为理想星载sar图像中p×n个像素点各自的幅度。

上述的公式(7)是表示的一个近似关系，而且并没有考虑到雷达的工作频率、极化方式、地物类型等因素对后向散射系数的影响；要想得出更为准确的反演sar图像，需要用到ulaby地物散射系数经验模型。

4.4ulaby地物散射系数经验模型是美国密歇根大学f.t.ulaby和m.c.dobson提出的一个简单灵活的模型，他们二人综合了20实际60年代至80年代公开发表的大量地面后向散射系数测量数据，提出了ulaby模型经验公式，来更加准确的由电磁波入射角计算场景的后向散射系数，进而计算得到基于ulaby地物散射系数经验模型的理想星载sar图像中第m行、第n列处像素点的幅度σmn：

σmn(db)＝p1+p2exp(-p3θmn)+p4cos(p5θmn+p6)(8)

其中，db是分贝的单位，cos表示求余弦操作，exp表示指数函数，p1表示ulaby地物散射系数经验模型的常数项，p2表示ulaby地物散射系数经验模型的指数系数项，p3表示ulaby地物散射系数经验模型的指数相位项，p4表示ulaby地物散射系数经验模型的三角系数项，p5表示ulaby地物散射系数经验模型的三角相位系数项，p6表示ulaby地物散射系数经验模型的三角相位常数项；p1、p2、p3、p4、p5和p6都是考虑星载雷达工作波段、星载雷达发射信号的极化方式、观测星载雷达所在场景的地物类型等因素，并经过大量实验而得出的经验值，详细数值可以参见雷达手册。

4.5令m不变，且n分别取1至n，重复执行4.4，进而分别得到基于ulaby地物散射系数经验模型的理想星载sar图像中第m行、第1列处像素点的幅度σm1至基于ulaby地物散射系数经验模型的理想星载sar图像中第m行、第n列处像素点的幅度σmn，记为基于ulaby地物散射系数经验模型的理想星载sar图像中第m行n个像素点的幅度σm。

4.6令m分别取1至p，重复执行4.5，进而分别得到基于ulaby地物散射系数经验模型的理想星载sar图像中第1行n个像素点的幅度σ1至基于ulaby地物散射系数经验模型的理想星载sar图像中第p行n个像素点的幅度σp，记为基于ulaby地物散射系数经验模型的理想星载sar图像中p×n个像素点各自的幅度，并将所述基于ulaby地物散射系数经验模型的理想星载sar图像中p×n个像素点各自的幅度作为经过反演得到的理想星载sar图像，简记为最终反演的理想星载sar图像sima。

步骤5，分别计算待校正的实际星载sar图像相对于最终反演的理想星载sar图像sima在方位向的像素偏差量na，以及待校正的实际星载sar图像相对于最终反演的理想星载sar图像sima在距离向的像素偏差量nr。

步骤5的子步骤为：

5.1对待校正的实际星载sar图像和最终反演的理想星载sar图像sima进行配准可以使用很多方法，本发明采用相关函数法在频域实现，将待校正的实际星载sar图像和最终反演的理想星载sar图像sima分别变换到二维频域，分别得到待校正的实际星载sar图像的二维频谱s2f，以及最终反演的理想星载sar图像的二维频谱σ2f，其表达式分别为：

s2f＝fft(fft(sreal))(9)

σ2f＝fft(fft(σmn))(10)

其中，fft表示做快速傅立叶变换操作，sreal表示待校正的实际星载sar图像，σmn表示最终反演的理想星载sar图像sima，下标2f表示二维频域。

使用max函数分别确定待校正的实际星载sar图像的二维频谱s2f的峰值位置和最终反演的理想星载sar图像的二维频谱σ2f的峰值位置，然后比较s2f的峰值位置和σ2f的峰值位置的偏差，分别将s2f的峰值位置记作(samax，srmax)，将σ2f的峰值位置记作(σamax，σrmax)，samax表示待校正的实际星载sar图像的二维频谱最大值所在位置在方位向的坐标值，srmax表示待校正的实际星载sar图像的二维频谱最大值所在位置在距离向的坐标值，σamax表示最终反演的理想星载sar图像的二维频谱最大值所在位置在方位向的坐标值，σrmax表示最终反演的理想星载sar图像的二维频谱最大值所在位置在距离向的坐标值。

5.2进而分别计算得到待校正的实际星载sar图像相对于最终反演的理想星载sar图像sima在方位向的像素偏差量na，以及待校正的实际星载sar图像相对于最终反演的理想星载sar图像sima在距离向的像素偏差量nr，其表达式分别为：

na＝samax-σamax(11)

nr＝srmax-σrmax(12)。

步骤6，根据待校正的实际星载sar图像相对于最终反演的理想星载sar图像sima在方位向的像素偏差量na，以及待校正的实际星载sar图像相对于最终反演的理想星载sar图像sima在距离向的像素偏差量nr，使用matlab中二维插值函数interp2对待校正的实际星载sar图像进行几何精校正，得到几何精校正后的实际星载sar图像。

具体地，根据待校正的实际星载sar图像相对于最终反演的理想星载sar图像在方位向的像素偏差量na，以及待校正的实际星载sar图像相对于最终反演的理想星载sar图像距离向的像素偏差量nr，使用matlab中二维插值函数interp2对几何精校正后的实际星载sar图像进行几何精校正，得到几何精校正后的实际星载sar图像，即分别令x为几何精校正后的实际星载sar图像中第x行,x＝1,2,3,…,p；y为几何精校正后的实际星载sar图像中第y列,y＝1,2,3,…,n，p表示待处理的实际星载sar图像的方位向包含的像素数个数，与几何精校正后的实际星载sar图像的列数取值相等；n表示待校正的实际星载sar图像的距离向包含的像素点个数，与几何精校正后的实际星载sar图像的行数取值相等。

将几何精校正后的实际星载sar图像中第x行、第y列处的像素点记为重采样后的实际星载sar图像srxy，srxy＝interp2(x，y，sreal，x+nr，y+na)，sreal表示待校正的实际星载sar图像，interp2表示matlab中二维插值函数。

x＝1,2,3,…,p，y＝1,2,3,…,n，遍历待校正的实际星载sar图像中p×n个像素点，进而得到几何精校正后的实际星载sar图像。

步骤7，对几何精校正后的实际星载sar图像进行地理编码，得到地理编码后的几何精校正星载sar图像。

具体地，几何精校正后的实际星载sar图像中p×n个像素点的坐标都是wgs84坐标系中的，现在将几何精校正后的实际星载sar图像中p×n个像素点分别转换到经纬度坐标系中；使用matlab中自带的一个函数ecef2lla将几何精校正后的实际星载sar图像中p×n个像素点分别转换到经纬度坐标系中，其转换过程如图5所示。

图5中的左图是几何精校正后的实际星载sar图像，图4中的右图是由经度轴和纬度轴组成的经纬度坐标系；当知道几何精校正后的实际星载sar图像中每一个像素点的wgs84坐标后，通过ecef2lla函数将几何精校正后的实际星载sar图像中每一个像素点的wgs84坐标分别转换到经纬度坐标系中，具体为：

7.1令几何精校正后的实际星载sar图像中第m行、第n列处的像素点为t'mn，t'mn的wgs84坐标为(x'mn,y'mn,z'mn)，x'mn表示几何精校正后的实际星载sar图像中第m行、第n列处的像素点t'mn在wgs84坐标系中x轴的坐标值，y'mn表示几何精校正后的实际星载sar图像中第m行、第n列处的像素点t'mn在wgs84坐标系中y轴的坐标值，z'mn表示几何精校正后的实际星载sar图像中第m行、第n列处的像素点t'mn在wgs84坐标系中z轴的坐标值。

将几何精校正后的实际星载sar图像中第m行、第n列处的像素点t'mn的wgs84坐标对应转换到经纬度坐标系中，得到几何精校正后的实际星载sar图像中第m行、第n列处的像素点经纬度坐标(lmn,bmn,hmn)，lmn表示几何精校正后的实际星载sar图像中第m行、第n列处像素点的经度值，bmn表示几何精校正后的实际星载sar图像中第m行、第n列处像素点的纬度值，hmn表示几何精校正后的实际星载sar图像中第m行、第n列处像素点的高度值。

7.2令m不变，且n分别取1至n，重复执行7.1，进而分别得到几何精校正后的实际星载sar图像中第m行、第1列处的像素点经纬度坐标(lm1,bm1,hm1)至几何精校正后的实际星载sar图像中第m行、第n列处的像素点经纬度坐标(lmn,bmn,hmn)，记为几何精校正后的实际星载sar图像中第m行n个像素点的经纬度坐标。

7.3令m分别取1至p，重复执行7.2，进而分别得到几何精校正后的实际星载sar图像中第1行n个像素点的经纬度坐标至几何精校正后的实际星载sar图像中第p行n个像素点的经纬度坐标，记为几何精校正后的实际星载sar图像中p×n个像素点各自的经纬度坐标，并将所述几何精校正后的实际星载sar图像中p×n个像素点各自的经纬度坐标作为地理编码后的几何精校正星载sar图像。

结合实测数据对本发明效果作进一步验证说明。

参照图6，为高分三号星载sar图像示意图，图6中每个像素点的经纬度和实际经纬度都是存在偏差的，这些偏差是由于星载雷达的系统误差造成的。

图7a为图6中地形起伏较大的山地部分示意图；用所述地形起伏较大的山地部分确定星载sar图像的系统误差，从而完成几何精校正。

图7b为根据本发明方法对几何精校正后的图7a中每一个像素点分别进行反演后得到的理想星载sar图像示意图；

图8为对图6进行几何精校正后得到的星载sar图像示意图；

选取图8所在场景中的一个像素点来验证几何精校正的效果，图8上半部分是城区，在密集的城区中有一块矩形的建筑区，如图8中的十字线处，选取图8矩形建筑区的右下角这个像素点进行验证。

用googleearth准确找出这个像素点的经纬度高程坐标为：(34.896336，112.816433，113)，第一个参数是纬度，第二个参数是经度，第三个参数是高程，该经纬度高程坐标是该像素点的实际坐标。

在未进行几何精校正之前，经过定位这个点的经纬度高程为：(34.897016，112.817699，113)，这个坐标和实际坐标的距离误差为138.1368m。

进行几何精校正之后，这个点的经纬度高程为：(34.896536，112.816483，113)，这个坐标和实际坐标的距离误差为22.6720m。

可以看出，经过本专利的几何精校正技术，sar图像中像素点的位置更加准确，sar图像的利用价值更高；图6和图8整体看上去差距不大，但是选择其中的一个点去进行评估，就会发现经过几何精校正后的图像比校正之前，每个像素点的位置有了更为精确的改善；图9为图8经过地理编码后得到的结果示意图；图10为对图9中的每一个像素点进行评估后的结果示意图。

综上所述，仿真实验验证了本发明的正确性，有效性和可靠性。

显然，本领域的技术人员可以对本发明进行各种改动和变型而不脱离本发明的精神和范围；这样，倘若本发明的这些修改和变型属于本发明权利要求及其等同技术的范围之内，则本发明也意图包含这些改动和变型在内。