MATLAB 基于LMS算法的通道失配校正算法仿真

MATLAB 基于LMS算法的通道失配校正算法仿真

上一章节介绍了LMS算法LMS算法_小小低头哥的博客-CSDN博客，这次介绍基于LMS算法的通道失配校正校正结构

1.利用线性滤波器实现通道失配校正

  在数字信号处理系统中，有时候需要在接收端同时接收多个通道信号，由于每个通道包含多个模拟处理单元的组成不完全相同，导致通道之间出现频率响应的幅度与相位的不一致性，即通道失配。例如，在天线阵系统中，如果多个接收通道之间的频率响应不一致，将严重影响天线阵列的性能。因此，在进行阵列数字信号处理前，必须对通道间的失配进行校正，以保证经系统后续处理后获得较高的性能。

  所有的校正方法可分别应用于下面两种情况:随频率变化的通道失配情况，以及与频率响应无关的恒定幅相误差校正问题。本节的实例正是针对独立于频率响应的恒定通道失配问题采用的一种十分易于工程实现的校正方法。

假设第k个通道的频率响应为$H_k(jw)$，以参考通道KaTeX parse error: Expected group as argument to '\H' at position 3: \H_̲{ref}(jw)为基准，第k个通道相对于它的通道失配特性定义为
$$D_k(jw)=H_k(jw)/H_{ref}(jw)\text{\hspace{1em}(1)}$$
考虑一个频率响应为$H_k(jw)$的通道，令
$$D_k(jw)=H(jw)/H_{ref}(jw)=[1+\Delta H(jw)]e^{j\Delta \phi (w)}\text{\hspace{1em}(2)}$$
对于与频率无关的恒定幅相误差校正来说，$\Delta H(jw)$及$\Delta \phi (jw)$均为常数，不妨设$\Delta H(jw)$为$\delta$，$\Delta \phi (jw)$为$\theta$，则
$$D_k(jw)=H(jw)/H_{ref}(jw)=[1+\delta ]e^{j\theta }\text{\hspace{1em}(3)}$$
对于一般通信系统来讲，接收机输入的信号可以看成平稳信号，因此在一次通信过程中可以认为天线阵元通道的特性是不变的。这样，每次通信之前对阵列通道进行校正，校正完成后再开始正常通信，即从时间上把校正和正常通信分开，避免注入信号对后续的阵列处理产生影响。不失一般性，本次以两路通道为例对通道失配校正进行分析。

  对于与频率响应无关的恒定幅相误差校正来说，由以上分析可知，通道失配校正就是要找出校正因子，即
$$C(jw)=[1/(1+\delta )]e^{-j\theta }\text{\hspace{1em}(4)}$$
  基于LMS算法的校正结构，如图1所示。

图1通道失配校正模型

  在校正时，注入信号同时送至参考通道和待校正通道，当LMS算法收敛后，参考信号与经校正后的通道输出信号之间的均方差值趋向最小，从而迫使$D(jw)C(jw)$的乘积近似为1，也即待校正通道得到了校正。图7-10中参考信号通道中的延时时间等于通道中与校正因子C(jw)相乘所需要的时间。由于需要同时校正幅度及相位因子，因此LMS算法中的运算为复数运算，且滤波器长度为1，也即只需一级复加权就可以实现幅度与相位校正。

2.校正算法的MATLAB仿真

  对于通道失配校正算法的仿真，有两种仿真模式:一是仿真算法的性能、找出校正算法所能校正的幅度及相位失配度，以及校正后的精度;二是针对确定的失配度仿真校正前后的输入输出数据，并与FPGA实现后的测试仿真数据进行对比，并作为判断FPGA实现是否正确的依据。

2.1对基于LMS算法的通道失配校正模型进行仿真

  对基于LMS 算法的通道失配校正模型进行仿真。要求分别进行算法的性能仿真、生成FPGA测试用的输入信号、仿真权值在运算过程中的数据范围。对于性能仿真，要求仿真不同信噪比条件下校正前后的通道失配度改善情况，参考信号为叠加有噪声的随机方波信号，步长因子为1/64。

function output=My_E7_2_lms_mismatch_in(SNR)
%对基于LMS算法的通道失配校正模型进行仿真。要求分别进行算法的性能仿真、
%生成FPGA测试用的输入信号、仿真权值在运算过程中数据范围。对于性能仿真，要
%求仿真不同信噪比条件下校正前后的通道失配度改善情况，参考信号为叠加有噪声的随机方波信号。

