CN101915920A - 一种地球同步轨道合成孔径雷达卫星的高分辨率成像方法 - Google Patents

Abstract

本发明公开了一种地球同步轨道合成孔径雷达卫星的高分辨率成像方法，包括以下几个步骤，步骤一：子孔径划分；步骤二：子孔径数据距离向压缩；步骤三：子孔径数据距离向插值；步骤四：地面网格划分；步骤五：确定距离徙动曲线；步骤六：进行方位压缩处理；步骤七：子孔径图像叠加，得到最终的图像。本发明提出的成像方法可处理地球同步轨道合成孔径雷达卫星的全孔径回波数据，实现高分辨率成像；可实现多点目标成像，即对不同距离向和方位向目标的不同距离徙动曲线均可实现方位压缩；还可实现广域场景范围内的成像，具有成像区域自由化的特点，即可根据需求对不同位置的场景实现成像。

Description

一种地球同步轨道合成孔径雷达卫星的高分辨率成像方法

技术领域

本发明涉及一种地球同步轨道(GEO)合成孔径雷达(SAR)卫星的高分辨率成像方法，属于信号处理技术领域。

背景技术

星载合成孔径雷达具备全天时、全天候的工作能力，是目前备受关注的对地观测系统。遥感应用中，重复观测周期是非常重要的技术指标。地形动态监测、地表现象理解等领域都要求持续观测或频繁的重访。实验表明，对于土壤湿度检测，具有1天重访时间的观测系统能获取55％有效信息，而35天重访时间的系统仅能获取20％的有效信息。然而，现有的SAR卫星都运行在低轨，重复周期较长，难以满足高频重访的用户需求。解决此问题的一种有效途径就是将卫星轨道升高至地球同步轨道，它的重访周期短，仅为1天，并且观测范围大，具有很好的应用前景。

1978年Kiyo Tomiyasu首先在《Synthetic Aperture Radar in GeosynchronousOrbit》一文中提出了GEO SAR的概念，并初步设计了系统参数。2000年以后，GEO SAR成为国内外研究的热点问题。S.N.Madsen等提出了四种工作模式和相应的技术指标，指出天线技术是GEO SAR发展的关键。Wendy Edelstein等分析了不同轨道高度对瞬时视场范围的影响和对天线尺寸的要求，并在考虑高轨卫星系统复杂性的基础上绘制了技术路线图。S.E.Hobbs等研究了L波段GEO SAR的成像技术，重点分析了在超长积累时间中对流层和电离层扰动以及地球潮汐对雷达成像的影响。Davide Bruno等指出GEO SAR具有良好的覆盖和重访性能，在干涉处理等方面能够获得比LEO系统更好的效果，有助于理解地表运动现象。在国内，北京航空航天大学于泽等从成像处理的角度出发研究了GEO SAR的一些总体特性，包括方位分辨率、距离徙动量等等，其中方位分辨率通过理论和仿真分析证明可达到3米以下。

GEO SAR的成像技术一直是一个难点。现今，合成孔径雷达的载荷一般是飞机和低轨道卫星，其合成孔径时间是秒量级，在孔径时间内雷达与地面目标的相对运动轨迹近似为直线，回波信号的多普勒历程近似为二次多项式，传统的SAR成像算法均可适用。而GEO SAR具有超高的轨道高度(约36000公里)，其合成孔径时间达到千秒量级，在孔径时间内它的星地相对运动轨迹已经无法再用直线近似，其回波信号的多普勒历程需要用高阶多项式进行描述。同时波束地面移动速度受地球自转影响较大，这都使得广域场景目标的多普勒特性具有强空变特征，给成像处理带来很大难度。传统的SAR成像方法(RD、CS、ωk算法)都是基于直线运动展开推导的，并且它们的方位压缩处理都是在方位频域进行，对多普勒历程沿方位向快变的情况，无法实现方位参考信号的频繁更新，不能实现对方位向不同点目标的不同距离徙动曲线的压缩处理，从而造成目标散焦，因此不适用于GEO SAR。到目前为止，只有北航于泽等人提出过GEO SAR成像方法，但该方法只能用于子孔径数据处理，且只能实现单点目标成像，方位分辨率仅达到15米，存在一定的局限性。至今，国内外还未有学者提出用于GEO SAR全孔径数据处理的高分辨率成像方法。

发明内容

本发明的目的是为了解决上述问题，提出一种地球同步轨道合成孔径雷达卫星的高分辨率成像方法，该方法能够处理地球同步轨道(GEO)合成孔径雷达(SAR)卫星全孔径回波数据，实现广域场景内多点目标的高分辨率成像，方位分辨率可达到2米以下。

一种地球同步轨道合成孔径雷达卫星的高分辨率成像方法，包括以下几个步骤：

步骤一：子孔径划分；

全孔径回波数据是一个二维数组，大小为X×Y，一维是方位向，有X个采样点，另一维是距离向，有Y个采样点，表示雷达获取了X个方位时刻的一维脉冲回波数据，每个时刻的一维脉冲回波数据有Y个采样点，Y为2的整数次幂；

对全孔径回波数据在方位向进行分块，将其分成N等份，即子孔径划分，得到N块回波子孔径数据A₁、A₂、…、A_N，每块数据大小为(X/N)×Y；

步骤二：子孔径数据距离向压缩；

对每块回波子孔径数据A₁、A₂、…、A_N进行距离向匹配滤波，即距离向压缩，得到距离向压缩后的子孔径数据B₁、B₂、…、B_N；

步骤三：子孔径数据距离向插值；

对距离向压缩后的子孔径数据B₁、B₂、…、B_N做距离向快速傅里叶变换，变换到距离向频域，在频域的高频部分补零，补零的个数为(M-1)·Y，即在每个方位时刻的距离向第个和第

个采样点之间增加(M-1)·Y个采样点，每个采样点的值为0，之后做快速傅里叶逆变换，实现距离向压缩后数据的M倍插值，得到距离向插值后的子孔径数据C₁、C₂、…、C_N，每块数据大小为(X/N)×M·Y，其中，M取1、2、4、8、16、32、64、128、256或512，当M＝1时，为不插值；

步骤四：地面网格划分；

成像区域为一个大小为P×Q的矩阵，其中P＜X，Q＜Y，方位向有P个采样点，设其P个方位时刻为

对应全孔径回波数据方位向最中心的P个时刻，即

其中，k＝1、2、…、P，t₁、t₂、…、t_X为全孔径回波数据的X个方位时刻，

为全孔径回波数据的方位中心时刻，当X、P为奇数时，

和按四舍五入取整；成像区域距离向有Q个采样点，设其Q个斜距分别为

对应全孔径回波数据距离向最中心的Q个斜距，即

其中，l＝1、2、…、Q，R₁、R₂、…、R_Y为全孔径回波数据的Y个斜距，

为全孔径回波数据距离向的中心斜距，当Q为奇数时，按四舍五入取整；对于成像区域的每个方位时刻，得到成像区域上该时刻距离向每个采样点按照斜距对应地面上目标点的位置，方位向第k个、距离向第l个采样点对应地面目标点的位置坐标表示为(x_gt(k，l)，y_gt(k，l)，z_gt(k，l))，每个网格为地面上一个目标点，对应成像区域上一个采样点，最后得到所有P×Q个目标点的坐标，完成地面网格划分；

