基于DEM和DLCD的机载PD雷达高保真非均匀地杂波仿真方法与流程

本发明属于机载脉冲多普勒(pd)雷达非均匀地杂波仿真技术领域，特别是涉及一种基于数字高程模型(dem)和数字地表分类数据(dlcd)的机载pd雷达高保真非均匀地杂波仿真方法。

背景技术

与地基雷达相比，机载雷达具有良好的视野和机动性，较强的生存能力，能够发现超低空飞行目标和地面运动目标。当机载pd雷达工作在下视模式时，有用信号会被强地杂波覆盖，因此地杂波的抑制效果将直接影响目标检测和参数估计的性能。当前，机载雷达地杂波抑制主要采用自适应处理技术，建立在均匀杂波基础上的自适应处理技术需要满足rmb准则来估计杂波协方差矩阵，但机载雷达往往工作在非均匀杂波环境中，导致无法获得足够的独立同分布(independentandidenticallydistributed，iid)训练样本来估计杂波协方差矩阵，从而严重影响了机载雷达的杂波抑制性能。在这种情况下，2002年美国国防高级研究计划署(defenseadvancedresearchprojectsagency,darpa)启动了“知识辅助的传感器信号处理和专家咨询计划(knowledgeaidedsensorsignalprocessingandexpertreasoning,kassper)”，研究利用先验地理信息提升雷达目标检测和参数估计的性能，然而在实际测量中，为了获得含有丰富先验信息的实采地杂波数据，需要耗费大量的人力和物力，所需成本过高。因此，构建符合真实地理环境的高保真非均匀地杂波模型，对于研究提升机载pd雷达在非均匀地杂波环境下的目标检测和参数估计性能具有十分重要的意义。

根据地杂波仿真过程中是否引入实际场景的地理信息，可以将地杂波仿真分为基于统计模型的地杂波仿真和地理信息辅助的地杂波仿真。其中，基于统计模型的地杂波仿真假设雷达作用区域为平面，不考虑实际复杂多变的地形和地理地貌特征，多用于实现满足一定概率分布的均匀地杂波仿真。地理信息辅助的地杂波仿真利用实际的地面场景信息(如dem数据，dlcd等)计算杂波单元的后向散射系数，用于仿真非均匀地杂波，可以为机载pd雷达提供有效的先验地理信息，进而更好地验证地杂波抑制算法的性能。

围绕地理信息辅助的非均匀地杂波建模，研究者开展了大量的研究。范国忠等人提出了一种自然场景的机载pd雷达地杂波仿真方法，从dem数据出发，对各单元进行遮蔽判断，并利用pd雷达的回波信号模型仿真地杂波。洪丽娜同样利用dem数据确定杂波散射单元后向散射系数，并计算雷达散射截面积，实现了单一地貌的地杂波仿真。dlhellard等人利用数字地形和地物信息计算镜面反射系数，并根据镜面反射系数计算特定场景的雷达散射截面，最后结合雷达方程实现地杂波仿真。kurekina等人提出利用多模式遥感影像和数字地形数据生成陆基雷达特定杂波图的新方法，该方法基于陆基雷达构建后向散射系数模型，不适合描述机载雷达地杂波特性。饶妮妮等人针对机载认知雷达，提出了一种基于数字高程模型数据的机载认知雷达地杂波建模方法，该方法充分利用dem数据提取地形因子，并结合地理地貌学理论构建地杂波，具有一定的实际价值。上述方法将先验地理信息数据融合到地杂波仿真中，取得了较为真实的地杂波仿真结果，但上述方法采用的后向散射系数模型较为简单，并且仅从dem数据中计算地形因子确定地貌类型，数据分辨率不高，地貌分类比较粗糙，不能很好地对非均匀的真实地形进行高保真拟合。

技术实现要素：

为了解决上述问题，本发明的目的在于提供一种基于dem和dlcd的机载pd雷达高保真非均匀地杂波仿真方法。

为了达到上述目的，本发明提供的基于dem和dlcd的机载pd雷达高保真非均匀地杂波仿真方法包括按顺序进行的下列步骤：

1)截取仿真所用的dem数据和dlcd；

2)利用dem数据提取高程值，并计算每个地表散射点相对于载机的斜距、方位角、俯仰角、擦地角在内的信息；

3)根据步骤2)获得的擦地角和斜距信息对地表散射点进行遮蔽处理，去除对机载pd雷达回波无贡献的地表散射点：

