CN101667855B - 多普勒频率的估测系统及方法 - Google Patents

Abstract

本发明公开了一种多普勒频率的估测系统及方法，该估测系统包含一基础基底投影装置，接收多个通道取样信号，并将多个通道取样信号投影至一组正交基底，以产生多个通道相关向量、一多项式产生装置连接至基础基底投影装置，依据多个通道相关向量、一估测频道包络至噪声干扰功率比率、及一频道包络功率，以产生一个目标多项式，以及一极值决定装置连接至多项式产生装置，以决定目标多项式的极值，并输出该极值对应的一频率作为一估测多普勒频率。应用本发明，能估测一移动装置的多普勒频率，而无需大量地计算。

Description

多普勒频率的估测系统及方法

技术领域

本发明涉及无线传输技术领域，特别涉及一种用于多普勒频率的估测系统及方法。

背景技术

由于传送端和接收端的相对运动，使得无线传输通道会动态地改变，无线传输通道的时变特性的统计量与多普勒频谱有着密切的关系。该多普勒频谱通过对一通道脉冲回应的自相关函数(autocorrelation function)执行傅立叶变换而产生，亦称为多普勒扩展(Doppler spread)，其中，多普勒扩展与传送端和接收端的相对运动速度成正比。

由于多普勒效应会使通道的特性随着时间而改变，进而增加信号品质的不确定性。同时，多普勒扩展会在接收端产生频率偏移效应，增加接收端的误比特率(Bit Error Rate，BER)。故估测多普勒扩展及接收端的运动速度会直接影响到移动式无线通讯系统的性能。例如：在IS-05、WCDMA、及CDMA2000等码分多址系统(CDMA)中，接收端的运动速度作为切换手机制式参考的重要参数。在OFDM系统中，准确的估测多普勒扩展及接收端的运动速度有助于接收端的同步估测及时变通道的估测。

美国专利USP6,563,861号使用傅立叶变换以直接估测一多普勒频谱的最大频宽。图1为现有技术使用傅立叶变换以直接估测一多普勒频谱的最大频宽的结构示意图。如图1所示，乘法器26接收采样输入信号r_n，并乘上传输符元

以产生对应的自相关函数，再通过一低通滤波器28去除高频噪声。第二处理模组30执行快速傅立叶变换，以将自相关函数转换至频率域，从而产生一信号频谱。第三处理模组32具有多个频率域的相关器(Correlator)，其中，每一个相关器系对该信号频谱的频带(Band)与多普勒频谱的频带执行相关运算(Correlation)，以产生一相似比率计量(likelihood ratio metric)。最大值模组34则选择最大的相似比率计量对应的频率，作为一估测多普勒频率。

上述技术虽可直接估测出一估测多普勒频率，然而，由于接收端运动速度的限制，一般多普勒频率约为几十赫兹至1.5K赫兹(Hz)，为能准确估测多普勒频率，第二处理模组30执行快速傅立叶变换时，其傅立叶变换的解析度需求较大，亦即在频率域需要处理的资料量变得很多，进而增加整个系统成本。然而，若将第二处理模组30的傅立叶变换的解析度降低，虽可减少资料量，但是却无法准确地估测接收端的多普勒频率，进而影响系统的性能。

发明内容

有鉴于此，本发明的主要目的在于提供一种多普勒频率的估测系统，降低处理的资料量、节省系统成本，同时能准确且快速地估测多普勒频率。

本发明的另一目的在于提供一种多普勒频率的估测方法，降低处理的资料量、节省系统成本，同时能准确且快速地估测多普勒频率。

为实现本发明目的，本发明提出一种多普勒频率的估测系统，包含一基础基底投影装置、一多项式产生装置、及一极值决定装置。基础基底投影装置接收多个通道取样信号，并将多个通道取样信号投影至一组正交基底，以产生多个通道相关向量。多项式产生装置连接至基础基底投影装置，依据多个通道相关向量、一估测频道包络至噪声干扰功率比率、及一频道包络功率，以产生一个目标多项式。极值决定装置连接至多项式产生装置，以决定目标多项式的极值，并输出该极值对应的一频率作为一估测多普勒频率。

为实现本发明的另一目的，本发明提出一种多普勒频率的估测方法，包含：一基础基底投影步骤，用以接收多个通道取样信号，并将多个通道取样信号投影至一组正交基底，以产生多个通道相关向量；一多项式产生步骤，依据多个通道相关向量、一估测频道包络至噪声干扰功率比率、及一频道包络功率，以产生一个目标多项式；一极值决定步骤，决定目标多项式的极值，并输出该极值对应的一频率作为一估测多普勒频率。

由上述技术方案可见，本发明的多普勒频率的估测系统及方法，可以降低处理的资料量、节省系统成本，同时能准确且快速地估测多普勒频率，能估测一移动装置的多普勒频率，而无需大量地计算。

附图说明

图1为现有技术使用傅立叶变换以直接估测一多普勒频谱的最大频宽的结构示意图。

图2为本发明多普勒频率的估测系统的结构示意图。

图3为本发明基础基底投影装置的结构示意图。

图4为本发明区域极值测试装置工作的流程示意图。

图5为本发明全面极值决定装置工作的流程示意图。

图6至图11为本发明f_p ^(cand)分别为1Hz、41Hz、81Hz、121Hz、161Hz、201Hz时w1、w2、w3的值的示意图。

图12至图17为本发明f_p ^(cand)分别为1Hz、41Hz、81Hz、121Hz、161Hz、201Hz时z1、z2、z3、z4的值的示意图。

图18为本发明估测出

值的电路结构示意图。

图19为本发明估测出L_p(υ，γ，Ω)值的电路结构示意图。

图20为一多普勒频率估测的比较示意图。

图21为一正规划均方根误差与估测频道包络至噪声干扰功率比率γ的比较示意图。