步骤五：确定距离徙动曲线；

步骤四得到成像区域每个采样点对应地面目标点的位置，对于每个采样点，确定其在距离向插值后的每个子孔径数据C₁、C₂、…、C_N中的距离徙动曲线，即求出在子孔径数据中每个方位时刻雷达天线相位中心与该采样点对应地面目标点的斜距，确定斜距所在的子孔径数据中的距离向采样点，确定的每个时刻的距离向采样点构成距离徙动曲线，当成像区域所有采样点的距离徙动曲线都确定后，完成确定距离徙动曲线；

步骤六：进行方位压缩处理；

对成像区域每个采样点进行方位压缩处理；

具体为：对当前成像区域的采样点在每个子孔径数据C₁、C₂、…、C_N中的距离徙动曲线上的采样点进行相位补偿，然后相加，将相加后的数值作为当前成像区域采样点的值，表示为：

$f_i=\underset{j=1}{\overset{X/N}{\mathrm{\Σ}}}S\left(t_{\mathrm{ij}}\right)\·\exp \left\{j\frac{4\mathrm{\π}}{\mathrm{\λ}}R\left(t_{\mathrm{ij}}\right)\right\}$

式中，S(t_ij)表示当前成像区域的采样点在第i个子孔径数据上的距离徙动曲线上时刻为t_ij的采样点的值，

表示参考相位，λ表示雷达波长，f_i表示通过第i个子孔径数据得到的当前成像区域采样点的值，i＝1、2、…、N；

当成像区域的所有采样点都进行完压缩处理，完成方位压缩处理，得到N个子孔径图像I₁、I₂、…、I_N，I_i表示通过第i个子孔径数据得到的成像区域数据，i＝1、2、…、N；

步骤七：子孔径图像叠加，得到最终的图像；

将步骤六中得到的N个子孔径图像相加，即N个子孔径图像上对应采样点的值相加，实现全孔径数据压缩处理，得到最终的图像。

本发明的优点在于：

(1)本发明通过处理全孔径回波数据，可实现方位分辨率优于2米的成像；

(2)本发明提出的成像方法可实现多点目标成像，即对不同距离向和方位向目标的不同距离徙动曲线均可实现方位压缩；

(3)本发明提出的成像方法可实现广域场景范围内的成像，具有成像区域自由化的特点，即可根据需求对不同位置的场景实现成像。

附图说明

图1是本发明的为成像方法流程图；

图2是本发明的为成像区域网格划分的示意图；

图3是本发明的为点阵目标成像结果。

具体实施方式

下面将结合附图和实施例对本发明作进一步的详细说明。

本发明是一种地球同步轨道合成孔径雷达卫星的高分辨率成像方法，处理的对象是雷达的全孔径回波数据，得到的结果是一幅高分辨率图像。全孔径回波数据是一个二维数组(矩阵)，设其大小为X×Y，一维是方位向，有X个采样点，另一维是距离向，有Y个采样点，表示雷达获取了X个方位时刻的一维脉冲回波数据，每个时刻的一维脉冲回波数据有Y个采样点。Y为2的整数次幂，这样保证距离向数据可直接进行快速傅里叶变换(FFT)。不同的方位向采样点对应不同的方位时刻，两相邻采样点相差方位时刻1/prf，prf表示脉冲重复频率；不同的距离向采样点对应不同的斜距(雷达天线相位中心到地面目标点的距离)，两相邻采样点相差斜距c/2f_s，c表示光速，f_s表示采样率。

方法的流程如图1所示，包括以下几个步骤：

步骤一：子孔径划分；

由于全孔径回波数据量极大，现有计算机不能一次全部处理，并且也为提高处理速度，对全孔径回波数据在方位向进行分块，将其分成N等份，即子孔径划分，得到N块回波子孔径数据A₁、A₂、…、A_N，每块数据大小为(X/N)×Y；划分的子孔径数据的数据量要小于数据处理程序可开辟的最大内存：若采用处理器是32位的计算机，或者计算机的操作系统是32位的操作系统，或者采用Microsoft Visual C++6.0软件编程，一个程序一次开辟的内存最多为2GB，子孔径数据量要小于2GB；当采用处理器是64位的计算机，并且计算机的操作系统是64位的操作系统，并且采用Microsoft Visual Studio软件编程，一个程序一次可开辟大于2GB的内存，这由操作系统决定，子孔径数据量小于可开辟的最大内存即可。

步骤二：子孔径数据距离向压缩；

对每块回波子孔径数据A₁、A₂、…、A_N进行距离向匹配滤波，即距离向压缩，得到距离向压缩后的子孔径数据B₁、B₂、…、B_N。

具体为：对回波子孔径数据A₁、A₂、…、A_N沿距离向做快速傅里叶变换(FFT)，即对回波子孔径数据的每个方位时刻的距离向数据做快速傅里叶变换，将其变换到距离向频域，构造雷达发射的线性调频信号作为参考信号，参考信号为一维离散信号，其采样点个数为n₀，n₀＜Y且取n₀为偶数，在参考信号的前后各增加

个值为0的采样点，使得参考信号的采样点数为Y，对其做快速傅里叶变换将其变换到频域，在频域内每个方位时刻的距离向数据与参考信号的共轭相乘，将得到的结果进行快速傅里叶逆变换(IFFT)，再将距离向数据的前

个点和后

个点对换。最后得到距离向压缩后的子孔径数据B₁、B₂、…、B_N。

步骤三：子孔径数据距离向插值；

对距离向压缩后的子孔径数据B₁、B₂、…、B_N做距离向快速傅里叶变换，变换到距离向频域，在频域的高频部分补零，补零的个数为(M-1)·Y，即在每个方位时刻的距离向第

个和第

个采样点之间增加(M-1)·Y个采样点，每个采样点的值为0，之后做快速傅里叶逆变换，实现距离向压缩后数据的M倍插值，得到距离向插值后的子孔径数据C₁、C₂、…、C_N，每块数据大小为(X/N)×M·Y。其中，M取1、2、4、8、16、32、64、128、256或512，当M＝1时，为不插值。

步骤四：地面网格划分；

如图2所示，设成像区域为一个大小为P×Q(P＜X，Q＜Y)的矩阵，方位向有P个采样点，设其P个方位时刻为

对应全孔径回波数据方位向最中心的P个时刻，即

其中，k＝1、2、…、P，t₁、t₂、…、t_X为全孔径回波数据的X个方位时刻，

为全孔径回波数据的方位中心时刻，当X、P为奇数时，

和

按四舍五入取整。成像区域距离向有Q个采样点，设其Q个斜距分别为

对应全孔径回波数据距离向最中心的Q个斜距，即

其中，l＝1、2、…、Q，R₁、R₂、…、R_Y为全孔径回波数据的Y个斜距，

为全孔径回波数据距离向的中心斜距，当Q为奇数时，