4)通过dlcd读取各地表散射点的地表地貌分类码，确定其地貌类型，并分类匹配相应的后向散射系数模型：

5)根据步骤4)确定的后向散射系数模型以及步骤2)获得的各地表散射点的擦地角计算出机载pd雷达作用范围内各地表散射点的后向散射系数，并将步骤2)获得的各地表散射点相对于载机的斜距、方位角和俯仰角以及后向散射系数代入ward模型，由此实现高保真的非均匀地杂波仿真。

在步骤1)中，所述的dem数据采用先进星载热发射和反射辐射仪全球数字高程模型数据；dlcd采用美国地质勘探局的国家土地覆盖数据。

在步骤2)中，所述的利用dem数据提取高程值，并计算每个地表散射点相对于载机的斜距、方位角、俯仰角、擦地角在内的信息的方法是：通过设定载机位置坐标，并根据机载pd雷达作用范围内各地表散射点与载机的位置关系计算出斜距、方位角、俯仰角和擦地角在内的信息。

在步骤3)中，所述的对地表散射点进行遮蔽处理的方法是：将机载pd雷达工作在下视模式时，在其探测范围内由于地形起伏而无法被照射到的地表散射点的擦地角置零，由此剔除掉存在遮蔽的地表散射点。

在步骤4)中，所述的通过dlcd读取各地表散射点的地表地貌分类码，确定其地貌类型，并分类匹配相应的后向散射系数模型的方法是：针对dlcd中多样的地貌划分类型，采用修正的morchin模型、morchin模型和γ-f模型分别匹配路面自然地貌、水面地貌和城市地貌。

在步骤5)中，所述的根据步骤4)确定的后向散射系数模型以及步骤2)获得的各地表散射点的擦地角计算出机载pd雷达作用范围内各地表散射点的后向散射系数，并将步骤2)获得的各地表散射点相对于载机的斜距、方位角和俯仰角以及上述后向散射系数代入ward模型，由此实现高保真的非均匀地杂波仿真的方法是：根据步骤4)确定的各地表散射点的地貌类型，将步骤2)中获得的各地表散射点的擦地角分别带入到步骤4)中匹配的各模型公式中，进而计算得到各地表散射点的后向散射系数；然后将各地表散射点的斜距、方位角、俯仰角和上述后向散射系数带入ward模型中，得到高保真的非均匀地杂波仿真。

本发明提供的基于dem和dlcd的机载pd雷达高保真非均匀地杂波仿真方法是在dem数据的基础上引入了dlcd，利用dlcd中的地貌类型信息匹配不同的后向散射模型计算后向散射系数，并结合ward地杂波模型实现了机载pd雷达高保真非均匀地杂波仿真。本发明所提方法由于融合了dlcd，使地貌类型划分更加丰富，同时基于地貌类型分类匹配不同的后向散射系数模型，大大提高了仿真地杂波的准确性和真实性。

附图说明

图1为本发明提供的基于dem和dlcd的机载pd雷达高保真非均匀地杂波仿真方法的实现框图；

图2为机载pd雷达地杂波模型；

图3为基于dem和dlcd的地理信息处理流程图；

图4为截取的dem和dlcd数据的原始图；

图5为载机位置设定及机载pd雷达作用范围内高程值显示图；

图6为地表散射点相关信息；

图7为擦地角计算及遮蔽处理；

图8为后向散射系数对比图；

图9为谷歌地图匹配结果；

图10为杂波功率谱对比图；

图11为距离多普勒谱对比图；

图12为地杂波空时二维谱。

具体实施方式

下面结合附图和具体实施例对本发明提供的基于dem和dlcd的机载pd雷达高保真非均匀地杂波仿真方法进行详细说明。

如图1所示，本发明提供的基于dem和dlcd的机载pd雷达高保真非均匀地杂波仿真方法包括按顺序进行的下列步骤：

1)截取仿真所用的dem数据和dlcd；

本发明选用的dem数据是先进星载热发射和反射辐射仪全球数字高程模型(advancedspacebornethermalemissionandreflectionradiometerglobaldigitalelevationmodel,astergdem)数据，该模型数据分辨率为30米。dlcd选取美国地质勘探局(unitedstatesgeologicalsurvey，usgs)的国家土地覆盖数据(nationallandcoverdata，nlcd)，该数据将地貌类型划分为9大类，21小类，数据分辨率与gdem数据一致，可以很好地实现对同一地区的数据匹配。仿真中计算机载pd雷达作用范围，通过globalmapper软件直接截取所需经纬度的dem数据和dlcd。