附图中的标号说明

乘法器26                       低通滤波器28

第二处理模组30                 第三处理模组32

最大值模组34

基础基底投影装置210            多项式产生装置220

极值决定装置230          区域极值测试装置231

全面极值决定装置233      关连器310

步骤S405～S450           步骤S510～S530

具体实施方式

为使本发明的目的、技术方案及优点更加清楚明白，以下参照附图并举具体实施例，对本发明作进一步详细说明。

请参见图2，图2为本发明多普勒频率的估测系统的结构示意图，应用于一无线传输系统的接收端。该估测系统200包含一基础基底投影装置210、一多项式产生装置220及一极值决定装置230。

基础基底投影装置210用以接收N个通道取样信号h＝[h₁，h₂，...，h_N，]^T，并将该N个通道取样信号h＝[h₁，h₂，...，h_N，]^T投影至一组正交基底V_K，进而产生多个通道相关向量υ，亦即 $\mathrm{\υ}=V_K^{T}h,$ 当中，K为近似的阶数。

多项式产生装置220连接至基础基底投影装置210，依据多个通道相关向量υ、一估测频道包络至噪声干扰功率比率γ、及一频道包络功率Ω(channel-envelop power)，以产生一个目标多项式。

极值决定装置230连接至多项式产生装置220，以决定该目标多项式的极值，并输出该极值对应的频率作为估测多普勒频率f_m ^*。

估测多普勒频率f_m ^*依据最大相似估测(Maximum Likelihood Estimation)法获得。亦即，将该估测多普勒频率f_m ^*代入一费用函数(Cost Function)时，可获得该费用函数的极值。该费用函数为瑞利衰落通道(Rayleigh FadingChannel)相关函数，且依据最大相似估测法，该费用函数可为一相似函数Λ(f_m)或一对数相似函数L(f_m)，其中L(f_m)＝-ln(Λ(f_m))。

当该费用函数为相似函数Λ(f_m)时，对应的目标多项式在估测多普勒频率f_m ^*时，会获得一极大值。当该费用函数为对数相似函数L(f_m)时，对应的目标多项式在估测多普勒频率f_m ^*时，会获得一极小值。本实施例中，该费用函数为对数相似函数L(f_m)。依据最大相似估测法及瑞利衰落通道，该对数相似函数L(f_m)可用公式(1)表示：

$L\left(f_m\right)=\ln \left(\det \left(C\right)\right)+\frac{1}{\mathrm{\Ω}}{h}^T{C}^{-1}h,---\left(1\right)$

当中，Ω为频道包络功率，h为N个通道取样信号，h^T为h的转置矩阵(Transpose Matrix)，C为拓普利兹(Toeplitz)对称矩阵，其对应于N个通道取样信号h的共变异数矩阵(Covariance Matrix)。该矩阵C可用公式(2)表示：

γ为估测频道包络至噪声干扰功率比率，I_N为NxN的单位矩阵，J₀(.)为第一类零阶贝索函数(Zero-order Bessel Function of the First Kind)，f_m为最大的多普勒扩展(Doppler spread)，T_s为一采样间隔，估测频道包络至噪声干扰功率比率γ可表示为

Ω为频道包络功率，σ_n+i ²为噪声(n)及干扰(i)的标准差。故本发明多普勒频率的估测系统所估测的多普勒频率f_m ^*的方法，可转化为找寻一多普勒频率f_m ^*，使得对数相似函数L(f_m)为最小，亦即用数学式可表示为：

${\hat{f}}_m^{\left(\mathrm{ML}\right)}=\arg \underset{f_m}{\min }L\left(f_m\right)---\left(3\right)$

公式(2)中的第一类零阶贝索函数J₀(.)可表示成一具有偶次方的多项式，亦即， $J_0\left(x\right)=\underset{k=0}{\overset{\∞}{\mathrm{\Σ}}}{g}_k{x}^{2k}.$ 故公式(2)中的第一类零阶贝索函数J₀(nπT_sf_m)可近似为公式(4)所示：

$J_0\left(\mathrm{n\π}{T}_s{f}_m\right)=\underset{k=\mathrm{0,1},...}{\overset{K}{\mathrm{\Σ}}}\left(g_k{\left(\mathrm{n\π}{T}_s\right)}^{2k}\right)f_m^{2k},---\left(4\right)$

当中{g_k}_{k＝0，1，...}为常数，k表示索引。将公式(2)中矩阵C的各项以公式(4)表示，矩阵C可近似为公式(5)：

$C\≅{\tilde{C}}_K\left(f_m\right)={\mathrm{\γ}}^{-1}{I}_n+\underset{k=0}{\overset{K}{\mathrm{\Σ}}}{g}_k{\left(\mathrm{\π}{T}_s{f}_m\right)}^{2k}{A}_k\≡{\mathrm{\γ}}^{-1}{I}_n+B_K\left(f_m\right),---\left(5\right)$

当中，K为近似的阶数，B_K(f_m)为{A_k}_k＝0 ^K的权重组合。{A_k}_k＝0 ^K为Toeplitz对称矩阵，其可表示为公式(6)：

其中，0≤k≤K

由公式(5)可知，该B_K(f_m)为对称矩阵。利用一正规正交矩阵V_K及一对称矩阵M_K(f_m)，矩阵B_K(f_m)可近似为公式(7)：

$B_K\left(f_m\right)\≅V_K{M}_K\left(f_m\right)V_K^T,---\left(7\right)$

V_K其大小为N×ρ_K，亦即V_K由ρ_K个正规正交列向量构成，每一正规正交列向量为N维。V_K的每一元素与f_m无关。

由公式(7)，对称矩阵M_K(f_m)可由以下公式得到：

$M_K\left(f_m\right)=V_K^T{B}_K\left(f_m\right)V_K,---\left(8\right)$