按四舍五入取整。对于成像区域的每个方位时刻，得到成像区域上该时刻距离向每个采样点按照斜距对应地面上目标点的位置(目标点在转动地心坐标系下的坐标)，方位向第k个、距离向第l个采样点对应地面目标点的位置坐标表示为(x_gt(k，l)，y_gt(k，l)，z_gt(k，l))。上述过程也就是按成像区域的采样点对地面场景进行网格划分，每个网格为地面上一个目标点，对应成像区域上一个采样点，最后得到所有P×Q个目标点的坐标，完成地面网格划分。

获取目标点位置坐标(x_gt(k，l)，y_gt(k，l)，z_gt(k，l))的具体方法为：

表示成像区域的一个方位时刻，将该时刻卫星在不转动地心坐标系中的坐标

经矩阵相乘转换到转动地心坐标系下，得到卫星位置矢量如式(1)所示，

$\stackrel{\→}{A_s}=A_{\mathrm{go}}\left(\begin{array}{c}{x}_{\mathrm{os}}\left({t_C}_k\right)\\ y_{\mathrm{os}}\left({t_C}_k\right)\\ z_{\mathrm{os}}\left({t_C}_k\right)\end{array}\right)---\left(1\right)$

$A_{\mathrm{go}}=\left(\begin{array}{ccc}\cos H_G& \sin H_G& 0\\ -\sin H_G& \cos H_G& 0\\ 0& 0& 1\end{array}\right)---\left(2\right)$

式(1)中，A_go表示由不转动地心坐标系到转动地心坐标系的转换矩阵，式(2)中，

HG表示格林威治时角，ω_e是地球自转角速度，t₀是格林威治时间。上述涉及的两个坐标系定义如下：不转动地心坐标系原点位于地心，X轴在赤道平面内指向春分点，Z轴沿地球的自转轴指向正北极，Y轴在赤道平面内使该坐标系构成右手直角坐标系。转动地心坐标系由不转动地心坐标系绕Z轴逆时针旋转一个角H_G得到。

将天线相位中心在卫星星体坐标系下的坐标(x_e，y_e，z_e)经矩阵相乘转换到转动地心坐标系下得到天线偏置矢量

如式(3)所示，

$\stackrel{\→}{A_e}=A_{\mathrm{go}}{A}_{\mathrm{ov}}{A}_{\mathrm{vr}}{A}_{\mathrm{re}}\left(\begin{array}{c}{X}_e\\ y_e\\ z_e\end{array}\right)---\left(3\right)$

$A_{\mathrm{re}}=\left(\begin{array}{ccc}\cos {\mathrm{\θ}}_y& 0& -\sin {\mathrm{\θ}}_y\\ 0& 1& 0\\ \sin {\mathrm{\θ}}_y& 0& \cos {\mathrm{\θ}}_y\end{array}\right)\left(\begin{array}{ccc}\cos {\mathrm{\θ}}_p& -\sin {\mathrm{\θ}}_p& 0\\ \sin {\mathrm{\θ}}_p& \cos {\mathrm{\θ}}_p& 0\\ 0& 0& 1\end{array}\right)\left(\begin{array}{ccc}1& 0& 0\\ 0& \cos {\mathrm{\θ}}_r& -\sin {\mathrm{\θ}}_r\\ 0& \sin {\mathrm{\θ}}_r& \cos {\mathrm{\θ}}_r\end{array}\right)---\left(4\right)$

${\mathrm{\θ}}_y=a\tan \left[\mathrm{\ϵ}\frac{\cos (\mathrm{\θ}+\mathrm{\ω})\sin \mathrm{\η}}{n/{\mathrm{\ω}}_e-\cos \mathrm{\η}}\right]$ 其中， $n=\sqrt{\frac{\mathrm{\μ}}{a^3}}---\left(5\right)$

$A_{\mathrm{vr}}=\left(\begin{array}{ccc}-\sin (\mathrm{\θ}-\mathrm{\γ})& -\cos (\mathrm{\θ}-\mathrm{\γ})& 0\\ \cos (\mathrm{\θ}-\mathrm{\γ})& -\sin (\mathrm{\θ}-\mathrm{\γ})& 0\\ 0& 0& 1\end{array}\right)---\left(6\right)$

$\mathrm{\γ}=a\tan \left(\frac{e\sin \mathrm{\θ}}{1+e\cos \mathrm{\θ}}\right)---\left(7\right)$

$A_{\mathrm{ov}}=\left(\begin{array}{ccc}\cos \mathrm{\Ω}& -\sin \mathrm{\Ω}& 0\\ \sin \mathrm{\Ω}& \cos \mathrm{\Ω}& 0\\ 0& 0& 1\end{array}\right)\left(\begin{array}{ccc}1& 0& 0\\ 0& \cos \mathrm{\η}& -\sin \mathrm{\η}\\ 0& \sin \mathrm{\η}& \cos \mathrm{\η}\end{array}\right)\left(\begin{array}{ccc}\cos \mathrm{\ω}& -\sin \mathrm{\ω}& 0\\ \sin \mathrm{\ω}& \cos \mathrm{\ω}& 0\\ 0& 0& 1\end{array}\right)---\left(8\right)$

式(3)中，A_re表示由卫星星体坐标系到卫星平台坐标系的转换矩阵，如式(4)所示，θ_y、θ_p和θ_r分别表示卫星姿态偏航角、俯仰角和横滚角，其中偏航角由式(5)给出，表明卫星采用偏航控制，式(5)中，ε为卫星左视或右视的系数，ε＝1为右视，ε＝-1为左视，n为卫星平均角速率，μ为引力常数，a为卫星轨道半长轴，ω_e为地球自转角速度，ω为近地点幅角，η为轨道倾角，θ是当前方位时刻

卫星的真近心角；A_vr表示由卫星平台坐标系到轨道平面坐标系的转换矩阵，如式(6)所示，γ是当前方位时刻

卫星的航迹角，γ由式(7)给出，e表示偏心率；A_ov表示由轨道平面坐标系到不转动地心坐标系的转换矩阵，如式(8)所示，Ω为升交点赤经。上述涉及的三个坐标系定义如下：轨道平面坐标系由不转动地心坐标系经三次旋转得到，第一次将不转动地心坐标系绕Z轴逆时针旋转一个角Ω，第二次将得到的坐标系绕X轴逆时针旋转一个角η，第三次再将得到的坐标系绕Z轴逆时针旋转一个角ω，最后得到轨道平面坐标系。卫星平台坐标系由轨道平面坐标系绕Z轴顺时针旋转一个角度270°-θ+γ得到，且原点由地心平移至卫星平台质心。卫星星体坐标系由卫星平台坐标系经三次旋转得到，第一次将卫星平台坐标系绕Y轴顺时针旋转一个角度θ_y，第二次将得到的坐标系绕Z轴逆时针旋转一个角度θ_p，第三次将得到的坐标系绕X轴逆时针旋转一个角度θ_r。

设雷达天线相位中心到天线波束中心线指在地面上目标点(天线波束中心线与地球表面的交点)的斜距(距离)为R，则在天线坐标系下斜距矢量为(0，R，0)，经矩阵相乘，在转动地心坐标系下斜距矢量表示为：

