CN1248370A - 利用奇异值分解进行信道脉冲响应估算 - Google Patents

Abstract

本发明涉及CDMA接收机中CIR(信道脉冲响应)和SINR(信号对干扰加噪声比)的估算。具体地说,本发明涉及一种确定通信系统的信道脉冲响应(CIR)比如数字移动无线网(GSM网)的无线信道的CIR的方法。本发明尤其涉及确定这样的CIR,它基于已知训练序列的接收。本发明还可应用于干扰消除和用于CDMA接收机中。本发明基于这样的事实:通过先将网络分析所得到的矩阵变换为“对角矩阵”,再对所得到的对角矩阵求逆,可以减轻现有技术中存在的问题。进一步的发明是将逆阵中本征值的倒数乘以一个递减的数,直至这些旁瓣被去掉或很小以致没有什么影响。

Description

利用奇异值分解进行信道脉冲响应估算

本发明涉及CDMA接收机中CIR(信道脉冲响应)和SINR(信号对干扰加噪声比)的估算。具体地说，本发明涉及一种确定通信系统的信道脉冲响应(CIR)比如数字移动无线网(GSM网)的无线信道的CIR的方法。本发明尤其涉及确定这样的CIR，它基于已知训练序列的接收。本发明还可应用于干扰消除和用于CDMA接收机中。

为了确定CIR，必须知道一部分发射信号。对于GSM网，同步脉冲串(SB)是信号的有用部分。这些SB从每个基站以至少一个信道发送，并以规则码型发送。没有必要对GSM协议进行解码。无论SB中的数据还是与其共存的码型是不变的，并且对所有基站而言基本相同。

用这种SB来确定CIR的优点在于，它们可以表示一个较长的类噪声的预定发射信号。通常，在一段237μs的时间内发射64比特。因此，CIR的确定要求对接收信号中的脉冲串进行充分同步，以便能提取SB，然后对这些SB进行处理以确定CIR。

该处理采用估算技术来实现。为了估算CIR，发送一种已知的训练序列S_tx(t)，在产生接收信号S_rx(t)之前，该序列在通信信道中会受影响而恶化。估算CIR时的问题是要确定FIR滤波器(该滤波器用来近似CIR)的各权分量{α}，以便已知的S_tx(t)经该滤波器后尽可能接近于所接收信号S_rx(t)。

假定，发射信号S_tx(t)和接收信号S_rx(t)已知，则根据下述量来估算CIR：

1.已知的发射信号的抽样T_k＝S_tx(t₀+kτ)，-N_c≤k＜N(注意，N+N_c个抽样表示整个训练序列，该编号已被安排从已知的接收信号的干净(clean)抽样开始，在训练序列开始后它们有N_c个抽样)，和

2.测量到的接收信号的抽样R_k＝S_rx(t₀+kτ)，0≤k＜N(训练序列的头N_c个抽样假定受影响而恶化并被忽略)。

权分量{α_k}用简单的相关确定如下： ${\mathrm{\α}}_K=\underset{j=0}{\overset{N-1}{\mathrm{\Σ}}}{\stackrel{\‾}{T}}_{j-k}{R}_{j}0\≤k\≤N_c---\left(2.1\right)$

这种算法取决于发射信号的类噪声特性，据此，其自相关函数应具有低时间旁瓣。然而，由于训练序列的“干净”部分的相关特性和由于对后面的权只完成了局部的相关，因此这种算法的旁瓣特性受限。

US 5,473,632采用“最小二乘方方法”来估算CIR，其中复脉冲响应(CIR)通过使信道的实际输出与所估算的预测模型的输出之间的误差的平方最小化来推算。这种方法给出了US 5,473,632中第6专栏第35行的方程，这种方程可从各种课本中看到。

应考虑到US 5,473,632方法存在一个问题。US 5,473,632大体上提出了一种矩阵，但实际上，问题是如何用这种矩阵来计算CIR。

在US 5,473,632中，参照第5专栏第27至61行及其权利要求2，在矩阵求逆前先将“因数α”加到该矩阵，加进“因数α”后再对矩阵求逆。实际上，在计算CIR时，可以看到矩阵求逆会带来一些问题。

本发明的目的是要减轻现有技术中确定CIR和SINR所遇到的问题。

本发明大体上基于这样的事实：通过先将网络分析所得到的矩阵变换为“对角矩阵”，再对所得到的对角矩阵求逆，可以减轻现有技术中存在的问题。因此不必加“因数α”。

利用对角矩阵的一个优点在于，可以去除或忽略一些本征值(如一些小的有问题的本征值)，这样，在矩阵求逆时，可避免US5,473,632中放大的误差和噪声。