2)利用dem数据提取高程值，并计算每个地表散射点相对于载机的斜距、方位角、俯仰角、擦地角在内的信息；

dem数据是实际地形空间分布特征的离散化模型，主要描述的是实际地形的高程信息。dem数据采用等间隔的正方形网格来表示地表散射点的高程值(即海拔)，并且以矩阵的形式进行数据存储，通过载机位置坐标即可从dem数据中提取出各地表散射点的高程值及其相对于载机的方位角、俯仰角、斜距、擦地角在内的信息。

如图2所示，设载机位置坐标为(xp,yp,h)，dem数据中任意地表散射点m相对于载机的斜距、方位角、俯仰角分别用下式进行计算：

式中，(xm,ym,hm)为地表散射点m的坐标，hm为地表散射点m的高程值，可以通过截取的dem数据直接获取。

擦地角是机载pd雷达与地表散射点的连线和地表散射点所在平面之间的夹角，在计算后向散射系数时必须要确定擦地角的大小。与均匀杂波仿真中擦地角为地表散射点的俯仰角不同，当机载pd雷达工作在复杂多变的起伏地形中时，擦地角并不等于该地表散射点的俯仰角。此时，计算擦地角首先需要计算待测地表散射点到载机的视线向量，然后计算地表散射点所在平面的法向量，并根据两向量之间的关系来确定地表散射点的擦地角。

如图2所示，dem数据中三个相邻的地表散射点a(xa,ya,ha)，b(xb,yb,hb)，c(xc,yc,hc)组成了一个真实地形曲面(当数据分辨率较高时，网格间距较小，可认为该曲面为平面)，对于待测地表散射点a，u为地表散射点a到载机的视线向量，w为地表散射点a所在三角平面的法向量，则地表散射点a的擦地角为地表散射点a到载机的视线向量u和地表散射点a所在三角平面的法向量w的夹角的余角。

根据上述载机位置坐标可得地表散射点a到载机的视线向量为：

u＝(xp-xa,yp-ya,h-ha)(2)

同时由图2可得，向量和分别为：

假设地表散射点a所在三角平面的法向量w＝(x,y,z)，由于w为地表散射点a所在平面的法向量，因此地表散射点a所在三角平面的法向量w与向量和垂直，可得：

为了保证地表散射点a所在三角平面的法向量w指向z轴正方向，取z＝1。解上述方程组可得：

w＝((ha-hb)/ds,(ha-hc)/ds,1)(5)

式中，ds为数据分辨率(即网格间距)。

此时可得地表散射点a的擦地角为

φa＝π/2-arccos(w·u/(||w||||u||))(6)

以此类推可以得到机载pd雷达作用范围内所有地表散射点的擦地角。

3)根据步骤2)获得的擦地角和斜距信息对地表散射点进行遮蔽处理，去除对机载pd雷达回波无贡献的地表散射点：

当机载pd雷达工作在下视模式时，在其探测范围内会存在由于地形起伏而无法被照射到的地表散射点，这些地表散射点对地杂波仿真无贡献，因此在地杂波仿真过程中，必须考虑遮蔽处理。遮蔽主要包括视线遮蔽和自身遮蔽，如图2所示，以机载pd雷达同一照射方向上的两个地表散射点d和e为例，对遮蔽处理进行详细说明。假设地表散射点e为待测地表散射点，根据dem数据计算得到两地表散射点相对于载机的斜距分别为rd和re，擦地角分别为φd和φe，则：

a)若rd＜re，并且φd＜φe，则低海拔的待测地表散射点e被高海拔的地表散射点d遮蔽，此时待测地表散射点e无法被机载pd雷达照射到(即地表散射点e对机载pd雷达不可见)，这种情况即为视线遮蔽；

b)若地表散射点e所在平面的法向量和视线向量之间的夹角大于90°(即擦地角φe小于0°)，这种情况即为自身遮蔽。

当待测地表散射点e存在上述视线遮蔽或者自身遮蔽时，将擦地角φe置为0°即可完成对待测地表散射点e的遮蔽处理。

4)通过dlcd读取各地表散射点的地表地貌分类码，确定其地貌类型，并分类匹配相应的后向散射系数模型：