%定义参数
N = 4000;   %数据长度4000
output = zeros(3,N);    %初始化输出信号

%计算信噪比所对应的增益比值 噪声、信号功率分别为pn、ps,
%则增益比分别为pn/(pn+ps) ps/(pn+ps) 化简为1/(1+ps/pn) 1/(pn/ps+1)
psn = 10^(SNR/20);  %ps/pn
normal_gs = 1/(1/psn + 1);  %信号增益比
normal_gn = 1/(1 + psn);    %噪声增益比

%产生信号
rs = randn(1,N)>0;  %生成随机方波信号
noise1 = randn(1,N);%产生2路互不相关的高斯噪声序列，仿真通过模拟信道的噪声,
noise2 = randn(1,N);%如果是L路通道，则应有L个噪声，每个通道一个

%产生混有声的正交参考注入信号 两个通道传输 乘以i的代表第二通道（正交）
%虽然第一通道才有原始信号，但是每个通道都难以避免有噪声，且是独立的
rn = (normal_gs*rs+normal_gn*noise1)+normal_gn*noise2*i;

%产生失配信号作为FPGA测试输入信号
mis_amp = 1.5;  %失配幅度
mis_ang = pi/4; %失配角度
xn = rn.*(1+mis_amp)*exp(j*mis_ang);    %D(jw)=[1+a]exp(j*b);

%对参考信号及测试输入信号进行归一化处理
rr = real(rs);  mrr = max(abs(rr)); %得到实部最大值
ri = imag(rs);  mri = max(abs(ri)); %得到虚部最大值
xr = real(xn);  mxr = max(abs(xr)); %得到实部最大值
xi = imag(xn);  mxi = max(abs(xi)); %得到虚部最大值
m = max([mrr,mri,mxr,mxi])
rr=rr/m;ri=ri/m;xr=xr/m;xi=xi/m;