利用对角矩阵的另一个优点在于对角矩阵的性质，也就是说，可以通过只求对角元素的倒数来求得矩阵的逆阵。这暴露了小本征值所引起的困难——这些小本征值会导致逆阵中有很大的元素。

然而，忽略或去除某些本征值，会带来又一问题，如图1所示，即产生一些大旁瓣。

因此，本申请公开了进一步的发明，该发明减轻了产生大旁瓣的问题。

该进一步的发明对上面所述的发明作了改进，并且基于，不是弃掉小本征值(这是上面所述的方法)，而是代之以将逆阵中本征值的倒数乘以一个递减的数，直至这些旁瓣被去掉或很小以致没有什么影响。

下面，将参照附图来描述以上所公开的发明的优选实施方式，其中：

图1示出了采用US 5,473,632中公开的算法并取ε＝200、S/N＝0dB时所估算的CIR，

图2示出了采用US 5,473,632中公开的算法并取ε＝400、无噪声时所估算的CIR，

图3示出了利用SVD并计及10个最大本征值所估算的CIR，

图4示出了利用SVD并计及16个最大本征值所估算的CIR，

图5示出了利用SVD并计及25个最大本征值所估算的CIR，

图6示出了利用SVD并计及36个最大本征值所估算的CIR，

图7示出了利用SVD并利用n₁＝4、n₂＝25的斜坡函数且无附加噪声时所估算的CIR，

图8示出了利用SVD并利用n₁＝4，n₂＝25的斜坡函数且S/N-0dB时所估算的CIR，

图9示出了利用SVD并利用n₁＝4、n₂＝20的斜坡函数且无噪声时所估算的CIR，

图10示出了利用SVD并利用n₁＝4、n₂＝20的斜坡函数且S/N＝0dB时所估算的CIR，

图11示出了利用SVD并利用n₁＝1、n₂＝25斜坡函数且无噪声时所估算的CIR，

图12示出了利用SVD并利用n₁＝1、n₂＝25的斜坡函数且S/N＝0dB时所估算的CIR，

图13示出了利用SVD并利用n₁＝1、n₂＝64的斜坡函数且无噪声时所估算的CIR(注意dB标到-60dB)，

图14示出了利用SVD并利用n₁＝1、n₂＝64的斜坡函数且S/N-40dB时所估算的CIR，

图15示出了利用SVD并利用n₁＝1、n₂＝25的斜坡函数且S/N＝40dB时所估算的CIR，

图16示出了利用SVD并利用n₁＝1、n₂＝25的斜坡函数且无噪声时所估算的CIR，轨迹1在6.9μs 0dB处，轨迹2在18.5μs-20dB处，轨迹3在29.5μs-10dB处，

图17示出了利用SVD并利用n₁＝1、n₂＝25的斜坡函数且S/N-20dB时所估算的CIR，轨迹1在6.9μs-10dB处，轨迹2在18.5μs0dB处，轨迹3在29.5μs 0dB处，轨迹4在40.4μs-6dB处，和

图18示出了利用SVD并利用n₁＝1、n₂＝25的斜坡函数且无噪声时所估算的CIR，轨迹1在6.9μs 0dB处，轨迹2在12.2μs 0dB处。

估算权分量的一种方法是要确定这样一些权，它们可以预测最精确地与所测量的信号相符的接收信号。

假定发射信号S_tx(t)具有已知抽样：

T_k＝S_tx(t₀+kτ)，-M≤k＜N

(M为CIR的虚构的拓延——“脏(dirty)”抽样数)

而接收信号S_rx(t)具有所测量抽样：

R_k＝S_rx(t₀+kτ)，0≤k＜N

所要求的FIR权分量(α_k)可这样来估算：要求使测量到的接收信号与接收信号的估算(即通过近似信道的FIR的已知发射信号)之间的累积平方误差最小化，

即，使下式最小： ${\mathrm{\ϵ}}^2=\underset{K=0}{\overset{N-1}{\mathrm{\Σ}}}{|R_k-\underset{j=0}{\overset{M-1}{\mathrm{\Σ}}}{\mathrm{\α}}_j{T}_{k-j}|}^2---\left(2.2\right)$

这可以简化为下列方程组： $\underset{j=0}{\overset{M-1}{\mathrm{\Σ}}}\left({\mathrm{\α}}_i\underset{k=0}{\overset{N-1}{\mathrm{\Σ}}}{T}_{n-j}{\stackrel{\‾}{T}}_{n-i}\right)=\underset{k=0}{\overset{N-1}{\mathrm{\Σ}}}{R}_k{\stackrel{\‾}{T}}_{k-i}0\≤i\≤M--\left(2.3\right)$

或采用矩阵形式：

Ac＝TR(2.4)