$\stackrel{\→}{A_r}=A_{\mathrm{go}}{A}_{\mathrm{ov}}{A}_{\mathrm{vr}}{A}_{\mathrm{re}}{A}_{\mathrm{ea}}\left(\begin{array}{c}0\\ R\\ 0\end{array}\right)---\left(9\right)$

$A_{\mathrm{ea}}=\left(\begin{array}{ccc}1& 0& 0\\ 0& \cos {\mathrm{\θ}}_L& \sin {\mathrm{\θ}}_L\\ 0& -\sin {\mathrm{\θ}}_L& \cos {\mathrm{\θ}}_L\end{array}\right)---\left(10\right)$

式(9)中，A_ea表示由天线坐标系到卫星星体坐标系的转换矩阵，如式(10)所示，θ_L是雷达波束中心视角。天线坐标系由卫星星体坐标系绕X轴顺时针旋转一个角度θ_L得到，且原点由卫星平台质心平移到天线相位中心。

根据矢量和原理，得到天线波束中心线指在地面的目标点在转动地心坐标系下坐标为

$\left(\begin{array}{c}{x}_{\mathrm{gt}}\\ y_{\mathrm{gt}}\\ z_{\mathrm{gt}}\end{array}\right)=\stackrel{\→}{A_s}+\stackrel{\→}{A_e}+\stackrel{\→}{A_r}---\left(11\right)$

(x_gt，y_gt，z_gt)满足地球椭球方程，如式(12)所示，E_a表示地球半长轴，E_b表示地球半短轴。

$\frac{x^2}{E_a^2}+\frac{y^2}{E_a^2}+\frac{z^2}{E_b^2}=1---\left(12\right)$

式(11)代入式(12)，解方程求得斜距R，将R代入式(11)求得坐标(x_gt，y_gt，z_gt)。雷达波束中心视角θ_L不同，求出的R及坐标就不同，改变不同的雷达波束中心视角θ_L，求出一组不同的斜距R时目标点的坐标，在其中选出斜距为

时的目标点的坐标作为成像区域上

时刻距离向每个采样点对应地面目标点的坐标，记为(x_gt(k，l)，y_gt(k，l)，z_gt(k，l))，l＝1、2、…、Q。

获取卫星在每个

时刻斜距为的目标点的坐标，最后得到成像区域内P×Q个采样点对应地面目标点的坐标(x_gt(k，l)，y_gt(k，l)，z_gt(k，l))，k＝1、2、…、P，l＝1、2、…、Q。

步骤五：确定距离徙动曲线；

步骤四得到了成像区域每个采样点对应地面目标点的位置，对于每个采样点，确定其在距离向插值后的每个子孔径数据C₁、C₂、…、C_N中的距离徙动曲线，也就是求出在子孔径数据中每个方位时刻雷达天线相位中心与该采样点对应地面目标点的斜距(距离)，找到斜距所在的子孔径数据中的距离向采样点，找到的每个时刻的距离向采样点构成距离徙动曲线。当成像区域所有采样点的距离徙动曲线都确定后，完成此步骤。

确定成像区域上方位向第k个、距离向第l个采样点在每个子孔径数据C₁、C₂、…、C_N中的距离徙动曲线具体为：

用t_ij表示第i个子孔径数据的第j个方位采样时刻，i＝1，2，…，N，j＝1，2，…，X/N。如式(13)所示，首先，此时刻卫星在不转动地心坐标系下的坐标(x_os(t_ij)，y_os(t_ij)，z_os(t_ij))经矩阵相乘得到转动地心坐标系下的位置矢量，转换矩阵A_go如式(2)所示，其中H_G＝ω_e(t_ij-t₀)；然后，天线相位中心在卫星星体坐标系下的坐标(x_e，y_e，z_e)经矩阵相乘得到转动地心坐标系下的天线偏置矢量，转换矩阵如式(4)～式(8)所示，其中转换矩阵A_vr中θ为当前方位时刻t_ij卫星的真近心角，γ为当前方位时刻t_ij卫星的航迹角；最后，将卫星位置矢量与天线偏置矢量求和，再减去当前成像区域的采样点对应地面目标点的坐标(x_gt(k，l)，y_gt(k，l)，z_gt(k，l))，得到斜距矢量

$\stackrel{\→}{A_r}=A_{\mathrm{go}}\left(\begin{array}{c}{x}_{\mathrm{os}}\left(t_{\mathrm{ij}}\right)\\ y_{\mathrm{os}}\left(t_{\mathrm{ij}}\right)\\ z_{\mathrm{os}}\left(t_{\mathrm{ij}}\right)\end{array}\right)+A_{\mathrm{go}}{A}_{\mathrm{ov}}{A}_{\mathrm{vr}}{A}_{\mathrm{re}}\left(\begin{array}{c}{x}_e\\ y_e\\ z_e\end{array}\right)-\left(\begin{array}{c}{x}_{\mathrm{gt}}(k,l)\\ y_{\mathrm{gt}}(k,l)\\ z_{\mathrm{gt}}(k,l)\end{array}\right)---\left(13\right)$

斜距为

取模，即

R(t_ij)表示时刻t_ij的斜距，找到此方位时刻斜距为R(t_ij)的距离向采样点。将j取遍所有方位时刻，找到的这些采样点构成距离徙动曲线，确定了当前成像区域的采样点在第i个子孔径数据中的距离徙动曲线。将i取遍所有子孔径数据，就确定了当前成像区域的采样点在每个子孔径数据中的距离徙动曲线。

对成像区域的每个采样点都确定其在每个子孔径数据C₁、C₂、…、C_N中的距离徙动曲线，就完成了此步骤。

步骤六：进行方位压缩处理；

对成像区域每个采样点进行方位压缩处理。

具体为：对当前成像区域的采样点在每个子孔径数据C₁、C₂、…、C_N中的距离徙动曲线上的采样点进行相位补偿(采样点的值与参考相位相乘)，然后相加，将相加后的数值作为当前成像区域采样点的值。表示为：

$f_i=\underset{j=1}{\overset{X/N}{\mathrm{\Σ}}}S\left(t_{\mathrm{ij}}\right)\·\exp \left\{j\frac{4\mathrm{\π}}{\mathrm{\λ}}R\left(t_{\mathrm{ij}}\right)\right\}---\left(14\right)$

式(14)中，S(t_ij)表示当前成像区域的采样点在第i个子孔径数据上的距离徙动曲线上时刻为t_ij的采样点的值，

表示参考相位，λ表示雷达波长，f_i表示通过第i个子孔径数据得到的当前成像区域采样点的值，i＝1、2、…、N。

当成像区域的所有采样点都进行完压缩处理，完成此步骤，得到N个子孔径图像I₁、I₂、…、I_N，I_i(i＝1、2、…、N)表示通过第i个子孔径数据得到的成像区域数据。

步骤七：子孔径图像叠加，得到最终的图像。

将步骤六中得到的N个子孔径图像相加，即N个子孔径图像上对应采样点的值相加，实现全孔径数据压缩处理，得到最终的图像。

