本发明提出了一种基于降阶核函数的高阶多项式相位信号参数估计方法，用于解决现有高阶多项式相位信号的参数估计方法存在的估计精度低的技术问题，实现过程为：对混合高斯白噪声的高阶多项式相位信号均匀采样，得到相位信号序列；构造相位信号序列的降阶核函数；对相位信号序列进行样条插值，得到非一致间隔序列；计算降阶信号序列；对降阶信号序列进行快速傅立叶变换；计算相位信号序列中待估计参数的估计量并输出；对相位信号序列进行解调，得到降阶解调序列；更新相位信号序列，循环构造降阶核函数，依次计算更新后的相位信号序列待估计参数的估计量并输出；更新降阶信号序列，得到一阶相位信号序列，并计算一阶参数估计量并输出。
1.一种基于降阶核函数的高阶多项式相位信号参数估计方法，包括如下步骤：(1)对混合高斯白噪声的高阶多项式相位信号x(t)进行均匀采样，得到高阶多项式相位信号序列x(n)：其中，n为相位信号序列x(n)的离散时间变量，w(n)为高斯白噪声序列，Δ为相位信号x(t)的采样间隔，A0为相位信号x(t)的幅度，{a1,...,ar,...,aP}为阶数r对应的所需估计的相位参数，P为相位信号序列x(n)的最高阶数，且P≥3；
(2)利用高阶多项式相位信号序列x(n)的最高阶数P，构造高阶多项式相位信号序列x(n)的降阶核函数ker[x(n)]，实现步骤为：
2a)定义构造高阶多项式相位信号序列x(n)的降阶核函数ker[x(n)]所需的参数，并对这些参数进行初始化：定义待估计相位参数aP对应的阶数变量P′、奇降阶迭代次数i、偶降阶迭代次数v和标记向量l(i+v)，对这些参数进行初始化，得到初始阶数变量P′＝P，初始奇降阶迭代次数i＝0，初始偶降阶迭代次数v＝0，初始标记向量l(0)＝[]；
2b)对阶数变量P′进行判定；若P′≠1且为偶数，执行步骤2c)，若P′≠1且为奇数，执行步骤2d)，若P′＝1，输出标记向量l(M+H)，并执行步骤2e)，其中l(M+H)＝[l1,...,lk,...,lM+H]，编号k＝1,...,M+H，M为标记向量l(M+H)中0元素个数，H为标记向量l(M+H)中1元素个数；
2c)对阶数变量P′进行降阶，得到P′＝P′/2，令偶降阶迭代次数v＝v+1，在标记向量l(i+v-1)末端增添元素0v，得到标记向量l(i+v)＝[l(i+v-1),0v]，并执行步骤2b)，其中0v为标记向量l(i+v)中偶降阶迭代次数为v时得到的0元素；
2d)对阶数变量P′进行降阶，得到P′＝P′-1，令奇降阶迭代次数i＝i+1，奇降阶阶数Pi＝P′，在标记向量l(i+v-1)末端增添元素1i，得到标记向量l(i+v)＝[l(i+v-1),1i]，并执行步骤2b)，其中，1i为标记向量l(i+v)中奇降阶迭代次数为i时得到的1元素，Pi为奇降阶迭代次数为i时的阶数；
2e)根据标记向量l 中元素的取值顺序，构造高阶多项式相位信号序列x(n)的降阶核函数ker[x(n)]；
(3)根据降阶核函数ker[x(n)]的需要，对相位信号序列x(n)进行样条插值，得到相位信号序列x(n)的非一致间隔序列(4)计算降阶信号序列xp：将非一致间隔序列 代入降阶核函数ker[x(n)]，得到降阶信号序列xp：
(5)对降阶信号序列xp进行快速傅立叶变换，得到降阶信号序列xp的频域函数f(Ω)；
(6)计算相位信号序列x(n)中待估计参数的估计量 并输出：采用一维搜索方法，计算频域函数f(Ω)的幅度函数|f(Ω)|的最大值对应的频率点 并根据频率点计算相位信号序列x(n)中待估计参数的估计量 并输出；
(7)对相位信号序列x(n)进行解调，得到降阶解调序列x′(n)：(8)更新相位信号序列，令x(n)＝x′(n)，得到降阶后的相位信号序列x(n)，循环构造降阶核函数，依次计算更新后的相位信号序列待估计参数的估计量 并输出：令最高阶数减1，执行步骤(2)～(7)，直到P＝1；
(9)令降阶信号序列xp＝x(n)，执行步骤(4)～(6)，并计算一阶参数估计量 并输出。
2.根据权利要求1所述的基于降阶核函数的高阶多项式相位信号参数估计方法，步骤
2e)中所述降阶核函数ker[x(n)]，其表达式为：
其中，a0,a0+1,...b0为相位信号序列x(n)中离散时间变量n的取值范围， 为序列算子，包含 和 两种运算，lk＝0时，对信号进行 运算，并更新序列的离散时间变量m的取值范围，lk＝1时，对信号进行 运算，并更新序列 的离散时间变量n的取值范围：
n＝ak,ak+1,...,bk ak＝ak-1+τi,bk＝bk-1-τi，j＝1,2,...,M，i＝1,2,...,H，k＝1,...,M+H，c1＝c2＝…＝cM-1＝1，
其中，ak,ak+1,...,bk为 运算后序列的离散时间变量的取值范围， 为向上取整，为向下取整， 为计算过程引入的非一致间隔序列。
3.根据权利要求2所述的基于降阶核函数的高阶多项式相位信号参数估计方法，步骤(6)中所述的计算相位信号序列x(n)中待估计参数的估计量 其计算公式为：其中，NP为降阶信号序列xP的长度， 为幅度函数|f(Ω)|的最大值对应的频率点，Ω为频率变量。
4.根据权利要求1所述的基于降阶核函数的高阶多项式相位信号参数估计方法，步骤(3)中所述的对相位信号序列x(n)进行样条插值，其插值公式为：其中，x0(n′)为所需的插值点，n′为所需离散时间，ceil(n′)为取大于n′且紧邻的离散时间，floor(n′)为取小于n′且紧邻的离散时间。
5.根据权利要求1所述的基于降阶核函数的高阶多项式相位信号参数估计方法，步骤(5)中所述的降阶信号序列xp的频域函数f(Ω)，其表达式为：其中，NP为降阶信号序列xP的长度，m为离散时间变量，Ω为频率变量。
6.根据权利要求1所述的基于降阶核函数的高阶多项式相位信号参数估计方法，步骤(9)中所述的计算一阶参数估计量 其计算公式为：其中，NP为降阶信号序列xP的长度， 为幅度函数|f(Ω)|的最大值对应的频率点，Ω为频率变量。
基于降阶核函数的高阶多项式相位信号参数估计方法
技术领域
[0001] 本发明属于非平稳信号分析技术领域，涉及一种非平稳信号的参数化估计方法，具体涉及一种基于降阶核函数的高阶多项式相位信号参数估计方法，可用于机动目标雷达回波信号的检测和估计以及ISAR成像技术等领域。
背景技术
[0002] 多项式相位信号(PPS)己被广泛用于非平稳信号如雷达(脉冲多普勒雷达、SAR和ISAR)、声纳、通信、生物医药、地震分析和动物声学应用建模中。线性、高斯和平稳是传统信号处理的主要方面。而现代信号处理的特点是非线性、非高斯和非平稳信号。非平稳信号的研究中，多项式相位信号是比较常见的一种。噪声背景下(特别是色噪声、杂波等复杂环境)多项式相位信号的时频分析和参数估计问题在雷达、通信、声纳、生物医学、自然界、地震信号分析等领域函待解决。例如，在移动通信领域，发射机与接收机存在相对运动，会造成发射信号的相位随着时间变化。瞬时相位会随着收发间距离的连续变化而变化，此时接收信号就是多项式相位信号的近似表示。通信中时变衰落信道估计采用的信号是调频信号(二阶PPS即是线性调频信号，高阶PPS即是非线性调频信号)。
[0003] 作为一类重要的信号模型，多项式相位信号可对许多实际信号进行建模，多项式相位信号中的相位系数在不同的应用中都具有重要的意义，能提供丰富的信息，因此，多项式相位信号的参数估计具有广阔的应用前景和重要的研究意义。通常情况下，一个白噪声中的离散P阶多项式相位信号可以表示为：
[0004] x(n)＝s(n)+w(n),-N/2≤n≤N/2，
[0005]
[0006] 其中w(n)为高斯白噪声，Δ为信号采样间隔，信号的采样点数为N+1，b0为信号幅度，为一个常数，{a1,a2,...,aP}为信号相位的1到P阶参数。需要用混合噪声的多项式相位信号的N+1个采样点值估计参数{a1,a2,...,aP}。多项式相位信号的参数估计有多种方法，最大似然估计可以得到最小的均方误差(MSE)，然而由于其采用了多维搜索，这种方法的计算量巨大。为了减少计算量，高阶模糊函数(HAF)和立方相位函数(CPF)的方法以MSE为代价，采用了一维搜索，大大减少了计算复杂度。HAF反复利用相位差分(PD)算子，对P阶多项式相位信号采用P-1次PD算子降阶，将信号变为单频信号，采用快速傅立叶变化方法估计信号最高阶参数，然后降阶解调相位信号，重复操作循环计算低阶参数估计量；然而HAF方法采用P-1次PD算子降阶，非线性度较高，使得估计的门限很高，算法在信噪比(SNR)低于门限的时候会彻底失效，在针对于高阶多项式相位信号时，由于误差增殖效应，HAF方法的高门限会严重影响估计性能，进而导致估计的精度低。而立方相位函数和混合高阶模糊函数(CPF-HAF)方法对HAF方法进行了改进，先通过相位差分将高阶多项式相位信号降阶到3阶，利用CPF方法估计3阶多项式相位信号的相位参数，通过设置不同的时间点，进行谱峰搜索，建立方程组，计算参数估计量。虽然CPF-HAF方法在低SNR时的性能比HAF要好，但非线性程度仍旧较高，噪声分量较多，对参数估计影响较大，导致估计方法的精度较低。
发明内容
[0007] 本发明的目的在于克服上述已有技术存在的缺陷，提出了一种基于降阶核函数的高阶多项式相位信号参数估计方法，用于解决现有高阶多项式相位信号的参数估计方法存在的估计精度低的技术问题。
[0008] 为实现上述目的，本发明采取的技术方案包括如下步骤：
[0009] (1)对混合高斯白噪声的高阶多项式相位信号x(t)进行均匀采样，得到高阶多项式相位信号序列x(n)：
[0010]
[0011] 其中，n为相位信号序列x(n)的离散时间变量，w(n)为高斯白噪声序列，Δ为相位信号x(t)的采样间隔，A0为相位信号x(t)的幅度，{a1,...,ar,...,aP}为阶数r对应的所需估计的相位参数，P为相位信号序列x(n)的最高阶数；
[0012] (2)利用高阶多项式相位信号序列x(n)的最高阶数P，构造高阶多项式相位信号序列x(n)的降阶核函数ker[x(n)]，实现步骤为：
[0013] 2a)定义构造高阶多项式相位信号序列x(n)的降阶核函数ker[x(n)]所需的参数，并对这些参数进行初始化：定义待估计相位参数aP对应的阶数变量P′、奇降阶迭代次数i、偶降阶迭代次数v和标记向量l(i+v)，对这些参数进行初始化，得到初始阶数变量P′＝P，初始(0)奇降阶迭代次数i＝0，初始偶降阶迭代次数v＝0，初始标记向量l ＝[]；
[0014] 2b)对阶数变量P′进行判定；若P′≠1且为偶数，执行步骤2c)，若P′≠1且为奇数，执行步骤2d)，若P′＝1，输出标记向量l(M+H)，并执行步骤2e)，其中l(M+H)＝[l1,...,lk,...,lM+H]，编号k＝1,...,M+H，M为标记向量l(M+H)中0元素个数，H为标记向量l(M+H)中1元素个数；
[0015] 2c)对阶数变量P′进行降阶，得到P′＝P′/2，令偶降阶迭代次数v＝v+1，在标记向量l(i+v-1)末端增添元素0v，得到标记向量l(i+v)＝[l(i+v-1),0v]，并执行步骤2b)，其中0v为标记向量l(i+v)中偶降阶迭代次数为v时得到的0元素；
[0016] 2d)对阶数变量P′进行降阶，得到P′＝P′-1，令奇降阶迭代次数i＝i+1，奇降阶阶数Pi＝P′，在标记向量l(i+v-1)末端增添元素1i，得到标记向量l(i+v)＝[l(i+v-1),1i]，并执行步骤2b)，其中，1i为标记向量l(i+v)中奇降阶迭代次数为i时得到的1元素，Pi为奇降阶迭代次数为i时的阶数；
[0017] 2e)根据标记向量l(M+H)中元素的取值顺序，构造高阶多项式相位信号序列x(n)的降阶核函数ker[x(n)]；
[0018] (3)根据降阶核函数ker[x(n)]的需要，对相位信号序列x(n)进行样条插值，得到相位信号序列x(n)的非一致间隔序列
[0019] (4)计算降阶信号序列xp：将非一致间隔序列 代入降阶核函数ker[x(n)]，得到降阶信号序列xp：
[0020]
[0021] (5)对降阶信号序列xp进行快速傅立叶变换，得到降阶信号序列xp的频域函数f(Ω)；
[0022] (6)计算相位信号序列x(n)中待估计参数的估计量 并输出：采用一维搜索方法，计算频域函数f(Ω)的幅度函数|f(Ω)|的最大值对应的频率点 并根据频率点 计算相位信号序列x(n)中待估计参数的估计量 并输出；
[0023] (7)对相位信号序列x(n)进行解调，得到降阶解调序列x′(n)：
[0024]
[0025] (8)更新相位信号序列，令x(n)＝x′(n)，得到降阶后的相位信号序列x(n)，循环构造降阶核函数，依次计算更新后的相位信号序列待估计参数的估计量 并输出：令最高阶数减1，执行步骤(2)～(6)，直到P＝1；
[0026] (9)令降阶信号序列xp＝x(n)，执行步骤(4)～(6)，并计算一阶参数估计量 并输出。
[0027] 本发明与现有技术相比，具有如下优点：
[0028] 1)本发明在对相位参数进行估计时，采用了降阶核函数，对高阶多项式相位信号进行奇偶降阶，降低了核函数的非线性程度，消除了现有技术估计门限高，且均方误差较大的缺陷，有效地提高了参数估计的精度。
[0029] 2)本发明在计算相位参数估计量时，采用快速傅立叶变换，以及一维搜索方法，降低了参数估计的复杂度。
附图说明
[0030] 图1是本发明的实现流程图；
[0031] 图2是采用本发明、HAF方法和CPF-HAF方法，分别对8阶相位多项式信号的8阶、7阶、6阶和5阶参数a8,a7,a6,a5估计的最小均方误差对比的仿真结果图。
具体实施方式
[0032] 以下结合附图和具体实施例，对本发明作进一步详细说明。
[0033] 参照图1，基于降阶核函数的高阶多项式相位信号参数估计方法，包括如下步骤：
[0034] 步骤1，对混合高斯白噪声的高阶多项式相位信号x(t)进行均匀采样，得到高阶多项式相位信号序列x(n)：
[0035]
[0036] 其中，n为相位信号序列x(n)的离散时间变量，w(n)为高斯白噪声序列，Δ为相位信号x(t)的采样间隔，A0为相位信号x(t)的幅度，a1,...,ar,...,aP为阶数r对应的待估计相位参数，P＝8为相位信号序列x(n)的最高阶数；
[0037] 步骤2，利用高阶多项式相位信号序列x(n)的最高阶数P，构造高阶多项式相位信号序列x(n)的降阶核函数ker[x(n)]，实现步骤为：
[0038] 2a)定义构造高阶多项式相位信号序列x(n)的降阶核函数ker[x(n)]所需的参数，并对这些参数进行初始化：定义待估计相位参数aP对应的阶数变量P′、奇降阶迭代次数i、偶降阶迭代次数v和标记向量l(i+v)，对这些参数进行初始化，得到初始阶数变量P′＝P，初始奇降阶迭代次数i＝0，初始偶降阶迭代次数v＝0，初始标记向量l(0)＝[]；
[0039] 2b)对阶数变量P′进行判定；若P′≠1且为偶数，执行步骤2c)，若P′≠1且为奇数，(M+H) (M+H)执行步骤2d)，若P′＝1，输出标记向量l ，并执行步骤2e)，其中l ＝[l1,...,lk,...,lM+H]，编号k＝1,...,M+H，M为标记向量l(M+H)中0元素个数，H为标记向量l(M+H)中1元素个数；
[0040] 2c)对阶数变量P′进行降阶，得到P′＝P′/2，令偶降阶迭代次数v＝v+1，在标记向量l(i+v-1)末端增添元素0v，得到标记向量l(i+v)＝[l(i+v-1),0v]，并执行步骤2b)，其中0v为标(i+v)记向量l 中偶降阶迭代次数为v时得到的0元素；
[0041] 2d)对阶数变量P′进行降阶，得到P′＝P′-1，令奇降阶迭代次数i＝i+1，奇降阶阶数Pi＝P′，在标记向量l(i+v-1)末端增添元素1i，得到标记向量l(i+v)＝[l(i+v-1),1i]，并执行步骤2b)，其中，1i为标记向量l(i+v)中奇降阶迭代次数为i时得到的1元素，Pi为奇降阶迭代次数为i时的阶数；
[0042] 2e)根据标记向量l(M+H)中元素的取值顺序，构造高阶多项式相位信号序列x(n)的降阶核函数ker[x(n)]其表达式为：
[0043]
[0044] 其中，a0,a0+1,...b0为相位信号序列x(n)中离散时间变量n的取值范围， 为序列算子，包含 和 两种运算，lk＝0时，对信号进行 运算，并更新序列的离散时间变量m的取值范围，lk＝1时，对信号进行 运算，并更新序列
的离散时间变量n的取值范围：
[0045]
[0046]
[0047] n＝ak,ak+1,...,bk ak＝ak-1+τi,bk＝bk-1-τi，
[0048] j＝1,2,...,M，i＝1,2,...,H，k＝1,...,M+H，
[0049]
[0050]
[0051] c1＝c2＝…＝cM-1＝1，
[0052]
[0053] 其中，ak,ak+1,...,bk为 运算后序列的离散时间变量的取值范围， 为向上取整， 为向下取整， 为计算过程引入的非一致间隔序列。
[0054] 步骤3，根据降阶核函数ker[x(n)]的需要，对相位信号序列x(n)进行样条插值，得到相位信号序列x(n)的非一致间隔序列 最终的降阶核函数ker[x(n)]为2M+H个非一致间隔序列相乘形式，而这些序列需要插值得到，在插值计算时，插值点利用与其紧邻的信号序列x(n)采样点插值得到，插值公式为：
[0055]
[0056] 其中，x0(n′)为所需的插值点，n′为所需离散时间，ceil(n′)为取大于n′且紧邻的离散时间，floor(n′)为取小于n′且紧邻的离散时间。
[0057] 步骤4，计算降阶信号序列xp：将非一致间隔序列 代入降阶核函数ker[x(n)]，得到降阶信号序列xp：
[0058]
[0059] 步骤5，对降阶信号序列xp进行快速傅立叶变换，得到降阶信号序列xp的频域函数f(Ω)，其表达式为：
[0060]
[0061] 其中,NP为降阶信号序列xP的长度，m为离散时间变量，Ω为频率变量。
[0062] 步骤6，计算相位信号序列x(n)中待估计参数的估计量 并输出：采用一维搜索方法，计算频域函数f(Ω)的幅度函数|f(Ω)|的最大值对应的频率点 并根据频率点 计算相位信号序列x(n)中待估计参数的估计量 并输出，其计算公
式为：
[0063]
[0064] 其中，NP为降阶信号序列xP的长度， 为幅度函数|f(Ω)|的最大值对应的频率点，Ω为频率变量。
[0065] 步骤7，对相位信号序列x(n)进行解调，得到降阶解调序列x′(n)：
[0066]
[0067] 步骤8，更新相位信号序列，令x(n)＝x′(n)，得到降阶后的相位信号序列x(n)，此时，相位信号序列x(n)的最高阶数变为P-1，即待估计参数对应阶数减1；循环构造降阶核函数，并更新相位信号序列，依次计算更新后的相位信号序列的待估计参数的估计量并输出：令最高阶数减1，执行步骤(2)～(7)，直到P＝1；
[0068] 步骤9，令降阶信号序列xp＝x(n)，执行步骤(4)～(6)，并计算一阶参数估计量 并输出，其计算公式为：
[0069]
[0070] 其中，NP为降阶信号序列xP的长度， 为幅度函数|f(Ω)|的最大值对应的频率点，Ω为频率变量。
[0071] 以下结合仿真实验，对本发明技术效果进行进一步说明：
[0072] 1.仿真条件和内容：
[0073] 仿真条件：MATLAB 7.5.0,Intel(R)Pentium(R)2CPU 3.0GHz，Window 7Professional。
[0074] 仿真内容：计算机仿真本发明、HAF方法和CPF-HAF方法分别对8阶相位多项式信号的8阶、7阶、6阶和5阶参数a8,a7,a6,a5估计的最小均方误差随信噪比变化，其结果如图2所示。
[0075] 2.仿真结果分析：
[0076] 参照图2，图(a)、图(b)、图(c)和图(d)分别是本发明、HAF方法和CPF-HAF方法对8阶相位多项式信号的8阶、7阶、6阶和5阶参数a8,a7,a6,a5估计的最小均方误差对比的仿真结果图，横坐标为信噪比，纵坐标为参数估计的最小均方误差；其中，CRLB为参数估计的克拉美罗下界，即为参数估计的最小均方误差下限。
[0077] 从图2可以看出本发明与HAF方法和CPF-HAF方法相比，在较低信噪比时，本发明估计的最小均方误差已经接近克拉美罗下界，即具有更低的估计门限；同时，本发明估计的最小均方误差更接近克拉美罗下界，具有更高的估计精度。