不同的后向散射系数模型匹配的地貌类型不同。修正的morchin模型可以很好地匹配陆面自然地貌；morchin模型仅与水面地貌相匹配；γ-f模型对城市地貌进行了建模，可以很好地用于城市地杂波仿真。因此，本发明针对dlcd中多样的地貌划分类型，分别采用修正的morchin模型、morchin模型和γ-f模型对不同地貌类型的地杂波进行仿真。下面就水面、陆面和城市模型匹配及其后向散射系数计算分别进行论述。

a)水面地貌是机载pd雷达实际工作中较为常见的地貌类型，也是dlcd中明确标注的地貌类型。由c.m.william提出的morchin模型可以很好地匹配水面地貌，进而得到较准确的水面后向散射系数，此时第m(m＝1,2,…,m)个地表散射点的水面后向散射系数可用下式进行计算：

其中，m为仿真区域内地表散射点个数，φm为擦地角；ss为海情的级数(也即水面的波动起伏情况)；β＝[2.44(ss+1)^1.08]/57.29，单位为rad；σc为初始后向散射系数，可以表示为：

式中，φc＝arcsin(λ/4πhe)，单位为rad；he≈0.025+0.046ss^1.72，单位为米，表示水面的粗糙程度；建议k＝1.9。由上述内容可以看出，擦地角和海情的级数是计算水面后向散射系数的关键。

b)在高分辨率dlcd中，陆面自然地貌被细致地划分为多种地貌类型，并对各地貌类型进行了详细的编码。修正的morchin模型适应的陆面地貌类型较为丰富，可以很好地拟合dlcd中给出的各陆面地貌。因此本发明结合dlcd中丰富的地貌类型选用修正的morchin模型计算各陆面自然地貌的后向散射系数，此时第m(m＝1,2,…,m)个地表散射点的陆面后向散射系数可用下式进行计算：

其中，m为仿真区域内地表散射点个数；φm为擦地角；f为雷达工作频率，单位为ghz；下标lax,(x＝1,2,3…)表示dlcd中划分的各陆面自然地貌类型，σc为初始后向散射系数，ρ、χ、β0是由地貌类型确定的常量，其取值如表1所示

表1各陆面地貌类型的参数值

c)在高分辨率dlcd中，城市等人造地貌类型被明确的地表分类码标注出来，因此必须确定合适的后向散射系数模型来拟合城市地貌，进而准确地计算出城市地貌的后向散射系数。γ-f模型定义了城市的后向散射系数，可以很好地用于城市地杂波仿真，此时第m(m＝1,2,…,m)个地表散射点的城市后向散射系数可用下式进行计算：

σci(m)＝γaf^bsin(φm+c)(10)

其中，m为仿真区域内地表散射点个数；γ、a、b、c是与地貌类型有关的参数，针对城市地貌，γ＝0.316，a＝0.36，b＝0.18，c＝0.7；f为雷达工作频率，单位为ghz。

5)根据步骤4)确定的后向散射系数模型以及步骤2)获得的各地表散射点的擦地角计算出机载pd雷达作用范围内各地表散射点的后向散射系数，并将步骤2)获得的各地表散射点相对于载机的斜距、方位角和俯仰角以及后向散射系数代入ward模型，由此实现高保真的非均匀地杂波仿真：

图3为基于dem和dlcd的地理信息处理流程图。从图中可以看出，前文中由dem数据确定的各地表散射点相对于载机的斜距、方位角、俯仰角以及由dlcd的地表分类码匹配后向散射系数模型计算得到的后向散射系数σq(m)(下标q＝wa,lax或ci代表dlcd中第m个地表散射点对应的地貌类型)是计算各地表散射点回波强度的关键。假设载机做匀速直线运动，速度大小为v，在机载平台上沿航向方向均匀放置n元线阵，阵元间隔d＝0.5λ，其中λ为机载pd雷达波长，脉冲重复频率为fr，并且假设第l个距离门包含m个地表散射点，根据ward模型，该距离门的第m个地表散射点在机载pd雷达第n个阵元第k个脉冲采样下的数据可以表示为：

式中，表示第l个距离门内第m个地表散射点的回波强度，pt为机载pd雷达发射峰值功率，gt为阵列天线发射增益，gt为阵列天线接收增益，ls为系统损耗，为天线方向图，可以表示为：

式中，θm为第l个距离门内第m个地表散射点的方位角，为第l个距离门内第m个地表散射点的俯仰角，in为子阵加权，θ0和分别表示天线主瓣方位角和俯仰角。式(11)中，ωs,lm和ωt,lm分别表示第l个距离门内第m个地表散射点相对于载机的空间角频率和时间角频率，其表达式为：