实施例：

卫星轨道参数、卫星姿态参数和雷达参数如表1、表2和表3所示。

表1卫星轨道参数

表2卫星姿态参数

表3雷达参数

按照表1至表3的参数分别对位于赤道附近的地面上500km×500km范围内五个场景的3×3点阵目标进行了回波仿真，得到五个全孔径回波数据，利用表1至表3的参数分别对每个数据按以下步骤做成像处理。

步骤一：子孔径划分；

仿真中的一个全孔径回波数据大小为450000×4096，方位向450000个点，距离向4096(2¹²)个点，数据量为13.733GB。采用装有windows xp操作系统(32位系统)的计算机处理，将全孔径回波数据分成30块子孔径数据A₁、A₂、…、A₃₀，每块数据大小为15000×4096，数据量为0.458GB。

步骤二：子孔径数据距离向压缩；

对每块回波子孔径数据A₁、A₂、…、A₃₀进行距离向匹配滤波，即距离向压缩，得到距离向压缩后的子孔径数据B₁、B₂、…、B₃₀。

具体为：对回波子孔径数据沿距离向做快速傅里叶变换，即对回波子孔径数据的每个方位时刻的距离向数据做快速傅里叶变换，将其变换到距离向频域，构造仿真中雷达发射的线性调频信号作为参考信号，如式(15)所示，

s＝exp(-jπbτ²)    (15)

式(15)中，s表示参考信号，b表示调频率，为b＝Bw/τ₀，τ是离散时间变量，取值范围为

步长为

表示τ的取值从

开始，每次增长

直到

共200个采样点，在其前后各增加1948个值为0的采样点，使得参考信号的采样点数为4096，对其做快速傅里叶变换将其变换到频域，在频域内每个方位时刻的距离向数据与参考信号的共轭相乘，将得到的结果进行快速傅里叶逆变换，再将距离向数据的前2048个点和后2048个点对换。最后得到距离向压缩后的子孔径数据。

步骤三：子孔径数据距离向插值；

对距离向压缩后的子孔径数据B₁、B₂、…、B₃₀做距离向快速傅里叶变换，变换到距离向频域，在频域的高频部分补零，补零的个数为3×4096，即在每个方位时刻的距离向第2048个和第2049个采样点之间增加3×4096个采样点，每个采样点的值为0，之后做快速傅里叶逆变换，实现距离向压缩后数据的4倍插值，得到距离向插值后的子孔径数据C₁、C₂、…、C₃₀，每块数据大小为15000×16384。

步骤四：地面网格划分；

定义一个大小为512×256的矩阵，作为成像区域，方位向有512个采样点，其512个方位时刻为对应全孔径回波数据方位向最中心的512个时刻，即

其中，k＝1、2、…、512，t₁、t₂、…、t₄₅₀₀₀₀为全孔径回波数据的450000个方位时刻，距离向有256个采样点，其256个斜距分别为

对应全孔径回波数据距离向最中心的256个斜距，即

其中，l＝1、2、…、256，R₁、R₂、…、R₄₀₉₆为全孔径回波数据的4096个斜距。对于成像区域的每个方位时刻，得到成像区域上该时刻距离向每个采样点按照斜距对应地面上目标点的位置(目标点在转动地心坐标系下的坐标)，方位向第k个、距离向第l个采样点对应地面目标点的位置坐标表示为(x_gt(k，l)，y_gt(k，l)，z_gt(k，l))。求出所有512×256个目标点的坐标，即完成地面网格划分。

获取目标点位置坐标的具体方法如下：

(k＝1，2，…，512)表示成像区域的一个方位时刻，将该时刻卫星在不转动地心坐标系中的坐标

经矩阵相乘转换到转动地心坐标系下，得到卫星位置矢量

如式(1)所示，转换矩阵A_go如式(2)所示，其中，

将天线相位中心在卫星星体坐标系下的坐标(x_e，y_e，z_e)经矩阵相乘转换到转动地心坐标系下得到天线偏置矢量如式(3)所示，转换矩阵如式(4)～式(8)所示，其中，θ是当前方位时刻

卫星的真近心角，γ是当前方位时刻

卫星的航迹角，卫星为右视ε＝1。

设雷达天线相位中心到天线波束中心线指在地面上目标点的斜距(距离)为R，则在天线坐标系下斜距矢量为(0，R，0)，经矩阵相乘，得到在转动地心坐标系下斜距矢量

如式(9)所示。

根据矢量和原理，得到转动地心坐标系下目标点坐标，如式(11)所示，式(11)代入式(12)，解方程求得斜距R，将R代入式(11)求得坐标(x_gt，y_gt，z_gt)。改变不同的雷达波束中心视角，求出一组不同的斜距R时目标点的坐标，在其中选出斜距为

时的目标点的坐标作为成像区域上

时刻距离向每个采样点对应地面目标点的坐标，记为(x_gt(k，l)，y_gt(k，l)，z_gt(k，l))，l＝1、2、…、256。

获取卫星在每个

时刻斜距为

的目标点的坐标，最后得到成像区域内512×256个采样点对应地面目标点的坐标(x_gt(k，l)，y_gt(k，l)，z_gt(k，l))，k＝1、2、…、512，l＝1、2、…、256。

步骤五：确定距离徙动曲线；

步骤四已求出了成像区域每个采样点对应地面目标点的位置，对于每个采样点，确定其在距离向插值后的每个子孔径数据C₁、C₂、…、C₃₀中的距离徙动曲线。当成像区域所有采样点的距离徙动曲线都确定后，完成此步骤。

确定成像区域上方位向第k个、距离向第l个采样点在每个子孔径数据C₁、C₂、…、C₃₀中的距离徙动曲线具体为：

t_ij表示第i个子孔径数据的第j个方位采样时刻，i＝1，2，…，30，j＝1，2，…，15000。如式(13)所示，首先，此时刻卫星在不转动地心坐标系下的坐标(x_os(t_ij)，y_os(t_ij)，z_os(t_ij))经矩阵相乘得到转动地心坐标系下的位置矢量，转换矩阵A_go如式(2)所示，其中H_G＝ω_e(t_ij-t₀)；然后，天线相位中心在卫星星体坐标系下的坐标(x_e，y_e，z_e)经矩阵相乘得到转动地心坐标系下的天线偏置矢量，转换矩阵如式(4)～式(8)所示，其中θ为当前方位时刻t_ij卫星的真近心角，γ为当前方位时刻t_ij卫星的航迹角；最后，将卫星位置矢量与天线偏置矢量求和，再减去当前成像区域的采样点对应地面目标点的坐标(x_gt(k，l)，y_gt(k，l)，z_gt(k，l))，得到斜距矢量

。斜距等于

取模，即

R(t_ij)表示时刻t_ij的斜距，找到此方位时刻斜距为R(t_ij)的距离向采样点。将j取遍所有方位时刻，找到的这些采样点构成距离徙动曲线，即确定了当前成像区域的采样点在第i个子孔径数据中的距离徙动曲线。将i取遍所有子孔径数据，就确定了当前成像区域的采样点在每个子孔径数据中的距离徙动曲线。