其中，矩阵A是M×M相关矩阵，由下式给出： ${\mathrm{\α}}_{i,j}=\underset{n=0}{\overset{N-1}{\mathrm{\Σ}}}{T}_{n-j}{\stackrel{\‾}{T}}_{n-i}---\left(2.5\right)$

T是M×N的共轭延时信号的矩阵：

t_i，j＝T_ij 0≤i≤M  0≤j≤N(2.6)

而矢量c是所要求的CIR(α_k)

注意，简单相关可表示为：

c＝TR(2.7)

矩阵方程(2.4)有下列形式解：

c＝(A^-1T)R(2.8)

这表明，采用简单的矩阵乘法运算，便可根据接收信号的抽样来确定CIR。

采用US 5,473,632进行最小二乘方估算

方程(2.8)等同于US 5,473,632中第6专栏第35行所公开的方程。

如果采用这种方程，在直接解这些方程时会带来复杂问题，这是因为它们通常是病态的。事实上，可以看到，对于SB中的GSM训练序列而言，矩阵的本征值在超出了10¹⁰的范围内变化。应当注意，矩阵求逆时极易受到舍入误差的影响，并且结果也极易受到测量信号中噪声的影响。

可以看到，如果采用常规技术来求这种矩阵的逆阵或解这些方程，即便采用IEEE双精度运算，也得不到好的结果。

US 5,473,632公开了一种算法，这种算法可使矩阵的逆阵稳定并得到不易受噪声影响的结果。US 5,473,632用下式来估算CIR：

c＝((A+εl)^-1T)R(2.9)

其中，l是单位矩阵，而ε是一个常数(它们称为噪声项)。这一方程等同于US 5,473,632中第6专栏第44行所公开的方程。采用这种算法的结果如下：

这种算法的特性关键取决于变量ε的选择。对于具有192个对角元素的矩阵A，后面的图说明了US 5,473,632算法的特性，其中在确定含有单个轨迹的CIR时采用不同的ε值。

例如参照图1和2，可以看到，在分辨力(主瓣宽度)与噪声特性之间有一个折衷。约为-18dB的旁瓣电平看来是这种算法的特征。

取ε≡100，我们得到最类似的特性，据此与US 5,473,632进行比较。

采用SVD进行最小二乘方估算

在本发明中，我们已选定将奇异值分解(SVD)运用于病态矩阵，以实现最小二乘方估算。这样，可以看到，我们把要被求逆的矩阵变换为对角矩阵，再对新产生的对角矩阵求逆。

采用该SVD，矩阵A(正数限定)可以表示为：

A＝USV^t(2.10)

其中，U和V是正交矩阵，而S是对角矩阵。事实上，S的对角元素是A的本征值(均为正)并按递减次序排列：

s_1，1≥s_2，2≥..≥s_n，n(2.11)

如果该矩阵是奇异的，则某些对角元素为零。那么，A的逆阵为：

A^-1＝VS^-1U^t(2.12)

其中，S^-1是通过求S的(对角)元素的倒数得到。这里，病态矩阵的影响很明显，(几乎为零)的本征值的倒数对逆阵有很大的影响。它们会放大噪声。

传统的SVD技术仅求一些本征值的倒数，到某处便停止。本发明不是将1/λ₁置入S的逆阵中，而是用零来取而代之。这些技术的特性及其他情况综述如下。

方程(2.8)中矩阵中的系数纯粹是已知训练序列的函数，因此，可以预先计算出。对这一矩阵可求一次逆，此后，求解中随时都要进行矩阵乘法运算。因此，为了根据一组信号抽样得到CIR估算，要进行M^*N次复数乘法和M^*(N-1)次复数加法，总共要进行6^*M^*N+2^*M^*(N-1)步运算。

SVD在超定的最小二乘方问题中的应用，只要求涉及该矩阵的伪逆阵中N个最大的本征值的作用。这样，矩阵被模拟时，它有64维，因此，N在1至64范围内。后面的图说明了N＝10，16，25，36时SVD算法的特性。增加本征值个数的作用是为了提高估算量的分辨力，即减小主瓣宽度。

尽管没有示出，进行一些模拟可以看到，在10至20个本征值的范围内，噪声特性迅速下降。

在有噪声情况下，当使伪逆阵的秩最小化时，本SVD方法将取得良好的特性。结果中最为麻烦的问题是持久的约为-13dB的时间旁瓣。以傅里叶原理进行类推，其中截短的时间波形将给出-13dB的频率旁瓣，可以发现，不是在某一数目后进行简单的截断而是斜坡降低连续本征值的影响，将很有利于在时间旁瓣上进行控制。