且对称矩阵M_K(f_m)为一ρ_K×ρ_K全秩(full rank)对称矩阵，且其元素为f_m的多项式。公式(7)中使用近似符号但是当V_K可展开矩阵B_K(f_m)的向量空间时，公式(7)中矩阵B_K(f_m)则等于V_KM_K(f_m)V_K ^T。亦即，ρ_K代表矩阵B_K(f_m)的全秩(full rank)或减秩(reduced rank)，因此1≤ρ_K≤N。

经过适当挑选的V_K，使其所展开的空间等于矩阵B_K(f_m)的ρ_K个较大特征值(eigenvalue)所对应的特征向量所展开的空间。也就是说，通过适当选择V_K，可以使对称矩阵M_K(f_m)的ρ_K个特征值等于矩阵B_K(f_m)的ρ_K个较大特征值。

将

取代C，公式(1)中的对数相似函数L(f_m)可改写成：

$L\left(f_m\right)=\ln \left(\det \left({\tilde{C}}_K\left(f_m\right)\right)\right)+\frac{1}{\mathrm{\Ω}}{h}^T{\tilde{C}}_K{\left(f_m\right)}^{-1}h\≡L_i\left(f_m\right)+L_d\left(f_m\right)---\left(9\right)$

该对数相似函数L(f_m)可分成两部分，一部分是与输入无关的矩阵行列式 $L_i\left(f_m\right)=\ln \left(\det \left({\tilde{C}}_K\left(f_m\right)\right)\right),$ 另一部分是与输入相关的二次式乘积项 $L_d\left(f_m\right)={\mathrm{\Ω}}^{-1}{h}^T{\tilde{C}}_K{\left(f_m\right)}^{-1}h,$ 其中，矩阵行列式

可通过离线计算而获得，二次式乘积项L_d(f_m)则需线上依据N个通道取样信号h＝[h₁，h₂，...，h_N]^T计算而产生。

依据逆矩阵定理，二次式乘积项L_d(f_m)中的

的逆矩阵

可由V_K及对称矩阵M_K(f_m)近似为：

${\tilde{C}}_K{\left(f_m\right)}^{-1}\≅\mathrm{\γ}(I_N-V_K{V}_K^T)+V_K{({\mathrm{\γ}}^{-1}{I}_{{\mathrm{\ρ}}_K}+M_K\left(f_m\right))}^{-1}{V}_K^T,---\left(10\right)$

矩阵行列式

可由V_K及对称矩阵M_K(f_m)近似为：

$\det \left({\tilde{C}}_K\left(f_m\right)\right)\≅{\mathrm{\γ}}^{-(N-{\mathrm{\ρ}}_K)}\det ({\mathrm{\γ}}^{-1}{I}_{{\mathrm{\ρ}}_K}+M_K\left(f_m\right)),---\left(11\right)$

对于特定的N及K，公式(10)中的逆矩阵

的元素为f_m的多项式，可使用高斯消去法予以计算。公式(11)中的行列式

亦为f_m的多项式，其值可使用数值方法获得。由前述可知，矩阵

的行列式及逆矩阵

与N个通道取样信号h＝[h₁，h₂，...，h_N]^T无关，且可事先计算。

行列式

为f_m的多项式，故可用一r_i次的f_m多项式表示，亦即行列式的自然对数L_i(f_m)可改写为：

$L_i\left(f_m\right)=\ln (c_0\left(\mathrm{\γ}\right)+c_1\left(\mathrm{\γ}\right)f_m+...+c_{r_i}\left(\mathrm{\γ}\right)f_m^{r_i})=\ln \left(c{\left(\mathrm{\γ}\right)}^T{f}^{r_i}\right),---\left(12\right)$

当中，

为多项式的系数，且为估测频道包络至噪声干扰功率比率γ的函数， $c\left(\mathrm{\γ}\right)={[c_0\left(\mathrm{\γ}\right),c_1\left(\mathrm{\γ}\right),...,c_{r_i}\left(\mathrm{\γ}\right)]}^T$ 及 $f^{r_i}={[1,f_m,...,f_m^{r_i}]}^T$ 为相对应的向量表示形式。

由于V_K ^Th可以使用υ表示，故二次式乘积项L_d(f_m)可改写为：

$L_d\left(f_m\right)=\frac{1}{\mathrm{\Ω}}{\mathrm{\υ}}^T{({\mathrm{\γ}}^{-1}{I}_{{\mathrm{\ρ}}_K}+M_K\left(f_m\right))}^{-1}\mathrm{\υ}+\mathrm{\η},---\left(13\right)$

公式(13)中右边第一项为

对应于该多个通道相关向量υ的权重和，公式(13)中右边第二项η为Ω^-1γh^T(I_N-V_KV_K ^T)h，η为一与f_m无关的偏移(offset)，η不影响

的计算值，为简化起见，可将η忽略。因此L_d(f_m)可视为f_m的一分式多项式，且可改写为：

$L_d\left(f_m\right)=\frac{1}{\mathrm{\Ω}}\frac{b_0(\mathrm{\υ},\mathrm{\γ})+b_1(\mathrm{\υ},\mathrm{\γ})f_m+...+b_{r_{d2}}(\mathrm{\υ},\mathrm{\γ})f_m^{r_{d2}}}{a_0\left(\mathrm{\γ}\right)+a_1\left(\mathrm{\γ}\right)f_m+...+a_{r_{d1}}\left(\mathrm{\γ}\right)f_m^{r_{d1}}}=\frac{1}{\mathrm{\Ω}}\frac{b{(\mathrm{\υ},\mathrm{\γ})}^T{f}^{r_{d2}}}{a{\left(\mathrm{\γ}\right)}^T{f}^{r_{d1}}}---\left(14\right)$

当中，b(υ，γ)及a(γ)分别为分子多项式及分母多项式的系数向量。亦即公式(9)中的对数相似函数L(f_m)可改写成：