对成像区域的每个采样点都确定其在每个子孔径数据C₁、C₂、…、C₃₀中的距离徙动曲线，就完成了此步骤。

步骤六：进行方位压缩处理；

对成像区域每个采样点按照下述具体过程进行方位压缩处理。

具体为：对当前成像区域的采样点在每个子孔径数据C₁、C₂、…、C₃₀中的距离徙动曲线上的采样点进行相位补偿(采样点的值与参考相位相乘)并相加，将相加后的数值作为当前成像区域采样点的值。表示为：

$f_i=\underset{j=1}{\overset{15000}{\mathrm{\Σ}}}S\left(t_{\mathrm{ij}}\right)\·\exp \left\{j\frac{4\mathrm{\π}}{\mathrm{\λ}}R\left(t_{\mathrm{ij}}\right)\right\}---\left(16\right)$

式(16)中，S(t_ij)表示当前成像区域的采样点在第i个子孔径数据上的距离徙动曲线上时刻为t_ij的采样点的值，

表示参考相位，λ表示雷达波长，f_i表示通过第i个子孔径数据得到的当前成像区域采样点的值，i＝1，2，…，30。

当成像区域的所有采样点都进行完压缩处理，完成此步骤，得到30个子孔径图像I₁、I₂、…、I₃₀，I_i(i＝1，2，…，30)表示通过第i个子孔径数据得到的成像区域数据。

步骤七：子孔径图像叠加，得到最终的图像。

将步骤六中得到的30个子孔径图像相加，即30个子孔径图像上对应采样点的值相加，实现全孔径数据压缩处理，得到最终的图像。

经过以上步骤的成像处理，分别对得到的五个最终图像上的采样点的值取模并以bmp格式输出，得到图3所示的结果。图3为在地面上500km×500km范围内五个场景的3×3点阵目标成像结果。其中，(c)为场景中心成像结果，(a)、(b)、(d)、(e)分别为处于成像范围四个角的成像结果，每个场景中心点在地面上的坐标标于图下方。每个场景内横向为距离向，点地面间距500m，纵向为方位向，点地面间距250m。

表3和表4是图3中点目标成像质量评估结果。表3给出了每个场景中9个点的分辨率、峰值旁瓣比和积分旁瓣比的均值，表4给出了标准差。其中，距离分辨率为斜距分辨率。

评估结果表明距离分辨率为3m，方位分辨率为1.3m。点目标峰值旁瓣比和积分旁瓣比都接近理论值。场景内指标标准差和场景间指标均值的差别都不大，说明成像一致性较好。

因此，本发明提出的GEO SAR成像方法是可以实现广域场景高分辨率成像的。

表3点目标质量评估结果——指标均值

表4点目标质量评估结果——指标标准差

Claims (10)

1.一种地球同步轨道合成孔径雷达卫星的高分辨率成像方法，其特征在于，包括以下几个步骤：

步骤一：子孔径划分；

全孔径回波数据是一个二维数组，大小为X×Y，一维是方位向，有X个采样点，另一维是距离向，有Y个采样点，表示雷达获取了X个方位时刻的一维脉冲回波数据，每个时刻的一维脉冲回波数据有Y个采样点，Y为2的整数次幂；

对全孔径回波数据在方位向进行分块，将其分成N等份，即子孔径划分，得到N块回波子孔径数据A₁、A₂、…、A_N，每块数据大小为(X/N)×Y；

步骤二：子孔径数据距离向压缩；

对每块回波子孔径数据A₁、A₂、…、A_N进行距离向匹配滤波，即距离向压缩，得到距离向压缩后的子孔径数据B₁、B₂、…、B_N；

步骤三：子孔径数据距离向插值；

对距离向压缩后的子孔径数据B₁、B₂、…、B_N做距离向快速傅里叶变换，变换到距离向频域，在频域的高频部分补零，补零的个数为(M-1)·Y，即在每个方位时刻的距离向第

个和第

个采样点之间增加(M-1)·Y个采样点，每个采样点的值为0，之后做快速傅里叶逆变换，实现距离向压缩后数据的M倍插值，得到距离向插值后的子孔径数据C₁、C₂、…、C_N，每块数据大小为(X/N)×M·Y，其中，M取1、2、4、8、16、32、64、128、256或512，当M＝1时，为不插值；

步骤四：地面网格划分；

成像区域为一个大小为P×Q的矩阵，其中P＜X，Q＜Y，方位向有P个采样点，设其P个方位时刻为

对应全孔径回波数据方位向最中心的P个时刻，即

其中，k＝1、2、…、P，t₁、t₂、…、t_X为全孔径回波数据的X个方位时刻，

为全孔径回波数据的方位中心时刻，当X、P为奇数时，和

按四舍五入取整；成像区域距离向有Q个采样点，设其Q个斜距分别为

对应全孔径回波数据距离向最中心的Q个斜距，即

其中，l＝1、2、…、Q，R₁、R₂、…、R_Y为全孔径回波数据的Y个斜距，

为全孔径回波数据距离向的中心斜距，当Q为奇数时，

按四舍五入取整；对于成像区域的每个方位时刻，得到成像区域上该时刻距离向每个采样点按照斜距对应地面上目标点的位置，方位向第k个、距离向第l个采样点对应地面目标点的位置坐标表示为(x_gt(k，l)，y_gt(k，l)，z_gt(k，l))，每个网格为地面上一个目标点，对应成像区域上一个采样点，最后得到所有P×Q个目标点的坐标，完成地面网格划分；

步骤五：确定距离徙动曲线；

步骤四得到成像区域每个采样点对应地面目标点的位置，对于每个采样点，确定其在距离向插值后的每个子孔径数据C₁、C₂、…、C_N中的距离徙动曲线，即求出在子孔径数据中每个方位时刻雷达天线相位中心与该采样点对应地面目标点的斜距，确定斜距所在的子孔径数据中的距离向采样点，确定的每个时刻的距离向采样点构成距离徙动曲线，当成像区域所有采样点的距离徙动曲线都确定后，完成确定距离徙动曲线；

步骤六：进行方位压缩处理；

对成像区域每个采样点进行方位压缩处理；

具体为：对当前成像区域的采样点在每个子孔径数据C₁、C₂、…、C_N中的距离徙动曲线上的采样点进行相位补偿，然后相加，将相加后的数值作为当前成像区域采样点的值，表示为：

$f_i=\underset{j=1}{\overset{X/N}{\mathrm{\Σ}}}S\left(t_{\mathrm{ij}}\right)\·\exp \left\{j\frac{4\mathrm{\π}}{\mathrm{\λ}}R\left(t_{\mathrm{ij}}\right)\right\}---\left(14\right)$

式(14)中，S(t_ij)表示当前成像区域的采样点在第i个子孔径数据上的距离徙动曲线上时刻为t_ij的采样点的值，表示参考相位，λ表示雷达波长，f_i表示通过第i个子孔径数据得到的当前成像区域采样点的值，i＝1、2、…、N；