有若干种斜坡函数可以采用，例如Bartlett、Blackman-Harris、Poisson、Riesz、Hamming、Kaiser-Bessel、Lanczos、Tukey和其他众所周知的窗口函数。本实施方法中所选的函数从Hanning窗导出，尽管下面所定义的升余弦斜坡被证明最为成功：

其中，W_i权对应于第i个本征值的作用。作为例子，图7示出了n₁＝4、n₂＝25的结果。

注意，-13dB的旁瓣已被降至约-22dB(大大好于US 5,473,632的特性)。主瓣宽度也比ε＝100时US 5,473,632算法得到的窄，并具有较低的旁瓣。图8中示出了S/N值为0dB的特性。

这一特性类似于ε＝100时US 5,473,632方案所得到特性，不过注意其主瓣更窄，这将降低噪声特性。

以波瓣宽度为代价可改善噪声特性，如图9和10所示，其中n₁＝4、n₂＝20。

经过大量的试验，得出了令人满意的方案。目标是：

^*大于20dB的动态范围(最好25dB)，以得到超过US 5,473,632的特性的有效边沿

^*类似于与S/N-0dB相应的噪声特性。(也就是说，不想使该方案毫无必要地放大噪声的影响)。

以上所确定的结果利用了n₁＝1、n₂＝25时所得到的加权因数，其特性如图11和12所示。

下面马上将看到在更复杂的CIR情况下本算法的特性。所关心的是，斜坡算法能提供多大的动态范围，如图13中的例子所示。

这一例子对噪声极为敏感。图14示出了S/N＝40dB时的特性。可看到，甚至这样小量的噪声也被大大地放大了。同样大小的噪声在优选的如以上所确定的结果中产生的噪声不明显，如图15所示。

图16和17示出了一些更复杂的CIR情况的结果。

还可以看到，本实施方式能使该算法进行峰间辨别，而图18示出了近的峰峰间隔，时间间隙大于约5μs。

在干扰消除的内容中，当估算信道脉冲响应矢量h时，也需要上述矩阵的逆阵：

Φh＝ψ    (2.14)

其中，Φ为本地产生的信号之间的互相关矩阵，而ψ为本地产生的信号与所接收信号之间的互相关矩阵。在这一方案中，同样可以采用上述原理和方法。

类似地，在干扰受限的环境中，通常利用从多个传感器(天线阵)接收到的信号的最佳组合来应用CDMA。最佳组合使信号对干扰加噪声比(SINR)最大。对于天线阵，为了使SINR最大，其权矢量例如可用下式得到：

w＝αR^-1U_d ^*(2.15)

其中，w是天线元的权矢量，α是一个常数，R是所接收的干扰加噪声相关矩阵，而U_d ^*是所需信号矢量的共轭。我们已经知道通过先计算CIR来得到权矢量的有效方法。一些最佳组合技术先计算出CIR(和噪声协方差矩阵Q)，然后将其应用于多维的MLSE，以对来自天线阵的信号进行组合。通常，在移动无线信道中，矩阵求逆会成问题。上述实施方式可直接应用于矩阵求逆。为得到权矢量或要利用Q的逆阵，就要对信号进行最佳组合。

总之，当本申请要完成一些选择和变形时，所公开的本发明的优选实施方法是，采用方程(2.8)所指定的矩阵乘法来得到CIR，然后用SVD来确定矩阵A的逆阵。根据n₁＝1、n₂＝25时的方程(2.13)，使逆阵中所含的本征值成斜坡(被加权)。

Claims (6)

1.一种估算通信系统的信道脉冲响应(CIR)的方法，通过通信系统的估算预测模型来应用，该方法包括以下步骤：

通过使所测量的接收信号与所述接收信号的预估算之间的误差的平方最小化来推算CIR，所述估算基于无线信道的近似；

产生所述CIR的第一矩阵形式；

将所述第一矩阵变换为对角矩阵；和

对所述对角矩阵求逆。

2.如权利要求1所述的方法，其特征在于，将奇异值分解(SVD)运用于第一矩阵，来得到对角矩阵。

3.如权利要求1或权利要求2任一所述的方法，在紧随矩阵求逆步骤之后，还包括以下步骤：斜坡降低求逆后的对角矩阵中所含的过大的本征值的影响。

4.如权利要求3所述的方法，其特征在于，用窗口函数来斜坡减小过大的本征值。

5.如权利要求4所述的方法，其特征在于，用Bartlett、Blackman、Harris、Gaussian、Poisson、Riesz、Hanning(升余弦)、Kaiser-Bessel、Lanczos、Tukey和其他类似的众所周知的窗口函数的任意之一来斜坡减小所述过大的本征值。

6.如上述权利要求任一所述的方法，可用于GSM网、干扰消除或CDMA接收机中。

Images