$L\left(f_m\right)=\{\ln \left(c{\left(\mathrm{\γ}\right)}^T{f}^{r_i}\right)+\frac{1}{\mathrm{\Ω}}\frac{b{(\mathrm{\υ},\mathrm{\γ})}^T{f}^{r_{d2}}}{a{\left(\mathrm{\γ}\right)}^T{f}^{r_{d1}}}\},---\left(15\right)$

非线性的ML最佳化可转换为：

${\hat{f}}_m^{\left(\mathrm{ML}\right)}=\arg \underset{f_m}{\min }\{\ln \left(c{\left(\mathrm{\γ}\right)}^T{f}^{r_i}\right)+\frac{1}{\mathrm{\Ω}}\frac{b{(\mathrm{\υ},\mathrm{\γ})}^T{f}^{r_{d2}}}{a{\left(\mathrm{\γ}\right)}^T{f}^{r_{d1}}}\}.---\left(16\right)$

为获得分子多项式的系数向量b(υ，γ)，需先决定一组向量 $\mathrm{\υ}={[{\mathrm{\υ}}_1,{\mathrm{\υ}}_2,...,{\mathrm{\υ}}_{{\mathrm{\ρ}}_K}]}^T,$ 该向量υ为N个通道取样信号h＝[h₁，h₂，...，h_N]^T投影至一正规正交基底V_K而获得。

基础基底投影装置210接收N个通道取样信号h＝[h₁，h₂，...，h_N]^T，并将该N个通道取样信号h＝[h₁，h₂，...，h_N]^T投影至一正交基底V_K，以产生多个通道相关向量υ。图3本发明基础基底投影装置210的结构示意图。如图3所示，基础基底投影装置210包含ρ_K个关连器(correlator)310，第i个(1≤i≤ρ_K)关连器的系数为[v_1，i，v_2，i，...，v_N，i]^T。亦即，{v_i，j}为正交基底V_K的元素，故正交基底V_K可表示为：

$V_K=\left[\begin{array}{cccc}{v}_{\mathrm{1,1}}& v_{\mathrm{1,2}}& \·\·\·& v_{1,{\mathrm{\ρ}}_k}\\ v_{\mathrm{2,1}}& v_{\mathrm{2,2}}& \·\·\·& v_{2,{\mathrm{\ρ}}_k}\\ \·& \·& & \·\\ \·& \·& \·\·\·& \·\\ \·& \·& & \·\\ v_{N,1}& v_{N,2}& \·\·\·& v_{N,{\mathrm{\ρ}}_k}\end{array}\right].---\left(17\right)$

多项式产生装置220连接至基础基底投影装置210，依据多个通道相关向量υ、一估测频道包络至噪声干扰功率比率γ、以及一频道包络功率Ω，以产生一个目标多项式。多项式产生装置220产生目标多项式的系数，以供极值决定装置230决定目标多项式的极值。目标多项式相关于对数相似函数L(f_m)。亦即，将非线性的ML估测转换为多项式寻根程式(polynomial root-searchprocedure)。

本发明主要关键在于：(1)基础基底投影装置210中的正规正交矩阵V_K的设计；(2)极值决定装置230中找寻对数相似函数L(f_m)极值的演算法，或找寻对数相似函数L(f_m)极小值的演算法。

在图3中的基础基底投影装置210选择一组与Toeplitz对称矩阵A_k′相关的正规正交基底V_K，以便适当地分解B_K(f_m)，其中

通过对该Toeplitz对称矩阵A_k′的特征值分解，可据以选择对称矩阵A_k′的ρ_K个特征向量。该ρ_K个特征向量对应对称矩阵A_k′的ρ_K个绝对值较大值的特征值，并构成V_K。具体来说，V_K的产生方式，是挑选某一个矩阵A_k，其中，k＝0～K，由该矩阵对应前ρ_K个较大特征值的特征向量构成V_K，实际应用中，由于矩阵A_k(k＝0～K)的秩是2k+1，且其秩随k值的增大而增加，只要满足选择的k符合2k+1≥ρ_K的条件的矩阵(A₀～A_K)即可，较佳地，选择符合条件的最小值，该最小值

另一种取得V_K的方法如下。定义矩阵如下

$T_{{\mathrm{\ρ}}_K}=\left[\begin{array}{cccc}{0}^0& 0^1& \·\·\·& 0^{{\mathrm{\ρ}}_K}\\ 1^0& 1^1& \·\·\·& 1^{{\mathrm{\ρ}}_K}\\ \·& \·& & \·\\ \·& \·& \·\·\·& \·\\ \·& \·& & \·\\ {(N-1)}^0& {(N-1)}^1& \·\·\·& {(N-1)}^{{\mathrm{\ρ}}_K}\end{array}\right].---\left(18\right)$

依序对的每一列做Gram-Schmidt正交化，所得ρ_K个N×1向量即可构成另一组V_K。

通过上述两种方式确定的V_K能够对B_K(f_m)进行适当的分解，且V_K保留与B_k(f_m)相同的空间，使对称矩阵M_K(f_m)的ρ_K个特征值等于矩阵B_K(f_m)的ρ_K个较大特征值，从而满足前述对于V_K所提的要求。

极值决定装置230包含一区域极值测试装置231及一全面极值决定装置233。区域极值测试装置231连接至多项式产生装置220，将最大的多普勒扩展(估测的频率范围)f_m区分为多个子频率频段，并计算每一个子频率频段中是否具有目标多项式的极值，当目标多项式在该子频率频段中有极值存在，则输出该子频率频段中的序号(index)。

全面极值决定装置233连接至区域极值测试装置231，使用插值方法以获得一子频率频段区域最小值所对应的区域最小值频率，再由该区域最小值频率使用插值方法获得一个子频率频段中目标多项式的极值，再从多个目标多项式的极值中选择最小值，用以作为目标多项式的最小值，并选择对应的频率作为估测多普勒频率。