当成像区域的所有采样点都进行完压缩处理，完成方位压缩处理，得到N个子孔径图像I₁、I₂、…、I_N，I_i表示通过第i个子孔径数据得到的成像区域数据，i＝1、2、…、N；

步骤七：子孔径图像叠加，得到最终的图像；

将步骤六中得到的N个子孔径图像相加，即N个子孔径图像上对应采样点的值相加，实现全孔径数据压缩处理，得到最终的图像。

2.根据权利要求1所述的一种地球同步轨道合成孔径雷达卫星的高分辨率成像方法，其特征在于，步骤一中所述子孔径划分，其子孔径数据的数据量要小于数据处理程序可开辟的最大内存。

3.根据权利要求1所述的一种地球同步轨道合成孔径雷达卫星的高分辨率成像方法，其特征在于，步骤二中所述距离向压缩具体为：对回波子孔径数据A₁、A₂、…、A_N沿距离向做快速傅里叶变换，即对回波子孔径数据的每个方位时刻的距离向数据做快速傅里叶变换，将其变换到距离向频域，构造雷达发射的线性调频信号作为参考信号，参考信号为一维离散信号，其采样点个数为n₀，n₀＜Y且取n₀为偶数，在参考信号的前后各增加个值为0的采样点，使得参考信号的采样点数为Y，对其做快速傅里叶变换将其变换到频域，在频域内每个方位时刻的距离向数据与参考信号的共轭相乘，将得到的结果进行快速傅里叶逆变换，再将距离向数据的前个点和后

个点对换，最后得到距离向压缩后的子孔径数据B₁、B₂、…、B_N。

4.根据权利要求1所述的一种地球同步轨道合成孔径雷达卫星的高分辨率成像方法，其特征在于，步骤四中所述获取目标点位置坐标(x_gt(k，l)，y_gt(k，l)，z_gt(k，l))的具体方法为：

表示成像区域的一个方位时刻，将该时刻卫星在不转动地心坐标系中的坐标

经矩阵相乘转换到转动地心坐标系下，得到卫星位置矢量

$\stackrel{\→}{A_s}=A_{\mathrm{go}}\left(\begin{array}{c}{x}_{\mathrm{os}}\left({t_C}_k\right)\\ y_{\mathrm{os}}\left({t_C}_k\right)\\ z_{\mathrm{os}}\left({t_C}_k\right)\end{array}\right)---\left(1\right)$

$A_{\mathrm{go}}=\left(\begin{array}{ccc}\cos H_G& \sin H_G& 0\\ -\sin H_G& \cos H_G& 0\\ 0& 0& 1\end{array}\right)---\left(2\right)$

式(1)中，A_go表示由不转动地心坐标系到转动地心坐标系的转换矩阵，式(2)中，

H_G表示格林威治时角，ω_e是地球自转角速度，t₀是格林威治时间；

将天线相位中心在卫星星体坐标系下的坐标(x_e，y_e，z_e)经矩阵相乘转换到转动地心坐标系下得到天线偏置矢量

$\stackrel{\→}{A_e}=A_{\mathrm{go}}{A}_{\mathrm{ov}}{A}_{\mathrm{vr}}{A}_{\mathrm{re}}\left(\begin{array}{c}{X}_e\\ y_e\\ z_e\end{array}\right)---\left(3\right)$

$A_{\mathrm{re}}=\left(\begin{array}{ccc}\cos {\mathrm{\θ}}_y& 0& -\sin {\mathrm{\θ}}_y\\ 0& 1& 0\\ \sin {\mathrm{\θ}}_y& 0& \cos {\mathrm{\θ}}_y\end{array}\right)\left(\begin{array}{ccc}\cos {\mathrm{\θ}}_p& -\sin {\mathrm{\θ}}_p& 0\\ \sin {\mathrm{\θ}}_p& \cos {\mathrm{\θ}}_p& 0\\ 0& 0& 1\end{array}\right)\left(\begin{array}{ccc}1& 0& 0\\ 0& \cos {\mathrm{\θ}}_r& -\sin {\mathrm{\θ}}_r\\ 0& \sin {\mathrm{\θ}}_r& \cos {\mathrm{\θ}}_r\end{array}\right)---\left(4\right)$

${\mathrm{\θ}}_y=a\tan \left[\mathrm{\ϵ}\frac{\cos (\mathrm{\θ}+\mathrm{\ω})\sin \mathrm{\η}}{n/{\mathrm{\ω}}_e-\cos \mathrm{\η}}\right]$ 其中， $n=\sqrt{\frac{\mathrm{\μ}}{a^3}}---\left(5\right)$

$A_{\mathrm{vr}}=\left(\begin{array}{ccc}-\sin (\mathrm{\θ}-\mathrm{\γ})& -\cos (\mathrm{\θ}-\mathrm{\γ})& 0\\ \cos (\mathrm{\θ}-\mathrm{\γ})& -\sin (\mathrm{\θ}-\mathrm{\γ})& 0\\ 0& 0& 1\end{array}\right)---\left(6\right)$

$\mathrm{\γ}=a\tan \left(\frac{e\sin \mathrm{\θ}}{1+e\cos \mathrm{\θ}}\right)---\left(7\right)$

$A_{\mathrm{ov}}=\left(\begin{array}{ccc}\cos \mathrm{\Ω}& -\sin \mathrm{\Ω}& 0\\ \sin \mathrm{\Ω}& \cos \mathrm{\Ω}& 0\\ 0& 0& 1\end{array}\right)\left(\begin{array}{ccc}1& 0& 0\\ 0& \cos \mathrm{\η}& -\sin \mathrm{\η}\\ 0& \sin \mathrm{\η}& \cos \mathrm{\η}\end{array}\right)\left(\begin{array}{ccc}\cos \mathrm{\ω}& -\sin \mathrm{\ω}& 0\\ \sin \mathrm{\ω}& \cos \mathrm{\ω}& 0\\ 0& 0& 1\end{array}\right)---\left(8\right)$

式(3)中，A_re表示由卫星星体坐标系到卫星平台坐标系的转换矩阵，如式(4)所示，θ_y、θ_p和θ_r分别表示卫星姿态偏航角、俯仰角和横滚角，其中偏航角由式(5)得到，表明卫星采用偏航控制，式(5)中，ε为卫星左视或右视的系数，ε＝1为右视，ε＝-1为左视，n为卫星平均角速率，μ为引力常数，a为卫星轨道半长轴，ω_e为地球自转角速度，ω为近地点幅角，η为轨道倾角，θ是当前方位时刻

卫星的真近心角；A_vr表示由卫星平台坐标系到轨道平面坐标系的转换矩阵，如式(6)所示，γ是当前方位时刻

卫星的航迹角，γ由式(7)给出，e表示偏心率；A_ov表示由轨道平面坐标系到不转动地心坐标系的转换矩阵，如式(8)所示，Ω为升交点赤经；

设雷达天线相位中心到天线波束中心线指在地面上目标点的斜距为R，所述的目标点为天线波束中心线与地球表面的交点，则在天线坐标系下斜距矢量为(0，R，0)，经矩阵相乘，在转动地心坐标系下斜距矢量表示为：