%量化
Q_rr=round(rr*(2^(16-1)-1));%16比特量化
Q_ri=round(ri*(2^(16-1)-1));%16比特量化
Q_xr=round(xr*(2^(16-1)-1));%16比特量化
Q_xi=round(xi*(2^(16-1)-1));%16比特量化
%fid=fopen('\\\','w'); 写入文件中供FPGA读取
%for k=1:length(Q_xi)
%    B_si=dec2bin(Q_xi(k)+(Q_xi(k)<0)*2^Q,Q);
%    for q=1:Q
%       if B_si(q)=='1'
%           tb=1;
%       else
%           tb=0;
%       end
%       fprintf(fid,'%d',tb);  
%    end
%    fprintf(fid,'\r\n');
%end
%fprintf(fid,';'); 
%fclose(fid);
output(1,:)=rn;    %原始参考信号
output(2,:)=rr+ri*j;
output(3,:)=xr+xi*j;

clear 
u = 1/128;
SNR = -10;  %信噪比
output = My_E7_2_lms_mismatch_in(SNR);
x = output(1,:);    %得到混有声的正交参考注入信号
N = 15; %运行15次
for i=1:N   
    mis_amp = i*0.1;    %步长为0.1 范围0.1-1.5
    mis_ang = pi/4;     %失配角度
    mis = (1+mis_amp)*exp(j*mis_ang);   %失配值
    xinput = x*mis;
    %w=zeros(1,length(x)+1);w(1)=1;
    w = 1;  %LMS滤波器系数长度为1 
    e = zeros(1,length(x));
    for n=1:length(x)
        y = w'* xinput(n);   %别忘了w需要共轭转置'
        e(n) = x(n) - y;
        %w(n+1)=w(n)+2*u*e(n)'*xinput(n);%别忘了e也需要共轭转置'
        w = w + 2*u*e(n)'*xinput(n);
    end     
    %看abs(w(w可能为复数))*(1+mis_amp)*rn是否会最终等于rn
    %即(1+mis_amp)*abs(w)是否为1，偏差范围多少
    amp_10(i)=abs( (1+mis_amp) * abs(w) -1 );
    %amp_10(i)=abs( (1+mis_amp) * abs(w(length(x))) -1 );
end

SNR = 0;  %信噪比
output = My_E7_2_lms_mismatch_in(SNR);
x = output(1,:);    %得到混有声的正交参考注入信号
N = 15; %运行15次
for i=1:N   
    mis_amp = i*0.1;    %步长为0.1 范围0.1-1.5
    mis_ang = pi/4;     %失配角度
    mis = (1+mis_amp)*exp(j*mis_ang);   %失配值
    xinput = x*mis;
    %w=zeros(1,length(x)+1);w(1)=1;
    w = 1;  %LMS滤波器系数长度为1 
    e = zeros(1,length(x));
    for n=1:length(x)
        y = w'* xinput(n);
        e(n) = x(n) - y;
        %w(n+1)=w(n)+2*u*e(n)'*xinput(n);
        w = w + 2*u*e(n)'*xinput(n);
    end     
    %看abs(w(w可能为复数))*(1+mis_amp)*rn是否会最终等于rn
    %即(1+mis_amp)*abs(w)是否为1，偏差范围多少
    amp0(i)=abs( (1+mis_amp) * abs(w) -1 );
    %amp0(i)=abs( (1+mis_amp) * abs(w(length(x))) -1 );
end  


SNR = 10;  %信噪比
output = My_E7_2_lms_mismatch_in(SNR);
x = output(1,:);    
N = 15; 
for i=1:N   
    mis_amp = i*0.1;    
    mis_ang = pi/4;     
    mis = (1+mis_amp)*exp(j*mis_ang);  
    xinput = x*mis;     
    w = 1;  
    e = zeros(1,length(x));
    for n=1:length(x)
        y = w' * xinput(n);
        e(n) = x(n) - y;
        w = w + 2*u*e(n)'*xinput(n);
    end     
    amp10(i)=abs( (1+mis_amp) * abs(w) -1 );
end

snr=-10;
output=My_E7_2_lms_mismatch_in(snr);
x=output(1,:);
for i=-6:6
    mis_amp = 0.5;    %失配幅度不变
    mis_ang = i*0.1*pi;   %失配角度变化
    mis = (1+mis_amp)*exp(j*mis_ang);
    xinput = x*mis;
    w = 1;  
    e = zeros(1,length(x));
    for n=1:length(x)
       y = w' * xinput(n);
       e(n) = x(n)-y;
       w = w + 2*u*e(n)'*xinput(n);
    end
    %angle(mis) 表示失配系数 mis 的幅角，angle(w) 表示校正后的
    % 权值系数 w 的幅角，二者的差值表示实际失配角度和校正后的
    % 权值系数之间的角度误差
    ang_10(i+7)=abs(angle(mis) - angle(w))*180/pi;
end

snr=0;
output=My_E7_2_lms_mismatch_in(snr);
x=output(1,:);
for i=-6:6
    mis_amp = 0.5;    
    mis_ang = i*0.1*pi;   
    mis = (1+mis_amp)*exp(j*mis_ang);
    xinput = x*mis;
    w = 1;  
    e = zeros(1,length(x));
    for n=1:length(x)
       y = w' * xinput(n);
       e(n) = x(n)-y;
       w = w + 2*u*e(n)'*xinput(n);
    end
   ang0(i+7)=abs(angle(mis) - angle(w))*180/pi;
end

snr=10;
output=My_E7_2_lms_mismatch_in(snr);
x=output(1,:);
for i=-6:6
    mis_amp = 0.5;    
    mis_ang = i*0.1*pi;   
    mis = (1+mis_amp)*exp(j*mis_ang);
    xinput = x*mis;
    w = 1;  
    e = zeros(1,length(x));
    for n=1:length(x)
       y = w' * xinput(n);
       e(n) = x(n)-y;
       w = w + 2*u*e(n)'*xinput(n);
    end
   ang10(i+7)=abs(angle(mis) - angle(w))*180/pi;
end

figure(2)
len=1:15;
len=len*0.1;
subplot(211);plot(len,amp_10,'+',len,amp0,'-',len,amp10,'-.');legend('snr=-10dB','snr=0dB','snr=10dB');
xlabel('校正前幅度失配度');ylabel('校正后幅度失配度');
grid;
len=-6:6;
len=len*18;
subplot(212);plot(len,ang_10,'+',len,ang0,'-',len,ang10,'-.');legend('snr=-10dB','snr=0dB','snr=10dB');
xlabel('校正前相位失配度');ylabel('校正后相位失配度');
grid;

2.2 实验结果及分析

图2 通道失配矫正算法性能仿真

由图2可知，当幅度失配度为0 ~ 1.5时，基本能够完全矫正到一致；当相位失配度在-150°~150°时，也可以完全矫正到一致。输入信噪比越大，信号质量就越好，校正效果就越好。通过仿真可知，该矫正方法可以较好的校正与频率响应无关的恒定幅相误差。

3.相关疑问解答（ChatGPT）

请你以数学表达式的形式描述整个过程，比如一开始失配的数学表达式是什么，失配过程中数学表达式又是怎么变的？最后得到校正后的数学表达式又是什么？

失配校正的目的是在失配的情况下，通过调整参考信号和权值系数之间的关系，使得输出信号能够逆向抵消失配引入的相位偏差。具体来说，失配校正的过程可以用以下数学表达式表示：

假设原始信号为 $x[n]$，失配后的参考信号为 $r[n]$，失配系数为 $m$，校正后的权值系数为 $w[n]$，则失配校正的过程可以用以下步骤表示：

1. 原始信号 $x[n]$ 经过混合得到失配后的参考信号 $r[n]$：

$$r[n]=m\cdot x[n]$$

其中，$m$ 表示失配系数，通常由失配角度和失配幅度确定。

1. 使用LMS算法对失配后的参考信号 $r[n]$ 进行校正，得到校正后的权值系数 $w[n]$：

$$w[n]=w[n-1]+2\mu e[n]x[n]$$

其中，$w[n]$ 表示第 $n$ 个时刻的权值系数，$w[n-1]$ 表示第 $n-1$ 个时刻的权值系数，$\mu$ 表示学习率，$e[n]$ 表示当前时刻的误差，$x[n]$ 表示当前时刻的输入信号。

1. 使用校正后的权值系数 $w[n]$ 对输入信号 $x[n]$ 进行预测，得到输出信号 $y[n]$：

$$y[n]=w^H[n]\cdot x[n]$$

其中，$w^H[n]$ 表示权值系数的共轭转置。

1. 计算误差 $e[n]$：

$$e[n]=x[n]-y[n]$$

1. 重复步骤2到步骤4，直到收敛或达到预设的迭代次数。

在失配校正的过程中，失配系数 $m$ 和校正后的权值系数 $w[n]$ 都是复数，其幅角分别表示失配角度和校正后的相位偏差。通过调整校正后的权值系数 $w[n]$，可以使其幅角与失配系数 $m$ 的幅角相等，从而达到逆向抵消失配引入的相位偏差的目的

为什么幅角与失配系数 $m$ 的幅角相等，就能达到逆向抵消失配引入的相位偏差的目的？

当失配系数 $m$ 和校正后的权值系数 $w[n]$ 的幅角相等时，它们的幅角差值为0，即 $|angle(m)-angle(w[n])|=0$，说明它们的相位偏差相等，但方向相反。因此，如果将校正后的权值系数 $w[n]$ 乘以失配系数 $m$ 的共轭，即 $w^{\ast }[n]m$，则可以得到一个复数，其幅角为0，表示没有相位偏差。因此，将输出信号 $y[n]$ 乘以 $w^{\ast }[n]m$，即可将失配引入的相位偏差逆向抵消，得到一个相位正确的输出信号：

$$y_c[n]=y[n]\cdot w^{\ast }[n]m$$

其中，$w^{\ast }[n]$ 表示 $w[n]$ 的共轭，$m$ 表示失配系数。由于 $w^{\ast }[n]m$ 的幅角为0，因此 $y_c[n]$ 的相位与参考信号 $r[n]$ 的相位一致，即消除了失配引入的相位偏差。因此，当幅角与失配系数 $m$ 的幅角相等时，就能达到逆向抵消失配引入的相位偏差的目的。

为什么方向相反，是因为求输出的时候权值系数$w[n]$需要共轭转置吗？

是的，当使用复数权值系数进行信号处理时，一般需要对其进行共轭转置操作，以实现正确的信号处理。具体来说，对于复数权值系数 $w[n]$，其共轭转置表示为 $w^{\ast }[n]$，即将 $w[n]$ 中的所有元素取共轭并转置得到的矩阵。

在失配校正中，输出信号的计算公式为 $y[n]=w^H[n]\cdot x[n]$，其中 $w^H[n]$ 表示权值系数的共轭转置。因此，如果失配系数 $m$ 和校正后的权值系数 $w[n]$ 的幅角相等，即 $|angle(m)-angle(w[n])|=0$，则可以将 $w[n]$ 的共轭转置 $w^{\ast }[n]$ 作为逆向抵消失配引入的相位偏差的系数，即：

$$y_c[n]=y[n]\cdot w^{\ast }[n]m$$

其中，$w^{\ast }[n]m$ 表示逆向抵消失配引入的相位偏差的系数，$m$ 表示失配系数。由于 $w^{\ast }[n]m$ 的幅角为0，因此实现了逆向抵消失配引入的相位偏差。