在区域极值测试装置231中，将估测的频率范围f_m分成P个子频率频段，其中区分可为均等区分，亦可为非均等区分。均等区分表示每个子频率频段的大小均相同，非均等区分表示每个子频率频段的大小无需均相同。这些子频率频段的边界以{f₁，f₂，...，f_P+1}表示，其中，f₁及f_P+1分别为频率范围f_m的左右边界。选择这些子频率频段大小时，该子频率频段较佳仅包含一个区域极小值。

区域极值测试装置231用于在每一子频率频段中找出区域极小值的位置。在每一子频率频段中，区域极小值会发生在对数相似函数L(f_m)的微分由负转正的位置，亦即当存在区域极小值时，

等于0。由于选取子频率频段大小时，该子频率频段仅包含一个区域极小值，故对p＝1，2，...，P+1而言，可在f_p处定义一微分对数相似值 ${\stackrel{\·}{L}}_p=\frac{\∂L\left(f_m\right)}{\∂f_m}{|}_{f_m=f_p}.$ 因此当左边边界

小于0、且右边边界大于0时，第p个子频率频段为具有区域极小值的子频率频段。同理，当

大于0时，该频率范围f_m的左边边界具有区域极小值，当

小于0时，该频率范围f_m的右边边界具有区域极小值。为清楚起见，将频率范围f_m的左边边界、右边边界分别标示为第0个子频率频段及第P+1个子频率频段。需注意的是，

的值依据υ、γ、Ω而决定。

图4为本发明区域极值测试装置231工作的流程示意图。在步骤S405中，读取

的值，并设定一指标p为1。在步骤S410中，判断该指标p是否大于或等于1、且小于或等于一预设值P。若是，执行步骤S415，若否，执行步骤S430。

在步骤S415中，判断

是否小于或等于0、且

是否大于或等于0。若是，表示第p个子频率频段中有区域极小值，故在步骤S420中储存第p个子频率频段的序号p。若否，表示第p个子频率频段中没有区域极小值，故执行步骤S425，以判断下一个子频率频段中是否有区域极小值。

在步骤S430中，判断是否大于等于0。若是，表示该频率范围f_m的左边边界具有区域极小值，故在步骤S435中储存第0个子频率频段的序号0。若否，表示频率范围f_m的左边边界没有区域极小值，并执行步骤S440，以判断

是否小于0。若是，表示频率范围f_m的右边边界有区域极小值，故在步骤S445中储存第P+1个子频率频段的序号P+1。若否，表示频率范围f_m的右边边界没有区域极小值，执行步骤S450，输出具有极小值子频率频段的序号。

全面极值决定装置233连接至区域极值测试装置231，并接收区域极值测试装置231输出的子频率频段的序号，全面极值决定装置233依据子频率频段的序号，使用插值方法以获得一子频率频段区域最小值所对应的区域最小值频率，再通过区域最小值频率使用插值方法获得一个子频率频段中目标多项式的极值，再从多个目标多项式的极值中选择最小值，以作为目标多项式的最小值，并选择对应的频率作为估测多普勒频率。

区域极值测试装置231输出p的序号，以代表该p子频率频段具有一区域极小值。全面极值决定装置233则依据公式(19)求出在该p子频率频段中的区域最小值对应的区域最小值频率f_p ^(cand)：

$f_p^{\left(\mathrm{cand}\right)}=\frac{{\stackrel{\·}{L}}_{P+1}{f}_p-{\stackrel{\·}{L}}_P{f}_{p+1}}{{\stackrel{\·}{L}}_{P+1}-{\stackrel{\·}{L}}_P},---\left(19\right)$

其中，公式(19)中的p取值范围为区域极值测试装置输出的子频率频段序号。对于p＝0及p＝P+1，分别对应频率范围f_m的左边边界f₁、及右边边界f_P+1。

全面极值决定装置233依据公式(19)求出区域最小值频率f_p ^(cand)，再利用式(20)计算该区域最小值频率f_p ^(cand)对应的目标多项式的极值L(f_p ^(cand))：

$L\left(f_p^{\left(\mathrm{cand}\right)}\right)=L_p+\frac{{\stackrel{\·}{L}}_P}{2}(f_p^{\left(\mathrm{cand}\right)}-f_p),---\left(20\right)$

其中， $L_p=L\left(f_m\right){|}_{f_m=f_p}$ 为目标多项式在f_m＝f_p的值，p为区域极值测试装置输出的子频率频段序号。如果在一子频率频段中有一个以上的目标多项式的极值时，这些极值也加入全面极值决定装置233的选择中。全面极值决定装置233依据公式(20)计算各个区域最小值频率对应的目标多项式的极值，并从多个目标多项式的极值中选择最小值，以作为目标多项式的最小值，并选择对应的频率作为估测多普勒频率

${\hat{f}}_m^{\left(\mathrm{ML}\right)}\≅\arg \min \left\{L\left(f_p^{\left(\mathrm{cand}\right)}\right)\right\}.---\left(21\right)$

图5为本发明全面极值决定装置233工作的流程示意图。在步骤S510中，全面极值决定装置233使用插值方法求出在p子频率频段中的区域最小值对应的区域最小值频率f_p ^(cand)，其中，p为区域极值测试装置输出的子频率频段序号。

在步骤S520中，全面极值决定装置233使用插值方法求出区域最小值频率f_p ^(cand)对应的目标多项式的极值L(f_p ^(cand))。

在步骤S530中，全面极值决定装置233接收步骤S520中所产生的至少一个目标多项式的极值L(f_p ^(cand))，并从至少一个目标多项式的极值中选择最小值，以作为目标多项式的最小值，并选择其相对应的频率作为估测多普勒频率