将该距离门上所有地表散射点的回波采样数据进行叠加，可得到机载pd雷达第n个阵元第k个脉冲接收的第l个距离门的杂波采样数据：

将接收到的杂波采样数据cl按照阵元-脉冲的规则排列成一列向量，即可得到第l个距离门地杂波的空时快拍数据cl：

由此可以得到第m个地表散射点的空时二维杂波的功率谱为：

其中，rrc＝cu·cu^h表示杂波协方差矩阵，cu＝[c1c2…cl]nk×l。

表示第m个地表散射点的空时二维导向矢量，at(v)和分别表示时间和空间导向矢量，可以表示为：

本发明提供的基于dem和dlcd的机载pd雷达高保真非均匀地杂波仿真方法的效果可以通过以下仿真结果进一步说明。

仿真条件描述：机载pd雷达的天线阵列为阵元数n＝8的均匀线阵，阵元间距d＝λ/2，主瓣波束水平方位角为90°，俯仰角为15°，机载pd雷达工作波长为0.05m，脉冲重复频率为7000hz，最小可分辨距离150m，相干处理脉冲数k＝64，杂噪比100db，载机速度为87.5m/s，飞行高度为600m。

图4为截取的dem和dlcd的原始图。其中图4(a)为dem数据原始图，图4(b)为dlcd原始图。选定的经纬度范围为北纬36°45’～37°15’，西经112°50’～113°20’。从dlcd图中可以看出，所选区域包含城市、高山、水面以及庄稼地等地形，具有很好的普适性。

图5为载机位置设定及机载pd雷达作用范围内高程值显示图，从图中可以看出，选取的仿真区域，地形起伏剧烈，非均匀特性明显。

图6为机载pd雷达作用范围内各地面散射点相对于载机的斜距、方位角、俯仰角、多普勒信息图，从图中可以清楚地看出地杂波的非均匀特性非常明显，尤其是俯仰角对地形的变化尤为剧烈。

图7为擦地角计算及遮蔽处理结果图。其中，图7(a)为机载pd雷达作用范围内的地表散射点的擦地角，图7(b)为进行遮蔽处理以后的擦地角。两图对比可以看出，所选地形中存在大量自身遮蔽和视线遮蔽的地表散射点，尤其是在山体附近，存在明显的视线遮蔽，这些遮蔽点对于仿真地杂波是没有任何贡献的，因此被置零处理了。

图8为不同方法得到的后向散射系数图，其中图8(a)为基于传统统计模型得到的高斯分布后向散射系数图，图8(b)为仅利用dem数据，通过计算各地形因子(包括坡度、地面曲率、地形起伏、地表粗糙度)，判断得到的后向散射系数图。从图8(a)和图8(b)中可以看出，图8(b)中引入dem数据，并通过计算地形因子得到的后向散射系数有明显的非均匀性。图8(c)为基于本发明方法计算得到的各地表散射点的后向散射系数图。它是根据nlcd的地貌类型划分结果匹配后向散射系数模型并根据dem数据计算的擦地角得到的，其中划圈处为美国hurricane市。

图9是在谷歌地图上截取的相应区域的卫星实测图。经过对比，可以清晰地看出，本发明方法计算的后向散射系数能够更好地匹配实际地形。

图10为杂波功率谱对比图，其中，图10(a)和图10(b)分别为仅利用dem数据得到的杂波功率谱以及本发明方法得到的杂波功率谱，两图进行对比可以看出，本发明方法由于引入了dlcd，可以很好地拟合城市(划圈内为美国hurricane市的回波功率谱)等地貌特征，使得复杂的城市杂波功率谱也能很好地体现。

图11为基于高斯统计模型均匀地杂波仿真方法与本发明方法得到的距离多普勒谱对比图。其中图11(a)为均匀地杂波的距离多普勒谱，图11(b)为本发明方法的距离多普勒谱，经过对比可以得到，通过引入dem和dlcd仿真实现的地杂波有明显的非均匀性和非平稳性。

图12为仿真的机载pd雷达地杂波的空时二维谱，从图中可以看出，由于载机运动引起的频率扩展，机载pd雷达侧视阵地杂波的空时二维谱沿对角线呈直线分布，这是机载pd雷达侧视阵地杂波的主要特征。