$\stackrel{\→}{A_r}=A_{\mathrm{go}}{A}_{\mathrm{ov}}{A}_{\mathrm{vr}}{A}_{\mathrm{re}}{A}_{\mathrm{ea}}\left(\begin{array}{c}0\\ R\\ 0\end{array}\right)---\left(9\right)$

$A_{\mathrm{ea}}=\left(\begin{array}{ccc}1& 0& 0\\ 0& \cos {\mathrm{\θ}}_L& \sin {\mathrm{\θ}}_L\\ 0& -\sin {\mathrm{\θ}}_L& \cos {\mathrm{\θ}}_L\end{array}\right)---\left(10\right)$

式(9)中，A_ea表示由天线坐标系到卫星星体坐标系的转换矩阵，如式(10)所示，θ_L是雷达波束中心视角；

得到天线波束中心线指在地面的目标点在转动地心坐标系下坐标为：

$\left(\begin{array}{c}{x}_{\mathrm{gt}}\\ y_{\mathrm{gt}}\\ z_{\mathrm{gt}}\end{array}\right)=\stackrel{\→}{A_s}+\stackrel{\→}{A_e}+\stackrel{\→}{A_r}---\left(11\right)$

(x_gt，y_gt，z_gt)满足地球椭球方程，如式(12)所示，E_a表示地球半长轴，E_b表示地球半短轴；

$\frac{x^2}{E_a^2}+\frac{y^2}{E_a^2}+\frac{z^2}{E_b^2}=1---\left(12\right)$

式(11)代入式(12)，解方程求得斜距R，将R代入式(11)求得坐标(x_gt，y_gt，z_gt)；

改变不同的雷达波束中心视角θ_L，求出一组不同的斜距R时目标点的坐标，在其中选出斜距为

时的目标点的坐标作为成像区域上

时刻距离向每个采样点的对应地面目标点的坐标，记为(x_gt(k，l)，y_gt(k，l)，z_gt(k，l))，l＝1、2、…、Q；

获取卫星在每个

时刻斜距为

的目标点的坐标，最后得到成像区域内P×Q个采样点对应地面目标点的坐标(x_gt(k，l)，y_gt(k，l)，z_gt(k，l))，k＝1、2、…、P，l＝1、2、…、Q。

5.根据权利要求1所述的一种地球同步轨道合成孔径雷达卫星的高分辨率成像方法，其特征在于，步骤四中所述的不转动地心坐标系，其原点位于地心，X轴在赤道平面内指向春分点，Z轴沿地球的自转轴指向正北极，Y轴在赤道平面内使该坐标系构成右手直角坐标系；所述的转动地心坐标系，是由不转动地心坐标系绕Z轴逆时针旋转一个角HG得到。

6.根据权利要求1所述的一种地球同步轨道合成孔径雷达卫星的高分辨率成像方法，其特征在于，步骤四中所述的轨道平面坐标系，由不转动地心坐标系经三次旋转得到，第一次将不转动地心坐标系绕Z轴逆时针旋转一个角Ω，第二次将得到的坐标系绕X轴逆时针旋转一个角η，第三次再将得到的坐标系绕Z轴逆时针旋转一个角ω，最后得到轨道平面坐标系。

7.根据权利要求1所述的一种地球同步轨道合成孔径雷达卫星的高分辨率成像方法，其特征在于，步骤四中所述的卫星平台坐标系，是由轨道平面坐标系绕Z轴顺时针旋转一个角度270°-θ+γ得到，且原点由地心平移至卫星平台质心。

8.根据权利要求1所述的一种地球同步轨道合成孔径雷达卫星的高分辨率成像方法，其特征在于，步骤四中所述的卫星星体坐标系，是由卫星平台坐标系经三次旋转得到，第一次将卫星平台坐标系绕Y轴顺时针旋转一个角度θ_y，第二次将得到的坐标系绕Z轴逆时针旋转一个角度θ_p，第三次将得到的坐标系绕X轴逆时针旋转一个角度θ_r；

所述的天线坐标系由卫星星体坐标系绕X轴顺时针旋转一个角度θ_L得到，且原点由卫星平台质心平移到天线相位中心。

9.根据权利要求1所述的一种地球同步轨道合成孔径雷达卫星的高分辨率成像方法，其特征在于，步骤五中所述的确定成像区域上方位向第k个、距离向第l个采样点在每个子孔径数据C₁、C₂、…、C_N中的距离徙动曲线具体方法为：

用t_ij表示第i个子孔径数据的第j个方位采样时刻，i＝1，2，…，N，j＝1，2，…，X/N；首先，此时刻卫星在不转动地心坐标系下的坐标(x_os(t_ij)，y_os(t_ij)，z_os(t_ij))经矩阵相乘得到转动地心坐标系下的位置矢量，转换矩阵A_go如式(2)所示，其中H_G＝ω_e(t_ij-t₀)；然后，天线相位中心在卫星星体坐标系下的坐标(x_e，y_e，z_e)经矩阵相乘得到转动地心坐标系下的天线偏置矢量，转换矩阵如式(4)～式(8)所示，其中转换矩阵A_vr中θ为当前方位时刻t_ij卫星的真近心角，γ为当前方位时刻t_ij卫星的航迹角；最后，将卫星位置矢量与天线偏置矢量求和，再减去当前成像区域的采样点对应地面目标点的坐标(x_gt(k，l)，y_gt(k，l)，z_gt(k，l))，得到斜距矢量

$\stackrel{\→}{A_r}=A_{\mathrm{go}}\left(\begin{array}{c}{x}_{\mathrm{os}}\left(t_{\mathrm{ij}}\right)\\ y_{\mathrm{os}}\left(t_{\mathrm{ij}}\right)\\ z_{\mathrm{os}}\left(t_{\mathrm{ij}}\right)\end{array}\right)+A_{\mathrm{go}}{A}_{\mathrm{ov}}{A}_{\mathrm{vr}}{A}_{\mathrm{re}}\left(\begin{array}{c}{x}_e\\ y_e\\ z_e\end{array}\right)-\left(\begin{array}{c}{x}_{\mathrm{gt}}(k,l)\\ y_{\mathrm{gt}}(k,l)\\ z_{\mathrm{gt}}(k,l)\end{array}\right)---\left(13\right)$

斜距为

取模，即R(t_ij)表示时刻t_ij的斜距，找到此方位时刻斜距为R(t_ij)的距离向采样点；将j取遍所有方位时刻，确定的采样点构成距离徙动曲线，确定了当前成像区域的采样点在第i个子孔径数据中的距离徙动曲线；将i取遍所有子孔径数据，就确定了当前成像区域的采样点在每个子孔径数据中的距离徙动曲线；

对成像区域的每个采样点都确定其在每个子孔径数据C₁、C₂、…、C_N中的距离徙动曲线，就完成了确定距离徙动曲线。

10.根据权利要求1所述的一种地球同步轨道合成孔径雷达卫星的高分辨率成像方法，其特征在于，步骤六中所述的采样点进行相位补偿，是采样点的值与参考相位相乘。

Images