由前述说明可知，全面极值决定装置233需利用{L_p}_p＝1 ^P+1及的值。目标多项式的微分可表示成公式(22)：

$\stackrel{\·}{L}\left(f_m\right)=\frac{{\stackrel{\·}{b}}_i^T\left(\mathrm{\γ}\right)f^{r_i-1}}{{\stackrel{\·}{a}}_i^T\left(\mathrm{\γ}\right)f^{r_i}}+\frac{1}{\mathrm{\Ω}}\frac{{\stackrel{\·}{b}}_d^T(\mathrm{\υ},\mathrm{\γ})f^{r_{d1}+r_{d2}-1}}{{\stackrel{\·}{a}}_d^T\left(\mathrm{\γ}\right)f^{2r_{d1}}}.---\left(22\right)$

故

的值为：

$\stackrel{\·}{L}\left(f_p\right)=\frac{{\stackrel{\·}{b}}_i^T\left(\mathrm{\γ}\right)f^{r_i-1}}{{\stackrel{\·}{a}}_i^T\left(\mathrm{\γ}\right)f^{r_i}}+\frac{1}{\mathrm{\Ω}}\frac{{\stackrel{\·}{b}}_d^T(\mathrm{\υ},\mathrm{\γ})f^{r_{d1}+r_{d2}-1}}{{\stackrel{\·}{a}}_d^T\left(\mathrm{\γ}\right)f^{2r_{d1}}}{|}_{f_m=f_p}.---\left(23\right)$

亦即，对一给定的f_p，该 ${\stackrel{\·}{b}}_d^T(\mathrm{\υ},\mathrm{\γ})f^{r_{d1}+r_{d2}-1}{|}_{f_m=f_p}$ 具有下列的形式：

$\underset{s=0}{\overset{r_{d1}+r_{d2}-1}{\mathrm{\Σ}}}\underset{l=1}{\overset{{\mathrm{\ρ}}_K}{\mathrm{\Σ}}}\underset{j=1}{\overset{{\mathrm{\ρ}}_K}{\mathrm{\Σ}}}c(j,l,s){\mathrm{\υ}}_j{\mathrm{\υ}}_l{\mathrm{\γ}}^s,---\left(24\right)$

而 ${\stackrel{\·}{a}}_i^T\left(\mathrm{\γ}\right)f^{r_i}{|}_{f_m=f_p},{\stackrel{\·}{b}}_i^T\left(\mathrm{\γ}\right)f^{r_i-1}{|}_{f_m=f_p},{\stackrel{\·}{a}}_d^T\left(\mathrm{\γ}\right)f^{2r_{d1}{r}_i}{|}_{f_m=f_p}$ 则为下列的形式：

$\underset{s}{\mathrm{\Σ}}c\left(s\right){\mathrm{\γ}}^s,---\left(25\right)$

其中，c(s)与c(j，l，s)为纯量系数。可见，对一给定的f_p，

可以由三个以γ为变量的一元函数以及一个以υ、γ为变量的二元函数构成。同理，对一给定的f_p，L(f_p)可以由两个以γ为变量的一元函数以及一个以υ、γ为变量的二元函数构成。其中，各个函数中变量的系数可以根据前述公式(11)和(13)确定。

为方便对L(f_p)和

的计算，基于上述结论，对于一给定的f_p，可将公式(15)中的 $\{\ln \left(c{\left(\mathrm{\γ}\right)}^T{f}^{r_i}\right)+\frac{1}{\mathrm{\Ω}}\frac{b{(\mathrm{\υ},\mathrm{\γ})}^T{f}^{r_{d2}}}{a{\left(\mathrm{\γ}\right)}^T{f}^{r_{d1}}}\}$ 改写为 $L\left(f_p\right)=w1\left(\mathrm{\γ}\right)+\frac{w2(\mathrm{\γ},\mathrm{\υ})}{w3\left(\mathrm{\γ}\right)}\frac{1}{\mathrm{\Ω}},$ 并将 $\{\ln \left(c{\left(\mathrm{\γ}\right)}^T{f}^{r_i}\right)+\frac{1}{\mathrm{\Ω}}\frac{b{(\mathrm{\υ},\mathrm{\γ})}^T{f}^{r_{d2}}}{a{\left(\mathrm{\γ}\right)}^T{f}^{r_{d1}}}\}$ 的微分改写为 $\stackrel{\·}{L}\left(f_p\right)=\frac{z1\left(\mathrm{\γ}\right)}{z2\left(\mathrm{\γ}\right)}+\frac{z3(\mathrm{\γ},\mathrm{\υ})}{d_h\left(\mathrm{\γ}\right)}\frac{1}{\mathrm{\Ω}}.$ 同时，在确定f_p的取值后，根据公式(11)与公式(13)并利用数值方法得到w1(γ)、w2(γ，υ)和w3(γ)中的各个系数(以下简称为w1、w2、w3的值)以及z1(γ)、z2(γ)、z3(γ，υ)和d_h(γ)中的各个系数(以下简称为z1、z2、z3、d_h的值)，并通过表格方式记录，从而可以简化L(f_p)和

的计算。

当N为40、K为2、ρ_K为5时，图6至图11为本发明显示当f_p ^(cand)(即图6至图11表格中的频率fp，也就是预先设定的各个子频率频段的边界值，p＝1，2，...，P)分别为1Hz、41Hz、81Hz、121Hz、161Hz、201Hz时，w1、w2、w3的值的示意图。其中，w1、w3两列的系数表示的是γ的相应幂次下对应的系数，w3列的系数表示的是γ相应幂次与υ相应幂次的乘积对应的系数。依据图6至图11所示的w1、w2、w3，即可算出L(f_m)在频率f_p ^(cand)分别为1Hz、41Hz、81Hz、121Hz、161Hz、201Hz时的值。

图12至图17所示的分别为第一至第六表格。图12至图17为本发明显示当f_p ^(cand)分别为1Hz、41Hz、81Hz、121Hz、161Hz、201Hz时，z1、z2、z3、d_h的值。依据图12至图17所示的z1、z2、z3、d_h，即可算出

在频率f_p ^(cand)分别为1Hz、41Hz、81Hz、121Hz、161Hz、201Hz时的值。

在图12至图17中，每一行表示对应至{γ^s}_s＝0 ^Sp的系数，S_p为

中的最大幂次(order)。如图12至图17所示，其中包含一些非常小的系数，可省略并不影响精确度，故S_p在不同图中并不相同。同时，z3栏包含子栏，其与输入参数的关系以下举例说明。在K为2、ρ_K为5时，z3栏包含9个子栏(υ₁ ²，υ₂ ²，υ₃ ²，υ₄ ²，υ₅ ²，υ₁υ₂，υ₃υ₄，υ₃υ₅，υ₄υ₅)。例如：第3行第1子栏则代表为γ²υ₁ ²的系数。通过图12至图17的表格及υ，γ，Ω即可估测

在各个子频率频段的边界值。图18为本发明估测出

值的电路结构示意图。其中，t_i，j ^(p)为第p个表格中第i行第j列的系数。通过查表及图18所示的电路即可估测出

的值。

同理，图6至图11分别列出公式(26)中w1、w2、w3的值，分为第一至第六表格。通过图6至图11的表格及υ，γ，Ω即可估测L_p(υ，γ，Ω)在各个子频率频段的边界值。图19为本发明估测出L_p(υ，γ，Ω)值的电路结构示意图。其中，u_i，j ^(p)为第p个表格中第i行第j列的系数。通过查表及图19所示的电路即可估测出L_p(υ，γ，Ω)值。

图20为一多普勒频率估测的比较示意图。与一正确的多普勒频率做比较。其中横轴为多普勒测试频率的正确值，纵轴为一正规划均方根误差，其定义如下：

${\mathrm{RMSE}}_{\mathrm{Norm}}=\frac{{\left(E\right[{(f_m-{\hat{f}}_m)}^2\left]\right)}^{\frac{1}{2}}}{f_m}.---\left(28\right)$

其中，该比较图在WCDMA通讯系统下进行类比。相关参数中，符号时间(Symbol duration)设定为66.67μs，通道为瑞利衰落通道且估测频道包络至噪声干扰功率比率γ为7dB，观测向量的长度N为50，同时K为3。在图20中，本发明的多普勒频率估测系统所产生的多普勒频率估测值使用实心线绘制，而由理论推导得出的无偏估测子的理论下界-正规划均方根克拉美-罗下界(Cramer-Rao lower bound，CRLB)则使用虚线表示。由图20所示，本发明的多普勒频率估测准确度，可在多普勒频率测试范围内逼近理论的克拉美-罗下界。

图21为一正规划均方根误差与估测频道包络至噪声干扰功率比率γ的比较示意图。在一极低的多普勒频率f_m为20Hz的环境下进行类比，类比相关参数同图20的设定。由图21所示，本发明的多普勒频率的正规划均方根误差，在估测频道包络至噪声干扰功率比率测试范围内，均与理论的克拉美-罗下界差距有限，可提供准确的多普勒频率估测值。

在本实施例中，图6至图11的表格及图12至图17的表格由多项式产生装置220来实现，亦即多项式产生装置220产生目标多项式(15)与(22)，并根据上述表格中的系数和p，υ，γ，Ω输出对数相似函数L(f_m)及其微分的值。

由上述说明可知，本发明可解决现有技术中，为应付资料量增多而增加整个系统成本的问题，或是现有技术中在低成本架构下无法准确地估测接收端的多普勒频率的问题。相较于直接使用最大相似估测法，先将N个通道取样讯号投影成为ρ_K个值(ρ_K≤N)再进行多普勒频率估测，避免直接使用N个原始讯号做估测运算，省去大量地计算。本发明能降低处理的资料量，以节省系统成本，同时本发明通过查表方法而能准确且快速地估测多普勒频率。

上述实施例仅仅是为了方便说明而举例而已，本发明所主张的权利范围应以发明申请范围所述为准，而非仅限于上述实施例。

Claims (14)

1.一种多普勒频率的估测系统，其特征在于，该估测系统包含：

一基础基底投影装置，用以接收多个通道取样信号，并将所述多个通道取样信号投影至一组正交基底，从而产生多个通道相关向量；

一多项式产生装置，连接至所述基础基底投影装置，依据所述多个通道相关向量、一估测频道包络与噪声干扰功率比率、及一频道包络功率，以产生一个目标多项式；以及

一极值决定装置，连接至所述多项式产生装置，用以决定所述目标多项式的极值及其相对应的频率，并输出一估测多普勒频率，所述极值决定装置包含：

一区域极值测试装置，连接至所述多项式产生装置，用以将多普勒频率区分为多个子频率频段，并分别依据所述多个子频率频段计算所述目标多项式是否有极值，若是，则输出该子频率频段相对应的序号；以及

一全面极值决定装置，连接至所述区域极值测试装置，用以获得前述子频率频段中区域最小值及其相对应的区域最小值频率，再根据所述区域最小值频率获得该子频率频段中所述目标多项式的极值，进而从所述目标多项式的极值中选择最小值，作为所述目标多项式的最小值，并选择其相对应的频率作为所述估测多普勒频率。

2.如权利要求1所述的估测系统，其特征在于，所述目标多项式对应于一费用函数，所述费用函数相关于瑞利衰落通道，且所述费用函数为多普勒频率的函数。

3.如权利要求2所述的估测系统，其特征在于，当所述费用函数为一相似函数时，对应的所述目标多项式在所述估测多普勒频率时获得一极大值。

4.如权利要求2所述的估测系统，其特征在于，当所述费用函数为一对数相似函数时，对应的所述目标多项式在所述估测多普勒频率时获得一极小值。

5.如权利要求1所述的估测系统，其特征在于，所述极值决定装置利用查表来计算每一个所述子频率频段中所述目标多项式的极值。

6.如权利要求5所述的估测系统，其特征在于，所述区域极值测试装置执行下列步骤：

(A)读取

的值，并设定一指标即子频率频段p为1，当中，f_m为估测的频率范围，P为子频率频段数目，L(f_m)为对数相似函数；

(B)判断该指标即子频率频段p是否大于等于1、且小于等于一预设值P，若是，执行步骤(C)，若否，执行步骤(F)；

(C)判断

是否小于等于0、且

是否大于等于0，若是，表示第p个子频率频段中有区域极小值，执行步骤(D)，若否，表示第p个子频率频段中没有区域极小值，执行步骤(E)；

(D)储存第p个子频率频段的序号；

(E)判断下一个子频率频段中是否有区域极小值，即p＝p+1，然后返回执行步骤B；

(F)判断

是否大于等于0，若是，表示该频率范围f_m的左边边界具有区域极小值，执行步骤(G)，若否，表示该频率范围f_m的左边边界没有区域极小值，并执行步骤(H)；

(G)储存第0个子频率频段的序号0；

(H)判断

是否小于0，若是，表示该频率范围f_m的右边边界有区域极小值，执行步骤(I)，若否，表示该频率范围f_m的右边边界没有区域极小值，执行步骤(J)；

(I)储存第p+1个子频率频段的序号；以及

(J)输出具有极小值子频率频段的序号。

7.如权利要求6所述的估测系统，其特征在于，所述全面极值决定装置执行下列步骤：

(K)使用插值方法求出所述子频率频段中的区域最小值对应的区域最小值频率；

(L)使用插值方法求出所述区域最小值频率对应的所述目标多项式的极值；以及

(M)接收步骤(L)中所产生的至少一个目标多项式的极值，并从所述至少一个目标多项式的极值中选择最小值，以作为所述目标多项式的最小值，并选择其相对应的频率作为所述估测多普勒频率。

8.一种多普勒频率的估测方法，适用于无线通讯系统，其特征在于，该方法包含：

一基础基底投影步骤，用以接收多个通道取样信号，并将所述多个通道取样信号投影至一组正交基底，从而产生多个通道相关向量；

一多项式产生步骤，依据所述多个通道相关向量、一估测频道包络与噪声干扰功率比率、及一频道包络功率，以产生一个目标多项式；以及

一极值决定步骤，用以决定所述目标多项式的极值及其相对应的频率，进而输出一估测多普勒频率，包含：

一区域极值测试步骤，用以将多普勒频率区分为多个子频率频段，并分别依据所述多个子频率频段计算所述目标多项式是否有极值，若是，则输出该子频率频段相对应的序号；以及

一全面极值决定步骤，用以获得前述子频率频段中区域最小值及其相对应的区域最小值频率，再根据所述区域最小值频率获得所述子频率频段中所述目标多项式的极值，进而从所述目标多项式的极值中选择最小值，作为所述目标多项式的最小值，并选择其相对应的频率作为所述估测多普勒频率。

9.如权利要求8所述的估测方法，其特征在于，所述目标多项式对应于一费用函数，所述费用函数相关于瑞利衰落通道，且所述费用函数为多普勒频率的函数。

10.如权利要求9所述的估测方法，其特征在于，当所述费用函数为一相似函数时，对应的所述目标多项式在所述估测多普勒频率时获得一极大值。

11.如权利要求9所述的估测方法，其特征在于，当所述费用函数为一对数相似函数时，对应的所述目标多项式在所述估测多普勒频率时获得一极小值。

12.如权利要求8所述的估测方法，其特征在于，所述极值决定步骤利用查表来计算每一个所述子频率频段中所述目标多项式的极值。

13.如权利要求12所述的估测方法，其特征在于，所述区域极值测试步骤进一步包含：

(A)读取

的值，并设定一指标即子频率频段p为1，当中，

f_m为估测的频率范围，P为子频率频段数目，L(f_m)为对数相似函数；

(B)判断所述指标即子频率频段p是否大于等于1、且小于等于一预设值P，若是，执行步骤(C)，若否，执行步骤(F)；

(C)判断

是否小于等于0、且

是否大于等于0，若是，表示第p个子频率频段中有区域极小值，执行步骤(D)，若否，表示第p个子频率频段中没有区域极小值，故执行步骤(E)；

(D)储存第p个子频率频段的序号；

(E)判断下一个子频率频段中是否有区域极小值，即p＝p+1，然后返回执行步骤B；

(F)判断

是否大于等于0，若是，表示所述频率范围f_m的左边边界具有区域极小值，执行步骤(G)，若否，表示所述频率范围f_m的左边边界没有区域极小值，并执行步骤(H)；

(G)储存第0个子频率频段的序号0；

(H)判断是否小于0，若是，表示所述频率范围f_m的右边边界有区域极小值，执行步骤(I)，若否，表示所述频率范围f_m的右边边界没有区域极小值，执行步骤(J)；

(I)储存第p+1个子频率频段的序号；以及

(J)输出具有极小值子频率频段的序号。

14.如权利要求13所述的估测方法，其特征在于，所述全面极值决定步骤进一步包含：

(K)使用插值方法求出子频率频段中的区域最小值对应的区域最小值频率；

(L)使用插值方法求出所述区域最小值频率对应的该目标多项式的极值；以及

(M)接收步骤(L)中所产生至少一个所述目标多项式的极值，并从所述至少一个该目标多项式的极值中选择最小值，以作为所述目标多项式的最小值，并选择其相对应的频率作为所述估测多普勒频率。

Images